We often process LiDAR in tiles for two reasons: first, to keep the number of points per file low and use main memory efficient, and second, to speed up the computation with parallel tile processing and keep all cores of a modern CPU busy. However, it is very (!!!) important to take the necessary precautions to avoid „edge artifacts“ when processing LiDAR in tiles. We have to include points from neighboring tiles during certain LAStools processing steps to avoid edge artifacts. Why? Here is an illustration from our PHIL LiDAR tour earlier this year:
What you see is the temporary TIN of ground points created (internally) by las2dem or blast2dem that is then rastered at the user-specified step size onto a grid. Without a buffer (right side) there will not always be a triangle to cover every pixel. Especially in the corners of the tile you will often find empty pixels. Furthermore the poorly shaped „sliver triangles“ along the boundary of the TIN do not interpolate the ground elevations properly. In contrast, with buffer (left side) the TIN generously covers the entire area that is to be rastered with nicely shaped triangles.
Here the christmas cookies analogy: You need to roll out the dough larger than the cookies you will cut to make sure your cookies will have nice edges. Think of the TIN as the dough and the square tile as your cookie cutter. You need to use a sufficiently large buffer when you roll out your TIN to assure an edge without crumbles when you cut out the tile … (-: … otherwise you are pretty much guaranteed to get results that – upon closer inspection – have these kind of artifacts:
Without buffers processing artifacts also happen when classifying points with lasground or lasclassify, when calculating height above ground or height-normalizing LiDAR tiles with lasheight, when removing noise with lasnoise, when creationg contours with las2iso or blast2iso, or any other operation where an incomplete neighborhood of points can affect the results. Hence, we need to surround each tile with a temporary buffer of points. Currently there are two ways of working with buffers with LAStools:
- creating buffered tiles during the initial tiling step with the ‚-buffer 25‘ option of lastile, maintaining buffered tiles throughout processing and finally using the ‚-use_tile_bb‘ option of lasgrid, las2dem, blast2dem, or lascanopy to raster the tiles without the temporary buffer.
- creating buffered tiles from non-overlapping (= unbuffered) tiles with „on-the-fly“ buffering using the ‚-buffered 25‘ option of most LAStools such as lasground, lasheight, or las2dem. For some workflows it is useful to also add ‚-remain_buffered‘ if buffers are needed again in the next step. Finally, we use the ‚-use_orig_bb‘ option of lasgrid, las2dem, blast2dem, or lascanopy to raster the tiles without the temporary buffer.
In the following three (tiny) examples using the venerable ‚fusa.laz‘ sample that is distributed with LAStools to illustrate the two types of buffering as well as to show what happens when no buffers are used. In each example we will first cut the small ‚fusa.laz‘ sample into nine smaller tiles and then process these separately on 4 cores in parallel.
1. Initial buffer creation with lastile
This is what most of my tutorials teach. It assumes you are the one creating the tiling in the first place. If you do it with lastile and add a buffer right from the start things are pretty easy.
lastile -i ..\data\fusa.laz ^
-set_classification 0 -set_user_data 0 ^
-tile_size 100 -buffer 20 ^
-odir 1_raw -o futi.laz
We cut the input into 100 meter by 100 meter tiles but add a 20 meter buffer around each tile. That means that each tile on disk will contain the points for an area of up to 140 meter by 140 meter. The GUI for LAStools shows the overlap and if you scrutinize the bounding box values that the cursor points to you notice the extra 20 meters in each direction.
Now we can forget about the buffers and run the standard workflow consiting of lasground, lasheight, and lasclassify to distinguish ground, vegetation, and building points in the LiDAR tiles.
lasground -i 1_raw\futi*.laz ^ -city ^ -odir 1_ground -olaz ^ -cores 4 lasheight -i 1_ground\futi*.laz ^ -drop_above 50 ^ -odir 1_height -olaz ^ -cores 4 lasclassify -i 1_height\futi*.laz ^ -odir 1_classify -olaz ^ -cores 4
At the end – when we generate raster products – we have to remember that the tiles were buffered by lastile and cut off the buffers when we raster the TIN with option ‚-use_tile_bb‘ of las2dem.
las2dem -i 1_classify\futi*.laz ^
-keep_class 2 6 ^
-step 0.25 -use_tile_bb ^
-odir 1_dbm -obil ^
-cores 4
We created a digital terrain model with buildings (DBM) by keeping the points with classification 2 (ground) and 6 (building). After loading the resulting 9 tiles into QGIS and generating a virtual raster we see a nice seamless DBM without any edge artifacts.
If you need to deliver the LiDAR files you should remove the buffers with lastile and option ‚-remove_buffer‘.
lastile -i 1_classify\futi*.laz ^
-remove_buffer ^
-odir 1_final -olaz ^
-cores 4
2. On-the-fly buffering
Now assume you are given LiDAR tiles without buffers. We generate them here with lastile.
lastile -i ..\data\fusa.laz ^ -set_classification 0 -set_user_data 0 ^ -tile_size 100 ^ -odir 2_raw -o futi.laz
The only difference is that we do not request the 20 meter buffer and the result is a typical tiling as you may receive it from a vendor or download it from a LiDAR portal. The GUI for LAStools shows that there is no overlap and if you scrutinize the bounding box values that the cursor points to, you see that the tiles is exactly 100 meters by 100 meters.
Now we have to think about buffers a lot. When using on-the-fly buffering we should first spatially index the tiles with lasindex for faster access to the points from neighbouring tiles.
lasindex -i 1_raw\futi*.laz -cores 4
Below in red are the modifications for on-the-fly buffering to the standard workflow of lasground, lasheight, and lasclassify. The first lasground run uses ‚-buffered 20‘ to add buffers to each tile and ‚-remain_buffered‘ to write those buffers to disk. This way they do not have to created again by lasheight and lasclassify.
lasground -i 2_raw\futi*.laz ^ -buffered 20 -remain_buffered ^ -city ^ -odir 2_ground -olaz ^ -cores 4 lasheight -i 2_ground\futi*.laz ^ -remain_buffered ^ -drop_above 50 ^ -odir 2_height -olaz ^ -cores 4 lasclassify -i 2_height\futi*.laz ^ -remain_buffered ^ -odir 2_classify -olaz ^ -cores 4
At the end we have to remember that the tiles still have on-the-fly buffers and them cut off with option ‚-use_orig_bb‘ of las2dem.
las2dem -i 2_classify\futi*.laz ^
-keep_class 2 6 ^
-step 0.25 -use_orig_bb ^
-odir 2_dbm -obil ^
-cores 4
Again, we created a digital terrain model with buildings (DBM) by keeping the points with classification 2 (ground) and 6 (building). The resulting hillshade computed from a virtual raster that combines the 9 BIL rastera into one looks perfectly smooth in QGIS.
If you need to deliver the LiDAR files you should probably remove the buffers first … (-:
lastile -i 2_classify\futi*.laz ^ -remove_buffer ^ -odir 2_final -olaz ^ -cores 4
3. Bad: No buffering
Here what you are *not* supposed to do. Assuming you get unbuffered tiles.
lastile -i ..\data\fusa.laz ^ -set_classification 0 -set_user_data 0 ^ -tile_size 100 ^ -odir 3_raw -o futi.laz
Bad. You do not take care about buffering when processing the tiles.
lasground -i 3_raw\futi*.laz ^ -city ^ -odir 3_ground -olaz ^ -cores 4 lasheight -i 3_ground\futi*.laz ^ -drop_above 50 ^ -odir 3_height -olaz ^ -cores 4 lasclassify -i 3_height\futi*.laz ^ -odir 3_classify -olaz ^ -cores 4
Bad. You do not take care about buffering when generating the DBM.
las2dem -i 3_classify\futi*.laz ^ -keep_class 2 6 ^ -step 0.25 ^ -odir 3_dbm -obil ^ -cores 4
Bad. You get crappy results with edge artifacts clearly visible in the hillshade.
Bad. If you zoom in on a corner where 4 tiles meet you find missing pixels and incorrect elevation values. Bad. Bad. Bad. So please folks. Try this on your own data. Notice the horrible edge artifacts. Then always use buffers … (-:
PS: Usually no buffers are needed for running lasgrid, lasoverlap, or lascanopy as they perform simple binning operations that do not make use of neighbour information.
Nice post something I can use definitely on flood plains right now but tried using blast2dem can’t seem to manage to work with the option -use_orig_bb when used instead of las2dem? But somehow works with -use_tile_bb though
Nice, thanks a lot.
set PATH=%PATH%;C:\LAStools\bin;
rmdir .\fusa\1_raw /s /q
mkdir .\fusa\1_raw
lastile -i c:\lastools\data\fusa.laz ^
-set_classification 0 -set_user_data 0 ^
-tile_size 100 -buffer 20 ^
-odir fusa\1_raw -o futi.laz
rmdir .\fusa\1_ground /s /q
mkdir .\fusa\1_ground
lasground -i fusa\1_raw\futi*.laz ^
-city ^
-odir fusa\1_ground -olaz
rmdir .\fusa\1_height /s /q
mkdir .\fusa\1_height
lasheight -i fusa\1_ground\futi*.laz ^
-drop_above 50 ^
-odir fusa\1_height -olaz
rmdir .\fusa\1_classify /s /q
mkdir .\fusa\1_classify
lasclassify -i fusa\1_height\futi*.laz ^
-odir fusa\1_classify -olaz
rmdir .\fusa\1_dbm /s /q
mkdir .\fusa\1_dbm
las2dem -i fusa\1_classify\futi*.laz ^
-keep_class 2 6 ^
-step 0.25 -use_tile_bb ^
-odir fusa\1_dbm -obil
rmdir .\fusa\1_final /s /q
mkdir .\fusa\1_final
lastile -i fusa\1_classify\futi*.laz ^
-remove_buffer ^
-odir fusa\1_final -olaz
Pingback: Complete LiDAR Processing Pipeline: from raw Flightlines to final Products | rapidlasso GmbH
Pingback: National Open LiDAR Strategy of Latvia humiliates Germany, Austria, and other European “Closed Data” States | rapidlasso GmbH
Pingback: Generating DTM for fluffy Livox MID-40 LiDAR via “median ground” points | rapidlasso GmbH