add grid spacing and the ratio of grid spacing to chipsize0 so that adjacent chips can be dependent (finer sampling grid)

geo_autoRIFT/autoRIFT/autoRIFT.py

@@ -24,9 +24,9 @@
 # authority as may be required before exporting this software to any 'EAR99'
 # embargoed foreign country or citizen of those countries.
 #
-# Author: Yang Lei
+# Author: Yang Lei, Alex S. Gardner
 #
-# Note: this is based on the MATLAB code, "auto-RIFT", written by Alex Gardner,
+# Note: this is based on the MATLAB code, "auto-RIFT", written by Alex S. Gardner,
 #       and has been translated to Python and further optimized.
 #~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
@@ -324,6 +324,33 @@ def autorift(self):
         Dy.fill(np.nan)
 
         Flag = 3
+
+
+        if self.ChipSize0X > self.GridSpacingX:
+            if np.mod(self.ChipSize0X,self.GridSpacingX) != 0:
+                sys.exit('when GridSpacing < smallest allowable chip size (ChipSize0), ChipSize0 must be integer multiples of GridSpacing')
+            else:
+                ChipSize0_GridSpacing_oversample_ratio = int(self.ChipSize0X / self.GridSpacingX)
+        else:
+            ChipSize0_GridSpacing_oversample_ratio = 1
+
+
+        DispFiltC = DISP_FILT()
+        overlap_c = np.max((1 - self.sparseSearchSampleRate / ChipSize0_GridSpacing_oversample_ratio,0))
+        DispFiltC.FracValid = self.FracValid * (1 - overlap_c) + overlap_c**2
+        DispFiltC.FracSearch = self.FracSearch
+        DispFiltC.FiltWidth = (self.FiltWidth - 1) * ChipSize0_GridSpacing_oversample_ratio + 1
+        DispFiltC.Iter = self.Iter - 1
+        DispFiltC.MadScalar = self.MadScalar
+
+        DispFiltF = DISP_FILT()
+        overlap_f = 1 - 1 / ChipSize0_GridSpacing_oversample_ratio
+        DispFiltF.FracValid = self.FracValid * (1 - overlap_f) + overlap_f**2
+        DispFiltF.FracSearch = self.FracSearch
+        DispFiltF.FiltWidth = (self.FiltWidth - 1) * ChipSize0_GridSpacing_oversample_ratio + 1
+        DispFiltF.Iter = self.Iter
+        DispFiltF.MadScalar = self.MadScalar
+
 
         for i in range(ChipSizeUniX.__len__()):
 
@@ -377,22 +404,22 @@ def autorift(self):
             SearchLimitX0[(np.logical_not(idxZero)) & (SearchLimitX0 < self.minSearch)] = self.minSearch
             SearchLimitY0[(np.logical_not(idxZero)) & (SearchLimitY0 < self.minSearch)] = self.minSearch
 
-            if ((xGrid0.shape[0] - 2)/self.sparseSearchSampleRate < 5) | ((xGrid0.shape[1] - 2)/self.sparseSearchSampleRate < 5):
+            if ((xGrid0.shape[0] - 2)/(self.sparseSearchSampleRate * ChipSize0_GridSpacing_oversample_ratio) < 5) | ((xGrid0.shape[1] - 2)/(self.sparseSearchSampleRate * ChipSize0_GridSpacing_oversample_ratio) < 5):
                 Flag = 2
                 return Flag
 
             # Setup for coarse search: sparse sampling / resize
