From fe6f84bc1ea63caed314f94d2ad8e433095f45fa Mon Sep 17 00:00:00 2001
From: Zhang Yunjun <yunjunz@outlook.com>
Date: Thu, 14 Mar 2024 11:00:15 +0800
Subject: [PATCH] utils0: add `epsg_code2utm_zone()` (#1159)

+ utils.utils0.py: add `epsg_code2utm_zone()` to convert the EPSG code to the UTM zone string.

---------

Co-authored-by: Scott Staniewicz <scott.stanie@gmail.com>
---
 docs/api/attributes.md     |  4 ++--
 src/mintpy/utils/utils0.py | 24 ++++++++++++++++++++++++
 2 files changed, 26 insertions(+), 2 deletions(-)

diff --git a/docs/api/attributes.md b/docs/api/attributes.md
index b034b6584..d611e6918 100644
--- a/docs/api/attributes.md
+++ b/docs/api/attributes.md
@@ -51,8 +51,8 @@ The following attributes vary for each interferogram:
 +  SUBSET_XMIN/XMAX/YMIN/YMAX = start/end column/row number of subset in the original coverage
 +  MODIFICATION_TIME = dataset modification time, exists in ifgramStack.h5 file for 3D dataset, used for "--update" option of unwrap error corrections.
 +  NCORRLOOKS = number of independent looks, as explained in [SNAPHU](https://web.stanford.edu/group/radar/softwareandlinks/sw/snaphu/snaphu.conf.full)
-+  UTM_ZONE = UTM zone, comprises a zone number and a hemisphere, e.g. 11N, 60S, for geocoded file with UTM projection only.
-+  EPSG = EPSG code for coordinate systems, for geocoded files only.
++  UTM_ZONE = [UTM zone](https://docs.up42.com/data/reference/utm#utm-wgs84), comprises a zone number and a hemisphere, e.g. 11N, 60S, for geocoded file with UTM projection only.
++  EPSG = EPSG code for coordinate systems, for geocoded files only. Check [here] for its relationship with UTM zone.
 +  CENTER_INCIDENCE_ANGLE = incidence angle in degrees at the scene center, read from the 2D incidence angle matrix, for isce2 files only.
 
 ### Reference ###
diff --git a/src/mintpy/utils/utils0.py b/src/mintpy/utils/utils0.py
index 1eb0ebf46..7cbd70905 100644
--- a/src/mintpy/utils/utils0.py
+++ b/src/mintpy/utils/utils0.py
@@ -270,6 +270,8 @@ def touch(fname_list, times=None):
 def utm_zone2epsg_code(utm_zone):
     """Convert UTM Zone string to EPSG code.
 
+    Reference: https://docs.up42.com/data/reference/utm#utm-wgs84
+
     Parameters: utm_zone  - str, atr['UTM_ZONE'], comprises a zone number
                             and a hemisphere, e.g. 11N, 36S, etc.
     Returns:    epsg_code - str, EPSG code
@@ -285,6 +287,28 @@ def utm_zone2epsg_code(utm_zone):
     return epsg_code
 
 
+def epsg_code2utm_zone(epsg_code):
+    """Convert EPSG code to UTM Zone string.
+
+    Reference: https://docs.up42.com/data/reference/utm#utm-wgs84
+
+    Parameters: epsg_code - str / int, EPSG code
+    Returns:    utm_zone  - str, atr['UTM_ZONE'], comprises a zone number
+                            and a hemisphere, e.g. 11N, 36S, etc. None for
+                            a EPSG code not in a UTM coordnate system
+    Examples:   utm_zone = epsg_code2utm_zone('32736')
+    """
+    from pyproj import CRS
+
+    # link: https://pyproj4.github.io/pyproj/stable/examples.html#initializing-crs
+    # alternatively, use CRS.from_user_input(epsg_code) or CRS.from_string(f'EPSG:{epsg_code}')
+    crs = CRS.from_epsg(epsg_code)
+    utm_zone = crs.utm_zone
+    if not utm_zone:
+        print(f'WARNING: input EPSG code ({epsg_code}) is NOT a UTM zone, return None and continue.')
+    return utm_zone
+
+
 def to_latlon(infile, x, y):
     """Convert x, y in the projection coordinates of the file to lat/lon in degree.