距离最近的点对

韩冰

<>给定n 个点（xi，yi）（1≤i≤n），要求找出其中距离最近的两个点。
& Z' b7 p5 p2 P
例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进行递归求解。
% k) W3 E1 x& m3 |! b  K* V9 I

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

R e t u r n ; }

$ j# }0 ?# w& [6 t' I& q1 U/ F3 d// n较大

9 P3 P- Q/ i' I" ~5 m7 V( Q将点集分成大致相等的两个部分A和B

( x% |8 V( F; w& }! C1 @- I7 ?5 a7 F确定A和B中的最近点对

确定一点在A中、另一点在B中的最近点对

从上面得到的三对点中，找出距离最小的一对点

* G% E& \, j1 A( L7 w( U图14-13 寻找最近的点对
# W( M6 F, O1 _% x! {" t1 c
, _: t. r6 {8 h+ W, s/ l1 U- Z
为了确定第三种情况下的最近点对，需要采用一种不同的方法。这种方法取决于点集是如何被划分成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。

设是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）。
% q9 m, e) t. F% a
例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) 是最近的点对。

8 ^/ l1 v7 {7 |  P* _/ B为了确定一个距离更小的第三类点，RA 中的每个点最多只需和RB 中的6个点比较，如图1 4 - 1 6所示。

2 k8 H# V3 T. n& j, W( Q& ~5 K4 ?/ W1. 选择数据结构

为了实现图1 4 - 1 3的分而治之算法，需要确定什么是“小问题”以及如何表示点。由于集合中少于两点时不存在最近点对，因此必须保证分解过程不会产生少于两点的点集。如果将少于四点的点集做为“小问题”，就可以避免产生少于两点的点集。
  |: V/ p, `6 l5 b
每个点可有三个参数：标号， x 坐标，y 坐标。假设标号为整数，每个点可用P o i n t l类（见程序1 4 - 8）来表示。为了便于按x 坐标对各个点排序，可重载操作符＜=。归并排序程序如1 4 -3所示。
& P9 L; z$ H: b/ d( ^
% d: c4 A9 \! ~; P1 N' \程序14-8 点类

class Point1 {

' A+ s/ S: S8 h' L, efriend float dist(const Point1&amp;, const Point1&amp;);

; M+ y: l7 h# ]. I# zfriend void close(Point1 *, Point2 *, Point2 *, int, int, Point1&amp;, Point1&amp;, float&amp;);

9 S  ~& `# c$ ~friend bool closest(Point1 *, int, Point1&amp;, Point1&amp;,float&amp;);

friend void main();

) s8 G' j0 U* o3 L, W% V/ Pp u b l i c :

int operator&lt;=(Point1 a) const

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

p r i v a t e :

. W! y9 ~9 K# N: t2 Y" ?" z) B4 n; ]int ID; // 点的编号
1 x. _/ E% T, T9 I; y1 O" w$ p- b
float x, y; // 点坐标
# k6 }! [6 k, M. a( \
+ j3 T1 P) d+ R: T: M4 B( I& p2 N} ;

# q  G  W4 a0 ?9 e' U; l; {class Point2 {
) L2 ?1 |$ L, B, L1 |
friend float dist(const Point2&amp;, const Point2&amp;);
8 d4 ?% n  k7 r5 t5 D
9 D: h' A+ k7 B) V* A" y. Ifriend void close(Point1 *, Point2 *, Point2 *, int, int, Point1&amp;, Point1&amp;, float&amp;);

0 Y& T3 g$ Y$ M, @1 H' Sfriend bool closest(Point1 *, int, Point1&amp;, Point1&amp;, float&amp;);

friend void main();
0 z' C" D# r$ l$ @, O7 X5 a/ D
p u b l i c :
4 s# B7 K1 Q) e
int operator&lt;=(Point2 a) const

{return (y &lt;= a.y);}
& u8 P$ k4 a- L5 }7 X( A% a
p r i v a t e :

% q$ i6 i5 u  Q, K! O+ q0 \3 B$ Rint p; // 数组X中相同点的索引

float x, y; // 点坐标

} ;

7 X" H; M/ s1 D5 T# Q所输入的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中的对应点。
9 ?1 g% R9 f+ n; [" C& 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类的友元。

程序14-9 计算两点距离

template<CLASS T>
: D) v9 a/ T7 E9 K9 c) T9 h
inline float dist(const T&amp; u, const T&amp; v)

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

7 t% d; b" R9 M9 l: S- Hfloat dx = u.x-v. x ;
, P+ _* S8 C' n  Z2 N
% @+ i2 S9 C- P5 z5 e. a( Pfloat dy = u.y-v. y ;
. |! b! F/ w9 m% |9 U
3 ^2 j: r8 n! Lreturn sqrt(dx * dx + dy * dy);

}
( e1 F5 a1 T1 H1 u6 e: l
4 e. E& l7 \0 i5 J如果点的数目少于两个，则函数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 )，该函数实际求解最近点对。

程序14-10 预处理及调用c l o s e

bool closest(Point1 X[], int n, Point1&amp; a, Point1&amp; b, float&amp; d)

