# -*- coding: utf-8 -*-
import math
import pyproj


def lonlat2utm(lon,lat):
    z=int(lon/6+31)
    proj = pyproj.Proj(proj='utm',zone=z,ellps='WGS84')
    return proj(lon, lat),z

def utm2lonlat(x,y,z):
    proj = pyproj.Proj(proj='utm',zone=z,ellps='WGS84')
    return proj(x, y,inverse=True)

def insec(p1,r1,p2,r2):
    x = p1[0]
    y = p1[1]
    R = r1
    a = p2[0]
    b = p2[1]
    S = r2
    d = math.sqrt((abs(a-x))**2 + (abs(b-y))**2)
    if d > (R+S) or d < (abs(R-S)):
#        print ("没有公共点")
        return 
    elif d == 0 and R==S :
#        print ("两个圆同心")
        return
    else:
        A = (R**2 - S**2 + d**2) / (2 * d)
        h = math.sqrt(R**2 - A**2)
        x2 = x + A * (a-x)/d
        y2 = y + A * (b-y)/d
        x3 = x2 - h * (b - y) / d
        y3 = y2 + h * (a - x) / d
        x4 = x2 + h * (b - y) / d
        y4 = y2 - h * (a - x) / d
        c1=[x3, y3]
        c2=[x4, y4]
        return c1,c2

def location_trans(p1,r1,p2,r2):
    z1=lonlat2utm(p1[0],p1[1])
    z2=lonlat2utm(p2[0],p2[1])
    z=int((z1[1]+z2[1])/2)    
    C=insec(z1[0],r1,z2[0],r2)
    if C:
        a=utm2lonlat(C[0][0],C[0][1],z)
        b=utm2lonlat(C[1][0],C[1][1],z)
        return a,b
    else:
        return None,None

def location_min(p1,p2,p,r):
    d1=math.fabs(r-math.sqrt((p[0]-p1[0])**2+(p[1]-p1[1])**2))
    d2=math.fabs(r-math.sqrt((p[0]-p2[0])**2+(p[1]-p2[1])**2))
    if d1<d2:
        return p1
    else:
        return p2
    
def location_judg(p1,r1,p2,r2,p3,r3):
    li=[]
    
    z1=lonlat2utm(p1[0],p1[1])
    z2=lonlat2utm(p2[0],p2[1])
    z3=lonlat2utm(p3[0],p3[1])
    
    z12=int((z1[1]+z2[1])/2)
    z13=int((z1[1]+z3[1])/2)
    z23=int((z2[1]+z3[1])/2)
    z=int((z12+z13+z23)/3)
    
    C12=insec(z1[0],r1,z2[0],r2)
    C13=insec(z1[0],r1,z3[0],r3)
    C23=insec(z2[0],r2,z3[0],r3)
    
    if C12:
        m12=location_min(C12[0],C12[1],z3[0],r3)
        li.append(utm2lonlat(m12[0],m12[1],z12))
    else:
        li.append(None)
    if C13:
        m13=location_min(C13[0],C13[1],z2[0],r2)
        li.append(utm2lonlat(m13[0],m13[1],z13))
    else:
        li.append(None)
    if C23:
        m23=location_min(C23[0],C23[1],z1[0],r1)
        li.append(utm2lonlat(m23[0],m23[1],z23))
    else:
        li.append(None)
        
    if C12 and C13 and C23:
#        print("三个坐标作的圆都有公共点")
        m=[(m12[0]+m13[0]+m23[0])/3,(m12[1]+m13[1]+m23[1])/3]
        li.append(utm2lonlat(m[0],m[1],z))
        return li
    elif C12 or C13 or C23:
#        print("三个坐标作的圆不全有公共点")
        li.append(None)
        return li
    else:
#        print("三个坐标作的圆都没有公共点")
        return 
    
    
if __name__ == "__main__":
    a=[[114.304569,30.593354],300000]
    b=[[115.857972,28.682976],400000]
    c=[[116.378517,39.865246],900000]
    
    print(location_trans(b[0],b[1],c[0],c[1]))    
    print(location_judg(a[0],a[1],b[0],b[1],c[0],c[1]))