Reclassifying Noise and Overlap Points with PDAL

This page covers the two ASPRS classes that are routinely mishandled and quietly corrupt every downstream statistic — class 7 (low/high noise) and class 12 (overlap) — and the withheld bit that is not a class at all. Getting these right is what makes a density figure mean something, stops a stray return two metres underground from dragging the terrain surface with it, and prevents flight-line overlap from reporting twice the coverage a survey actually achieved.

Why you hit this

Noise and overlap are the classes nobody specifies and everybody assumes. A delivery may arrive with noise already marked, marked as class 1 (unclassified), silently deleted, or flagged with the withheld bit and left in class 2. Each of those needs different handling, and treating them the same means either discarding real measurements or keeping returns that do not correspond to any surface. The overlap case is worse because it is invisible: a cloud with unmarked flight-line overlap looks like a denser survey, and a density acceptance check run over it passes a specification the survey did not meet.

Prerequisites

  • PDAL 2.6+ with Python bindings, plus laspy>=2.5 and numpy>=1.24.
  • LAS 1.4 point format 6 or 7. Formats 0–5 store classification in a 5-bit field with the synthetic, key-point and withheld flags packed into the same byte, which is the origin of most of the confusion below.
  • A cloud in a projected metric CRS, with PointSourceId populated — it identifies the flight line and is what makes overlap detectable.

Step-by-Step

1. Find out what the delivery actually did

Before changing anything, count what is there. The withheld bit and the classification field are separate, and a file can use either, both or neither.

python
import laspy
import numpy as np

las = laspy.read("survey_delivery.laz")
cls = np.asarray(las.classification)
codes, counts = np.unique(cls, return_counts=True)

print("point format:", las.header.point_format.id)
for c, n in zip(codes, counts):
    print(f"  class {c:>3}: {n:>12,}")

# The withheld bit lives outside the classification field.
if hasattr(las, "withheld"):
    w = np.asarray(las.withheld).astype(bool)
    print(f"withheld flag set on {w.sum():,} points "
          f"({100 * w.mean():.2f}%), of which class 7: {(cls[w] == 7).sum():,}")

Three outcomes are common and each means something different. Class 7 present with no withheld bits means the producer classified noise and left it in — good, and you decide whether to drop it. Withheld bits set with everything in class 1 means the producer flagged noise without classifying it — the flag is authoritative and the class field is not. Neither present, in a delivery that has clearly been cleaned, means the noise was deleted, which is a data-loss decision somebody made on your behalf.

Classification and the withheld bit are different fields The classification field says what a point is. The withheld bit says whether it should be used. A point can be class 2 ground and withheld, class 7 noise and not withheld, or any other combination, and code that reads only one of the two fields will act on a subset of what the producer meant. classificationwhat the point IS — ground, building, noise withheld bitwhether the point SHOULD BE USED class 2, not withheld — use it class 2, withheld — do not use it class 7, not withheld — noise, still present class 7, withheld — noise, excluded All four combinations occur in real deliveries, and only reading both fields tells you which one you have In LAS 1.4 formats 6+ they are separate bytes; in formats 0–5 they share one, which is where the confusion started
Two independent statements about the same point. A filter that reads classification alone will happily consume points the producer marked as unusable.

2. Classify noise by height relative to the ground surface

Absolute elevation thresholds do not survive a hilly site. Height above ground does.

python
import json
import pdal

noise = {
    "pipeline": [
        "classified_smrf.laz",
        {"type": "filters.hag_nn"},
        # Below the ground surface by more than the survey's vertical accuracy: noise.
        {"type": "filters.assign",
         "value": ["Classification = 7 WHERE HeightAboveGround < -0.30"]},
        # Absurdly high returns — birds, cloud, multipath off water.
        {"type": "filters.assign",
         "value": ["Classification = 7 WHERE HeightAboveGround > 120.0"]},
        {"type": "writers.las", "filename": "noise_marked.laz", "forward": "all"},
    ]
}
print(pdal.Pipeline(json.dumps(noise)).execute(), "points processed")

