use isce3::geocode::geocodeSlc array mode instead of raster mode in…

… `s1_geocode_slc` (#127)

* array mode in s1_geocode_slc to replace raster mode
* default geo block to NaN
* replace has rdr bbox logic with raster fill value
* replace block processing with processing of entire array
* add runtime convenience function
* updated to write geocoded phase to HDF5
* remove block proc from schema
* add missing imports
* add output options for compression and different complex types
* wrap phases
* fix unit runconfig template
* pin isce3 version

requirements.txt

@@ -2,7 +2,7 @@ python>=3.9     # sentinel1-reader requirement
 numpy           # sentinel1-reader requirement
 lxml            # sentinel1-reader requirement
 gdal>=3
-#isce3        # since the conda-installed isce3 is not the most updated version, installing isce3 from stratch is recommended, to stay in sync with isce3 development.
+isce3>=0.13.0   # since the conda-installed isce3 is not the most updated version, installing isce3 from stratch is recommended, to stay in sync with isce3 development.
 #journal      # as of Mar 2022, journal from conda does not support python3.9; since it is included during isce3 installation above, comment this out temporarily.
 pandas
 pyproj

src/compass/defaults/s1_cslc_geo.yaml

@@ -46,7 +46,6 @@ runconfig:
           polarization: co-pol
           geocoding:
               flatten: True
-              lines_per_block: 1000
               x_posting:
               y_posting:
               x_snap:
@@ -123,3 +122,16 @@ runconfig:
               equalize: False
           perform_qa: True
           output_to_json: False
+
+      # Output options for CSLC raster (e.g. compression, chunking, shuffle filter)
+      output:
+          # Data type of CSLC raster
+          cslc_data_type: complex64_zero_mantissa
+          # Enable gzip compression of rasters
+          compression_enabled: True
+          # Level of compression applied to rasters
+          compression_level: 4
+          # Chunk size of rasters
+          chunk_size: [128, 128]
+          # Enable shuffle filtering of rasters
+          shuffle: True

src/compass/s1_geo2rdr.py

@@ -2,14 +2,13 @@
 
 '''wrapper for geo2rdr'''
 
-from datetime import timedelta
 import os
 import time
 
 import isce3
 import journal
 
-from compass.utils.helpers import get_module_name
+from compass.utils.helpers import get_module_name, get_time_delta_str
 from compass.utils.runconfig import RunConfig
 from compass.utils.yaml_argparse import YamlArgparse
 
@@ -85,7 +84,7 @@ def run(cfg: dict):
         # Execute geo2rdr
         geo2rdr_obj.geo2rdr(topo_raster, out_paths.output_directory)
 
-    dt = str(timedelta(seconds=time.time() - t_start)).split(".")[0]
+    dt = get_time_delta_str(t_start)
     info_channel.log(f"{module_name} burst successfully ran in {dt} (hr:min:sec)")
 
 

src/compass/s1_geocode_metadata.py

@@ -2,7 +2,6 @@
 
 '''wrapper to geocode metadata layers'''
 
-from datetime import timedelta
 import time
 
 import h5py
@@ -20,12 +19,14 @@
                                       init_geocoded_dataset,
                                       metadata_to_h5group, DATA_PATH,
                                       ROOT_PATH)
-from compass.utils.helpers import bursts_grouping_generator, get_module_name
+from compass.utils.helpers import (bursts_grouping_generator, get_module_name,
+                                   get_time_delta_str)
 from compass.utils.yaml_argparse import YamlArgparse
 from compass.utils.radar_grid import get_decimated_rdr_grd
 
 
