距离最近的点对

韩冰

<>给定n 个点（xi，yi）（1≤i≤n），要求找出其中距离最近的两个点。

例14-7 假设在一片金属上钻n 个大小一样的洞，如果洞太近，金属可能会断。若知道任意两个洞的最小距离，可估计金属断裂的概率。这种最小距离问题实际上也就是距离最近的点对问题。
2 t5 O9 j# P; r
2 K8 A. ~% \$ {) n8 F3 G! |  R3 d通过检查所有的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进行递归求解。
9 D" f% }, ?9 ~1 l% g

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

R e t u r n ; }
/ b4 v6 r' }: b' F6 A$ f7 w
// n较大
5 T5 X8 W7 W0 U, b$ a
2 R3 b/ m5 P+ u) u0 g将点集分成大致相等的两个部分A和B

确定A和B中的最近点对
& h& k' y1 H7 \5 ^" G# ~' y) z
确定一点在A中、另一点在B中的最近点对

, |4 h9 i, `. K" N从上面得到的三对点中，找出距离最小的一对点

图14-13 寻找最近的点对
4 D4 {, w( g3 P$ J- b3 `
0 {# g4 T1 A0 G9 b! C
为了确定第三种情况下的最近点对，需要采用一种不同的方法。这种方法取决于点集是如何被划分成A、B的。一个合理的划分方法是从xi（中间值）处划一条垂线，线左边的点属于A，线右边的点属于B。位于垂线上的点可在A和B之间分配，以便满足A、B的大小。
6 ~5 @8 P1 e! J) m
例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。

9 ~7 \9 [; t9 \设是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）。
" P) x# g$ U) {* f/ q- ]3 ?, W
# _( O& p9 |* N5 F0 d; u例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) 是最近的点对。
3 ]' Q# r8 U- ]
/ x' n9 E- Y( M5 Q为了确定一个距离更小的第三类点，RA 中的每个点最多只需和RB 中的6个点比较，如图1 4 - 1 6所示。

1. 选择数据结构
- ~& E4 @' V3 _+ J4 I
0 O6 a4 ?0 {3 }6 g为了实现图1 4 - 1 3的分而治之算法，需要确定什么是“小问题”以及如何表示点。由于集合中少于两点时不存在最近点对，因此必须保证分解过程不会产生少于两点的点集。如果将少于四点的点集做为“小问题”，就可以避免产生少于两点的点集。

