构建的总体函数
void buildIndexImpl()
{
// Create a permutable array of indices to the input vectors.
vind_.resize(size_); //size是点的个数参数 vind_为存放所有点索引的std::vector<int>变量 (kdtree_single_index.h 中)
for (size_t i = 0; i < size_; i++) {
vind_[i] = i; //vind_向量初始化
}
computeBoundingBox(root_bbox_); //根部包围盒root_bbox_ 是一个vector vector中 每个元素是struct struct中的两个参数为 T类型的 low 和high 两个变量
root_node_ = divideTree(0, size_, root_bbox_ ); // construct the tree
if (reorder_) {
data_ = flann::Matrix<ElementType>(new ElementType[size_*veclen_], size_, veclen_);
for (size_t i=0; i<size_; ++i) {
std::copy(points_[vind_[i]], points_[vind_[i]]+veclen_, data_[i]);
}
}
}
重点函数computeBoundingBox(root_bbox_);分析:
该函数是计算根部包围盒的
我的理解 veclen_应该是要构建成kd树的数据的维数 比如三维点云数据 该参数应该是3 下面以点云数据为例 分析代码
points_存储所有点云的数据 例如points_[0][0]表示点云数据中第一个点的x坐标数据,points_[0][1]表示点云数据中第一个点的y坐标数据,points_[0][2]表示点云数据中第一个点的z坐标数据
bbox[i].low和bbox[i].high 表示点云最外层包围盒的第i维的最小值和最大值
void computeBoundingBox(BoundingBox& bbox)
{
bbox.resize(veclen_);
//首先对vector<struct>bbox 进行初始化
for (size_t i=0; i<veclen_; ++i) {
bbox[i].low = (DistanceType)points_[0][i];
bbox[i].high = (DistanceType)points_[0][i];
}
//循环整个点云数据找到 能够完整包围整个点云的最小的外层包围长方体,每一个维度包含一个最大值一个最小值
for (size_t k=1; k<size_; ++k) {
for (size_t i=0; i<veclen_; ++i) {
if (points_[k][i]<bbox[i].low) bbox[i].low = (DistanceType)points_[k][i]; //k为点云中点的索引值 i为该点坐标维度索引
if (points_[k][i]>bbox[i].high) bbox[i].high = (DistanceType)points_[k][i];
}
}
}
root_bbox_参数分析:
kdtree_single_index.h 中
678行
/**
* Root bounding box根部包围盒
*/
BoundingBox root_bbox_;
334行
typedef std::vector<Interval> BoundingBox;
320行
struct Interval
{
DistanceType low, high;
private:
template <typename Archive>
void serialize(Archive& ar)
{
ar & low;
ar & high;
}
friend struct serialization::access;
};
flann.hpp
80行
typedef typename Distance::ResultType DistanceType;
dist.h
79行
typedef typename Accumulator<T>::Type ResultType;
50-63行
template<typename T>
struct Accumulator { typedef T Type; };
template<>
struct Accumulator<unsigned char> { typedef float Type; };
template<>
struct Accumulator<unsigned short> { typedef float Type; };
template<>
struct Accumulator<unsigned int> { typedef float Type; };
template<>
struct Accumulator<char> { typedef float Type; };
template<>
struct Accumulator<short> { typedef float Type; };
template<>
struct Accumulator<int> { typedef float Type; };
重点函数root_node_ = divideTree(0, size_, root_bbox_ ); // construct the tree分析:
kdtree_single_index.h 中
673行:
/**
* Array of k-d trees used to find neighbours.
*/
NodePtr root_node_;
317行:
typedef Node* NodePtr;
262行:
struct Node
{
/**
* Indices of points in leaf node 叶子节点中点的索引
*/
int left, right;
/**
* Dimension used for subdivision. 用于划分的尺寸
*/
int divfeat;
/**
* The values used for subdivision. 用于划分的值
*/
DistanceType divlow, divhigh;
/**
* The child nodes.
*/
Node* child1, * child2;
~Node()
{
if (child1) child1->~Node();
if (child2) child2->~Node();
}
private:
template<typename Archive>
void serialize(Archive& ar)
{
typedef KDTreeSingleIndex<Distance> Index;
Index* obj = static_cast<Index*>(ar.getObject());
ar & left;
ar & right;
ar & divfeat;
ar & divlow;
ar & divhigh;
bool leaf_node = false;
if (Archive::is_saving::value) {
leaf_node = ((child1==NULL) && (child2==NULL));
}
ar & leaf_node;
if (!leaf_node) {
if (Archive::is_loading::value) {
child1 = new(obj->pool_) Node();
child2 = new(obj->pool_) Node();
}
ar & *child1;
ar & *child2;
}
}
friend struct serialization::access;
};
/**创建一个树状节点,将vecs的列表从vind[first]细分到vind[last]。 该例程在每个子列表上都被递归调用。在pTree位置放置一个指向这个新树节点的指针。
* Create a tree node that subdivides the list of vecs from vind[first]
* to vind[last]. The routine is called recursively on each sublist.
* Place a pointer to this new tree node in the location pTree.
*
* Params: pTree = the new node to create
* first = index of the first vector
* last = index of the last vector
*/
NodePtr divideTree(int left, int right, BoundingBox& bbox)
{
NodePtr node = new (pool_) Node(); // allocate memory
/* If too few exemplars remain, then make this a leaf node. */
//当前待分配的区域的点数的大小小于叶节点可容纳的最大点数时,认为当前节点为叶节点,
//注意此处的left是和right 是索引vector(vind_)中点云索引的索引值(这个索引vector中放置的是原始点云中点云数据在点云中的索引值)points_为原始点云数据
// 则 points_[vind_[left]][i] 表示points_中第n个点的第i维的值 这个点是索引值n 这个索引值为vind_向量中第left+1个值的内容
if ( (right-left) <= leaf_max_size_) {
//因为是个叶节点 ,所以该叶节点的左孩子和有孩子的树 都是空的 并且这个叶节点的取值范围是left~right
//然后这个叶节点及所有的他的子节点 所在区域的包围盒 是两个for循环计算得到的。
node->child1 = node->child2 = NULL; /* Mark as leaf node. */
node->left = left;
node->right = right;
// compute bounding-box of leaf points
for (size_t i=0; i<veclen_; ++i) {
bbox[i].low = (DistanceType)points_[vind_[left]][i];
bbox[i].high = (DistanceType)points_[vind_[left]][i];
}
for (int k=left+1; k<right; ++k) {
for (size_t i=0; i<veclen_; ++i) {
if (bbox[i].low>points_[vind_[k]][i]) bbox[i].low=(DistanceType)points_[vind_[k]][i];
if (bbox[i].high<points_[vind_[k]][i]) bbox[i].high=(DistanceType)points_[vind_[k]][i];
}
}
}
else {
int idx;
int cutfeat;
DistanceType cutval;
middleSplit(&vind_[0]+left, right-left, idx, cutfeat, cutval, bbox);
node->divfeat = cutfeat;
BoundingBox left_bbox(bbox);
left_bbox[cutfeat].high = cutval;
node->child1 = divideTree(left, left+idx, left_bbox);
BoundingBox right_bbox(bbox);
right_bbox[cutfeat].low = cutval;
node->child2 = divideTree(left+idx, right, right_bbox);
node->divlow = left_bbox[cutfeat].high;
node->divhigh = right_bbox[cutfeat].low;
for (size_t i=0; i<veclen_; ++i) {
bbox[i].low = std::min(left_bbox[i].low, right_bbox[i].low);
bbox[i].high = std::max(left_bbox[i].high, right_bbox[i].high);
}
}
return node;
}
pool_
kdtree_single_index.h 中
690行:
/** 池化内存分配器。当有大量小的内存分配时,使用池式内存分配器比直接分配内存更有效率。
* Pooled memory allocator.
*
* Using a pooled memory allocator is more efficient
* than allocating memory directly when there is a large
* number small of memory allocations.
*/
PooledAllocator pool_;
PooledAllocator类的具体实现与解释在 allocator.h中
NodePtr node = new (pool_) Node(); // allocate memory
上一句是定位new重载的一个用法 (X)
不是定位new 是operator new 的一种带参重载 带的参数是PooledAllocator类对象而已
PooledAllocator类实现的功能有些类似定位new
其实new的全写是:new (expression-list)_optional Type (constructor arguments)_optional
pMemoryBlock = new (m_uiGrowNum, m_uiUnitSize) MemoryBlock(m_uiGrowNum, m_uiUnitSize);
就相当于调用 operator new(sizeof(tagMemoryBlock), m_uiGrowNum, m_uiUnitSize)
然后再调用constructor tagMemoryBlock(m_uiGrowNum, m_uiUnitSize)
具体功能及实现:
在allocator.h中的 197行
inline void* operator new (std::size_t size, flann::PooledAllocator& allocator)
{
return allocator.allocateMemory(size) ;
}
void middleSplit(int* ind, int count, int& index, int& cutfeat, DistanceType& cutval, const BoundingBox& bbox)
{
// find the largest span from the approximate bounding box
//1.从近似的边界盒中找出最大的跨度。
ElementType max_span = bbox[0].high-bbox[0].low;
cutfeat = 0;
cutval = (bbox[0].high+bbox[0].low)/2;
for (size_t i=1; i<veclen_; ++i) {
ElementType span = bbox[i].high-bbox[i].low;
if (span>max_span) {
max_span = span;
cutfeat = i;
cutval = (bbox[i].high+bbox[i].low)/2;
}
}
// compute exact span on the found dimension
//2.在找到的维度上计算精确的跨度
ElementType min_elem, max_elem;
computeMinMax(ind, count, cutfeat, min_elem, max_elem);
cutval = (min_elem+max_elem)/2;
max_span = max_elem - min_elem;
// check if a dimension of a largest span exists
//3.检查是否存在一个最大跨度的尺寸
size_t k = cutfeat;
for (size_t i=0; i<veclen_; ++i) {
if (i==k) continue;
ElementType span = bbox[i].high-bbox[i].low;
if (span>max_span) {
computeMinMax(ind, count, i, min_elem, max_elem);
span = max_elem - min_elem;
if (span>max_span) {
max_span = span;
cutfeat = i;
cutval = (min_elem+max_elem)/2;
}
}
}
//以上三个步骤不断检测是否找到的维度是否正确,以及精确的跨度。
//上述步骤解算出来的值为cutfeat为找到的哪个维度所包含的点在该维度上跨度最大。
//min_elem表示该维度上的最小值
//max_elem表示该维度上的最大值
//cutval表示该维度上的最大值和最小值的平均值
int lim1, lim2;
planeSplit(ind, count, cutfeat, cutval, lim1, lim2);
if (lim1>count/2) index = lim1;
else if (lim2<count/2) index = lim2;
else index = count/2;
assert(index > 0 && index < count);
}
//计算一堆点的数据,在dim这个维度上的最大值和最小值
//这一堆点数据 数量为count个,在原始点云中的索引存在了ind所指向空间向后数count个空间内(类型为int)
void computeMinMax(int* ind, int count, int dim, ElementType& min_elem, ElementType& max_elem)
{
min_elem = points_[ind[0]][dim];
max_elem = points_[ind[0]][dim];
for (int i=1; i<count; ++i) {
ElementType val = points_[ind[i]][dim];
if (val<min_elem) min_elem = val;
if (val>max_elem) max_elem = val;
}
}
/**在'cutval'位置用垂直于'cutfeat'尺寸的轴线对点列表进行细分。
//列表细分 调整索引位置的主要调整函数 其实就是找中位数的函数
重要 要仔细的读清楚
* Subdivide the list of points by a plane perpendicular on axe corresponding
* to the 'cutfeat' dimension at 'cutval' position.
*
* On return:
* dataset[ind[0..lim1-1]][cutfeat]<cutval
* dataset[ind[lim1..lim2-1]][cutfeat]==cutval
* dataset[ind[lim2..count]][cutfeat]>cutval
*/函数就是将索引列表中该包围盒中的索引调整至以cutval为分界线 左侧索引的点在该维的数据小于cutval 右侧索引的点在该维的数据大于cutval 该函数功能做了一个循环调整位置
void planeSplit(int* ind, int count, int cutfeat, DistanceType cutval, int& lim1, int& lim2)
{
int left = 0;
int right = count-1;
for (;; ) {
while (left<=right && points_[ind[left]][cutfeat]<cutval) ++left;
while (left<=right && points_[ind[right]][cutfeat]>=cutval) --right;
if (left>right) break;
std::swap(ind[left], ind[right]); ++left; --right;
}
lim1 = left;
right = count-1;
for (;; ) {
while (left<=right && points_[ind[left]][cutfeat]<=cutval) ++left;
while (left<=right && points_[ind[right]][cutfeat]>cutval) --right;
if (left>right) break;
std::swap(ind[left], ind[right]); ++left; --right;
}
lim2 = left;
}