根据两站点的经纬度求两站点间的距离
/**** 根据两站点的经纬度求两站点间的距离 ****/ double D_jw(double wd1,double jd1,double wd2,double jd2) { double x,y,out; double PI=3.14159265; double R=6.371229*1e6;x=(jd2-jd1)*PI*R*cos( ((wd1 wd2)/2) *PI/180)/180;
y=(wd2-wd1)*PI*R/180; out=hypot(x,y); return out/1000; }==
一个经纬度相关计算的C 类 写了一个经纬度距离计算的类
--------------CJWD.h--------------
#ifndef __JWD_AND_HELPER_20051005
#define __JWD_AND_HELPER_20051005#include "stdafx.h"
#include <math.h> #include <iostream> using namespace std;#ifndef PI
#define PI 3.14159265;
#endif
static double Rc = 6378137; // 赤道半径static double Rj = 6356725; // 极半径
namespace CDYW{
class JWD { public: double m_LoDeg, m_LoMin, m_LoSec; // longtitude 经度 double m_LaDeg, m_LaMin, m_LaSec; double m_Longitude, m_Latitude; double m_RadLo, m_RadLa; double Ec; double Ed; public: // 构造函数, 经度: loDeg 度, loMin 分, loSec 秒; 纬度: laDeg 度, laMin 分, laSec秒 JWD(double loDeg, double loMin, double loSec, double laDeg, double laMin, double laSec) { m_LoDeg=loDeg; m_LoMin=loMin; m_LoSec=loSec; m_LaDeg=laDeg; m_LaMin=laMin; m_LaSec=laSec; m_Longitude = m_LoDeg m_LoMin / 60 m_LoSec / 3600; m_Latitude = m_LaDeg m_LaMin / 60 m_LaSec / 3600; m_RadLo = m_Longitude * PI / 180.; m_RadLa = m_Latitude * PI / 180.; Ec = Rj (Rc - Rj) * (90.- m_Latitude) / 90.; Ed = Ec * cos(m_RadLa); }
//！
JWD(double longitude, double latitude) { m_LoDeg = int(longitude); m_LoMin = int((longitude - m_LoDeg)*60); m_LoSec = (longitude - m_LoDeg - m_LoMin/60.)*3600; m_LaDeg = int(latitude); m_LaMin = int((latitude - m_LaDeg)*60); m_LaSec = (latitude - m_LaDeg - m_LaMin/60.)*3600; m_Longitude = longitude; m_Latitude = latitude; m_RadLo = longitude * PI/180.; m_RadLa = latitude * PI/180.; Ec = Rj (Rc - Rj) * (90.-m_Latitude) / 90.; Ed = Ec * cos(m_RadLa); } };
class CJWDHelper
{ public: CJWDHelper() {}; ~CJWDHelper() {};
//! 计算点A 和 点B的经纬度，求他们的距离和点B相对于点A的方位
/*! * \param A A点经纬度 * \param B B点经纬度 * \param angle B相对于A的方位, 不需要返回该值，则将其设为空 * \return A点B点的距离 */ static double distance(JWD A, JWD B, double *angle) { double dx = (B.m_RadLo - A.m_RadLo) * A.Ed; double dy = (B.m_RadLa - A.m_RadLa) * A.Ec; double out = sqrt(dx * dx dy * dy); if( angle != NULL) { *angle = atan(fabs(dx/dy))*180./PI; // 判断象限 double dLo = B.m_Longitude - A.m_Longitude; double dLa = B.m_Latitude - A.m_Latitude; if(dLo > 0 && dLa <= 0) { *angle = (90. - *angle) 90.; } else if(dLo <= 0 && dLa < 0) { *angle = *angle 180.; } else if(dLo < 0 && dLa >= 0) { *angle = (90. - *angle) 270; } }
return out/1000;
}
//! 计算点A 和 点B的经纬度，求他们的距离和点B相对于点A的方位
/*! * \param longitude1 A点经度 * \param latitude1 A点纬度 * \param longitude2 B点经度 * \param latitude2 B点纬度 * \param angle B相对于A的方位, 不需要返回该值，则将其设为空 * \return A点B点的距离 */ static double distance( double longitude1, double latitude1, double longitude2, double latitude2, double *angle) { JWD A(longitude1,latitude1); JWD B(longitude2,latitude2);
return distance(A, B, angle);
}
//! 已知点A经纬度，根据B点据A点的距离，和方位，求B点的经纬度
/*! * \param A 已知点A * \param distance B点到A点的距离 * \param angle B点相对于A点的方位 * \return B点的经纬度坐标 */ static JWD GetJWDB(JWD A, double distance, double angle) { double dx = distance*1000 * sin(angle * PI /180.); double dy = distance*1000 * cos(angle * PI /180.); //double dx = (B.m_RadLo - A.m_RadLo) * A.Ed; //double dy = (B.m_RadLa - A.m_RadLa) * A.Ec;
double BJD = (dx/A.Ed A.m_RadLo) * 180./PI;
double BWD = (dy/A.Ec A.m_RadLa) * 180./PI; JWD B(BJD, BWD); return B; }
//! 已知点A经纬度，根据B点据A点的距离，和方位，求B点的经纬度
/*! * \param longitude 已知点A经度 * \param latitude 已知点A纬度 * \param distance B点到A点的距离 * \param angle B点相对于A点的方位 * \return B点的经纬度坐标 */ static JWD GetJWDB(double longitude, double latitude, double distance, double angle) { JWD A(longitude,latitude); return GetJWDB(A, distance, angle); }};
} #endif=========== 测试程序==========
#include "stdafx.h"
#include <math.h> #include <iostream>#include "CJWD.h" using namespace std;using namespace CDYW; double Rc = 6378137; // 赤道半径 double Rj = 6356725; // 极半径// 绵阳 double jd1 = 104.740999999; double wd1 = 31.4337;// 成都 double jd2 = 104.01; double wd2 = 30.40; int main(int argc, char* argv[]) { double angle = 0; cout << "A(绵阳): JD = " << jd1 << " WD = " << wd1 << endl; cout << "B(成都): JD = " << jd2 << " WD = " << wd2 << endl; cout << "--------------------" << endl; cout << D_jw(wd1,jd1,wd2,jd2, angle) << endl; cout << "angle: " << angle <<endl; cout << "==============" <<endl; JWD A(jd1,wd1),B(jd2,wd2); double distance = CJWDHelper::distance(jd1,wd1,jd2,wd2, &angle); //cout << CJWDHelper::distance(A,B, &angle) << endl; cout << distance << endl; cout << "angle: " << angle <<endl; cout << "==============" <<endl; JWD C = CJWDHelper::GetJWDB(A, distance, angle); cout << "JD = " << C.m_Longitude << " WD = " << C.m_Latitude << endl; cout << "==============" <<endl; cout << A.m_LoDeg << " " << A.m_LoMin << " " << A.m_LoSec << endl; return 0; }=====
通过两个点的经纬度计算距离 关键词： gis从google maps的脚本里扒了段代码，没准啥时会用上。大家一块看看是怎么算的。
private const double EARTH_RADIUS = 6378.137;
private static double rad(double d) { return d * Math.PI / 180.0; }public static double GetDistance(double lat1, double lng1, double lat2, double lng2)
{ double radLat1 = rad(lat1); double radLat2 = rad(lat2); double a = radLat1 - radLat2; double b = rad(lng1) - rad(lng2); double s = 2 * Math.Asin(Math.Sqrt(Math.Pow(Math.Sin(a/2),2) Math.Cos(radLat1)*Math.Cos(radLat2)*Math.Pow(Math.Sin(b/2),2))); s = s * EARTH_RADIUS; s = Math.Round(s * 10000) / 10000; return s; } 非常感谢，帮了我大忙了：） 虽然我也没看明白到底原理是什么，但验算了A（60，30），B（60，90）两点之间，此段代码和我用余弦定理算出来的结果很一致。 余弦定理的步骤是：1、算A、B弦长：地球半径R*cos（经度差60）=R/2； 2、算角AOB,O为地球圆心，利用余弦定理， cosAOB=(2R*R-(R/2)^2) /2*R*R=7/8; 3、弧AB的长为：R*arc cos(7/8);求毕 。