Attention

You are viewing an older version of the documentation. The latest version is 2.2.

Spatial Partitioning and Search Operations with Octrees

An octree is a tree data structure in which each internal node has exactly eight children. Octrees are most often used to partition a three-dimensional space by recursively subdividing it into eight octants. Octrees are the three-dimensional analog of quadtrees. The word is derived from oct (Greek root meaning “eight”) + tree. Octrees are often used in 3D graphic and 3D game engines.

In this tutorial we will learn how to use the octree for spatial partitioning and neighbor search within the point cloud data. We will explain how to perform “Neighbors within Radius Search”, “Approximate Nearest Neighbor (ANN) Search” and “K-Nearest Neighbors (KNN) Search”.

Note

This tutorial can be run both inside and outside a Docker* image. We assume that the pcl-oneapi-tutorial Deb package has been installed, and the user has copied the tutorial directory from /opt/intel/pcl/oneapi/tutorials/ to a user-writable directory.

  1. Prepare the environment:

    cd <path-to-oneapi-tutorials>/octree
    Copy to clipboard

  2. oneapi_octree_search.cpp should be in the directory with following content:

      1#include <iostream>
  2#include <fstream>
  3#include <numeric>
  4#include <pcl/oneapi/octree/octree.hpp>
  5#include <pcl/oneapi/containers/device_array.h>
  6#include <pcl/point_cloud.h>
  7
  8using namespace pcl::oneapi;
  9
 10float dist(Octree::PointType p, Octree::PointType q) {
 11    return std::sqrt((p.x-q.x)*(p.x-q.x) + (p.y-q.y)*(p.y-q.y) + (p.z-q.z)*(p.z-q.z));
 12}
 13
 14int main (int argc, char** argv)
 15{
 16    std::cout << "Running on device: " << dpct::get_default_queue().get_device().get_info<sycl::info::device::name>() << "\n";
 17
 18    std::size_t data_size = 871000;
 19    std::size_t query_size = 10000;
 20    float cube_size = 1024.f;
 21    float max_radius    = cube_size / 30.f;
 22    float shared_radius = cube_size / 30.f;
 23    const int max_answers = 5;
 24    const int k = 5;
 25    std::size_t i;
 26    std::vector<Octree::PointType> points;
 27    std::vector<Octree::PointType> queries;
 28    std::vector<float> radiuses;
 29    std::vector<int> indices;
 30
 31    //Generate point cloud data, queries, radiuses, indices
 32    srand (0);
 33    points.resize(data_size);
 34    for(i = 0; i < data_size; ++i)
 35    {
 36        points[i].x = ((float)rand())/(float)RAND_MAX * cube_size;
 37        points[i].y = ((float)rand())/(float)RAND_MAX * cube_size;
 38        points[i].z = ((float)rand())/(float)RAND_MAX * cube_size;
 39    }
 40
 41    queries.resize(query_size);
 42    radiuses.resize(query_size);
 43    for (i = 0; i < query_size; ++i)
 44    {
 45        queries[i].x = ((float)rand())/(float)RAND_MAX * cube_size;
 46        queries[i].y = ((float)rand())/(float)RAND_MAX * cube_size;
 47        queries[i].z = ((float)rand())/(float)RAND_MAX * cube_size;
 48        radiuses[i]  = ((float)rand())/(float)RAND_MAX * max_radius;
 49    };
 50
 51    indices.resize(query_size / 2);
 52    for(i = 0; i < query_size / 2; ++i)
 53    {
 54        indices[i] = i * 2;
 55    }
 56
 57    //Prepare oneAPI cloud
 58    pcl::oneapi::Octree::PointCloud cloud_device;
 59    cloud_device.upload(points);
 60
 61    //oneAPI build 
 62    pcl::oneapi::Octree octree_device;
 63    octree_device.setCloud(cloud_device);
 64    octree_device.build();
 65
 66    //Upload queries and radiuses
 67    pcl::oneapi::Octree::Queries queries_device;
 68    pcl::oneapi::Octree::Radiuses radiuses_device;
 69    queries_device.upload(queries);
 70    radiuses_device.upload(radiuses);
 71
 72    //Prepare output buffers on device
 73    pcl::oneapi::NeighborIndices result_device1(queries_device.size(), max_answers);
 74    pcl::oneapi::NeighborIndices result_device2(queries_device.size(), max_answers);
 75    pcl::oneapi::NeighborIndices result_device3(indices.size(), max_answers);
 76    pcl::oneapi::NeighborIndices result_device_ann(queries_device.size(), 1);
 77    pcl::oneapi::Octree::ResultSqrDists dists_device_ann;
 78    pcl::oneapi::NeighborIndices result_device_knn(queries_device.size(), k);
 79    pcl::oneapi::Octree::ResultSqrDists dists_device_knn;
 80
 81    //oneAPI octree radius search with shared radius
 82    octree_device.radiusSearch(queries_device, shared_radius, max_answers, result_device1);
 83
 84    //oneAPI octree radius search with individual radius
 85    octree_device.radiusSearch(queries_device, radiuses_device, max_answers, result_device2);
 86
 87    //oneAPI octree radius search with shared radius using indices to specify 
 88    //the queries.
 89    pcl::oneapi::Octree::Indices cloud_indices;
 90    cloud_indices.upload(indices);
 91    octree_device.radiusSearch(queries_device, cloud_indices, shared_radius, max_answers, result_device3);
 92
 93    //oneAPI octree ANN search
 94    //if neighbor points distances results are not required, can just call
 95    //octree_device.approxNearestSearch(queries_device, result_device_ann)
 96    octree_device.approxNearestSearch(queries_device, result_device_ann, dists_device_ann);
 97
 98    //oneAPI octree KNN search
 99    //if neighbor points distances results are not required, can just call
100    //octree_device.nearestKSearchBatch(queries_device, k, result_device_knn)
101    octree_device.nearestKSearchBatch(queries_device, k, result_device_knn, dists_device_knn);
102
103    //Download results
104    std::vector<int> sizes1;
105    std::vector<int> sizes2;
106    std::vector<int> sizes3;
107    result_device1.sizes.download(sizes1);
108    result_device2.sizes.download(sizes2);
109    result_device3.sizes.download(sizes3);
110
111    std::vector<int> downloaded_buffer1, downloaded_buffer2, downloaded_buffer3, results_batch;
112    result_device1.data.download(downloaded_buffer1);
113    result_device2.data.download(downloaded_buffer2);
114    result_device3.data.download(downloaded_buffer3);
115
116    int query_idx = 2;
117    std::cout << "Neighbors within shared radius search at ("
118              << queries[query_idx].x << " "
119              << queries[query_idx].y << " "
120              << queries[query_idx].z << ") with radius=" << shared_radius << std::endl;
121    for (i = 0; i < sizes1[query_idx]; ++i)
122    {
123        std::cout << "    "  << points[downloaded_buffer1[max_answers * query_idx + i]].x
124                  << " "     << points[downloaded_buffer1[max_answers * query_idx + i]].y
125                  << " "     << points[downloaded_buffer1[max_answers * query_idx + i]].z
126                  << " (distance: " << dist(points[downloaded_buffer1[max_answers * query_idx + i]], queries[query_idx])  << ")" << std::endl;
127    }
128
129    std::cout << "Neighbors within individual radius search at ("
130              << queries[query_idx].x << " "
131              << queries[query_idx].y << " "
132              << queries[query_idx].z << ") with radius=" << radiuses[query_idx] << std::endl;
133    for (i = 0; i < sizes2[query_idx]; ++i)
134    {
135        std::cout << "    "  << points[downloaded_buffer2[max_answers * query_idx + i]].x
136                  << " "     << points[downloaded_buffer2[max_answers * query_idx + i]].y
137                  << " "     << points[downloaded_buffer2[max_answers * query_idx + i]].z
138                  << " (distance: " << dist(points[downloaded_buffer2[max_answers * query_idx + i]], queries[query_idx])  << ")" << std::endl;
139    }
140
141    std::cout << "Neighbors within indices radius search at ("
142              << queries[query_idx].x << " "
143              << queries[query_idx].y << " "
144              << queries[query_idx].z << ") with radius=" << shared_radius << std::endl;
145    for (i = 0; i < sizes3[query_idx/2]; ++i)
146    {
147        std::cout << "    "  << points[downloaded_buffer3[max_answers * query_idx / 2 + i]].x
148                  << " "     << points[downloaded_buffer3[max_answers * query_idx / 2 + i]].y
149                  << " "     << points[downloaded_buffer3[max_answers * query_idx / 2 + i]].z
150                  << " (distance: " << dist(points[downloaded_buffer3[max_answers * query_idx / 2 + i]], queries[2])  << ")" << std::endl;
151    }
152
153    std::cout << "Approximate nearest neighbor at ("
154              << queries[query_idx].x << " "
155              << queries[query_idx].y << " "
156              << queries[query_idx].z << ")" << std::endl;
157    std::cout << "    "  << points[result_device_ann.data[query_idx]].x
158              << " "     << points[result_device_ann.data[query_idx]].y
159              << " "     << points[result_device_ann.data[query_idx]].z
160              << " (distance: " << std::sqrt(dists_device_ann[query_idx])  << ")" << std::endl;
161
162    std::cout << "K-nearest neighbors (k = " << k << ") at ("
163              << queries[query_idx].x << " "
164              << queries[query_idx].y << " "
165              << queries[query_idx].z << ")" << std::endl;
166    for (i = query_idx * k; i < (query_idx + 1) * k; ++i)
167    {
168        std::cout << "    "  << points[result_device_knn.data[i]].x
169                  << " "     << points[result_device_knn.data[i]].y
170                  << " "     << points[result_device_knn.data[i]].z
171                  << " (distance: " << std::sqrt(dists_device_knn[i])  << ")" << std::endl;
172    }
173}
    Copy to clipboard

  3. Source the Intel® oneAPI Base Toolkit environment:

    source /opt/intel/oneapi/setvars.sh
    Copy to clipboard

  4. (Optional) Setup proxy setting to download test data:

    export http_proxy="http://<http_proxy>:port"
