Add query_geodetic_crs_from_datum

Author: jjimenezshaw

pyproj/database.pyi

@@ -1,5 +1,6 @@
 from typing import NamedTuple
 
+from pyproj._crs import _CRS
 from pyproj.aoi import AreaOfInterest, AreaOfUse
 from pyproj.enums import PJType
 
@@ -44,3 +45,9 @@ def query_utm_crs_info(
     contains: bool = False,
 ) -> list[CRSInfo]: ...
 def get_database_metadata(key: str) -> str | None: ...
+def query_geodetic_crs_from_datum(
+    crs_auth_name: str | None,
+    datum_auth_name: str,
+    datum_code: str,
+    pj_type: PJType | None = None,
+) -> list[_CRS]: ...

pyproj/database.pyx

@@ -6,8 +6,13 @@ from collections import namedtuple
 from libc.stdlib cimport free, malloc
 
 from pyproj._compat cimport cstrdecode, cstrencode
-from pyproj._datadir cimport pyproj_context_create, pyproj_context_destroy
+from pyproj._datadir cimport (
+    _clear_proj_error,
+    pyproj_context_create,
+    pyproj_context_destroy,
+)
 
+from pyproj._crs import _CRS
 from pyproj.aoi import AreaOfUse
 from pyproj.enums import PJType
 
@@ -483,3 +488,105 @@ def get_database_metadata(str key not None):
         return metadata
     finally:
         pyproj_context_destroy(context)
+
+
+def query_geodetic_crs_from_datum(
+    str crs_auth_name,
+    str datum_auth_name not None,
+    str datum_code not None,
+    pj_type=None
+    ):
+    """
+    .. versionadded:: 3.7.0
+
+    Return GeodeticCRS that use the specified datum
+
+    See: :c:func:`proj_query_geodetic_crs_from_datum`
+
+    Parameters
+    ----------
+    crs_auth_name: str | None
+        The authority name to filter by (e.g. EPSG, ESRI). None is all.
+    datum_auth_name: str
+        The authority of the datum
+    datum_code: str
+        Datum code
+    pj_type: pyproj.enums.PJType | None, optional
+        The type of object to get the CRSs. Can be PJType.GEOCENTRIC_CRS,
+        PJType.GEOGRAPHIC_3D_CRS, PJType.GEOGRAPHIC_2D_CRS or None for all.
+
+    Returns
+    -------
+    list[_CRS]
+    """
+
+    cdef const char* c_crs_type = NULL
+    cdef bytes b_crs_type
+    cdef str crs_type = ""
+    if pj_type is None:
+        pass
+    elif pj_type == PJType.GEOCENTRIC_CRS:
+        crs_type = "geocentric"
+    elif pj_type == PJType.GEOGRAPHIC_2D_CRS:
+        crs_type = "geographic 2D"
+    elif pj_type == PJType.GEOGRAPHIC_3D_CRS:
+        crs_type = "geographic 3D"
+    else:
+        raise ValueError("type must be GEOCENTRIC_CRS, GEOGRAPHIC_2D_CRS, GEOGRAPHIC_3D_CRS or None")
+
+    if pj_type is not None:
+        b_crs_type = cstrencode(crs_type)
+        c_crs_type = b_crs_type
+
+    cdef const char* c_crs_auth_name = NULL
+    cdef const char* c_datum_auth_name = NULL
+    cdef const char* c_datum_code = NULL
+    cdef bytes b_crs_auth_name
+    cdef bytes b_datum_auth_name
+    cdef bytes b_datum_code
+
+    if crs_auth_name is not None:
+        b_crs_auth_name = cstrencode(crs_auth_name)
+        c_crs_auth_name = b_crs_auth_name
+
+    if datum_auth_name is not None:
+        b_datum_auth_name = cstrencode(datum_auth_name)
+        c_datum_auth_name = b_datum_auth_name
+
+    if datum_code is not None:
+        b_datum_code = cstrencode(datum_code)
+        c_datum_code = b_datum_code
+
+    ret_list = []
+
+    cdef PJ_OBJ_LIST *proj_list = NULL
+    cdef int num_proj_objects = 0
+
+    cdef PJ_CONTEXT* context = pyproj_context_create()
+    proj_list = proj_query_geodetic_crs_from_datum(
+        context,
+        c_crs_auth_name,
+        c_datum_auth_name,
+        c_datum_code,
+        c_crs_type
+    )
+
+    if proj_list != NULL:
+        num_proj_objects = proj_list_get_count(proj_list)
+
+    cdef PJ* proj = NULL
+    try:
+        for iii in range(num_proj_objects):
+            proj = proj_list_get(context, proj_list, iii)
+            ret_list.append(_CRS(proj_as_wkt(context, proj, PJ_WKT2_2019, NULL)))
+            proj_destroy(proj)
+            proj = NULL
+    finally:
+        # If there was an error we have to call proj_destroy
+        # If there was none, calling it on NULL does nothing
+        proj_destroy(proj)
+        proj_list_destroy(proj_list)
+        pyproj_context_destroy(context)
+        _clear_proj_error()
+
+    return ret_list

