#!/usr/local/bin/python3
# Author: Dr. Robert Heckendorn, Computer Science Department, University of Idaho, 2013
#
# This finds the distance between the closest set of points
# in a set of points.  This is done by divide and conquer to get
# an O(n log(n)) running time

from random import *
from math import *

## ## ## ## ## ## ## ## ## ## ## ## ## ## ## ## 
# simple cp

# NOTE: this is the SQUARE of distance between two points
def dist(p1, p2) : return (p1[0] - p2[0])*(p1[0] - p2[0]) + (p1[1] - p2[1])*(p1[1] - p2[1])

# for a small number of points find the minimum dist
# This is O(len(pts)^2)
def simplecp(pts, start, end) :
    # WARNING: assumes there are at least 2 points
    dmin = dist(pts[start], pts[start+1])
    for i in range(start, end) :
        for j in range(i+1, end) :
            newdist = dist(pts[i], pts[j])
            if newdist<dmin : 
                dmin = newdist
    return dmin



## ## ## ## ## ## ## ## ## ## ## ## ## ## ## ## 
# various sorting by x and y

# point compare first on x then y
def lessXY(a, b) :
    return (a[0]<b[0]) or (a[0]==b[0] and a[1]<b[1])

# point compare first on y then x
def lessYX(a, b) :
    return (a[1]<b[1]) or (a[1]==b[1] and a[0]<b[0])


# This is O(len(pts))
# simple partition function
def partitionaux(pts, left, right, compare) :
    mid = (left+right)//2   # guard against presorted lists
    (pts[mid], pts[left]) = (pts[left], pts[mid])
    s = left
    pivot = pts[s]
    for i in range(s+1, right) :
        if compare(pts[i], pivot) :    # if  pts[i]<pivot then move it to left side
            s += 1
            (pts[s], pts[i]) = (pts[i], pts[s])

    (pts[s], pts[left]) = (pts[left], pts[s])
    return s

# sort first on x and then y 
def sortXY(pts) :
    sortaux(pts, 0, len(pts), lessXY)
    return pts

# sort first on y and then x
def sortYX(pts) :
    sortaux(pts, 0, len(pts), lessYX)
    return pts

# This is O(n log(n)) where n = len(pts)
def sortaux(pts, left, right, compare) :
    if right-left > 1 :
        s = partitionaux(pts, left, right, compare)
        sortaux(pts, left, s, compare)
        sortaux(pts, s+1, right, compare)


# split ypts into two halves around the point pivot
# each half is SELECTED BY X < PIVOT POINT
# the points in each half are SORTED BY Y
# This is O(len(ypts))
def splitsortYX(ypts, pivot) :
    left = []
    right = []
    for i in range(0, len(ypts)) :
#        if ypts[i][0]<pivot[0] :  # yes you sort by x first to do the split
        if lessXY(ypts[i], pivot) :  # yes you sort by x first to do the split
            left.append(ypts[i])
        else: 
            right.append(ypts[i])
    return (left, right)


# takes list of tuples representing points
# This is O(n log(n)) where n is len(pts)
def cp(pts) :
    ypts = sortYX(pts[0:len(pts)])   # makes copy and sorts by y
    sortXY(pts)                      # force x to be sorted by x coordinate

    # pts is sorted by x and ypts is sorted by y
    return cpaux(pts, ypts, 0, len(pts))


    
# look up to end but not including end
# This is O(n log(n)) where n is len(pts)
def cpaux(pts, ypts, start, end) :
    if end-start<4 :
        return simplecp(pts, start, end)

    split = (start + end)//2
    (ypts1, ypts2) = splitsortYX(ypts, pts[split])  # pivot goes on the right

    # pts[start:split] is sorted by x
    # ypts1 is the same set of points but sorted by y
    d1 = cpaux(pts, ypts1, start, split)     # upto split but not including split
    d2 = cpaux(pts, ypts2, split, end)

    return cpcrossedge(ypts, pts[split], min(d1, d2))


# get the array of points within delta of x-coordinate pivotx SORTED BY Y.
# NOTE: that delta is the square of the distance!
# This is O(len(ypts))
def selectY(ypts, pivotx, delta) :
    ans = []
    low = pivotx - sqrt(delta)
    high = pivotx + sqrt(delta)
    for i in range(0, len(ypts)) :
        if  low <= ypts[i][0] <= high :
            ans.append(ypts[i])
    return ans


# Need consider no more than 8 points before moving on to next point.
# NOTE: that delta is the square of the distance!
# This is O(len(ypts))
def cpcrossedge(ypts, pivot, delta) :
    pts = selectY(ypts, pivot[0], delta)
    sqrdelta = sqrt(delta)

    for i in range(0, len(pts)) :
        j = i+1
        while j<len(pts) and pts[i][1] - pts[j][1] < sqrdelta :
            if dist(pts[i], pts[j]) < delta :
                delta = dist(pts[i], pts[j])  # could update sqrdelta but not necessary
            j += 1
    return delta


## ## ## ## ## ## ## ## ## ## ## ## ## ## ## ## 
# various testing


def main() :
    seed(3333331)

    pts = []
    for i in range(0, 10000) :
        pts.append( (random(), random()) )

    print("List size:", len(pts))

    # ATTENTION: uncomment one statement or the other
    # WARNING: brute force way can take a long long time
#    print(cp(pts))                        # our new code
    print(simplecp(pts, 0, len(pts)))   # old brute force way 

main()