寻找最近点对问题（FCPP find the closest pair of point）

最新推荐文章于 2024-05-22 09:32:36 发布

WsHunTer

最新推荐文章于 2024-05-22 09:32:36 发布

阅读量6.3k

点赞数

分类专栏：算法-分治文章标签： pair float c 算法 ini n2

本文链接：https://blog.csdn.net/WsHunTer/article/details/6777834

版权

算法-分治专栏收录该内容

1 篇文章 0 订阅

订阅专栏

寻找最近点对问题（FCPP）

问题描述：

给定平面上n个点，找出距离最近的两个点。

思考过程：

   1）对于这种问题，我们首先想到的求解方法就是求出所有点对的距离，并找出最近的那个，当然这个是个显而

易见的方法，具体过程大体可以警醒如下描述。

① 定义变量

N (点的数目)，X[N] (点的x坐标值)，Y[N]（点的y坐标值），closest（最近的值），P（最近点对的一个点），Q（最近点对的另一个点）。

② 依次扫描并比较出最近点对

初始化closest=get_len(X[1],Y[1],X[2],Y[2]);

P=(X[1],Y[1]);Q=(X[2],Y[2]);

For i=1:N

For j=i+1:N

Tmp=get_len(X[i],Y[i],X[j],Y[j]);

If Tmp < closest

Closest=Tmp;

P=(X[i],Y[i]);

Q= (X[j],Y[j]);

End If

End For

这样运行后就可以求出所有的点中的最近的点对值，然而这是一种Brute-force（蛮力法），当然时间复杂度也是非常不能令人满意的Q(n2)。

2）对于1）中的算法基本上很简单而且很容易想到，的确，对于一些问题，人们往往想的简单做的难，想的难就做的简单啦(呵呵)。在此，可以利用分治法来对问题求解。具体算法步骤如下：

① 将所有的点按X坐标值排序，然后平均等分成两个部分PL与PR。在左右分开处画个分割线part_line。（如下图）

② 找出PL中的最近值与PR中的最近值分别为dL，dR。

③ 令d=min（dL，dR）；现在，问题来了，d一定最小么？

回答是相当确定的，就是肯定不是最小的。原因也显而易见，就是最近点对可能出现在PL与PR的分割处。

④ 根据③中的问题，可以以part_line为中心左右|d|的距离里的点考虑（为什么？超过这个距离显然不需要考虑啦，因为距离一定会大于d啦~~~）。现在，问题又来了，如何取得这个范围内的最近点呢？蛮力？全搜索计算？对！不行的这是，可能所有的点都落在这个范围里，就此一项我的时间复杂度就又会回到Q(n2)。怎么办？

⑤ 对于④中的问题，解决方法是这样的，将左右|d|范围里的点按Y值排序，然后依次从每个点出发做水平线，并在其之上d的距离做水平线，这样得到了一个2d*d的矩形，显然，在此矩形之外的点无需开率啦，可以证明，矩形内最多容下8个点（包括自己本身）。所以意味着在排好序的Y坐标点中，只需要考虑其后的7个点就行了，此外已有人证明了事实上只需要考虑其后的4个点就行啦。可见，这些的时间复杂度最多Q（N）。

⑥ 纵观整个过程，首先要对X坐标排序时间为O（NlogN），这是分治之前的预处理。然后问题分成了N/2规模，然后用Q（N）的复杂度得到中间的最近点对，然后可以在常数时间内得到最终的点对与最近值。So。。。T（N）=2T（N/2）+O（N），可以得到该算法的时间为O（NlogN）（即使加上一开始的排序预处理的时间复杂度）。

算法具体实现过程：

①定义函数 FCPP（st，ed，X，Y，IND，d，p，q）；

其中：

st，ed分别是X数组的最小和最大下标

X[]和Y[] 这个不用说，就是存放点的X和Y坐标的数组

IND[]是存放Y坐标对应X坐标的对应关系，比如(X[IND[i]],Y[i])表示一个点。（以上是输入的变量）

δ表示最近的值,p,q分别是最近点对的点（以上是输出变量）

②具体运行过程为：

FCPP(j,k,X,Ypres,INDpres,d,p,q)