-            rIdxC = slice(self.sparseSearchSampleRate-1,xGrid0.shape[0],self.sparseSearchSampleRate)
-            cIdxC = slice(self.sparseSearchSampleRate-1,xGrid0.shape[1],self.sparseSearchSampleRate)
+            rIdxC = slice((self.sparseSearchSampleRate * ChipSize0_GridSpacing_oversample_ratio)-1,xGrid0.shape[0],(self.sparseSearchSampleRate * ChipSize0_GridSpacing_oversample_ratio))
+            cIdxC = slice((self.sparseSearchSampleRate * ChipSize0_GridSpacing_oversample_ratio)-1,xGrid0.shape[1],(self.sparseSearchSampleRate * ChipSize0_GridSpacing_oversample_ratio))
             xGrid0C = xGrid0[rIdxC,cIdxC]
             yGrid0C = yGrid0[rIdxC,cIdxC]
 
 #            pdb.set_trace()
 
-            if np.remainder(self.sparseSearchSampleRate,2) == 0:
-                filtWidth = self.sparseSearchSampleRate + 1
+            if np.remainder((self.sparseSearchSampleRate * ChipSize0_GridSpacing_oversample_ratio),2) == 0:
+                filtWidth = (self.sparseSearchSampleRate * ChipSize0_GridSpacing_oversample_ratio) + 1
             else:
-                filtWidth = self.sparseSearchSampleRate
+                filtWidth = (self.sparseSearchSampleRate * ChipSize0_GridSpacing_oversample_ratio)
 
             SearchLimitX0C = colfilt(SearchLimitX0.copy(), (int(filtWidth), int(filtWidth)), 0)
             SearchLimitY0C = colfilt(SearchLimitY0.copy(), (int(filtWidth), int(filtWidth)), 0)
@@ -426,14 +453,6 @@ def autorift(self):
             # M0C is the mask for reliable estimates after coarse search, MC is the mask after disparity filtering, MC2 is the mask after area closing for fine search
             M0C = np.logical_not(np.isnan(DxC))
 
-            DispFiltC = DISP_FILT()
-            if ChipSizeUniX[i] == ChipSizeUniX[0]:
-                DispFiltC.FracValid = self.FracValid
-                DispFiltC.FracSearch = self.FracSearch
-                DispFiltC.FiltWidth = self.FiltWidth
-                DispFiltC.Iter = self.Iter
-                DispFiltC.MadScalar = self.MadScalar
-            DispFiltC.Iter = DispFiltC.Iter - 1
             MC = DispFiltC.filtDisp(DxC.copy(), DyC.copy(), SearchLimitX0C.copy(), SearchLimitY0C.copy(), M0C.copy(), overSampleRatio)
 
             MC[np.logical_not(M0C)] = False
@@ -444,7 +463,7 @@ def autorift(self):
                 continue
 
             MC2 = ndimage.distance_transform_edt(np.logical_not(MC)) < self.BuffDistanceC
-            dstShape = (int(MC2.shape[0]*self.sparseSearchSampleRate),int(MC2.shape[1]*self.sparseSearchSampleRate))
+            dstShape = (int(MC2.shape[0]*(self.sparseSearchSampleRate*ChipSize0_GridSpacing_oversample_ratio)),int(MC2.shape[1]*(self.sparseSearchSampleRate*ChipSize0_GridSpacing_oversample_ratio)))
 
             MC2 = cv2.resize(MC2.astype(np.uint8),dstShape[::-1],interpolation=cv2.INTER_NEAREST).astype(np.bool)
 #            pdb.set_trace()
@@ -473,14 +492,6 @@ def autorift(self):
 
 #            pdb.set_trace()
 
-            DispFiltF = DISP_FILT()
-            if ChipSizeUniX[i] == ChipSizeUniX[0]:
-                DispFiltF.FracValid = self.FracValid
-                DispFiltF.FracSearch = self.FracSearch
-                DispFiltF.FiltWidth = self.FiltWidth
-                DispFiltF.Iter = self.Iter
-                DispFiltF.MadScalar = self.MadScalar
-
 
             M0 = DispFiltF.filtDisp(DxF.copy(), DyF.copy(), SearchLimitX0.copy(), SearchLimitY0.copy(), np.logical_not(np.isnan(DxF)), overSampleRatio)
 #            pdb.set_trace()
