Data Access Viewer: Changes

We’ve made some changes to the back-end processing in the Data Access Viewer (DAV) and they will affect the outputs a bit. I’ll do a review of what’s different here and how the new processing works. Overall, the primary shift has been to a greater reliance on the open source PDAL and GDAL packages. That means you can replicate a lot of the same processing yourself if you want something slightly different. These changes are only to the processing of point clouds. Changes went into affect October 22, 2021.

Contours

Perhaps the biggest change is with contours. DAV has been using a commercial tool that creates the contours from a triangulated irregular network (TIN) and then tries to simplify the lines and make them a little smoother. However, they were still pretty jagged and the simplify flag could result in output that was not topologically correct. There was the potential for contours to cross, which is generally considered “not good”. Although extracting the contours from a TIN will make them truer to the data, you don’t necessarily want every little dip and rise in someone’s lawn represented. Nor do you want the inherent noise in any point cloud to add irrelevant detail.

We’ve taken the approach that contours are meant for human interpretation and smoothness is desirable. Instead of a TIN, we’re going to generate a raster and then extract the contours from it. There are multiple methods of interpolation to create that raster, but we’ve noticed that the rasters from the ANUDEM program, using splines, are a bit smoother than some of the other methods we have (more on those methods below). We’re going to use that method for now, but it may change as we learn more. It’s possible that the interpolation method may become an option for contours, although it will probably need to be an advanced option.

To give you some idea of how different the contours are, below are two sets of contours over the same raster. The raster was generated as just the average value per cell with filling for small gaps. The first set of contours are from a TIN with simplification turned on. The second set is using splines to make a raster and contouring with gdal_contour.

Contours from a TIN with some simplification.

Contours from a raster created with spline interpolation.

Just the reduction of small polygons that probably have little meaning is enormous. The contours from the raster have roughly a third as many features as the ones directly from the TIN. It’s entirely possible that some of the features that were smoothed out of existence had significance for someone’s project. Unfortunately, there is no perfect answer for contours. They are already a simplification of the data, so if you need all the detail, you’re likely better off going to a raster or TIN version of the data.

Raster Interpolation

The DAV systems offers several options to generate a raster from a point cloud. Several of them used commercial tools requiring regular acquisition of a new license, which put the system at risk if we couldn’t get the paperwork done in time. These methods included TIN, average/mean, minimum, and maximum. Luckily, the open source PDAL and GDAL packages were able to step in for all of them. PDAL was already in use for inverse distance weighting (IDW).

Min, Max, and Mean

The page for writers.gdal in the PDAL documentation gives a good example of how PDAL and GDAL are combined to generate a raster from a lidar point cloud and that’s the pattern we followed. For the case where you want to fill small gaps, we use the GDAL command gdal_fillnodata.py with a maximum distance of 5 pixels. If you would like to fill further than that, just take the outputs and run gdal_fillnodata.py with a larger maximum distance.

TIN

The triangulated irregular network was a bit trickier. There is an example in the PDAL documentation, but it isn’t obvious. If you look at the documentation for the writers.raster, the first example shows how to do a TIN. However, it uses a filter called filters.faceraster that doesn’t currently show up in the list of filters on the PDAL page. That filter does have some needed options for setting the bounds of the output. However, we still have a problem with the TIN output, especially when we’re working at the coastline. We want to be able to limit how far it extrapolates so we aren’t making triangles out into a harbor if we don’t have bathy data. At the coast, we’re very likely to have a jagged edge of valid data and the seaward or lakeward areas should be empty instead of having a TIN interpolated value.

I don’t believe there is a simple automated solution to removing the areas you don’t want interpolated while retaining the ones you do using only the data itself. We do want to fill in areas like lakes and the holes left by buildings, but don’t want to extrapolate out to sea. Our previous solution would remove triangles that had a leg longer than a threshold. This kept it from going out to sea, but the threshold small enough to do that also left a lot of holes.

PDAL does not appear to have a feature to remove long triangles, so a different approach was needed. Our current solution is to take the point cloud and derive the concave boundary with PDAL’s hexbin filter. This boundary has polygons describing the areas with data, including having all the holes. A python script using GDAL’s OGR routines deletes interior holes in the hexbin output as we want to allow interpolation over the holes left by buildings. We use the boundary as a cutline for the raster in GDAL’s gdalwarp command, thus clipping out all the areas where we didn’t want interpolation.

Raster interpolation using a TIN in Winyah Bay, SC. A clipping polygon was created based on the points for the entire request. This dataset was made as six output tiles, each of which was clipped with the polygon.

The result is by no means perfect, but it’s pretty good for only using the data to drive it. The example above shows a large job in the Winyah Bay, SC area. Little lakes and building footprints have been filled in, but the bay has not, though it retains the islands. This is what we were looking for. However, you’ll note at the bottom that there are a couple of rivers that aren’t filled for a little bit because their gaps extend to the edge of the area I selected. That’s the trade-off. Gaps that extend to the edge will have their interpolated values clipped out.

I think the next step that would be required to improve the output is to use breaklines to determine the clipping or filling. While those exist for many datasets, the variations in their construction make using them in an automated system challenging. A lot of manual labor will need to go into cataloging and normalizing what we have.

Other Changes

Tile Sizes

Some of the changes above ended up enabling other beneficial changes. Previous products that used a TIN, including the contours, were using a 32-bit program and that had us limited to a tile layout that would keep the number of input points below about 15 million per tile. Using 64-bit PDAL, we aren’t hitting that limit and can make significantly bigger tiles. We haven’t tested the limits, but 120 million points is not a problem and we may go larger than that. That means a lot less tiles for you and fewer potential issues at tile boundaries. Tiles are now targeted to have 10,000 pixels per side.

DXF Format

We originally created DXF format from the shapefiles, created by the contour program, and then modified to add the z-component in each vertex. The converter was old, home-grown, circa 2000 software that never seemed to work well for AutoCAD users. While gdal_contour is still outputting shapefiles (maybe geopackage in the future), we’re now using ogr2ogr to create the DXF files and I’m hoping that will be a better product. Time, and users, will tell.

Future

Incorporating PDAL as more fundamental to the processing stream could allow some additional products to be generated. Already in the plans are potential products that can be used without a GIS system, though they will lack the ability to extract data values. PDAL can also output a mesh in PLY format and that could be translated to formats used by 3D printers, making it an interesting candidate. If you know of other outputs you’d like to see from the system, please let us know. Just recognize that it is automated, so anything that requires someone to put eyes on isn’t feasible.