Merge pull request #29 from leiyangleon/upgrade_to_v1p3p0

upgrade autoRIFT and Geogrid to v1.3.0 for memory and runtime improvement
nasa-jpl · May 24, 2021 · e940875 · e940875
2 parents b973c1b 1557a51
commit e940875
Show file tree

Hide file tree

Showing 20 changed files with 2,528 additions and 686 deletions.
diff --git a/geo_autoRIFT/SConscript b/geo_autoRIFT/SConscript
@@ -34,4 34,4 @@ autorift='autoRIFT/SConscript'
 SConscript(autorift)
 
 geogrid='geogrid/SConscript'
-SConscript(geogrid)
 SConscript(geogrid)
diff --git a/geo_autoRIFT/autoRIFT/__init__.py b/geo_autoRIFT/autoRIFT/__init__.py
@@ -13,4 13,4 @@
  # this means ISCE support not available. Don't raise error. Allow standalone use
  pass
 
-__version__ = '1.2.0'
 __version__ = '1.3.0'
diff --git a/geo_autoRIFT/autoRIFT/autoRIFT.py b/geo_autoRIFT/autoRIFT/autoRIFT.py
@@ -105,7 105,7 @@ def preprocess_filt_wal(self):
 #
 # self.I2 = (self.I2 - lp)
 