@@ -636,6 +647,7 @@ def __init__(self):
         self.ChipSizeMinX = 32
         self.ChipSizeMaxX = 64
         self.ChipSize0X = 32
+        self.GridSpacingX = 32
         self.ScaleChipSizeY = 1
         self.SearchLimitX = 25
         self.SearchLimitY = 25
@@ -645,7 +657,7 @@ def __init__(self):
         self.minSearch = 6
         self.sparseSearchSampleRate = 4
         self.FracValid = 8/25
-        self.FracSearch = 0.25
+        self.FracSearch = 0.20
         self.FiltWidth = 5
         self.Iter = 3
         self.MadScalar = 4
@@ -1381,7 +1393,7 @@ def __init__(self):
         ##filter parameters; try different parameters to decide how much fine-resolution estimates we keep, which can make the final images smoother
 
         self.FracValid = 8/25
-        self.FracSearch = 0.25
+        self.FracSearch = 0.20
         self.FiltWidth = 5
         self.Iter = 3
         self.MadScalar = 4
@@ -1390,6 +1402,9 @@ def filtDisp(self, Dx, Dy, SearchLimitX, SearchLimitY, M, OverSampleRatio):
 
         import numpy as np
 
+        if np.mod(self.FiltWidth,2) == 0:
+            sys.exit('NDC filter width must be an odd number')
+
         dToleranceX = self.FracValid * self.FiltWidth**2
         dToleranceY = self.FracValid * self.FiltWidth**2
 #        pdb.set_trace()

geo_autoRIFT/autoRIFT/autoRIFT_ISCE.py

@@ -65,7 +65,14 @@
         default = 32,
         type = int,
         mandatory = False,
-        doc = 'Minimum acceptable size (in X direction) of the reference data window to be used for correlation')
+        doc = 'Minimum acceptable size (in X direction) of the reference data window to be used for correlation without resampling the grid; if a chip size greater than this value is provided, need to resize the sampling grid')
+
+GRID_SPACING_X = Component.Parameter('GridSpacingX',
+        public_name='GRID_SPACING_X',
+        default = 32,
+        type = int,
+        mandatory = False,
+        doc = 'Spacing (in X direction) of the sampling grid')
 
 SCALE_CHIP_SIZE_Y = Component.Parameter('ScaleChipSizeY',
         public_name='SCALE_CHIP_SIZE_Y',
@@ -132,7 +139,7 @@
 
 FRAC_SEARCH = Component.Parameter('FracSearch',
         public_name = 'FRAC_SEARCH',
-        default = 0.25,
+        default = 0.20,
         type = float,
         mandatory=False,
         doc = 'Fraction of search')
@@ -214,6 +221,7 @@ class autoRIFT_ISCE(autoRIFT, Component):
                       CHIP_SIZE_MIN_X,
                       CHIP_SIZE_MAX_X,
                       CHIP_SIZE_0X,
+                      GRID_SPACING_X,
                       SCALE_CHIP_SIZE_Y,
                       SEARCH_LIMIT_X,
                       SEARCH_LIMIT_Y,

geo_autoRIFT/geogrid/Geogrid.py

@@ -252,6 +252,7 @@ def setState(self):
         geogrid.setEPSG_Py(self._geogrid, self.epsg)
         geogrid.setIncidenceAngle_Py(self._geogrid, self.incidenceAngle)
         geogrid.setChipSizeX0_Py(self._geogrid, self.chipSizeX0)
+        geogrid.setGridSpacingX_Py(self._geogrid, self.gridSpacingX)
 
         geogrid.setXLimits_Py(self._geogrid, self._xlim[0], self._xlim[1])
         geogrid.setYLimits_Py(self._geogrid, self._ylim[0], self._ylim[1])
@@ -319,6 +320,7 @@ def __init__(self):
         self.repeatTime = None
         self.incidenceAngle = None
         self.chipSizeX0 = None
+        self.gridSpacingX = None
 
         ##Input related parameters
         self.demname = None

geo_autoRIFT/geogrid/GeogridOptical.py

@@ -163,6 +163,7 @@ def geogrid(self):
         print("Map inputs: ")
         print("EPSG: " + str(self.epsgDem))
         print("Smallest Allowable Chip Size in m: " + str(self.chipSizeX0))
+        print("Grid spacing in m: " + str(self.gridSpacingX))
         print("Repeat Time: " + str(self.repeatTime))
         print("XLimits: " + str(self._xlim[0]) + "  " + str(self._xlim[1]))
         print("YLimits: " + str(self._ylim[0]) + "  " + str(self._ylim[1]))
@@ -811,6 +812,7 @@ def __init__(self):
         self.numberOfLines = None
         self.repeatTime = None
         self.chipSizeX0 = None
+        self.gridSpacingX = None
 
         ##Input related parameters
         self.dat1name = None

geo_autoRIFT/geogrid/bindings/geogridmodule.cpp

@@ -116,6 +116,18 @@ PyObject* setChipSizeX0(PyObject *self, PyObject *args)
     return Py_BuildValue("i", 0);
 }
 
