• kd树
2022-03-25 15:49:14

## 1、三维空间点代码：

1、kdtree.h

/**

* @file   kdtree.h
* @brief Thisis a brief description.
* @author dongjian
* @date   2018:04:24
*  @note   mattersneeding attention
*/
#include <vector>
#include <queue>
#include <iostream>

namespace kdtree
{
class Point //空间中的三维点
{
public:
double x;
double y;
double z;
// 末尾const的作用是 不会改变本类对象中的实例成员变量
double Distance(const Point& p) const;
// 有参构造函数。
Point(double ox, double oy, double oz);
void printPoint(){
std::cout << "最近点：("<< this->x << "," << this->y << "," << this->z << ")"<<std::endl;

}
};

class BBox
{
public:
double xmin, ymin, zmin;
double xmax, ymax, zmax;

bool Contains(const Point& p) const;
bool Contains(const BBox& bbox) const;
// 相交的意思
bool Intersects(const BBox& bbox) const;

BBox(double x0, double x1, double y0, double y1, double z0, double z1);

// 静态函数，返回值是BBox
static BBox UniverseBBox();
double ShortestDistance(const Point& p) const;
};

enum Dimension
{
X = 0,
Y = 1,
Z = 2
};

class Node
{
public:
Node* left;			//left child
Node* right;		//right child
int begin;			//start index [close
int end;				//end index  (open
Dimension dimension;	//cut dimension
double pivot;		//cut value

Node(int b, int e, Dimension dim, double piv);
bool IsLeaf() const;
// 左子叶信封 什么鬼？
BBox LeftSubTreeEnvelope(const BBox& current_entext) const;
BBox RightSubTreeEnvelope(const BBox& current_entext) const;
};

class KDTree
{
public:
// 设置输入，构建树
void SetInput(const std::vector<Point>& inputs);

// 搜索
void RangeSearch(const BBox& bbox, std::vector<int>& indices) const ;

int NearestNeighbourSearch(const Point& searchpoint) const ;

void  NearestNeighbourSearch(const Point& searchpoint, int k, std::vector<int>& indices) const;

static int _leaf_max_size;
private:
const std::vector<Point>* pData; // 2维数组
std::vector<int> indices;
Node* pRoot;

private:
struct comparepair
{
bool operator()(const std::pair<int, double>&, const std::pair<int, double>&);
};

typedef std::priority_queue<std::pair<int, double>, std::vector<std::pair<int, double>>, comparepair> index_distance;

Node* DivideTree(const std::vector<Point>& inputs,int start,int end);
Dimension SelectDimension(const std::vector<Point>& inputs, int start, int end);
int Partition(int start, int end, Dimension dim);

double GetAt(int index, Dimension dim);

void _range_search(const BBox& search_bbox, Node* pNode, const BBox& current_extent,std::vector<int>& indices) const ;

int _nearest_neighbour(const Point& searchPoint, Node* pNode, const BBox& current_extent, double& current_shorest_dis) const ;

void _nearest_neighbour(const Point& searchPoint, int k, index_distance& ins, Node* pNode, const BBox& current_extent, double& current_shorest_dis) const ;
};

}


2、kdtree.cpp 含demo实现

#include "kdtree.h"
#include <algorithm>
#include <iostream>
#include <math.h>
using namespace std;
using namespace kdtree;

Point::Point(double ox, double oy, double oz) :
x(ox),
y(oy),
z(oz)
{}

double Point::Distance(const Point& p) const
{
double dx = p.x - x;
double  dy = p.y - y;
double dz = p.z - z;
return sqrt(dx*dx + dy*dy + dz*dz);
}

BBox Node::LeftSubTreeEnvelope(const BBox& current_entext) const
{
BBox leftRegion(current_entext);
switch (dimension)
{
case X:
leftRegion.xmax = pivot;
break;
case Y:
leftRegion.ymax = pivot;
break;
case Z:
leftRegion.zmax = pivot;
break;
}
return leftRegion;
}

BBox Node::RightSubTreeEnvelope(const BBox& current_entext) const
{
BBox rightRegion(current_entext);
switch (dimension)
{
case X:
rightRegion.xmin = pivot;
break;
case Y:
rightRegion.ymin = pivot;
break;
case Z:
rightRegion.zmin = pivot;
break;
}
return rightRegion;
}

bool BBox::Contains(const BBox& bbox) const
{
if (bbox.xmin < xmin)		return false;
if (bbox.xmax > xmax)		return false;
if (bbox.ymin < ymin)		return false;
if (bbox.ymax < ymax)		return false;
if (bbox.zmin < zmin)		return false;
if (bbox.zmax > zmax)		return false;
return true;
}

bool BBox::Contains(const Point& p) const
{
return p.x >= xmin && p.x <= xmax && p.y >= ymin && p.y <= ymax
&& p.z >= zmin && p.z <= zmax;
}

bool BBox::Intersects(const BBox& bbox) const
{
if (bbox.xmin > xmax || bbox.xmax < xmin) return false;
if (bbox.ymin > xmax || bbox.ymax < ymin) return false;
if (bbox.zmin > zmax || bbox.zmax < zmin) return false;
return true;
}

BBox BBox::UniverseBBox()
{
double DOUBLE_MAX = std::numeric_limits<double>::max();
return BBox(-DOUBLE_MAX, DOUBLE_MAX, -DOUBLE_MAX,
DOUBLE_MAX, -DOUBLE_MAX, DOUBLE_MAX);
}

BBox::BBox(double x0, double x1, double y0, double y1, double z0, double z1) :
xmin(x0),
xmax(x1),
ymin(y0),
ymax(y1),
zmin(z0),
zmax(z1)
{}

double BBox::ShortestDistance(const Point& p) const
{
double dx = xmin - p.x > 0 ? xmin - p.x : (p.x - xmax > 0 ? p.x - xmax : 0);
double dy = ymin - p.y > 0 ? ymin - p.y : (p.y - ymax > 0 ? p.y - ymax : 0);
double dz = zmin - p.z > 0 ? zmin - p.z : (p.z - zmax > 0 ? p.z - zmax : 0);
return sqrt(dx*dx + dy*dy + dz*dz);
}

