Add PointCloud subclass of Vector

dev-environment.yml

@@ -19,6 +19,7 @@ dependencies:
   - pip
 
   # Optional dependencies
+  - laspy
   - scikit-image
 
   # Test dependencies
@@ -49,4 +50,4 @@ dependencies:
   - typing-extensions
 
   - pip:
-    - -e ./
+      - -e ./

geoutils/__init__.py

@@ -2,10 +2,12 @@
 GeoUtils is a Python package for the analysis of geospatial data.
 """
 
-from geoutils import examples, pointcloud, projtools, raster, vector  # noqa
+from geoutils import examples, projtools  # noqa
 from geoutils._config import config  # noqa
-from geoutils.raster import Mask, Raster  # noqa
-from geoutils.vector import Vector  # noqa
+
+from geoutils.raster import Mask, Raster  # noqa isort:skip
+from geoutils.vector import Vector  # noqa isort:skip
+from geoutils.pointcloud import PointCloud  # noqa isort:skip
 
 try:
     from geoutils._version import __version__ as __version__  # noqa

geoutils/interface/gridding.py

@@ -3,8 +3,10 @@
 import warnings
 from typing import Literal
 
+import affine
 import geopandas as gpd
 import numpy as np
+import rasterio as rio
 from scipy.interpolate import griddata
 
 from geoutils._typing import NDArrayNum
@@ -16,7 +18,7 @@ def _grid_pointcloud(
     data_column_name: str = "b1",
     resampling: Literal["nearest", "linear", "cubic"] = "linear",
     dist_nodata_pixel: float = 1.0,
-) -> NDArrayNum:
+) -> tuple[NDArrayNum, affine.Affine]:
     """
     Grid point cloud (possibly irregular coordinates) to raster (regular grid) using delaunay triangles interpolation.
 
@@ -88,4 +90,7 @@ def _grid_pointcloud(
     # Flip Y axis of grid
     aligned_dem = np.flip(aligned_dem, axis=0)
 
-    return aligned_dem
+    # 3/ Derive output transform from input grid
+    transform_from_coords = rio.transform.from_origin(min(grid_coords[0]), max(grid_coords[1]), res_x, res_y)
+
+    return aligned_dem, transform_from_coords

geoutils/interface/raster_point.py

@@ -100,7 +100,7 @@ def _raster_to_pointcloud(
     as_array: bool,
     random_state: int | np.random.Generator | None,
     force_pixel_offset: Literal["center", "ul", "ur", "ll", "lr"],
-) -> NDArrayNum | gu.Vector:
+) -> NDArrayNum | gu.PointCloud:
     """
     Convert a raster to a point cloud. See Raster.to_pointcloud() for details.
     """
@@ -227,15 +227,13 @@ def _raster_to_pointcloud(
     )
 
     if not as_array:
-        points = gu.Vector(
-            gpd.GeoDataFrame(
-                pixel_data.T,
-                columns=all_column_names,
-                geometry=gpd.points_from_xy(x_coords_2, y_coords_2),
-                crs=source_raster.crs,
-            )
+        pc = gpd.GeoDataFrame(
+            pixel_data.T,
+            columns=all_column_names,
+            geometry=gpd.points_from_xy(x_coords_2, y_coords_2),
+            crs=source_raster.crs,
         )
-        return points
+        return gu.PointCloud(pc, data_column=data_column_name)
     else:
         # Merge the coordinates and pixel data an array of N x K
         # This has the downside of converting all the data to the same data type

geoutils/pointcloud/__init__.py

@@ -1 +1 @@
-from geoutils.pointcloud.pointcloud import *  # noqa
+from geoutils.pointcloud.pointcloud import PointCloud  # noqa