【转】通过经纬度坐标计算距离的方法（经纬度距离计算）

最近在网上搜索“通过经纬度坐标计算距离的方法”，发现网上大部分都是如下的代码：

#define PI 3.14159265

static double Rc = 6378137;  // 赤道半径

static double Rj = 6356725;  // 极半径

class JWD

{

public:

double m_Longitude, m_Latitude;

double m_RadLo, m_RadLa;

double Ec;

double Ed;

public:

JWD(double longitude, double latitude)

{

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);

}

~JWD() {};

};

static JWD GetJWDB(JWD A, double x,double y)

{

double dx=x;

double dy=y;

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;

}

void main()

{

double referla=30.0;

double referlo=60.0;

double dx=500.0;

double dy=60.0;

JWD A(referla,referlo),B(0.0,0.0);

B=GetJWDB(A,dx,dy);

cout < < "  LA = " < < B.m_Latitude < < "  LO= " < < B.m_Longitude < < endl;

}

上面这段与之类似的代码是最容易通过搜索引擎找到的。大部分都是相互的转载，从来没有说明和注释。给初学者造成了很大的困扰。特别是：

Ec = Rj + (Rc - Rj) * (90.-m_Latitude) / 90.;

Ed = Ec * cos(m_RadLa);

Ec、Ed这2个参数，有人还在CSDN上发帖询问“函数中Ec，Ed是什么意思”。但从来没有见到有人回答。

提问的链接：

http://www.gisforum.net/bbs/TopicOther.asp?t=5&BoardID=33&id=155609

http://bbs.csdn.net/topics/320024634

我刚开始接触这个问题，在中文世界中也只能搜到这些Ctrl+C 到Ctrl+V的东西。全球最大的中文互联网上也没有任何解答。已经明白的人知道很简单，但初学者肯定很疑惑。更何况现在LBS应用已经非常多了。肯定有非常多的人已经找到了更好的解答方式。 但对于我来说，用最常用的关键词，最容易的找到了这些复制的答案，但疑点重重。本着好奇的心，经过一番了解，我把结果给大家分析一下。下面是C#的代码：

        public const double Ea = 6378137;     //   赤道半径  

        public const double Eb = 6356725;     //   极半径 

        private static void GetJWDB(double GLAT, double GLON, double distance, double angle, out double BJD, out double BWD)

        {

            double dx = distance * 1000 * Math.Sin(angle * Math.PI / 180.0);

            double dy = distance * 1000 * Math.Cos(angle * Math.PI / 180.0);

            //double ec = 6356725 + 21412 * (90.0 - GLAT) / 90.0;

            // 21412 是赤道半径与极半径的差

            double ec = Eb + (Ea-Eb) * (90.0 - GLAT) / 90.0;

            double ed = ec * Math.Cos(GLAT * Math.PI / 180);

            BJD = (dx / ed + GLON * Math.PI / 180.0) * 180.0 / Math.PI;

            BWD = (dy / ec + GLAT * Math.PI / 180.0) * 180.0 / Math.PI;

        }

上面这个函数一看就是懂中文的人搞的，名字叫GetJWDB，取得经纬度。他根据输入的经度、纬度、距离和一个角度，得到另外一个经纬度坐标，输出参数为BJD、BWD。

1）angle * Math.PI / 180.0 作用是将角度转换为弧度，经纬度坐标是角度值，计算时需要换为弧度。这里所有的计算都是用弧度。

2）函数以正北方(due north) 也就是指南针的方向为0度，顺时针方向增加。如下图，Distance距离如果是d的话，dx就是x轴方向的长度，即longitude经度方向的长度；dy就是y轴方向的长度，即latitude纬度方向的长度。

dx、dy的计算方式也可以是以正东(due east)方向为0度。那：dx=distance* Cos(θ)，dy=distance*Sin(θ)。区别是Cos 与Sin互换。（图上未标）

图1

3）Ea 表示赤道半径，Eb表示极半径，地球是一个近似球体，Ea与Eb稍微有点差距。ec的作用就是修正因为纬度不断变化的球半径长度。如果在GLAT=0，即在赤道上的时候，ec=Eb+(Ea-Eb)*(90-0)/90=Ea，那ec就刚好是赤道半径Ea；如果在极点GLAT=90，ec=Eb+(Ea-Eb)*(90-90)/90=Eb，那ec 就刚好是极半径Eb。