Node::Node(int b, int e, Dimension dim, double piv) :
begin(b),
end(e),
dimension(dim),
pivot(piv),
left(nullptr),
right(nullptr)
{

}

bool Node::IsLeaf() const
{
return left == nullptr && right == nullptr;
}
// 选择维度？
/**
* @brief
*
* @param inputs 输入原始数据
* @param start 起始索引
* @param end
* @return Dimension  输入的数据中，在x y z轴上哪个跨度大，就输出哪个枚举类型
*/
Dimension KDTree::SelectDimension(const std::vector<Point>& inputs, int start, int end)
{
struct cmpobj
{
Dimension dim;
cmpobj(Dimension d) :
dim(d){}

bool operator()(const Point& p0, const Point p1) const
{
if (dim == X)
return p0.x < p1.x;
if (dim == Y)
return p0.y < p1.y;
if (dim == Z)
return p0.z < p1.z;
return false;
}
};
// span存的是最大－最小值
double span[3];
// 该pair 的first指向[first,last)范围内最小值元素，second指向最大值元素。
auto pair =  std::minmax_element(inputs.begin() + start, inputs.begin() + end,cmpobj(X));
span[0] = pair.second->x - pair.first->x;

pair = std::minmax_element(inputs.begin() + start, inputs.begin() + end, cmpobj(Y));
span[1] = pair.second->y - pair.first->y;

pair = std::minmax_element(inputs.begin() + start, inputs.begin() + end, cmpobj(Z));
span[2] = pair.second->z - pair.first->z;

// 获得的索引是 span里面最大元素与第一个元素之间的距离，也就是获得最大元素的索引。
auto index = std::distance(span, std::max_element(span, span + 3));

return static_cast<Dimension>(index);
}

double KDTree::GetAt(int index, Dimension dim)
{
auto p = (*pData)[indices[index]];
return  dim == X ? p.x : (dim == Y ? p.y : p.z);
}

int KDTree::Partition(int start, int end, Dimension dimension)
{
int size = end - start;
if (size <= 0)
{
cout << "a serious error occurs " << start << "\t" << endl;
return -1;
}

struct cmpobj
{
Dimension dim;
const std::vector<Point>* pData;
bool operator()(int i, int j) const
{
if (X == dim)
return (*pData)[i].x < (*pData)[j].x;
if (Y == dim)
return (*pData)[i].y < (*pData)[j].y;
if (Z == dim)
return (*pData)[i].z < (*pData)[j].z;

return true;
}
cmpobj(Dimension dimension, const std::vector<Point>* pInputData):
dim(dimension),
pData(pInputData)
{}
};

int median = start +size / 2;
std::nth_element(indices.begin() + start, indices.begin()+median , indices.begin() + end, cmpobj(dimension, pData));
return median;
}

int KDTree::NearestNeighbourSearch(const Point& searchpoint) const
{
double shortestDistance = std::numeric_limits<double>::max();
return _nearest_neighbour(searchpoint, pRoot, BBox::UniverseBBox(), shortestDistance);
}

void KDTree::NearestNeighbourSearch(const Point& searchpoint, int k, std::vector<int>& indices) const
{
double shortestDistance = std::numeric_limits<double>::max();
indices.clear();

index_distance neighbours;
_nearest_neighbour(searchpoint, k, neighbours, pRoot, BBox::UniverseBBox(), shortestDistance);

while (!neighbours.empty())
{
indices.push_back(neighbours.top().first);
neighbours.pop();
}
}

bool KDTree::comparepair::operator()(const std::pair<int, double>& p0, const std::pair<int, double>& p1)
{
return p0.second < p1.second;
}

void KDTree::_nearest_neighbour(const Point& searchPoint, int k, index_distance& ins,
Node* pNode, const BBox& current_extent, double& current_shorest_dis) const
{
double min_shortest_distance = current_extent.ShortestDistance(searchPoint);
if (min_shortest_distance >= current_shorest_dis)
return;

if (pNode->IsLeaf())
{
for (auto i = pNode->begin; i < pNode->end; ++i)
{
double distance = (*pData)[indices[i]].Distance(searchPoint);
if (ins.size() < k)
{
if (k == ins.size())
current_shorest_dis = ins.top().second;
}
else
{
if (distance < current_shorest_dis)
{
ins.pop();
current_shorest_dis = ins.top().second;
}
}
}
}
else
{
BBox leftRegion = pNode->LeftSubTreeEnvelope(current_extent);
BBox rightRegion = pNode->RightSubTreeEnvelope(current_extent);

double dis_to_left = leftRegion.ShortestDistance(searchPoint);
double dis_to_right = rightRegion.ShortestDistance(searchPoint);

if (dis_to_left < dis_to_right)
{
_nearest_neighbour(searchPoint,k,ins,pNode->left,leftRegion,current_shorest_dis);
_nearest_neighbour(searchPoint, k,ins,pNode->right, rightRegion, current_shorest_dis);
}
else
{
_nearest_neighbour(searchPoint,k,ins, pNode->right, rightRegion, current_shorest_dis);
_nearest_neighbour(searchPoint,k,ins, pNode->left, leftRegion, current_shorest_dis);
}
}
}

int KDTree::_nearest_neighbour(const Point& searchPoint, Node* pNode, const BBox& current_extent, double& current_shorest_dis) const
{
double min_shortest_distance = current_extent.ShortestDistance(searchPoint);
if (min_shortest_distance >= current_shorest_dis)
return -1;

if (pNode->IsLeaf())
{
int shortestindex =-1;
for (auto i = pNode->begin; i < pNode->end;++i)
{
double distance = (*pData)[indices[i]].Distance(searchPoint);
if (distance < current_shorest_dis)
{
shortestindex = indices[i];
current_shorest_dis = distance;
}
}
return shortestindex;
}
else
{
BBox leftRegion = pNode->LeftSubTreeEnvelope(current_extent);
BBox rightRegion = pNode->RightSubTreeEnvelope(current_extent);

double dis_to_left = leftRegion.ShortestDistance(searchPoint);
double dis_to_right = rightRegion.ShortestDistance(searchPoint);

if (dis_to_left < dis_to_right)
{
int left = _nearest_neighbour(searchPoint, pNode->left, leftRegion, current_shorest_dis);
int right = _nearest_neighbour(searchPoint, pNode->right, rightRegion, current_shorest_dis);

return right == -1 ? left : right;
}
else
{
int right = _nearest_neighbour(searchPoint, pNode->right, rightRegion, current_shorest_dis);
int left = _nearest_neighbour(searchPoint, pNode->left, leftRegion, current_shorest_dis);
return left == -1 ? right : left;
}

return -1;
}
}

