refactor: move grounding zone fits to separate module (#39)

* feat: add ICESat GLA12 GZ calculator * feat: add GLA12 grounding zone mask program * refactor: moved DEM and symbolic tools from IS2tk to GZ * refactor: moved subsetting programs from IS2tk to GZ * docs: add subset group for moved tools
tsutterley · May 24, 2024 · f5b3ced · f5b3ced
1 parent 3bcef03
commit f5b3ced
Show file tree

Hide file tree

Showing 108 changed files with 14,409 additions and 767 deletions.
diff --git a/DAC/interp_DAC_ICESat_GLA12.py b/DAC/interp_DAC_ICESat_GLA12.py
@@ -44,6 +44,7 @@
 UPDATE HISTORY:
     Updated 05/2024: use wrapper to importlib for optional dependencies
         fix memory allocation for output 40HZ data
+        use ellipsoid transformation function from pyTMD
     Updated 04/2024: use timescale for temporal operations
     Updated 08/2023: create s3 filesystem when using s3 urls as input
     Updated 12/2022: single implicit import of grounding zone tools
@@ -69,13 +70,12 @@
 
 # attempt imports
 h5py = gz.utilities.import_dependency('h5py')
-is2tk = gz.utilities.import_dependency('icesat2_toolkit')
 netCDF4 = gz.utilities.import_dependency('netCDF4')
 pyproj = gz.utilities.import_dependency('pyproj')
 pyTMD = gz.utilities.import_dependency('pyTMD')
 timescale = gz.utilities.import_dependency('timescale')
 
-# PURPOSE: read ICESat ice sheet HDF5 elevation data (GLAH12) from NSIDC
+# PURPOSE: read ICESat ice sheet HDF5 elevation data (GLAH12)
 # calculate and interpolate the dynamic atmospheric correction
 def interp_DAC_ICESat_GLA12(base_dir, INPUT_FILE,
     OUTPUT_DIRECTORY=None,
@@ -164,8 +164,9 @@ def interp_DAC_ICESat_GLA12(base_dir, INPUT_FILE,
     topex = pyTMD.datum(ellipsoid='TOPEX', units='MKS')
     wgs84 = pyTMD.datum(ellipsoid='WGS84', units='MKS')
     # convert from Topex/Poseidon to WGS84 Ellipsoids
-    lat_40HZ,elev_40HZ = is2tk.spatial.convert_ellipsoid(lat_TPX, elev_TPX,
-        topex.a_axis, topex.flat, wgs84.a_axis, wgs84.flat, eps=1e-12, itmax=10)
+    lat_40HZ,elev_40HZ = pyTMD.spatial.convert_ellipsoid(lat_TPX, elev_TPX,
+        topex.a_axis, topex.flat, wgs84.a_axis, wgs84.flat,
+        eps=1e-12, itmax=10)
 
     # pyproj transformer for converting from input coordinates (EPSG)
     # to model coordinates

diff --git a/DAC/interp_IB_response_ICESat2_ATL06.py b/DAC/interp_IB_response_ICESat2_ATL06.py
@@ -210,7 +210,7 @@ def ncdf_pressure(FILENAMES,VARNAME,TIMENAME,LATNAME,MEAN,OCEAN,AREA):
     # return the sea level pressure anomalies and times
     return (SLP, TPX, latitude, MJD)
 
-# PURPOSE: read ICESat-2 land ice data (ATL06) from NSIDC
+# PURPOSE: read ICESat-2 land ice data (ATL06)
 # calculate and interpolate the instantaneous inverse barometer response
 def interp_IB_response_ICESat2(base_dir, INPUT_FILE, MODEL,
     OUTPUT_DIRECTORY=None,

diff --git a/DAC/interp_IB_response_ICESat2_ATL07.py b/DAC/interp_IB_response_ICESat2_ATL07.py
@@ -204,7 +204,7 @@ def ncdf_pressure(FILENAMES,VARNAME,TIMENAME,LATNAME,MEAN,OCEAN,AREA):
     # return the sea level pressure anomalies and times
     return (SLP, TPX, latitude, MJD)
 
-# PURPOSE: read ICESat-2 sea ice height (ATL07) from NSIDC
+# PURPOSE: read ICESat-2 sea ice height (ATL07)
 # calculate and interpolate the instantaneous inverse barometer response
 def interp_IB_response_ICESat2(base_dir, INPUT_FILE, MODEL,
     OUTPUT_DIRECTORY=None,

diff --git a/DAC/interp_IB_response_ICESat2_ATL11.py b/DAC/interp_IB_response_ICESat2_ATL11.py
@@ -222,7 +222,7 @@ def ncdf_pressure(FILENAMES,VARNAME,TIMENAME,LATNAME,MEAN,OCEAN,INDICES,AREA):
     # return the sea level pressure anomalies, latitudes and times
     return (SLP, TPX, latitude, MJD)
 
-# PURPOSE: read ICESat-2 annual land ice height data (ATL11) from NSIDC
+# PURPOSE: read ICESat-2 annual land ice height data (ATL11)
 # calculate and interpolate the instantaneous inverse barometer response
 def interp_IB_response_ICESat2(base_dir, INPUT_FILE, MODEL,
     OUTPUT_DIRECTORY=None,

diff --git a/DAC/interp_IB_response_ICESat_GLA12.py b/DAC/interp_IB_response_ICESat_GLA12.py
@@ -52,6 +52,7 @@
 UPDATE HISTORY:
     Updated 05/2024: use wrapper to importlib for optional dependencies
         fix memory allocation for output 40HZ data
+        use ellipsoid transformation function from pyTMD
     Updated 04/2024: use timescale for temporal operations
     Updated 08/2023: create s3 filesystem when using s3 urls as input
         use time functions from pyTMD.time
@@ -79,7 +80,6 @@
 
 # attempt imports
 h5py = gz.utilities.import_dependency('h5py')
-is2tk = gz.utilities.import_dependency('icesat2_toolkit')
 netCDF4 = gz.utilities.import_dependency('netCDF4')
 pyproj = gz.utilities.import_dependency('pyproj')
 pyTMD = gz.utilities.import_dependency('pyTMD')
@@ -210,7 +210,7 @@ def ncdf_pressure(FILENAMES,VARNAME,TIMENAME,LATNAME,MEAN,OCEAN,AREA):
     # return the sea level pressure anomalies and times
     return (SLP, TPX, latitude, MJD)
 
-# PURPOSE: read ICESat ice sheet HDF5 elevation data (GLAH12) from NSIDC
+# PURPOSE: read ICESat ice sheet HDF5 elevation data (GLAH12)
 # calculate and interpolate the instantaneous inverse barometer response
 def interp_IB_response_ICESat(base_dir, INPUT_FILE, MODEL,
     OUTPUT_DIRECTORY=None,
@@ -341,8 +341,9 @@ def interp_IB_response_ICESat(base_dir, INPUT_FILE, MODEL,
     topex = pyTMD.datum(ellipsoid='TOPEX', units='MKS')
     wgs84 = pyTMD.datum(ellipsoid='WGS84', units='MKS')
     # convert from Topex/Poseidon to WGS84 Ellipsoids
-    lat_40HZ,elev_40HZ = is2tk.spatial.convert_ellipsoid(lat_TPX, elev_TPX,
-        topex.a_axis, topex.flat, wgs84.a_axis, wgs84.flat, eps=1e-12, itmax=10)
+    lat_40HZ,elev_40HZ = pyTMD.spatial.convert_ellipsoid(lat_TPX, elev_TPX,
+        topex.a_axis, topex.flat, wgs84.a_axis, wgs84.flat,
+        eps=1e-12, itmax=10)
     # colatitude in radians
     theta_40HZ = (90.0 - lat_40HZ)*np.pi/180.0