Add geocode entrypoint

CHANGELOG.md

@@ -8,12 +8,16 @@ and uses [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
 ## [0.5.0]
 
+### Added
+* Geocode-only entrypoint.
+
 ### Changed
 * Update readme with docker instructions and layer info.
 
 ### Removed
 * Output of area projection metadata.
 * Prototype polar grid support in favor of full support via the new docker image.
+* Old per-file entrypoints.
 
 ### Fixed
 * Plotting issues with the cal/val submodule.

README.md

@@ -45,6 +45,14 @@ Where `PLATFORM` is the name of the satellite platform (currently `S1`, `CAPELLA
 
 Output RTC pixel values represent gamma0 power.
 
+To create an image that is geocoded but not radiometricly corrected, use the `geocoded` flag instead:
+
+```bash
+multirtc geocode PLATFORM SLC-GRANULE --resolution RESOLUTION --work-dir WORK-DIR
+```
+
+Output geocoded pixel values represent sigma0 power.
+
 ### Running via Docker
 In addition to the main python interface, I've also provided an experimental docker container that contains full support for polar grid format SICD data. Encapsulating this functionality in a docker container is ncessary for now because it requires re-compiling a development version of ISCE3. The docker container can be run using a similar interface, with exception of needing to pass your EarthData credentials and the need to pass a mounted volume with an `input` and `output` directory inside:
 
@@ -78,8 +86,10 @@ MultiRTC outputs one main RTC image and seven metadata images as GeoTIFFs. All l
 
 More information on the metadata images can be found in the OPERA RTC Static Product guide on the [OPERA RTC Product website](https://www.jpl.nasa.gov/go/opera/products/rtc-product/).
 
+All metadata images other than `FILEID_mask.tif`, and `FILEID_number_of_looks.tif` are omitted for geocode-only products.
+
 ### DEM options
-Currently, only the OPERA DEM is supported. This is a global Height Above Ellipsoid DEM sourced from the [COP-30 DEM](https://portal.opentopography.org/raster?opentopoID=OTSDEM.032021.4326.3). In the future, we hope to support a wider variety of automatically retrieved and user provided DEMs.
+Currently, only the OPERA DEM is supported. This is a global Height Above Ellipsoid DEM sourced from the [COP-30 DEM](https://portal.opentopography.org/raster?opentopoID=OTSDEM.032021.4326.3). In the future, we hope to support a wider variety of automatically retrieved and user provided DEMs. If the low resolution of the default DEM is causing radiometry issues, try using the `geocode` instead of `rtc` workflow.
 
 ## Calibration & Validation Subcommands
 

src/multirtc/__main__.py

@@ -1,6 +1,6 @@
 import argparse
 
-from multirtc import multirtc
+from multirtc import geocode, multirtc
 from multirtc.multimetric import ale, point_target, rle
 
 
@@ -12,8 +12,11 @@ def main():
     )
     subparsers = global_parser.add_subparsers(title='command', help='MultiRTC sub-commands')
 
-    multirtc_parser = multirtc.create_parser(subparsers.add_parser('rtc', help=multirtc.__doc__))
-    multirtc_parser.set_defaults(func=multirtc.run)
+    rtc_parser = multirtc.create_parser(subparsers.add_parser('rtc', help=multirtc.__doc__))
+    rtc_parser.set_defaults(func=multirtc.run)
+
+    geocode_parser = multirtc.create_parser(subparsers.add_parser('geocode', help=geocode.__doc__))
+    geocode_parser.set_defaults(func=geocode.run)
 
     ale_parser = ale.create_parser(subparsers.add_parser('ale', help=ale.__doc__))
     ale_parser.set_defaults(func=ale.run)

src/multirtc/create_rtc.py

@@ -518,11 +518,16 @@ def rtc(slc, geogrid, opts):
         geocode_kwargs['input_layover_shadow_mask_raster'] = slantrange_layover_shadow_mask_raster
 
     out_geo_nlooks_obj = isce3.io.Raster(nlooks_file, geogrid.width, geogrid.length, 1, gdal.GDT_Float32, raster_format)
-    out_geo_rtc_obj = isce3.io.Raster(rtc_anf_file, geogrid.width, geogrid.length, 1, gdal.GDT_Float32, raster_format)
-    out_geo_rtc_gamma0_to_sigma0_obj = isce3.io.Raster(
-        rtc_anf_gamma0_to_sigma0_file, geogrid.width, geogrid.length, 1, gdal.GDT_Float32, raster_format
-    )
-    geocode_kwargs['out_geo_rtc_gamma0_to_sigma0'] = out_geo_rtc_gamma0_to_sigma0_obj
+    out_geo_rtc_obj = None
+    if opts.apply_rtc:
+        out_geo_rtc_obj = isce3.io.Raster(
+            rtc_anf_file, geogrid.width, geogrid.length, 1, gdal.GDT_Float32, raster_format
+        )
+        out_geo_rtc_gamma0_to_sigma0_obj = isce3.io.Raster(
+            rtc_anf_gamma0_to_sigma0_file, geogrid.width, geogrid.length, 1, gdal.GDT_Float32, raster_format
+        )
+        geocode_kwargs['out_geo_rtc_gamma0_to_sigma0'] = out_geo_rtc_gamma0_to_sigma0_obj
+
     if opts.apply_bistatic_delay or opts.apply_static_tropo:
         rg_lut, az_lut = compute_correction_lut(
             slc.source,
@@ -611,26 +616,27 @@ def rtc(slc, geogrid, opts):
     out_geo_nlooks_obj.close_dataset()
     del out_geo_nlooks_obj
 
-    out_geo_rtc_obj.close_dataset()
-    del out_geo_rtc_obj
+    if opts.apply_rtc:
+        out_geo_rtc_gamma0_to_sigma0_obj.close_dataset()
+        del out_geo_rtc_gamma0_to_sigma0_obj
 
-    out_geo_rtc_gamma0_to_sigma0_obj.close_dataset()
-    del out_geo_rtc_gamma0_to_sigma0_obj
+        out_geo_rtc_obj.close_dataset()
+        del out_geo_rtc_obj
 
-    radar_grid_file_dict = {}
-    save_intermediate_geocode_files(
-        geogrid,
-        opts.dem_interpolation_method_isce3,
-        product_id,
-        output_dir,
-        raster_extension,
-        dem_raster,
-        radar_grid_file_dict,
-        lookside,
-        wavelength,
-        orbit,
-        doppler=doppler,
-    )
+        radar_grid_file_dict = {}
+        save_intermediate_geocode_files(
+            geogrid,
+            opts.dem_interpolation_method_isce3,
+            product_id,
+            output_dir,
+            raster_extension,
+            dem_raster,
+            radar_grid_file_dict,
+            lookside,
+            wavelength,
+            orbit,
+            doppler=doppler,
+        )
     t_end = time.time()
     logger.info(f'elapsed time: {t_end - t_start}')
 

src/multirtc/geocode.py

@@ -0,0 +1,19 @@
+"""Create a geocoded dataset for a multiple satellite platforms"""
+
+from pathlib import Path
+
+from multirtc.multirtc import SUPPORTED, run_multirtc
+
+
+def create_parser(parser):
+    parser.add_argument('platform', choices=SUPPORTED, help='Platform to create geocoded dataset for')
+    parser.add_argument('granule', help='Data granule to create geocoded dataset for.')
+    parser.add_argument('--resolution', type=float, help='Resolution of the output dataset (m)')
+    parser.add_argument('--work-dir', type=Path, default=None, help='Working directory for processing')
+    return parser
+
+
+def run(args):
+    if args.work_dir is None:
+        args.work_dir = Path.cwd()
+    run_multirtc(args.platform, args.granule, args.resolution, args.work_dir, apply_rtc=False)

src/multirtc/multimetric/ale.py

@@ -1,6 +1,5 @@
 """Absolute Location Error (ALE) analysis"""
 
-from argparse import ArgumentParser
 from datetime import datetime
 from pathlib import Path
 
@@ -243,14 +242,3 @@ def run(args):
     assert args.filepath.exists(), f'File {args.filepath} does not exist.'
 
     ale(args.filepath, args.date, args.azmangle, args.project, basedir=args.basedir)
-
-
-def main():
-    parser = ArgumentParser(description=__doc__)
-    parser = create_parser(parser)
-    args = parser.parse_args()
-    run(args)
-
-
-if __name__ == '__main__':
-    main()

src/multirtc/multimetric/point_target.py

@@ -1,6 +1,5 @@
 """Point target (resolution, PSLR, and ISLR) analysis"""
 
-from argparse import ArgumentParser
 from pathlib import Path
 
 import matplotlib.pyplot as plt
@@ -111,14 +110,3 @@ def run(args):
     assert args.filepath.exists()
     args.basedir = Path(args.basedir).expanduser()
     analyze_point_targets(args.platform, args.filepath, args.project, basedir=args.basedir)
-
-
-def main():
-    parser = ArgumentParser(description='Point target (resolution, PSLR, and ISLR) analysis')
-    parser = create_parser(parser)
-    args = parser.parse_args()
-    run(args)
-
-
-if __name__ == '__main__':
-    main()

src/multirtc/multimetric/rle.py

@@ -1,6 +1,5 @@
 """Relative Location Error (RLE) analysis"""
 
-from argparse import ArgumentParser
 from dataclasses import dataclass
 from pathlib import Path
 from tempfile import NamedTemporaryFile
@@ -196,14 +195,3 @@ def run(args):
     args.basedir = Path(args.basedir)
     assert args.basedir.exists(), f'Base directory {args.basedir} does not exist'
     rle(args.reference, args.secondary, project=args.project, basedir=args.basedir)
-
-
-def main():
-    parser = ArgumentParser(description=__doc__)
-    parser = create_parser(parser)
-    args = parser.parse_args()
-    run(args)
-
-
-if __name__ == '__main__':
-    main()

src/multirtc/multirtc.py

@@ -1,6 +1,5 @@
 """Create an RTC dataset for a multiple satellite platforms"""
 
-import argparse
 from pathlib import Path
 
 from burst2safe.burst2safe import burst2safe
@@ -61,14 +60,15 @@ def get_slc(platform: str, granule: str, input_dir: Path) -> Slc:
     return slc
 
 
-def run_multirtc(platform: str, granule: str, resolution: int, work_dir: Path) -> None:
+def run_multirtc(platform: str, granule: str, resolution: int, work_dir: Path, apply_rtc=True) -> None:
     """Create an RTC or Geocoded dataset using the OPERA algorithm.
 
     Args:
         platform: Platform type (e.g., 'UMBRA').
         granule: Granule name if data is available in ASF archive, or filename if granule is already downloaded.
         resolution: Resolution of the output RTC (in meters).
         work_dir: Working directory for processing.
+        apply_rtc: If True perform radiometric correction; if False, only geocode.
     """
     input_dir, output_dir = prep_dirs(work_dir)
     slc = get_slc(platform, granule, input_dir)
@@ -79,6 +79,7 @@ def run_multirtc(platform: str, granule: str, resolution: int, work_dir: Path) -
         opts = RtcOptions(
             dem_path=str(dem_path),
             output_dir=str(output_dir),
+            apply_rtc=apply_rtc,
             resolution=resolution,
             apply_bistatic_delay=slc.supports_bistatic_delay,
             apply_static_tropo=slc.supports_static_tropo,
@@ -102,20 +103,4 @@ def create_parser(parser):
 def run(args):
     if args.work_dir is None:
         args.work_dir = Path.cwd()
-    run_multirtc(args.platform, args.granule, args.resolution, args.work_dir)
-
-
-def main():
-    """Create a RTC or geocoded dataset for a multiple satellite platforms
-
-    Example command:
-    multirtc UMBRA umbra_image.ntif --resolution 40
-    """
-    parser = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.ArgumentDefaultsHelpFormatter)
-    parser = create_parser(parser)
-    args = parser.parse_args()
-    run(args)
-
-
-if __name__ == '__main__':
-    main()
+    run_multirtc(args.platform, args.granule, args.resolution, args.work_dir, apply_rtc=True)

src/multirtc/rtc_options.py

@@ -10,11 +10,13 @@ class RtcOptions:
 
     output_dir: str
     dem_path: str
-    apply_rtc: bool = True
+    apply_rtc: bool
+    resolution: float
     apply_thermal_noise: bool = True
     apply_abs_rad: bool = True
     apply_bistatic_delay: bool = True
     apply_static_tropo: bool = True
+    terrain_radiometry: str = 'gamma0'  # 'gamma0' or 'sigma0'
     apply_valid_samples_sub_swath_masking: bool = True
     apply_shadow_masking: bool = True
     dem_interpolation_method: str = 'biquintic'
@@ -28,33 +30,22 @@ class RtcOptions:
     clip_min: float = np.nan
     clip_max: float = np.nan
     upsample_radar_grid: bool = False
-    terrain_radiometry: str = 'gamma0'  # 'gamma0' or 'sigma0'
     rtc_algorithm_type: str = 'area_projection'  # 'area_projection' or 'bilinear_distribution'
     input_terrain_radiometry: str = 'beta0'
-    rtc_min_value_db: int = -30.0
+    rtc_min_value_db: float = -30.0
     rtc_upsampling: int = 2
     rtc_area_beta_mode: str = 'auto'
     geo2rdr_threshold: float = 1.0e-7
     geo2rdr_numiter: int = 50
     rdr2geo_threshold: float = 1.0e-7
     rdr2geo_numiter: int = 25
-    output_epsg: int = None
-    resolution: int = 30
+    output_epsg: int | None = None
 
     def __post_init__(self):
         if not self.apply_rtc:
-            if self.save_rtc_anf:
-                raise ValueError('RTC ANF flags are only available with RTC enabled')
-            if self.save_rtc_anf_gamma0_to_sigma0:
-                raise ValueError('RTC ANF gamma0 to sigma0 flags are only available with RTC enabled')
-
-        if self.terrain_radiometry == 'sigma0' and self.save_rtc_anf_gamma0_to_sigma0:
-            raise ValueError('RTC ANF gamma0 to sigma0 flags are only available with output type set to gamma0')
+            self.terrain_radiometry: str = 'sigma0'
 
-        if self.apply_rtc:
-            self.layer_name_rtc_anf = f'rtc_anf_{self.terrain_radiometry}_to_{self.input_terrain_radiometry}'
-        else:
-            self.layer_name_rtc_anf = ''
+        self.layer_name_rtc_anf = f'rtc_anf_{self.terrain_radiometry}_to_{self.input_terrain_radiometry}'
 
         if self.dem_interpolation_method == 'biquintic':
             self.dem_interpolation_method_isce3 = isce3.core.DataInterpMethod.BIQUINTIC
@@ -94,11 +85,6 @@ def __post_init__(self):
         else:
             raise ValueError(f'Invalid terrain radiometry: {self.terrain_radiometry}')
 
-        if self.apply_rtc:
-            self.layer_name_rtc_anf = f'rtc_anf_{self.terrain_radiometry}_to_{self.input_terrain_radiometry}'
-        else:
-            self.layer_name_rtc_anf = ''
-
         if self.rtc_area_beta_mode == 'pixel_area':
             self.rtc_area_beta_mode_isce3 = isce3.geometry.RtcAreaBetaMode.PIXEL_AREA
         elif self.rtc_area_beta_mode == 'projection_angle':