[Done] Add pygeos.empty (#381)

CHANGELOG.rst

@@ -9,6 +9,8 @@ Version 0.11 (unreleased)
 
 * Optionally output to a user-specified array (``out`` keyword argument) when constructing
   geometries from ``indices`` (#380).
+* Added ``pygeos.empty`` to create a geometry array pre-filled with None or
+  with empty geometries (#381).
 
 **API Changes**
 

pygeos/creation.py

@@ -16,6 +16,7 @@
     "box",
     "prepare",
     "destroy_prepared",
+    "empty",
 ]
 
 
@@ -503,3 +504,36 @@ def destroy_prepared(geometry, **kwargs):
     prepare
     """
     lib.destroy_prepared(geometry, **kwargs)
+
+
+def empty(shape, geom_type=None, order="C"):
+    """Create a geometry array prefilled with None or with empty geometries.
+
+    Parameters
+    ----------
+    shape : int or tuple of int
+        Shape of the empty array, e.g., ``(2, 3)`` or ``2``.
+    geom_type : pygeos.GeometryType, optional
+        The desired geometry type in case the array should be prefilled
+        with empty geometries. Default ``None``.
+    order : {'C', 'F'}, optional, default: 'C'
+        Whether to store multi-dimensional data in row-major
+        (C-style) or column-major (Fortran-style) order in
+        memory.
+
+    Examples
+    --------
+    >>> empty((2, 3)).tolist()
+    [[None, None, None], [None, None, None]]
+    >>> empty(2, geom_type=GeometryType.POINT).tolist()
+    [<pygeos.Geometry POINT EMPTY>, <pygeos.Geometry POINT EMPTY>]
+    """
+    if geom_type is None:
+        return np.empty(shape, dtype=object, order=order)
+
+    geom_type = GeometryType(geom_type)  # cast int to GeometryType
+    if geom_type is GeometryType.MISSING:
+        return np.empty(shape, dtype=object, order=order)
+
+    fill_value = Geometry(geom_type.name + " EMPTY")
+    return np.full(shape, fill_value, dtype=object, order=order)

pygeos/tests/test_creation.py

@@ -2,6 +2,7 @@
 import pytest
 
 import pygeos
+from pygeos.geometry import GeometryType
 
 from .common import (
     empty_polygon,
@@ -468,3 +469,17 @@ def test_subclass_is_geometry(with_point_in_registry):
 
 def test_subclass_is_valid_input(with_point_in_registry):
     assert pygeos.is_valid_input(Point("POINT (1 1)"))
+
+
+@pytest.mark.parametrize("geom_type", [None, GeometryType.MISSING, -1])
+def test_empty_missing(geom_type):
+    actual = pygeos.empty((2,), geom_type=geom_type)
+    assert pygeos.is_missing(actual).all()
+
+
+@pytest.mark.parametrize("geom_type", range(8))
+def test_empty(geom_type):
+    actual = pygeos.empty((2,), geom_type=geom_type)
+    assert (~pygeos.is_missing(actual)).all()
+    assert pygeos.is_empty(actual).all()
+    assert (pygeos.get_type_id(actual) == geom_type).all()