diff --git a/src/hyp3_isce2/insar_stripmap.py b/src/hyp3_isce2/insar_stripmap.py index b0cc8412..86015c3b 100644 --- a/src/hyp3_isce2/insar_stripmap.py +++ b/src/hyp3_isce2/insar_stripmap.py @@ -23,7 +23,7 @@ log = logging.getLogger(__name__) -def insar_stripmap(user: str, password: str, reference_scene: str, secondary_scene: str) -> Path: +def insar_stripmap(reference_scene: str, secondary_scene: str) -> Path: """Create a Stripmap interferogram Args: @@ -35,12 +35,22 @@ def insar_stripmap(user: str, password: str, reference_scene: str, secondary_sce Returns: Path to the output files """ - session = asf_search.ASFSession().auth_with_creds(user, password) - - reference_product, secondary_product = asf_search.search( + scenes = sorted([reference_scene, secondary_scene]) + print(scenes) + reference_scene = scenes[0] + secondary_scene = scenes[1] + products = asf_search.search( granule_list=[reference_scene, secondary_scene], - processingLevel=asf_search.L1_0, + processingLevel="L1.0", ) + + if products[0].properties['sceneName']==reference_scene: + reference_product = products[0] + secondary_product = products[1] + else: + reference_product = products[1] + secondary_product = products[0] + assert reference_product.properties['sceneName'] == reference_scene assert secondary_product.properties['sceneName'] == secondary_scene products = (reference_product, secondary_product) @@ -51,7 +61,7 @@ def insar_stripmap(user: str, password: str, reference_scene: str, secondary_sce dem_path = download_dem_for_isce2(insar_roi, dem_name='glo_30', dem_dir=Path('dem'), buffer=0) urls = [product.properties['url'] for product in products] - asf_search.download_urls(urls=urls, path=os.getcwd(), session=session, processes=2) + asf_search.download_urls(urls=urls, path=os.getcwd(), processes=2) zip_paths = [product.properties['fileName'] for product in products] for zip_path in zip_paths: @@ -93,7 +103,7 @@ def insar_stripmap(user: str, password: str, reference_scene: str, secondary_sce def get_product_file(product: asf_search.ASFProduct, file_prefix: str) -> str: paths = glob.glob(str(Path(product.properties['fileID']) / f'{file_prefix}*')) - assert len(paths) == 1 + assert len(paths) > 0 return paths[0] @@ -104,10 +114,8 @@ def main(): parser.add_argument('--bucket', help='AWS S3 bucket HyP3 for upload the final product(s)') parser.add_argument('--bucket-prefix', default='', help='Add a bucket prefix to product(s)') - parser.add_argument('--username', type=str, required=True) - parser.add_argument('--password', type=str, required=True) - parser.add_argument('--reference-scene', type=str, required=True) - parser.add_argument('--secondary-scene', type=str, required=True) + parser.add_argument('--reference', type=str, required=True) + parser.add_argument('--secondary', type=str, required=True) args = parser.parse_args() @@ -117,10 +125,8 @@ def main(): log.info('Begin InSAR Stripmap run') product_dir = insar_stripmap( - user=args.username, - password=args.password, - reference_scene=args.reference_scene, - secondary_scene=args.secondary_scene, + reference_scene=args.reference, + secondary_scene=args.secondary, ) log.info('InSAR Stripmap run completed successfully') diff --git a/src/hyp3_isce2/insar_tops.py b/src/hyp3_isce2/insar_tops.py index 0e1fe26e..dd836e4b 100644 --- a/src/hyp3_isce2/insar_tops.py +++ b/src/hyp3_isce2/insar_tops.py @@ -2,6 +2,7 @@ import argparse import logging +import os import sys from pathlib import Path from shutil import copyfile, make_archive @@ -13,6 +14,7 @@ from hyp3_isce2.dem import download_dem_for_isce2 from hyp3_isce2.logger import configure_root_logger from hyp3_isce2.s1_auxcal import download_aux_cal +from hyp3_isce2.utils import make_browse_image log = logging.getLogger(__name__) @@ -114,12 +116,16 @@ def insar_tops_packaged( Path to the output files """ pixel_size = packaging.get_pixel_size(f'{range_looks}x{azimuth_looks}') - product_name = packaging.get_product_name(reference, secondary, pixel_spacing=int(pixel_size)) log.info('Begin ISCE2 TopsApp run') - insar_tops(reference, secondary, download=False) + if os.path.exists(f'{reference}.SAFE') and os.path.exists(f'{reference}.SAFE'): + insar_tops(reference, secondary, download=False) + else: + insar_tops(reference, secondary) log.info('ISCE2 TopsApp run completed successfully') + product_name = packaging.get_product_name(reference, secondary, pixel_spacing=int(pixel_size)) + product_dir = Path(product_name) product_dir.mkdir(parents=True, exist_ok=True) @@ -129,7 +135,7 @@ def insar_tops_packaged( if apply_water_mask: packaging.water_mask(unwrapped_phase, f'{product_name}/{product_name}_water_mask.tif') - packaging.make_browse_image(unwrapped_phase, f'{product_name}/{product_name}_unw_phase.png') + make_browse_image(unwrapped_phase, f'{product_name}/{product_name}_unw_phase.png') packaging.make_readme( product_dir=product_dir, product_name=product_name, @@ -157,7 +163,7 @@ def main(): parser = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.ArgumentDefaultsHelpFormatter) parser.add_argument('--reference', type=str, help='Reference granule') parser.add_argument('--secondary', type=str, help='Secondary granule') - parser.add_argument('--polarization', type=str, defualt='VV', help='Polarization to use') + parser.add_argument('--polarization', type=str, default='VV', help='Polarization to use') parser.add_argument( '--looks', choices=['20x4', '10x2', '5x1'], default='20x4', help='Number of looks to take in range and azimuth' ) @@ -179,8 +185,8 @@ def main(): raise ValueError('Polarization must be one of VV, VH, HV, or HH') insar_tops_packaged( - reference_scene=args.reference, - secondary_scene=args.secondary, + reference=args.reference, + secondary=args.secondary, polarization=args.polarization, azimuth_looks=azimuth_looks, range_looks=range_looks, diff --git a/src/hyp3_isce2/insar_tops_burst.py b/src/hyp3_isce2/insar_tops_burst.py index 21ee2e1e..544bf645 100644 --- a/src/hyp3_isce2/insar_tops_burst.py +++ b/src/hyp3_isce2/insar_tops_burst.py @@ -30,6 +30,7 @@ from hyp3_isce2.utils import ( image_math, isce2_copy, + make_browse_image, oldest_granule_first, resample_to_radar_io, ) @@ -176,7 +177,7 @@ def insar_tops_single_burst( if apply_water_mask: packaging.water_mask(unwrapped_phase, f'{product_name}/{product_name}_water_mask.tif') - packaging.make_browse_image(unwrapped_phase, f'{product_name}/{product_name}_unw_phase.png') + make_browse_image(unwrapped_phase, f'{product_name}/{product_name}_unw_phase.png') packaging.make_readme( product_dir=product_dir, diff --git a/src/hyp3_isce2/merge_tops_bursts.py b/src/hyp3_isce2/merge_tops_bursts.py index b3c1f667..4ccc327a 100644 --- a/src/hyp3_isce2/merge_tops_bursts.py +++ b/src/hyp3_isce2/merge_tops_bursts.py @@ -39,6 +39,7 @@ import hyp3_isce2 import hyp3_isce2.burst as burst_utils from hyp3_isce2.dem import download_dem_for_isce2 +from hyp3_isce2.packaging import get_pixel_size, translate_outputs from hyp3_isce2.utils import ( ParameterFile, create_image, @@ -56,8 +57,6 @@ logging.basicConfig(stream=sys.stdout, format=log_format, level=logging.INFO, force=True) log = logging.getLogger(__name__) -from hyp3_isce2.insar_tops_burst import get_pixel_size, translate_outputs # noqa - BURST_IFG_DIR = 'fine_interferogram' BURST_GEOM_DIR = 'geom_reference' diff --git a/src/hyp3_isce2/packaging.py b/src/hyp3_isce2/packaging.py index dc0fcf73..2b79f324 100644 --- a/src/hyp3_isce2/packaging.py +++ b/src/hyp3_isce2/packaging.py @@ -19,7 +19,7 @@ import hyp3_isce2.metadata.util from hyp3_isce2.burst import BurstPosition from hyp3_isce2.slc import get_geometry_from_manifest -from hyp3_isce2.utils import get_projection, utm_from_lon_lat +from hyp3_isce2.utils import get_projection, ParameterFile, utm_from_lon_lat @dataclass @@ -284,43 +284,6 @@ def water_mask(unwrapped_phase: str, water_mask: str) -> None: subprocess.run(cmd.split(' '), check=True) -class GDALConfigManager: - """Context manager for setting GDAL config options temporarily""" - - def __init__(self, **options): - """ - Args: - **options: GDAL Config `option=value` keyword arguments. - """ - self.options = options.copy() - self._previous_options = {} - - def __enter__(self): - for key in self.options: - self._previous_options[key] = gdal.GetConfigOption(key) - - for key, value in self.options.items(): - gdal.SetConfigOption(key, value) - - def __exit__(self, exc_type, exc_val, exc_tb): - for key, value in self._previous_options.items(): - gdal.SetConfigOption(key, value) - - -def make_browse_image(input_tif: str, output_png: str) -> None: - with GDALConfigManager(GDAL_PAM_ENABLED='NO'): - stats = gdal.Info(input_tif, format='json', stats=True)['stac']['raster:bands'][0]['stats'] - gdal.Translate( - destName=output_png, - srcDS=input_tif, - format='png', - outputType=gdal.GDT_Byte, - width=2048, - strict=True, - scaleParams=[[stats['minimum'], stats['maximum']]], - ) - - def make_readme( product_dir: Path, product_name: str, @@ -361,96 +324,6 @@ def make_readme( f.write(content) -@dataclass -class ParameterFile: - reference_granule: str - secondary_granule: str - reference_orbit_direction: str - reference_orbit_number: str - secondary_orbit_direction: str - secondary_orbit_number: str - baseline: float - utc_time: float - heading: float - spacecraft_height: float - earth_radius_at_nadir: float - slant_range_near: float - slant_range_center: float - slant_range_far: float - range_looks: int - azimuth_looks: int - insar_phase_filter: bool - phase_filter_parameter: float - range_bandpass_filter: bool - azimuth_bandpass_filter: bool - dem_source: str - dem_resolution: int - unwrapping_type: str - speckle_filter: bool - water_mask: bool - radar_n_lines: Optional[int] = None - radar_n_samples: Optional[int] = None - radar_first_valid_line: Optional[int] = None - radar_n_valid_lines: Optional[int] = None - radar_first_valid_sample: Optional[int] = None - radar_n_valid_samples: Optional[int] = None - multilook_azimuth_time_interval: Optional[float] = None - multilook_range_pixel_size: Optional[float] = None - radar_sensing_stop: Optional[datetime] = None - - def __str__(self): - output_strings = [ - f'Reference Granule: {self.reference_granule}\n', - f'Secondary Granule: {self.secondary_granule}\n', - f'Reference Pass Direction: {self.reference_orbit_direction}\n', - f'Reference Orbit Number: {self.reference_orbit_number}\n', - f'Secondary Pass Direction: {self.secondary_orbit_direction}\n', - f'Secondary Orbit Number: {self.secondary_orbit_number}\n', - f'Baseline: {self.baseline}\n', - f'UTC time: {self.utc_time}\n', - f'Heading: {self.heading}\n', - f'Spacecraft height: {self.spacecraft_height}\n', - f'Earth radius at nadir: {self.earth_radius_at_nadir}\n', - f'Slant range near: {self.slant_range_near}\n', - f'Slant range center: {self.slant_range_center}\n', - f'Slant range far: {self.slant_range_far}\n', - f'Range looks: {self.range_looks}\n', - f'Azimuth looks: {self.azimuth_looks}\n', - f'INSAR phase filter: {"yes" if self.insar_phase_filter else "no"}\n', - f'Phase filter parameter: {self.phase_filter_parameter}\n', - f'Range bandpass filter: {"yes" if self.range_bandpass_filter else "no"}\n', - f'Azimuth bandpass filter: {"yes" if self.azimuth_bandpass_filter else "no"}\n', - f'DEM source: {self.dem_source}\n', - f'DEM resolution (m): {self.dem_resolution}\n', - f'Unwrapping type: {self.unwrapping_type}\n', - f'Speckle filter: {"yes" if self.speckle_filter else "no"}\n', - f'Water mask: {"yes" if self.water_mask else "no"}\n', - ] - - # TODO could use a more robust way to check if radar data is present - if self.radar_n_lines: - radar_data = [ - f'Radar n lines: {self.radar_n_lines}\n', - f'Radar n samples: {self.radar_n_samples}\n', - f'Radar first valid line: {self.radar_first_valid_line}\n', - f'Radar n valid lines: {self.radar_n_valid_lines}\n', - f'Radar first valid sample: {self.radar_first_valid_sample}\n', - f'Radar n valid samples: {self.radar_n_valid_samples}\n', - f'Multilook azimuth time interval: {self.multilook_azimuth_time_interval}\n', - f'Multilook range pixel size: {self.multilook_range_pixel_size}\n', - f'Radar sensing stop: {datetime.strftime(self.radar_sensing_stop, "%Y-%m-%dT%H:%M:%S.%f")}\n', - ] - output_strings += radar_data - - return ''.join(output_strings) - - def __repr__(self): - return self.__str__() - - def write(self, out_path: Path): - out_path.write_text(self.__str__()) - - def find_available_swaths(base_dir: Path | str) -> list[str]: """Find the available swaths in the given directory diff --git a/src/hyp3_isce2/slc.py b/src/hyp3_isce2/slc.py index b26893f4..02be9468 100644 --- a/src/hyp3_isce2/slc.py +++ b/src/hyp3_isce2/slc.py @@ -33,7 +33,7 @@ def get_granule(granule: str) -> Path: def unzip_granule(zip_file: Path, remove: bool = False) -> Path: with ZipFile(zip_file) as z: z.extractall() - safe_dir = next(item.filename for item in z.infolist() if item.is_dir() and item.filename.endswith('.SAFE/')) + safe_dir = zip_file.split('.')[0]+'.SAFE/' if remove: os.remove(zip_file) return safe_dir.strip('/') diff --git a/src/hyp3_isce2/utils.py b/src/hyp3_isce2/utils.py index cce0af7c..b1d7d191 100644 --- a/src/hyp3_isce2/utils.py +++ b/src/hyp3_isce2/utils.py @@ -1,5 +1,8 @@ import shutil import subprocess +from dataclasses import dataclass +from datetime import datetime +from pathlib import Path from typing import Optional import isceobj @@ -12,6 +15,119 @@ gdal.UseExceptions() +class GDALConfigManager: + """Context manager for setting GDAL config options temporarily""" + + def __init__(self, **options): + """ + Args: + **options: GDAL Config `option=value` keyword arguments. + """ + self.options = options.copy() + self._previous_options = {} + + def __enter__(self): + for key in self.options: + self._previous_options[key] = gdal.GetConfigOption(key) + + for key, value in self.options.items(): + gdal.SetConfigOption(key, value) + + def __exit__(self, exc_type, exc_val, exc_tb): + for key, value in self._previous_options.items(): + gdal.SetConfigOption(key, value) + + +@dataclass +class ParameterFile: + reference_granule: str + secondary_granule: str + reference_orbit_direction: str + reference_orbit_number: str + secondary_orbit_direction: str + secondary_orbit_number: str + baseline: float + utc_time: float + heading: float + spacecraft_height: float + earth_radius_at_nadir: float + slant_range_near: float + slant_range_center: float + slant_range_far: float + range_looks: int + azimuth_looks: int + insar_phase_filter: bool + phase_filter_parameter: float + range_bandpass_filter: bool + azimuth_bandpass_filter: bool + dem_source: str + dem_resolution: int + unwrapping_type: str + speckle_filter: bool + water_mask: bool + radar_n_lines: Optional[int] = None + radar_n_samples: Optional[int] = None + radar_first_valid_line: Optional[int] = None + radar_n_valid_lines: Optional[int] = None + radar_first_valid_sample: Optional[int] = None + radar_n_valid_samples: Optional[int] = None + multilook_azimuth_time_interval: Optional[float] = None + multilook_range_pixel_size: Optional[float] = None + radar_sensing_stop: Optional[datetime] = None + + def __str__(self): + output_strings = [ + f'Reference Granule: {self.reference_granule}\n', + f'Secondary Granule: {self.secondary_granule}\n', + f'Reference Pass Direction: {self.reference_orbit_direction}\n', + f'Reference Orbit Number: {self.reference_orbit_number}\n', + f'Secondary Pass Direction: {self.secondary_orbit_direction}\n', + f'Secondary Orbit Number: {self.secondary_orbit_number}\n', + f'Baseline: {self.baseline}\n', + f'UTC time: {self.utc_time}\n', + f'Heading: {self.heading}\n', + f'Spacecraft height: {self.spacecraft_height}\n', + f'Earth radius at nadir: {self.earth_radius_at_nadir}\n', + f'Slant range near: {self.slant_range_near}\n', + f'Slant range center: {self.slant_range_center}\n', + f'Slant range far: {self.slant_range_far}\n', + f'Range looks: {self.range_looks}\n', + f'Azimuth looks: {self.azimuth_looks}\n', + f'INSAR phase filter: {"yes" if self.insar_phase_filter else "no"}\n', + f'Phase filter parameter: {self.phase_filter_parameter}\n', + f'Range bandpass filter: {"yes" if self.range_bandpass_filter else "no"}\n', + f'Azimuth bandpass filter: {"yes" if self.azimuth_bandpass_filter else "no"}\n', + f'DEM source: {self.dem_source}\n', + f'DEM resolution (m): {self.dem_resolution}\n', + f'Unwrapping type: {self.unwrapping_type}\n', + f'Speckle filter: {"yes" if self.speckle_filter else "no"}\n', + f'Water mask: {"yes" if self.water_mask else "no"}\n', + ] + + # TODO could use a more robust way to check if radar data is present + if self.radar_n_lines: + radar_data = [ + f'Radar n lines: {self.radar_n_lines}\n', + f'Radar n samples: {self.radar_n_samples}\n', + f'Radar first valid line: {self.radar_first_valid_line}\n', + f'Radar n valid lines: {self.radar_n_valid_lines}\n', + f'Radar first valid sample: {self.radar_first_valid_sample}\n', + f'Radar n valid samples: {self.radar_n_valid_samples}\n', + f'Multilook azimuth time interval: {self.multilook_azimuth_time_interval}\n', + f'Multilook range pixel size: {self.multilook_range_pixel_size}\n', + f'Radar sensing stop: {datetime.strftime(self.radar_sensing_stop, "%Y-%m-%dT%H:%M:%S.%f")}\n', + ] + output_strings += radar_data + + return ''.join(output_strings) + + def __repr__(self): + return self.__str__() + + def write(self, out_path: Path): + out_path.write_text(self.__str__()) + + def utm_from_lon_lat(lon: float, lat: float) -> int: """Get the UTM zone EPSG code from a longitude and latitude. See https://en.wikipedia.org/wiki/Universal_Transverse_Mercator_coordinate_system @@ -49,6 +165,20 @@ def extent_from_geotransform(geotransform: tuple, x_size: int, y_size: int) -> t return extent +def make_browse_image(input_tif: str, output_png: str) -> None: + with GDALConfigManager(GDAL_PAM_ENABLED='NO'): + stats = gdal.Info(input_tif, format='json', stats=True)['stac']['raster:bands'][0]['stats'] + gdal.Translate( + destName=output_png, + srcDS=input_tif, + format='png', + outputType=gdal.GDT_Byte, + width=2048, + strict=True, + scaleParams=[[stats['minimum'], stats['maximum']]], + ) + + def oldest_granule_first(g1, g2): if g1[14:29] <= g2[14:29]: return g1, g2