ENH: Added pyproj.crs.CRS.to_3d

docs/advanced_examples.rst

@@ -137,6 +137,31 @@ which transformation operation is selected in the transformation.
     - bounds: (-136.46, 49.0, -60.72, 83.17)
 
 
+Promote CRS to 3D
+-------------------
+
+.. versionadded:: 3.1
+
+
+In PROJ 6+ you neeed to explictly change your CRS to 3D if you have
+2D CRS and you want the ellipsoidal height taken into account.
+
+
+.. code-block:: python
+
+    >>> from pyproj import CRS, Transformer
+    >>> transformer = Transformer.from_crs("EPSG:4326", "EPSG:2056", always_xy=True)
+    >>> transformer.transform(8.37909, 47.01987, 1000)
+    (2671499.8913080636, 1208075.1135782297, 1000.0)
+    >>> transformer_3d = Transformer.from_crs(
+    ...     CRS("EPSG:4326").to_3d(),
+    ...     CRS("EPSG:2056").to_3d(),
+    ...     always_xy=True,
+    ...)
+    >>> transformer_3d.transform(8.37909, 47.01987, 1000)
+    (2671499.8913080636, 1208075.1135782297, 951.4265527743846)
+
+
 Multithreading
 --------------
 

docs/history.rst

@@ -10,7 +10,7 @@ Change Log
 * ENH: Added :meth:`pyproj.transformer.Transformer.to_proj4` (pull #798)
 * ENH: Added authority, accuracy, and allow_ballpark kwargs to :meth:`pyproj.transformer.Transformer.from_crs` (issue #754)
 * ENH: Added support for "AUTH:CODE" input to :meth:`pyproj.transformer.Transformer.from_pipeline` (issue #755)
-
+* ENH: Added :meth:`pyproj.crs.CRS.to_3d` (pull #808)
 
 3.0.1
 -----

pyproj/_crs.pxd

@@ -1,5 +1,10 @@
 include "proj.pxi"
 
+cdef extern from "proj_experimental.h":
+    PJ *proj_crs_promote_to_3D(PJ_CONTEXT *ctx,
+                               const char* crs_3D_name,
+                               const PJ* crs_2D)
+
 
 cdef _get_concatenated_operations(PJ_CONTEXT*context, PJ*concatenated_operation)
 cdef _to_proj4(

pyproj/_crs.pyi

@@ -218,6 +218,7 @@ class _CRS(Base):
     def to_authority(
         self, auth_name: Optional[str] = None, min_confidence: int = 70
     ): ...
+    def to_3d(self, name: Optional[str] = None) -> "_CRS": ...
     @property
     def is_geographic(self) -> bool: ...
     @property

pyproj/_crs.pyx

@@ -2782,6 +2782,43 @@ cdef class _CRS(Base):
 
         return None
 
+    def to_3d(self, name=None):
+        """
+        .. versionadded:: 3.1
+
+        Convert the current CRS to the 3D version if it makes sense.
+
+        New vertical axis attributes:
+          - ellipsoidal height
+          - oriented upwards
+          - metre units
+
+        Parameters
+        ----------
+        name: str, optional
+            CRS name. Defaults to use the name of the original CRS.
+
+        Returns
+        -------
+        _CRS
+        """
+        cdef char* name_3D = NULL
+        if name is not None:
+            b_name = cstrencode(name)
+            name_3D = b_name
+
+        cdef PJ * projobj = proj_crs_promote_to_3D(
+            self.context, name_3D, self.projobj
+        )
+        CRSError.clear()
+        if projobj == NULL:
+            return self
+        try:
+            crs_3d = _CRS(_to_wkt(self.context, projobj))
+        finally:
+            proj_destroy(projobj)
+        return crs_3d
+
     def _is_crs_property(self, property_name, property_types, sub_crs_index=0):
         """
         .. versionadded:: 2.2.0

pyproj/crs/crs.py

@@ -1294,6 +1294,28 @@ def to_authority(self, auth_name: Optional[str] = None, min_confidence: int = 70
             auth_name=auth_name, min_confidence=min_confidence
         )
 
+    def to_3d(self, name: Optional[str] = None) -> "CRS":
+        """
+        .. versionadded:: 3.1
+
+        Convert the current CRS to the 3D version if it makes sense.
+
+        New vertical axis attributes:
+          - ellipsoidal height
+          - oriented upwards
+          - metre units
+
+        Parameters
+        ----------
+        name: str, optional
+            CRS name. Defaults to use the name of the original CRS.
+
+        Returns
+        -------
+        CRS
+        """
+        return CRS(self._crs.to_3d(name=name))
+
     @property
     def is_geographic(self) -> bool:
         """

test/crs/test_crs.py

@@ -1445,3 +1445,28 @@ def test_coordinate_operation__to_proj4__pretty():
     proj_string = operation.to_proj4(pretty=True)
     assert "+proj=pipeline" in proj_string
     assert "\n" in proj_string
+
+
+@pytest.mark.parametrize(
+    "crs_input",
+    [
+        "EPSG:4326",
+        "EPSG:2056",
+    ],
+)
+def test_to_3d(crs_input):
+    crs = CRS(crs_input)
+    assert len(crs.axis_info) == 2
+    crs_3d = crs.to_3d()
+    assert len(crs_3d.axis_info) == 3
+    vert_axis = crs_3d.axis_info[-1]
+    assert vert_axis.name == "Ellipsoidal height"
+    assert vert_axis.unit_name == "metre"
+    assert vert_axis.direction == "up"
+    assert crs_3d.to_3d() == crs_3d
+    assert crs_3d.name == crs.name
+
+
+def test_to_3d__name():
+    crs_3d = CRS("EPSG:2056").to_3d(name="TEST")
+    assert crs_3d.name == "TEST"