距离最近的点对

韩冰

<>给定n 个点（xi，yi）（1≤i≤n），要求找出其中距离最近的两个点。
! d9 H$ Y8 p9 S) B6 C1 \1 M+ u
例14-7 假设在一片金属上钻n 个大小一样的洞，如果洞太近，金属可能会断。若知道任意两个洞的最小距离，可估计金属断裂的概率。这种最小距离问题实际上也就是距离最近的点对问题。

通过检查所有的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进行递归求解。


if (n较小) {用直接法寻找最近点对

R e t u r n ; }
; X: f8 g; S! S+ _7 g. |
6 U% o+ p: E# Y  r// n较大

将点集分成大致相等的两个部分A和B

确定A和B中的最近点对
' W% p- ?/ n- J* A. W4 f
确定一点在A中、另一点在B中的最近点对
# w0 ]1 f' ~' Z& [% o) _
+ Z9 A$ i$ x. H8 `6 v( X3 n从上面得到的三对点中，找出距离最小的一对点

图14-13 寻找最近的点对


为了确定第三种情况下的最近点对，需要采用一种不同的方法。这种方法取决于点集是如何被划分成A、B的。一个合理的划分方法是从xi（中间值）处划一条垂线，线左边的点属于A，线右边的点属于B。位于垂线上的点可在A和B之间分配，以便满足A、B的大小。
6 A+ I$ c  q0 g. F' y8 B
* q6 M+ E* ^" K. W3 H" z例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。
3 M3 p; Y5 x  h6 }' L
3 ]9 s# {, U$ j0 r设是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）。

: r  Q" S! Q- [5 z1 D. H例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) 是最近的点对。

为了确定一个距离更小的第三类点，RA 中的每个点最多只需和RB 中的6个点比较，如图1 4 - 1 6所示。
% V( W/ q  }2 i7 w1 e3 K( G
1. 选择数据结构
/ }2 N5 n% a! U* f1 l; ]- v
为了实现图1 4 - 1 3的分而治之算法，需要确定什么是“小问题”以及如何表示点。由于集合中少于两点时不存在最近点对，因此必须保证分解过程不会产生少于两点的点集。如果将少于四点的点集做为“小问题”，就可以避免产生少于两点的点集。
: r; ]% U4 W0 X; c/ c( Q, h- Y
9 g" P9 A- k0 m! m每个点可有三个参数：标号， x 坐标，y 坐标。假设标号为整数，每个点可用P o i n t l类（见程序1 4 - 8）来表示。为了便于按x 坐标对各个点排序，可重载操作符＜=。归并排序程序如1 4 -3所示。

' J) e: h/ _. q  W9 F6 S# T) I6 B3 v程序14-8 点类

: ?2 ?# M$ Z0 a7 @9 D6 N: f1 {class Point1 {
7 g3 b0 ?5 a. ~2 ?' g  F, J
- V  H: ^) e" X8 Bfriend float dist(const Point1&amp;, const Point1&amp;);
: h  s, s7 `$ V" P6 e  P
friend void close(Point1 *, Point2 *, Point2 *, int, int, Point1&amp;, Point1&amp;, float&amp;);
& U/ H! y( `- e4 Y7 t
friend bool closest(Point1 *, int, Point1&amp;, Point1&amp;,float&amp;);
( F! h  W% Y( y3 K3 w
7 N  a' Y, X( Y/ Nfriend void main();

6 L1 t1 p7 z# r/ V: P9 F; W4 Xp u b l i c :
- }+ ^3 F  X4 B
$ ~$ _4 L& O% rint operator&lt;=(Point1 a) const

{return (x &lt;= a.x);}
; j# a# N5 h% g, P
3 }" j( }6 s- q7 t8 bp r i v a t e :

int ID; // 点的编号
2 |& b# X- c3 U! I& @% ~2 F
: |% C: f" q" ?float x, y; // 点坐标

7 m: j8 q. s% m" Z2 o, ~- @} ;
- z. |! b6 Z7 a
! p( Z2 M  g3 Hclass Point2 {
0 s# K. Y3 T4 ]; K! {
friend float dist(const Point2&amp;, const Point2&amp;);

$ c9 u% H! d0 ], ?+ Yfriend void close(Point1 *, Point2 *, Point2 *, int, int, Point1&amp;, Point1&amp;, float&amp;);

friend bool closest(Point1 *, int, Point1&amp;, Point1&amp;, float&amp;);
! J* T! ?3 w3 e& L% _
9 }' L) |, ~& n$ g7 wfriend void main();

5 j0 ~7 C* S, L/ E& ]p u b l i c :
5 w/ j5 \$ N: X9 i0 F3 V! `3 e
int operator&lt;=(Point2 a) const
9 I2 W, R1 m2 L# C3 r9 ~
{return (y &lt;= a.y);}

p r i v a t e :
1 @8 N8 ~  \# C
5 O) W2 u9 T4 K5 {/ xint p; // 数组X中相同点的索引
3 ]% s0 i" D3 j7 U6 w
7 Q) I* }5 C1 |5 s4 u2 l1 F, u; V8 Jfloat x, y; // 点坐标