export https_proxy="http://<https_proxy>:port"
    Copy to clipboard

  5. Build the code:

    mkdir build && cd build
cmake ../
make -j
    Copy to clipboard

  6. Run the binary:

    ./oneapi_octree_search
    Copy to clipboard

  7. Expected results example:

The search only finds the first five neighbors (as specified by max_answers), so a different radius finds different points.

Neighbors within shared radius search at (671.675 733.78 466.178) with radius=34.1333
  660.296 725.957 439.677 (distance: 29.8829)
  665.768 721.884 442.919 (distance: 26.7846)
  683.988 714.608 445.164 (distance: 30.9962)
  677.927 725.08 446.531 (distance: 22.3788)
  695.066 723.509 445.762 (distance: 32.7028)
Neighbors within individual radius search at (671.675 733.78 466.178) with radius=19.3623
  672.71 736.679 447.835 (distance: 18.6)
  664.46 731.504 452.074 (distance: 16.0048)
  671.238 725.881 461.408 (distance: 9.23819)
  667.707 718.527 466.622 (distance: 15.7669)
  654.552 733.636 467.795 (distance: 17.1993)
Neighbors within indices radius search at (671.675 733.78 466.178) with radius=34.1333
  660.296 725.957 439.677 (distance: 29.8829)
  665.768 721.884 442.919 (distance: 26.7846)
  683.988 714.608 445.164 (distance: 30.9962)
  677.927 725.08 446.531 (distance: 22.3788)
  695.066 723.509 445.762 (distance: 32.7028)
