Kd-Tree讲解--（HDUOJ-4347：The Closest M Points）

 

The Closest M Points Time Limit: 16000/8000 MS (Java/Others)    Memory Limit: 98304/98304 K (Java/Others) Total Submission(s): 8330    Accepted Submission(s): 2704   Problem Description The course of Software Design and Development Practice is objectionable. ZLC is facing a serious problem .There are many points in K-dimensional space .Given a point. ZLC need to find out the closest m points. Euclidean distance is used as the distance metric between two points. The Euclidean distance between points p and q is the length of the line segment connecting them.In Cartesian coordinates, if p = (p1, p2,..., pn) and q = (q1, q2,..., qn) are two points in Euclidean n-space, then the distance from p to q, or from q to p is given by: Can you help him solve this problem?     Input In the first line of the text file .there are two non-negative integers n and K. They denote respectively: the number of points, 1 <= n <= 50000, and the number of Dimensions,1 <= K <= 5. In each of the following n lines there is written k integers, representing the coordinates of a point. This followed by a line with one positive integer t, representing the number of queries,1 <= t <=10000.each query contains two lines. The k integers in the first line represent the given point. In the second line, there is one integer m, the number of closest points you should find,1 <= m <=10. The absolute value of all the coordinates will not be more than 10000. There are multiple test cases. Process to end of file.     Output For each query, output m+1 lines: The first line saying :”the closest m points are:” where m is the number of the points. The following m lines representing m points ,in accordance with the order from near to far It is guaranteed that the answer can only be formed in one ways. The distances from the given point to all the nearest m+1 points are different. That means input like this: 2 2 1 1 3 3 1 2 2 1 will not exist.   Sample Input 3 2 1 1 1 3 3 4 2 2 3 2 2 3 1 Sample Output the closest 2 points are: 1 3 3 4 the closest 1 points are: 1 3

题意：给你一个n维的点m，求图中所有已知点中离m最近的k个点？

题解：kd-tree模板题啊，但是不了解kd树的孩子可能只能想到暴力遍历的方式，但是这道题图中有50000的点，并且查询次数可能有10000次这么多，暴力怎么可能过的了。。。

kd-Tree是本科阶段就应该学习的，但是因为懒就一直晾着，最近在搞机器学习，碰到了kNN算法，这里不再科普这个东西了，反正就是想解决这个一直没有学的算法。。。。

本讲解参考：https://www.jianshu.com/p/abcaaf754f92（我觉得他讲解的很棒！）

好了言归正传，我先介绍完kdtree相关知识，再去贴这道题的代码，相信你很快就能读懂代码。

什么是Kd-Tree？

Kd-Tree，即K-dimensional tree，是一种高维索引树形数据结构，常用于在大规模的高维数据空间进行最近邻查找和近似最近邻查找，例如：图像检索识别中的高维图像特征向量的K近邻查找与匹配；还有我们众所周知的机器学习里的kNN算法，一开始要找最近的k个点，相信封装好的kNN的优化一定是可以用kd-Tree的，因此，学习kd-Tree真的很有必要哦~

首先我们要明确，kd-Tree本质上是一颗二叉树，书中存储的是一些K维数据。那么我们在一个K维数据集合上构建一棵kd-Tree代表了对该K维数据集合构成的K维空间的一个划分，即树中的每个结点就对应了一个K维的超矩形区域。

二叉树的概念我们不再回顾，相信学过数据结构的童鞋这些知识都不在话下（其实是我真的懒得讲了哈哈哈），我们脑海里过一遍二叉查找树，然后再往下看，你一定觉得so easy！

此时我们假设我们要处理的对象集合是一个K维空间中的数据集，那么我们是否可以构建一棵类似于1维空间中的二叉查找树呢，答案是肯定的，知识推广到K维空间后，创建二叉树和查询二叉树的算法会有一些变化。

如何构造一棵Kd-Tree？

对于Kd-Tree这样一棵二叉树，我们首先需要确定怎样划分左子树和右子树，即一个K维数据是依据什么被划分到左子树或右子树的。

在构造1维BST树时，一个1维数据根据其与树的根结点和中间结点进行大小比较的结果来决定是划分到左子树还是右子树，同理，我们也可以按照这样的方式，将一个K维数据与Kd-tree的根结点和中间结点进行比较，只不过不是对K维数据进行整体的比较，而是选择某一个维度Di，然后比较两个K维数在该维度Di上的大小关系，即每次选择一个维度Di来对K维数据进行划分，相当于用一个垂直于该维度Di的超平面将K维数据空间一分为二，平面一边的所有K维数据在Di维度上的值小于平面另一边的所有K维数据对应维度上的值。也就是说，我们每选择一个维度进行如上的划分，就会将K维数据空间划分为两个部分，如果我们继续分别对这两个子K维空间进行如上的划分，又会得到新的子空间，对新的子空间又继续划分，重复以上过程直到每个子空间都不能再划分为止。以上就是构造Kd-Tree的过程，上述过程中涉及到两个重要的问题：

