Counting Points per Class with PDAL
TL;DR: For a quick look, pdal info --stats --enumerate Classification tile.laz lists the classes present. For counts, run a reader-only pipeline and np.bincount(arr["Classification"], minlength=256); for files that do not fit in memory, sum bincount over pipeline.iterator(chunk_size=...). Join the counts to ASPRS class names and compute shares for a QA table.
# Context and Motivation
This guide is part of ASPRS Classification Codes. A per-class point count is the most basic classification QA there is, and it answers more questions than it looks like it should. Is anything left unclassified? Did the vendor use classes the specification does not allow? Did ground classification fail on this tile — ground share of 8 percent where neighbours show 45? Did a reclassification step move the points it was meant to, and only those? Every one of those is a glance at a histogram, provided the histogram is easy to produce for one tile or a thousand.
# Prerequisites and Assumptions
- PDAL 2.x with Python bindings; NumPy and pandas.
- LAS files with classification populated.
- A list of classes your specification allows, to flag unexpected ones.
# Step-by-Step Implementation
# Step 1 — List classes present
pdal info --stats --enumerate Classification computes statistics and reports the distinct classification values. It reads every point but needs no code.
# Step 2 — Count with bincount
A reader-only pipeline returns an array; np.bincount(a["Classification"], minlength=256) counts every class in one vectorized call.
# Step 3 — Stream for large files
pipeline.iterator(chunk_size=2_000_000) yields chunks; summing bincount across them gives exact counts in bounded memory.
# Step 4 — Add names and shares
Join counts to a dictionary of ASPRS class names and divide by the total.
# Step 5 — Aggregate a project
Run over every tile, keep per-tile rows for QA and sum for the project total.
# Complete Working Example
"""Per-class counts for tiles and a project, with ASPRS names and anomaly flags."""
from __future__ import annotations
import json
from pathlib import Path
import numpy as np
import pandas as pd
import pdal
ASPRS = {0: "never classified", 1: "unclassified", 2: "ground", 3: "low vegetation",
4: "medium vegetation", 5: "high vegetation", 6: "building", 7: "low noise",
8: "model key point (legacy)", 9: "water", 10: "rail", 11: "road surface",
12: "overlap (legacy)", 13: "wire guard", 14: "wire conductor",
15: "transmission tower", 16: "wire connector", 17: "bridge deck",
18: "high noise", 19: "overhead structure", 20: "ignored ground",
21: "snow", 22: "temporal exclusion"}
ALLOWED = {1, 2, 3, 4, 5, 6, 7, 9, 17, 18}
def class_counts(path: Path, chunk: int = 2_000_000) -> np.ndarray:
p = pdal.Pipeline(json.dumps({"pipeline": [str(path)]}))
counts = np.zeros(256, dtype=np.int64)
for arr in p.iterator(chunk_size=chunk):
counts += np.bincount(arr["Classification"], minlength=256)
return counts
def project_table(tiles: list[Path]) -> pd.DataFrame:
rows = []
for t in tiles:
c = class_counts(t)
total = c.sum()
for cls in np.nonzero(c)[0]:
rows.append({"tile": t.stem, "class": int(cls), "name": ASPRS.get(int(cls), "user-defined"),
"points": int(c[cls]), "share": c[cls] / total,
"allowed": int(cls) in ALLOWED})
return pd.DataFrame(rows)
if __name__ == "__main__":
df = project_table(sorted(Path("tiles").glob("*.laz")))
project = (df.groupby(["class", "name"]).points.sum().reset_index()
.assign(share=lambda d: (d.points / d.points.sum()).round(4)))
print(project.to_string(index=False))
bad = df[~df.allowed].groupby("tile").points.sum()
if len(bad):
print("tiles with disallowed classes:\n", bad)
ground = df[df["class"] == 2].set_index("tile").share
median = ground.median()
suspicious = ground[ground < 0.5 * median]
print(f"median ground share {median:.1%}; suspicious tiles: {list(suspicious.index)}")The streaming iterator requires a streamable pipeline — a bare reader always is — so the same function handles 10-million and 500-million-point files.
# Turning Counts into a QA Report
A project-wide table of counts is useful; a short report built from it is what reviewers read. Four derived views cover most needs.
Project totals by class, with names and shares, as the example prints. This is the headline table for a delivery report and the first thing a client’s QA team compares with their expectations.
Tiles with disallowed classes. Any class outside the specification — class 0, legacy class 12, an undocumented vendor code — listed per tile with counts. These are usually quick to fix with conditional filters.assign or a vendor-code remapping, but only if someone sees them.
Outlier tiles per class. For ground, vegetation and buildings, tiles whose share is far from the median of their neighbours. A tile with 9 percent ground among neighbours at 40 percent almost always has a failed ground run, a large water body or a data gap, and each of those needs a different response.
Change between versions. When a project is reprocessed, a table of per-class differences between the old and new delivery shows at a glance what the reprocessing actually changed — and whether it changed anything it was not supposed to.
Keep the raw per-tile counts as a CSV beside the report so any figure can be recomputed.
# Key Parameter Table
| Method | Reads | Memory | Output |
|---|---|---|---|
pdal info --stats --enumerate Classification |
all points | low | distinct values in JSON |
np.bincount on pipeline.arrays |
all points | whole tile | exact counts |
bincount over iterator(chunk_size) |
all points | one chunk | exact counts |
filters.stats with count |
all points | low | per-value counts in metadata |
laspy chunk_iterator |
all points | one chunk | exact counts without PDAL |
# Verification
- Totals match the header. The sum of all class counts equals the header point count.
- Before/after diffs. After any reclassification step, only the classes you meant to change should differ.
- Cross-tool agreement. On one tile, compare with laspy:
np.bincount(laspy.read(p).classification)must be identical.
# Gotchas and Edge Cases
Legacy formats and 5-bit classes. In PDRF 0–5, PDAL exposes the 5-bit class; the withheld, synthetic and key-point bits are separate dimensions. Counts above class 31 are impossible there.
Class 0 versus class 1. Class 0 (never classified) means no classification was ever attempted; class 1 (unclassified) means a classifier looked and assigned nothing. Many specifications forbid class 0 in deliveries.
Shares mislead across landscapes. A tile over a lake has little ground; a tile over a city has little high vegetation. Compare shares against neighbouring tiles, not against a fixed number.
Withheld points. Decide whether withheld points belong in the counts. For QA of the delivery, count everything; for product statistics, exclude them.
# Frequently Asked Questions
How do I count points by classification in a LAS file?
Read the file with a PDAL pipeline and apply NumPy’s bincount to the Classification array, with minlength 256 so every code has a slot. For large files, sum bincount over chunks from the pipeline iterator.
How do I list the classes present with the PDAL command line?
Run pdal info with the stats flag and enumerate Classification. The output includes the distinct classification values found in the file.
What share of points should be ground?
It depends on land cover. Open agricultural land can exceed 60 percent ground; dense forest or urban centres can be below 20 percent. Compare each tile with its neighbours rather than with a fixed value.
Why does my file contain class 0?
Class 0 means the points were never classified, typically because a processing step was skipped or a tool wrote default values. Most delivery specifications do not allow it.
# Related
- ASPRS Classification Codes — what each code means
- Understanding ASPRS Classification Codes — the full code table
- Flagging Overlap and Withheld Points — flags alongside classes
- Checking Pipeline Output with pdal info --stats — automated output checks
- Evaluating Point Classification with a Confusion Matrix — beyond counts to correctness