Refs #35058 -- Added support for measured geometries to GDAL Point.

Co-authored-by: Nick Pope <nick@nickpope.me.uk>

django/contrib/gis/gdal/geometries.py

@@ -287,6 +287,22 @@ class OGRGeometry(GDALBase):
         else:
             raise ValueError(f"Input to 'set_3d' must be a boolean, got '{value!r}'.")
 
+    @property
+    def is_measured(self):
+        """Return True if the geometry has M coordinates."""
+        return capi.is_measured(self.ptr)
+
+    def set_measured(self, value):
+        """Set if this geometry has M coordinates."""
+        if value is True:
+            capi.set_measured(self.ptr, 1)
+        elif value is False:
+            capi.set_measured(self.ptr, 0)
+        else:
+            raise ValueError(
+                f"Input to 'set_measured' must be a boolean, got '{value!r}'."
+            )
+
     # #### SpatialReference-related Properties ####
 
     # The SRS property
@@ -386,14 +402,22 @@ class OGRGeometry(GDALBase):
         sz = self.wkb_size
         # Creating the unsigned character buffer, and passing it in by reference.
         buf = (c_ubyte * sz)()
-        capi.to_wkb(self.ptr, byteorder, byref(buf))
+        # For backward compatibility, export old-style 99-402 extended
+        # dimension types when geometry does not have an M dimension.
+        # https://gdal.org/api/vector_c_api.html#_CPPv417OGR_G_ExportToWkb12OGRGeometryH15OGRwkbByteOrderPh
+        to_wkb = capi.to_iso_wkb if self.is_measured else capi.to_wkb
+        to_wkb(self.ptr, byteorder, byref(buf))
         # Returning a buffer of the string at the pointer.
         return memoryview(string_at(buf, sz))
 
     @property
     def wkt(self):
         "Return the WKT representation of the Geometry."
-        return capi.to_wkt(self.ptr, byref(c_char_p()))
+        # For backward compatibility, export old-style 99-402 extended
+        # dimension types when geometry does not have an M dimension.
+        # https://gdal.org/api/vector_c_api.html#_CPPv417OGR_G_ExportToWkt12OGRGeometryHPPc
+        to_wkt = capi.to_iso_wkt if self.is_measured else capi.to_wkt
+        return to_wkt(self.ptr, byref(c_char_p()))
 
     @property
     def ewkt(self):
@@ -568,11 +592,21 @@ class Point(OGRGeometry):
         if self.is_3d:
             return capi.getz(self.ptr, 0)
 
+    @property
+    def m(self):
+        """Return the M coordinate for this Point."""
+        if self.is_measured:
+            return capi.getm(self.ptr, 0)
+
     @property
     def tuple(self):
         "Return the tuple of this point."
+        if self.is_3d and self.is_measured:
+            return self.x, self.y, self.z, self.m
         if self.is_3d:
-            return (self.x, self.y, self.z)
+            return self.x, self.y, self.z
+        if self.is_measured:
+            return self.x, self.y, self.m
         return self.x, self.y
 
     coords = tuple
@@ -753,7 +787,9 @@ GEO_CLASSES = {
     6: MultiPolygon,
     7: GeometryCollection,
     101: LinearRing,
-    1 + OGRGeomType.wkb25bit: Point,
+    2001: Point,  # POINT M
+    3001: Point,  # POINT ZM
+    1 + OGRGeomType.wkb25bit: Point,  # POINT Z
     2 + OGRGeomType.wkb25bit: LineString,
     3 + OGRGeomType.wkb25bit: Polygon,
     4 + OGRGeomType.wkb25bit: MultiPoint,

django/contrib/gis/gdal/prototypes/geom.py

@@ -51,6 +51,7 @@ to_kml = string_output(
 getx = pnt_func(lgdal.OGR_G_GetX)
 gety = pnt_func(lgdal.OGR_G_GetY)
 getz = pnt_func(lgdal.OGR_G_GetZ)
+getm = pnt_func(lgdal.OGR_G_GetM)
 
 # Geometry creation routines.
 if GDAL_VERSION >= (3, 3):
@@ -82,6 +83,8 @@ geom_sym_diff = geom_output(lgdal.OGR_G_SymmetricDifference, [c_void_p, c_void_p
 geom_union = geom_output(lgdal.OGR_G_Union, [c_void_p, c_void_p])
 is_3d = bool_output(lgdal.OGR_G_Is3D, [c_void_p])
 set_3d = void_output(lgdal.OGR_G_Set3D, [c_void_p, c_int], errcheck=False)
+is_measured = bool_output(lgdal.OGR_G_IsMeasured, [c_void_p])
+set_measured = void_output(lgdal.OGR_G_SetMeasured, [c_void_p, c_int], errcheck=False)
 
 # Geometry modification routines.
 add_geom = void_output(lgdal.OGR_G_AddGeometry, [c_void_p, c_void_p])
@@ -94,9 +97,13 @@ destroy_geom = void_output(lgdal.OGR_G_DestroyGeometry, [c_void_p], errcheck=Fal
 to_wkb = void_output(
     lgdal.OGR_G_ExportToWkb, None, errcheck=True
 )  # special handling for WKB.
+to_iso_wkb = void_output(lgdal.OGR_G_ExportToIsoWkb, None, errcheck=True)
 to_wkt = string_output(
     lgdal.OGR_G_ExportToWkt, [c_void_p, POINTER(c_char_p)], decoding="ascii"
 )
+to_iso_wkt = string_output(
+    lgdal.OGR_G_ExportToIsoWkt, [c_void_p, POINTER(c_char_p)], decoding="ascii"
+)
 to_gml = string_output(
     lgdal.OGR_G_ExportToGML, [c_void_p], str_result=True, decoding="ascii"
 )
@@ -115,7 +122,7 @@ get_geom_srs = srs_output(lgdal.OGR_G_GetSpatialReference, [c_void_p])
 get_area = double_output(lgdal.OGR_G_GetArea, [c_void_p])
 get_centroid = void_output(lgdal.OGR_G_Centroid, [c_void_p, c_void_p])
 get_dims = int_output(lgdal.OGR_G_GetDimension, [c_void_p])
-get_coord_dim = int_output(lgdal.OGR_G_GetCoordinateDimension, [c_void_p])
+get_coord_dim = int_output(lgdal.OGR_G_CoordinateDimension, [c_void_p])
 set_coord_dim = void_output(
     lgdal.OGR_G_SetCoordinateDimension, [c_void_p, c_int], errcheck=False
 )

django/contrib/gis/utils/layermapping.py

@@ -488,6 +488,10 @@ class LayerMapping:
         if necessary (for example if the model field is MultiPolygonField while
         the mapped shapefile only contains Polygons).
         """
+        # Measured geometries are not yet supported by GeoDjango models.
+        if geom.is_measured:
+            geom.set_measured(False)
+
         # Downgrade a 3D geom to a 2D one, if necessary.
         if self.coord_dim == 2 and geom.is_3d:
             geom.set_3d(False)