- 
 
  def preprocess_filt_hps(self):
  '''
  Do the pre processing using (orig - low-pass filter) = high-pass filter filter (3.9/5.3 min).
@@ -205,7 205,7 @@ def preprocess_filt_lap(self):
 
 
 
-  
 
  def uniform_data_type(self):
 
  import numpy as np
@@ -342,6 342,7 @@ def autorift(self):
  DispFiltC.FiltWidth = (self.FiltWidth - 1) * ChipSize0_GridSpacing_oversample_ratio   1
  DispFiltC.Iter = self.Iter - 1
  DispFiltC.MadScalar = self.MadScalar
  DispFiltC.colfiltChunkSize = self.colfiltChunkSize
 
  DispFiltF = DISP_FILT()
  overlap_f = 1 - 1 / ChipSize0_GridSpacing_oversample_ratio
@@ -350,6 351,7 @@ def autorift(self):
  DispFiltF.FiltWidth = (self.FiltWidth - 1) * ChipSize0_GridSpacing_oversample_ratio   1
  DispFiltF.Iter = self.Iter
  DispFiltF.MadScalar = self.MadScalar
  DispFiltF.colfiltChunkSize = self.colfiltChunkSize
 
 
  for i in range(ChipSizeUniX.__len__()):
@@ -369,13 371,13 @@ def autorift(self):
  yGrid0 = np.round(yGrid0)
 
  M0 = (ChipSizeX == 0) & (self.ChipSizeMinX <= ChipSizeUniX[i]) & (self.ChipSizeMaxX >= ChipSizeUniX[i])
- M0 = colfilt(M0.copy(), (int(1/Scale*6), int(1/Scale*6)), 0)
  M0 = colfilt(M0.copy(), (int(1/Scale*6), int(1/Scale*6)), 0, self.colfiltChunkSize)
  M0 = cv2.resize(np.logical_not(M0).astype(np.uint8),dstShape[::-1],interpolation=cv2.INTER_NEAREST).astype(np.bool)
 
- SearchLimitX0 = colfilt(self.SearchLimitX.copy(), (int(1/Scale), int(1/Scale)), 0)   colfilt(self.Dx0.copy(), (int(1/Scale), int(1/Scale)), 4)
- SearchLimitY0 = colfilt(self.SearchLimitY.copy(), (int(1/Scale), int(1/Scale)), 0)   colfilt(self.Dy0.copy(), (int(1/Scale), int(1/Scale)), 4)
- Dx00 = colfilt(self.Dx0.copy(), (int(1/Scale), int(1/Scale)), 2)
- Dy00 = colfilt(self.Dy0.copy(), (int(1/Scale), int(1/Scale)), 2)
  SearchLimitX0 = colfilt(self.SearchLimitX.copy(), (int(1/Scale), int(1/Scale)), 0, self.colfiltChunkSize)   colfilt(self.Dx0.copy(), (int(1/Scale), int(1/Scale)), 4, self.colfiltChunkSize)
  SearchLimitY0 = colfilt(self.SearchLimitY.copy(), (int(1/Scale), int(1/Scale)), 0, self.colfiltChunkSize)   colfilt(self.Dy0.copy(), (int(1/Scale), int(1/Scale)), 4, self.colfiltChunkSize)
  Dx00 = colfilt(self.Dx0.copy(), (int(1/Scale), int(1/Scale)), 2, self.colfiltChunkSize)
  Dy00 = colfilt(self.Dy0.copy(), (int(1/Scale), int(1/Scale)), 2, self.colfiltChunkSize)
 
  SearchLimitX0 = np.ceil(cv2.resize(SearchLimitX0,dstShape[::-1]))
  SearchLimitY0 = np.ceil(cv2.resize(SearchLimitY0,dstShape[::-1]))
@@ -421,8 423,8 @@ def autorift(self):
  else:
  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)
  SearchLimitX0C = colfilt(SearchLimitX0.copy(), (int(filtWidth), int(filtWidth)), 0, self.colfiltChunkSize)
  SearchLimitY0C = colfilt(SearchLimitY0.copy(), (int(filtWidth), int(filtWidth)), 0, self.colfiltChunkSize)
  SearchLimitX0C = SearchLimitX0C[rIdxC,cIdxC]
  SearchLimitY0C = SearchLimitY0C[rIdxC,cIdxC]
 
@@ -499,8 501,8 @@ def autorift(self):
  DyF[np.logical_not(M0)] = np.nan
 
  # Light interpolation with median filtered values: DxFM (filtered) and DxF (unfiltered)
- DxFM = colfilt(DxF.copy(), (self.fillFiltWidth, self.fillFiltWidth), 3)
- DyFM = colfilt(DyF.copy(), (self.fillFiltWidth, self.fillFiltWidth), 3)
  DxFM = colfilt(DxF.copy(), (self.fillFiltWidth, self.fillFiltWidth), 3, self.colfiltChunkSize)
  DyFM = colfilt(DyF.copy(), (self.fillFiltWidth, self.fillFiltWidth), 3, self.colfiltChunkSize)
 
  # M0 is mask for original valid estimates, MF is mask for filled ones, MM is mask where filtered ones exist for filling
  MF = np.zeros(M0.shape, dtype=np.bool)
@@ -529,16 531,16 @@ def autorift(self):
 
  # DxF0 (filtered) / Dx (unfiltered) is the result from earlier iterations, DxFM (filtered) / DxF (unfiltered) is that of the current iteration
  # first colfilt nans within 2-by-2 area (otherwise 1 nan will contaminate all 4 points)
- DxF0 = colfilt(Dx.copy(),(int(Scale 1),int(Scale 1)),2)
  DxF0 = colfilt(Dx.copy(),(int(Scale 1),int(Scale 1)),2, self.colfiltChunkSize)
  # then resize to half size using area (similar to averaging) to match the current iteration
  DxF0 = cv2.resize(DxF0,dstShape[::-1],interpolation=cv2.INTER_AREA)
- DyF0 = colfilt(Dy.copy(),(int(Scale 1),int(Scale 1)),2)
  DyF0 = colfilt(Dy.copy(),(int(Scale 1),int(Scale 1)),2, self.colfiltChunkSize)
  DyF0 = cv2.resize(DyF0,dstShape[::-1],interpolation=cv2.INTER_AREA)
 
  # Note this DxFM is almost the same as DxFM (same variable) in the light interpolation (only slightly better); however, only small portion of it will be used later at locations specified by M0 and MF that are determined in the light interpolation. So even without the following two lines, the final Dx and Dy result is still the same.
  # to fill out all of the missing values in DxF
- DxFM = colfilt(DxF.copy(), (5,5), 3)
- DyFM = colfilt(DyF.copy(), (5,5), 3)
  DxFM = colfilt(DxF.copy(), (5,5), 3, self.colfiltChunkSize)
  DyFM = colfilt(DyF.copy(), (5,5), 3, self.colfiltChunkSize)
 
  # fill the current-iteration result with previously determined reliable estimates that are not searched in the current iteration
  idx = np.isnan(DxF) & np.logical_not(np.isnan(DxF0))
@@ -578,7 580,7 @@ def autorift(self):
 
 
 
-
  
  def runAutorift(self):
  '''
  quick processing routine which calls autorift main function (user can define their own way by mimicing the workflow here).
@@ -661,6 663,7 @@ def __init__(self):
  self.FiltWidth = 5
  self.Iter = 3
  self.MadScalar = 4
  self.colfiltChunkSize = 4
  self.BuffDistanceC = 8
  self.CoarseCorCutoff = 0.01
  self.OverSampleRatio = 16