The two thresholds are asymmetric on purpose. Below the surface, anything past the survey’s stated vertical accuracy is physically impossible and can be marked confidently. Above it, the ceiling has to clear the tallest real thing in the extent — a 120 m cutoff is right for a low-rise city and wrong for one with towers, so it belongs in the site manifest rather than hard-coded.

Why one sub-surface return costs a whole neighbourhood of terrain A ground filter estimates the terrain as the lower envelope of the cloud, so a single multipath return two metres beneath the true surface becomes the local minimum. The estimated ground is pulled down toward it across the filter's whole window, not just at the offending point. one multipath return, 2 m below the surface true ground, from the real returns Which is why noise has to be marked before the ground filter runs, not cleaned out of its output afterwards Dashed red is what the filter returns once the low point is in scope
The estimated surface bends over the filter's whole window, not at one cell, and removing the point afterwards does not undo it: the surrounding ground points were already rejected for sitting too far above a surface that should never have dipped.

3. Mark flight-line overlap as class 12, do not delete it

Overlap is real measurement collected twice. It should be identifiable and excludable, not discarded.

python
overlap = {
    "pipeline": [
        "noise_marked.laz",
        # Points seen by more than one flight line within a small radius.
        {"type": "filters.radialdensity", "radius": 1.0},
        {"type": "filters.assign",
         "value": ["Classification = 12 WHERE RadialDensity > 18.0 && Classification == 2"]},
        {"type": "writers.las", "filename": "overlap_marked.laz", "forward": "all"},
    ]
}
print(pdal.Pipeline(json.dumps(overlap)).execute(), "points processed")

A density threshold is the crude version and it works where the flight plan is regular. The precise version uses PointSourceId directly: group points into cells, and where a cell contains returns from two or more source ids, the returns from all but the nearest-nadir line are overlap. That is more code and it is the right implementation for a survey whose lines vary in altitude or speed.

python
import numpy as np
import laspy

las = laspy.read("noise_marked.laz")
x, y = np.asarray(las.x), np.asarray(las.y)
src = np.asarray(las.point_source_id)
cls = np.asarray(las.classification)