{// 在n &gt;= 2 个点中寻找最近点对
5 ~* Q0 h. K1 M. E$ J. R7 l
% D# C: f5 P: i0 c6 P4 A// 如果少于2个点，则返回f a l s e

/ V- w/ U* W2 b" _: J* l! G/ G4 I# B// 否则，在a 和b中返回距离最近的两个点
, M4 a9 `; L2 }: j( F- G: B" k
! P0 E6 R$ z/ C/ e: l& iif (n &lt; 2) return false;
8 j9 V9 H, S  W; g
// 按x坐标排序
2 S, @6 K- v" I! {8 p, D( T
M e r g e S o r t ( X , n ) ;

// 创建一个按y坐标排序的点数组

Point2 *Y = new Point2 [n];
0 S) U* F$ m4 C2 D% ]
. E; F4 I3 w- t# _5 C& V% y. g! U6 afor (int i = 0; i &lt; n; i++) {

' u# @! }& ~/ p' n. m// 将点i 从X 复制到Y
3 Z. M! m4 I. D5 }6 [" A6 b
Y.p = i;
7 s+ s7 T/ Z$ H0 A6 G) y
0 V# i% l5 d$ u  `Y.x = X.x;
# i/ h7 S5 R8 M! z7 y! C
( S, F! c& v7 U. I& VY.y = X.y;
2 h; t4 {- N4 [
0 r7 }9 z& v9 C}
2 k& ^( {8 ~& ^" Q7 L% g" O
M e r g e S o r t ( Y,n); // 按y坐标排序
6 J. a& U9 Z' R9 \( w9 O2 u
// 创建临时数组

5 I8 S9 Z1 r$ @' |& WPoint2 *Z = new Point2 [n];

// 寻找最近点对
9 r' S7 i$ c7 G
c l o s e ( X , Y, Z , 0 , n - 1 , a , b , d ) ;

2 L2 A: w+ ~! m1 _; a( k6 C// 删除数组并返回
7 \! I9 C4 M. z
delete [] Y;

1 m# ^& |# u' o3 \delete [] Z;
5 q* e- J3 ~) G7 j
return true;
$ ^2 e* S0 B# F( Z, @
( C& L6 d, v: v% f! |}
) Y' t# x' L8 s( n$ {
" {, }5 ~2 |; w8 K4 A1 ?& [( I' ^程序1 4 - 11 计算最近点对
  }+ m) A6 Z! c0 K! r0 t! q
void close(Point1 X[], Point2 Y[], Point2 Z[], int l, int r, Point1&amp; a, Point1&amp; b, float&amp; d)
  q* O; ?7 r" z
{//X[l:r] 按x坐标排序
( \6 a1 C# c* B. F
( V) b/ |$ q0 B8 n//Y[l:r] 按y坐标排序

( ^0 V# Q# p4 rif (r-l == 1) {// 两个点

" w/ w( L( s. T' K& a8 ^, ^a = X[l];

b = X[r];

6 Y+ T% D7 g9 q' U8 Z: T3 |d = dist(X[l], X[r]);
# I% [+ ~$ K2 J5 j/ j6 B8 X
% [7 V$ w5 y& H% a' lr e t u r n ; }
( I* ], U# O# J5 W: t. |# `
if (r-l == 2) {// 三个点
" ^) K( n1 c6 i$ n( {; I
; @. A2 m4 Z3 I0 f0 }// 计算所有点对之间的距离

# ?% `- C6 T7 ^& k; \float d1 = dist(X[l], X[l+1]);

: h, x# m# Y# ?( u+ \float d2 = dist(X[l+1], X[r]);
  B( y% V* }  P# U8 u
float d3 = dist(X[l], X[r]);
8 ~- X2 m+ B" X+ h1 G3 v
// 寻找最近点对

+ w  G% ^5 }3 g  `, qif (d1 &lt;= d2 &amp;&amp; d1 &lt;= d3) {
, H, s% Q) L" \. M: ?+ F
. x/ ]9 y* e$ V' @+ `5 ta = X[l];
3 `* y7 E" S" K) C4 H
b = X[l+1];
# I/ P  p' ?1 x- a
7 A& w( e2 y+ g7 ^( y/ W, ?5 ^d = d1;

+ @2 j! X8 Z9 u  ir e t u r n ; }

) O) ?5 d& l3 x5 o1 u, x8 M4 B" Lif (d2 &lt;= d3) {a = X[l+1];
$ z! \$ d) c$ }+ z* V0 |9 s
. m$ V4 r6 X' F# |% c( F  ?3 e% Q6 T. nb = X[r];
: w! \+ c5 i% F" v
; m: i* z9 z% j8 Nd = d2;}
- L' x! A+ ], h* m
) o, J  |5 B( j1 `( y/ oelse {a = X[l];
. q0 x: J/ N8 Q0 i& P$ c% _
b = X[r];
0 N! e& i* d/ S. e2 k5 l4 c8 n
d = d3;}
' F+ |- z3 D7 d6 O) [/ r
: W" z' W6 V' N* W1 b; Ur e t u r n ; }

* W0 U  a8 N0 Y  n  v: v/ e/ /多于三个点，划分为两部分

* i, {$ Z& C; G! Q4 I; E: h+ pint m = (l+r)/2; // X[l:m] 在A中，余下的在B中

' s! L; V( R+ v( F9 T& k( c" v// 在Z[l:m] 和Z [ m + 1 : r ]中创建按y排序的表
% i; n7 n% m- |9 k) n1 Z" Q7 \" f
% i$ V* G2 i: V$ {" V; F' w, S) |int f = l, // Z[l:m]的游标

g = m+1; // Z[m+1:r]的游标
- n' C' T* z# C/ E! m5 w0 p+ B
2 h2 Z& i/ k# b! s9 y7 S% e8 Xfor (int i = l; i &lt;= r; i++)
% {6 |2 p- J- B0 ]0 R9 c/ y
if (Y.p &gt; m) Z[g++] = Y;

else Z[f++] = Y;
2 ]4 C2 k6 `4 q' z9 T! ?5 ]. _# |
  J: t9 }$ I& t" n  ^// 对以上两个部分进行求解
5 ^/ V: Q4 n  k7 [, K" @
: i  ^& y& W! s! F; Z5 Bc l o s e ( X , Z , Y, l , m , a , b , d ) ;

float dr;
1 h; J; `7 M5 [) d$ n5 k# M! z4 C, W
: ]# ?7 ]* C$ {, QPoint1 ar, br;
8 V; ]* I8 |5 t$ e/ f
c l o s e ( X , Z , Y, m + 1 , r, a r, b r, d r ) ;
- \6 e9 @" Y8 ^( f) U
// (a,b) 是两者中较近的点对

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

$ }( d( [& {1 Q3 W) Z" lb = br;
. |2 ]$ ?2 s9 g; f6 e+ p
8 M$ T1 {/ J& t5 }& W. Bd = dr;}

M e r g e ( Z , Y,l,m,r);// 重构Y
7 R2 b( u& i3 m! r+ {- L- k
/ /距离小于d的点放入Z

int k = l; // Z的游标
1 j8 I; w# x3 O, s! `
for (i = l; i &lt;= r; i++)
2 s5 W$ {5 `$ A! o- Q; w
+ f5 |5 ]6 s$ uif (fabs(Y[m].x - Y.x) &lt; d) Z[k++] = Y;
* H! _8 F4 [4 r
7 {/ r( g  E2 E& b// 通过检查Z [ l : k - 1 ]中的所有点对，寻找较近的点对
  U# R3 C: _4 A' m* k# m# |
2 M( F6 K: ]# x, P" v- k) q* ^for (i = l; i &lt; k; i++){
1 ^; U1 }# X0 F" ]
for (int j = i+1; j &lt; k &amp;&amp; Z[j].y - Z.y &lt; d;
# j: s& _; y1 i6 ?3 n  I
9 B, c& O/ s( }) _# n3 ^! jj + + ) {

% G- d! Q% `7 O5 s& M0 s( cfloat dp = dist(Z, Z[j]);
3 A* |( B  q, e. T
if (dp &lt; d) {// 较近的点对
# |: G( p& `: C  e! _7 Y
# e: L5 [+ m$ n4 ]4 |% Dd = dp;
) }2 [& S8 s. Y
a = X[Z.p];
& P7 x. R( j5 B" k
! A7 O- M3 ?9 u- u, Ib = X[Z[j].p];}
) X4 y3 L( t1 W9 x4 E. ?
$ m, n$ i5 `$ R2 s}

}
9 i( m( O! }& O% M' K$ u2 L
, v3 m6 _) P* K}

函数c l o s e（见程序1 4 - 11）用来确定X[1:r] 中的最近点对。假定这些点按x 坐标排序。在Y [ 1 : r ]中对这些点按y 坐标排序。Z[ 1 : r ]用来存放中间结果。找到最近点对以后，将在a, b中返回最近点对，在d 中返回距离，数组Y被恢复为输入状态。函数并未修改数组X。
% [: Y, H2 C- m, m) o- F
$ H  W3 k' i4 e+ A) h$ I首先考察“小问题”，即少于四个点的点集。因为分割过程不会产生少于两点的数组，因此只需要处理两点和三点的情形。对于这两种情形，可以尝试所有的可能性。当点数超过三个时，通过计算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 ]。
( X0 Q/ g0 `0 G% R! p+ E
9 l/ `$ |+ T# v+ |8 l/ J为实现图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 ]最多有七个。
0 i! H9 }% H5 F6 z; X3 g' G) E" r# }
2. 复杂性分析
# U6 V8 ]# l6 z
令t (n) 代表处理n 个点时，函数close 所需要的时间。当n＜4时，t (n) 等于某个常数d。当n≥4时，需花费(n) 时间来完成以下工作：将点集划分为两个部分，两次递归调用后重构Y，淘汰距分割线很远的点，寻找更好的第三类点对。两次递归调用需分别耗时t (「n /2ù」和t (?n /2?).

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