from PADS.SVG import SVG,colors
import sys

pattern = [0] + [2**(i/2)*1j**(k/2) for i in range(-2,3) for k in range(8)]

def tanrad(A,p):
    """Radius of circle at p tangent to circle A"""
    Actr,Arad = A
    return abs(Actr-p)-Arad

def discrepancy(A,B,C,p):
    """How far are we from being triply tangent?"""
    rads = [tanrad(X,p) for X in (A,B,C)]
    return max(rads) - min(rads)

def positive(x,y,imag):
    if imag:
        x,y = x.imag,y.imag
    else:
        x,y = x.real,y.real
    return y > x

def inscribe(A,B,C,imag):
    """Find triply tangent circle in given direction from B.
    imag=0 for positive real dir, 1 for positive imaginary dir"""
    p = B[0] + 1.5*(1j)**imag*B[1]
    tries = 0
    d = discrepancy(A,B,C,p)
    while tries < 200 and d > 0.0001:
        candidates = [p+d*x for x in pattern if positive(B[0],p+d*x,imag)]
        best = min ([(discrepancy(A,B,C,q),q) for q in candidates])
        d,p = best
        tries += 1
    return (p,tanrad(A,p))

bbox = 480+360j

fills = [colors.blue,colors.yellow,colors.red]
fillidx = 0

svg = SVG(bbox,sys.stdout)
svg.group(stroke=colors.black)
svg.group(fill="#666")

def pack(slope,limit):
    Acenter = 240 - slope*240*1j
    Bcenter = 180j - 180/slope
    Arad = abs(Acenter)
    Brad = abs(Bcenter)
    A = (Acenter,Arad)
    B = (Bcenter,Brad)
    C = (bbox-Bcenter,Brad)
    D = (bbox-Acenter,Arad)
    circles = [A,B,C,D,None]
    offidx = 0
    while len(circles) < limit:
        E = inscribe(A,B,D,offidx)
        Cdist = tanrad(C,E[0])-E[1]
        if Cdist < 0:  # flip!
            circles.append(None)
            A,B,C,D = B,A,D,C
            offidx ^= 1
        else:
            circles.append(E)
            B = E
            if Cdist < 0.001:
                break
    return circles

def confrac(circles):
    f = 0
    for i in range(len(circles)-1,4,-1):
        if circles[i]:
            f += 1
        elif f:
            f = 1/f
    return f

target = 31/9
lo = 0.5 # gives 4 + 1/(3 + 1/...)
hi = 1   # gives 1 + 1/(1 + 1/...)

for i in range(30):
    mid = (lo+hi)/2
    circles = pack(mid,30)
    brooks = confrac(circles)
    if brooks > target:
        lo = mid
    else:
        hi = mid

circles = pack(mid,30)
for C in circles:
    if C:
        svg.circle(*C)
    else:
        svg.ungroup()
        svg.group(fill=fills[fillidx])
        fillidx = (fillidx + 1)%len(fills)

svg.ungroup()
svg.ungroup()
svg.close()