How to Handle Efficiently the High Number of Points?

Point clouds are often huge and can contain a very high number of points data (several millions of points). Handling them efficiently in an application is not an easy task.

By default the point shape (RED::IPointShape) contains an array of unordered vertices. Searching specific points in such a data structure is particularly costly because we have to loop through the entire array. As we will have to do search operations (to find the nearest neighbours of a point for example), an efficient space-partitioning data structure must be used.

Using a k-d Tree

A k-d tree (k-dimensional) is a binary space partitioning tree. It is very useful for searches involving a multidimensional search key. It is perfect for our need of nearest-neighbours searching.

Note

We will not detail how we implemented the k-d tree here. The full implementation is available with the source code of the application.

The essence of a k-d tree is to be very efficient for queries (compared to an array) but it requires more time to construct (than an array).

Note

The time complexity of the search operation in a k-d tree is O(log n) at average whereas it is O(n) for a search in a basic array.

The k-d tree is a binary tree. Each node of the tree splits the space into two parts using a plane chosen along one of the k directions (3 in our case: x, y or z). The leaves of the tree contain the point cloud samples. All the other nodes contain the splitting informations.

An example of k-d tree of dimension 2

During the construction, the splitting stops when the node contains the desired number of points.

In order to be really effective during the search operation, the final tree needs to be well balanced, i.e. its branches need to be well-proportioned (all the branches have to contain approximately the same number of nodes from root to leaves).

Here is how we create our k-d tree in the application: First the vertex array is retrieved from the RED::IPointShape object:

// Get point positions array.
void* data;
RC_TEST( ipointshape->GetArray( data, RED::MCL_VERTEX, state ) );
float* vertex_in = (float*)data;

Then the structure is built from the point positions and balanced:

// Create the k-d tree:
// Split the space when the number of samples is reached in one node.
g_tree = new REDPointCloudKdTree( KD_SPLIT_ON_COUNT );

double bbmin[3], bbmax[3];
RED::Vector< REDPointCloudSample > samples;
samples.resize( count );

RFK::TutorialApplication::ShowProgressBar( true );

// Create point samples:
for( int i = 0; i < count; ++i )
{
    RFK::TutorialApplication::SetMessage( "Step 2: Creating the k-d tree.");
    RFK::TutorialApplication::UpdateProgressBar( (float)(i + 1) / (float)count );

    samples[i]._p[0] = vertex_in[ i * 3 + 0 ];
    samples[i]._p[1] = vertex_in[ i * 3 + 1 ];
    samples[i]._p[2] = vertex_in[ i * 3 + 2 ];

    samples[i]._index = i;
}

// Add samples to the tree and balance it:
g_tree->AddSamples( samples );
g_tree->Balance();

// Get the global bounding box to compute a shape offset.
// We will translate the points around world origin to minimize the float precision loss.
g_tree->GetBoundingBox( bbmin, bbmax );

g_mesh_offset[0] = ( bbmin[0] + bbmax[0] ) * 0.5;
g_mesh_offset[1] = ( bbmin[1] + bbmax[1] ) * 0.5;
g_mesh_offset[2] = ( bbmin[2] + bbmax[2] ) * 0.5;

The Find operation can then be very fast when we will query the n neighbours of a point:

// Get point.
p[0] = vertex_in[ i * 3 + 0 ];
p[1] = vertex_in[ i * 3 + 1 ];
p[2] = vertex_in[ i * 3 + 2 ];

// Retrieve the point nearest neighbours.
RC_TEST( g_tree->Find( neighbours, max_dist2, p, POINT_CLOUD_NEIGHBOURS_COUNT, workfind ) );

The Find function parameters are:

• The output list of neighbours

• The output maximum distance to the furthest neighbour

• The point to search around

• The number of neighbours we want

• A list of internal find information

Splitting the Heavy Point Shape

In addition to the enhancing of the processing phase with the k-d tree, we split the too heavy point shape into smaller ones. The objective here is to improve the rendering frame rate and allow the culling of the non-visible shapes when moving the camera. Moreover, the primitive index array can stay of type unsigned short integer (2 bytes) if the number of vertices per mesh does not exceed 65536.

