@@ -406,6 +406,76 @@ USE-CASE : Take a LiDAR map, create a mesh from the ground points, split into ti
406406 buffer.close()
407407
408408
409+ Digital Terrain Model (DTM) Creation Example
410+ ................................................................................
411+
412+ The following is a script sample that can be used to create a DTM from a PDAL-
413+ readable pointcloud.
414+
415+ Method:
416+ 1. read point cloud file
417+ 2. remove noise
418+ 3. clean up invalid values
419+ 4. classify ground points using `SMRF <https://pdal.io/en/2.9.2/stages/filters.smrf.html >`__
420+ 5. write with `GDAL writer <https://pdal.io/en/2.9.2/stages/writers.gdal.html >`__
421+
422+ .. note :: If your pointcloud already has ground classified, you can skip all but
423+ the reader and writer and achieve the same result.
424+
425+ .. code-block :: python
426+
427+ import pdal
428+
429+ pc_path = ' https://github.com/PDAL/data/raw/refs/heads/main/autzen/autzen.laz'
430+ out_file = ' autzen_dtm.tif'
431+
432+
433+ # read
434+ reader = pdal.Reader.las(pc_path)
435+
436+ # remove noisy points
437+ lownoise_filter = pdal.Filter.range(
438+ limits = ' Classification![7:7]' , tag = ' lownoise'
439+ )
440+ highnoise_filter = pdal.Filter.range(
441+ limits = ' Classification![18:]' , tag = ' highnoise'
442+ )
443+
444+ # saving incorrectly labeled returns here, some people want this, some don't
445+ prepare_ground = pdal.Filter.assign(
446+ value = [
447+ ' Classification=0' ,
448+ ' ReturnNumber=1 WHERE ReturnNumber < 1' ,
449+ ' NumberOfReturns=1 WHERE NumberOfReturns < 1' ,
450+ ],
451+ tag = ' prepare_ground_classifier' ,
452+ )
453+
454+ # classify ground
455+ smrf_classifier = pdal.Filter.smrf(tag = ' ground_classifier' )
456+
457+ # write with gdal, resolution in feet for autzen
458+ gdal_writer = pdal.Writer.gdal(
459+ filename = out_file,
460+ where = ' Classification == 2' ,
461+ data_type = ' float32' ,
462+ resolution = 10 ,
463+ output_type = ' idw' ,
464+ window_size = 3 ,
465+ pdal_metadata = True ,
466+ )
467+
468+ # collect pdal stages and execute pipline
469+ pipeline = (
470+ reader
471+ | lownoise_filter
472+ | highnoise_filter
473+ | prepare_ground
474+ | smrf_classifier
475+ | gdal_writer
476+ )
477+ pipeline.execute()
478+
409479
410480 .. _`Numpy` : http://www.numpy.org/
411481.. _`schema` : http://www.pdal.io/dimensions.html
0 commit comments