42 lines
1.6 KiB
Plaintext
42 lines
1.6 KiB
Plaintext
closestPair[ptsIn_] :=
|
|
Module[{xP, yP,
|
|
pts},(*Top level function.Sorts the pts by x and by y and then \
|
|
calls closestPairR[]*)pts = N[ptsIn];
|
|
xP = Sort[pts, #1[[1]] < #2[[1]] &];
|
|
yP = Sort[pts, #1[[2]] < #2[[2]] &];
|
|
closestPairR[xP, yP]]
|
|
|
|
closestPairR[xP_, yP_] :=
|
|
Module[{n, mid, xL, xR, xm, yL, yR, dL, pairL, dmin, pairMin, yS, nS,
|
|
closest, closestP, k,
|
|
cDist},(*where xP is P(1).. P(n) sorted by x coordinate,
|
|
and yP is P(1).. P(n) sorted by y coordinate (ascending order)*)
|
|
n = Length[xP];
|
|
If[n <= 3,(*Brute Force*)
|
|
Piecewise[{{{\[Infinity], {}},
|
|
n < 2}, {{EuclideanDistance[xP[[1]], xP[[2]]], {xP[[1]],
|
|
xP[[2]]}},
|
|
n == 2}, {Last@
|
|
MinimalBy[{{EuclideanDistance[xP[[1]], xP[[2]]], {xP[[1]],
|
|
xP[[2]]}}, {EuclideanDistance[xP[[1]], xP[[3]]], {xP[[1]],
|
|
xP[[3]]}}, {EuclideanDistance[xP[[3]], xP[[2]]], {xP[[3]],
|
|
xP[[2]]}}}, First], n == 3}}], mid = Ceiling[n/2];
|
|
xL = xP[[1 ;; mid]];
|
|
xR = xP[[mid + 1 ;; n]];
|
|
xm = xP[[mid]];
|
|
yL = Select[yP, #[[1]] <= xm[[1]] &];
|
|
yR = Select[yP, #[[1]] > xm[[1]] &];
|
|
{dL, pairL} = closestPairR[xL, yL];
|
|
{dmin, pairMin} = closestPairR[xR, yR];
|
|
If[dL < dmin, {dmin, pairMin} = {dL, pairL}];
|
|
yS = Select[yP, Abs[#[[1]] - xm[[1]]] <= dmin &];
|
|
nS = Length[yS];
|
|
{closest, closestP} = {dmin, pairMin};
|
|
Table[k = i + 1;
|
|
While[(k <= nS) && (yS[[k, 2]] - yS[[i, 2]] < dmin),
|
|
cDist = EuclideanDistance[yS[[k]], yS[[i]]];
|
|
If[cDist <
|
|
closest, {closest, closestP} = {cDist, {yS[[k]], yS[[i]]}}];
|
|
k = k + 1], {i, 1, nS - 1}];
|
|
{closest, closestP}](*end if*)]
|