int KDTree::_leaf_max_size = 2; //原来是15

void KDTree::SetInput(const std::vector<Point>& inputs)
{
pData = &inputs;// 指针的首地址 指向了inputs的地址
indices.resize(inputs.size());
for (int i = 0, n = inputs.size(); i < n; ++i)
indices[i] = i;// 1 2 3 4 5 6 ....
// Node类型的指针
pRoot = DivideTree(inputs, 0, inputs.size()); //划分树
}

Node* KDTree::DivideTree(const std::vector<Point>& inputs, int start, int end)
{
//cout << "build " << start << "\t" << end << endl;

int size = end - start;
if (size <= 0)
return nullptr;

Dimension dim = SelectDimension(inputs, start, end); // 选择划分维度，是按X轴划分还是y轴还是z轴
// 中位数
int median = Partition(start, end, dim); //获取中位数

Node* pNode = new Node(start, end, dim, GetAt(median, dim));
//递归终止条件是输入的点的个数小于叶子节点最大点个数
if (size > _leaf_max_size)
{
pNode->left = DivideTree(inputs, start, median);
pNode->right = DivideTree(inputs, median, end);
}

return pNode;
}

void KDTree::RangeSearch(const BBox& bbox, std::vector<int>& indices) const
{
BBox universe_bbox = BBox::UniverseBBox();
_range_search(bbox, pRoot, universe_bbox, indices);
}

void KDTree::_range_search(const BBox& search_bbox, Node* pNode, const BBox& current_extent, std::vector<int>& ins) const
{
if (nullptr == pNode)
return;

if (pNode->IsLeaf())
{
for (int i = pNode->begin; i < pNode->end; ++i)
{
const Point& p = (*pData)[indices[i]];
if (search_bbox.Contains(p))
ins.push_back(indices[i]);
}
}
else
{
//trim bounding box
BBox leftRegion = pNode->LeftSubTreeEnvelope(current_extent);

if (search_bbox.Contains(leftRegion))
{
for (int i = pNode->left->begin; i < pNode->left->end; ++i)
ins.push_back(indices[i]);
}
else if (search_bbox.Intersects(leftRegion))
{
_range_search(search_bbox, pNode->left, leftRegion, ins);
}

BBox rightRegion = pNode->RightSubTreeEnvelope(current_extent);

if (search_bbox.Contains(rightRegion))
{
for (int i = pNode->right->begin; i < pNode->right->end; ++i)
ins.push_back(indices[i]);
}
else if (search_bbox.Intersects(rightRegion))
{
_range_search(search_bbox, pNode->right, rightRegion, ins);
}
}
}
int main(int argc, char *argv[]){

std::vector<kdtree::Point> input;

// 1、输入数据
for (int  i = 0; i < 10; i++)
{
double x,y,z;
std::cin>>x>>y>>z;
std::cout<<x<<y<<z<<std::endl;
kdtree::Point p(x, y, z);

input.push_back(p);
}
//2、构建树
KDTree kdtree;
kdtree.SetInput(input);

//3、最近邻搜索
//创建寻找的点
int k =3;
Point p1(1.2 , 2.1, 3.2);
std::vector<int> indices;//索引
kdtree.NearestNeighbourSearch(p1, k, indices);

//4、打印K个最近点
for(int i = 0; i<k; i++ ){

input[indices[i]].printPoint();
std::cout << "距离是：" << p1.Distance(input[indices[i]]) << std::endl;
}

return 0;
}


结果：

## 2、二维空间代码

1、kdtree_new.h

//
//  KDTree.h
//  Test
//
//  Created by xiuzhu on 2021/7/13.
//

#ifndef kdtree_new_h
#define kdtree_new_h

#include <vector>

using namespace std;
//template <class T>
class KDTree {
private:
int key;
vector<double> root;//树及其子树的根节点
KDTree *parent;
KDTree *left_child;
KDTree *right_child;

public:
KDTree();
~KDTree();
bool is_empty();//判断KD树是否为空
bool is_leaf();//判断树是否只有一个叶子节点
bool is_root();//判断是否是树的根节点
bool is_left();//判断该子kd树的根结点是否是其父kd树的左结点
bool is_right();//判断该子kd树的根结点是否是其父kd树的右结点
vector<vector<double> > Transpose(const vector<vector<double> > &Matrix);//坐标转换
double findMiddleValue(vector<double> vec);//查找中指
void buildKdTree(KDTree* tree, vector<vector<double> > data, unsigned depth);//构建kd树
void printKdTree(KDTree* tree, unsigned int depth);//逐层打印KD树
double measureDistance(vector<double> point1, vector<double> point2, unsigned method);//计算空间中两个点的距离
vector<double> searchNearestNeighbor(vector<double> goal, KDTree *tree);在kd树tree中搜索目标点goal的最近邻

};

#endif /* KDTree_h */



2、kdtree_new.cpp

//
//  KDTree.cpp
//  Test
//
//  Created by xiuzhu on 2021/7/13.
//
#include "kdtree_new.h"
#include <stdio.h>
#include <iostream>
#include <vector>
#include <algorithm>
#include <string>
#include <cmath>

using namespace std;

