【算法】经纬度常用计算

最近工作中遇到经纬度搜索的需求，初步想法是计算所有目标城市距该点的距离，然后进行筛选，但头疼的是，没有所有产品的缓存，计算距离的操作只能放到DB端，这样是不可接受的；所以打算先将所有产品放到缓存中，再进行计算。可这么做的话，一方面改造工时比较长，另一方面目前的缓存系统不是很稳定，几番思考征得产品经理同意后得出一个不精确的方形搜索方案。

即以目标点为中心，画一个正方型，在应用端根据目标点经纬度、范围距离、角度算出正方型左下点和右上点的经纬度，然后以此去DB里between。恩，在要求不精确且没有缓存的情况下这是一个较好的折中方案。

于是接下来就开始考虑算法，参考了博客园的帖子（http://www.cnblogs.com/hellofox2000/archive/2010/07/13/1776159.html#2042746），试验后发现计算两点间距离的方法偏差有点大，于是对其做了一些修改，作为工具类收藏起来，代码如下：

实体类：

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace Xinwei.Test.Geography
{
    public class GeographicPoint
    {
        /// <summary>
        /// 赤道半径
        /// </summary>
        public const double EARTH_RADIUS = 6378137;

        /// <summary>
        /// 极半径
        /// </summary>
        public const double POLAR_RADIUS = 6356725;

        /// <summary>
        /// 
        /// </summary>
        public GeographicPoint()
        { }

        /// <summary>
        /// 构造函数
        /// </summary>
        /// <param name="lat">维度</param>
        /// <param name="lon">经度</param>
        public GeographicPoint(double lat, double lon)
        {
            this.Latitude = lat;
            this.Longitude = lon;
        }

        /// <summary>
        /// 纬度
        /// </summary>
        public double Latitude
        { get; set; }

        /// <summary>
        /// 经度
        /// </summary>
        public double Longitude
        { get; set; }
        
        /// <summary>
        /// 纬度的弧度
        /// </summary>
        public double RadianOfLatitude
        {
            get
            {
                return Latitude * Math.PI / 180;
            }
        }
        
        /// <summary>
        /// 经度的弧度
        /// </summary>
        public double RadianOfLongitude
        {
            get
            {
                return Longitude * Math.PI / 180;
            }
        }

        /// <summary>
        /// 暂时不知意义，请大神们帮助
        /// </summary>
        public double Ec
        {
            get
            {
                return POLAR_RADIUS + (EARTH_RADIUS - POLAR_RADIUS) * (90 - Latitude) / 90;
            }
        }

        /// <summary>
        /// 暂时不知意义，请大神们帮助
        /// </summary>
        public double Ed
        {
            get
            {
                return Ec * Math.Cos(RadianOfLatitude);
            }
        }
    }
}

Helper类：

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace Xinwei.Test.Geography
{
    public static class GeographyHelper
    {
        /// <summary>
        /// 根据两点的经纬度计算两点距离
        /// </summary>
        /// <param name="sourcePoint">A点维度</param>        
        /// <param name="destinationPoint">B点经度</param>
        /// <returns></returns>
        public static double GetDistance(GeographicPoint sourcePoint, GeographicPoint destinationPoint)
        {
            if (Math.Abs(sourcePoint.Latitude) > 90 || Math.Abs(destinationPoint.Latitude) > 90 || Math.Abs(sourcePoint.Longitude) > 180 || Math.Abs(destinationPoint.Longitude) > 180)
                throw new ArgumentException("经纬度信息不正确！");

            double distance = GeographicPoint.EARTH_RADIUS / 1000 * Math.Acos(Math.Cos(sourcePoint.RadianOfLatitude)
                * Math.Cos(destinationPoint.RadianOfLatitude) * Math.Cos(destinationPoint.RadianOfLongitude - sourcePoint.RadianOfLongitude)
                + Math.Sin(sourcePoint.RadianOfLatitude) * Math.Sin(destinationPoint.RadianOfLatitude));

            return distance;
        }

        /// <summary>
        /// 已知点A经纬度，根据B点据A点的距离，和方位，求B点的经纬度
        /// </summary>
        /// <param name="sourcePoint">已知点A</param>
        /// <param name="distance">B点到A点的距离 </param>
        /// <param name="angle">B点相对于A点的方位，12点钟方向为零度，角度顺时针增加</param>
        /// <returns>B点的经纬度坐标</returns>
        public static GeographicPoint GetGeographicPoint(GeographicPoint sourcePoint, double distance, double angle)
        {
            double dx = distance * 1000 * Math.Sin(angle * Math.PI / 180);
            double dy = distance * 1000 * Math.Cos(angle * Math.PI / 180);

            double longitude = (dx / sourcePoint.Ed + sourcePoint.RadianOfLongitude) * 180 / Math.PI;
            double latitude = (dy / sourcePoint.Ec + sourcePoint.RadianOfLatitude) * 180 / Math.PI;

            GeographicPoint destinationPoint = new GeographicPoint(latitude, longitude);
            return destinationPoint;
        }
    }
}