LiDAR DTM and DSM LiDAR DTM separated

LiDAR data processing

GEOMATRIX UAB has developed an original automated work-flow for processing of "raw" air-born LiDAR datasets (point clouds in LAS format) sampled with ~2 points/m2 density and delivered in low level of thematic classification (only "land", "water" and "other" classes). The solution is designed for automated parallel processing on a cluster of servers and is completely based on open source software (GRASS GIS version 7.0 with libLAS version 1.7). Except of the original LAS datasets, the solution does not require any other external data inputs for masking and clean-up of the results. LiDAR data processing work-flow currently is capable of performing the following operations:

  • Analysis of the original LAS points cloud and extraction of solid ground, water and surface of the objects by using a more flexible classification approach (instead of just taking "first" and "last" returns) which gives better results while separating DTM from DSM;
  • Creation of "land" and "water" binary masks and iterative flattening of water surfaces to cover large gaps in the original data sample and get rid of sampling/extrapolation errors. Elevation gradients along the streams are well preserved.
  • Filtering and cleaning of the point datasets by standard statistical methods and re-combination of LiDAR points into DTM and DSM surfaces (our original approach gives excellent quality);
  • Stitching of neighboring LiDAR point tiles into seamless surfaces (e.g. corresponding to the aerial photo images) for the further interpolation and processing. This is very useful, as the total number of data "pieces" is reduced considerably, thus reducing the number of interpolation errors along the internal edges of tiles;
  • Interpolation of DTM and DSM raster layers at a given resolution (suggested optimum resolution for typical LiDAR samples is 0.5 m);
  • Production and cleaning of DTM vector iso-lines (e.g. with vertical steps of 1 m);
  • Subtraction of DSM offsets from the DTM base, producing layer of surface objects volume;
  • Filtering out low grass vegetation and other "noise" (in the following examples DSM_base pixels with elevation less than 1m are filtered out);
  • Extraction of tree cover maps at different elevation levels (e.g. 1-10 m - "low trees", 10-20 m - "medium-size trees", more than 20 m - "tall trees"). These are flexible thresholds, which can be defined separately for each sample/project, as well as defined dynamically by a given function;
  • Computing statistical products of DSM and DSM_change within a regular grid or a layer of statistical units.

Images on the left side show examples of our LiDAR processing work-flow. Top row - examples of LiDAR point cloud separation and surface interpolation: left - DTM+DSM raster surface (0.5 m resolution) of absolute elevation values and right - DTM raster surface (0.5 m resolution) produced from filtered LiDAR sampling points. Bottom row - DSM-DTM offsets raster surface (volume of surface objects) product displayed on top of a panchromatic aerial photo image. Pixels below 1m (low vegetation, errors, "noise") are removed from the DSM. Middle row - examples of masking and flattening of water surfaces: left - panchromatic aerial photo image of a small lake in Sweden with poor coverage of LiDAR points, right - corrected DEM with masked and flattened water surface.

LiDAR uneven distribution of points on lake surface LiDAR flattened water surface
LiDAR DSM on top pf earial photo

Technology behind the LiDAR processing system is rather simple: it does not include any proprietary software code or custom processing functions, except those provided by a standard GRASS GIS software package. External libLAS software library is used ony for reading and filtering of the original LiDAR points cloud provided in LAS file format. Filtering and processing of LiDAR data is based on standard masking, smoothing, interpolation functions, as well as a powerful raster algebra system provided by the GRASS GIS software. Due to extremely high density and measurement precision of the original sample of points, LiDAR DEM raster datasets by default are produced in 64-bit float data type and have very large file sizes. It is possible to generate datasets of lower data depth on a special request.

Processing algorithm is optimised for speed, completely automated and very flexible. It allows easy manipulation of processing masks, resampling parameters and methods, pixel resolution and data type of the output raster products, as well as intervals and clean-up parameters of vector isolines. Parallel processing of multiple LiDAR sampling units makes the production relatively fast and efficient. It also enables step-wise processing approach, as well as distributed processing of large datasets on distributed productioon systems - even installed as robotic "boxes" within the production centers of our clients - in order to avoid delays and save money on download/delivery of large LiDAR datasets. This could be a valuable asset for the companies dealing with large amounts of LiDAR data on a regular basis, or providing emergency services.

LiDAR forest DSM above 1 m height LiDaR forest DSM above 1 m in 10x10 grid LiDAR forest DSM above 1 m coverage (%) in 100x100 grid

Analysis of interpolated LiDAR raster datasets is based on a broad range of raster algebra functions, resampling statistics, moving window operations, buffering, cross products and other functions of raster manipulations available in the GRASS GIS software. Statistical analysis layers are stored in a GRASS database, allowing easy automated post-processing, combining with external data layers, automated production of digital maps, and even direct broadcasting as OGC-compliant WMS and CSW services. Due to a wide range of possible data analysis options, we did not develop a strictly defined statistical analysis work-flow for the interpolated LiDAR raster datasets. However, algorithms for automated computing of some "standard" derivative products and the main statistical analysis functions were implemented in the main LiDAR data processing work-flow.

A pilot  testing use-case related to extraction of certain ranges of the forest cover, statistical resampling  and assessmentt of the forest cover show some of the typical post-processing and statistical analysis capabilities of our automated LiDAR processing system. Left column of images show complete forest cover displayed on top of panchromatic aerial photo image: top - forest cover above 1 m, middle - forest cover above 10 m and bottom - forest cover above 20 m height. Middle column of images show forest cover resampled into a 10m resolution grid (specifically for fusion with Sentinel 1/2 products) and displayed on top of panchromatic aerial photo image: top - forest cover above 1 m, middle - forest cover above 10 m and bottom - forest cover above 20 m height. Right column of images show forest cover resampled into a 100m resolution grid (specifically for a European-scale statistics): top - forest cover above 1 m, middle - forest cover above 10 m and bottom - forest cover above 20 m height.

LiDAR forest DSM above 10 m height LiDaR forest DSM above 10 m in 10x10 grid LiDAR forest DSM above 10 m coverage (%) in 100x100 grid
LiDAR forest DSM above 20 m height LiDAR forest DSM above 20 m resampled into 10x10 grid LiDAR forest DSM above 20 m coverage (%) in 100x100 grid