So recently, I was asked to design a few applications which processed large data structures such as point clouds. The algorithms which I needed to implement were mostly smoothing algorithms where for each point or data in the cloud, you had to look at all the other data around it and perform some computations with that data. Obviously, if you are using a nearest neighbor approach, this is not that bad. It still suffers some drawbacks, such as when you are not sure of the area to process on or if you simply can't use a nearest neighbors approach to reduce the computational overhead. An obvious choice would be to create some type of partitioning structure which will allow each data point to only take into account those other points within its geometric range ( an individual cell ). So when dealing with arbitrary data points where you don't know how many cells you will need or even if the data would do well with a uniform grid partition (perhaps most points lie along the same X, Y planes and don't span the Z plane that much), then a different type of grid partitioning system would need to be used in these situations. The method that I chose to use and aim to discuss here is the geometricaly hashed grid.
In this formulation, each grid cell has a coordinate set of 3 integers. When looking for a cell to add the point to, the point's position is quickly and effiently converted to the grid cell's coordinate system and added to the bucket represented by that grid cell containing all the other poitns in the cell. The grid cell's integer coordinates are used as the key in a has table linking the cell's key to a bucket containing the points within that cell. The one caveat here is that cells are added IF AND ONLY IF there are points for this cell. This means that unneeded memory (and processing) overhead are reduced by making sure that empty cells are not stored in memory. It turns out that this structure is easy to code, efficient to process, and works well with testing algorithms which require some kind of partitioning structure to run in a reasonable amount of time, but don't need a complex data structure such as an octree quite yet.
So how did I want to interface this structure? It turns out that the amount of coding for this structure can be minimized using some of C++'s STL containers/algorithms, as well as templating it for easy expansion. When I show this code, please be aware that it is not a true hash algorithm in the sense that the coordinates used to index the table or not hashed to an index value for constant time access. Instead, the hash table remains sorted in memory, so whenever a grid's cell coordinate is used to index the hash table, there is a time complexity of Log N to perform the binary search. Some hash schemes could be extended to hash the 3 integer coordinates of the grid cells into a single integer for constant time access, but I'll leave that up to you :-D
Throughout this, Point<3, CT> == a 3-dimensional point where dimensions are type CT
First is the overall structure for the grid class:
template <typename CT>
class Grid {
public:
//
//Constructor/Desctructors
//
Grid ( const std::vector< Point< 3, CT > > &, CT );
~Grid();
//
//Gets a vector a points from each cell and those within (dist)cells away from it
//for example, if dist = 1, then all cells immediately touching the cell accessed
//would contribute their points. If dist = 2, then cells touching the immediate cells would contribute.
//etc. This is nessacary for operations needing points near the grid cells.
//
void GetNeighborCells( int, int, std::vector< Point< 3, CT > > & );
//
//Insert single point into grid cell.
//
void InsertIntoGrid( const Point< 3, CT > &);
//
//Retrieve a grid cell
//
std::vector< Point< 3, CT > > GetCellElements( int, int, int );
//
//Gets number of cells in grid
//
int CurrNumberOfCells() { return mNCells; }
//
//Gets number of elements in the grid
//
int CurrNumberOfElements() {
int retVal = 0;
for (std::map<Point<3, float>, Bucket >::const_iterator cit = mTable.begin();
cit != mTable.end();
++cit)
retVal += cit->second.data.size();
return retVal;
}
//
//Prints all of the index values for each grid cell (integers)
//
void PrintCellIndexValues(){
for (std::map<Point<3, float>, Bucket >::const_iterator cit = mTable.begin();
cit != mTable.end();
++cit)
}
//
//Individually access each cell with a single index value.
//
std::vector< Point< 3, CT> > operator[](int i){
assert ( (i+1) <= mTable.size() && (i >= 0) );
//need to get iterator since map[] does not do what you would expect
std::map< Point<3, float>, Bucket >::const_iterator cit = mTable.begin();
for (int j = 0; j < i; j++)
++cit; //increment to the proper position
std::vector < Point<3, CT> > retVal = cit->second.data;
return (retVal);
}
protected:
private:
//Methods
Grid ();
//Used to hold all the point elements for each cell.
typedef struct _Bucket{
std::vector< Point< 3, CT > > data;
} Bucket;
//Members
CT mEdgeLength; //Edge length of each cell
int mNCells; //number of cells total
std::map< Point< 3, float >, Bucket > mTable; //Actual Hash Table
};
While I'm not going to show each and ever method, I will show how the constructor sets up the initial grid, how the accessors work using the stl algorithms.
//
//Constructor
//
template <typename CT>
Grid<CT>::Grid( const std::vector< Point<3, CT> > & inPts, CT edgeLength){
//Set # cells to 0
mNCells = 0; //Initialize number of cells to 0
mEdgeLength = edgeLength; //Initialize edgeLength
Point< 3, float > tempPoint; //Used for intermediate data
Bucket tempBucket; //Used for intermediate data
int i, j, k; //Cell index values.
std::pair< std::map< Point<3, float>, Bucket>::iterator, bool > check; //Temporary individual elements of the hash table.
//insert the data into the cells
if (inPts.size() == 0){
std::cout << "Input points size == 0" << std::endl;
}
else{
//Insert to see if data is in the hash table
for (std::vector < Point< 3, CT> >::const_iterator cit = inPts.begin();
cit != inPts.end();
++cit){
//Take the 3-d coordinates of the point, and
//convert into cell coordinates. This is the
//simple way to go about it.
i = floor(cit->X() / mEdgeLength);
j = floor(cit->Y() / mEdgeLength);
k = floor(cit->Z() / mEdgeLength);
tempBucket.data.clear(); //clear temporary data from previous iteration.
tempBucket.data.push_back( (*cit) ); //Add the point to the bucket
tempPoint.X() = i;
tempPoint.Y() = j;
tempPoint.Z() = k;
//Since this is a map, if this insert operation fails due to
//The element already being in the hash table, then check.second will == false.
check = mTable.insert( std::pair< Point<3, float>, Bucket > (tempPoint, tempBucket ) );
//Use the return value to see if a new cell was added (previously empty)
//or if the insertion failed and the point just needs to be added to the current
//bucket
if (check.second == false){ //element already in hash table,
//So simply add point to current bucket for cell's grid coordinate.
(((check.first)->second).data).push_back(*cit); //Put point in bucket.
}
//else element successfully inserted.
else{
++mNCells; //A new cell was added.
}
}
}
}
And that's how the grid is set up. Basically, it converts each dimensional coordinate of each point into a 3-d grid integar coordinate. This coordinate is used to hash the grid cells. While iterating through the points, if the cells exist (the hash value already has an existing entry/bucket in the table) then the point is added to that entry's bucket. If not, the a new entry/buck are created and inserted into the std::map hash table.
To access each cell, I have created a method which will take the coordinate or index value and return a vector of points associated with that index value. Although this is simple and can be extended in many ways, it gets the job done.
//
// Get all poitns associated with hash index
//
template <typename CT>
std::vector< Point< 3, CT > > Grid<CT>::GetCellElements( int i, int j, int k ){
Point<3, float> tempPoint;
std::vector< Point<3, CT> > retVal;
std::map< Point<3, float>, Bucket >::const_iterator cit;
tempPoint.X() = i;
tempPoint.Y() = j;
tempPoint.Z() = k;
cit = this->mTable.find( tempPoint ); //looks for value
if ( cit != mTable.end() ){ //found the key
// assert (cit->second.data.size() != 0); //if 0, something went horribly wrong
retVal.insert( retVal.begin(), cit->second.data.begin(), cit->second.data.end() );
}
return retVal;
}
Also, the operator[] which instead takes an index value (great for iteration algorithms where the cell coordinates are not known and you simply want to iterate through cells [0] to [N] where N == number of cells in grid table) is inlined and shown up in the first code posted in this blog.
A proper destructor should be used to deallocate any memory, and if storage is a concern, this table's buckets could instead store references or pointers to an externally allocated array of points and dereference these when sorting/checking for grid coordinates. In this case, a custom sort algorithm must be provided to std::map for it to sort properly (it uses the default std::less algorithm which uses the templated classe's < operator to keep the map sorted, so make sure your point class has these operators overloaded!!!!! I hope this helps, and let me know if you have any questions or comments.
