From cc8da781486190b30a92d783c21ea8acd73b80de Mon Sep 17 00:00:00 2001
From: Simran S Sangha <sssangha@ucla.edu>
Date: Tue, 25 Oct 2022 12:30:12 -0700
Subject: [PATCH] Pass output projection and spacing flexibility

Pushing following fixes/enhancements:
1. By default, the output project was set, but not applied to the output raster. This lead to issues when trying to interact with/visualize products in e.g. a GIS interface. Now, the output project is passed and set.
2. By default, spacing is set to '0.001' degrees, but one could not toggle full-resolution geocoding in a straightforward way. Now if 'None' is specified for '-y' and '-x', a user may geocode in full-resolution.
---
 contrib/stack/topsStack/geocodeGdal.py | 36 ++++++++++++++++++--------
 1 file changed, 25 insertions(+), 11 deletions(-)

diff --git a/contrib/stack/topsStack/geocodeGdal.py b/contrib/stack/topsStack/geocodeGdal.py
index 919369ed..0af5fa96 100755
--- a/contrib/stack/topsStack/geocodeGdal.py
+++ b/contrib/stack/topsStack/geocodeGdal.py
@@ -7,6 +7,7 @@
 import isce
 import isceobj
 import os
+import re
 from osgeo import gdal
 import numpy as np
 import xml.etree.ElementTree as ET
@@ -146,6 +147,14 @@ def runGeo(inps):
        #cmd = 'isce2gis.py envi -i ' + file
        #os.system(cmd)
 
+    # set format and file extenstion
+    if inps.istiff:
+        ext = '.tif'
+        outformat = 'GTiff'
+    else:
+        ext = ''
+        outformat = 'ENVI'
+
 
     if not inps.isAlexGrid:
         WSEN = str(inps.bbox[2]) + ' ' + str(inps.bbox[0]) + ' ' + str(inps.bbox[3]) + ' ' + str(inps.bbox[1])
@@ -162,11 +171,21 @@ def runGeo(inps):
             #os.system(cmd)
             writeVRT(infile, latFile, lonFile)
 
-            cmd = 'gdalwarp -of ENVI -geoloc  -te '+ WSEN + ' -tr ' + str(inps.latStep) + ' ' + str(inps.lonStep) + ' -srcnodata 0 -dstnodata 0 ' + ' -r ' + inps.resamplingMethod + ' ' + infile +'.vrt '+ outFile
+            step_cmd = ' '
+            if re.match(r'^-?\d+(?:\.\d+)$', str(inps.lonStep)) and \
+                 re.match(r'^-?\d+(?:\.\d+)$', str(inps.latStep)):
+                step_cmd = ' -tr ' + str(inps.latStep) + \
+                           ' ' + str(inps.lonStep)
+
+            cmd = 'gdalwarp -of ' + outformat + ' -t_srs ' + inps.outproj + \
+                  ' -geoloc  -te '+ WSEN + step_cmd + \
+                  ' -srcnodata 0 -dstnodata 0 ' + \
+                  ' -r ' + inps.resamplingMethod + ' ' + infile + '.vrt ' + \
+                  outFile + ext
             print (cmd)
             os.system(cmd)
 
-            write_xml(outFile)
+            write_xml(outFile+ext)
 
 
     else:
@@ -174,19 +193,16 @@ def runGeo(inps):
         ylims, xlims = getGridLimits(latfile=latFile, lonfile=lonFile)
 
         WSEN = str(xlim[0]) + ' ' + str(ylim[0]) + ' ' + str(xlim[1]) + ' ' + str(ylim[1])
-        if inps.istiff:
-            ext = '.tif'
-            outformat = 'GTiff'
-        else:
-            ext = '.ant'
-            outformat = 'ENVI'
 
         for infile in inps.prodlist:
             print('geocoding: ' + infile)
 
             writeVRT(infile, latFile, lonFile)
 
-            cmd = 'gdalwarp -of ' + outformat + ' -t_srs ' + inps.outproj + ' -geoloc -te ' + WSEN + ' -tr ' + str(inps.lonStep) + ' ' + str(inps.latStep) + ' -srcnodata 0 -dstnodata 0 -r ' + inps.resamplingMethod + ' ' + infile + '.vrt ' + infile+'ext'
+            cmd = 'gdalwarp -of ' + outformat + ' -t_srs ' + inps.outproj + \
+                  ' -geoloc -te ' + WSEN + ' -tr ' + str(inps.lonStep) + \
+                  ' ' + str(inps.latStep) + ' -srcnodata 0 -dstnodata 0 -r ' + \
+                  inps.resamplingMethod + ' ' + infile + '.vrt ' + infile + ext
             status = os.system(cmd)
             if status:
                 raise Exception('Command {0} Failed'.format(cmd))
@@ -266,5 +282,3 @@ def main(iargs=None):
    
 if __name__ == '__main__':
     main()
-
-