Merge d03542b into a2ebd88

pyroSAR/_dev_config.py

@@ -91,14 +91,15 @@ def __dir__(self):
                'Calibration': 'Cal',
                'Cross-Correlation': '',
                'LinearToFromdB': 'dB',
+               'Multilook': 'ML',
+               'Read': '',
                'Remove-GRD-Border-Noise': 'bnr',
-               'ThermalNoiseRemoval': 'tnr',
                'SAR-Simulation': 'Sim',
                'SARSim-Terrain-Correction': 'TC',
                'Subset': '',
                'Terrain-Correction': 'TC',
                'Terrain-Flattening': 'TF',
-               'Read': '',
+               'ThermalNoiseRemoval': 'tnr',
                'Write': '',
                'Write (2)': ''}
 

pyroSAR/ancillary.py

@@ -7,6 +7,7 @@
 """
 import os
 import re
+import math
 from datetime import datetime
 from ._dev_config import product_pattern
 
@@ -80,6 +81,47 @@ def groupbyTime(images, function, time):
     return [x[0] if len(x) == 1 else x for x in groups]
 
 
+def multilook_factors(sp_rg, sp_az, tr_rg, tr_az, geometry, incidence):
+    """
+    compute multilooking factors to approximate a target ground range resolution
+    
+    Parameters
+    ----------
+    sp_rg: int or float
+        the range pixel spacing
+    sp_az: int or float
+        the azimuth pixel spacing
+    tr_rg: int or float
+        the range target resolution
+    tr_az: int or float
+        the azimuth target resolution
+    geometry: str
+        the imaging geometry; either 'SLANT_RANGE' or 'GROUND_RANGE'
+    incidence: int or float
+        the angle of incidence
+
+    Returns
+    -------
+    tuple
+        the multilooking factors as (rangelooks, azimuthlooks)
+    """
+    if geometry == 'SLANT_RANGE':
+        # compute the ground range resolution
+        groundRangePS = sp_rg / (math.sin(math.radians(incidence)))
+        rlks = int(math.floor(float(tr_rg) / groundRangePS))
+    elif geometry == 'GROUND_RANGE':
+        rlks = int(math.floor(float(tr_rg) / sp_rg))
+    else:
+        raise ValueError("parameter 'geometry' must be either 'SLANT_RANGE' or 'GROUND_RANGE'")
+    # compute the azimuth looks
+    azlks = int(math.floor(float(tr_az) / sp_az))
+
+    # set the look factors to 1 if they were computed to be 0
+    rlks = rlks if rlks > 0 else 1
+    azlks = azlks if azlks > 0 else 1
+    return rlks, azlks
+
+
 def seconds(filename):
     """
     function to extract time in seconds from a file name.

pyroSAR/drivers.py

@@ -1297,16 +1297,17 @@ def scanMetadata(self):
         meta['totalSlices'] = int(tree.find('.//s1sarl1:totalSlices', namespaces).text)
 
         annotations = self.findfiles(self.pattern_ds)
-        ann_xml = self.getFileObj(annotations[0])
-        ann_tree = ET.fromstring(ann_xml.read())
-        ann_xml.close()
+        with self.getFileObj(annotations[0]) as ann_xml:
+            ann_tree = ET.fromstring(ann_xml.read())
+
         meta['spacing'] = tuple([float(ann_tree.find('.//{}PixelSpacing'.format(dim)).text)
                                  for dim in ['range', 'azimuth']])
         meta['samples'] = int(ann_tree.find('.//imageAnnotation/imageInformation/numberOfSamples').text)
         meta['lines'] = int(ann_tree.find('.//imageAnnotation/imageInformation/numberOfLines').text)
         heading = float(ann_tree.find('.//platformHeading').text)
         meta['heading'] = heading if heading > 0 else heading + 360
         meta['incidence'] = float(ann_tree.find('.//incidenceAngleMidSwath').text)
+        meta['image_geometry'] = ann_tree.find('.//projection').text.replace(' ', '_').upper()
 
         return meta
 

pyroSAR/gamma/util.py

@@ -30,6 +30,7 @@
 from ..S1 import OSV
 from ..drivers import ID, CEOS_ERS, CEOS_PSR, ESA, SAFE, TSX, identify
 from . import ISPPar, Namespace, par2hdr