Copy to clipboard

Code Explanation

Generate point cloud data, queries, radiuses, indices with a random number.

    //Generate point cloud data, queries, radiuses, indices
    srand (0);
    points.resize(data_size);
    for(i = 0; i < data_size; ++i)
    {
        points[i].x = ((float)rand())/(float)RAND_MAX * cube_size;
        points[i].y = ((float)rand())/(float)RAND_MAX * cube_size;
        points[i].z = ((float)rand())/(float)RAND_MAX * cube_size;
    }

    queries.resize(query_size);
    radiuses.resize(query_size);
    for (i = 0; i < query_size; ++i)
    {
        queries[i].x = ((float)rand())/(float)RAND_MAX * cube_size;
        queries[i].y = ((float)rand())/(float)RAND_MAX * cube_size;
        queries[i].z = ((float)rand())/(float)RAND_MAX * cube_size;
        radiuses[i]  = ((float)rand())/(float)RAND_MAX * max_radius;
    };

    indices.resize(query_size / 2);
    for(i = 0; i < query_size / 2; ++i)
    {
        indices[i] = i * 2;
    }
Copy to clipboard

Create and build the oneAPI™ point cloud; then upload the queries and radiuses to a oneAPI™ device.

    //Prepare oneAPI cloud
    pcl::oneapi::Octree::PointCloud cloud_device;
    cloud_device.upload(points);

    //oneAPI build 
    pcl::oneapi::Octree octree_device;
    octree_device.setCloud(cloud_device);
    octree_device.build();

    //Upload queries and radiuses
    pcl::oneapi::Octree::Queries queries_device;
    pcl::oneapi::Octree::Radiuses radiuses_device;
    queries_device.upload(queries);
    radiuses_device.upload(radiuses);