+PyObject* setGridSpacingX(PyObject *self, PyObject *args)
+{
+    uint64_t ptr;
+    double gridSpacingX;
+    if (!PyArg_ParseTuple(args, "Kd", &ptr, &gridSpacingX))
+    {
+        return NULL;
+    }
+    ((geoGrid*)(ptr))->gridSpacingX = gridSpacingX;
+    return Py_BuildValue("i", 0);
+}
+
 PyObject* setRepeatTime(PyObject *self, PyObject *args)
 {
     uint64_t ptr;

geo_autoRIFT/geogrid/include/geogrid.h

@@ -53,6 +53,7 @@ struct geoGrid
     std::string ssmname;     //Stable surface mask
     int epsgcode;
     double chipSizeX0;
+    double gridSpacingX;
 
     //Bounding box related
     double xmin, xmax;

geo_autoRIFT/geogrid/include/geogridmodule.h

@@ -65,6 +65,7 @@ extern "C"
         PyObject * setEPSG(PyObject *, PyObject *);
         PyObject * setIncidenceAngle(PyObject *, PyObject *);
         PyObject * setChipSizeX0(PyObject *, PyObject *);
+        PyObject * setGridSpacingX(PyObject *, PyObject *);
         PyObject * setXLimits(PyObject *, PyObject *);
         PyObject * setYLimits(PyObject *, PyObject *);
         PyObject * getXPixelSize(PyObject *, PyObject *);
@@ -88,6 +89,7 @@ static PyMethodDef geogrid_methods[] =
         {"setEPSG_Py", setEPSG, METH_VARARGS, " "},
         {"setIncidenceAngle_Py", setIncidenceAngle, METH_VARARGS, " "},
         {"setChipSizeX0_Py", setChipSizeX0, METH_VARARGS, " "},
+        {"setGridSpacingX_Py", setGridSpacingX, METH_VARARGS, " "},
         {"setVelocities_Py", setVelocities, METH_VARARGS, " "},
         {"setSearchRange_Py", setSearchRange, METH_VARARGS, " "},
         {"setChipSizeMin_Py", setChipSizeMin, METH_VARARGS, " "},

geo_autoRIFT/geogrid/src/geogrid.cpp

@@ -56,6 +56,7 @@ void geoGrid::geogrid()
     std::cout << "\nMap inputs: \n";
     std::cout << "EPSG: " << epsgcode << "\n";
     std::cout << "Smallest Allowable Chip Size in m: " << chipSizeX0 << "\n";
+    std::cout << "Grid spacing in m: " << gridSpacingX << "\n";
     std::cout << "Repeat Time: " << dt << "\n";
     std::cout << "XLimits: " << xmin << "  " << xmax << "\n";
     std::cout << "YLimits: " << ymin << "  " << ymax << "\n";