KDTree::KDTree():parent(nullptr),left_child(nullptr),right_child(nullptr){}
KDTree::~KDTree(){}
bool KDTree::is_empty()
{
return root.empty();
}
bool KDTree::is_leaf()
{
return (!root.empty()) && right_child == nullptr && left_child == nullptr;
}
bool KDTree::is_root()
{
return (!is_empty()) && parent == nullptr;
}
bool KDTree::is_left()
{
return parent->left_child->root == root;
}
bool KDTree::is_right()
{
return parent->right_child->root == root;
}
//用于转换坐标
vector<vector<double> > KDTree::Transpose(const vector<vector<double>> &Matrix)
{
unsigned row = Matrix.size();
unsigned col = Matrix[0].size();
vector<vector<double> > Trans(col,vector<double>(row,0));
for (unsigned i = 0; i < col; ++i)
{
for (unsigned j = 0; j < row; ++j)
{
Trans[i][j] = Matrix[j][i];
}
}
return Trans;
}
//在不同的坐标轴上寻找中值

double KDTree::findMiddleValue(vector<double> vec)
{
sort(vec.begin(),vec.end());
auto pos = vec.size() / 2;
return vec[pos];
}
void KDTree::buildKdTree(KDTree *tree, vector<vector<double>> data, unsigned int depth)
{
//样本的数量
unsigned long samplesNum = data.size();
//终止条件
if (samplesNum == 0)
{
return;
}
if (samplesNum == 1)
{
tree->root = data[0];
return;
}
//样本的维度
unsigned long k = data[0].size();//坐标轴个数
vector<vector<double> > transData = Transpose(data);
//选择切分属性
unsigned splitAttribute = depth % k;
vector<double> splitAttributeValues = transData[splitAttribute];
//选择切分值
double splitValue = findMiddleValue(splitAttributeValues);
//cout << "splitValue" << splitValue  << endl;

// 根据选定的切分属性和切分值，将数据集分为两个子集
vector<vector<double> > subset1;
vector<vector<double> > subset2;
for (unsigned i = 0; i < samplesNum; ++i)
{
if (splitAttributeValues[i] == splitValue && tree->root.empty())
tree->root = data[i];
else
{
if (splitAttributeValues[i] < splitValue)
subset1.push_back(data[i]);
else
subset2.push_back(data[i]);
}
}

//子集递归调用buildKdTree函数
tree->left_child = new KDTree;
tree->left_child->parent = tree;
tree->right_child = new KDTree;
tree->right_child->parent = tree;
buildKdTree(tree->left_child, subset1, depth + 1);
buildKdTree(tree->right_child, subset2, depth + 1);
}
void KDTree::printKdTree(KDTree *tree, unsigned int depth)//打印
{
for (unsigned i = 0; i < depth; ++i)
cout << "\t";

for (vector<double>::size_type j = 0; j < tree->root.size(); ++j)
cout << tree->root[j] << ",";
cout << endl;
if (tree->left_child == nullptr && tree->right_child == nullptr )//叶子节点
return;
else //非叶子节点
{
if (tree->left_child != nullptr)
{
for (unsigned i = 0; i < depth + 1; ++i)
cout << "\t";
cout << " left:";
printKdTree(tree->left_child, depth + 1);
}

cout << endl;
if (tree->right_child != nullptr)
{
for (unsigned i = 0; i < depth + 1; ++i)
cout << "\t";
cout << "right:";
printKdTree(tree->right_child, depth + 1);
}
cout << endl;
}
}
//计算空间中两个点的距离
double KDTree::measureDistance(vector<double> point1, vector<double> point2, unsigned int method){
if (point1.size() != point2.size())
{
cerr << "Dimensions don't match！！" ;
exit(1);
}
switch (method)
{
case 0://欧氏距离
{
double res = 0;
for (vector<double>::size_type i = 0; i < point1.size(); ++i)
{
res += pow((point1[i] - point2[i]), 2);
}
return sqrt(res);
}
case 1://曼哈顿距离
{
double res = 0;
for (vector<double>::size_type i = 0; i < point1.size(); ++i)
{
res += abs(point1[i] - point2[i]);
}
return res;
}
default:
{
cerr << "Invalid method!!" << endl;
return -1;
}
}

}
//在kd树tree中搜索目标点goal的最近邻
//输入：目标点；已构造的kd树
//输出：目标点的最近邻
vector<double> KDTree::searchNearestNeighbor(vector<double> goal, KDTree *tree)
{
/*第一步：在kd树中找出包含目标点的叶子结点：从根结点出发，
递归的向下访问kd树，若目标点的当前维的坐标小于切分点的
坐标，则移动到左子结点，否则移动到右子结点，直到子结点为
叶结点为止,以此叶子结点为“当前最近点”
*/
unsigned long k = tree->root.size();//计算出数据的维数
unsigned d = 0;//维度初始化为0，即从第1维开始
KDTree* currentTree = tree;
vector<double> currentNearest = currentTree->root;
while(!currentTree->is_leaf())
{
unsigned index = d % k;//计算当前维
if (currentTree->right_child->is_empty() || goal[index] < currentNearest[index])
{
currentTree = currentTree->left_child;
}
else
{
currentTree = currentTree->right_child;
}
++d;
}
currentNearest = currentTree->root;

/*第二步：递归地向上回退， 在每个结点进行如下操作：
(a)如果该结点保存的实例比当前最近点距离目标点更近，则以该例点为“当前最近点”
(b)当前最近点一定存在于某结点一个子结点对应的区域，检查该子结点的父结点的另
一子结点对应区域是否有更近的点（即检查另一子结点对应的区域是否与以目标点为球
心、以目标点与“当前最近点”间的距离为半径的球体相交）；如果相交，可能在另一
个子结点对应的区域内存在距目标点更近的点，移动到另一个子结点，接着递归进行最
近邻搜索；如果不相交，向上回退*/

//当前最近邻与目标点的距离
double currentDistance = measureDistance(goal, currentNearest, 0);

//如果当前子kd树的根结点是其父结点的左孩子，则搜索其父结点的右孩子结点所代表
//的区域，反之亦反
KDTree* searchDistrict;
if (currentTree->is_left())
{
if (currentTree->parent->right_child == nullptr)
searchDistrict = currentTree;
else
searchDistrict = currentTree->parent->right_child;
}
else
{
searchDistrict = currentTree->parent->left_child;
}

//如果搜索区域对应的子kd树的根结点不是整个kd树的根结点，继续回退搜索
while (searchDistrict->parent != nullptr)
{
//搜索区域与目标点的最近距离
double districtDistance = abs(goal[(d+1)%k] - searchDistrict->parent->root[(d+1)%k]);

//如果“搜索区域与目标点的最近距离”比“当前最近邻与目标点的距离”短，表明搜索
//区域内可能存在距离目标点更近的点
if (districtDistance < currentDistance )//&& !searchDistrict->isEmpty()
{

double parentDistance = measureDistance(goal, searchDistrict->parent->root, 0);

if (parentDistance < currentDistance)
{
currentDistance = parentDistance;
currentTree = searchDistrict->parent;
currentNearest = currentTree->root;
}
if (!searchDistrict->is_empty())
{
double rootDistance = measureDistance(goal, searchDistrict->root, 0);
if (rootDistance < currentDistance)
{
currentDistance = rootDistance;
currentTree = searchDistrict;
currentNearest = currentTree->root;
}
}
if (searchDistrict->left_child != nullptr)
{
double leftDistance = measureDistance(goal, searchDistrict->left_child->root, 0);
if (leftDistance < currentDistance)
{
currentDistance = leftDistance;
currentTree = searchDistrict;
currentNearest = currentTree->root;
}
}
if (searchDistrict->right_child != nullptr)
{
double rightDistance = measureDistance(goal, searchDistrict->right_child->root, 0);
if (rightDistance < currentDistance)
{
currentDistance = rightDistance;
currentTree = searchDistrict;
currentNearest = currentTree->root;
}
}
}//end if

if (searchDistrict->parent->parent != nullptr)
{
searchDistrict = searchDistrict->parent->is_left()?
searchDistrict->parent->parent->right_child:
searchDistrict->parent->parent->left_child;
}
else
{
searchDistrict = searchDistrict->parent;
}
++d;
}//end while
return currentNearest;
}

