Using PDAL to clean legacy LiDAR data

The CRREL/UCSV Energy Site (CUES) is a snow research site on Mammoth Mountain, CA:

In February 2011 a Riegl LMS-Z390 LiDAR scanner came online at the site, mounted in a glass enclosure on the central structure. This scanner regularly captures snow depths over a portion of the study area, and stores these snow depths on the CUES web server as .3dd files:

While these data have been used with some research, in general they have been under-utilized, partially due to their file format. .3dd is an old Riegl file format, so old that the company no longer produces binaries to read from the format — the only binaries are available in riscanlib for Windows 32 systems. And because Riegl doesn't release the layout of their file formats to the public, we can't¹ write our own multi-platform reader.

riscanlib, which requires a Riegl website account to download, comes with an example executable named conv2asc.exe that reads a .3dd file and prints the X, Y, Z, and Intensity coordinates to standard output. For now, this is what we're using to extract information from .3dd files. However, there's a couple of hangups:

1. Our trusty ol' Z390 scanner produces a lot of bogus data at (0., 0., 0.) that should be filtered out.
2. .3dd intensities are floats in [0, 1], but other point cloud formats (e.g. las) work best when intensity is scaled from [0, 65536].
3. Our scanner shifts and moves through time, both intentionally when someone makes and adjustment and naturally due to shifts in the ground. We need to apply a rigid transform to our data to ensure all of our scans are in the same reference frame.

We'll tackle parts 1 and 2 in this post, and save part 3 for a later date. If you want to follow along at home with the same data, download it from the CUES website. We'll assume you've already used conv2asc.exe on the 3dd file and redirected the output to txt/original/20161229-1315-30.Z390.frame.txt.

PDAL (with the Python plugin enabled) is the only required software for this exercise (assuming you have a unix-y system that includes sed).

Filtering and scaling

PDAL's text reader expects reasonably well-formed text files. Our conv2asc.exe output looks like this, which isn't good enough:

$ head txt/original/20161229-1315-30.Z390.frame.txt
0.0000, 0.0000, 0.0000, 0.0000
0.0000, 0.0000, 0.0000, 0.0000
0.0000, 0.0000, 0.0000, 0.0000
0.0000, 0.0000, 0.0000, 0.0000
0.0000, 0.0000, 0.0000, 0.0000
0.0000, 0.0000, 0.0000, 0.0000
0.0000, 0.0000, 0.0000, 0.0000
0.0000, 0.0000, 0.0000, 0.0000
0.0000, 0.0000, 0.0000, 0.0000
0.0000, 0.0000, 0.0000, 0.0000

We need to remove those spaces between the comma and the next number, and add a header. This is a simple two-step process:

mkdir -p txt/cleaned
echo "X,Y,Z,Intensity" > txt/cleaned/20161229-1315-30.Z390.frame.txt
sed 's/, /,/g' txt/original/20161229-1315-30.Z390.frame.txt >> txt/cleaned/20161229-1315-30.Z390.frame.txt

We'll use PDAL to do the (0., 0., 0.) filtering and intensity scaling. The filtering can be done out of the box with range filter, but we need to write some code to do the intensity scaling. We'll use the programmable filter to do the scaling work.

Create a file called scripts/scale_intensity.py with the following code:

def scale_intensity(input, output):
    intensity = input["Intensity"]
    intensity = intensity * 65536
    output["Intensity"] = intensity
    return True

This takes our [0, 1] intensity values and scales them to [0, 65536], which happens to be the maximum range of intensity values supported by the las format.²

We then can created a "cleaned" laz file with no origin points³ and scaled intensity values with this one-liner:

mkdir -p laz/cleaned
pdal translate txt/cleaned/20161229-1315-30.Z390.frame.txt laz/cleaned/20161229-1315-30.Z390.frame.laz \
    -f range --filters.range.limits='X![0:0],Y![0:0],Z![0:0]' \
    -f programmable --filters.programmable.script=scripts/scale_intensity.py \
        --filters.programmable.module=scale_intensity \
        --filters.programmable.function=scale_intensity

Conclusion

By combining simple external tools (e.g. sed) with PDAL's powerful filters, we can transform uncomfortable data formats into more usable layouts. I've gisted a Makefile and the scale intensity script that are useful for batch processing:

In a later post I'll discuss how to use transformation matrices to bring all of these data into the same coordinate system.

Footnotes

Back to top