; q1 B- ]7 E5 f1 B3 [} ;
9 o' t: {9 r2 J6 p# m3 J/ m
, x: G* Y0 |1 H8 j所输入的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中的对应点。

1 D3 y, O7 [( C8 e. R/ d确定了必要的数据结构之后，再来看看所要产生的代码。首先定义一个模板函数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类的友元。

2 c, N$ \6 M6 J5 P: d程序14-9 计算两点距离
: u; f- }6 w$ @; a
) |- b+ c+ l# otemplate<CLASS T>

inline float dist(const T&amp; u, const T&amp; v)

{ / /计算点u 和v之间的距离

# z" a! }" _* \& Y; Q# ~$ w7 Y. n: [float dx = u.x-v. x ;

float dy = u.y-v. y ;
7 N; b% i7 C$ c) C6 q0 \
return sqrt(dx * dx + dy * dy);

% J/ e, @1 [! T}

  u2 p2 }, K9 }$ O/ ?* O- c$ U2 q如果点的数目少于两个，则函数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 )，该函数实际求解最近点对。
( @$ v' @4 C! n8 M9 X2 r. z$ U
5 p1 a$ w$ A4 e2 w6 {" k程序14-10 预处理及调用c l o s e
2 [) R: p: t: w! L- e' z) V
  l8 W" W& g' [" ?0 q$ n8 g& Tbool closest(Point1 X[], int n, Point1&amp; a, Point1&amp; b, float&amp; d)

( M$ [5 @( S+ l6 H{// 在n &gt;= 2 个点中寻找最近点对

// 如果少于2个点，则返回f a l s e

% |' w) @* t! G, N# {8 J' _// 否则，在a 和b中返回距离最近的两个点
5 I7 A. {2 O2 U- q7 T0 u
& y9 u0 G7 t. k2 H9 h. X5 N' Jif (n &lt; 2) return false;

4 s/ B( p6 f2 N// 按x坐标排序
4 B- X: r3 ?/ |) J4 k) n
; r& _  P+ j. R+ CM e r g e S o r t ( X , n ) ;

# s$ @( C! y6 J7 A// 创建一个按y坐标排序的点数组

1 N8 G& j) J2 t8 x' ~Point2 *Y = new Point2 [n];

. p7 Z7 l6 n5 }- H$ l! e6 bfor (int i = 0; i &lt; n; i++) {

// 将点i 从X 复制到Y
- g1 `8 P9 j; x/ A' `, o3 l- e
# ~, \) U/ D1 KY.p = i;
9 {2 `; L+ W6 X1 }7 L; h' C
Y.x = X.x;
! n$ K+ O+ M, K2 r7 C: A5 i3 ^7 S
Y.y = X.y;
! k" F5 {( U8 |! G' r" L7 p* u
; s9 `; O' g7 I. F1 K( B$ c}
3 g- q2 Y! O  d& ^! m' j
) x( F- Y+ w) KM e r g e S o r t ( Y,n); // 按y坐标排序

; g2 I. g3 O. f) w8 x$ e// 创建临时数组

# Z! B/ @9 i1 \- uPoint2 *Z = new Point2 [n];
* X. m/ M4 n$ n. E6 S
' `7 ~7 g8 t' v! W6 ]+ ]// 寻找最近点对

) g" _  N9 A7 ^1 f, F+ j9 C3 A& Oc l o s e ( X , Y, Z , 0 , n - 1 , a , b , d ) ;
/ e1 j$ V6 J- B( [2 {. ?
" U9 t, \1 s  Y9 C" o; T// 删除数组并返回
9 i3 K( M, I+ P( `" f
delete [] Y;

; V* g( q4 J$ Ndelete [] Z;

6 B+ k8 {6 P1 b# }6 U& Z$ k( Lreturn true;
9 _& r6 _; N3 m* i
}

3 @* C) T) T0 [程序1 4 - 11 计算最近点对
: h& P9 n$ u, ?2 N7 x6 _: F
( O. S/ B9 w( ?/ S6 avoid close(Point1 X[], Point2 Y[], Point2 Z[], int l, int r, Point1&amp; a, Point1&amp; b, float&amp; d)
- ~' x+ ^/ u  f& m! \; e- R
{//X[l:r] 按x坐标排序

1 {# u/ }) z6 q2 [5 \: Y//Y[l:r] 按y坐标排序

if (r-l == 1) {// 两个点
- e7 ~9 [3 r2 b* }+ ]' T
8 V- i1 G# W5 J( m' u$ K9 e% a: aa = X[l];

b = X[r];
8 I3 N! u! N. I$ l4 i3 I
: u$ Q  b# K; f# @d = dist(X[l], X[r]);
  f( H- r  m$ W: q1 ~, s$ f8 j
& T8 F( z6 M5 m5 I1 [/ [r e t u r n ; }
  R6 |. \: {2 p2 ^
if (r-l == 2) {// 三个点

* }; V: d' |4 y& W' c' Y' i; L/ z& R// 计算所有点对之间的距离
7 E5 s* i4 l7 E" ~6 x& e- A  W
7 J/ x  c4 v# x9 @3 u& u: Lfloat d1 = dist(X[l], X[l+1]);
! _9 n! j6 d' B: `& i
5 Y. A2 F* {7 g  R  m; o+ L% Z: kfloat d2 = dist(X[l+1], X[r]);
& s! Y* @, \7 P6 P7 ]
float d3 = dist(X[l], X[r]);
& D3 M8 r% c2 N+ [5 w( m
# g3 b3 l) ^- i% W// 寻找最近点对
' z; l; I; l. A( J; @: L
/ u6 T8 Q8 V, z: X/ L2 Dif (d1 &lt;= d2 &amp;&amp; d1 &lt;= d3) {

8 y0 f4 `# I; g: N, H/ Xa = X[l];

b = X[l+1];

/ r6 c9 @) v2 q6 R$ l7 X: R7 Yd = d1;

r e t u r n ; }

if (d2 &lt;= d3) {a = X[l+1];

/ R- a1 Z1 Q/ a, h& f: |# F1 V9 mb = X[r];
) ?% w' P5 e" p& W: g5 n
d = d2;}

5 @/ ~% @9 X/ u: |else {a = X[l];

9 j! z1 \1 K4 {0 k* v  v& [- ob = X[r];
9 }# X7 E- `5 j8 ^1 j! c
* x. S- p, N9 j+ s' g- md = d3;}
( u" p1 s: q& K* o6 b& K
6 q! ]( d/ O0 G# Ir e t u r n ; }

/ /多于三个点，划分为两部分
" U% k/ Q) S. g
3 T) M  j, A' k. @- Yint m = (l+r)/2; // X[l:m] 在A中，余下的在B中

' a5 `# g, v7 z* n1 X// 在Z[l:m] 和Z [ m + 1 : r ]中创建按y排序的表

int f = l, // Z[l:m]的游标

- Q! ^, D* [4 Q; L! y. j) Mg = m+1; // Z[m+1:r]的游标

' t6 `5 L- f0 p, S: }! M" l; L1 R' Ffor (int i = l; i &lt;= r; i++)
; t1 c* V1 L' R2 s
if (Y.p &gt; m) Z[g++] = Y;

else Z[f++] = Y;

// 对以上两个部分进行求解
1 w) t9 [, K" [5 c( K
+ e. I% Y$ Y! S" yc l o s e ( X , Z , Y, l , m , a , b , d ) ;
9 Y- G. {- h& ?' q" Y' w
0 r4 L! I% U# V6 [& h. Tfloat dr;
9 z1 v  v  O, K" l* W1 \8 n' b
Point1 ar, br;

c l o s e ( X , Z , Y, m + 1 , r, a r, b r, d r ) ;
0 x: K! A& W# `; _' G
// (a,b) 是两者中较近的点对
2 F* v" z! o) j' h
if (dr &lt; d) {a = ar;

+ I, H' U7 x; ]* {9 ib = br;
' y7 E  f* {& }. `9 `
& ]# k. Y% N% B+ w# Dd = dr;}

" T9 h1 @  s- c3 {6 nM e r g e ( Z , Y,l,m,r);// 重构Y

$ u, C  f" O5 ]/ /距离小于d的点放入Z
  o# p) T4 a$ N# m  i
% P: a; v# D) v. P' fint k = l; // Z的游标

for (i = l; i &lt;= r; i++)
0 k$ T/ ^! F4 l$ M, [6 c  e! m
+ Q# J) n$ Z( `) \* Xif (fabs(Y[m].x - Y.x) &lt; d) Z[k++] = Y;

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

for (i = l; i &lt; k; i++){

for (int j = i+1; j &lt; k &amp;&amp; Z[j].y - Z.y &lt; d;
/ V$ C" a  `) N1 m5 a% }
j + + ) {

3 Q. h5 G: E  ]/ b' W7 pfloat dp = dist(Z, Z[j]);
- Q0 Z3 i8 v/ n
& @( `" V0 d5 r% _% H: `if (dp &lt; d) {// 较近的点对

d = dp;

! @3 M( F/ A. ~% C; D9 Ka = X[Z.p];

b = X[Z[j].p];}

}

}
8 ]- ]* B) X8 g! |# V4 S
  p8 f! C. a  J  x: a6 m" s, ]}

) P. J- X/ k. ]) _$ z3 E函数c l o s e（见程序1 4 - 11）用来确定X[1:r] 中的最近点对。假定这些点按x 坐标排序。在Y [ 1 : r ]中对这些点按y 坐标排序。Z[ 1 : r ]用来存放中间结果。找到最近点对以后，将在a, b中返回最近点对，在d 中返回距离，数组Y被恢复为输入状态。函数并未修改数组X。

8 G3 j6 J( g9 u9 X6 H$ e首先考察“小问题”，即少于四个点的点集。因为分割过程不会产生少于两点的数组，因此只需要处理两点和三点的情形。对于这两种情形，可以尝试所有的可能性。当点数超过三个时，通过计算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 ]。

1 @: H7 E8 U, V% ~4 z为实现图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
$ A, p7 N. W0 z6 i  k
; d2 m. G' f, X: T2 I, z# s0 C3 }还是在RB中）与Z[j] 配对来实现，其中i＜j 且Z [ j ] . y - Z [ i ] . y＜ 。对每一个Z [ i ]，在2 × 区域内所检查的点如图1 4 - 1 7所示。由于在每个2 × 子区域内的点至少相距。因此每一个子区域中的点数不会超过四个，所以与Z [ i ]配对的点Z [ j ]最多有七个。
. F& k) u. h$ Q
/ q9 h; g& |' m! K- \1 P2 a" g2. 复杂性分析
8 P: y; q& E8 ?# S; n
令t (n) 代表处理n 个点时，函数close 所需要的时间。当n＜4时，t (n) 等于某个常数d。当n≥4时，需花费(n) 时间来完成以下工作：将点集划分为两个部分，两次递归调用后重构Y，淘汰距分割线很远的点，寻找更好的第三类点对。两次递归调用需分别耗时t (「n /2ù」和t (?n /2?).

4 |8 g' ~) f1 x! V: d  J这个递归式与归并排序的递归式完全一样，其结果为t (n) = (nl o gn)。另外，函数c l o s e s t还需耗时(nl o gn)来完成如下额外工作：对X进行排序，创建Y和Z，对Y进行排序。因此分而治之最近点对求解算法的时间复杂性为(nl o gn)。</P>