您现在的位置是：首页 >其他 >C#中高斯坐标转经纬度坐标网站首页其他

学习目标：

开发项目中，会遇到高斯坐标转经纬度坐标的情况，下面记录一个转换的工具方法

学习时间：

2025年2月10号星期一

学习产出：

代码逻辑如下

                /// <summary>
        /// 高斯转经纬度
        /// </summary>
        /// <param name="X">坐标X</param>
        /// <param name="Y">坐标Y</param>
        /// <returns>返回纬度X，经度Y</returns>
        public static double2 GasstoLatLon(double X, double Y)
        {
            int l = 0;
            double f, a, b, e2, e3;
            double K0, B0, Bf, Mf, Nf, Uf, y, Tf;
            double XX = X;
            double YY = Y;
            // 84坐标系统参数
            a = 6378137.0;
            f = 298.257223563;
            b = (f - 1.0) * a / f;
            e2 = (a * a - b * b) / (a * a);
            e3 = (a * a - b * b) / (b * b);
            int daihao = (int)(YY / 1000000);
            y = YY - daihao * 1000000 - 500000;
            // 带的中心经度
            l = daihao * 6 - 3;
            K0 = 1 + 0.75 * e2 + 45 * e2 * e2 / 64 + 350 * e2 * e2 * e2 / 512 + 11025 * e2 * e2 * e2 * e2 / 16834;
            B0 = XX / (a * (1 - e2) * K0);
            K0 = -(26328 * e2 * e2 * e2 * e2 / 16384) / 3;
            Bf = K0 * Math.Sin(B0);
            K0 = (604 * e2 * e2 * e2 / 512 + 68484 * e2 * e2 * e2 * e2 / 16384) / 3;
            Bf = (K0 + Bf) * Math.Sin(B0) * Math.Sin(B0);
            K0 = -(63 * e2 * e2 / 64 + 1108 * e2 * e2 * e2 / 512 + 58239 * e2 * e2 * e2 * e2 / 16384) / 3;
            Bf = (K0 + Bf) * Math.Sin(B0) * Math.Sin(B0);
            K0 = 0.5 * (0.75 * e2 + 45 * e2 * e2 / 64 + 350 * e2 * e2 * e2 / 512 + 11025 * e2 * e2 * e2 * e2 / 16384);
            Bf = B0 + Math.Sin(2 * B0) * (K0 + Bf);

            Mf = a * (1 - e2) / Math.Pow((1 - e2 * Math.Sin(Bf) * Math.Sin(Bf)), 1.5);
            Nf = a / Math.Sqrt(1 - e2 * Math.Sin(Bf) * Math.Sin(Bf));
            Tf = Math.Tan(Bf);
            Uf = Math.Sqrt(e3) * Math.Cos(Bf);
            K0 = Tf * (61 + 90 * Tf * Tf + 45 * Math.Pow(Tf, 4)) / (-720.0) / Mf / Math.Pow(Nf, 5);
            K0 = K0 + Tf * Math.Pow(y, 4) * (5 + 3 * Tf * Tf + Uf * Uf - 9 * Uf * Uf * Tf * Tf) / 24 / Mf / Math.Pow(Nf, 3);

            X = (Bf - Tf * y * y / (2 * Mf * Nf) + K0) / CD;//  B纬度

            K0 = 5 + 28 * Math.Tan(Bf) * Math.Tan(Bf) + 24 * Math.Pow(Math.Tan(Bf), 4) + 6 * e3 * Math.Cos(Bf) * Math.Cos(Bf) + 8 * e3 * Math.Sin(Bf) * Math.Sin(Bf);
            K0 = K0 * Math.Pow(y, 5) / (120 * Math.Pow(Nf, 5) * Math.Cos(Bf));
            K0 = K0 - Math.Pow(y, 3) * (1 + 2 * Math.Tan(Bf) * Math.Tan(Bf) + e3 * Math.Cos(Bf) * Math.Cos(Bf)) / (6 * Math.Pow(Nf, 3) * Math.Cos(Bf));
            K0 = K0 + y / (Nf * Math.Cos(Bf));

            Y = K0 / CD + l;//L

            return new double2(Y, X);
        }