// Compute the number of sub-meshes and the number of remaining points in the last shape.
// 'count' is the total number of points.
g_submesh_count = ( count / POINT_CLOUD_SUBMESH_SIZE ) + 1;
const int last_points_count = count % POINT_CLOUD_SUBMESH_SIZE;

// Define the arrays of new point shapes:
RED::Vector< RED::IPointShape* > inewpointshape;
inewpointshape.resize( g_submesh_count );
// A point data contains vertices, ...
RED::Vector< float* > vertex_out;
vertex_out.resize( g_submesh_count );
// ... color, ...
RED::Vector< unsigned char* > color_out;
color_out.resize( g_submesh_count );
// ... and custom data stored in texture coordinates.
RED::Vector< float* > texcoord_out;
texcoord_out.resize( g_submesh_count );

int vertcount;
for( int i = 0; i < g_submesh_count; ++i )
{
    // Create a new mesh shape.
    RED::Object* newshape = RED::Factory::CreateInstance( CID_REDPointShape );
    if( newshape == NULL )
        return RED_ALLOC_FAILURE;

    newshape->SetID( RED::String( "pointshape%1" ).Arg( i ).GetIDFromString() );
    RED::IShape* inewshape = newshape->As< RED::IShape >();
    inewpointshape[i] = newshape->As< RED::IPointShape >();

    // Set the material.
    RC_TEST( inewshape->SetMaterial( matr, state ) );

    // Attach the shape to the camera.
    RC_TEST( icamera->AddShape( newshape, state ) );

    // Define the number of points contained in the new shape.
    // The number of points in each shape is constant and the last shape stores the remaining points.
    vertcount = i < g_submesh_count - 1 ? POINT_CLOUD_SUBMESH_SIZE : last_points_count;

    // Create the vertex position array and get a pointer to it.
    RC_TEST( inewpointshape[i]->SetArray( RED::MCL_VERTEX, NULL, vertcount, 3, RED::MFT_FLOAT, state ) );
    RC_TEST( inewpointshape[i]->GetArray( data, RED::MCL_VERTEX, state ) );
    vertex_out[i] = (float*)data;

    // Create the color array and get a pointer to it.
    RC_TEST( inewpointshape[i]->SetArray( RED::MCL_COLOR, NULL, vertcount, 4, RED::MFT_UBYTE, state ) );
    RC_TEST( inewpointshape[i]->GetArray( data, RED::MCL_COLOR, state ) );
    color_out[i] = (unsigned char*)data;

    // Create the texture coordinates array and get a pointer to it.
    RC_TEST( inewpointshape[i]->SetArray( RED::MCL_TEX0, NULL, vertcount, 2, RED::MFT_FLOAT, state ) );
    RC_TEST( inewpointshape[i]->GetArray( data, RED::MCL_TEX0, state ) );
    texcoord_out[i] = (float*)data;

    // Create points in the new point shape.
    RC_TEST( inewpointshape[i]->AddPoints( NULL, vertcount, state ) );
}

At this point, the new point shapes are initialized but empty. They do not contain the correct points data yet. The description of this operation is done in the next page: Preparing the Point Shape Data.

The original point shape can then be removed from the scene graph (RED::ITransformShape::RemoveChild). Only the new ones will be seen as they were added in the camera scene graph (RED::IViewpoint::AddShape).

// Get root.
RED::Object *root;
RC_TEST( icamera->GetRootShape( root ) );
RED::IShape* iroot = root->As< RED::IShape >();
RED::ITransformShape* irootshape = root->As< RED::ITransformShape >();

// Get point shape and remove it from the scene graph.
RED::Object *shape;
RC_TEST( iroot->GetChild( shape, 0 ) );
RED::IShape* ishape = shape->As< RED::IShape >();
RC_TEST( irootshape->RemoveChild( shape, RED_SHP_DAG_NO_UPDATE, state ) );