grid_terrain crashes without error message

[andrew-plowright]

Hi Jean-Romain,

I've encountered an explained bug when using grid_terrain. Running the lines of code below will crash R.

inLAS <- readLAS("dem_bug.laz")
DEM <- grid_terrain(inLAS, res = 1, algorithm = tin(), keep_lowest = FALSE)

I've uploaded the point cloud subset that is causing the issue here.

Some observations:

• Only this small section of the point cloud is causing the problem. It is part of a larger dataset that processes without an issue
• This only happens when the algorithm is tin() (i.e.: knnidw() works)

I've made sure that my lidR and geometry libraries are up-to-date.

Thanks!

[Jean-Romain]

I fixed this issue. Because I always use lasnormalize and never grid_terrain I missed the fact that I forget to filter degenerated ground points in grid_terrain. This is my bad and I fixed it however you should check your data also on your side with lascheck. I do my best to internally manage "special case" but it is better when the user is aware of potential problem before to run any function. Here 30% of your ground points are degenerated. Why? Is it a problem? What to do with these points?

The triangulation crashes with degenerated ground points so internally they are filtered however duplicated points are not a problem for function like grid_metrics. Your metrics will simply have no meaning. For tree_detection it is not a problem but you will find duplicated trees. So my advice is to always check your data to ensure to get meaningful outputs.

lascheck(inLAS)

#> Checking the data
#>  - Checking coordinates... ok
#>  - Checking coordinates type... ok
#>  - Checking attributes type... ok
#>  - Checking ReturnNumber validity... ok
#>  - Checking NumberOfReturns validity... ok
#>  - Checking ReturnNumber vs. NumberOfReturns... ok
#>  - Checking RGB validity... ok
#>  - Checking absence of NAs... ok
#>  - Checking duplicated points...
#>    1955 points are duplicated and share XYZ coordinates with other points
#>  - Checking degenerated ground points...
#>    There were 1955 degenerated ground points. Some X Y Z coordinates were repeated.
#>    There were 4335 degenerated ground points. Some X Y coordinates were repeated but with different Z coordinates.
#>  - Checking attribute population...
#>    'PointSourceID' attribute is not populated.
#>    'ScanDirectionFlag' attribute is not populated.
#>    'EdgeOfFlightline' attribute is not populated.
#> Checking the header
#>  - Checking header completeness... ok
#>  - Checking scale factor validity... ok
#>  - Checking Point Data Format ID validity... ok
#>  - Checking extra bytes attributes validity... ok
#>  - Checking coordinate reference sytem... ok
#> Checking header vs data adequacy
#>  - Checking attributes vs. point format... ok
#>  - Checking header bbox vs. actual content... ok
#>  - Checking header number of points vs. actual content... ok
#>  - Checking header return number vs. actual content... ok
#> Checking preprocessing already done 
#>  - Checking ground classification... yes
#>  - Checking normalization... no
#>  - Checking negative outliers... ok
#>  - Checking flightline classification... no

[andrew-plowright]

Thanks JR. Using lasfilterduplicates solved the problem.

Interestingly, instead of removing duplicates, ADDING points seems to also circumvent the bug for some reason.

I tried adding 10 manually created points just to see what happens.

library(lidR)
library(raster)

# Read LAS dataset
inLAS <- readLAS("dem_bug.laz")

# Add extra points
extraLAS <- readLAS("extra_pts.las")
crs(extraLAS) <- crs(inLAS)
inLAS <- rbind(inLAS, extraLAS)

# Try grid_terrain
DEM <- grid_terrain(inLAS, res = 1, algorithm = tin(), keep_lowest = FALSE)

Oddly enough, the bug doesn't occur now. Even if I didn't use lasfilterduplicates.

I don't know if this is useful information for you or not, but I thought I'd mention it.

[Jean-Romain]

Well actually at this stage it relies on the geometry package that itself relies on the Qhull library. I passed some options to the triangulation methods that may explain that it does not crash but do not trust too much the output of a triangulation of degenerated points. I acknowledge that I don't know the behavior of Qhull. Anyway, I now remove the duplicates before to triangulate