Skip to content

Commit

Permalink
Speed up lmf with fixed windows size
Browse files Browse the repository at this point in the history
  • Loading branch information
Jean-Romain committed Mar 5, 2024
1 parent d7f27b4 commit 029220e
Show file tree
Hide file tree
Showing 2 changed files with 21 additions and 3 deletions.
1 change: 1 addition & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@ Fix: regression of the `stars` package makes `rasterize_terrain()` extremely slo
New: `catalog_intersects()` support a `SpatExtent`
Fix: `lidR` can fully works without `raster` ans `sp`
Fix: [#742](https://github.com/r-lidar/lidR/discussions/742)
Fix: local maximum filter `lmf()` with a fixed windows is now 20 times faster.

## lidR v4.1.0 (Release date: 2024-01-31)

Expand Down
23 changes: 20 additions & 3 deletions src/LAS.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -358,12 +358,20 @@ NumericVector LAS::compute_range(DataFrame flightlines)

void LAS::filter_local_maxima(NumericVector ws, double min_height, bool circular)
{
char UKN = 0;
char NLM = 1;
char LMX = 2;

bool abort = false;
bool vws = ws.length() > 1;

SpatialIndex tree(las);
Progress pb(npoints, "Local maximum filter: ");

std::vector<char> state(npoints);
std::fill(state.begin(), state.end(), 0);
for (int i = 0 ; i < npoints ; i++) { if (Z[i] < min_height) state[i] = NLM; }

#pragma omp parallel for num_threads(ncpu)
for (unsigned int i = 0 ; i < npoints ; i++)
{
Expand All @@ -373,8 +381,7 @@ void LAS::filter_local_maxima(NumericVector ws, double min_height, bool circular

double hws = (vws) ? ws[i]/2 : ws[0]/2;

if (Z[i] < min_height)
continue;
if (state[i] == NLM) continue;

// Get the points within a windows centered on the current point
std::vector<PointXYZ> pts;
Expand Down Expand Up @@ -403,7 +410,13 @@ void LAS::filter_local_maxima(NumericVector ws, double min_height, bool circular
if(z > zmax)
{
is_lm = false;
break;
}

// If the windows size is fixed, then, if a point is below the central point
// we are sure it is not a LM
if (!vws && z < Z[i])
{
state[pt.id] = NLM;
}

// Found one equal. If this one was already tagged LM we can't have two lm
Expand Down Expand Up @@ -456,6 +469,10 @@ void LAS::filter_local_maxima(NumericVector ws)
else
Rcpp::stop("C++ unexpected internal error in 'filter_local_maxima': invalid windows."); // # nocov

std::vector<char> state(npoints);
std::fill(state.begin(), state.end(), 0);
for (int i = 0 ; i < npoints ; i++) { if (Z[i] < min_height) state[i] = NLM; }

SpatialIndex tree(las, skip);
Progress pb(npoints, "Local maximum filter: ");

Expand Down

0 comments on commit 029220e

Please sign in to comment.