-def _fix_layover_shadow_mask(static_layers_dict, h5_root, geo_grid):
+def _fix_layover_shadow_mask(static_layers_dict, h5_root, geo_grid,
+                             output_params):
     '''
     kludge correctly mask invalid pixel in geocoded layover shadow to address
     isce3::geocode::geocodeCov's inability to take in an user defined invalid
@@ -59,7 +60,7 @@ def _fix_layover_shadow_mask(static_layers_dict, h5_root, geo_grid):
         _ = init_geocoded_dataset(h5_root[DATA_PATH], dst_ds_name, geo_grid,
                                   dtype=None,
                                   description=np.string_(dst_ds_name),
-                                  data=temp_arr)
+                                  data=temp_arr, output_cfg=output_params)
 
 
 def run(cfg, burst, fetch_from_scratch=False):
@@ -162,17 +163,21 @@ def run(cfg, burst, fetch_from_scratch=False):
             if not enabled:
                 continue
 
+            # init value is invalid value for the single/float32
             dtype = np.single
             # layoverShadowMask is last option, no need to change data type
             # and interpolator afterwards
             if dataset_name == 'layover_shadow_mask':
                 geocode_obj.data_interpolator = 'NEAREST'
                 dtype = np.byte
+                # layover shadow is a char (no NaN char, 0 represents unmasked
+                # value)
 
             # Create dataset with x/y coords/spacing and projection
             topo_ds = init_geocoded_dataset(static_layer_data_group,
                                             dataset_name, geo_grid, dtype,
-                                            np.string_(dataset_name))
+                                            np.string_(dataset_name),
+                                            output_cfg=cfg.output_params)
 
             # Init output and input isce3.io.Raster objects for geocoding
             output_raster = isce3.io.Raster(f"IH5:::ID={topo_ds.id.id}".encode("utf-8"),
@@ -193,7 +198,8 @@ def run(cfg, burst, fetch_from_scratch=False):
             del output_raster
 
             if dataset_name == 'layover_shadow_mask':
-                _fix_layover_shadow_mask(static_layers, h5_root, geo_grid)
+                _fix_layover_shadow_mask(static_layers, h5_root, geo_grid,
+                                         cfg.output_params)
 
     if cfg.quality_assurance_params.perform_qa:
         cslc_qa = QualityAssuranceCSLC()
@@ -204,13 +210,13 @@ def run(cfg, burst, fetch_from_scratch=False):
             if cfg.quality_assurance_params.output_to_json:
                 cslc_qa.write_qa_dicts_to_json(out_paths.stats_json_path)
 
-    dt = str(timedelta(seconds=time.time() - t_start)).split(".")[0]
+    dt = get_time_delta_str(t_start)
     info_channel.log(
         f"{module_name} burst successfully ran in {dt} (hr:min:sec)")
 
 
 def geocode_luts(geo_burst_h5, burst, cfg, dst_group_path, item_dict,
-                 dec_factor_x_rng=20, dec_factor_y_az=5):
+                 output_params, dec_factor_x_rng=20, dec_factor_y_az=5):
     '''
     Geocode the radiometric calibration parameters,
     and write them into output HDF5.
@@ -295,7 +301,8 @@ def geocode_luts(geo_burst_h5, burst, cfg, dst_group_path, item_dict,
                               item_name,
                               decimated_geogrid,
                               'float32',
-                              f'geocoded {item_name}')
+                              f'geocoded {item_name}',
+                              output_cfg=cfg.output_params)
 
         dst_dataset = geo_burst_h5[f'{dst_group_path}/{item_name}']
 
@@ -398,8 +405,7 @@ def geocode_calibration_luts(geo_burst_h5, burst, cfg,
               None]
         }
     geocode_luts(geo_burst_h5, burst, cfg, dst_group_path, item_dict_calibration,
-                 dec_factor_x_rng,
-                 dec_factor_y_az)
+                 cfg.output_params, dec_factor_x_rng, dec_factor_y_az)
 
 
 def geocode_noise_luts(geo_burst_h5, burst, cfg,
@@ -430,8 +436,7 @@ def geocode_noise_luts(geo_burst_h5, burst, cfg,
                                        burst.burst_noise.azimuth_lut]
                                        }
     geocode_luts(geo_burst_h5, burst, cfg, dst_group_path, item_dict_noise,
-                 dec_factor_x_rng,
-                 dec_factor_y_az)
+                 cfg.output_params, dec_factor_x_rng, dec_factor_y_az)
 
 
 if __name__ == "__main__":