Last Updated on 2022年1月4日
KD树与SKD树
首先给出两个搜到的有点内容的KD树文章,论述的比我说的更完整(更冗长),可以先看,也可以看完本文再看.
* https://www.zybuluo.com/l1ll5/note/967681
* http://www.whudj.cn/?p=920
主体思想
- KD树和SKD树都使用坐标轴对齐的最小包围盒来描述空间.
- 例如,平面内,一堆点的点集对应的空间可以用点
A
=(min(all_x),min(all_y))
,B
=(max(all_x),max(all_y))
对应的矩形空间来描述. 当点的个数变为1时,这个矩形空间也会自然地退化为一个点.
- 例如,平面内,一堆点的点集对应的空间可以用点
- 构建时的主要思想: 每个节点
node
都对应了一个空间box(node)
,node->plane
用一个超平面把这个空间box(node)
一分为二,记为sub_left
,sub_right
- 左子树
node->left
对应的空间box(node->left)
包含于sub_left
- 右子树
node->left
对应的空间box(node->right)
包含于sub_right
- 左子树
- 访问时的主要思想: 由于我们知道
box(node->left)
和box(node->right)
,那么我们就可以利用这些box
进行快捷的加速.- 例如, 如果我们需要找到距离
A
距离为5的点, 且我们已知A
到box(node->left)
的最近距离为10,那么,距离为5的点就一定不在node->left
中.
- 例如, 如果我们需要找到距离
注意: KD树和SKD的主体思想基本都是这样,但是具体的实现/应用中会有一些细节上的变化.
补充: 算法中使用的超平面都是垂直于坐标轴的超平面,很容易用数学语言描述,
x_i=K
就是穿过x_i=K
处且与x_i
轴垂直的超平面,例如三维空间中的z=8
KD树
直接上代码, It tells everything. 注意,这份代码主要关注的是可读性,而非性能.
#define DIM_NUM 3
struct Point {
double d[DIM_NUM];
double &operator[](int i) {
return d[i];
}
const double &operator[](int i) const {
return d[i];
}
};
struct Box {
Point lo;
Point hi;
Box() {
for (int i = 0; i < DIM_NUM; ++i) {
lo[i] = INFINITY;
hi[i] = -1.0 * INFINITY;
}
}
};
static std::function<bool(const Point &, const Point &)> GetCmpAtAxis(int axis_ind) {
return [axis_ind](const Point &l, const Point &r) -> bool { return l[axis_ind] < r[axis_ind]; };
}
struct Node {
Node *left;
Node *right;
Point p; // 既用于存储p,又用于作为分割平面
int axis_id;
Box node_box;
};
Box BuildNodeBoxFromPoints(const std::vector<Point> &points,int left,int right) {
Box ret;
for (int j = left; j <= right; ++j>) {
auto & p = points[j];
for (int i = 0; i < DIM_NUM; ++i) {
if (p[i] < ret.lo[i]) {
ret.lo[i] = p[i];
}
if (p[i] > ret.hi[i]) {
ret.hi[i] = p[i];
}
}
}
return ret;
}
int FindPlaneAxis(std::vector<Point> &points, int left, int right) {
// There are many way to determine the axis, here just shows a simple one.
static int i = -1;
++i;
return i % DIM_NUM;
}
// Note: points[right] is access-able
Node *BuildKDTree(std::vector<Point> &points, int left, int right) {
if (left > right) {
return nullptr;
}
int axis_id = FindPlaneAxis(points, left, right);
auto cmp = GetCmpAtAxis(axis_id);
int mid = left + (right - mid) / 2;
std::nth_element(points.begin() + left, points.begin() + mid, points.begin() + right + 1, cmp);
Node *node = new Node;
node->p = points[mid];
node->axis_id = axis_id;
node->node_box = BuildNodeBoxFromPoints(points,left,right);
node->left = BuildKDTree(points, left, mid - 1);
node->right = BuildKDTree(points, mid + 1, right);
return node;
}
- 补充:
- 显然,KD树是一个完全二叉树
- 因此,它是平衡树
- 必要时,可以使用数组来存储节点(第i个节点的子节点分别位于2i和2i+1处)
FindPlaneAxis
的较好实现是选择点集中标准差最大的方向.
- 寻找距离
target
最近点的一种方案:- 先设置
d=inf
作为初始距离. - 如果
PointToBoxDist(target,root->right) <= d
,那么最近点可能会在左子树中,否则一定不在左子树中 - 如果
PointToBoxDist(target,root->right) <= d
,那么最近点可能会在右子树中.否则一定不在右子树中 - 求解过程中不断更新
d
即可,很容易用递归实现.
- 先设置
SKD树: KD树的简单变种.
- 在KD树中,
Point
是树中存储的元素,而在SKD树中,变为了Box
, 在很多问题中, 最小的元素是存在体积的,例如三维重建/有限元仿真等使用的三维模型.- SKD树中,
Point
可以用退化的Box
表示. - 最小的体素是Box,划分空间时会用新的大Box来包络这些有体积的体素Box
- SKD树中,
- SKD建树过程中,一般用
box
的中心作为比较的基准.- 在确定了分割面所在的坐标
axis_point
后,一般认为,只要box的中心cent <= axis_point
,那么都划分到左子树,否则为右子树. - SKD树中,左子树和右子树的Box一般总会有重叠的部分,这和KD树是不同的.但这并不影响效率.
- SKD树中,非叶节点一般仅存储划分点的坐标,划分点并不是一个已经存在的体素,而是一个抽象的点. 也就是说,非叶节点并不会对应到某个体素.
- 在确定了分割面所在的坐标
附上: SKD 树的原始论文
Spatial KD tree A Data Structure for Geographic Database ( 访问密码:gZixPB )