From 8127ecf6c1db12e42cbe9d7ae78773a469cffe8d Mon Sep 17 00:00:00 2001 From: rexionmars Date: Tue, 24 Oct 2023 10:56:11 -0300 Subject: [PATCH] Geogrid update for post-processing correction --- testGeogrid_ISCE.py | 83 ++++++++++++++++++++++++++++++++++++++++++--- 1 file changed, 79 insertions(+), 4 deletions(-) diff --git a/testGeogrid_ISCE.py b/testGeogrid_ISCE.py index 22ecdcf..0617d22 100755 --- a/testGeogrid_ISCE.py +++ b/testGeogrid_ISCE.py @@ -67,6 +67,12 @@ def cmdLineParse(): help='Input stable surface mask') parser.add_argument('-fo', '--flag_optical', dest='optical_flag', type=bool, required=False, default=0, help='flag for reading optical data (e.g. Landsat): use 1 for on and 0 (default) for off') + parser.add_argument('-b', '--buffer', dest='buffer', type=bool, required=False, default=0, + help='buffer to add to the starting/end range accounting for all passes from the same relative orbit') + parser.add_argument('-p', '--parse', dest='parse', action='store_true', + default=False, help='Parse the SAFE zip file to get radar image and orbit metadata; no need to run ISCE') + parser.add_argument('--orbit-dir', hep='Directory of Sentinel-1 orbit files') + parser.add_argument('--aux-dir', hep='Directory of Sentinel-1 aux files') return parser.parse_args() @@ -113,7 +119,7 @@ def getMergedOrbit(product): return orb -def loadMetadata(indir): +def loadMetadata(indir,buffer=0): ''' Input file. ''' @@ -135,13 +141,78 @@ def loadMetadata(indir): info.prf = 1.0 / frames[0].bursts[0].azimuthTimeInterval info.rangePixelSize = frames[0].bursts[0].rangePixelSize info.lookSide = -1 + + info.startingRange -= buffer * info.rangePixelSize + info.farRange += buffer * info.rangePixelSize + info.numberOfLines = int( np.round( (info.sensingStop - info.sensingStart).total_seconds() * info.prf)) + 1 - info.numberOfSamples = int( np.round( (info.farRange - info.startingRange)/info.rangePixelSize)) + 1 + info.numberOfSamples = int( np.round( (info.farRange - info.startingRange)/info.rangePixelSize)) + 1 + 2 * buffer info.orbit = getMergedOrbit(frames) return info +def get_polarizations(s1_safe: str): + ''' + Return: + Tuple[str] + ''' + from typing import Tuple + from pathlib import Path + + mapping = { + 'SH': ('hh',), + 'SV': ('vv',), + 'DH': ('hh', 'hv'), + 'DV': ('vv', 'vh'), + } + key = Path(s1_safe).name[14:16] + + return mapping[key] + + +def loadParsedata(indir, orbit_dir, aux_dir, buffer=0): + ''' + Input file. + ''' + import os + import numpy as np + import isce + from isceobj.Sensor.TOPS.Sentinel1 import Sentinel1 + + + frames = [] + for swath in range(1,4): + rdr=Sentinel1() + rdr.configure() +# rdr.safe=['./S1A_IW_SLC__1SDH_20180401T100057_20180401T100124_021272_024972_8CAF.zip'] + rdr.safe=[indir] + rdr.output='reference' + rdr.orbitDir=orbit_dir + rdr.auxDir=aux_dir + rdr.swathNumber=swath + rdr.polarization=get_polarizations(indir)[0] + rdr.parse() + frames.append(rdr.product) + + info = Dummy() + info.sensingStart = min([x.sensingStart for x in frames]) + info.sensingStop = max([x.sensingStop for x in frames]) + info.startingRange = min([x.startingRange for x in frames]) + info.farRange = max([x.farRange for x in frames]) + info.prf = 1.0 / frames[0].bursts[0].azimuthTimeInterval + info.rangePixelSize = frames[0].bursts[0].rangePixelSize + info.lookSide = -1 + + info.startingRange -= buffer * info.rangePixelSize + info.farRange += buffer * info.rangePixelSize + + info.numberOfLines = int( np.round( (info.sensingStop - info.sensingStart).total_seconds() * info.prf)) + 1 + info.numberOfSamples = int( np.round( (info.farRange - info.startingRange)/info.rangePixelSize)) + 1 + 2 * buffer + info.orbit = getMergedOrbit(frames) + + return info + def coregisterLoadMetadataOptical(indir_m, indir_s): ''' Input file. @@ -383,8 +454,12 @@ def main(): metadata_m, metadata_s = coregisterLoadMetadataOptical(inps.indir_m, inps.indir_s) runGeogridOptical(metadata_m, metadata_s, inps.demfile, inps.dhdxfile, inps.dhdyfile, inps.vxfile, inps.vyfile, inps.srxfile, inps.sryfile, inps.csminxfile, inps.csminyfile, inps.csmaxxfile, inps.csmaxyfile, inps.ssmfile) else: - metadata_m = loadMetadata(inps.indir_m) - metadata_s = loadMetadata(inps.indir_s) + if inps.parse: + metadata_m = loadParsedata(inps.indir_m, inps.orbit_dir, inputs.aux_dir, inps.buffer) + metadata_s = loadParsedata(inps.indir_s, inps.orbit_dir, inputs.aux_dir, inps.buffer) + else: + metadata_m = loadMetadata(inps.indir_m,inps.buffer) + metadata_s = loadMetadata(inps.indir_s,inps.buffer) runGeogrid(metadata_m, metadata_s, inps.demfile, inps.dhdxfile, inps.dhdyfile, inps.vxfile, inps.vyfile, inps.srxfile, inps.sryfile, inps.csminxfile, inps.csminyfile, inps.csmaxxfile, inps.csmaxyfile, inps.ssmfile)