@@ -1346,44 1349,93 @@ def arImgDisp_s(I1, I2, xGrid, yGrid, ChipSizeX, ChipSizeY, SearchLimitX, Search
 
 
 
-def colfilt(A, kernelSize, option):
 
 ################## Chunked version of column filter
 def colfilt(A, kernelSize, option, chunkSize=4):
 
  from skimage.util import view_as_windows as viewW
  import numpy as np
 
- A = np.lib.pad(A,((int((kernelSize[0]-1)/2),int((kernelSize[0]-1)/2)),(int((kernelSize[1]-1)/2),int((kernelSize[1]-1)/2))),mode='constant',constant_values=np.nan)
-
- B = viewW(A, kernelSize).reshape(-1,kernelSize[0]*kernelSize[1]).T[:,::1]
-
- output_size = (A.shape[0]-kernelSize[0] 1,A.shape[1]-kernelSize[1] 1)
- C = np.zeros(output_size,dtype=A.dtype)
- if option == 0:# max
- C = np.nanmax(B,axis=0).reshape(output_size)
- elif option == 1:# min
- C = np.nanmin(B,axis=0).reshape(output_size)
- elif option == 2:# mean
- C = np.nanmean(B,axis=0).reshape(output_size)
- elif option == 3:# median
- C = np.nanmedian(B,axis=0).reshape(output_size)
- elif option == 4:# range
- C = np.nanmax(B,axis=0).reshape(output_size) - np.nanmin(B,axis=0).reshape(output_size)
- elif option == 6:# MAD (Median Absolute Deviation)
- m = B.shape[0]
- D = np.abs(B - np.dot(np.ones((m,1),dtype=A.dtype), np.array([np.nanmedian(B,axis=0)])))
- C = np.nanmedian(D,axis=0).reshape(output_size)
- elif option[0] == 5:# displacement distance count with option[1] being the threshold
- m = B.shape[0]
- c = int(np.round((m   1) / 2)-1)
-# c = 0
- D = np.abs(B - np.dot(np.ones((m,1),dtype=A.dtype), np.array([B[c,:]])))
- C = np.sum(D<option[1],axis=0).reshape(output_size)
- else:
- sys.exit('invalid option for columnwise neighborhood filtering')
-
- C = C.astype(A.dtype)
  chunkInds = int(A.shape[1]/chunkSize)
  chunkRem = A.shape[1] - chunkSize * chunkInds
 
  O = 0
 
- return C
  for ii in range(chunkSize):
  startInds = ii*chunkInds
  if ii == chunkSize-1:
  endInds = (ii 1)*chunkInds   chunkRem
  else:
  endInds = (ii 1)*chunkInds
 
  if (ii == 0)&(ii == chunkSize-1):
  A1 = np.lib.pad(A[:,startInds:endInds],((int((kernelSize[0]-1)/2),int((kernelSize[0]-1)/2)),(int((kernelSize[1]-1)/2),int((kernelSize[1]-1)/2))),mode='constant',constant_values=np.nan)
  else:
  if ii == 0:
  A1 = np.lib.pad(A[:,startInds:endInds int((kernelSize[1]-1)/2)],((int((kernelSize[0]-1)/2),int((kernelSize[0]-1)/2)),(int((kernelSize[1]-1)/2),0)),mode='constant',constant_values=np.nan)
  elif ii == chunkSize-1:
  A1 = np.lib.pad(A[:,startInds-int((kernelSize[1]-1)/2):endInds],((int((kernelSize[0]-1)/2),int((kernelSize[0]-1)/2)),(0,int((kernelSize[1]-1)/2))),mode='constant',constant_values=np.nan)
  else:
  A1 = np.lib.pad(A[:,startInds-int((kernelSize[1]-1)/2):endInds int((kernelSize[1]-1)/2)],((int((kernelSize[0]-1)/2),int((kernelSize[0]-1)/2)),(0,0)),mode='constant',constant_values=np.nan)
 
  B = viewW(A1, kernelSize).reshape(-1,kernelSize[0]*kernelSize[1]).T[:,::1]
 
  Adtype = A1.dtype
  Ashape = A1.shape
  del A1
 
  output_size = (Ashape[0]-kernelSize[0] 1,Ashape[1]-kernelSize[1] 1)
  C = np.zeros((B.shape[1],),dtype=Adtype)
 
  if option == 0:# max
  C = np.nanmax(B,axis=0)
  del B
  C = C.reshape(output_size)
  elif option == 1:# min
  C = np.nanmin(B,axis=0)
  del B
  C = C.reshape(output_size)
  elif option == 2:# mean
  C = np.nanmean(B,axis=0)
  del B
  C = C.reshape(output_size)
  elif option == 3:# median
  C = np.nanmedian(B,axis=0, overwrite_input=True)
  del B
  C = C.reshape(output_size)
  elif option == 4:# range
  C = np.nanmax(B,axis=0) - np.nanmin(B,axis=0)
  del B
  C = C.reshape(output_size)
  elif option == 6:# MAD (Median Absolute Deviation)
  m = B.shape[0]
  D = np.zeros((B.shape[1],),dtype=Adtype)
  D = np.nanmedian(B,axis=0)
  D = np.abs(B - np.dot(np.ones((m,1),dtype=Adtype), np.array([D])))
  del B
  C = np.nanmedian(D,axis=0, overwrite_input=True)
  del D
  C = C.reshape(output_size)
  elif option[0] == 5:# displacement distance count with option[1] being the threshold
  m = B.shape[0]
  c = int(np.round((m   1) / 2)-1)
  # c = 0
  D = np.abs(B - np.dot(np.ones((m,1),dtype=Adtype), np.array([B[c,:]])))
  del B
  C = np.sum(D<option[1],axis=0)
  del D
  C = C.reshape(output_size)
  else:
  sys.exit('invalid option for columnwise neighborhood filtering')
 
  C = C.astype(Adtype)
 
  if np.isscalar(O):
  O = C.copy()
  else:
  O = np.append(O,C,axis=1)
 
  return O
 
 
 
@@ -1397,6 1449,8 @@ def __init__(self):
  self.FiltWidth = 5
  self.Iter = 3
  self.MadScalar = 4
  self.colfiltChunkSize = 4
 
 
  def filtDisp(self, Dx, Dy, SearchLimitX, SearchLimitY, M, OverSampleRatio):
 
@@ -1419,21 1473,22 @@ def filtDisp(self, Dx, Dy, SearchLimitX, SearchLimitY, M, OverSampleRatio):
  for i in range(self.Iter):
  Dx[np.logical_not(M)] = np.nan
  Dy[np.logical_not(M)] = np.nan
- M = (colfilt(Dx.copy(), (self.FiltWidth, self.FiltWidth), (5,self.FracSearch)) >= dToleranceX) & (colfilt(Dy.copy(), (self.FiltWidth, self.FiltWidth), (5,self.FracSearch)) >= dToleranceY)
- 
  M = (colfilt(Dx.copy(), (self.FiltWidth, self.FiltWidth), (5,self.FracSearch), self.colfiltChunkSize) >= dToleranceX) & (colfilt(Dy.copy(), (self.FiltWidth, self.FiltWidth), (5,self.FracSearch), self.colfiltChunkSize) >= dToleranceY)
 
 # if self.Iter == 3:
 # pdb.set_trace()
 
  for i in range(np.max([self.Iter-1,1])):
  Dx[np.logical_not(M)] = np.nan
  Dy[np.logical_not(M)] = np.nan
 
- DxMad = colfilt(Dx.copy(), (self.FiltWidth, self.FiltWidth), 6)
- DyMad = colfilt(Dy.copy(), (self.FiltWidth, self.FiltWidth), 6)
-
- DxM = colfilt(Dx.copy(), (self.FiltWidth, self.FiltWidth), 3)
- DyM = colfilt(Dy.copy(), (self.FiltWidth, self.FiltWidth), 3)
  DxMad = colfilt(Dx.copy(), (self.FiltWidth, self.FiltWidth), 6, self.colfiltChunkSize)
  DyMad = colfilt(Dy.copy(), (self.FiltWidth, self.FiltWidth), 6, self.colfiltChunkSize)
 
  DxM = colfilt(Dx.copy(), (self.FiltWidth, self.FiltWidth), 3, self.colfiltChunkSize)
  DyM = colfilt(Dy.copy(), (self.FiltWidth, self.FiltWidth), 3, self.colfiltChunkSize)
 
 
  M = (np.abs(Dx - DxM) <= np.maximum(self.MadScalar * DxMad, DxMadmin)) & (np.abs(Dy - DyM) <= np.maximum(self.MadScalar * DyMad, DyMadmin)) & M
 
  return M
@@ -1458,5 1513,3 @@ def bwareaopen(image,size1):
  return image
 
 
-
-
diff --git a/geo_autoRIFT/autoRIFT/autoRIFT_ISCE.py b/geo_autoRIFT/autoRIFT/autoRIFT_ISCE.py
@@ -165,6 165,13 @@
  mandatory=False,
  doc = 'Mad Scalar')
 
 COLFILT_CHUNK_SIZE = Component.Parameter('colfiltChunkSize',
  public_name = 'COLFILT_CHUNK_SIZE',
  default = 4,
  type = int,
  mandatory=False,
  doc = 'column filter chunk size')
 
 BUFF_DISTANCE_C = Component.Parameter('BuffDistanceC',
  public_name = 'BUFF_DISTANCE_C',
  default = 8,
@@ -235,6 242,7 @@ class autoRIFT_ISCE(autoRIFT, Component):
  FILT_WIDTH,
  ITER,
  MAD_SCALAR,
  COLFILT_CHUNK_SIZE,
  BUFF_DISTANCE_C,
  COARSE_COR_CUTOFF,
  OVER_SAMPLE_RATIO,