const int data[6][2]={{2,3},{5,4},{9,6},{4,7},{8,1},{7,2}};
int main()
{
vector<vector<double> > train(6, vector<double>(2, 0)); // 初始化为6个 值为（2,0）的坐标
for (unsigned i = 0; i < 6; ++i)
for (unsigned j = 0; j < 2; ++j)
train[i][j] = data[i][j];
KDTree* kdTree = new KDTree();
kdTree->buildKdTree(kdTree, train, 0);
kdTree->printKdTree(kdTree, 0);

vector<double> goal;
goal.push_back(4);
goal.push_back(6);
vector<double> nearestNeighbor = kdTree->searchNearestNeighbor(goal, kdTree);
vector<double>::iterator beg = nearestNeighbor.begin();
cout << "The nearest neighbor is: ";
while(beg != nearestNeighbor.end()) cout << *beg++ << ",";
cout << endl;
return 0;
}


更多相关内容
• 通过在REFERENCE中从数据点构建kD树，并在MODEL中为每个数据点查询树，以有效的方式执行搜索。 PTS = KDRANGEQUERY( ROOT, QUERYPT, DISTLIM ) 查找存储在 kD 树 ROOT 中且在 QUERYPT 的 DISTLIM 单元内的所有点。...
• python K近近邻邻算算法法的的kd树树实实现现 这篇文章主要介绍了python K近邻算法的kd树实现小编觉得挺不错的现在分享给大家也给大家做个参考 一起跟随小编过 看看吧 k近近邻算算法法的的介介绍 k近邻算法是一种...
• 主要的方法就是，将数据存储在kd树这种空间数据结构中，树的思想类似于二叉搜索树、红黑树等，真的很强大（相信童靴们明白的）。改进的k-means方法主要考虑了合理的初始候选质心的方式，并且利用kd树的特性，减少了...
• 主要介绍了python K近邻算法的kd树实现，小编觉得挺不错的，现在分享给大家，也给大家做个参考。一起跟随小编过来看看吧
• 一种用于求解振动问题非线性微分方程的计算方法。一般来说，求解振动问题微分方程的计算方法主要有：时域法、频域法和时频变换方法
• 利用Matlab对三维点云建立KD树，搜索一点或多点的柱状邻域、球状邻域与KNN点。 其中，柱状邻域、球状邻域搜索半径为r内的三维点； KNN搜索最邻近k个三维点。
• 用C++实现KD树的构建以及递归输出。。。
• ikd-Tree是为机器人应用程序设计的增量kd树。 ikd-Tree仅用新的到来点来增量更新kd树，因此与现有的静态kd树相比，计算时间要短得多。 除了逐点操作外，ikd-Tree还支持盒式操作和下采样等多项功能，这些功能在机器人...
• matlab 的大部分 kdtree 代码都是通过 mex 文件实现的。 我决定提出一个纯粹基于 matlab 的实现，所以这里是......代码显然预计会比那里的一些 c/c++ 实现慢，但它在 matlab 中实现的事实可能会使它在某些情况下很...
• 基于KD树实现的python算法，用于KD树的生成及K近邻算法
• matlab实现kd树的创建和配套的搜索程序，注释详细，还附有算法思路讲解。欢迎下载
• 经典的kd树数据结构C#代码实现，代码通俗易懂，适合初学者学习
• kd树(K-dimension tree)是一种对k维空间中的实例点进行存储以便对其进行快速检索的树形数据结构。kd树是是一种二叉树，表示对k维空间的一个划分，构造kd树相当于不断地用垂直于坐标轴的超平面将K维空间切分，构成一...
• 问题, 提出使用KD树作为模式库的存储结构, 能够有效提高搜索速度, 并且能够在实际运行中不断将新发现的交通 流模式实时地加入模式库. 提出使用遗传算法对非参数回归中的重要参数进行优化, 实验表明能够得到相对...
• kd_filter 使用Kd树过滤点噪声 KD-树代码： :
• 利用Java语言实现kd树。利用快速排序找到每一维度的中位数构建Kd树。搜索Kd树时候利用K大堆维护K个最近邻的点。
• 基于Cuda的光线跟踪器。 在GPU上使用kd树以加快渲染速度。 在cuda v6.0上测试。 包括一些样本obj文件。
• 主要介绍了Python语言描述KNN算法与Kd树，具有一定借鉴价值，需要的朋友可以参考下。
• VL库测试程序文件，安装开源的vl库，可以免费使用kd树的matlab函数，很方便。
• 基于Kd树改进的高效K-means聚类算法.pdf
• 基于KD树划分的云计算DBSCAN优化算法.pdf
• 包含kd树的建立以及搜索两部分程序，具体算法过程如下： 给定一个目标点，搜索其最近邻，首先找到包含目标点的叶节点，然后从该叶节点出发，依次退回到其父节点，不断查找是否存在比当前最近点更近的点，直到退回到...
• 用matlab编写的关于Kd树算法，很实用的
• kd树(K-dimension tree)是一种对k维空间中的实例点进行存储以便对其进行快速检索的树形数据结构。kd树是是一种二叉树，表示对k维空间的一个划分，构造kd树相当于不断地用垂直于坐标轴的超平面将K维空间切分，构成一...
• ## KD树