Copy to clipboard

Create output buffers where we can download output from the oneAPI™ device.

    //Prepare output buffers on device
    pcl::oneapi::NeighborIndices result_device1(queries_device.size(), max_answers);
    pcl::oneapi::NeighborIndices result_device2(queries_device.size(), max_answers);
    pcl::oneapi::NeighborIndices result_device3(indices.size(), max_answers);
    pcl::oneapi::NeighborIndices result_device_ann(queries_device.size(), 1);
    pcl::oneapi::Octree::ResultSqrDists dists_device_ann;
    pcl::oneapi::NeighborIndices result_device_knn(queries_device.size(), k);
    pcl::oneapi::Octree::ResultSqrDists dists_device_knn;
Copy to clipboard

The fist radius search method is “search with shared radius”. In this search method, all queries use the same radius to find the neighbors.

    //oneAPI octree radius search with shared radius
    octree_device.radiusSearch(queries_device, shared_radius, max_answers, result_device1);

    //oneAPI octree radius search with individual radius
    octree_device.radiusSearch(queries_device, radiuses_device, max_answers, result_device2);
Copy to clipboard

The second radius search method is “search with individual radius”. In this search method, each query uses its own specific radius to find the neighbors.

    //oneAPI octree radius search with individual radius
    octree_device.radiusSearch(queries_device, radiuses_device, max_answers, result_device2);
Copy to clipboard

The third radius search method is “search with shared radius using indices”. In this search method, all queries use the same radius, and indices specify the queries.

    //oneAPI octree radius search with shared radius using indices to specify 
    //the queries.
    pcl::oneapi::Octree::Indices cloud_indices;
    cloud_indices.upload(indices);
    octree_device.radiusSearch(queries_device, cloud_indices, shared_radius, max_answers, result_device3);
Copy to clipboard

Perform ANN search.

    //oneAPI octree ANN search
    //if neighbor points distances results are not required, can just call
    //octree_device.approxNearestSearch(queries_device, result_device_ann)
    octree_device.approxNearestSearch(queries_device, result_device_ann, dists_device_ann);
Copy to clipboard

Perform KNN search.

    //oneAPI octree KNN search
    //if neighbor points distances results are not required, can just call
    //octree_device.nearestKSearchBatch(queries_device, k, result_device_knn)
    octree_device.nearestKSearchBatch(queries_device, k, result_device_knn, dists_device_knn);
Copy to clipboard

