diff --git a/python/packages/CMakeLists.txt b/python/packages/CMakeLists.txt
index 9b2f62f6a..17d660f33 100644
--- a/python/packages/CMakeLists.txt
+++ b/python/packages/CMakeLists.txt
@@ -39,6 +39,7 @@ set (list_of_exe
nisar/workflows/gslc.py
nisar/workflows/gslc_point_target_analysis.py
nisar/workflows/insar.py
+ nisar/workflows/insar_kz.py
nisar/workflows/nisar_l0b_dm1_to_science.py
nisar/workflows/nisar_l0b_dm2_to_dbf.py
nisar/workflows/noise_estimator.py
diff --git a/python/packages/nisar/workflows/insar.py b/python/packages/nisar/workflows/insar.py
index 91ee8c862..08d34f997 100644
--- a/python/packages/nisar/workflows/insar.py
+++ b/python/packages/nisar/workflows/insar.py
@@ -2,10 +2,12 @@
import time
import journal
+import isce3
from nisar.workflows import (bandpass_insar, baseline, crossmul, dense_offsets,
filter_interferogram, geo2rdr, geocode_insar,
h5_prep, ionosphere, offsets_product,
- prepare_insar_hdf5, rdr2geo, resample_slc_v2,
+ prepare_insar_hdf5, rdr2geo,
+ resample_slc, resample_slc_v2,
rubbersheet, solid_earth_tides, split_spectrum,
troposphere, unwrap)
from nisar.workflows.geocode_insar import InputProduct
@@ -36,7 +38,13 @@ def run(cfg: dict, out_paths: dict, run_steps: dict):
prepare_insar_hdf5.run(cfg)
if run_steps['coarse_resample']:
- resample_slc_v2.run(cfg, 'coarse')
+ # Resample SLC V2 doesn't have GPU implemented, so run the old version if
+ # GPU was requested.
+ if isce3.core.gpu_check.use_gpu(cfg['worker']['gpu_enabled'],
+ cfg['worker']['gpu_id']):
+ resample_slc.run(cfg, 'coarse')
+ else:
+ resample_slc_v2.run(cfg, 'coarse')
if (run_steps['dense_offsets']) and \
(cfg['processing']['dense_offsets']['enabled']):
@@ -57,7 +65,13 @@ def run(cfg: dict, out_paths: dict, run_steps: dict):
and cfg['processing']['fine_resample']['enabled']
and 'RIFG' in out_paths
):
- resample_slc_v2.run(cfg, 'fine')
+ # Resample SLC V2 doesn't have GPU implemented, so run the old version if
+ # GPU was requested.
+ if isce3.core.gpu_check.use_gpu(cfg['worker']['gpu_enabled'],
+ cfg['worker']['gpu_id']):
+ resample_slc.run(cfg, 'fine')
+ else:
+ resample_slc_v2.run(cfg, 'fine')
# If fine_resampling is enabled, use fine-coregistered SLC
# to run crossmul
diff --git a/python/packages/nisar/workflows/insar_kz.py b/python/packages/nisar/workflows/insar_kz.py
new file mode 100644
index 000000000..4260c4363
--- /dev/null
+++ b/python/packages/nisar/workflows/insar_kz.py
@@ -0,0 +1,173 @@
+#!/usr/bin/env python3
+import time
+
+import journal
+import tempfile
+from osgeo import gdal
+import numpy as np
+import plant
+
+from nisar.workflows import (bandpass_insar, crossmul, dense_offsets, geo2rdr,
+ geocode_insar, h5_prep, filter_interferogram,
+ offsets_product, prepare_insar_hdf5, rdr2geo,
+ resample_slc, rubbersheet,
+ split_spectrum, unwrap, ionosphere, baseline,
+ troposphere, solid_earth_tides)
+
+from nisar.workflows.geocode_insar import InputProduct
+from nisar.workflows.insar_runconfig import InsarRunConfig
+from nisar.workflows.persistence import Persistence
+from nisar.workflows.yaml_argparse import YamlArgparse
+
+
+def _generate_dem_with_offset(dem_file, dem_offset, scratch_dir):
+ '''
+ Generate DEM with offset. It may require:
+
+ >> export GDAL_VRT_ENABLE_PYTHON=YES
+ '''
+
+ dem_with_offset = tempfile.NamedTemporaryFile(
+ dir=scratch_dir, suffix='.vrt').name
+
+ print('building DEM VRT with height offset [m]:', dem_offset)
+ gdal.BuildVRT(dem_with_offset, dem_file)
+
+ pixel_function_str = (
+ ' offset\n'
+ ' Python\n'
+ f' \n'
+ ' \n'
+ ' \n')
+
+ with open(dem_with_offset, "r") as f:
+ lines = f.readlines()
+
+ with open(dem_with_offset, "w") as f:
+ for line in lines:
+ if '', ' subClass="VRTDerivedRasterBand">')
+ if '' in line or '' in line:
+ f.write(pixel_function_str)
+ f.write(line)
+
+ return dem_with_offset
+
+
+def run(cfg: dict, out_paths: dict, run_steps: dict):
+ '''
+ Run INSAR workflow with parameters in cfg dictionary
+ '''
+ info_channel = journal.info("insar.run")
+ info_channel.log("starting INSAR")
+
+ t_all = time.time()
+
+ # if run_steps['bandpass_insar']:
+ # bandpass_insar.run(cfg)
+ scratch_dir = cfg['product_path_group']['scratch_path']
+
+ # save original parameters
+ product_type_orig = cfg['primary_executable']['product_type']
+ print('product_type_orig:', product_type_orig)
+ # scratch_dir
+ # rdr2geo, geo2rdr, coarse_resample enabled
+
+ # k_z is computed as the InSAR phase difference using dem + 0.5 and
+ # dem - 0.5
+
+ dem_file = cfg['dynamic_ancillary_file_group']['dem_file']
+
+ for dem_offset in [0.5, 0, -0.5]:
+
+ if dem_offset == 0.5:
+ dem_offset_str = '_p0.5'
+ elif dem_offset == -0.5:
+ dem_offset_str = '_m0.5'
+ elif dem_offset == 0:
+ dem_offset_str = '_zero'
+ else:
+ raise NotImplementedError("Invalid option")
+
+ # udpate dem with dem_offset, create VRT and update dem_file ?
+
+ # write new parameters
+ cfg['primary_executable']['product_type'] = 'RIFG'
+
+ dem_with_offset = _generate_dem_with_offset(dem_file, dem_offset,
+ scratch_dir)
+
+ cfg['dynamic_ancillary_file_group']['dem_file'] = dem_with_offset
+
+ # create RFIG writer
+ if run_steps['prepare_insar_hdf5']:
+ prepare_insar_hdf5.run(cfg)
+
+ if run_steps['rdr2geo']:
+ rdr2geo.run(cfg)
+
+ if run_steps['geo2rdr']:
+ geo2rdr.run(cfg)
+
+ if run_steps['coarse_resample']:
+ resample_files_dict = resample_slc.run(cfg, 'coarse', flatten=True,
+ flag_constant_value=True,
+ suffix=dem_offset_str)
+
+ if dem_offset == 0.5:
+ resample_files_dict_p_05 = resample_files_dict
+ elif dem_offset == -0.5:
+ resample_files_dict_m_05 = resample_files_dict
+
+ for (_, pol_dict_p_05), (_, pol_dict_m_05) in \
+ zip(resample_files_dict_p_05.items(),
+ resample_files_dict_m_05.items()):
+
+ for (_, f_p_05), (_, f_m_05) \
+ in zip(pol_dict_p_05.items(), pol_dict_m_05.items()):
+
+ print(f'DEM +0.5 file: {f_p_05}')
+ print(f'DEM -0.5 file: {f_m_05}')
+ data_p_05_band = gdal.Open(f_p_05, gdal.GA_ReadOnly)
+ data_p_05 = data_p_05_band.GetRasterBand(1).ReadAsArray()
+ data_m_05_band = gdal.Open(f_m_05, gdal.GA_ReadOnly)
+ data_m_05 = data_m_05_band.GetRasterBand(1).ReadAsArray()
+ data_diff = np.angle(data_p_05 * np.conj(data_m_05))
+ output_filename = f_p_05.replace('_p0.5', '_kz')
+ plant.save_image(data_diff, output_filename, force=True)
+ print(f'file saved: {output_filename}')
+
+ print('*** resample_files_dict:', resample_files_dict)
+
+ # if run_steps['baseline']:
+ # baseline.run(cfg, out_paths)
+
+ t_all_elapsed = time.time() - t_all
+ info_channel.log(f"successfully ran INSAR in {t_all_elapsed:.3f} seconds")
+
+
+if __name__ == "__main__":
+ # parse CLI input
+ yaml_parser = YamlArgparse()
+ args = yaml_parser.parse()
+
+ # convert CLI input to run configuration
+ insar_runcfg = InsarRunConfig(args)
+
+ # To allow persistence, a logfile is required. Raise exception
+ # if logfile is None and persistence is requested
+ logfile_path = insar_runcfg.cfg['logging']['path']
+ if (logfile_path is None) and insar_runcfg.args.restart:
+ raise ValueError('InSAR workflow persistence requires to specify a logfile')
+ persist = Persistence(logfile_path, insar_runcfg.args.restart)
+
+ # run InSAR workflow
+ if persist.run:
+ _, out_paths = h5_prep.get_products_and_paths(insar_runcfg.cfg)
+
+ run(insar_runcfg.cfg, out_paths, persist.run_steps)
diff --git a/python/packages/nisar/workflows/resample_slc.py b/python/packages/nisar/workflows/resample_slc.py
index 7d3fff109..c6aa91eb0 100644
--- a/python/packages/nisar/workflows/resample_slc.py
+++ b/python/packages/nisar/workflows/resample_slc.py
@@ -7,17 +7,54 @@
import os
import pathlib
import time
+import tempfile
import journal
import isce3
from osgeo import gdal
from nisar.products.readers import SLC
-from nisar.workflows.helpers import copy_raster
+from nisar.workflows.helpers import copy_raster, complex_raster_path_from_h5
from nisar.workflows.resample_slc_runconfig import ResampleSlcRunConfig
from nisar.workflows.yaml_argparse import YamlArgparse
-def run(cfg, resample_type):
+def build_empty_vrt(filename, length, width, fill_value, dtype='CFloat32'):
+ """Build an empty VRT file, i.e, not pointing to any rasters,
+ with given input dimensions (length and width), data type, and
+ fill value.
+ Parameters
+ ----------
+ filename: str
+ VRT file name
+ length: int
+ VRT data length
+ width: int
+ VRT data width
+ fill_value: scalar
+ VRT data fill value
+ dtype: str
+ VRT data type
+ """
+
+ vrt_contents = f' \n'
+
+ vrt_contents += (
+ f' \n'
+ f' {fill_value} \n'
+ f' {fill_value} \n'
+ f' \n'
+ f' \n')
+
+ with open(filename, 'w') as out:
+ out.write(vrt_contents)
+
+ if os.path.isfile(filename):
+ print('file saved:', filename)
+
+
+def run(cfg, resample_type, flatten=False, flag_constant_value=False,
+ suffix=''):
'''
run resample_slc
'''
@@ -30,6 +67,10 @@ def run(cfg, resample_type):
# Get SLC parameters
slc = SLC(hdf5file=input_hdf5)
+ # If flattening is enabled, load the reference SLC
+ if flatten:
+ reference_hdf5 = cfg['input_file_group']['reference_rslc_file']
+ ref_slc = SLC(hdf5file=reference_hdf5)
info_channel = journal.info('resample_slc.run')
info_channel.log('starting resampling SLC')
@@ -45,7 +86,12 @@ def run(cfg, resample_type):
t_all = time.time()
+ resample_files_dict = {}
+
for freq in freq_pols.keys():
+
+ resample_files_dict[freq] = {}
+
# Get frequency specific parameters
radar_grid = slc.getRadarGrid(frequency=freq)
native_doppler = slc.getDopplerCentroid(frequency=freq)
@@ -79,7 +125,14 @@ def run(cfg, resample_type):
else:
Resamp = isce3.image.ResampSlc
- resamp_obj = Resamp(radar_grid, native_doppler)
+ # If flattening is enabled, add the reference SLC radar grid to
+ # the call to the constructor of the ResampSlc module
+ if flatten:
+ ref_radar_grid = ref_slc.getRadarGrid(frequency=freq)
+ resamp_obj = Resamp(radar_grid, native_doppler,
+ ref_rdr_grid=ref_radar_grid)
+ else:
+ resamp_obj = Resamp(radar_grid, native_doppler)
# If lines per tile is > 0, assign it to resamp_obj
if resamp_args['lines_per_tile']:
@@ -92,19 +145,46 @@ def run(cfg, resample_type):
# Create directory for each polarization
out_dir = resample_slc_scratch_path / pol
out_dir.mkdir(parents=True, exist_ok=True)
- out_path = out_dir / 'coregistered_secondary.slc'
+
+ out_path = out_dir / f'coregistered_secondary{suffix}.slc'
+
+ # If necessary, perform complex32 to complex64 conversion on input
+ input_as_c32_path = str(out_dir/'secondary.slc')
+ input_raster_path, _ = complex_raster_path_from_h5(
+ slc, freq, pol, input_hdf5, resamp_args['lines_per_tile'],
+ input_as_c32_path)
+
+ input_raster = isce3.io.Raster(input_raster_path)
# Dump secondary RSLC on disk
- raster_path = str(out_dir/'secondary.slc')
- copy_raster(input_hdf5, freq, pol, 1000,
- raster_path, file_type='ENVI')
- input_raster = isce3.io.Raster(raster_path)
+ # raster_path = str(out_dir/'secondary.slc')
+ # copy_raster(input_hdf5, freq, pol, 1000,
+ # raster_path, file_type='ENVI')
+ # input_raster = isce3.io.Raster(raster_path)
+
+ # Dump secondary RSLC on disk
+ # raster_path = str(out_dir/'secondary.slc')
+ # copy_raster(input_hdf5, freq, pol, 1000,
+ # raster_path, file_type='ENVI')
+ # input_raster = isce3.io.Raster(raster_path)
+
+ if flag_constant_value:
+ print('*** constant value!')
+ fill_value = 1
+ input_raster_vrt = tempfile.NamedTemporaryFile(
+ dir=scratch_path, suffix='.vrt').name
+ build_empty_vrt(input_raster_vrt, input_raster.length,
+ input_raster.width, fill_value)
+ input_raster = isce3.io.Raster(input_raster_vrt)
# Create output raster
resamp_slc = isce3.io.Raster(str(out_path), rg_off.width,
rg_off.length, rg_off.num_bands,
gdal.GDT_CFloat32, 'ENVI')
- resamp_obj.resamp(input_raster, resamp_slc, rg_off, az_off)
+ resamp_obj.resamp(input_raster, resamp_slc, rg_off, az_off,
+ flatten=flatten)
+
+ resample_files_dict[freq][pol] = str(out_path)
t_all_elapsed = time.time() - t_all
info_channel.log(f"successfully ran resample in {t_all_elapsed:.3f} seconds")
diff --git a/python/packages/nisar/workflows/resample_slc_v2.py b/python/packages/nisar/workflows/resample_slc_v2.py
index dac720ddb..a9b03da19 100644
--- a/python/packages/nisar/workflows/resample_slc_v2.py
+++ b/python/packages/nisar/workflows/resample_slc_v2.py
@@ -53,6 +53,8 @@ def run(cfg: dict, resample_type: str) -> None:
t_all = perf_counter()
+ resample_files_dict = {}
+
for freq in freq_pols.keys():
# Open offsets
@@ -73,7 +75,7 @@ def run(cfg: dict, resample_type: str) -> None:
raise ValueError(
"resample_type must be 'coarse' or 'fine', instead got "
f"{resample_type!r}")
-
+
az_off_path = offsets_path / "azimuth.off"
rg_off_path = offsets_path / "range.off"
@@ -91,11 +93,11 @@ def run(cfg: dict, resample_type: str) -> None:
# Get polarization list for which resample SLCs
pol_list = freq_pols[freq]
-
+
info_channel.log(f"Resampling SLC for frequency {freq}.")
t_freq_elapsed = -perf_counter()
- resample_secondary_rslc_onto_reference(
+ freq_output_files_dict = resample_secondary_rslc_onto_reference(
ref_file_path=ref_file_path,
sec_file_path=sec_file_path,
out_path=resample_slc_scratch_path,
@@ -107,6 +109,7 @@ def run(cfg: dict, resample_type: str) -> None:
block_size_rg=block_width,
with_gpu=with_gpu,
)
+ resample_files_dict[freq] = freq_output_files_dict
t_freq_elapsed += perf_counter()
info_channel.log(f"successfully ran resample for frequency {freq} in "
@@ -114,6 +117,7 @@ def run(cfg: dict, resample_type: str) -> None:
t_all_elapsed = perf_counter() - t_all
info_channel.log(f"successfully ran resample in {t_all_elapsed:.3f} seconds")
+ return resample_files_dict
def resample_secondary_rslc_onto_reference(
@@ -189,6 +193,8 @@ def resample_secondary_rslc_onto_reference(
mode='r+',
)
+ output_files_dict = {}
+
# For each polarization being output, create a GDALRaster to write to it.
out_writers: list[GDALRaster] = []
for pol in pols:
@@ -228,6 +234,7 @@ def resample_secondary_rslc_onto_reference(
fill_value=0.0 + 0.0j,
with_gpu=with_gpu,
)
+ return output_files_dict
if __name__ == "__main__":