GIS 纠偏，点距

妈的，这么个技术点都成香馍馍了，人家老早把算法私下颁布了啊（呵呵，我们公司好像10多年也没搞定，也是我来搞定的，不过不是我写的，我抄来的，写过游戏的直觉是对的，估计也是sin一类的）；

百度的二次偏移不支持，需要再次运算，以下算法为一次偏移纠偏运算（比如bing，google都只是一次偏移）；

上个图证明是OK的哈，图中的是silverlight bingmap ——

代码继续往下——

 

 

    internal class EvilTransform
    {
        private const double pi = 3.14159265358979324;

        //
        // Krasovsky 1940
        //
        // a = 6378245.0, 1/f = 298.3
        // b = a * (1 - f)
        // ee = (a^2 - b^2) / a^2;
        private const double a = 6378245.0;
        private const double ee = 0.00669342162296594323;

        //
        // World Geodetic System ==> Mars Geodetic System
        public static void Transform(double wgLat, double wgLon, out double mgLat, out double mgLon)
        {
            if (OutOfChina(wgLat, wgLon))
            {
                mgLat = wgLat;
                mgLon = wgLon;
                return;
            }
            double dLat = TransformLat(wgLon - 105.0, wgLat - 35.0);
            double dLon = TransformLon(wgLon - 105.0, wgLat - 35.0);
            double radLat = wgLat/180.0*pi;
            double magic = Math.Sin(radLat);
            magic = 1 - ee*magic*magic;
            double sqrtMagic = Math.Sqrt(magic);
            dLat = (dLat*180.0)/((a*(1 - ee))/(magic*sqrtMagic)*pi);
            dLon = (dLon*180.0)/(a/sqrtMagic*Math.Cos(radLat)*pi);
            mgLat = wgLat + dLat;
            mgLon = wgLon + dLon;
        }

        private static bool OutOfChina(double lat, double lon)
        {
            if (lon < 72.004 || lon > 137.8347)
                return true;
            if (lat < 0.8293 || lat > 55.8271)
                return true;
            return false;
        }

        private static double TransformLat(double x, double y)
        {
            double ret = -100.0 + 2.0*x + 3.0*y + 0.2*y*y + 0.1*x*y + 0.2*Math.Sqrt(Math.Abs(x));
            ret += (20.0*Math.Sin(6.0*x*pi) + 20.0*Math.Sin(2.0*x*pi))*2.0/3.0;
            ret += (20.0*Math.Sin(y*pi) + 40.0*Math.Sin(y/3.0*pi))*2.0/3.0;
            ret += (160.0*Math.Sin(y/12.0*pi) + 320*Math.Sin(y*pi/30.0))*2.0/3.0;
            return ret;
        }

        private static double TransformLon(double x, double y)
        {
            double ret = 300.0 + x + 2.0*y + 0.1*x*x + 0.1*x*y + 0.1*Math.Sqrt(Math.Abs(x));
            ret += (20.0*Math.Sin(6.0*x*pi) + 20.0*Math.Sin(2.0*x*pi))*2.0/3.0;
            ret += (20.0*Math.Sin(x*pi) + 40.0*Math.Sin(x/3.0*pi))*2.0/3.0;
            ret += (150.0*Math.Sin(x/12.0*pi) + 300.0*Math.Sin(x/30.0*pi))*2.0/3.0;
            return ret;
        }
    }

在来个算两点距离的（二维平面，不算三维的z哈）

public struct EarthPoint
    {
        public const double Ea = 6378137; // 赤道半径 WGS84标准参考椭球中的地球长半径(单位:m)  
        public const double Eb = 6356725; // 极半径   
        public readonly double Ec;
        public readonly double Ed;
        public readonly double Jd;
        public readonly double Latidute;
        public readonly double Longitude;
        public readonly double Wd;

        public EarthPoint(double longitude, double latidute)
        {
            Longitude = longitude;
            Latidute = latidute;
            Jd = Longitude*Math.PI/180; //转换成角度
            Wd = Latidute*Math.PI/180; //转换成角度
            Ec = Eb + (Ea - Eb)*(90 - Latidute)/90;
            Ed = Ec*Math.Cos(Wd);
        }

        public double Distance(EarthPoint point)
        {
            double dx = (point.Jd - Jd)*Ed;
            double dy = (point.Wd - Wd)*Ec;
            return Math.Sqrt(dx*dx + dy*dy);
        }

        public static double GetDistance(
            double latidute1,
            double longitude1,
            double latidute2,
            double longitude2)
        {
            var p1 = new EarthPoint(longitude1, latidute1);
            var p2 = new EarthPoint(longitude2, latidute2);
            return p1.Distance(p2);
        }
    }