4）Ed是GLAT所在纬度的纬度圈的半径，如下图：

图2

截面过球心，此时截面的面积最大，此圆叫球的大圆(Great Cycle)，沿着经线进行截面，得到的都是大圆(Great Cycle)。球面被不经过球心. 的截面所截得的圆. 叫做小圆。纬度圈所在的圆是一个小圆。地球半径R，平均值R=6371.0Km， Ed如果用地球半径R表示，那就是Ed=R*Cos(θ)，可以参看

《根据2个经纬度点，计算这2个经纬度点之间的距离(通过经度纬度得到距离)》这里提到的“A、D所在纬度圆圈的半径（AO`）”。放到上面函数里，因为它不断修正地球半径ec，那就是ed = ec * Math.Cos(GLAT * Math.PI / 180);

5）按照地球经纬度的划分，如下图：

每等分的纬度之间，经度线段的长度是一定的。 A段，B段长度是一样的。

每一等分的经度之间，纬度线段的长度是从赤道向2极点减小的。C端，D段的长度是不一样的。

图3

参看上图，那dx / ed 就相当于是在GLAT这个纬度上dx长度与总长度的占比，算出来应该是个经度跨度。如果这个经度跨度加上起始给定的经度就是最终的经度。

同理 dy/R就是在GLON这个经度上的dy长度与地球平均半径R的占比，算出来应该是一个纬度跨度。如果这个纬度跨度加上起始给定的纬度就是最终的纬度。这里使用了R，取地球平均半径。

dy/ec 就是用不断修正的ec半径替代了平均半径R。

(dy / ec + GLAT * Math.PI / 180.0) 就是起始纬度加上distance距离的最终纬度，同时需要将该结果转换为角度。 转换角度方法是：弧度* 180.0 / Math.PI。

BWD = (dy / ec + GLAT * Math.PI / 180.0) * 180.0 / Math.PI;

(dx / ed + GLON * Math.PI / 180.0)就是起始经度加上distance距离的最终经度，同时需要将该结果转换为角度。

BJD = (dx / ed + GLON * Math.PI / 180.0) * 180.0 / Math.PI;

这个根据一个经纬度坐标、距离然后求另外一个经纬度坐标的作用，主要就是确定一个最小外包矩形(Minimum bounding rectangle，简称MBR)。例如，我要找一个坐标点(lat,lon)的5公里范围内的所有商户信息、景点信息等。这个MBR就是一个最大的范围，这个矩形是包含5公里范围内所有这些有效信息的一个最小矩形。利用公式，求出四个方向0度、90度、180度、270度方向上的四个坐标点就可以得到这个MBR。

        public const double Ea = 6378137;     //   赤道半径  

        public const double Eb = 6356725;     //   极半径 

        private static void GetlatLon(double LAT, double LON, double distance, double angle, out double newLon, out double newLat)

        {

            double dx = distance * 1000 * Math.Sin(angle * Math.PI / 180.0);

            double dy = distance * 1000 * Math.Cos(angle * Math.PI / 180.0);

            double ec = Eb + (Ea-Eb) * (90.0 - LAT) / 90.0;

            double ed = ec * Math.Cos(LAT * Math.PI / 180);

            newLon = (dx / ed + LON * Math.PI / 180.0) * 180.0 / Math.PI;

            newLat = (dy / ec + LAT * Math.PI / 180.0) * 180.0 / Math.PI;

        }

        public static void GetRectRange(double centorlatitude, double centorLogitude, double distance, out double maxLatitude, out double minLatitude, out double maxLongitude, out double minLongitude)

        {

            double temp = 0.0;

            GetlatLon(centorlatitude, centorLogitude, distance, 0, out temp, out maxLatitude);

            GetlatLon(centorlatitude, centorLogitude, distance, 180, out temp, out minLatitude);

            GetlatLon(centorlatitude, centorLogitude, distance, 90, out maxLongitude, out temp);

            GetlatLon(centorlatitude, centorLogitude, distance, 270, out minLongitude , out temp);

        }

这里得到的maxLatitude, minLatitude, maxLongitude, minLongitude就组成矩形的四个顶点。

网上另外的方法，

http://www.movable-type.co.uk/scripts/latlong.html

这里的“Destination point given distance and bearing from start point”一节介绍了。我直接把代码贴上来：

这里GetRectRange 这个函数 也是以正北方（due north）为起始角度，顺时针方向，求得maxLatitude, minLatitude, maxLongitude, minLongitude 这样一个MBR。2种方法的误差很小。我觉得都是可以用的公式。

         /// <summary>

        /// where    φ is latitude, λ is longitude, θ is the bearing (clockwise from north),

        /// δ is the angular distance d/R; d being the distance travelled, R the earth’s radius

        /// bearing 方位 0，90，180，270

        /// </summary>

        private static void GetNewLatLon(double lat, double lon, double d, double bearing, out double lat2, out double lon2)

        {

            lat2 = 0.0;

            lon2 = 0.0;

            double R = 6378.137;

            var φ1 = ConvertDegreesToRadians(lat);

            var λ1 = ConvertDegreesToRadians(lon);

            var θ = ConvertDegreesToRadians(bearing);

            var φ2 = Math.Asin(Math.Sin(φ1) * Math.Cos(d / R) +

                             Math.Cos(φ1) * Math.Sin(d / R) * Math.Cos(θ));

            var λ2 = λ1 + Math.Atan2(Math.Sin(θ) * Math.Sin(d / R) * Math.Cos(φ1),

                                     Math.Cos(d / R) - Math.Sin(φ1) * Math.Sin(φ2));

            λ2 = (λ2 + 3 * Math.PI) % (2 * Math.PI) - Math.PI; // normalise to -180..+180°

            lat2 = ConvertRadiansToDegrees(φ2);

            lon2 = ConvertRadiansToDegrees(λ2);

        }

如果有一个应用，表里存有100万的数据，数据包含一个lat、lon的经纬度信息。就可以先根据输入的经纬度和距离得到一个MBR，然后通过类似

SELECT Id

FROM IdInfoTable

WHERE latitude >= minLat AND latitude < maxLat

AND longitude >= minLon AND longitude < maxLon

的方式过滤掉大部分的数据，然后再通过《根据2个经纬度点，计算这2个经纬度点之间的距离(通过经度纬度得到距离)》这里提到的方法进行精细过滤。

完整代码：

using System;

using System.Collections.Generic;

using System.Linq;

using System.Text;

namespace GetJWD

{

    public class GetMBR

    {

        public double MaxLatitude;

        public double MinLatitude;

        public double MaxLongitude;

        public double MinLongitude;

        public double MaxLatitude2;

        public double MinLatitude2;

        public double MaxLongitude2;

        public double MinLongitude2;

        public GetMBR(double centorlatitude, double centorLogitude, double distance)

        {

            GetRectRange(centorlatitude, centorLogitude, distance, out MaxLatitude, out MinLatitude, out MaxLongitude, out MinLongitude);

            GetRectRange2(centorlatitude, centorLogitude, distance, out MaxLatitude2, out MinLatitude2, out MaxLongitude2, out MinLongitude2);

        }

        public const double Ea = 6378137;     //   赤道半径  

        public const double Eb = 6356725;     //   极半径 

        private static void GetlatLon(double LAT, double LON, double distance, double angle, out double newLon, out double newLat)

        {

            double dx = distance * 1000 * Math.Sin(angle * Math.PI / 180.0);

            double dy = distance * 1000 * Math.Cos(angle * Math.PI / 180.0);

            double ec = Eb + (Ea - Eb) * (90.0 - LAT) / 90.0;

            double ed = ec * Math.Cos(LAT * Math.PI / 180);

            newLon = (dx / ed + LON * Math.PI / 180.0) * 180.0 / Math.PI;

            newLat = (dy / ec + LAT * Math.PI / 180.0) * 180.0 / Math.PI;

        }

        public static void GetRectRange(double centorlatitude, double centorLogitude, double distance,

                                      out double maxLatitude, out double minLatitude, out double maxLongitude, out double minLongitude)

        {

            double temp = 0.0;

            GetlatLon(centorlatitude, centorLogitude, distance, 0, out temp, out maxLatitude);

            GetlatLon(centorlatitude, centorLogitude, distance, 180, out temp, out minLatitude);

            GetlatLon(centorlatitude, centorLogitude, distance, 90, out maxLongitude, out temp);

            GetlatLon(centorlatitude, centorLogitude, distance, 270, out minLongitude, out temp);

        }

        public static void GetRectRange2(double centorlatitude, double centorLogitude, double distance,

                                      out double maxLatitude, out double minLatitude, out double maxLongitude, out double minLongitude)

        {

            double temp = 0.0;

            GetNewLatLon(centorlatitude, centorLogitude, distance, 0, out maxLatitude, out temp);

            GetNewLatLon(centorlatitude, centorLogitude, distance, 180, out minLatitude, out temp);

            GetNewLatLon(centorlatitude, centorLogitude, distance, 90, out temp, out maxLongitude);

            GetNewLatLon(centorlatitude, centorLogitude, distance, 270, out temp, out minLongitude);

        }

        /// <summary>

        /// where    φ is latitude, λ is longitude, θ is the bearing (clockwise from north),

        /// δ is the angular distance d/R; d being the distance travelled, R the earth’s radius

        /// bearing 方位 0，90，180，270

        /// </summary>

        private static void GetNewLatLon(double lat, double lon, double d, double bearing, out double lat2, out double lon2)

        {

            lat2 = 0.0;

            lon2 = 0.0;

            double R = 6378.137;

            var φ1 = ConvertDegreesToRadians(lat);

            var λ1 = ConvertDegreesToRadians(lon);

            var θ = ConvertDegreesToRadians(bearing);

            var φ2 = Math.Asin(Math.Sin(φ1) * Math.Cos(d / R) +

                             Math.Cos(φ1) * Math.Sin(d / R) * Math.Cos(θ));

            var λ2 = λ1 + Math.Atan2(Math.Sin(θ) * Math.Sin(d / R) * Math.Cos(φ1),

                                     Math.Cos(d / R) - Math.Sin(φ1) * Math.Sin(φ2));

            λ2 = (λ2 + 3 * Math.PI) % (2 * Math.PI) - Math.PI; // normalise to -180..+180°

            lat2 = ConvertRadiansToDegrees(φ2);

            lon2 = ConvertRadiansToDegrees(λ2);

        }

        public static double ConvertDegreesToRadians(double degrees)

        {

            return degrees * Math.PI / 180;

        }

        public static double ConvertRadiansToDegrees(double radian)

        {

            return radian * 180.0 / Math.PI;

        }

    }

    class Test

    {

        static void Main(string[] args)

        {

            double latorg = 22.54587746, lonorg = 114.12873077;

            var gpsdis = new GetMBR(latorg, lonorg, 5);

            Console.WriteLine("maxlat:{0}, minlat:{1}, maxlon:{2}, minlon:{3}",

                gpsdis.MaxLatitude, gpsdis.MinLatitude, gpsdis.MaxLongitude, gpsdis.MinLongitude);

            Console.WriteLine("maxlat:{0}, minlat:{1}, maxlon:{2}, minlon:{3}",

                gpsdis.MaxLatitude2, gpsdis.MinLatitude2, gpsdis.MaxLongitude2, gpsdis.MinLongitude2);

        }

    }

}

本文转自:http://www.cnblogs.com/softfair/p/lat_lon_distance_bearing_new_lat_lon.html

时间： 2024-10-12 16:38:33

【转】通过经纬度坐标计算距离的方法（经纬度距离计算）的相关文章

php计算两个坐标（经度，纬度）之间距离的方法

本文实例讲述了php计算两个坐标(经度,纬度)之间距离的方法.分享给大家供大家参考.具体如下: 这里使用php计算两个坐标(经度,纬度)之间的距离,返回结果为米或者千米 function distance($lat1, $lng1, $lat2, $lng2, $miles = true) {  $pi80 = M_PI / 180;  $lat1 *= $pi80;  $lng1 *= $pi80;  $lat2 *= $pi80;  $lng2 *= $pi80;  $r = 6372.79

计算球面上经纬度坐标方法比较

计算球面上的两点(坐标为经纬度)之间的距离可以直接通过公式计算得到,也可以先将经纬度坐标转化为墨卡托投影坐标来,然后用平面中两点之间的距离公式来计算. 在网上找了一些代码,然后简单进行了测试,发现前者精度更高: 资料来源:http://0414.iteye.com/blog/2039199 http://blog.sina.com.cn/s/blog_8ab3b2cf0100xd69.html 1 package com.suncreate.spatialquery.web.utils; 2 3

Java 根据两个经纬度坐标计算距离

public class Distance{    private static final double EARTH_RADIUS = 6378137;    private static double rad(double d)    {       return d * Math.PI / 180.0;    }        /**     * 根据两点间经纬度坐标(double值),计算两点间距离,单位为米     * @param lng1     * @param lat1    

微信lbs---返回两个经纬度坐标点的距离

微信开发:lbs附近的商家,在数据库里记录商家的坐标,lbs设置里管理搜索半径,查询的时候,查询 客户当前坐标的半径内的所有商家列表.个人喜欢不一样,我选择了执行sql ,毕竟效果高点.微信开发必须得将就效率问题.不然等半天出不来,急死人,半天出不来结果,客户直接走人.不多说,直接上代码 \\返回两个经纬度坐标点的距离(单位:米) \\C#  方法 ///<summary>返回两个经纬度坐标点的距离(单位:米) by Alex.Y</summary>         ///<

根据2个经纬度点，计算这2个经纬度点之间的距离

原文  http://www.cnblogs.com/softfair/p/distance_of_two_latitude_and_longitude_points.html 根据 2 个经纬度点,计算这 2 个经纬度点之间的距离 球面上任意两点之间的距离计算公式可以参考维基百科上的下述文章. Great-circle distance Haversine formula 值得一提的是,维基百科推荐使用Haversine公式,理由是Great-circle distance公式用到了大量余弦函

使用ArcGIS实现WGS84经纬度坐标到北京54高斯投影坐标的转换

[摘 要] 本文针对从事测绘工作者普遍遇到的坐标转换问题,简要介绍ArcGIS实现WGS84经纬度坐标到北京54高斯投影坐标转换原理和步骤. [关键词] ArcGIS 坐标转换 投影变换 1 坐标转换简介 坐标系统之间的坐标转换既包括不同的参心坐标之间的转换,或者不同的地心坐标系之间的转换,也包括参心坐标系与地心坐标系之间的转换以及相同坐标系的 直角坐标与大地坐标之间的坐标转换,还有大地坐标与高斯平面坐标之间的转换.在两个空间角直坐标系中,假设其分别为O--XYZ和O--XYZ,如果两个 坐标系

计算GPS两点间的距离[单位为:米]

/**     * 计算GPS两点间的距离[单位为:米]     * @param center GPS当前数据(LonLat对象表示,LonLat.lon表示经度,LonLat.lat表示纬度)     * @param turnPoint 转向点经纬度对象     * @return     */    private double gpsDistance( LonLat center, LonLat turnPoint )    {        double distance = 0; 

地图经纬度坐标与屏幕坐标的转换(android版)

我们在开发GIS系统的时候,首先要解决的就是地图的可视化问题,这个问题的关键就在于如何把地图的坐标转换成屏幕坐标,然后才到渲染着色.标注等.以下以wgs84经纬度坐标为基准,介绍一下地图经纬度坐标与屏幕坐标转换的问题.   一.比例尺 首先地图本身就是表达地球上各种事物的空间分布.联系及时间中的发展变化状态绘制的图形.如何把地球上大量的空间信息显示在有限的地图上,这就是涉及到比例尺的问题.比例尺=图上距离与实际距离的比. 比例尺有resolution和scale两个概念 Resolution是指

由经纬度坐标得到腾讯地图的瓦片/切片行列号

参考:https://blog.csdn.net/mygisforum/article/details/22997879 https://blog.csdn.net/u013929284/article/details/53614281 https://blog.csdn.net/shaxiaozilove/article/details/54908569 项目需要根据数据的瓦片行列号建立文件目录,而已知的位置信息是经纬度,下面是解决方法: 一.经纬度坐标 => 投影坐标:腾讯地图使用的是Web