From 3bf2bf0b5feb7d43c252956fe8ae1c829de1b629 Mon Sep 17 00:00:00 2001 From: Jean-Romain Date: Wed, 22 Nov 2023 14:26:08 -0500 Subject: [PATCH] Fix #732 --- NEWS.md | 1 + inst/include/SpatialIndex.h | 4 ++-- inst/include/lidR/GridPartition.h | 11 +++++++-- inst/include/lidR/Octree.h | 10 +++++++-- tests/testthat/test-spatialindex.R | 36 +++++++++++++++--------------- 5 files changed, 38 insertions(+), 24 deletions(-) diff --git a/NEWS.md b/NEWS.md index bb99b150..2c29f770 100644 --- a/NEWS.md +++ b/NEWS.md @@ -3,6 +3,7 @@ If you are viewing this file on CRAN, please check [the latest news on GitHub](h ## lidR v4.0.5 (Release date: ) - Fix: [#726](https://github.com/r-lidar/lidR/issues/726) character palette causes error in plot. +- Fix: [#732](https://github.com/r-lidar/lidR/issues/732) octree spatial indexes is not working properly. lidR now use voxel partition and no longer support octree until the problem is fixed. ## lidR v4.0.4 (Release date: 2023-09-07) diff --git a/inst/include/SpatialIndex.h b/inst/include/SpatialIndex.h index a02e5dd4..dea720b7 100644 --- a/inst/include/SpatialIndex.h +++ b/inst/include/SpatialIndex.h @@ -54,7 +54,7 @@ inline SpatialIndex::SpatialIndex(const Rcpp::S4 las) case GRIDPARTITION: case VOXELPARTITION: grid = GridPartition(las); break; case QUADTREE: quadtree = QuadTree(las); break; - case OCTREE: octree = Octree(las); break; + case OCTREE: Rcpp::stop("Error: octree no longer supported."); break; // # nocov default: Rcpp::stop("Internal error: spatial index code inccorect."); break; // # nocov } } @@ -186,7 +186,7 @@ inline int SpatialIndex::index_selector(const Rcpp::S4 las) const if (sensor == UKN || sensor == ALS) code = GRIDPARTITION; else if (sensor == TLS || sensor == UAV || sensor == DAP) - code = OCTREE; + code = VOXELPARTITION; else code = QUADTREE; } diff --git a/inst/include/lidR/GridPartition.h b/inst/include/lidR/GridPartition.h index 668446d6..c2e5f923 100644 --- a/inst/include/lidR/GridPartition.h +++ b/inst/include/lidR/GridPartition.h @@ -184,8 +184,15 @@ template void GridPartition::lookup(T& shape, std::vector& int colmax = std::ceil((xmax - this->xmin) / xres); int rowmin = std::floor((this->ymax - ymax) / yres); int rowmax = std::ceil((this->ymax - ymin) / yres); - int laymin = std::floor((zmin - this->zmin) / zres); - int laymax = std::ceil((zmax - this->zmin) / zres); + + int laymin = 0; + int laymax = nlayers; + if (zmin > this->zmin && zmax < this->zmax) + { + laymin = std::floor((zmin - this->zmin) / zres); + laymax = std::ceil((zmax - this->zmin) / zres); + } + int cell; res.clear(); diff --git a/inst/include/lidR/Octree.h b/inst/include/lidR/Octree.h index eb2ecbbb..2cac16ec 100644 --- a/inst/include/lidR/Octree.h +++ b/inst/include/lidR/Octree.h @@ -390,8 +390,14 @@ template Node::Ocnode* Octree::locate_region(T shape) double bbxmax = (shape.xmax - xmin)/(xmax-xmin); double bbymin = (shape.ymin - ymin)/(ymax-ymin); double bbymax = (shape.ymax - ymin)/(ymax-ymin); - double bbzmin = (shape.zmin - zmin)/(zmax-zmin); - double bbzmax = (shape.zmax - zmin)/(zmax-zmin); + + double bbzmin = 0; + double bbzmax = 1; + if (shape.zmin > this->zmin && shape.zmax < this->zmax) + { + bbzmin = (shape.zmin - zmin)/(zmax-zmin); + bbzmax = (shape.zmax - zmin)/(zmax-zmin); + } if (bbxmax < 0 || bbxmin > 1 || bbymax < 0 || bbymin > 1 || bbzmax < 0 || bbzmin > 1) return 0; diff --git a/tests/testthat/test-spatialindex.R b/tests/testthat/test-spatialindex.R index c39c106c..018994ad 100644 --- a/tests/testthat/test-spatialindex.R +++ b/tests/testthat/test-spatialindex.R @@ -162,7 +162,7 @@ las = LAS(D) test_that("circle lookup works for each spatial index", { - for(index in 0:4) + for(index in 0:3) { index(las) <- index @@ -183,7 +183,7 @@ test_that("circle lookup works for each spatial index", { test_that("oriented rectangle lookup works for each spatial index", { - for(index in 0:4) + for(index in 0:3) { index(las) <- index @@ -202,7 +202,7 @@ test_that("oriented rectangle lookup works for each spatial index", { test_that("knn in 2d works for each spatial indexes", { - for(index in 0:4) + for(index in 0:3) { index(las) <- index @@ -219,7 +219,7 @@ test_that("knn in 2d works for each spatial indexes", { test_that("knn in 3d works for each spatial indexes", { - for(index in 0:4) + for(index in 0:3) { index(las) <- index @@ -248,8 +248,8 @@ test_that("Spatial indexes work with more points (coverage)", { y = runif(1, 10, 900) z = runif(1, 1, 19) - u = vector("list", 4) - for(index in 0:4) + u = vector("list", 3) + for(index in 0:3) { index(las) <- index id = lidR:::C_knn3d_lookup(las, x,y,z,10) @@ -258,8 +258,8 @@ test_that("Spatial indexes work with more points (coverage)", { expect_true(all(sapply(u, identical, u[[1]]))) - u = vector("list", 4) - for(index in 0:4) + u = vector("list", 3) + for(index in 0:3) { index(las) <- index id = lidR:::C_circle_lookup(las, x, y, 5) @@ -278,8 +278,8 @@ test_that("Spatial indexes work with 0 point", { y = runif(1, 10, 900) z = runif(1, 1, 19) - u = vector("list", 4) - for(index in 0:4) + u = vector("list", 3) + for(index in 0:3) { index(las) <- index id = lidR:::C_knn3d_lookup(las, x, y, z,10) @@ -288,8 +288,8 @@ test_that("Spatial indexes work with 0 point", { expect_true(all(sapply(u, identical, integer(0)))) - u = vector("list", 4) - for(index in 0:4) + u = vector("list", 3) + for(index in 0:3) { index(las) <- index id = lidR:::C_circle_lookup(las, x, y, 5) @@ -308,8 +308,8 @@ test_that("Spatial indexes work with 1 point", { y = runif(1, 10, 11) z = runif(1, 1, 19) - u = vector("list", 4) - for(index in 0:4) + u = vector("list", 3) + for(index in 0:3) { index(las) <- index id = lidR:::C_knn3d_lookup(las, x, y, z, 10) @@ -318,8 +318,8 @@ test_that("Spatial indexes work with 1 point", { expect_true(all(sapply(u, identical, 1L))) - u = vector("list", 4) - for(index in 0:4) + u = vector("list", 3) + for(index in 0:3) { index(las) <- index id = lidR:::C_circle_lookup(las, x, y, 5) @@ -338,8 +338,8 @@ test_that("Spatial indexes work no point to find", { y = -100 z = 0 - u = vector("list", 4) - for(index in 0:4) + u = vector("list", 3) + for(index in 0:3) { index(las) <- index id = lidR:::C_circle_lookup(las, x, y, 5)