Repository: DHI/mikeio
Pull request #550 (2023-03-06 06:55:54)

Well Known Text for Points
https://www.ogc.org/standard/sfa/

The Shapely library is the Python wrapper for the GEOS library, which is a standard in geospatial analysis.

At this point we don't want to depend upon the shapely library (and it's dependencies), but if the user has it installed, we can easily convert mikeio points to shapely points.


Patch:
diff --git a/mikeio/spatial/geometry.py b/mikeio/spatial/geometry.py
index d8cd6129..350dd1f7 100644
--- a/mikeio/spatial/geometry.py
+++ b/mikeio/spatial/geometry.py
@@ -4,8 +4,9 @@ from collections import namedtuple
 
 BoundingBox = namedtuple("BoundingBox", ["left", "bottom", "right", "top"])
 
+
 class _Geometry(ABC):
-    def __init__(self, projection:str = "NON-UTM") -> None:
+    def __init__(self, projection: str = "NON-UTM") -> None:
         self._projstr = projection
 
     @property
@@ -34,7 +35,6 @@ class _Geometry(ABC):
         pass
 
 
-
 class GeometryUndefined:
     def __repr__(self):
         return "GeometryUndefined()"
@@ -46,13 +46,27 @@ class GeometryPoint2D(_Geometry):
         self.x = x
         self.y = y
 
-    def __repr__(self):
+    def __repr__(self) -> str:
         return f"GeometryPoint2D(x={self.x}, y={self.y})"
 
+    # TODO should we use wkt here
+    # def __str__(self) -> str:
+    #    return self.wkt
+
     @property
-    def ndim(self):
+    def wkt(self) -> str:
+        return f"POINT ({self.x} {self.y})"
+
+    @property
+    def ndim(self) -> int:
+        """Geometry dimension"""
         return 0
 
+    def to_shapely(self):
+        from shapely.geometry import Point
+
+        return Point(self.x, self.y)
+
 
 class GeometryPoint3D(_Geometry):
     def __init__(self, x: float, y: float, z: float, projection=None):
@@ -66,5 +80,14 @@ class GeometryPoint3D(_Geometry):
         return f"GeometryPoint3D(x={self.x}, y={self.y}, z={self.z})"
 
     @property
-    def ndim(self):
+    def wkt(self) -> str:
+        return f"POINT Z ({self.x} {self.y} {self.z})"
+
+    @property
+    def ndim(self) -> int:
         return 0
+
+    def to_shapely(self):
+        from shapely.geometry import Point
+
+        return Point(self.x, self.y, self.z)
diff --git a/pyproject.toml b/pyproject.toml
index 1732ffe4..ed0e0f85 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -48,7 +48,7 @@ dev = ["pytest",
        "netcdf4",
        "matplotlib"]
 
-test = ["pytest", "pytest-cov", "matplotlib!=3.5.0", "xarray","mypy"]
+test = ["pytest", "pytest-cov", "matplotlib!=3.5.0", "xarray","mypy","shapely"]
 
 notebooks= [
             "nbformat",


Test patch:
diff --git a/tests/spatial/test_point.py b/tests/spatial/test_point.py
new file mode 100644
index 00000000..a62974d6
--- /dev/null
+++ b/tests/spatial/test_point.py
@@ -0,0 +1,33 @@
+from mikeio.spatial.geometry import GeometryPoint2D, GeometryPoint3D
+
+# https://www.ogc.org/standard/sfa/
+
+
+def test_point2d_wkt():
+    p = GeometryPoint2D(10, 20)
+    assert p.wkt == "POINT (10 20)"
+
+    p = GeometryPoint2D(x=-5642.5, y=120.1)
+    assert p.wkt == "POINT (-5642.5 120.1)"
+
+
+def test_point3d_wkt():
+    p = GeometryPoint3D(10, 20, 30)
+    assert p.wkt == "POINT Z (10 20 30)"
+
+
+def test_point2d_to_shapely():
+    p = GeometryPoint2D(10, 20)
+    sp = p.to_shapely()
+    assert sp.x == 10
+    assert sp.y == 20
+    assert sp.wkt == p.wkt
+
+
+def test_point3d_to_shapely():
+    p = GeometryPoint3D(10, 20, -1)
+    sp = p.to_shapely()
+    assert sp.x == 10
+    assert sp.y == 20
+    assert sp.z == -1
+    assert sp.wkt == p.wkt