java经纬度凸包graham_python 凸包(经纬度)+面积［近似］

最近需要用到凸包算法，在网上找了一些不太对，可能是直角坐标系的。

发现一篇文章，简单修改一下就可以用了，效果很不错，在这里分享给大家。

def cross(A,B):

return A[0] * B[1] - A[1] * B[0]

def vectorMinus( a , b):

return ( (a[0] - b[0] )*1000,(a[1] - b[1] )*1000)

def getLTDis( A, B ):

lon1, lat1, lon2, lat2 = map(radians, [A[0], A[1], B[0], B[1]])

dlon = lon2 - lon1

dlat = lat2 - lat1

a = sin(dlat/2)**2 + cos(lat1) * cos(lat2) * sin(dlon/2)**2

c = 2 * asin(sqrt(a))

r = 6371.393

#print A,B

return c * r * 1000.0

def triangleAre(A,B,C):

x,y,z = getLTDis(A,B),getLTDis(B,C),getLTDis(C,A)

c = (x + y + z) /2

return sqrt((c)*(c-y)*(c-z)*(c-x))

def grahamScanArea(data):

data.sort(key=lambda x:(x[0],x[1]),reverse=False)

ans = [ 0 ] * (len(data)*2)

m = 0

for item in data:

top = len(item)

while( m > 1 and cross( vectorMinus(ans[ m -1 ] , ans [ m - 2 ]), vectorMinus( item , ans [ m - 2 ] )) <= 0 ) : m = m -1

ans[m] = item

m = m + 1

k = m

flag = True

data.reverse()

for item in data:

if flag :

flag = False

continue

while( m > k and cross( vectorMinus(ans[ m -1 ] , ans [ m - 2 ]), vectorMinus( item , ans [ m - 2 ] )) <= 0) : m = m - 1

ans [m] = item

m = m + 1

m = m -1

b = [ ans[i] for i in range(0, m)]

if len(b) < 3 : return 0

#if DEBUG : print b

return AREA(b)

def AREA(b):

ans = 0.0

for i in range(len(b)):

if i == 0 or i + 1 >= len(b) : continue

x , y = b[i] , b[i + 1]

ans += triangleAre( b[0] , x , y )

return ans 转自：http://www.cnblogs.com/shuly/p/5810253.html