diff --git a/docs/api/database.rst b/docs/api/database.rst index 80cfab7e3..bce40a89b 100644 --- a/docs/api/database.rst +++ b/docs/api/database.rst @@ -49,3 +49,9 @@ pyproj.database.get_database_metadata --------------------------------------- .. autofunction:: pyproj.database.get_database_metadata + + +pyproj.database.query_geodetic_crs_from_datum +--------------------------------------- + +.. autofunction:: pyproj.database.query_geodetic_crs_from_datum diff --git a/docs/history.rst b/docs/history.rst index 828b6334e..5ce50e146 100644 --- a/docs/history.rst +++ b/docs/history.rst @@ -5,6 +5,7 @@ Latest ------ - DEP: Minimum PROJ version 9.4 (pull #1481) - DEP: Minimum supported Python version 3.11 (pull #1483) +- ENH: Add :meth:`database.query_geodetic_crs_from_datum` (pull #1390) 3.7.1 ------ diff --git a/pyproj/database.pyi b/pyproj/database.pyi index 61bd7c040..f49024047 100644 --- a/pyproj/database.pyi +++ b/pyproj/database.pyi @@ -1,6 +1,7 @@ from typing import NamedTuple from pyproj.aoi import AreaOfInterest, AreaOfUse +from pyproj.crs import CRS from pyproj.enums import PJType class Unit(NamedTuple): @@ -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]: ... diff --git a/pyproj/database.pyx b/pyproj/database.pyx index a4481de2f..e3b3e3839 100644 --- a/pyproj/database.pyx +++ b/pyproj/database.pyx @@ -6,9 +6,10 @@ from collections import namedtuple from libc.stdlib cimport free, malloc from pyproj._compat cimport cstrdecode, cstrencode -from pyproj._context cimport pyproj_context_create +from pyproj._context cimport _clear_proj_error, pyproj_context_create from pyproj.aoi import AreaOfUse +from pyproj.crs import CRS from pyproj.enums import PJType @@ -465,3 +466,101 @@ def get_database_metadata(str key not None): if metadata == NULL: return None return metadata + + +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.8.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] + """ + + if pj_type is not None and not isinstance(pj_type, PJType): + pj_type = PJType.create(pj_type) + + cdef const char* c_crs_type = NULL + if pj_type is None: + pass + elif pj_type is PJType.GEOCENTRIC_CRS: + c_crs_type = b"geocentric" + elif pj_type is PJType.GEOGRAPHIC_2D_CRS: + c_crs_type = b"geographic 2D" + elif pj_type is PJType.GEOGRAPHIC_3D_CRS: + c_crs_type = b"geographic 3D" + else: + raise ValueError("type must be GEOCENTRIC_CRS, GEOGRAPHIC_2D_CRS, GEOGRAPHIC_3D_CRS or None") + + 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) + _clear_proj_error() + + return ret_list diff --git a/pyproj/proj.pxi b/pyproj/proj.pxi index d4db59d80..ebf7ffb27 100644 --- a/pyproj/proj.pxi +++ b/pyproj/proj.pxi @@ -564,3 +564,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) diff --git a/test/test_database.py b/test/test_database.py index 6b2f1628e..3f00f5df0 100644 --- a/test/test_database.py +++ b/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 @@ -268,3 +269,50 @@ 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", "GEOCENTRIC_CRS") + assert len(crss) == 1 + assert crss[0].to_authority()[1] == "6317" + + 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(ValueError): + query_geodetic_crs_from_datum("EPSG", "EPSG", "1116", "invalid string") + + with pytest.raises(TypeError): + query_geodetic_crs_from_datum("EPSG", "EPSG", None) + + with pytest.raises(TypeError): + query_geodetic_crs_from_datum("EPSG", None, "1116")