距离最近的点对

韩冰

<>给定n 个点（xi，yi）（1≤i≤n），要求找出其中距离最近的两个点。
1 s3 h. j% ]+ A9 G& z' l
6 ^; D, j% R; ?( Y例14-7 假设在一片金属上钻n 个大小一样的洞，如果洞太近，金属可能会断。若知道任意两个洞的最小距离，可估计金属断裂的概率。这种最小距离问题实际上也就是距离最近的点对问题。
, p2 ~8 T. y" p
% g& R! y# ~6 X; J9 }通过检查所有的n(n- 1 ) / 2对点，并计算每一对点的距离，可以找出距离最近的一对点。这种方法所需要的时间为(n2 )。我们称这种方法为直接方法。图1 4 - 1 3中给出了分而治之求解算法的伪代码。该算法对于小的问题采用直接方法求解，而对于大的问题则首先把它划分为两个较小的问题，其中一个问题（称为A）的大小为「n /2ù，另一个问题（称为B）的大小为「n /2ù。初始时，最近的点对可能属于如下三种情形之一： 1) 两点都在A中（即最近的点对落在A中）；2) 两点都在B中；3) 一点在A，一点在B。假定根据这三种情况来确定最近点对，则最近点对是所有三种情况中距离最小的一对点。在第一种情况下可对A进行递归求解，而在第二种情况下可对B进行递归求解。
1 V4 P- {+ B9 u
/ i6 k$ n. @* F* [8 L
; T0 O/ Y$ I* X% A8 Yif (n较小) {用直接法寻找最近点对
2 Q6 f, T( }  Q/ {8 `0 z- h
7 k) D. _% T! }4 rR e t u r n ; }

- V/ I: t- O/ V  ]$ p// n较大

1 T: s3 W6 E- v! w* I将点集分成大致相等的两个部分A和B
3 v  `) X0 w1 X) c) N; p6 a
  u! q- D* z8 O& A确定A和B中的最近点对

4 B, [8 p5 o! @2 ^" P确定一点在A中、另一点在B中的最近点对

从上面得到的三对点中，找出距离最小的一对点
% C4 `& r5 l) _0 ?+ [4 I2 O
% o$ i& b8 ?# ]% n  d5 T+ e* @0 c图14-13 寻找最近的点对

0 R: q2 j, b2 e9 k+ y; K
6 E. l9 |7 v2 O- j& e$ u7 ?. x为了确定第三种情况下的最近点对，需要采用一种不同的方法。这种方法取决于点集是如何被划分成A、B的。一个合理的划分方法是从xi（中间值）处划一条垂线，线左边的点属于A，线右边的点属于B。位于垂线上的点可在A和B之间分配，以便满足A、B的大小。

例2-8 考察图14-14a 中从a到n的1 4个点。这些点标绘在图14-14b 中。中点xi = 1，垂线x = 1如图14-14b 中的虚线所示。虚线左边的点(如b, c, h, n, i)属于A，右边的点(如a, e, f, j, k, l) 属于B。d, g, m 落在垂线上，可将其中两个加入A, 另一个加入B，以便A、B中包含相同的点数。假设d ,m加入A，g加入B。

+ r4 v# U( o6 U8 S% K设是i 的最近点对和B的最近点对中距离较小的一对点。若第三种情况下的最近点对比小。则每一个点距垂线的距离必小于，这样，就可以淘汰那些距垂线距离≥ 的点。图1 4 - 1 5中的虚线是分割线。阴影部分以分割线为中线，宽为2 。边界线及其以外的点均被淘汰掉，只有阴影中的点被保留下来，以便确定是否存在第三类点对（对应于第三种情况）其距离小于。用RA、RB 分别表示A和B中剩下的点。如果存在点对(p,q)，p?A, q?B且p, q 的距离小于，则p?RA，q?RB。可以通过每次检查RA 中一个点来寻找这样的点对。假设考察RA 中的p 点，p的y 坐标为p.y，那么只需检查RB 中满足p.y- ＜q.y＜p.y+ 的q 点，看是否存在与p 间距小于的点。在图14-16a 中给出了包含这种q 点的RB 的范围。因此，只需将RB 中位于×2 阴影内的点逐个与p 配对，以判断p 是否是距离小于的第三类点。这个×2 区域被称为是p 的比较区（comparing region）。

, W% Y" h) p; ~+ q8 z5 E例2-9 考察例2 - 8中的1 4个点。A中的最近点对为(b,h)，其距离约为0 . 3 1 6。B中最近点对为(f, j)，其距离为0 . 3，因此= 0 . 3。当考察是否存在第三类点时，除d, g, i, l, m 以外的点均被淘汰，因为它们距分割线x= 1的距离≥ 。RA ={d, i, m}，RB= {g, l}，由于d 和m 的比较区中没有点，只需考察i即可。i 的比较区中仅含点l。计算i 和l的距离，发现它小于，因此(i, l) 是最近的点对。
1 o' Z' h) j1 N2 {) \! f% u
0 n& s$ }1 r! Y6 ^  Q" S为了确定一个距离更小的第三类点，RA 中的每个点最多只需和RB 中的6个点比较，如图1 4 - 1 6所示。
. K6 H4 k2 W6 p, G! g5 T8 ]8 _5 H
3 N# C9 y7 x9 y. e2 c1. 选择数据结构
( [( Y& v/ y- @( Z( [1 |
, q/ o2 G* g( }- H为了实现图1 4 - 1 3的分而治之算法，需要确定什么是“小问题”以及如何表示点。由于集合中少于两点时不存在最近点对，因此必须保证分解过程不会产生少于两点的点集。如果将少于四点的点集做为“小问题”，就可以避免产生少于两点的点集。
$ ^, z1 E' \' e: Z: ~0 g
2 E7 y0 K. U) x3 R每个点可有三个参数：标号， x 坐标，y 坐标。假设标号为整数，每个点可用P o i n t l类（见程序1 4 - 8）来表示。为了便于按x 坐标对各个点排序，可重载操作符＜=。归并排序程序如1 4 -3所示。

, U, d" ^  I' `4 ~7 k程序14-8 点类

class Point1 {
4 q) x* F3 X' ?" Y, J# K
friend float dist(const Point1&amp;, const Point1&amp;);
' g. B4 p+ G8 p! b; @( {
/ u" M3 P6 G: e+ u$ ^6 f  Cfriend void close(Point1 *, Point2 *, Point2 *, int, int, Point1&amp;, Point1&amp;, float&amp;);

friend bool closest(Point1 *, int, Point1&amp;, Point1&amp;,float&amp;);

friend void main();

- @" i) E4 H. l6 Pp u b l i c :

+ b, _7 l' \" V9 J/ |1 V2 a- p" Xint operator&lt;=(Point1 a) const
, K! V9 l3 r7 g
{return (x &lt;= a.x);}

& T% Q2 r/ S! X, h  @9 [p r i v a t e :
& m' N/ t* v+ |; Y  V+ n; @
int ID; // 点的编号

float x, y; // 点坐标

) H& C! ?, N' S, V! N3 y; v; e} ;

class Point2 {
- K1 Q9 _3 P+ [' s
friend float dist(const Point2&amp;, const Point2&amp;);

friend void close(Point1 *, Point2 *, Point2 *, int, int, Point1&amp;, Point1&amp;, float&amp;);

) h$ J# `' |! E: Y1 ffriend bool closest(Point1 *, int, Point1&amp;, Point1&amp;, float&amp;);
( C# j, k  r9 c6 P# u
friend void main();
! G2 O$ _3 J0 @$ ?- I
p u b l i c :
( U" P7 D6 a: ]
int operator&lt;=(Point2 a) const
$ D  X, ]! j& D4 D1 y
{return (y &lt;= a.y);}

& |# d7 m; V+ P# C* [. D( Tp r i v a t e :

; S& w4 N# n# l: v( I! p; cint p; // 数组X中相同点的索引
! }: s2 j+ d! p" r, `8 ~0 h) R
float x, y; // 点坐标

( ?7 x& T1 {" S) c  ^6 |} ;
5 j6 M) t( S( \) r, r+ o8 q8 L
所输入的n 个点可以用数组X来表示。假设X中的点已按照x 坐标排序，在分割过程中如果当前考察的点是X [l :r]，那么首先计算m= (l+r) / 2，X[ l:m]中的点属于A，剩下的点属于B。计算出A和B中的最近点对之后，还需要计算RA 和RB，然后确定是否存在更近的点对，其中一点属于RA，另一点属于RB。如果点已按y 坐标排序，那么可以用一种很简单的方式来测试图1 4 - 1 6。按y 坐标排序的点保存在另一个使用类P o i n t 2 (见程序14-8) 的数组中。注意到在P o i n t 2类中，为了便于y 坐标排序，已重载了操作符＜＝。成员p 用于指向X中的对应点。

确定了必要的数据结构之后，再来看看所要产生的代码。首先定义一个模板函数d i s t (见程序1 4 - 9 )来计算点a, b 之间的距离。T可能是P o i n t 1或P o i n t 2，因此d i s t必须是P o i n t 1和P o i n t 2类的友元。
4 d( i6 H. ~4 x) b& L- |  N
程序14-9 计算两点距离

& c) u6 y9 J6 L+ Ktemplate<CLASS T>

inline float dist(const T&amp; u, const T&amp; v)
+ Z8 `1 C$ L7 t9 B
{ / /计算点u 和v之间的距离
, _, y" e' j2 g+ u
% d! S3 ^! Y  S( Mfloat dx = u.x-v. x ;
! b0 W+ Y, C1 d& Q9 Y2 p
& U8 W/ ]6 h- _+ n+ C3 Jfloat dy = u.y-v. y ;

; l6 }! W3 D3 L( ~$ @0 k6 Zreturn sqrt(dx * dx + dy * dy);

5 c! L2 t6 |2 B$ G2 }}

9 O4 n6 h  i+ |- j9 g. J7 l如果点的数目少于两个，则函数c l o s e s t (见程序1 4 - 1 0 )返回f a l s e，如果成功时函数返回t r u e。当函数成功时，在参数a 和b 中返回距离最近的两个点，在参数d 中返回距离。代码首先验证至少存在两点，然后使用M e rg e S o r t函数(见程序14-3) 按x 坐标对X中的点排序。接下来把这些点复制到数组Y中并按y 坐标进行排序。排序完成时，对任一个i，有Y [i ] . y≤Y [i+ 1 ] . y，并且Y [i ] .p给出了点i 在X中的位置。上述准备工作做完以后，调用函数close (见程序1 4 - 11 )，该函数实际求解最近点对。

4 W6 }, m3 R% E9 `' A  _/ e7 [程序14-10 预处理及调用c l o s e
; p/ P# H% e/ F" U# m
7 n# p/ ^0 ]* Z, D$ G- W, f8 V5 abool closest(Point1 X[], int n, Point1&amp; a, Point1&amp; b, float&amp; d)
2 e) _# v" U( \5 D
! {% Y0 p8 {- r, ~" D4 }+ P& ^8 i{// 在n &gt;= 2 个点中寻找最近点对
/ H1 y5 {5 h+ Q9 D! X
/ C* d% L6 s6 M8 K) t// 如果少于2个点，则返回f a l s e

// 否则，在a 和b中返回距离最近的两个点

if (n &lt; 2) return false;

% F  ^9 W8 Z  [3 s3 b// 按x坐标排序
! ^( m* @3 p# m& ^" v
" \$ S  P  }0 @, ^& U7 c0 ]+ [8 {M e r g e S o r t ( X , n ) ;
( W6 H5 T* K; @
1 b* H8 e. b. q" i: x) n, J// 创建一个按y坐标排序的点数组
# c: G$ F4 C: `' o/ W; `
) V4 D( C. N! p0 y+ }" cPoint2 *Y = new Point2 [n];

for (int i = 0; i &lt; n; i++) {

// 将点i 从X 复制到Y

Y.p = i;
' J! c9 F0 c/ e. ^' L) R1 U! g
1 R* D; q8 Q6 p9 B1 [7 WY.x = X.x;

Y.y = X.y;
. W. @" \4 A4 R: Q
}
8 V* |$ P2 z" Q! Z+ K( y: C; V
M e r g e S o r t ( Y,n); // 按y坐标排序
! f' r4 [, e  @; |* Z
! M+ t/ J* \8 t// 创建临时数组

" h# k- R! t& R# S  yPoint2 *Z = new Point2 [n];
6 k) J0 o- }! a; D: A4 ]( U" i4 [
// 寻找最近点对
8 u3 O: J6 Q) F+ S" v" V8 w7 c( K/ R$ O
2 o; K# w( O6 x' W" _$ rc l o s e ( X , Y, Z , 0 , n - 1 , a , b , d ) ;
, X: {' j, {7 V6 x! D
8 L) d7 x2 g$ D1 h// 删除数组并返回

& b+ s. t8 Z" N- W4 }4 Edelete [] Y;

. D) q2 K3 F8 X3 B, Mdelete [] Z;
5 w9 `. d0 {6 E* A$ q7 K
2 `  h- Z0 R: ^return true;

}

( Y3 n5 `1 R, ]5 D3 T9 [程序1 4 - 11 计算最近点对
+ ?* y. v5 O+ K/ W1 W
1 e  t* ]6 j9 Bvoid close(Point1 X[], Point2 Y[], Point2 Z[], int l, int r, Point1&amp; a, Point1&amp; b, float&amp; d)

{//X[l:r] 按x坐标排序
3 k7 {6 U6 h! x$ B! m
//Y[l:r] 按y坐标排序

if (r-l == 1) {// 两个点
4 _$ Q1 W7 n4 ?& x. w- w0 z& g  a5 e
, s  X3 b1 S! b  b9 A' a' Ua = X[l];

b = X[r];

7 s$ y# s' u: B4 R- m3 P# C( d* dd = dist(X[l], X[r]);
7 `+ ^* n0 h$ @4 n
# z. B9 m) `. {& |; Vr e t u r n ; }

if (r-l == 2) {// 三个点

+ E& R# ]: N' |$ M' @) w! D& E0 n$ I// 计算所有点对之间的距离
. `2 L3 ]1 ~# N( U/ F! N& o; ~
float d1 = dist(X[l], X[l+1]);
/ D; r$ J  s5 `/ W. Z/ z
' {2 y! s' h- B% }2 hfloat d2 = dist(X[l+1], X[r]);
( _( @: \7 ~; ?
float d3 = dist(X[l], X[r]);

// 寻找最近点对

5 b' o# g$ I, v- o; ^4 qif (d1 &lt;= d2 &amp;&amp; d1 &lt;= d3) {
9 D( ?! ^; i1 F( }2 G
, h) a% J3 ?- sa = X[l];
, o' u% g, L- Y7 d: G2 y8 M- L. R( k
3 \# i7 B0 z4 _& ?b = X[l+1];

d = d1;

r e t u r n ; }
: N. d0 ?: f$ C. f6 ~
if (d2 &lt;= d3) {a = X[l+1];
. d  y( Z7 J. Q; t8 n" P
b = X[r];
" J; O3 V; \. H, l: |! [& O, J
d = d2;}
4 e( y! |, c7 u( U( W6 r
; K, k! Z7 t' \9 Y; S' M- \' Z2 Xelse {a = X[l];
' J6 J  D1 ]% u7 {/ L
' b4 J- T7 K! Y  l# q! M1 {b = X[r];

d = d3;}
7 q; k1 q4 h7 o
( ^' j. S) d+ l- v8 d- xr e t u r n ; }

3 f* K# Q+ V: ?) K+ M" h/ /多于三个点，划分为两部分
% X$ W! ]' ~. k8 |% h! e8 ]4 S  h
' H" B; s. E6 C. V' Oint m = (l+r)/2; // X[l:m] 在A中，余下的在B中

// 在Z[l:m] 和Z [ m + 1 : r ]中创建按y排序的表

4 w, A- N8 _+ S/ Sint f = l, // Z[l:m]的游标
+ k/ s1 s% p0 a9 [
g = m+1; // Z[m+1:r]的游标
6 Y- y& `1 q+ o/ m7 b
) n; y! J) R  g9 F' y1 X  S1 ufor (int i = l; i &lt;= r; i++)
7 R8 }# m2 B' D5 o# K8 ~
! o# O3 n# ~% g1 wif (Y.p &gt; m) Z[g++] = Y;

else Z[f++] = Y;
6 Q8 e" o2 [4 U$ {. v# _
( J. K& r1 ~& S6 R+ a1 [// 对以上两个部分进行求解

; a: [- C8 v. dc l o s e ( X , Z , Y, l , m , a , b , d ) ;

float dr;

: P, f& Z: h9 u. `6 \6 LPoint1 ar, br;

4 V3 O8 T) M8 z  a. Sc l o s e ( X , Z , Y, m + 1 , r, a r, b r, d r ) ;
: q9 Q8 ~8 G8 F4 Y
6 m9 J8 G8 V/ g// (a,b) 是两者中较近的点对

) Z2 F1 |0 W" @* u" Wif (dr &lt; d) {a = ar;
0 h. F5 Y3 }1 z( B1 k
. I- P; b8 i: e5 ]1 |7 A6 z+ e# x7 Bb = br;
+ a5 g! c# r! w' l8 s+ T5 M
0 e' D, \) d5 H* b8 bd = dr;}
' F4 X. x) |0 y5 A# T+ k% E+ K3 ?/ u
M e r g e ( Z , Y,l,m,r);// 重构Y
& b: L: R( H# s* `3 e* q
$ X* n& I5 v  o( D. d) W/ /距离小于d的点放入Z
* k( n& C& f* k$ f" B8 A/ |
" s6 I- ]* F7 Mint k = l; // Z的游标

0 \6 f6 T& z; s. b; e; sfor (i = l; i &lt;= r; i++)
5 W+ I9 z& Y. s$ K% V
) h- i7 B' p& _0 n* f; ~if (fabs(Y[m].x - Y.x) &lt; d) Z[k++] = Y;

// 通过检查Z [ l : k - 1 ]中的所有点对，寻找较近的点对

# C9 z$ B5 d( X5 q% q- wfor (i = l; i &lt; k; i++){
" d) o0 o) X5 B" l
+ u; V/ `& V0 z+ M( ]0 zfor (int j = i+1; j &lt; k &amp;&amp; Z[j].y - Z.y &lt; d;
/ Q' ~4 c( L0 E9 C; }9 T+ n
' O  w' W$ E4 X+ Sj + + ) {

float dp = dist(Z, Z[j]);

if (dp &lt; d) {// 较近的点对
4 D) E6 I( g% ^3 _8 r/ I
d = dp;

  W/ F7 W7 E" B- _a = X[Z.p];

" U+ j! Z3 e2 b# e( H4 r' e) Yb = X[Z[j].p];}

' ?( X2 K5 I/ Q: q6 o}

}
! k( D0 H+ U0 l9 P; S: a4 @
}
2 y8 H# @0 ]2 V+ w; s8 D2 K
+ |" c, B, {( l/ Z5 Q+ ?4 Z% ?/ {函数c l o s e（见程序1 4 - 11）用来确定X[1:r] 中的最近点对。假定这些点按x 坐标排序。在Y [ 1 : r ]中对这些点按y 坐标排序。Z[ 1 : r ]用来存放中间结果。找到最近点对以后，将在a, b中返回最近点对，在d 中返回距离，数组Y被恢复为输入状态。函数并未修改数组X。

, y2 [7 i7 j% i/ }) L; S2 @5 Q2 S  j首先考察“小问题”，即少于四个点的点集。因为分割过程不会产生少于两点的数组，因此只需要处理两点和三点的情形。对于这两种情形，可以尝试所有的可能性。当点数超过三个时，通过计算m = ( 1 + r ) / 2把点集分为两组A和B，X [ 1 : m ]属于A，X [ m + 1 : r ]属于B。通过从左至右扫描Y中的点以及确定哪些点属于A，哪些点属于B，可以创建分别与A组和B组对应的，按y 坐标排序的Z [ 1 : m ]和Z [ m + 1 : r ]。此时Y和Z的角色互相交换，依次执行两个递归调用来获取A和B中的最近点对。在两次递归调用返回后，必须保证Z不发生改变，但对Y则无此要求。不过，仅Y [ l : r ]可能会发生改变。通过合并操作（见程序1 4 - 5）可以以Z [ 1 : r ]重构Y [ 1 : r ]。
% y: U0 n! W" s
为实现图1 4 - 1 6的策略，首先扫描Y [ 1 : r ]，并收集距分割线小于的点，将这些点存放在Z [ 1 : k - 1 ]中。可按如下两种方式来把RA中点p 与p 的比较区内的所有点进行配对：1) 与RB 中y 坐标≥p.y 的点配对；2) 与y 坐标≤p.y 的点配对。这可以通过将每个点Z [ i ]（1≤i &lt; k，不管该点是在RA
+ t, j8 |- e3 `/ t/ I
/ r( w0 W6 B" f1 b  p5 h( i还是在RB中）与Z[j] 配对来实现，其中i＜j 且Z [ j ] . y - Z [ i ] . y＜ 。对每一个Z [ i ]，在2 × 区域内所检查的点如图1 4 - 1 7所示。由于在每个2 × 子区域内的点至少相距。因此每一个子区域中的点数不会超过四个，所以与Z [ i ]配对的点Z [ j ]最多有七个。
0 B/ z$ y# \" s5 W
2. 复杂性分析

' a1 f% s8 ]6 v* V" j  G令t (n) 代表处理n 个点时，函数close 所需要的时间。当n＜4时，t (n) 等于某个常数d。当n≥4时，需花费(n) 时间来完成以下工作：将点集划分为两个部分，两次递归调用后重构Y，淘汰距分割线很远的点，寻找更好的第三类点对。两次递归调用需分别耗时t (「n /2ù」和t (?n /2?).
8 y9 O6 b  b  c2 G; j3 a
这个递归式与归并排序的递归式完全一样，其结果为t (n) = (nl o gn)。另外，函数c l o s e s t还需耗时(nl o gn)来完成如下额外工作：对X进行排序，创建Y和Z，对Y进行排序。因此分而治之最近点对求解算法的时间复杂性为(nl o gn)。</P>