千次阅读 多人点赞 2021-03-06 23:48:47
这里介绍的就是KD树 简介： 为了避免每次都重新计算一遍距离，算法会把距离信息保存在一棵树里，这样在计算之前从树里查询距离信息，尽量避免重新计算。其基本原理是，如果A和B距离很远，B和C距离很近，那么A和C的...

## 问题导入：

k近邻法最简单的实现是线性扫描（穷举搜索），即要计算输入实例与每一个训练实例的距离。计算并存储好以后，再查找K近邻。当训练集很大时，计算非常耗时。
为了提高kNN搜索的效率，可以考虑使用特殊的结构存储训练数据，以减小计算距离的次数。这里介绍的就是KD树

## 简介：

为了避免每次都重新计算一遍距离，算法会把距离信息保存在一棵树里，这样在计算之前从树里查询距离信息，尽量避免重新计算。其基本原理是，如果A和B距离很远，B和C距离很近，那么A和C的距离也很远。有了这个信息，就可以在合适的时候跳过距离远的点。（减少距离值的计算）

## 原理：

构造kd树相当于不断地用垂直于坐标轴的超平面将K维空间切分，构成一系列的K维超矩形区域。kd树的每个结点对应于一个k维超矩形区域。利用kd树可以省去对大部分数据点的搜索，从而减少搜索的计算量。
其原理有点类似于“二分查找”：给出一组数据：[9 1 4 7 2 5 0 3 8]，要查找8。如果挨个查找（线性扫描），那么将会把数据集都遍历一遍。而如果排一下序那数据集就变成了：[0 1 2 3 4 5 6 7 8 9]，按前一种方式我们进行了很多没有必要的查找，现在如果我们以5为分界点，那么数据集就被划分为了左右两个“簇” [0 1 2 3 4]和[6 7 8 9]。

因此，根本就没有必要进入第一个簇，可以直接进入第二个簇进行查找。把二分查找中的数据点换成k维数据点，这样的划分就变成了用超平面对k维空间的划分。空间划分就是对数据点进行分类，“挨得近”的数据点就在一个空间里面。

## kd树的构造：

1.构造根节点
2.通过递归的方法，不断地对k维空间进行切分，生成子节点
3.重复第二步骤，直到子区域中没有示例时终止
需要关注细节：a.选择向量的哪一维进行划分；b.如何划分数据


案例：
给定一个二维空间数据集：T={(2,3),(5,4),(9,6),(4,7),(8,1),(7,2)}，构造一个平衡kd树。
根结点对应包含数据集T的矩形，选择x(1)轴，6个数据点的x(1)坐标中位数是6，这里选最接近的(7,2)点，以平面x(1)=7将空间分为左、右两个子矩形（子结点）；接着左矩形以x(2)=4分为两个子矩形（左矩形中{(2,3),(5,4),(4,7)}点的x(2)坐标中位数正好为4），右矩形以x(2)=6分为两个子矩形，如此递归，最后得到如下图所示的特征空间划分和kd树。

## kd树的查找：

1.二叉树搜索比较待查询节点和分裂节点的分裂维的值，（小于等于就进入左子树分支，大于就进入右子树分支直到叶子结点）
2.顺着“搜索路径”找到最近邻的近似点
3.回溯搜索路径，并判断搜索路径上的结点的其他子结点空间中是否可能有距离查询点更近的数据点，如果有可能，则需要跳到其他子结点空间中去搜索
4.重复这个过程直到搜索路径为空


案例1：查找点(2.1,3.1)的最近邻点
在(7,2)点测试到达(5,4)，在(5,4)点测试到达(2,3)，然后search_path中的结点为<(7,2),(5,4), (2,3)>，从search_path中取出(2,3)作为当前最佳结点nearest, dist为0.141；

然后回溯至(5,4)，以(2.1,3.1)为圆心，以dist=0.141为半径画一个圆，并不和超平面y=4相交，如上图，所以不必跳到结点(5,4)的右子空间去搜索，因为右子空间中不可能有更近样本点了。

于是再回溯至(7,2)，同理，以(2.1,3.1)为圆心，以dist=0.141为半径画一个圆并不和超平面x=7相交，所以也不用跳到结点(7,2)的右子空间去搜索

