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. This tutorial covers the process of performing “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-tutorials 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

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

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

  3. Source the Intel® oneAPI Base Toolkit environment:

    source /opt/intel/oneapi/setvars.sh

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

    export http_proxy="http://<http_proxy>:port"
export https_proxy="http://<https_proxy>:port"

  5. Build the code:

    mkdir build && cd build
cmake ../
make -j

  6. Run the binary:

    ./oneapi_octree_search

  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)

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;
}
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);

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;

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);

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);

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);

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);

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);

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);

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;
}