[pdal] PDAL Python + Geopandas

[James Ford]

Found what I was looking for in the current PDAL python bindings source
code.



Pipeline array -> pandas DataFrame:

Adapted from
https://github.com/PDAL/python/blob/main/src/pdal/pipeline.py#L48

pipeline = ...

_ = pipeline.execute()

idx = 0

arr = pipeline.arrays[idx]

df = pd.DataFrame(arr)



Pandas DataFrame / Geopandas GeoDataFrame -> pipeline:

gdf = ...

df = gdf.drop(columns=["geometry"])

arr = df.to_records()

pipeline = pdal.pipeline.Pipeline(arr)

On Thu, 27 Jun 2024 at 09:29, James Ford <irvine.ford at gmail.com> wrote:

> Hello,
>
> I am using the PDAL python bindings to read point clouds into geopandas
> GeoDataFrames and then write them back to disc.
>
> The approach I have works, but is slow. Is there a faster way of doing
> this?
>
> The slowest part is converting from a GeoDataFrame back to a structured
> numpy array.
>
> Code below.
>
> Thanks,
>
> James
>
> import pdal
> import numpy as np
> import pandas as pd
> import geopandas as gpd
>
> input_point_cloud_filepath = "..."
>
> output_point_cloud_filepath = "..."
>
> crs = "..."
>
> #################### Read PC into Memory and convert to GeoDataFrame
> ####################
>
> pipeline_stages = [
>     pdal.Reader.copc(input_point_cloud_filepath),
>     pdal.Filter.hag_nn()
> ]
>
> pipeline = pdal.Pipeline(pipeline_stages)
>
> _ = pipeline.execute()
>
> pointcloud_dtype = pipeline.arrays[0].dtype
>
> pointcloud_df = pipeline.get_dataframe(0)
>
> pointcloud_gdf = gpd.GeoDataFrame(
>     pointcloud_df,
>     geometry=gpd.points_from_xy(pointcloud_df["X"], pointcloud_df["Y"],
> pointcloud_df["Z"]),
>     crs= crs ,
> )
>
> _ = pointcloud_gdf.sindex
>
> ########################## Manipulate PC via geopandas & numpy
> ##########################
>
> # Do stuff here
>
>
> ########################### GeoDataFrame -> numpy -> pipeline
> ###########################
>
> pointcloud_arr = np.array(
>     (
>         pointcloud_gdf
>         .drop(columns=["geometry"])
>         .apply(tuple, axis=1)
>     ),
>     dtype=pointcloud_dtype
> )
>
> output_pipeline = pdal.Filter.stats().pipeline(pointcloud_arr)
>
> output_pipeline |= pdal.Writer.copc(
>     output_point_cloud_filepath,
>     forward="all",
>     a_srs=crs,
> )
>
> _ = output_pipeline.execute()
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/pdal/attachments/20240628/7b2eded7/attachment.htm>