#!/usr/bin/env python3 ############################################################ # Program is part of MintPy # # Copyright (c) 2013, Zhang Yunjun, Heresh Fattahi # # Author: Antonio Valentino, Zhang Yunjun, Aug 2022 # ############################################################ import math import os import sys from mintpy.defaults.template import get_template_content from mintpy.utils import arg_utils ###################################################################################### TEMPLATE = get_template_content('geocode') EXAMPLE = """example: geocode.py velocity.h5 geocode.py velocity.h5 -b -0.5 -0.25 -91.3 -91.1 geocode.py velocity.h5 timeseries.h5 -t smallbaselineApp.cfg --outdir ./geo --update # geocode file using ISCE-2 lat/lon.rdr file geocode.py filt_fine.int --lat-file ../../geom_reference/lat.rdr --lon-file ../../geom_reference/lon.rdr # radar-code file in geo coordinates geocode.py swbdLat_S02_N01_Lon_W092_W090.wbd -l geometryRadar.h5 -o waterMask.rdr --geo2radar geocode.py geo_velocity.h5 --geo2radar """ DEG2METER = """ degrees --> meters on equator 0.000925926 --> 100 0.000833334 --> 90 0.000555556 --> 60 0.000462963 --> 50 0.000370370 --> 40 0.000277778 --> 30 0.000185185 --> 20 0.000092593 --> 10 """ def create_parser(subparsers=None): synopsis = 'Resample radar-coded files into geo-coordinates or vice versa' epilog = TEMPLATE + '\n' + EXAMPLE name = __name__.split('.')[-1] parser = arg_utils.create_argument_parser( name, synopsis=synopsis, description=synopsis, epilog=epilog, subparsers=subparsers) parser.add_argument('file', nargs='+', help='File(s) to be geocoded') parser.add_argument('-d', '--dset', help='dataset to be geocoded, for example:\n' + 'height for geometryRadar.h5\n' + 'unwrapPhase-20100114_20101017 for ifgramStack.h5') parser.add_argument('-l', '--lookup', dest='lookupFile', help='Lookup table file generated by InSAR processors.') parser.add_argument('--lat-file', dest='latFile', help='lookup table file for latitude.') parser.add_argument('--lon-file', dest='lonFile', help='lookup table file for longitude.') parser.add_argument('--geo2radar', '--geo2rdr', dest='radar2geo', action='store_false', help='resample geocoded files into radar coordinates.\n' + 'ONLY for lookup table in radar-coord (ISCE, Doris).') parser.add_argument('-t', '--template', dest='templateFile', help="Template file with geocoding options.") # output grid / geometry out = parser.add_argument_group('grid in geo-coordinates') out.add_argument('-b', '--bbox', dest='SNWE', type=float, nargs=4, metavar=('S', 'N', 'W', 'E'), help='Bounding box for the area of interest.\n' 'using coordinates of the upper left corner of the first pixel\n' ' and the lower right corner of the last pixel\n' "for radar2geo, it's the output spatial extent\n" "for geo2radar, it's the input spatial extent") out.add_argument('--lalo', '--lalo-step', dest='laloStep', type=float, nargs=2, metavar=('LAT_STEP', 'LON_STEP'), help=f'output pixel size in degree in latitude / longitude.{DEG2METER}') # interpolation / resampling interp = parser.add_argument_group('interpolation') interp.add_argument('-i', '--interp', dest='interpMethod', default='nearest', choices={'nearest', 'linear'}, help='interpolation/resampling method (default: %(default)s).') interp.add_argument('--fill', dest='fillValue', type=float, default=math.nan, help='Fill value for extrapolation (default: %(default)s or 0 for *.int/unw files).') interp.add_argument('-n','--nprocs', dest='nprocs', type=int, default=1, help='number of processors to be used for calculation (default: %(default)s).\n' 'Note: Do not use more processes than available processor cores.') interp.add_argument('--software', dest='software', default='pyresample', choices={'pyresample', 'scipy'}, help='software/module used for interpolation (default: %(default)s)\n' 'Note: --bbox is not supported for -p scipy') parser.add_argument('--update', dest='updateMode', action='store_true', help='skip resampling if output file exists and newer than input file') parser.add_argument('-o', '--output', dest='outfile', help="output file name. Default: add prefix 'geo_'") parser.add_argument('--outdir', '--output-dir', dest='out_dir', help='output directory.') # computing parser = arg_utils.add_memory_argument(parser) return parser def cmd_line_parse(iargs=None): # parse parser = create_parser() inps = parser.parse_args(args=iargs) # import from mintpy.utils import readfile, utils as ut # save argv (to check the manually specified arguments) # use iargs for python call # use sys.argv[1:] for command line call inps.argv = iargs if iargs else sys.argv[1:] # check if inps.templateFile: inps = read_template2inps(inps.templateFile, inps) # check: input file(s) existence inps.file = ut.get_file_list(inps.file) if not inps.file: raise Exception('ERROR: no input file found!') elif len(inps.file) > 1: inps.outfile = None # check: --lookup (lookup table existence) if not inps.lookupFile: # use isce-2 lat/lon.rdr file if not inps.lookupFile and inps.latFile: inps.lookupFile = inps.latFile # grab default lookup table inps.lookupFile = ut.get_lookup_file(inps.lookupFile) # final check if not inps.lookupFile: raise FileNotFoundError('No lookup table found! Can not geocode without it.') # check: src file coordinate & radar2geo operatioin atr = readfile.read_attribute(inps.file[0]) if 'Y_FIRST' in atr.keys() and inps.radar2geo: print('input file is already geocoded') print('to resample geocoded files into radar coordinates, use --geo2radar option') print('exit without doing anything.') sys.exit(0) elif 'Y_FIRST' not in atr.keys() and not inps.radar2geo: print('input file is already in radar coordinates, exit without doing anything') sys.exit(0) # check: --lalo-step # valid only if: # 1. radar2geo = True AND # 2. lookup table is in radar coordinates if inps.laloStep: if not inps.radar2geo: print('ERROR: "--lalo-step" can NOT be used together with "--geo2radar"!') sys.exit(0) atr = readfile.read_attribute(inps.lookupFile) if 'Y_FIRST' in atr.keys(): print('ERROR: "--lalo-step" can NOT be used with lookup table file in geo-coordinates!') sys.exit(0) # check: --nprocs (number of processors for multiprocessing) inps.nprocs = check_num_processor(inps.nprocs) # check: --geo2radar if not inps.radar2geo: if inps.SNWE: print('ERROR: "--geo2radar" can NOT be used together with "--bbox"!') sys.exit(0) if inps.software == 'scipy': print('ERROR: "--geo2radar" is NOT supported for "--software scipy"!') sys.exit(0) # default: --fill (set default to zero for .int/unw file) fext = os.path.splitext(inps.file[0])[1] if '--fill' not in inps.argv and fext in ['.int', '.unw']: inps.fillValue = 0 return inps def read_template2inps(template_file, inps): """Read input template options into Namespace inps""" print('read input option from template file:', template_file) from mintpy.utils import readfile, utils as ut iDict = vars(inps) template = readfile.read_template(template_file, skip_chars=['[', ']']) template = ut.check_template_auto_value(template) key_prefix = 'mintpy.geocode.' key_list = [i for i in list(iDict.keys()) if key_prefix + i in template.keys()] for key in key_list: value = template[key_prefix + key] if value: if key in ['SNWE', 'laloStep']: iDict[key] = [float(i) for i in value.split(',')] elif key in ['interpMethod']: iDict[key] = value elif key == 'fillValue': if 'nan' in value.lower(): iDict[key] = math.nan else: iDict[key] = float(value) # ensure laloStep is a list of two items key = 'laloStep' if key in iDict.keys() and iDict[key]: if len(iDict[key]) == 1: lalo_step = iDict[key] iDict[key] = [-1 * abs(lalo_step[0]), abs(lalo_step[0])] print(f'single laloStep input {lalo_step} detected, convert into two as {iDict[key]}') elif len(iDict[key]) > 2: raise ValueError(f'laloStep input {iDict[key]} could NOT have >2 items!') # computing configurations key = 'mintpy.compute.maxMemory' if key in template.keys() and template[key]: iDict['maxMemory'] = float(template[key]) return inps def check_num_processor(nprocs): """Check number of processors Note by Yunjun, 2019-05-02: 1. conda install pyresample will install pykdtree and openmp, but it seems not working: geocode.py is getting slower with more processors Test on a TS HDF5 file in size of (241, 2267, 2390) Memory: up to 10GB Run time: 2.5 mins for nproc=1, 3 mins for nproc=4 2. macports seems to have minor speedup when more processors Thus, default number of processors is set to 1; although the capability of using multiple processors is written here. """ if not nprocs: # OMP_NUM_THREADS is defined in environment variable for OpenMP if 'OMP_NUM_THREADS' in os.environ: nprocs = int(os.getenv('OMP_NUM_THREADS')) else: nprocs = int(os.cpu_count() / 2) nprocs = min(os.cpu_count(), nprocs) print(f'number of processor to be used: {nprocs}') return nprocs ###################################################################################### def main(iargs=None): # parse inps = cmd_line_parse(iargs) # import from mintpy.geocode import run_geocode # run run_geocode(inps) ###################################################################################### if __name__ == '__main__': main(sys.argv[1:])