PDAL vs laspy for Point Cloud Automation

Every Python point-cloud pipeline eventually picks a primary library, and the choice shapes everything built on top of it. PDAL is a C++ processing engine driven by JSON pipeline declarations, with a Python binding. laspy is a pure-Python LAS/LAZ reader and writer that hands you NumPy arrays and has no opinions about what you do with them.

They are not competitors in the way the framing suggests. PDAL is a processing framework; laspy is a file library. The interesting question is not which is better but which belongs at each point in a pipeline — and the answer is usually both, at different stages.

This page compares them on the axes that decide real pipelines: memory behaviour on large clouds, what each can express, how header control works, performance, and the operational cost of each dependency. It assumes the workflows in classifying point clouds with PDAL and Python and point cloud formats and interchange in Python.

Audience and prerequisites. Python 3.10+, and a working knowledge of both libraries’ basic use. The comparison assumes site-scale clouds — tens to hundreds of millions of points — where the differences matter; below ten million points, either library does anything you ask without complaint.

# Prerequisites

# Conceptual architecture

PDAL’s model is a directed pipeline of stages: a reader, some filters, a writer, declared as JSON. The engine handles memory, streaming and parallelism, and many pipelines execute in streaming mode where points flow through in chunks and the whole cloud never exists in memory at once. That is the property that makes PDAL the right tool for a 400-million-point survey on a 16 GB machine.

laspy’s model is a file object and NumPy arrays. laspy.read() loads everything; laspy.open() gives a chunked iterator that lets you stream manually. There is no processing framework — you write the loops — and in exchange you have complete control and no abstraction to fight.

The practical division is that PDAL is better wherever the operation is one of its filters or wherever the data does not fit in memory, and laspy is better wherever you need to do something PDAL does not implement, or to manipulate the header directly.

Figure 1 — The architectural difference, and the one consequence that decides most choices.

# Where PDAL wins

Clouds larger than memory. A streaming pipeline processes a 400-million-point file in a few gigabytes. The equivalent laspy implementation requires writing the chunking, and any neighbourhood operation requires writing the halo handling too — which is where the bugs live.

import json
import subprocess


def classify_stream(src: str, dst: str) -> None:
    """A pipeline that streams: peak memory is one chunk, not the cloud."""
    pipeline = {"pipeline": [
        src,
        {"type": "filters.assign", "value": "Classification = 0"},
        {"type": "filters.smrf", "window": 18.0, "threshold": 0.45},
        {"type": "writers.las", "filename": dst, "compression": "laszip",
         "forward": "all"},
    ]}
    subprocess.run(["pdal", "pipeline", "--stdin", "--stream"],
                   input=json.dumps(pipeline), text=True, check=True)

Operations that are already implemented. Ground classification, outlier removal, height above ground, covariance features, gridding, reprojection, COPC writing — each is a stage and a few parameters. Reimplementing any of them correctly in NumPy is days of work with a worse result.

Format breadth. PDAL reads and writes far more than LAS: E57, PLY, text, databases, COPC, EPT. A conversion between two of them is a two-line pipeline.

# Where laspy wins

Header control. Setting scale, offset, version, point format and CRS precisely, or reading a header without touching the points, is direct in laspy and awkward through PDAL’s abstraction. Every check in fixing LAS scale and offset precision loss is a laspy routine for this reason.

import laspy


def header_only(path: str) -> dict:
    """Read header fields without decompressing a single point."""
    with laspy.open(path) as fh:            # opens, does not read
        h = fh.header
        return {"points": h.point_count, "version": str(h.version),
                "point_format": h.point_format.id,
                "scales": h.scales.tolist(), "offsets": h.offsets.tolist(),
                "mins": h.mins.tolist(), "maxs": h.maxs.tolist()}

On a 6 GB LAZ this returns in milliseconds. The equivalent through a processing framework typically involves at least parsing metadata through a subprocess.

Anything PDAL does not implement. A bespoke statistic, a per-point calculation involving an external dataset, an unusual filter — these are a NumPy expression in laspy and either a custom C++ stage or a Python filter stage in PDAL, both of which cost more than the operation does.

Dependency weight. laspy plus lazrs installs with pip into any environment. PDAL wants conda, brings a substantial C++ dependency chain, and is the component most likely to complicate a container image — a real consideration for the workers described in containerising photogrammetry workers with Docker.

