两点经纬度计算方位

//弧度 返回弧度 正北为零

double bearing(double lat1Rad, double lon1Rad,double lat2Rad, double lon2Rad)

{  

    double dLon = (lon2Rad-lon1Rad); 
    double y = sin(dLon) * cos(lat2Rad);
    double x = cos(lat1Rad)*sin(lat2Rad) - sin(lat1Rad)*cos(lat2Rad)*cos(dLon);
    double brng = atan2(y, x);
    return brng;
}