1.每次对子空间的划分时，怎样确定在哪个维度上进行划分。

2.在某个维度上进行划分时，怎样确保在这一维度上的划分得到的两个子集合的数量尽量相等，即左子树和右子树中的结点个数尽量相等。

对于第一个问题：引入解法有点啰嗦（感兴趣的可以看我前边贴的链接），我这边直接告诉你正确的划分方法：

如果一个K维数据集合的分布像木条一样，那就是说明这K维数据在木条较长方向代表的维度上，这些数据的分布散得比较开，数学上来说，就是这些数据在该维度上的方差（invariance）比较大，换句话说，正因为这些数据在该维度上分散的比较开，我们就更容易在这个维度上将它们划分开，因此，这就引出了我们选择维度的另一种方法：最大方差法（max invarince），即每次我们选择维度进行划分时，都选择具有最大方差维度。

对于第二个问题：

假设当前我们按照最大方差法选择了在维度i上进行K维数据集S的划分，此时我们需要在维度i上将K维数据集合S划分为两个子集合A和B，子集合A中的数据在维度i上的值都小于子集合B中。首先考虑最简单的划分法，即选择第一个数作为比较对象（即划分轴，pivot），S中剩余的其他所有K维数据都跟该pivot在维度i上进行比较，如果小于pivot则划A集合，大于则划入B集合。把A集合和B集合分别看做是左子树和右子树，那么我们在构造一个二叉树的时候，当然是希望它是一棵尽量平衡的树，即左右子树中的结点个数相差不大。而A集合和B集合中数据的个数显然跟pivot值有关，因为它们是跟pivot比较后才被划分到相应的集合中去的。好了，现在的问题就是确定pivot了。给定一个数组，怎样才能得到两个子数组，这两个数组包含的元素个数差不多且其中一个子数组中的元素值都小于另一个子数组呢？方法很简单，找到数组中的中值（即中位数，median），然后将数组中所有元素与中值进行比较，就可以得到上述两个子数组。同样，在维度i上进行划分时，pivot就选择该维度i上所有数据的中值，这样得到的两个子集合数据个数就基本相同了。

解决了上面两个重要的问题后，就得到了Kd-Tree的构造算法了。

Kd-Tree的构造算法

（1） 在K维数据集合中选择具有最大方差的维度k，然后在该维度上选择中值m为pivot对该数据集合进行划分，得到两个子集合；同时创建一个树结点node，用于存储；

（2）对两个子集合重复（1）步骤的过程，直至所有子集合都不能再划分为止；如果某个子集合不能再划分时，则将该子集合中的数据保存到叶子结点（leaf node）。

利用Kd-Tree进行最近邻查找的算法

（1）将查询数据Q从根结点开始，按照Q与各个结点的比较结果向下访问Kd-Tree，直至达到叶子结点。

其中Q与结点的比较指的是将Q对应于结点中的k维度上的值与m进行比较，若Q(k) < m，则访问左子树，否则访问右子树。达到叶子结点时，计算Q与叶子结点上保存的数据之间的距离，记录下最小距离对应的数据点，记为当前“最近邻点”Pcur和最小距离Dcur。

（2）进行回溯（Backtracking）操作，该操作是为了找到离Q更近的“最近邻点”。即判断未被访问过的分支里是否还有离Q更近的点，它们之间的距离小于Dcur。

如果Q与其父结点下的未被访问过的分支之间的距离小于Dcur，则认为该分支中存在离P更近的数据，进入该结点，进行（1）步骤一样的查找过程，如果找到更近的数据点，则更新为当前的“最近邻点”Pcur，并更新Dcur。

如果Q与其父结点下的未被访问过的分支之间的距离大于Dcur，则说明该分支内不存在与Q更近的点。

回溯的判断过程是从下往上进行的，直到回溯到根结点时已经不存在与P更近的分支为止。

怎样判断未被访问过的树分支Branch里是否还有离Q更近的点？

从几何空间上来看，就是判断以Q为中心center和以Dcur为半径Radius的超球面（Hypersphere）与树分支Branch代表的超矩形（Hyperrectangle）之间是否相交。

在实现中，我们可以有两种方式来求Q与树分支Branch之间的距离。第一种是在构造树的过程中，就记录下每个子树中包含的所有数据在该子树对应的维度k上的边界参数[min, max]；第二种是在构造树的过程中，记录下每个子树所在的分割维度k和分割值m，（k, m），Q与子树的距离则为|Q(k) - m|。