Download the search results from the oneAPI™ device. The size vector contains the size of found neighbors for each query. The downloaded_buffer vector contains the index of all found neighbors for each query.

    //Download results
    std::vector<int> sizes1;
    std::vector<int> sizes2;
    std::vector<int> sizes3;
    result_device1.sizes.download(sizes1);
    result_device2.sizes.download(sizes2);
    result_device3.sizes.download(sizes3);

    std::vector<int> downloaded_buffer1, downloaded_buffer2, downloaded_buffer3, results_batch;
    result_device1.data.download(downloaded_buffer1);
    result_device2.data.download(downloaded_buffer2);
    result_device3.data.download(downloaded_buffer3);
Copy to clipboard

Print the query, radius, and found neighbors to verify that the result is correct.

    int query_idx = 2;
    std::cout << "Neighbors within shared radius search at ("
              << queries[query_idx].x << " "
              << queries[query_idx].y << " "
              << queries[query_idx].z << ") with radius=" << shared_radius << std::endl;
    for (i = 0; i < sizes1[query_idx]; ++i)
    {
        std::cout << "    "  << points[downloaded_buffer1[max_answers * query_idx + i]].x
                  << " "     << points[downloaded_buffer1[max_answers * query_idx + i]].y
                  << " "     << points[downloaded_buffer1[max_answers * query_idx + i]].z
                  << " (distance: " << dist(points[downloaded_buffer1[max_answers * query_idx + i]], queries[query_idx])  << ")" << std::endl;
    }

    std::cout << "Neighbors within individual radius search at ("
              << queries[query_idx].x << " "
              << queries[query_idx].y << " "
              << queries[query_idx].z << ") with radius=" << radiuses[query_idx] << std::endl;
    for (i = 0; i < sizes2[query_idx]; ++i)
    {
        std::cout << "    "  << points[downloaded_buffer2[max_answers * query_idx + i]].x
                  << " "     << points[downloaded_buffer2[max_answers * query_idx + i]].y
                  << " "     << points[downloaded_buffer2[max_answers * query_idx + i]].z
                  << " (distance: " << dist(points[downloaded_buffer2[max_answers * query_idx + i]], queries[query_idx])  << ")" << std::endl;
    }

    std::cout << "Neighbors within indices radius search at ("
              << queries[query_idx].x << " "
              << queries[query_idx].y << " "
              << queries[query_idx].z << ") with radius=" << shared_radius << std::endl;
    for (i = 0; i < sizes3[query_idx/2]; ++i)
    {
        std::cout << "    "  << points[downloaded_buffer3[max_answers * query_idx / 2 + i]].x
                  << " "     << points[downloaded_buffer3[max_answers * query_idx / 2 + i]].y
                  << " "     << points[downloaded_buffer3[max_answers * query_idx / 2 + i]].z
                  << " (distance: " << dist(points[downloaded_buffer3[max_answers * query_idx / 2 + i]], queries[2])  << ")" << std::endl;
    }

    std::cout << "Approximate nearest neighbor at ("
              << queries[query_idx].x << " "
              << queries[query_idx].y << " "
              << queries[query_idx].z << ")" << std::endl;
    std::cout << "    "  << points[result_device_ann.data[query_idx]].x
              << " "     << points[result_device_ann.data[query_idx]].y
              << " "     << points[result_device_ann.data[query_idx]].z
              << " (distance: " << std::sqrt(dists_device_ann[query_idx])  << ")" << std::endl;

    std::cout << "K-nearest neighbors (k = " << k << ") at ("
              << queries[query_idx].x << " "
              << queries[query_idx].y << " "
              << queries[query_idx].z << ")" << std::endl;
    for (i = query_idx * k; i < (query_idx + 1) * k; ++i)
    {
        std::cout << "    "  << points[result_device_knn.data[i]].x
                  << " "     << points[result_device_knn.data[i]].y
                  << " "     << points[result_device_knn.data[i]].z
                  << " (distance: " << std::sqrt(dists_device_knn[i])  << ")" << std::endl;
    }
Copy to clipboard