#diverse utilitati pentru lucrul cu coordonatele geografice;
#in principiu, in afara de haversine() restul metodelor sunt translatate din
#PHP (http://www.imaginerc.com/software/GeoCalc/) in Python
from math import *

def haversine(p1, p2):
    # convert to radians
    lon1, lon2 = p1['longitude'], p2['longitude']
    lat1, lat2 = p1['latitude'], p2['latitude']
    lon1 = lon1 * pi / 180 
    lon2 = lon2 * pi / 180 
    lat1 = lat1 * pi / 180 
    lat2 = lat2 * pi / 180 
    # haversine formula 
    dlon = lon2 - lon1 
    dlat = lat2 - lat1 
    a = sin(dlat/2)**2 + cos(lat1) * cos(lat2) * sin(dlon/2)**2
    c = 2 * atan2(sqrt(a), sqrt(1-a)) 
    km = 6367 * c 
    return km 

def deg2rad(grad_distance):
    #transformare din grade in radiani
    return grad_distance * pi / 180

def getKmPerLonAtLat(latitude):
    #transformam latitudinea din grade in radiani
    latitude = deg2rad(latitude)
    return 111.321 * cos(latitude)

def getLonPerKmAtLat(latitude):
    return 1 / getKmPerLonAtLat(latitude)

def getKmPerLat():
    return 1.0 / 111 

def get_polygon_area(polygon):
    #intoarce aria unui poligon
    #vezi http://local.wasp.uwa.edu.au/~pbourke/geometry/polyarea/
    sum = 0
    for i in range(0, len(polygon), 1):
        if i != (len(polygon) - 1):
            sum += polygon[i][0]*polygon[i+1][1] - polygon[i+1][0]*polygon[i][1]
        else:
            sum += polygon[i][0]*polygon[0][1] - polygon[0][0]*polygon[i][1]
    return abs(0.5 * sum)

def get_polygon_centroid(polygon, area):
    #obtine centroidul unui poligon
    #vezi http://local.wasp.uwa.edu.au/~pbourke/geometry/polyarea/
    sum_x = 0
    sum_y = 0
    for i in range(0, len(polygon), 1):
        if i != (len(polygon) - 1):
            sum_x += (polygon[i][0] + polygon[i+1][0]) * (polygon[i][0]*polygon[i+1][1] - polygon[i+1][0]*polygon[i][1])
            sum_y += (polygon[i][1] + polygon[i+1][1]) * (polygon[i][0]*polygon[i+1][1] - polygon[i+1][0]*polygon[i][1])
        else:
            sum_x += (polygon[i][0] + polygon[0][0]) * (polygon[i][0]*polygon[0][1] - polygon[0][0]*polygon[i][1])
            sum_y += (polygon[i][1] + polygon[0][1]) * (polygon[i][0]*polygon[0][1] - polygon[0][0]*polygon[i][1])
    return [abs(sum_x/(6*area)), abs(sum_y/(6*area))] 
