GeoNet - Channel Extraction from LiDAR

This post is a summary of one of our graduate students, Terry Barrett, who did some excellent work getting GeoNet to run on a sizable (billions of points) LiDAR dataset.

---

GeoNet is an open-source Matlab code for stream-channel extraction of high resolution digital elevation data, developed by Passalacqua et al. at the National Center for Earth-Surface Dynamics at the University of Minnesota.  With the release of version 1.0.1 last September, the code has matured enough to be capable of adaptation to production use.  One of our projects had need of stream-channel extraction for its one billion point, one meter resolution LiDAR DEM, and I was assigned the task of scaling up GeoNet and producing the stream-channel network for the study watershed.

GeoNet is designed to extract the channel network for a single watershed of maximum extent of about 2000x2000 elevation points, which it can do in a reasonable time (a few minutes) if the code is compiled.  When I began using GeoNet v.1.0.1, there were no compiled versions offered for Windows operating systems, but I was able to compile it for both 32 and 64-bit using the C++ compiler included in Visual Studio 2010.  

Development of the production code was done using a 3.75x3.75km subset of the watershed DEM, which produced a code that can: tile a geotiff raster into buffered tiles, sized for GeoNet; process those tiles through GeoNet to determine the channel network; trim the resulting channel-network shapefile tiles to the original tile dimensions; and connect the tiled-channel shapefiles with channels that conform to stream-flow logic. 

Here is the output on the subset DEM (runtime on a 32-bit machine was about two hours; the blue dots are the connecting channel locations):

Here are some representative channel connection locations from that run that show the robustness of the connection logic (green lines are trimmed-tile channels, red dots are identified channel ends, and blue lines are the channels added to provide logical connections; tile boundaries are shown as dashed lines):

 The logic correctly connected the nearest independent channel ends for all locations except one, and the code was fast enough for the available time on the production machine, and so it was then used to run the entire watershed DEM.

Here is the output for the entire watershed (run in two halves, in parallel, on a 64-bit machine; runtime  ~ 36 hours):

 

The connector channels:

There are three reasons why the channel network has a few holes and doesn’t go to the edge of the DEM:

1) GeoNet’s curvature calculation fails on tiles with more than about 20% missing data, which affected both the edge of the DEM and some interior areas with missing data.

2) A few tiles crashed the Dinf mexfile (which calculates contributing area); I examined the tiles, and could find no significant differences between them and nearby tiles that did not crash.

3) A few channel shapefile tiles (all lowland areas, possibly covered by water) crashed the connecting code; they were about 3X the size of the next largest files, and so generated overly-large connection matrices.

The project was successful, insofar as it produced an automatically generated shapefile network over most of the watershed, but there are areas for further improvement and development:  

1) Prevent holes in the network – by debugging/improving the areas of code mentioned above; using a watershed DEM with a larger buffer (perhaps 2km instead of the 10m that the watershed DEM currently has); and/or performing the tiling with another program, such as ERDAS IMAGINE (I modified some share-code, geotiffwrite, to write the geotiff tiles; may still have some bugs).

2) Validate and calibrate the channel network through comparison to networks generated by other methods (manual and automatic) and to ortho-photos of the watershed.  Both GeoNet and my production code have a number of input parameters (all kept at default settings in the above analysis), and my preliminary investigations into the sensitivity of the channel network to changes in these parameters suggest that there is a lot of room for optimization.

3) Further improve the connecting-channel logic – the subset DEM used for development still had one location where the logic produced an unnecessary connection channel (this location has overlapping independent channels on the stitching boundary), and there are areas where unnatural stream loops are created.

4) The final step required to create a channel network most useful for further geoprocessing requires that the connection channels and the tile channels be merged into a single network of poly lines that properly conforms to sub-watershed regions.  This may require “smart” tiling, with considerably more overlap (and processing time) than used by the current production method.