本文轉載自http://blog.csdn.net/zhjchengfeng5/article/details/7855241#
首先來一個問題：
給定平面上一個點集 E ，還有一個定點 V ，怎么在一群點中找出一個點 U，使得 V 與 U 的距離最近（歐幾里得距離）？
當然，我們能夠想到一種做法：枚舉 E 中所有的點，找出它們中距離V 最近的點 U。
但是，假設現在有兩個點集 E1 與 E2 ，對於 E2 中每一個點 Vi ，找出一個在E1 中的一個點 Ui，使得 Vi 到 Ui 的距離最短，這怎么做？還是枚舉？
既然枚舉的復雜度很高 ( O(n) 的復雜度 )，那有沒有辦法把復雜度降下來呢？答案是肯定的，引入一種數據結構：K-D tree
一、何為 K-D tree？
二叉樹（有左兒子，右兒子的那種樹形結構）
二、能解決哪些問題？
K-D tree 可以在 log(n) （ 最壞是 sqrt(n) ）的時間復雜度內求出一個點集 E 中，距離一個定點 V 最近的點（最近鄰查詢），稍稍處理一下，我們還可以求出點集 E 中距離距離 V 最近的 k 個點（k鄰近查詢）
三、怎么利用 K-D tree 解決上面的問題？
將點集 E中的點按照某種規則建成一棵二叉樹，查詢的時候就在這顆建好的二叉樹上面用 log(n) （最壞是 sqrt(n)）的時間復雜度查詢出距離最近的點
四、既然是二叉樹，怎么建樹？
這是最關鍵的地方，因為不管是 划分樹 ， 線段樹 ， 字典樹 ，甚至是其他的數據結構或者算法（例如 KMP 之類的） ，之所以能夠高效的處理問題，主要就是預處理的好。 K-D tree 之所以高效，就是因為建樹很高明，高明之處體現在 “將點集 E中的點按照某種規則建成一棵二叉樹” 的這種規則上
在講這種規則之前，我們先來看看 K-D tree 這種數據結構為什么叫做 K-D tree
K：K鄰近查詢中的k
D：空間是D維空間（Demension）
tree：你可以理解為是二叉樹，也可以單純的看做是一顆 tree
好了， K 我們已經用到了，tree 我們也已經用到了，但是 D 呢？貌似這篇文章到現在為止還沒有提到過 D 吧？
這種規則，就是針對空間的“維”的
既然要建樹，那么樹上的節點肯定要定義一些狀態：
節點的狀態：
分裂點（split_point）
分裂方式（split_method）
左兒子（left_son）
右兒子（right_son）
我們建樹的規則就是節點的狀態中的：分裂方式（split_method）
想必讀者已經看見上面的關鍵字了：分裂點 分裂方式，為什么反復的出現分裂這兩個字呢？難道建一顆 K-D tree 還要分裂什么，分裂空間？
對，K-D tree的建立就是分裂空間的過程！
怎么建樹呢？
建樹依據：
先計算當前區間 [ L , R ] 中（這里的區間是點的序號區間，而不是我們實際上的坐標區間），每個點的坐標的每一維度上的方差，取方差最大的那一維，設為 d，作為我們的分裂方式（split_method ），把區間中的點按照在 d 上的大小，從小到大排序，取中間的點 sorted_mid 作為當前節點記錄的分裂點，然后，再以 [ L , sorted_mid-1 ] 為左子樹建樹 ， 以 [sorted_mid+1 , R ] 為右子樹建樹，這樣，當前節點的所有狀態我們便確定下來了：
split_point= sorted_mid
split_method= d
left_son = [ L , sorted_mid-1 ]
right_son = [ sorted_mid+1 , R ]
為了便於理解，我先舉個例子：
假設現在我們有平面上的點集 E ，其中有 5 個二維平面上的點 ： （1,4）（5,8） （4,2） （7,9） （10，11）
它們在平面上的分布如圖：
首先，我們對區間 [ 1 , 5 ] 建樹：
先計算區間中所有點在第一維（也就是 x 坐標）上的方差：
平均值 ： ave_1 =5.4
方差 ： varance_1 =9.04
再計算區間中所有點在第二維（也就是 y 坐標）上的方差：
平均值：ave_2 =6.8
方差：varance_2 =10.96
明顯看見，varance_2 > varance_1 ，那么我們在本次建樹中，分裂方式 ：split_method =2 ， 再將所有的點按照 第 2 維 的大小從小到大排序，得到了新的點的一個排列：
（4,2） （1,4）（5,8） （7,9） （10,11）
取中間的點作為分裂點 sorted_mid =（5，8）作為根節點，再把區間 [ 1 , 2] 建成左子樹 , [ 4 , 5] 建成右子樹，此時，直線 ： y = 8 將平面分裂成了兩半，前面一半給左兒子，后面一半給了右兒子，如圖：
建左子樹 [1 , 3 ] 的時候可以發現，這時候是 第一維 的方差大 ，分裂方式就是1 ，把區間 [ 1, 2 ] 中的點按照 第一維 的大小，從小到大排序 ，取中間點（1,4） 根節點，再以區間 [ 2, 2] 建立右子樹 得到節點 （4,2）
建右子樹 [4 , 5 ] 的時候可以發現，這時還是 第一維 的方差大， 於是，我們便得到了這樣的一顆二叉樹 也就是 K-D tree，它把平面分成了如下的小平面，使得每個小平面中最多有一個點：
可以看見，我們實際上在建樹的過程中，把整個平面分成了 4 個部分
樹是建了，那么查詢呢?
查詢過程：
查詢，其實相當於我們要將一個點“添加”到已經建好的 K-D tree 中，但並不是真的添加進去，只是找到他應該處於的子空間即可，所以查詢就顯得簡單的毒攻了
每次在一個區間中查詢的時候，先看這個區間的分裂方式是什么，也就是說，先看這個區間是按照哪一維來分裂的，這樣如果這個點對應的那一維上面的值比根節點的小，就在根節點的左子樹上進行查詢操作，如果是大的話，就在右子樹上進查詢操作
每次回溯到了根節點（也就是說，對他的一個子樹的查找已經完成了）的時候，判斷一下，以該點為圓心，目前找到的最小距離為半徑，看是否和分裂區間的那一維所構成的平面相交，要是相交的話，最近點可能還在另一個子樹上，所以還要再查詢另一個子樹，同時，還要看能否用根節點到該點的距離來更新我們的最近距離。為什么是這樣的，我們可以用一幅圖來說明：
在查詢到左兒子的時候，我們發現，現在最小的距離是 r = 10 ，當回溯到父親節點的時候，我們發現，以目標點（10,1）為圓心，現在的最小距離 r = 10 為半徑做圓，與分割平面 y = 8 相交，這時候，如果我們不在父親節點的右兒子進行一次查找的話，就會漏掉 （10,9） 這個點，實際上，這個點才是距離目標點 （10,1） 最近的點
由於每次查詢的時候可能會把左右兩邊的子樹都查詢完，所以，查詢並不是簡單的 log(n) 的，最壞的時候能夠達到 sqrt(n)
好了，到此，K-D tree 就差不多了，寫法上與很多值得優化的地方，至於怎么把最鄰近查詢變換到 K 鄰近查詢，我們用一個數組記錄一個點是否可以用來更新最近距離即可，下面貼上 K-D tree 一個模板
#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <vector>
#include <string>
#include <queue>
#include <stack>
#define INT_INF 0x3fffffff
#define LL_INF 0x3fffffffffffffff
#define EPS 1e-12
#define MOD 1000000007
#define PI 3.141592653579798
#define N 60000
using namespace std;
typedef long long LL;
typedef unsigned long long ULL;
typedef double DB;
struct data
{
LL pos[10];
int id;
} T[N] , op , point;
int split[N],now,n,demension;
bool use[N];
LL ans,id;
DB var[10];
bool cmp(data a,data b)
{
return a.pos[split[now]]<b.pos[split[now]];
}
void build(int L,int R)
{
if(L>R) return;
int mid=(L+R)>>1;
//求出 每一維 上面的方差
for(int pos=0;pos<demension;pos++)
{
DB ave=var[pos]=0.0;
for(int i=L;i<=R;i++)
ave+=T[i].pos[pos];
ave/=(R-L+1);
for(int i=L;i<=R;i++)
var[pos]+=(T[i].pos[pos]-ave)*(T[i].pos[pos]-ave);
var[pos]/=(R-L+1);
}
//找到方差最大的那一維，用它來作為當前區間的 split_method ，split_method保存在split[mid]中
split[now=mid]=0;
for(int i=1;i<demension;i++)
if(var[split[mid]]<var[i]) split[mid]=i;
//對區間排排序，找到中間點
nth_element(T+L,T+mid,T+R+1,cmp);
build(L,mid-1);
build(mid+1,R);
}
//bulid過后split{i]表示以i結點為中心點的分裂方式
void query(int L,int R)
{
if(L>R) return;
int mid=(L+R)>>1;
//求出目標點 op 到現在的根節點的距離
LL dis=0;
for(int i=0;i<demension;i++)
dis+=(op.pos[i]-T[mid].pos[i])*(op.pos[i]-T[mid].pos[i]);
//如果當前區間的根節點能夠用來更新最近距離，並且 dis 小於已經求得的 ans
if(!use[T[mid].id] && dis<ans)
{
ans=dis; //更新最近距離
point=T[mid]; //更新取得最近距離下的點
id=T[mid].id; //更新取得最近距離的點的 id
}
//計算 op 到分裂平面的距離
LL radius=(op.pos[split[mid]]-T[mid].pos[split[mid]])*(op.pos[split[mid]]-T[mid].pos[split[mid]]);
//對子區間進行查詢
if(op.pos[split[mid]]<T[mid].pos[split[mid]])
{
query(L,mid-1);
if(radius<=ans) query(mid+1,R);
}
else
{
query(mid+1,R);
if(radius<=ans) query(L,mid-1);
}
}
int main()
{
while(scanf("%d%d",&n,&demension)!=EOF)
{
//讀入 n 個點
for(int i=1;i<=n;i++)
{
for(int j=0;j<demension;j++)
scanf("%I64d",&T[i].pos[j]);
T[i].id=i;
}
build(1,n); //建樹
int m,q; scanf("%d",&q); // q 個詢問
while(q--)
{
memset(use,0,sizeof(use));
for(int i=0;i<demension;i++)
scanf("%I64d",&op.pos[i]);
scanf("%d",&m);
printf("the closest %d points are:\n",m);
while(m--)
{
ans=(((LL)INT_INF)*INT_INF);
query(1,n);
for(int i=0;i<demension;i++)
{
printf("%I64d",point.pos[i]);
if(i==demension-1) printf("\n");
else printf(" ");
}
use[id]=1;
}
}
}
return 0;
}