

 This page uses content from Wikipedia. The original article was at Closest pair of points problem. The list of authors can be seen in the page history. As with Rosetta Code, the text of Wikipedia is available under the GNU FDL. (See links for details on variance)



Task

Provide a function to find the closest two points among a set of given points in two dimensions,   i.e. to solve the   Closest pair of points problem   in the   planar   case.

The straightforward solution is a   O(n2)   algorithm   (which we can call brute-force algorithm);   the pseudo-code (using indexes) could be simply:

bruteForceClosestPair of P(1), P(2), ... P(N)
if N < 2 then
  return ∞
else
  minDistance ← |P(1) - P(2)|
  minPoints ← { P(1), P(2) }
  foreach i ∈ [1, N-1]
    foreach j ∈ [i+1, N]
      if |P(i) - P(j)| < minDistance then
        minDistance ← |P(i) - P(j)|
        minPoints ← { P(i), P(j) } 
      endif
    endfor
  endfor
  return minDistance, minPoints
 endif

A better algorithm is based on the recursive divide&conquer approach,   as explained also at   Wikipedia's Closest pair of points problem,   which is   O(n log n);   a pseudo-code could be:

closestPair of (xP, yP)
               where xP is P(1) .. P(N) sorted by x coordinate, and
                     yP is P(1) .. P(N) sorted by y coordinate (ascending order)
if N ≤ 3 then
  return closest points of xP using brute-force algorithm
else
  xL ← points of xP from 1 to ⌈N/2⌉
  xR ← points of xP from ⌈N/2⌉+1 to N
  xm ← xP(⌈N/2⌉)x
  yL ← { p ∈ yP : px ≤ xm }
  yR ← { p ∈ yP : px > xm }
  (dL, pairL) ← closestPair of (xL, yL)
  (dR, pairR) ← closestPair of (xR, yR)
  (dmin, pairMin) ← (dR, pairR)
  if dL < dR then
    (dmin, pairMin) ← (dL, pairL)
  endif
  yS ← { p ∈ yP : |xm - px| < dmin }
  nS ← number of points in yS
  (closest, closestPair) ← (dmin, pairMin)
  for i from 1 to nS - 1
    k ← i + 1
    while k ≤ nS and yS(k)y - yS(i)y < dmin
      if |yS(k) - yS(i)| < closest then
        (closest, closestPair) ← (|yS(k) - yS(i)|, {yS(k), yS(i)})
      endif
      k ← k + 1
    endwhile
  endfor
  return closest, closestPair
endif



References and further readings

   Closest pair of points problem
   Closest Pair (McGill)
   Closest Pair (UCSB)
   Closest pair (WUStL)
   Closest pair (IUPUI) 


import scala.collection.mutable.ListBuffer
import scala.util.Random
 
object ClosestPair {
  case class Point(x: Double, y: Double){
    def distance(p: Point) = math.hypot(x-p.x, y-p.y)
 
    override def toString = "(" + x + ", " + y + ")"
  }
 
  case class Pair(point1: Point, point2: Point) {
    val distance: Double = point1 distance point2
 
    override def toString = {
      point1 + "-" + point2 + " : " + distance
    }
  }
 
  def sortByX(points: List[Point]) = {
    points.sortBy(point => point.x)
  }
 
  def sortByY(points: List[Point]) = {
    points.sortBy(point => point.y)
  }
 
  def divideAndConquer(points: List[Point]): Pair = {
    val pointsSortedByX = sortByX(points)
    val pointsSortedByY = sortByY(points)
 
    divideAndConquer(pointsSortedByX, pointsSortedByY)
  }
 
  def bruteForce(points: List[Point]): Pair = {
    val numPoints = points.size
    if (numPoints < 2)
      return null
    var pair = Pair(points(0), points(1))
    if (numPoints > 2) {
      for (i <- 0 until numPoints - 1) {
        val point1 = points(i)
        for (j <- i + 1 until numPoints) {
          val point2 = points(j)
          val distance = point1 distance point2
          if (distance < pair.distance)
            pair = Pair(point1, point2)
        }
      }
    }
    return pair
  }
 
 
  private def divideAndConquer(pointsSortedByX: List[Point], pointsSortedByY: List[Point]): Pair = {
    val numPoints = pointsSortedByX.size
    if(numPoints <= 3) {
      return bruteForce(pointsSortedByX)
    }
 
    val dividingIndex = numPoints >>> 1
    val leftOfCenter = pointsSortedByX.slice(0, dividingIndex)
    val rightOfCenter = pointsSortedByX.slice(dividingIndex, numPoints)
 
    var tempList = leftOfCenter.map(x => x)
    //println(tempList)
    tempList = sortByY(tempList)
    var closestPair = divideAndConquer(leftOfCenter, tempList)
 
    tempList = rightOfCenter.map(x => x)
    tempList = sortByY(tempList)
 
    val closestPairRight = divideAndConquer(rightOfCenter, tempList)
 
    if (closestPairRight.distance < closestPair.distance)
      closestPair = closestPairRight
 
    tempList = List[Point]()
    val shortestDistance = closestPair.distance
    val centerX = rightOfCenter(0).x
 
    for (point <- pointsSortedByY) {
      if (Math.abs(centerX - point.x) < shortestDistance)
        tempList = tempList :+ point
    }
 
    closestPair = shortestDistanceF(tempList, shortestDistance, closestPair)
    closestPair
  }
 
  private def shortestDistanceF(tempList: List[Point], shortestDistance: Double, closestPair: Pair ): Pair = {
    var shortest = shortestDistance
    var bestResult = closestPair
    for (i <- 0 until tempList.size) {
      val point1 = tempList(i)
      for (j <- i + 1 until tempList.size) {
        val point2 = tempList(j)
        if ((point2.y - point1.y) >= shortestDistance)
          return closestPair
        val distance = point1 distance point2
        if (distance < closestPair.distance)
        {
          bestResult = Pair(point1, point2)
          shortest = distance
        }
      }
    }
 
    closestPair
  }
 
  def main(args: Array[String]) {
    val numPoints = if(args.length == 0) 1000 else args(0).toInt
 
    val points = ListBuffer[Point]()
    val r = new Random()
    for (i <- 0 until numPoints) {
      points.+=:(new Point(r.nextDouble(), r.nextDouble()))
    }
    println("Generated " + numPoints + " random points")
 
    var startTime = System.currentTimeMillis()
    val bruteForceClosestPair = bruteForce(points.toList)
    var elapsedTime = System.currentTimeMillis() - startTime
    println("Brute force (" + elapsedTime + " ms): " + bruteForceClosestPair)
 
    startTime = System.currentTimeMillis()
    val dqClosestPair = divideAndConquer(points.toList)
    elapsedTime = System.currentTimeMillis() - startTime
    println("Divide and conquer (" + elapsedTime + " ms): " + dqClosestPair)
    if (bruteForceClosestPair.distance != dqClosestPair.distance)
      println("MISMATCH")
  }
}
 