diff --git a/README.md b/README.md index b8984d9..412a0ac 100644 --- a/README.md +++ b/README.md @@ -623,7 +623,23 @@ between three files as follows: 20061231 ``` - +### rtcApp.xml +The inputs are Sentinel GRD zipfiles +```xml + + /Users/data/sentinel1 + + 20 + sentinel1 + + $dir$/rtcApp/data/S1A_IW_GRDH_1SDV_20181221T225104_20181221T225129_025130_02C664_B46C.zip + $dir$/orbits + $dir$/rtcApp/output + [VV, VH] + + + +``` ----- ## Component Configurability diff --git a/applications/imageMath.py b/applications/imageMath.py index 3977a89..432e262 100755 --- a/applications/imageMath.py +++ b/applications/imageMath.py @@ -293,6 +293,7 @@ def main(args, files): #######Determine number of input and output bands bandList = [] + iMath['equations'] = [] for ii,expr in enumerate(args.equation.split(';')): #####Now parse the equation to get the file names used @@ -319,7 +320,11 @@ def main(args, files): ######Create input memmaps for ii,infile in enumerate(fileList): - fstr, files = parseInputFile(infile, files) + if type(files) == list: + fstr, files = parseInputFile(infile, files) + else: + fstr = getattr(files, infile) + logger.debug('Input string for File %d: %s: %s'%(ii, infile, fstr)) if len(fstr.split(';')) > 1: @@ -341,8 +346,9 @@ def main(args, files): if bbox is not None: iMath['bboxes'].append(bbox) - if len(files): - raise IOError('Unused input variables set:\n'+ ' '.join(files)) + if type(files) == list: + if len(files): + raise IOError('Unused input variables set:\n'+ ' '.join(files)) #######Some debugging logger.debug('List of available bands: ' + str(iMath['inBands'].keys())) diff --git a/applications/rtcApp.py b/applications/rtcApp.py index ab45a7e..836151e 100755 --- a/applications/rtcApp.py +++ b/applications/rtcApp.py @@ -155,7 +155,7 @@ NUMBER_RANGE_LOOKS = Application.Parameter('numberRangeLooks', ) POSTING = Application.Parameter('posting', - public_name='azimuth looks', + public_name='posting', default = 20.0, type = float, mandatory = False, @@ -363,6 +363,7 @@ class GRDSAR(Application): self.verifyDEM = RtcProc.createVerifyDEM(self) self.multilook = RtcProc.createLooks(self) self.runTopo = RtcProc.createTopo(self) + self.runNormalize = RtcProc.createNormalize(self) # self.runGeocode = RtcProc.createGeocode(self) return None @@ -392,6 +393,9 @@ class GRDSAR(Application): ##Run topo for each bursts self.step('topo', func=self.runTopo) + ##Run normalize to get gamma0 + self.step('normalize', func=self.runNormalize) + # Geocode # self.step('geocode', func=self.runGeocode, # args=(self.geocode_list, self.do_unwrap, self.geocode_bbox)) @@ -416,6 +420,9 @@ class GRDSAR(Application): ##Run topo for each burst self.runTopo() + + ##Run normalize to get gamma0 + self.runNormalize() ###Compute covariance # self.runEstimateCovariance() diff --git a/components/isceobj/RtcProc/Factories.py b/components/isceobj/RtcProc/Factories.py index c779864..038f1f7 100644 --- a/components/isceobj/RtcProc/Factories.py +++ b/components/isceobj/RtcProc/Factories.py @@ -46,5 +46,6 @@ createPreprocessor = _factory("runPreprocessor") createVerifyDEM = _factory("runVerifyDEM") createLooks = _factory("runLooks") createTopo = _factory("runTopo") +createNormalize = _factory("runNormalize") #createGeocode = _factory("runGeocode") diff --git a/components/isceobj/RtcProc/RtcProc.py b/components/isceobj/RtcProc/RtcProc.py index f529106..aac8b9a 100644 --- a/components/isceobj/RtcProc/RtcProc.py +++ b/components/isceobj/RtcProc/RtcProc.py @@ -69,7 +69,7 @@ INC_FILENAME = Component.Parameter( GAMMA0_FILENAME = Component.Parameter( 'gamma0FileName', public_name='Gamma0 backscatter file', - default = 'gamma0.rdr', + default = 'gamma0.img', type = str, mandatory = False, doc = 'Unmasked gamma0 backscatter file') diff --git a/components/isceobj/RtcProc/runNormalize.py b/components/isceobj/RtcProc/runNormalize.py index ce1e6c9..70ce9cb 100644 --- a/components/isceobj/RtcProc/runNormalize.py +++ b/components/isceobj/RtcProc/runNormalize.py @@ -1,4 +1,4 @@ -# +#!/usr/bin/env python3 # Author: Piyush Agram # Copyright 2016 # @@ -6,19 +6,23 @@ import logging import isceobj import mroipac +from .runTopo import filenameWithLooks +from .runLooks import takeLooks import os +import itertools import numpy as np from isceobj.Util.decorators import use_api +from applications import imageMath logger = logging.getLogger('isce.grdsar.looks') - +class Dummy: + pass def runNormalize(self): ''' Make sure that a DEM is available for processing the given data. ''' - refPol = self._grd.polarizations[0] master = self._grd.loadProduct( os.path.join(self._grd.outputFolder, 'beta_{0}.xml'.format(refPol))) @@ -26,17 +30,31 @@ def runNormalize(self): azlooks, rglooks = self._grd.getLooks( self.posting, master.groundRangePixelSize, master.azimuthPixelSize, self.numberAzimuthLooks, self.numberRangeLooks) - if (azlooks == 1) and (rglooks == 1): - return - - slantRange = False for pol in self._grd.polarizations: - inname = os.path.join( self._grd.outputFolder, 'beta_{0}.img'.format(pol) ) - takeLooks(inname, azlooks, rglooks) + if (azlooks == 1) and (rglooks == 1): + inname = os.path.join( self._grd.outputFolder, 'beta_{0}.img'.format(pol)) + else: + inname = os.path.join( self._grd.outputFolder, filenameWithLooks('beta_{0}.img'.format(pol), azlooks, rglooks)) - if not slantRange: - inname = master.slantRangeImage.filename - takeLooks(inname, azlooks, rglooks) - slantRange = True + basefolder, output = os.path.split(self._grd.outputFolder) + incname = os.path.join(basefolder, self._grd.geometryFolder, self._grd.incFileName) + outname = os.path.join(self._grd.outputFolder, filenameWithLooks('gamma_{0}'.format(pol)+'.img', azlooks, rglooks)) + maskname = os.path.join(basefolder, self._grd.geometryFolder, self._grd.slMaskFileName) + + args = imageMath.createNamespace() + args.equation = 'a*cos(b_0*PI/180.)/cos(b_1*PI/180.) * (c==0)' + args.dtype = np.float32 + args.scheme = 'BIL' + args.out = outname + #args.debug = True + + files = Dummy() + files.a = inname + files.b = incname + files.c = maskname + + + + imageMath.main(args, files) return diff --git a/components/isceobj/Sensor/GRD/Sentinel1.py b/components/isceobj/Sensor/GRD/Sentinel1.py index 45e648b..b910885 100755 --- a/components/isceobj/Sensor/GRD/Sentinel1.py +++ b/components/isceobj/Sensor/GRD/Sentinel1.py @@ -261,10 +261,10 @@ class Sentinel1(Component): self.validateUserInputs() - if self.xml.startswith('/vsizip'): #Read from zip file + if '.zip' in self.xml: try: parts = self.xml.split(os.path.sep) - zipname = os.path.join(*(parts[2:-3])) + zipname = os.path.join('/',*(parts[:-3])) fname = os.path.join(*(parts[-3:])) with zipfile.ZipFile(zipname, 'r') as zf: @@ -283,23 +283,22 @@ class Sentinel1(Component): self.populateMetadata() self.populateBbox() - ####Tru and locate an orbit file + ####Try and locate an orbit file if self.orbitFile is None: if self.orbitDir is not None: self.orbitFile = self.findOrbitFile() + print('Found this orbitfile: %s' %self.orbitFile) - ####Read in the orbits - if self.orbitFile: + if '_POEORB_' in self.orbitFile: orb = self.extractPreciseOrbit() - else: + elif '_RESORB_' in self.orbitFile: orb = self.extractOrbit() - + self.product.orbit.setOrbitSource('Header') for sv in orb: self.product.orbit.addStateVector(sv) - self.populateIPFVersion() self.extractBetaLUT() self.extractNoiseLUT() @@ -423,10 +422,11 @@ class Sentinel1(Component): nsp = "{http://www.esa.int/safe/sentinel-1.0}" - if self.manifest.startswith('/vsizip'): + if '.zip' in self.manifest: + import zipfile parts = self.manifest.split(os.path.sep) - zipname = os.path.join(*(parts[2:-2])) + zipname = os.path.join('/',*(parts[:-2])) fname = os.path.join(*(parts[-2:])) try: @@ -462,38 +462,40 @@ class Sentinel1(Component): datefmt = "%Y%m%dT%H%M%S" types = ['POEORB', 'RESORB'] + filelist = [] match = [] - timeStamp = self.product.sensingMid - + timeStamp = self.product.sensingStart+(self.product.sensingStop - self.product.sensingStart)/2. + for orbType in types: files = glob.glob( os.path.join(self.orbitDir, 'S1A_OPER_AUX_' + orbType + '_OPOD*')) - + filelist.extend(files) ###List all orbit files - for result in files: - fields = result.split('_') - taft = datetime.datetime.strptime(fields[-1][0:15], datefmt) - tbef = datetime.datetime.strptime(fields[-2][1:16], datefmt) - - #####Get all files that span the acquisition - if (tbef <= timeStamp) and (taft >= timeStamp): - tmid = tbef + 0.5 * (taft - tbef) - match.append((result, abs((timeStamp-tmid).total_seconds()))) - #####Return the file with the image is aligned best to the middle of the file - if len(match) != 0: - bestmatch = min(match, key = lambda x: x[1]) - return bestmatch[0] + for result in filelist: + fields = result.split('_') + taft = datetime.datetime.strptime(fields[-1][0:15], datefmt) + tbef = datetime.datetime.strptime(fields[-2][1:16], datefmt) + print(taft, tbef) + + #####Get all files that span the acquisition + if (tbef <= timeStamp) and (taft >= timeStamp): + tmid = tbef + 0.5 * (taft - tbef) + match.append((result, abs((timeStamp-tmid).total_seconds()))) + #####Return the file with the image is aligned best to the middle of the file + if len(match) != 0: + bestmatch = min(match, key = lambda x: x[1]) + return bestmatch[0] - if len(match) == 0: - raise Exception('No suitable orbit file found. If you want to process anyway - unset the orbitdir parameter') + if len(match) == 0: + raise Exception('No suitable orbit file found. If you want to process anyway - unset the orbitdir parameter') def extractOrbit(self): ''' Extract orbit information from xml node. ''' node = self._xml_root.find('generalAnnotation/orbitList') - + print('Extracting orbit from annotation XML file') frameOrbit = Orbit() frameOrbit.configure() @@ -516,13 +518,7 @@ class Sentinel1(Component): vec.setVelocity(vel) frameOrbit.addStateVector(vec) - - orbExt = OrbitExtender(planet=Planet(pname='Earth')) - orbExt.configure() - newOrb = orbExt.extendOrbit(frameOrbit) - - - return newOrb + return frameOrbit def extractPreciseOrbit(self): ''' @@ -534,11 +530,10 @@ class Sentinel1(Component): print("IOError: %s" % strerr) return - _xml_root = ElementTree(file=fp).getroot() - + _xml_root = ElementTree.ElementTree(file=fp).getroot() + node = _xml_root.find('Data_Block/List_of_OSVs') - print('Extracting orbit from Orbit File: ', self.orbitFile) orb = Orbit() orb.configure() @@ -582,10 +577,10 @@ class Sentinel1(Component): if self.calibrationXml is None: raise Exception('No calibration file provided') - if self.calibrationXml.startswith('/vsizip'): + if '.zip' in self.calibrationXml: import zipfile parts = self.calibrationXml.split(os.path.sep) - zipname = os.path.join(*(parts[2:-4])) + zipname = os.path.join('/',*(parts[:-4])) fname = os.path.join(*(parts[-4:])) try: @@ -723,7 +718,7 @@ class Sentinel1(Component): print('Extracting normalized image ....') - src = gdal.Open(self.tiff.strip(), gdal.GA_ReadOnly) + src = gdal.Open('/vsizip//'+self.tiff.strip(), gdal.GA_ReadOnly) band = src.GetRasterBand(1) if self.product.numberOfSamples != src.RasterXSize: diff --git a/contrib/stack/stripmapStack/crossmul.py b/contrib/stack/stripmapStack/crossmul.py index 8ffeeef..7b07716 100755 --- a/contrib/stack/stripmapStack/crossmul.py +++ b/contrib/stack/stripmapStack/crossmul.py @@ -33,6 +33,8 @@ def cmdLineParse(iargs = None): def run(imageSlc1, imageSlc2, resampName, azLooks, rgLooks): objSlc1 = isceobj.createSlcImage() + #right now imageSlc1 and 2 are just text files, need to open them as image + IU.copyAttributes(imageSlc1, objSlc1) objSlc1.setAccessMode('read') objSlc1.createImage() @@ -81,7 +83,6 @@ def run(imageSlc1, imageSlc2, resampName, azLooks, rgLooks): def main(iargs=None): - inps = cmdLineParse(iargs) img1 = isceobj.createImage() @@ -96,9 +97,8 @@ def main(iargs=None): run(img1, img2, inps.prefix, inps.azlooks, inps.rglooks) if __name__ == '__main__': - + + main() ''' Main driver. ''' - -