cell = 2.0
key = ((x // cell).astype(np.int64) << 32) + (y // cell).astype(np.int64)
order = np.argsort(key, kind="stable")
key_s, src_s = key[order], src[order]

# Cells whose returns come from more than one flight line.
boundaries = np.flatnonzero(np.diff(key_s)) + 1
multi = np.zeros(len(key_s), dtype=bool)
for start, end in zip(np.r_[0, boundaries], np.r_[boundaries, len(key_s)]):
    if len(np.unique(src_s[start:end])) > 1:
        multi[start:end] = True

in_overlap = np.zeros(len(key), dtype=bool)
in_overlap[order] = multi
print(f"{in_overlap.sum():,} points in multi-line cells "
      f"({100 * in_overlap.mean():.1f}% of the cloud)")
Why unmarked overlap passes a density specification the survey missed Two adjacent flight lines overlap by about thirty per cent of their swath width. In the overlap strip the point count is the sum of both lines, so a density check over the whole tile reports a figure well above what either line achieved. Marking the overlap as class 12 lets the check measure single-coverage density instead. line 1 — 9 pts/m² line 2 — 9 pts/m² overlap — 18 pts/m² tile mean: 11.7 pts/m² — comfortably above a 10 pts/m² specification single-coverage density: 9 pts/m² — the specification was missed The same cloud, two defensible density figures Class 12 is what lets the acceptance check exclude the double-counted strip and measure what each line actually delivered
Overlap is not an error and it is not extra coverage. Marking it is what lets a density figure answer the question the specification was written about.

4. Set the withheld bit rather than deleting

Deleting is irreversible and the file loses the record that anything was removed. The withheld bit is reversible and self-documenting.

python
import laspy
import numpy as np

las = laspy.read("overlap_marked.laz")
cls = np.asarray(las.classification)

las.withheld = (cls == 7)          # noise excluded by default, still present
las.write("final_flagged.laz")

check = laspy.read("final_flagged.laz")
print("withheld:", int(np.asarray(check.withheld).sum()),
      "| class 7 retained:", int((np.asarray(check.classification) == 7).sum()))

Class 12 deliberately does not get the withheld bit here. Overlap points are valid measurements — a strip-adjustment or an accuracy assessment wants them — so the right default is “identifiable but included”, with the exclusion applied by whichever consumer needs single coverage.

5. Consume the flags correctly downstream

Every consumer has to state its own policy, and PDAL makes that explicit.

python
import json
import pdal

# Terrain: ground only, noise and overlap excluded.
terrain = {"pipeline": [
    "final_flagged.laz",
    {"type": "filters.range", "limits": "Classification[2:2]"},
    {"type": "filters.expression", "expression": "Withheld == 0"},
    {"type": "writers.gdal", "filename": "dtm.tif", "resolution": 1.0,
     "output_type": "idw", "gdaldriver": "GTiff"},
]}

# Density acceptance: single coverage, so overlap is dropped.
density = {"pipeline": [
    "final_flagged.laz",
    {"type": "filters.expression", "expression": "Classification != 12 && Withheld == 0"},
    {"type": "filters.hexbin", "edge_size": 10.0},
]}

for name, spec in (("terrain", terrain), ("density", density)):
    p = pdal.Pipeline(json.dumps(spec))
    print(name, p.execute(), "points")

Expected Output & Verification

A typical delivery after the full pass:

text
point format: 6
  class   1:    4,102,338
  class   2:   28,441,902
  class   5:   11,204,776
  class   6:    9,880,145
  class   7:      142,006
  class  12:    6,318,440
withheld flag set on 142,006 points (0.23%), of which class 7: 142,006
1,204,882 points in multi-line cells (2.0% of the cloud)

Three sanity checks on those numbers. Noise above roughly one per cent of the cloud is not noise, it is a mis-set threshold or an uncorrected trajectory. Overlap materially above the flight plan’s nominal side-lap — thirty per cent is typical — means lines were flown closer than planned or the detection cell is too large. And every class 7 point should carry the withheld bit if that is your policy; a mismatch between the two counts means the flagging step did not run.

python
import laspy
import numpy as np

las = laspy.read("final_flagged.laz")
cls, w = np.asarray(las.classification), np.asarray(las.withheld).astype(bool)

noise_frac = (cls == 7).mean()
assert noise_frac < 0.01, f"noise is {noise_frac:.3%} of the cloud — threshold is wrong"
assert (cls[w] == 7).all(), "withheld set on something that is not class 7"
assert (w[cls == 7]).all(), "class 7 points exist without the withheld bit"
print("noise and overlap flags are internally consistent")

Common Errors

A density check passes and the survey visibly missed its specification. Overlap was never marked, so the tile mean includes the double-counted strips. Mark class 12 and exclude it from the acceptance calculation, as in step 5.

The terrain surface dips sharply in a few places. Sub-surface noise was not marked, and both SMRF and CSF trust the lowest return. Run the height-above-ground noise pass before the ground filter, not after, so the filter never sees the offending points.

filters.assign silently changes nothing. The expression referenced HeightAboveGround without a preceding filters.hag_nn, so the dimension does not exist and the WHERE clause matches nothing. PDAL does not error on an unknown dimension in an assignment expression — check that the dimension is present with filters.info first.

Frequently Asked Questions

Should I delete noise or flag it?

Flag it. Deleting is irreversible, removes the evidence that anything was removed, and makes a re-run with a different threshold impossible without going back to the original delivery. The withheld bit costs nothing and every serious consumer honours it.

Is class 12 overlap or “reserved”?

In the ASPRS LAS 1.4 specification, class 12 is overlap for point formats 0–5 and is deprecated in favour of the dedicated overlap bit for formats 6–10. In practice both are in circulation; write class 12 for compatibility and also set the overlap bit where the format has one.

What noise fraction should I expect?

Well under one per cent for airborne LiDAR over land. Terrestrial scans in wet or reflective environments run higher. A figure above a few per cent almost always means the threshold, not the sensor.

Back to LiDAR Classification and Ground Extraction.