+from ..ancillary import multilook_factors
 
 try:
     from .api import diff, disp, isp, lat
@@ -496,7 +497,7 @@ def geocode(scene, dem, tempdir, outdir, targetres, scaling='linear', func_geoba
         scene = identify(scene)
     else:
         raise RuntimeError("'scene' must be of type str or pyroSAR.ID")
-        
+
     if scene.sensor not in ['S1A', 'S1B']:
         raise IOError('this method is currently only available for Sentinel-1. Please stay tuned...')
 
@@ -900,19 +901,12 @@ def multilook(infile, outfile, targetres, logpath=None, outdir=None, shellscript
     # read the input parameter file
     par = ISPPar(infile + '.par')
 
-    # compute the range looks
-    if par.image_geometry == 'SLANT_RANGE':
-        # compute the ground range resolution
-        groundRangePS = par.range_pixel_spacing / (math.sin(math.radians(par.incidence_angle)))
-        rlks = int(round(float(targetres) / groundRangePS))
-    else:
-        rlks = int(round(float(targetres) / par.range_pixel_spacing))
-    # compute the azimuth looks
-    azlks = int(round(float(targetres) / par.azimuth_pixel_spacing))
-
-    # set the look factors to 1 if they were computed to be 0
-    rlks = rlks if rlks > 0 else 1
-    azlks = azlks if azlks > 0 else 1
+    rlks, azlks = multilook_factors(sp_rg=par.range_pixel_spacing,
+                                    sp_az=par.azimuth_pixel_spacing,
+                                    tr_rg=targetres,
+                                    tr_az=targetres,
+                                    geometry=par.image_geometry,
+                                    incidence=par.incidence_angle)
 
     pars = {'rlks': rlks,
             'azlks': azlks,

pyroSAR/snap/recipes/nodes/Multilook.xml

@@ -0,0 +1,13 @@
+<node id="Multilook">
+    <operator>Multilook</operator>
+    <sources>
+        <sourceProduct refid="Calibration"/>
+    </sources>
+    <parameters class="com.bc.ceres.binding.dom.XppDomElement">
+        <sourceBands/>
+        <nRgLooks/>
+        <nAzLooks/>
+        <outputIntensity>true</outputIntensity>
+        <grSquarePixel>true</grSquarePixel>
+    </parameters>
+</node>

pyroSAR/snap/util.py

@@ -3,8 +3,8 @@
 # John Truckenbrodt, 2016-2019
 ####################################################################
 import os
-
 import pyroSAR
+from ..ancillary import multilook_factors
 from .auxil import parse_recipe, parse_suffix, write_recipe, parse_node, insert_node, gpt
 
 from spatialist import crsConvert, Vector, Raster, bbox, intersect
@@ -198,6 +198,28 @@ def geocode(infile, outdir, t_srs=4326, tr=20, polarizations='all', shapefile=No
         raise RuntimeError('geocode_type not recognized')
 
     ############################################
+    # Multilook node configuration
+
+    try:
+        image_geometry = id.meta['image_geometry']
+        incidence = id.meta['incidence']
+    except KeyError:
+        raise RuntimeError('This function does not yet support sensor {}'.format(id.sensor))
+
+    rlks, azlks = multilook_factors(sp_rg=id.spacing[0],
+                                    sp_az=id.spacing[1],
+                                    tr_rg=tr,
+                                    tr_az=tr,
+                                    geometry=image_geometry,
+                                    incidence=incidence)
+
+    if azlks > 1 or rlks > 1:
+        insert_node(workflow, parse_node('Multilook'), after=tf.attrib['id'])
+        ml = workflow.find('.//node[@id="Multilook"]')
+        ml.find('.//parameters/sourceBands').text = bands_beta
+        ml.find('.//parameters/nAzLooks').text = str(azlks)
+        ml.find('.//parameters/nRgLooks').text = str(rlks)
+    ############################################
     # specify spatial resolution and coordinate reference system of the output dataset
     # print('-- configuring CRS')
     tc.find('.//parameters/pixelSpacingInMeter').text = str(tr)