{

令m= k-j+1；//即m是X[j]到X[k]这一段中的数据元素的个数。

if(m≤3){

根据Ypres，INDpres和X数组找出全部3个（或2个）点，

求出其中的最近点对(p,q)及最小距离d，然后返回d,p,q。

}

if(m≥4){

进行分治前的准备工作：

引入4个数组YL，INDL（长度均为ém/2ù），

YR，INDR（长度均为ëm/2û），

YL存放左边点集中点的y坐标（按从小到大的次序），

INDL存放与之对应的x坐标在X数组中的下标，

YR和INDR的作用与YL和INDL类似。

LengthL←0；LengthR←0;/* LengthL记录放入左边数组YL中的点数*/

For i=1 to m do /*对m个点逐一进行检查，Q(m)时间*/

{

if INDpres[i]≤j+m/2-1 /*当前检查的点应分在左边的点集中*/

{ LengthL增1； /*左边的点数增1*/

YL[LengthL]←Ypres [i]；/*当前检查点的y坐标放在左边Y数组中

INDL[LengthL]←INDpres [i]；/*在左边索引数组中保持y坐标与x坐标的对应关系*/

} else { /*当前被检查的点应分在右边的点集中*/

LengthR增1； /*右边的点数增1*/

YR[LengthR]←Ypres[i]；/*当前检查点的y坐标放在右边Y数组中*/

INDR[LengthR]←INDpres [i]；/*在右边索引数组中保持y坐标与x坐标的对应关系*/

} /*当前检查点已放在左边或右边的Y数组中*/

} /*m个点检查完毕时，*/

/*YL和YR仍然保持了从小到大的次序*/

/*然后分治：*/

FCPP(j, j+m/2-1, X, YL, INDL,dL,pL,qL)；

FCPP(j+m/2, k, X, YR, INDR,dR,pR,qR)；

/*以上两个调用执行时间为2T(m/2)*/

令d=min{dL，dR }并记录下相应的(p,q)； /*O(1)时间*/

/*求分割线的x坐标值：*/

xdiv=(X[j+m/2-1]+ X[j+m/2])/2； /*O(1)时间*/

/*引入数组Yd和INDd，长度最大为m，*/

/*Yd放带状域内点的y坐标，INDd放对应的x的下标*/

LengthYd←0；

For i=1 to m do /*对m个点逐一进行检查，Q(m)时间*/

{ if (xdiv-d < X[ INDpres[i] ] < xdiv+d)

/*当前点落到2d带状区域内*/

{ LengthYd增1；Yd [LengthYd]←Ypres [i]；

/*记录下该点的y坐标*/

INDd [LengthYd]←INDpres [i]；}

/*保持y坐标与x坐标的对应关系*/

} /*m个点检查完毕时，*/

/*Yd中的y坐标仍然保持了从小到大的次序*/

if LengthYd≤1 then return；/*只有一个或没有点落到带状区域内，无须算*/

For i=1 to LengthYd-1 do /* 至少有两个点落到*/

/*带状区域内时对每个点要逐一检查，*/

/*最后一个点因其后无点，故该点无需检查*/

{ ypres←Yd [i]；xpres←X[INDd[i]]；

u←min{i+4, LengthYd}； /*u是本次检查的最大下标*/

v←i+1；/* v是本次检查的游动变量，从i+1到u为止*/

while (v≤u & (Yd[v] - ypres) < d) /*y坐标值之差大于等于d就停止检查，因其后不可能有*/

{

计算点(xpres,ypres)与点(X[INDd[v]], Yd[v])之间的距离d’；if d’< d then {d←d’}；

点p,q更新为(xpres,ypres)和(X[INDd[v]], Yd[v])；}

v增1； /*v所指的点检查完毕*/

}/*该while循环最多执行4次，故为O(1)时间*/

} /* LengthYd最大为m，故For循环执行时间最多为O(m)*/

}//end for

}//end if m>=4

}//end FCPP

实现的代码为（C语言）：

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define    MAXLEN	200


int init(float *x,float *y,int *ind)
{
	int n,i;
	while(1==scanf("%d",&n))
	{
		if(n > 1) break;
		printf("n invalid!\n");
	}

	for(i=0;i<n;i++)
	{
		scanf("%f%f",x+i,y+i);
		ind[i]=i;
	}

	return n;
}

void pre_sort(float *x,float *y,int *ind,int n)
{//must finish in O(nlogn) but here using normorl sort method
	//bubble sort
	int i,j;
	int tmp;
	float tmp_y,tmp_x;

	for(i=0;i<n;i++)
		for(j=n-1;j > i;j--)
			if(x[j-1] > x[j])
			{
				//swap x
				tmp_x=x[j-1];
				x[j-1]=x[j];
				x[j]=tmp_x;
				//swap y
				tmp_y=y[j-1];
				y[j-1]=y[j];
				y[j]=tmp_y;
				//swap ind no need now!
				/*tmp=ind[j-1];
				ind[j-1]=ind[j];
				ind[j]=tmp;*/
			}
	for(i=0;i<n;i++)
		for(j=n-1;j > i;j--)
			if(y[j-1] > y[j])
			{
				//swap y
				tmp_y=y[j-1];
				y[j-1]=y[j];
				y[j]=tmp_y;
				//swap ind
				tmp=ind[j-1];
				ind[j-1]=ind[j];
				ind[j]=tmp;
			}
}

float get_pair_len(float x1,float y1,float x2,float y2)
{
	return sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
}