pyproj/proj.pxi

@@ -549,3 +549,5 @@ cdef extern from "proj.h" nogil:
 
     int proj_is_deprecated(const PJ *obj)
     PJ_OBJ_LIST *proj_get_non_deprecated(PJ_CONTEXT *ctx, const PJ *obj)
+
+    PJ_OBJ_LIST *proj_query_geodetic_crs_from_datum(PJ_CONTEXT *ctx, const char *crs_auth_name, const char *datum_auth_name, const char *datum_code, const char *crs_type)

test/test_database.py

@@ -8,6 +8,7 @@
     get_database_metadata,
     get_units_map,
     query_crs_info,
+    query_geodetic_crs_from_datum,
     query_utm_crs_info,
 )
 from pyproj.enums import PJType
@@ -281,3 +282,43 @@ def test_get_database_metadata():
 
 def test_get_database_metadata__invalid():
     assert get_database_metadata("doesnotexist") is None
+
+
+def test_query_geodetic_crs_from_datum():
+    crss = query_geodetic_crs_from_datum("EPSG", "EPSG", "1116", PJType.GEOCENTRIC_CRS)
+    assert len(crss) == 1
+    assert crss[0].to_authority()[1] == "6317"
+
+    crss = query_geodetic_crs_from_datum(None, "EPSG", "1116")
+    assert len(crss) == 3
+    codes = [x.to_authority()[1] for x in crss]
+    assert "6317" in codes
+    assert "6318" in codes
+    assert "6319" in codes
+
+    crss = query_geodetic_crs_from_datum("EPSG", "EPSG", "6269", None)
+    assert len(crss) == 1
+    assert crss[0].to_authority()[1] == "4269"
+
+    crss = query_geodetic_crs_from_datum(None, "EPSG", "6269")
+    assert len(crss) == 3  # EPSG, ESRI, OGC
+
+
+def test_query_geodetic_crs_from_datum_invalid():
+    crss = query_geodetic_crs_from_datum(None, "EPSG", "11")
+    assert len(crss) == 0
+
+    crss = query_geodetic_crs_from_datum(None, "EPSG", "32632")
+    assert len(crss) == 0
+
+    crss = query_geodetic_crs_from_datum("foo-bar", "EPSG", "6269", None)
+    assert len(crss) == 0
+
+    with pytest.raises(ValueError):
+        query_geodetic_crs_from_datum("EPSG", "EPSG", "1116", PJType.PROJECTED_CRS)
+
+    with pytest.raises(TypeError):
+        query_geodetic_crs_from_datum("EPSG", "EPSG", None)
+
+    with pytest.raises(TypeError):
+        query_geodetic_crs_from_datum("EPSG", None, "1116")