Figure 2 — Six axes, and the split is close to even. Which is why most mature pipelines carry both.

# The combination most pipelines converge on

The pattern that works is PDAL for bulk processing, laspy for header inspection and final validation. PDAL classifies, filters, reprojects and writes; laspy opens the result, checks the header, verifies the round trip, and reads metadata for the run record. Neither is asked to do the other’s job.

import json
import subprocess

import laspy


def process_and_verify(src: str, dst: str, expect_epsg: int) -> dict:
    """PDAL does the work; laspy proves it was done correctly."""
    pipeline = {"pipeline": [
        src,
        {"type": "filters.outlier", "method": "statistical", "mean_k": 12},
        {"type": "filters.range", "limits": "Classification![7:7]"},
        {"type": "filters.smrf", "window": 18.0},
        {"type": "writers.las", "filename": dst, "compression": "laszip",
         "forward": "all"},
    ]}
    subprocess.run(["pdal", "pipeline", "--stdin"],
                   input=json.dumps(pipeline), text=True, check=True)

    with laspy.open(dst) as fh:
        h = fh.header
        crs = h.parse_crs()
        problems = []
        if crs is None or crs.to_epsg() != expect_epsg:
            problems.append(f"CRS is {crs}, expected EPSG:{expect_epsg}")
        if float(h.scales.max()) > 0.005:
            problems.append(f"scale {h.scales.tolist()} coarser than 5 mm")
        return {"points": int(h.point_count), "problems": problems}

# Performance, measured rather than assumed

The received wisdom is that PDAL is fast because it is C++ and laspy is slow because it is Python. That is true for per-point work and misleading overall, because most real operations are dominated by decompression and I/O, which both libraries delegate to the same compiled code.

Three regimes show up consistently on site-scale clouds.

Pure reading. Both spend their time in LAZ decompression. laspy with the lazrs backend and PDAL land within about twenty percent of each other, with the difference depending more on chunk size than on library.

Standard filters. PDAL wins decisively, and not by a small factor. A statistical outlier filter over a hundred million points is a matter of minutes in PDAL and — implemented as a KD-tree query per chunk in NumPy — closer to an hour. The gap is not the language; it is that PDAL’s implementation has been tuned for a decade and a hand-rolled one has not.

Bespoke arithmetic. laspy wins, because the PDAL equivalent is either a Python filter stage — which pays the same Python cost plus marshalling overhead — or a custom C++ stage nobody is going to write for one project.

import time
from contextlib import contextmanager


@contextmanager
def timed(label: str, results: dict):
    """Minimal timing helper for comparing two implementations honestly."""
    start = time.perf_counter()
    try:
        yield
    finally:
        results[label] = time.perf_counter() - start


def compare_read(path: str) -> dict:
    """Read the same file both ways and time it, including decompression."""
    import json
    import laspy
    import pdal

    out: dict = {}
    with timed("laspy", out):
        las = laspy.read(path)
        _ = las.x.mean()
    with timed("pdal", out):
        pipe = pdal.Pipeline(json.dumps({"pipeline": [path]}))
        pipe.execute()
        _ = pipe.arrays[0]["X"].mean()
    out["ratio"] = out["pdal"] / out["laspy"]
    return out

The point of measuring on your own data is that the answer depends on the file. A LAZ with many small chunks behaves differently from one with few large ones, and a cloud with twelve extra dimensions costs more to marshal than a bare one. A five-minute benchmark on a representative file settles a question that otherwise recurs in every design discussion.

# Operational considerations that outlast the technical ones

Library choice has consequences beyond performance, and on a long-lived pipeline they usually matter more.

Reproducibility. A PDAL pipeline is JSON — it can be stored in a run manifest, diffed between runs, and replayed exactly. A laspy loop is code, which is reproducible only to the extent that the code is versioned with the run. Teams that care about being able to explain a result a year later find the JSON pipeline a genuine advantage, and the same argument appears in writing a manifest-driven batch runner for ODM.

Debuggability. When something inexplicable happens, laspy is pure Python: put a breakpoint in, look at the array, step through. PDAL’s internals are compiled, and the debugging surface is log verbosity and pipeline serialisation. For a subtle data problem, laspy will find it faster even in a pipeline that otherwise uses PDAL.

