diff --git a/applications/topsApp.py b/applications/topsApp.py index e79d9ee..6dd2b65 100755 --- a/applications/topsApp.py +++ b/applications/topsApp.py @@ -438,6 +438,20 @@ ION_DO_ION = Application.Parameter('ION_doIon', mandatory = False, doc = '') +ION_APPLY_ION = Application.Parameter('ION_applyIon', + public_name = 'apply ionosphere correction', + default = False, + type = bool, + mandatory = False, + doc = '') + +ION_CONSIDER_BURST_PROPERTIES = Application.Parameter('ION_considerBurstProperties', + public_name = 'consider burst properties in ionosphere computation', + default = False, + type = bool, + mandatory = False, + doc = '') + ION_START_STEP = Application.Parameter( 'ION_startStep', public_name='start ionosphere step', @@ -649,6 +663,8 @@ class TopsInSAR(Application): ######################################################## #for ionospheric correction ION_DO_ION, + ION_APPLY_ION, + ION_CONSIDER_BURST_PROPERTIES, ION_START_STEP, ION_END_STEP, ION_ION_HEIGHT, @@ -734,6 +750,9 @@ class TopsInSAR(Application): if(self.geocode_list is None): #if not provided by the user use the list from InsarProc self.geocode_list = self.insar.geocode_list + #for ionosphere + if 'topophase.ion' not in self.geocode_list: + self.geocode_list.append('topophase.ion') else: #if geocode_list defined here, then give it to InsarProc #for consistency between insarApp and InsarProc and warn the user diff --git a/components/isceobj/TopsProc/runIon.py b/components/isceobj/TopsProc/runIon.py index c2a98aa..06b48a1 100644 --- a/components/isceobj/TopsProc/runIon.py +++ b/components/isceobj/TopsProc/runIon.py @@ -42,6 +42,7 @@ def setup(self): #SECTION 1. PROCESSING CONTROL PARAMETERS #1. suggested default values of the parameters ionParam.doIon = False + ionParam.considerBurstProperties = False ionParam.startStep = ionParam.allSteps[0] ionParam.endStep = ionParam.allSteps[-1] @@ -77,6 +78,7 @@ def setup(self): #2. accept the above parameters from topsApp.py ionParam.doIon = self.ION_doIon + ionParam.considerBurstProperties = self.ION_considerBurstProperties ionParam.startStep = self.ION_startStep ionParam.endStep = self.ION_endStep @@ -2637,16 +2639,17 @@ def runIon(self): if run_step('filt_gaussian', ionParam): filt_gaussian(self, ionParam) + #only do the following steps when considering burst properties #ionosphere shift - if run_step('ionosphere_shift', ionParam): + if run_step('ionosphere_shift', ionParam) and ionParam.considerBurstProperties: ionosphere_shift(self, ionParam) #resample from ionospheric layer to ground layer, get ionosphere for each burst - if run_step('ion2grd', ionParam): + if run_step('ion2grd', ionParam) and ionParam.considerBurstProperties: ion2grd(self, ionParam) #esd - if run_step('esd', ionParam): + if run_step('esd', ionParam) and ionParam.considerBurstProperties: esd(self, ionParam) #pure esd without applying ionospheric correction diff --git a/components/isceobj/TopsProc/runMergeBursts.py b/components/isceobj/TopsProc/runMergeBursts.py index 8e30543..8839ae0 100755 --- a/components/isceobj/TopsProc/runMergeBursts.py +++ b/components/isceobj/TopsProc/runMergeBursts.py @@ -20,6 +20,38 @@ import logging from isceobj.Util.ImageUtil import ImageLib as IML +def interpolateDifferentNumberOfLooks(inputData, lengtho, widtho, nrli, nali, nrlo, nalo): + ''' + inputData: input numpy 2-d array + lengtho: length of output array + widtho: width of output array + nrli: number of range looks input + nali: number of azimuth looks input + nrlo: number of range looks output + nalo: number of azimuth looks output + ''' + import numpy as np + from scipy.interpolate import interp1d + + (lengthi, widthi) = inputData.shape + + indexi = np.linspace(0, widthi-1, num=widthi, endpoint=True) + indexo = np.linspace(0, widtho-1, num=widtho, endpoint=True) * nrli/nrlo + (nrli-nrlo)/(2.0*nrlo) + outputData0 = np.zeros((lengthi, widtho), dtype=inputData.dtype) + for i in range(lengthi): + f = interp1d(indexi, inputData[i,:], kind='cubic', fill_value="extrapolate") + outputData0[i, :] = f(indexo) + + indexi = np.linspace(0, lengthi-1, num=lengthi, endpoint=True) + indexo = np.linspace(0, lengtho-1, num=lengtho, endpoint=True) * nali/nalo + (nali-nalo)/(2.0*nalo) + outputData = np.zeros((lengtho, widtho), dtype=inputData.dtype) + for j in range(widtho): + f = interp1d(indexi, outputData0[:, j], kind='cubic', fill_value="extrapolate") + outputData[:, j] = f(indexo) + + return outputData + + def mergeBox(frame): ''' Merging using VRTs. @@ -666,14 +698,16 @@ def runMergeBursts(self, adjust=1): #totalLooksThreshold = 9 totalLooksThreshold = 99999999999999 #if doing ionospheric correction - ionCorrection = self.ION_doIon + doIon = self.ION_doIon + applyIon = self.ION_applyIon + considerBurstProperties = self.ION_considerBurstProperties ionDirname = 'ion/ion_burst' mergedIonname = 'topophase.ion' originalIfgname = 'topophase_ori.flat' ######################################### # backing out the tigher constraints for ionosphere as it could itnroduce gabs between along track products produced seperately - if not ionCorrection: + if not (doIon and considerBurstProperties): adjust=0 ######################################### @@ -712,7 +746,7 @@ def runMergeBursts(self, adjust=1): #restore frames = frames_bak else: - validOnly==True + validOnly=True ######################################### @@ -738,7 +772,7 @@ def runMergeBursts(self, adjust=1): mergeBursts2(frames, os.path.join(self._insar.fineIfgDirname, 'IW%d', 'burst_%02d.int'), burstIndex, box, os.path.join(mergedir, self._insar.mergedIfgname+suffix), virtual=virtual, validOnly=True) if self.numberAzimuthLooks * self.numberRangeLooks < totalLooksThreshold: mergeBursts2(frames, os.path.join(self._insar.fineIfgDirname, 'IW%d', 'burst_%02d.cor'), burstIndex, box, os.path.join(mergedir, self._insar.correlationFilename+suffix), virtual=virtual, validOnly=True) - if ionCorrection == True: + if doIon and considerBurstProperties: mergeBursts2(frames, os.path.join(ionDirname, 'IW%d', 'burst_%02d.ion'), burstIndex, box, os.path.join(mergedir, mergedIonname+suffix), virtual=virtual, validOnly=True) @@ -782,13 +816,61 @@ def runMergeBursts(self, adjust=1): os.remove(os.path.join(mergedir, pwrfile+'.xml')) os.remove(os.path.join(mergedir, pwrfile+'.vrt')) - if ionCorrection: - multilook(os.path.join(mergedir, mergedIonname+suffix), - outname = os.path.join(mergedir, mergedIonname), - alks = self.numberAzimuthLooks, rlks=self.numberRangeLooks) + if doIon: + if considerBurstProperties: + multilook(os.path.join(mergedir, mergedIonname+suffix), + outname = os.path.join(mergedir, mergedIonname), + alks = self.numberAzimuthLooks, rlks=self.numberRangeLooks) + else: + ionFilt = 'ion/ion_cal/filt.ion' + img = isceobj.createImage() + img.load(ionFilt+'.xml') + ionFiltImage = (np.fromfile(ionFilt, dtype=np.float32).reshape(img.length*2, img.width))[1:img.length*2:2, :] + img = isceobj.createImage() + img.load(os.path.join(mergedir, self._insar.mergedIfgname)+'.xml') + + #interpolate original + ionFiltImageOut = interpolateDifferentNumberOfLooks(ionFiltImage, img.length, img.width, self.numberRangeLooks, self.numberAzimuthLooks, self.ION_numberRangeLooks, self.ION_numberAzimuthLooks) + ionFiltOut = os.path.join(mergedir, mergedIonname) + ionFiltImageOut.astype(np.float32).tofile(ionFiltOut) + + image = isceobj.createImage() + image.setDataType('FLOAT') + image.setFilename(ionFiltOut) + image.extraFilename = ionFiltOut + '.vrt' + image.setWidth(img.width) + image.setLength(img.length) + #image.setAccessMode('read') + #image.createImage() + image.renderHdr() + #image.finalizeImage() else: print('Skipping multi-looking ....') + if self.doInSAR and doIon and (not considerBurstProperties): + ionFilt = 'ion/ion_cal/filt.ion' + img = isceobj.createImage() + img.load(ionFilt+'.xml') + ionFiltImage = (np.fromfile(ionFilt, dtype=np.float32).reshape(img.length*2, img.width))[1:img.length*2:2, :] + img = isceobj.createImage() + img.load(os.path.join(mergedir, self._insar.mergedIfgname+suffix)+'.xml') + + #interpolate original + ionFiltImageOut = interpolateDifferentNumberOfLooks(ionFiltImage, img.length, img.width, self.numberRangeLooks, self.numberAzimuthLooks, self.ION_numberRangeLooks, self.ION_numberAzimuthLooks) + ionFiltOut = os.path.join(mergedir, mergedIonname) + ionFiltImageOut.astype(np.float32).tofile(ionFiltOut) + + image = isceobj.createImage() + image.setDataType('FLOAT') + image.setFilename(ionFiltOut) + image.extraFilename = ionFiltOut + '.vrt' + image.setWidth(img.width) + image.setLength(img.length) + #image.setAccessMode('read') + #image.createImage() + image.renderHdr() + #image.finalizeImage() + ######################################### # STEP 4. APPLY CORRECTIONS @@ -796,8 +878,8 @@ def runMergeBursts(self, adjust=1): #do ionospheric and other corrections here #should also consider suffix, but usually we use multiple looks, so I ignore it for now. if self.doInSAR: - if ionCorrection: - print('user choose to do ionospheric correction') + if doIon and applyIon: + print('user choose to apply ionospheric correction') #define file names interferogramFilename = os.path.join(mergedir, self._insar.mergedIfgname) diff --git a/examples/input_files/topsApp.xml b/examples/input_files/topsApp.xml index 62c39bd..ae3241d 100644 --- a/examples/input_files/topsApp.xml +++ b/examples/input_files/topsApp.xml @@ -183,6 +183,8 @@ By default, ESD is always performed. #ionospheric correction module parameters #the values below are the default values used by the module False +False +False #choose from: ['subband', 'rawion', 'grd2ion', 'filt_gaussian', 'ionosphere_shift', 'ion2grd', 'esd'] subband esd