以上就是Kd-tree的构造过程和基于Kd-Tree的最近邻查找过程。

好了，在我们了解了整个过程后，我们来看下这道题的代码，相信你一定觉得so easy！

#include<queue>
#include<stdio.h>
#include<algorithm>
using namespace std;
#define maxn 50085
#define ll long long
#define K 5
struct kdnode
{
	ll x[K];
	int div;
	bool lef;
}arr[12];
struct node
{
	kdnode a;
	ll dis;
	bool operator < (const node &a) const
	{
		return dis<a.dis;
	}
	node(){}
	node(kdnode &tmp,ll d)
	{
		a=tmp;
		dis=d;
	}
};
priority_queue<node> qq;
int num,nowNum,m;
kdnode p[maxn],q;
ll ans;

int cmpNo;
bool cmp(kdnode a,kdnode b)
{
	return a.x[cmpNo]<b.x[cmpNo];
}

inline ll max(ll a,ll b)
{
	return a>b?a:b;
}
ll dis(kdnode a,kdnode b,int k)
{
	ll res=0;
	for(int i=0;i<k;i++)
		res+=(a.x[i]-b.x[i])*(a.x[i]-b.x[i]);
	return res;
}

void buildKDTree(int l,int r,kdnode *p,int d,int k)
{
	if(l>r) return;
	int m=(l+r)/2;
	cmpNo=d;
	//nth_element算法仅排序第nth个元素（从0开始的索引）.对给定范围[first,last)内的元素进行重新布置.
	//方法是,nth位置的元素放置的值就是把所有元素排序后在nth位置的值.把所有不大于nth的值放到nth的前面,
	//把所有不小于nth的值放到nth后面.对给定范围内的元素"排序"
	nth_element(p+l,p+m,p+r+1,cmp);
	p[m].div=d;
	if(l==r)
	{
		p[m].lef=1;
		return;
	}
	buildKDTree(l,m-1,p,(d+1)%k,k);
	buildKDTree(m+1,r,p,(d+1)%k,k);
}
void findKD(int l,int r,kdnode &tar,kdnode *p,int k)
{
	if(l>r) return;
	int m=(l+r)/2;
	ll d=dis(p[m],tar,k);
	if(p[m].lef)
	{
		if(nowNum<num)
		{
			nowNum++;
			ans=max(ans,d);
			qq.push(node(p[m],d));
		}
		else if(ans>d)
		{
			qq.pop();
			qq.push(node(p[m],d));
			ans=qq.top().dis;
		}
		return;
	}
	ll t=tar.x[p[m].div]-p[m].x[p[m].div];
	if(t>0)
	{
		findKD(m+1,r,tar,p,k);
		if(nowNum<num)
		{
			qq.push(node(p[m],d));
			nowNum++;
			ans=qq.top().dis;
			findKD(l,m-1,tar,p,k);
		}
		else
		{
			if(ans>d)
			{
				qq.pop();
				qq.push(node(p[m],d));
				ans=qq.top().dis;
			}
			if(ans>t*t)
				findKD(l,m-1,tar,p,k);
		}
	}
	else
	{
		findKD(l,m-1,tar,p,k);
		if(nowNum<num)
		{
			qq.push(node(p[m],d));
			nowNum++;
			ans=qq.top().dis;
			findKD(m+1,r,tar,p,k);
		}
		else
		{
			if(ans>d)
			{
				qq.pop();
				qq.push(node(p[m],d));
				ans=qq.top().dis;
			}
			if(ans>t*t)
				findKD(m+1,r,tar,p,k);
		}
	}
}
int main(void)
{
	int n,k;
	while(scanf("%d%d",&n,&k)==2)
	{
		for(int i=0;i<n;i++)
		{
			for(int j=0;j<k;j++)
				scanf("%lld",&p[i].x[j]);
			p[i].lef=0;
		}
		int t;
		scanf("%d",&t);
		buildKDTree(0,n-1,p,k-1,k);
		while(t--)
		{
			ans=-1;
			nowNum=0;
			for(int i=0;i<k;i++)
				scanf("%lld",&q.x[i]);
			while(!qq.empty()) qq.pop();
			scanf("%d",&num);
			findKD(0,n-1,q,p,k);
			for(int i=1;i<=num;i++)
			{
				arr[i]=qq.top().a;
				qq.pop();
			}
			printf("the closest %d points are:\n",num);
			for(int i=num;i>=1;i--)
			{
				for(int j=0;j<k;j++)
				{
					if(j==0)
						printf("%lld",arr[i].x[j]);
					else
						printf(" %lld",arr[i].x[j]);
				}
				printf("\n");
			}
		}
	}
	return 0;
}