至此，search_path为空，结束整个搜索，返回nearest(2,3)作为(2.1,3.1)的最近邻点，最近距离为0.141。
案例2：查找点(2,4.5)的最近邻点
在(7,2)处测试到达(5,4)，在(5,4)处测试到达(4,7），然后search_path中的结点为<(7,2),(5,4), (4,7)>，从search_path中取出(4,7)作为当前最佳结点nearest, dist为3.202；

然后回溯至(5,4)，以(2,4.5)为圆心，以dist=3.202为半径画一个圆与超平面y=4相交，所以需要跳到(5,4)的左子空间去搜索。所以要将(2,3)加入到search_path中，现在search_path中的结点为<(7,2),(2, 3)>；另外，(5,4)与(2,4.5)的距离为3.04 < dist = 3.202，所以将(5,4)赋给nearest，并且dist=3.04。

回溯至(2,3)，(2,3)是叶子节点，直接平判断(2,3)是否离(2,4.5)更近，计算得到距离为1.5，所以nearest更新为(2,3)，dist更新为(1.5)

回溯至(7,2)，同理，以(2,4.5)为圆心，以dist=1.5为半径画一个圆并不和超平面x=7相交, 所以不用跳到结点(7,2)的右子空间去搜索

至此，search_path为空，结束整个搜索，返回nearest(2,3)作为(2,4.5)的最近邻点，最近距离为1.5。

展开全文
• kd树（k-dimensional树的简称），是一种分割k维数据空间的数据结构，主要应用于多维空间关键数据的近邻查找(Nearest Neighbor)和近似最近邻查找(Approximate Nearest Neighbor)。
• 本文将从如下几个方面阐述kd树kd树的概念 从例子理解kd树的构建及查找 ball tree 和其他树类型简介 一 kd树的概念 kd（k-dimensional）树的概念自1975年提出，试图解决的是在k维空间为数据集建立索引的问题。...

作者 周知航

审核 张翔

ARGO曾经推出了ARGO：KNN算法介绍​zhuanlan.zhihu.com

介绍了KNN算法的原理，其核心是寻找待测样本在训练样本中的k个近邻，从而从k个近邻中获取数量最多的类别作为待测样本的类别。

KNN算法是机器学习中最简单的算法之一，但是工程实现上，如果训练样本过大，则传统的遍历全样本寻找k近邻的方式将导致性能的急剧下降。 因此，为了优化效率，不同的训练数据存储结构被纳入到实现方式之中。在sikit-learn中的KNN算法参数也提供了'kd_tree'之类的算法选择项。本文将从如下几个方面阐述kd树：

kd树的概念

从例子理解kd树的构建及查找

ball tree 和其他树类型简介

一 kd树的概念

kd(k-dimensional)树的概念自1975年提出，试图解决的是在k维空间为数据集建立索引的问题。依上文所述，已知样本空间如何快速查询得到其近邻？唯有以空间换时间，建立索引便是计算机世界的解决之道。但是索引建立的方式各有不同，kd树只是是其中一种。它的思想如同分治法，即：利用已有数据对k维空间进行切分。

我们先回顾一下二叉查找树，这也是在一维空间中的kd树的情形：

由图可知，二叉树在时间复杂度上是O(logN)，远远优于全遍历算法。对于该树，可以在空间上理解：树的每个节点把对应父节点切成的空间再切分，从而形成各个不同的子空间。查找某点的所在位置时，就变成了查找点所在子空间。上图仅仅是一维，如果换到二维，或者是更高维度上，这棵树该怎么建立？ 又该怎么理解所切成的子空间呢？

我们可以用一个在网络上较为普遍的例子来演示一下二维的kd树构建及查找方式。

二 从例子理解kd树的构建及查找

假设有6个二维数据点{(2,3)，(5,4)，(9,6)，(4,7)，(8,1)，(7,2)}。

依据上文，对于一维空间的切分的最终结果实际上切出的是线段，同理，二维空间中，对于该6个点集最后的切分结果的空间应该是子平面集合。这里，kd树最终目标就是构建一棵树，能够将二维平面切分得到的最终的子平面集，如下图所示：

切分后的空间示意图

找到切分这些平面的直线便是建立kd树的重点。我们将先介绍这种构建方式，然后再介绍基于kd树的最近邻查找的查询方法。

构建

将二维的平面想像成一块方型蛋糕，kd树构建就是面点师要将蛋糕切成上面示意图的模样。先将平面上的六个点在蛋糕上做好标记。面点师举起刀，他该如何下刀呢？这就涉及到kd树在切分时要注意两个要点 ，切分域和切分点的选择。

只要明确了这两个要点的实施规则，切分蛋糕就变得十分容易。因为对于子平面的切分，依然可以按照这些规则递归进行。

a 切分域的选择

切分(split)域，即切分的维度方向，在这个例子中，也就是在二维平面上，面点师必须决定，我这第一刀是平行于x轴横着切还是平行于y轴竖着切？作为面点师，最害怕的情况就是切出的子平面分布不均匀，存在大量过大或者过小的子平面。所以，借助统计学的手段，kd树的每轮中本次的切分域的选择将计算待切分空间内的点此轮未被切分维度上的方差，找方差最大的维度作为此次的切分域。 方差较大，表明在该维度上的点的分散度较高，按该维度切分分辨率比较高。

在这里，一轮是指的必须依照k个维度，对每个维度对平面进行一次切分。比如面点师如果第一刀选择了切分X轴，则第二刀必须选择Y轴。如果多余两个维度，则接下来依然按照剩下维度的方差值进行计算来选择切分域。直至遍历所有维度，一轮结束。

计算第一次切分时两个维度上的方差，得到

在X轴的方向方差较大，所以面点师选择首先切分X轴(也就是在二维平面中沿着Y轴方向)。

b 切分点的选择

切分点的选择策略比较简单，是将带切分平面上的所有数据点按照切分域的维度进行大小排序，选择正中间的点作为切分点。如果正中间点是偶数个数，则任取中间左边或者右边的值作为切分点。所以，对上述X轴的点位进行排序，(5 , 4)，(7 , 2)是X轴方向的中间点。这里选择(7, 2)作为第一次的切分点。

根据已经选择的切分点和切分域，面点师终于畅快的下了第一刀，如下所示。

第一次切分后的结果

依照切分点进行切分后，当前切分点将从子平面的数据集中排除。

c 本轮第二次切分及第二轮切分

本轮并未就此结束，在X轴方向进行切分后将继续对Y方向进行切分。由于面点师第一刀把方形平面切成了两个子平面，对Y轴的切分自然需要对两个子平面各切一刀。这里再次强调，对于二维情况，Y是最后本轮最后维度，不需要再计算方差得到切分域，但是对于多维情况，此次计算是排除第一次切分域的维度后剩下维度在子平面中计算分列域下选择的方差最大的维度。

选择对Y方向上切分的两个子平面的切分点，分别选择(5,4)(也可以是(2,3))(9, 6)(也可以是8,1)。按照切分点和切分域下刀，如下图所示的两条蓝色线条，则是面点师的下刀痕迹。

已切好的平面

在第二轮开始后，由于右上角的平面数据点数为0，不用再切分。而剩下的平面中数据点数都为1。因为此时方差相同，计算方差选择切分域没有意义，使用对X轴进行切分，如上上图的红线所示。

以上只是面点师作为kd树构建时的指导思想，在计算机中该算法究竟该如何表示？接下来看看完成对上述分割时建立的kd树的数据结构。

d 树形数据结构表示

如图所示，这是对应切分方法应该建立的树。第一个根节点(7, 2)是第一次的切分点，切分的是X方向。也就是树的左子树的X轴方向都比7小，右子树X方向点都比7大。再看下一级(5,4)，(9,6)点切分的是Y轴方向，也就是对于(5,4)为根的子树，左子树Y轴方向都比4小，右子树都比4大。以此类推。

值得注意的是，在kd树的实现过程中，切分域的选择并不一定采用方差选择方式，而可能只是简单的以 d mod n 维度的方式，来决定该层树的切分域(d代表树节点的深度，n代表维度数量，求模来确定切分域)。这样实现简单，效果也一般不会太差。

总结一下kd树的算法流程如下

在构建完kd树之后，就是kd树的应用了。如何使用kd树找到目标点的最近邻，kd树是如何帮助实现效率的提升的呢，这时候就得观察一下kd树是如何搜索最近邻的了。

查找

回到面点师切好的糕点平面图，用目标数据在kd树中寻找最近邻时，最核心的两个部分是：

1 寻找近似点-寻找最近邻的叶子节点作为目标数据的近似最近点。

2 回溯-以目标数据和最近邻的近似点的距离沿树根部进行回溯和迭代。

以两个例子来描述这一过程。例子来源：

例子1 搜索点(2.1,3.1)

计算近似最近点--很显然，通过查找寻找kd树，很容易定位到该点将落在点(2,3)所在的叶子节点上。 此时，将(2,3)作为该点的近似最近点。计算出到该近似点的距离为0.141。

为了理解方便，将(2,3)到(2.1,3.1)的距离为直径，以(2,3)为圆心作圆。

回溯--在根据目标节点寻找近似最近点时，其路径是(7,2)->(5,4)->(2,3)。 先计算该点与近似最近点的父节点(5,4)的距离，看看是否该距离是否小于0.141，(5,4)点的距离大于0.141，说明被(5,4) 切分的另外一个子平面与(2,3)点为圆心的圆没有交集。但不能排除跟(7,2)所切开的右子平面没有交集，所以再往上回溯(7,2)点，直至确认跟右子平面没有交集，回溯结束，确认最近邻点位(2,3)。

例子2 搜索点(2, 4.5)

重复上述搜索过程。

计算近似最近点--其落在叶子节点(4,7)点上。近似距离为3.20。

回溯--对于该点的回溯过程就会比刚才复杂一些。其父节点(5,4)与目标点的距离为3.04，小于3.20。 这就意味着(4,7)点作为最近近似点的假设不成立。以(5,4)作为最近近似点。以3.04做圆，圆是跟(5,4)所切分的上下两个平面相交的，所以需要检查(5,4)的另外一个子树的叶子节点(2,3)是否可以作为近似最近点。发现(2，3)的距离为1.5小于3.04，更新(2,3)为近似最近点。最后回溯至(7,2)，确认跟(7,2)切分的右子平面无关。 回溯结束，(2,3)为其最近点。

以上两例说明了为什么kd树能够减少查询次数。但是，例子2中，目标数据点的位置有可能也会要求查询另外一颗子树，加多查询次数，降低查询性能。在维度较高时，这种性能降低十分明显，所以一般kd树要求数据维度在20维以内。通常原则是，如果维度是k, 数据点数是N， 需要满足N >> 2k。

为了应对这一问题，不但出现了对于kd树本身算法的改进如BBF算法，也出现了ball tree一类的针对解决性能下降的适应高维的查找型数据结构。

三 ball tree 和其他树类型介绍

在kd tree 中，我们看到一个导致性能下降的最核心因素是因为kd树的平面是一个个的方形，求最近邻时使用的是圆形。方形平面跟圆形相交的可能性是极高的，如果方形的交汇点多的话，圆形和几个平面相交的可能性就变得更大。这样，凡是相交的平面，都需要进行检查，大大的降低运行效率。

为了解决这一个问题，ball tree抛弃了kd树画的方形，而建立球形，去掉棱角。简而言之，就是使用超球面而不是超矩形划分区域。

它的构造方式是：所有的数据先构成一个最大的球体。从球中选择一个离球的中心最远的点，然后选择第二个点离第一个点最远。将圆中所有离这两个点最近的观测点都赋给这两个簇的中心，而后计算每一个簇的中心点和包含所有其所属观测点的最小半径。不断递归得到最终结果如上图C所示。在上图B中所示数据结构中A，B，C , D, E非叶子节点存储了球的中心点和半径。

它的查找方式是一样先寻找叶子节点，并计算叶子节点和目标点的距离，此值确定目标最小距离上限。根据构造方式可知，可以根据目标点计算出其叶子节点。此时，唯一需要验证的是其兄弟节点的距离，如果小于上限就更新最近目标点。

也能够看出ball tree的构建比kd tree复杂的多，也耗时的多。但是为了节省查询时间这种耗时也是值得的。

对于高维度索引树，还有VP树和MVP树，这里就不再多做解释了。

四 总结

本文重点讲述了kd树的构建和查询，也同时介绍了kd树的改进数据结构ball树。kd树是KNN算法的重要基础，KNN是最简单的机器学习算法之一，而kd树的复杂估计也是读者没有想到的，所谓细节之处皆是魔鬼，此话诚然。

参考文献

好文推荐

kd树就给大家介绍到这了，欢迎大家多多提建议并关注我们的公众号。

微信文章列表索引：往期文章​mp.weixin.qq.com

展开全文
• 之前blog内曾经介绍过SIFT特征匹配算法，特征点匹配和数据库查、图像检索本质上是同一个问题，都可以归结为一个通过距离函数在高维矢量之间进行相似性检索的问题，如何快速而准确地找到查询点的近邻，不少人提出了很...

...