Environment stability. PDAL pulls GDAL, PROJ, GEOS and a compiled LAZ implementation. Those are the same dependencies the rest of a geospatial stack needs, so on a conda-managed environment the marginal cost is small — but in a pip-only environment, or a slim container, it is the dominant packaging problem. laspy adds essentially nothing.

Team knowledge. A JSON pipeline is readable by someone who does not write Python; a NumPy loop is not. On a team where operations staff need to adjust a classification parameter without a developer, that asymmetry is worth more than any benchmark.

None of these points to one library. What they point to is being deliberate: choose PDAL as the processing spine because of streaming, breadth and reproducibility, and keep laspy in the toolkit for headers, verification and anything bespoke, rather than drifting into a pipeline where the choice was made by whichever import happened first.

# Parameter deep-dive

# Verification and output inspection

A useful exercise when adopting either library is to prove they agree, on a small file, before trusting a pipeline built on one of them.

import json
import subprocess

import laspy
import numpy as np
import pdal


def cross_check(path: str, tol: float = 1e-6) -> dict:
    """Read the same file both ways and require the arrays to match."""
    las = laspy.read(path)
    a = np.column_stack([las.x, las.y, las.z])

    pipe = pdal.Pipeline(json.dumps({"pipeline": [path]}))
    pipe.execute()
    arr = pipe.arrays[0]
    b = np.column_stack([arr["X"], arr["Y"], arr["Z"]])

    if a.shape != b.shape:
        return {"agree": False, "reason": f"counts differ: {a.shape} vs {b.shape}"}
    order = np.lexsort((a[:, 2], a[:, 1], a[:, 0]))
    order_b = np.lexsort((b[:, 2], b[:, 1], b[:, 0]))
    worst = float(np.abs(a[order] - b[order_b]).max())
    return {"agree": worst < tol, "worst_difference": worst}

Sorting both before comparing is necessary because the two libraries do not guarantee the same point order after a read, and a naive element-wise comparison fails on a file where nothing is wrong.

# A decision procedure

For a new pipeline, four questions settle the choice in about a minute.

Will any single input exceed available memory? If yes, PDAL is the spine, because hand-written chunking with correct halo handling for neighbourhood operations is a multi-week project with a long tail of edge cases. If no, either works and the remaining questions decide.

Is the operation already a PDAL stage? Classification, outlier removal, height above ground, reprojection, gridding, COPC writing and format conversion all are. Reimplementing any of them is a bad trade regardless of preference.

Does the job need precise header control or bespoke per-point arithmetic? Then laspy handles that part, whatever else the pipeline uses. Mixing is normal and costs nothing — they read and write the same files.

Who will maintain it? A JSON pipeline that an operations technician can adjust is worth more on a long-running service than a marginally faster Python loop only its author understands.

Written out as code, the policy is short enough to live in a project’s conventions document, which is where it belongs — the point is that the decision is made once and recorded, rather than re-argued on every new module.

Figure 3 — The comparison resolves into a division of labour.

# Troubleshooting

A PDAL pipeline that ran before now exhausts memory. A stage was added that cannot stream — filters.smrf and most neighbourhood filters are not streamable. Check with pdal info --pipeline-serialization and tile instead.

laspy reads a LAZ file but PDAL cannot. Usually a LAZ variant or a compressor version mismatch. Convert with laspy, then process with PDAL.

The two libraries report different point counts. One is respecting a withheld or classification flag the other is not. Compare the raw counts from the headers rather than from processed arrays.

PDAL will not install alongside an existing GDAL. This is the most common operational cost of PDAL. Use conda-forge for the whole geospatial stack, or containerise it separately.

A laspy loop is far slower than the equivalent PDAL pipeline. Expected, and usually acceptable for a one-off. If it is not, the operation probably exists as a PDAL stage already, and reaching for it is cheaper than optimising the loop.

A pipeline works on one machine and fails on another. Almost always a PROJ data directory difference rather than a library one. Both PDAL and laspy resolve coordinate reference systems through PROJ, and a missing or differently versioned grid file changes the result without changing the code.

# Related

• Benchmarking PDAL Pipelines Against laspy Loops
• Choosing a Library for Streaming Versus In-Memory Work
• Point Cloud Formats and Interchange in Python
• Resolving Memory Errors in PDAL and laspy
• GDAL vs Rasterio for DEM and Orthomosaic Export

← Point Cloud Processing & 3D Deliverables