每个点可有三个参数：标号， x 坐标，y 坐标。假设标号为整数，每个点可用P o i n t l类（见程序1 4 - 8）来表示。为了便于按x 坐标对各个点排序，可重载操作符＜=。归并排序程序如1 4 -3所示。
0 T+ J. O$ P+ @. @( _
程序14-8 点类

class Point1 {

: O+ j- Z8 y5 l& Wfriend float dist(const Point1&amp;, const Point1&amp;);
- x- T; F; S4 Q% v. Q9 K- P! _, W$ ~
friend 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();

( Y- D: a. M" x0 t* {p u b l i c :

& f" k; |+ I/ C# p0 hint operator&lt;=(Point1 a) const

{return (x &lt;= a.x);}

* l, Q! y$ M" l; ]) V) _p r i v a t e :
+ w( W9 z% ]9 Z, A$ b5 \7 {; i
int ID; // 点的编号
- U9 S( J/ ~2 R, Z$ U7 k
' T& f+ m- h& Rfloat x, y; // 点坐标
. V: E' p* Z1 ~+ r
5 w5 K! w. {8 j; X* i' @} ;

' w# F" x1 O! D5 h' sclass Point2 {
: q0 t9 q% i2 u
% ~, C! |% D" \: x# c2 e* s2 Q" f$ \- @friend float dist(const Point2&amp;, const Point2&amp;);
* ^) Q5 e  ^3 _/ E6 v
friend 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();

  M8 E/ P* B9 {" g5 N2 l$ `  }p u b l i c :

int operator&lt;=(Point2 a) const
+ b: t/ P( B7 J" S, p2 o+ C
4 A9 f% X; d# n  r2 q, Z{return (y &lt;= a.y);}

9 Y- v3 q' U3 O' Z$ lp r i v a t e :
* [7 ]  L3 w  v( I  t$ j$ P$ f
  A  r5 ], T4 Y( aint p; // 数组X中相同点的索引
* k- j, J0 W1 j4 j
float x, y; // 点坐标

} ;

: f5 U; f: h5 b( {4 @+ b所输入的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类的友元。

程序14-9 计算两点距离
  m) ~. l3 i3 d
3 l  t  @, e* g' M8 {& Q  L, ^' Gtemplate<CLASS T>
4 v% k+ [/ t" t6 c) f7 ~0 s" u
* h4 d, Y" E9 Y/ ?, E5 [. winline float dist(const T&amp; u, const T&amp; v)

: r0 I% p2 r# w% I2 e+ h3 l{ / /计算点u 和v之间的距离
) ?$ i; h; E2 W+ S
- k! n6 Y, D2 p) x3 Rfloat dx = u.x-v. x ;
, l- C4 K$ K/ V3 v
( f3 J& x* A( g* pfloat dy = u.y-v. y ;

* G% d( W4 h. A3 p  lreturn sqrt(dx * dx + dy * dy);

* z  P+ t1 `! d}
: a) m4 Z/ t* V  m9 z& g
7 D. e" `- T* w4 L  d7 z如果点的数目少于两个，则函数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 )，该函数实际求解最近点对。

: {) F' {, |2 Q程序14-10 预处理及调用c l o s e
- ^5 z0 F1 {. w( F4 e
- G& D) [2 `& Q9 pbool closest(Point1 X[], int n, Point1&amp; a, Point1&amp; b, float&amp; d)

{// 在n &gt;= 2 个点中寻找最近点对

! L5 z( e% b. Z7 P" v: ^// 如果少于2个点，则返回f a l s e
7 A& J7 Z& R. ~" R6 A
// 否则，在a 和b中返回距离最近的两个点
& ]% \% v* M' ^$ Y6 ~% ]
& v/ s1 y, h/ |# ~6 Hif (n &lt; 2) return false;
# S+ x' s% N" k) f5 c
# O1 f; D# ?9 y, _// 按x坐标排序
# T% o; _! P% F
! c" \6 F& `8 rM e r g e S o r t ( X , n ) ;
. `  v1 P: F* E2 \6 w
// 创建一个按y坐标排序的点数组

Point2 *Y = new Point2 [n];
/ p! W, k9 {9 J' Q+ Q1 C+ h# F
& T% v$ x" b1 e1 P8 Zfor (int i = 0; i &lt; n; i++) {
/ R0 C8 j. c) t0 ~7 n7 E0 {
// 将点i 从X 复制到Y
# V) L: D& m* P2 l) G: M0 `
Y.p = i;

9 c* k, @2 u; a1 A3 @: ZY.x = X.x;
' P6 r6 q! S: Z  Z8 u
Y.y = X.y;
/ p8 R5 g* J8 X  S$ ~: f6 L; D  A
}
& z/ @% u) o+ U- X9 u
8 ?; _! a$ i8 O, b/ V) S0 FM e r g e S o r t ( Y,n); // 按y坐标排序

// 创建临时数组
3 b: S8 o4 M2 S6 s4 l, |
Point2 *Z = new Point2 [n];

// 寻找最近点对
' f* A! m% J' Z  w/ f! q
c l o s e ( X , Y, Z , 0 , n - 1 , a , b , d ) ;
0 b/ g. |# D% ^' ^
# C6 q* H6 e$ R  l! \// 删除数组并返回

delete [] Y;
) E: h' j3 K+ K7 P& I7 G# O# `! ~
delete [] Z;
3 G7 Q6 t0 f7 E6 ]* \( X" U' g& @
return true;

, E9 _; h5 l! s( @" x9 y% l}
" p+ P0 \* D, v" a' z' N
程序1 4 - 11 计算最近点对
! g: r( P' z1 n
3 c% z- T. G" T& j* H6 ~. Ivoid close(Point1 X[], Point2 Y[], Point2 Z[], int l, int r, Point1&amp; a, Point1&amp; b, float&amp; d)
' N2 l0 V4 L+ w, g2 r# A
{//X[l:r] 按x坐标排序
% i6 l# q7 ?& y3 x) a
//Y[l:r] 按y坐标排序
4 `& ?) y: t) I6 I: D6 ^
if (r-l == 1) {// 两个点

4 N% ?/ t: J6 @9 |+ {; Y! |7 ia = X[l];
9 G. n6 V% R' a$ _9 P. \( y
b = X[r];
. _. X. Z9 c: Q" {0 M( a* I
d = dist(X[l], X[r]);
7 i3 a! D6 y2 J4 e0 E1 [. x7 t
4 Y9 \5 M, G) x  ^) cr e t u r n ; }

1 ~" p" Z/ D: D% Oif (r-l == 2) {// 三个点

/ C! Z6 {6 `2 A. y0 a5 C2 W8 w  C// 计算所有点对之间的距离
3 f+ A3 f' m7 ^
3 r  N/ A0 p' O1 i/ o5 a' B5 pfloat d1 = dist(X[l], X[l+1]);
- X' b7 C, c- _: G" j' F
9 r6 r8 C: H6 a" E4 S: pfloat d2 = dist(X[l+1], X[r]);
' M$ w. \0 N, l8 y% t
- L- t- ]. ~# Y- o& Kfloat d3 = dist(X[l], X[r]);

8 t7 v- S4 b4 i9 M9 N0 {) [; M+ W; S// 寻找最近点对

+ [2 Y6 @& f% v" r3 c' q, Nif (d1 &lt;= d2 &amp;&amp; d1 &lt;= d3) {

+ t9 i2 E3 d  \9 w# ]" Y; Sa = X[l];
* E- z: M, {  F! [/ p( P( A
b = X[l+1];
/ ^% S) u5 q( q3 T1 b3 Q
2 u, c+ e1 q. d; _- O; j/ N# n1 ]d = d1;

' h$ E8 {" w  d1 B3 }r e t u r n ; }
1 _" N/ L: q& Y3 M: W
if (d2 &lt;= d3) {a = X[l+1];

b = X[r];

d = d2;}

2 r* p; B$ H, W# \: ~5 l4 Z9 ^' helse {a = X[l];

: c; ]" c/ ]9 I% ]; I2 [b = X[r];

0 F& l1 D4 d' c" L) @( L  u) e4 vd = d3;}
$ O7 B5 D. o; |$ f
r e t u r n ; }

/ /多于三个点，划分为两部分

) ?0 Z# E, X  W) F/ ?2 u" Fint m = (l+r)/2; // X[l:m] 在A中，余下的在B中
5 u; t% a7 u9 m
// 在Z[l:m] 和Z [ m + 1 : r ]中创建按y排序的表
! O5 O" j" P6 Z6 P/ Y, ?# C+ {
int f = l, // Z[l:m]的游标

g = m+1; // Z[m+1:r]的游标
4 z) v" o0 m) u
0 k& j/ }+ E7 K5 }9 Z4 f8 Kfor (int i = l; i &lt;= r; i++)
; x% l* X! e3 E2 E% D
5 S- Q# @- M9 H; G) oif (Y.p &gt; m) Z[g++] = Y;

else Z[f++] = Y;

// 对以上两个部分进行求解

c l o s e ( X , Z , Y, l , m , a , b , d ) ;

5 I4 |4 e) w& T2 y9 X( M# hfloat dr;
* H4 H4 l  v) p0 H4 z$ ]
7 k8 J8 B" g" ]6 Y1 h' G4 RPoint1 ar, br;
  y6 n# `3 K2 ]- S) j# W( s
( Q4 Y$ J* `; [0 ^# F7 C" \& j+ Qc l o s e ( X , Z , Y, m + 1 , r, a r, b r, d r ) ;
# B% u* q+ i5 t9 ]4 o, H4 p
1 [4 }& Q1 m% K' \) l; W// (a,b) 是两者中较近的点对

if (dr &lt; d) {a = ar;

b = br;
: V; G0 m) v# R5 \# H
" n- C( [1 U% Sd = dr;}

' R; s# M- S4 k: n. @M e r g e ( Z , Y,l,m,r);// 重构Y
1 D& A3 R1 K. n
7 w8 Y" j+ y' l' _; m/ /距离小于d的点放入Z
# L( C: `$ s! M! L2 o$ f
1 i4 _9 g5 T7 z9 y9 U" Bint k = l; // Z的游标
; L, T# `. o2 t
for (i = l; i &lt;= r; i++)

if (fabs(Y[m].x - Y.x) &lt; d) Z[k++] = Y;
- f+ M3 T. |; F0 q* J
// 通过检查Z [ l : k - 1 ]中的所有点对，寻找较近的点对
* M: b  N$ D4 ?9 J
: b3 ^+ Q" ]  {, s+ W, N( Hfor (i = l; i &lt; k; i++){

+ V/ X% j( y; k; G* E4 Zfor (int j = i+1; j &lt; k &amp;&amp; Z[j].y - Z.y &lt; d;
; t$ g3 @; Z; P  G! f$ x
j + + ) {

! r: i4 b7 Q5 n8 ?float dp = dist(Z, Z[j]);

. a/ X& n8 X9 [- K2 Dif (dp &lt; d) {// 较近的点对

2 E6 {) _1 z1 a$ x- N6 Q& w( ?+ _d = dp;

a = X[Z.p];
4 |- T3 R# [6 v/ G) W4 y6 X
  u! K7 _+ K- s8 i, L. X2 c4 u% t3 s8 Xb = X[Z[j].p];}

8 H+ T+ z$ {$ ~4 A! N}
/ ]6 d& e/ ?! y" L2 t; j4 Q2 f
" u' |) x; O# C# p7 U1 T}

" ~& ?: r  b& k}
& v3 L3 b: y' {6 E
! B' G% P3 E, {+ d0 e* V, q函数c l o s e（见程序1 4 - 11）用来确定X[1:r] 中的最近点对。假定这些点按x 坐标排序。在Y [ 1 : r ]中对这些点按y 坐标排序。Z[ 1 : r ]用来存放中间结果。找到最近点对以后，将在a, b中返回最近点对，在d 中返回距离，数组Y被恢复为输入状态。函数并未修改数组X。

0 n: y. S' q7 y# H首先考察“小问题”，即少于四个点的点集。因为分割过程不会产生少于两点的数组，因此只需要处理两点和三点的情形。对于这两种情形，可以尝试所有的可能性。当点数超过三个时，通过计算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 ]。
# s+ J2 @# c5 K
为实现图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

还是在RB中）与Z[j] 配对来实现，其中i＜j 且Z [ j ] . y - Z [ i ] . y＜ 。对每一个Z [ i ]，在2 × 区域内所检查的点如图1 4 - 1 7所示。由于在每个2 × 子区域内的点至少相距。因此每一个子区域中的点数不会超过四个，所以与Z [ i ]配对的点Z [ j ]最多有七个。
. |8 }! \: Z# r5 |
2. 复杂性分析
4 `4 p2 ]6 I/ e
令t (n) 代表处理n 个点时，函数close 所需要的时间。当n＜4时，t (n) 等于某个常数d。当n≥4时，需花费(n) 时间来完成以下工作：将点集划分为两个部分，两次递归调用后重构Y，淘汰距分割线很远的点，寻找更好的第三类点对。两次递归调用需分别耗时t (「n /2ù」和t (?n /2?).
5 Y% l0 ^) b9 ^
+ U$ W3 Z& t- Q4 O/ _这个递归式与归并排序的递归式完全一样，其结果为t (n) = (nl o gn)。另外，函数c l o s e s t还需耗时(nl o gn)来完成如下额外工作：对X进行排序，创建Y和Z，对Y进行排序。因此分而治之最近点对求解算法的时间复杂性为(nl o gn)。</P>