void find_closest_pair(int st,int ed,float *x,float *y,int *ind,float *closest,float *x1,float *y1,float *x2,float *y2)
{//key code
	if(ed-st == 1){//only 2 point
		*closest=get_pair_len(x[ind[0]],y[0],x[ind[1]],y[1]);
		*x1=x[ind[0]];
		*y1=y[0];
		*x2=x[ind[1]];
		*y2=y[1];
	}
	else if(ed-st ==2)
	{//only 3 point
		float len;
		*closest=get_pair_len(x[ind[0]],y[0],x[ind[1]],y[1]);
		*x1=x[ind[0]];
		*y1=y[0];
		*x2=x[ind[1]];
		*y2=y[1];

		len=get_pair_len(x[ind[0]],y[0],x[ind[2]],y[2]);
		if(len < *closest){
			*closest=len;
			*x1=x[ind[0]];
			*y1=y[0];
			*x2=x[ind[2]];
			*y2=y[2];
		}

		len=get_pair_len(x[ind[1]],y[1],x[ind[2]],y[2]);
		if(len < *closest){
			*closest=len;
			*x1=x[ind[1]];
			*y1=y[1];
			*x2=x[ind[2]];
			*y2=y[2];
		}
	}else{//at least 4 points
		int i,cl,cr,ct;
		int mid;
		float t_c1,t_c2;
		float t_c1_x1,t_c1_y1,t_c1_x2,t_c1_y2;
		float t_c2_x1,t_c2_y1,t_c2_x2,t_c2_y2;
		float *yl,*yr,*yct;
		int *indl,*indr,*indct;

		mid=(st+ed)/2;
		yl=(float*)malloc((mid-st+1)*sizeof(float));
		indl=(int*)malloc((mid-st+1)*sizeof(int));
		yr=(float*)malloc((ed-mid)*sizeof(float));
		indr=(int*)malloc((ed-mid)*sizeof(int));

		if(!yl || !yr || !indl || !indr) exit(-2);//non enough mm

		for(i=0,cl=0,cr=0;i<=ed-st;i++)
		{//store new order y&ind
			if(ind[i] <= mid){
				yl[cl]=y[i];
				indl[cl]=ind[i];
				cl++;
			}
			else{
				yr[cr]=y[i];
				indr[cr]=ind[i];
				cr++;
			}
		}

		//divided
		find_closest_pair(st,mid,x,yl,indl,&t_c1,&t_c1_x1,&t_c1_y1,&t_c1_x2,&t_c1_y2);
		find_closest_pair(mid+1,ed,x,yr,indr,&t_c2,&t_c2_x1,&t_c2_y1,&t_c2_x2,&t_c2_y2);

		if(t_c1 < t_c2){
			*closest=t_c1;
			*x1=t_c1_x1;
			*y1=t_c1_y1;
			*x2=t_c1_x2;
			*y2=t_c1_y2;
		}
		else{
			*closest=t_c2;
			*x1=t_c2_x1;
			*y1=t_c2_y1;
			*x2=t_c2_x2;
			*y2=t_c2_y2;
		}//get the closest pair between the two part(L&R)

		int v,ve;
		float closest_tmp;
		float part_line=(x[mid]+x[mid+1])/2.0;

		yct=(float*)malloc((ed-st+1)*sizeof(float));
		indct=(int*)malloc((ed-st+1)*sizeof(int));

		//get all the y in the rect 
		for(i=0,ct=0;i<=ed-st;i++)
			if(x[ind[i]]>part_line-*closest && x[ind[i]]<part_line+*closest ){
				yct[ct]=y[i];
				indct[ct]=ind[i];
				ct++;
			}

		//check followed the 4 points 
		if(ct > 1){
			for(i=0;i< ct-1;i++){
				ve=(i+4 < ct-1 ? i+4:ct-1);
				v=i+1;
				while(v<=ve && yct[v]-yct[i] <=*closest){
					closest_tmp=get_pair_len(x[indct[i]],yct[i],x[indct[v]],yct[v]);
					if(closest_tmp < *closest){
						*closest=closest_tmp;
						*x1=x[indct[i]];
						*y1=yct[i];
						*x2=x[indct[v]];
						*y2=yct[v];
					}
					v++;
				}
			}
		}

		free(yl);
		free(yr);
		free(indl);
		free(indr);
		free(yct);
		free(indct);
		
	}
}

int main()
{
	float X[MAXLEN],Y[MAXLEN];
	int IND[MAXLEN],n;
	n=init(X,Y,IND);

	float closest,x1,y1,x2,y2;
	
	pre_sort(X,Y,IND,n);
	
	find_closest_pair(0,n-1,X,Y,IND,&closest,&x1,&y1,&x2,&y2);
	printf("the closest value is:%-4.2f (%-4.2f,%-4.2f)  (%-4.2f,%-4.2f)\n",closest,x1,y1,x2,y2);

	return 0;
}

对与这个问题还有个更好的解法，可以在O（N）的时间内搞定，用的随机算法，具体的在以后贴出来的，that's all~~~

WsHunTer

关注

0
点赞
踩
11

收藏

觉得还不错? 一键收藏
3
评论
寻找最近点对问题（FCPP find the closest pair of point）

寻找最近点对问题（FCPP）问题描述：给定平面上n个点，找出距离最近的两个点。思考过程： 1）对于这种问题，我们首先想到的求解方法就是求出所有点对的距离，并找出最近的那个，当然这个是个显而易见的方法，具体过程大体可以警醒如下描述。 ① 定义变量 N (点的数目)，X[N] (点的x坐标值)，Y[N]（点的y坐标值），close
复制链接

扫一扫