1. 程式人生 > >大地測量學高斯投影正反算

大地測量學高斯投影正反算

namespace 大地測量學2
{
    public partial class Form1 : Form
    {
        public Form1()
        {
            InitializeComponent();
        }

        private void button1_Click(object sender, EventArgs e)
        {
            double ρ = 206264.806247096355, 
                 a = 6378245, b = 6356863.0187730473
, c = 6399698.9017827110, α = 1 / 298.3, e2 = 0.006693421622966, e3 = 0.006738525414683;//先BLA是角度,BmLmAm是弧度直接進行運算,可以直接進行三角函式運算,感覺在最上邊宣告MATH後就不用再一個一個的代用 double B = 51.0 + 38.0 / 60.0 + 43.9023 / 3600.0, L = 111.0 + 2.0 / 60.0 + 13.136 / 3600.0;//角度 B = B / 180 * Math.PI; L = L / 180 * Math.PI;//弧度 double
W = Math.Sqrt(1 - e2 * Math.Sin(B) * Math.Sin(B)); double V = Math.Sqrt(1 +e3 * Math.Cos(B) * Math.Cos(B)); // double N = a / W; double M = a * (1 - e2) / (W * W * W); double t = Math.Tan(B); double η = e3 * Math.Cos(B); B = 51.0 + 38.0
/ 60.0 + 43.9023 / 3600.0; L = 111 + 2.0 / 60.0+ 13.136 / 3600;//重新化為角度 B= 51.012195083333332 int n6,n3; n6 = (int)(L / 6) + 1;//投影帶編號19 n3 = (int)(L / 3 + 0.5); double L06 = 6 * n6 - 3;//L06即使這樣存在,是否有需要強制轉換 此時為111 double l2 = L - L06;//此時是角度 double zc = 2 * Math.PI * b + 4 * (a - b);//橢圓的周長40026876.244215891 double X6=zc/60/2*(2*(n6-1)+1);//12341620.1752999 B = B / 180 * Math.PI; L = L / 180 * Math.PI;//弧度 l2 = l2 / 180 * Math.PI;//沒卵用 //double x = X6 + N / 2 / ρ / ρ * Math.Sin(B) * Math.Cos(B) * l2 * l2 + N / 24 / ρ / ρ / ρ / ρ * Math.Sin(B) * Math.Cos(B) * Math.Cos(B) * Math.Cos(B) * (5 - t * t + 9 * η * η + 4 * η * η * η * η)*l2*l2 // * l2 * l2 + N / 720 / ρ / ρ / ρ / ρ / ρ / ρ * Math.Sin(B) * Math.Cos(B) * Math.Cos(B) * Math.Cos(B) * Math.Cos(B) * Math.Cos(B) * (61 - 58 * t * t + t * t * t * t) * l2 * l2 * l2 * l2 * l2 * l2; //double y = N / ρ * Math.Cos(B) * l2 + N / 6 / ρ / ρ / ρ * Math.Cos(B) * Math.Cos(B) * Math.Cos(B) * (1 - t * t + η * η) * l2 * l2 * l2 + N / 120 / ρ / ρ / ρ / ρ / ρ * Math.Cos(B) * Math.Cos(B) //* Math.Cos(B) * Math.Cos(B) * Math.Cos(B) * (5 - 18 * t * t + t * t * t * t + 14 * η * η - 58 * η * η * t * t) * l2 * l2 * l2 * l2 * l2; double N=6399698.902-(21562.267-(108.973-0.62*Math.Cos(B) * Math.Cos(B)) * Math.Cos(B)* Math.Cos(B)) * Math.Cos(B)* Math.Cos(B) ; double B2=(51.0 + 38.0 / 60.0 + 43.9023 / 3600.0)*3600;//這個為秒的格式 double a0 = (32140.404 - (135.3302 - (0.7092 - 0.004 * Math.Cos(B) * Math.Cos(B)) * Math.Cos(B) * Math.Cos(B)) * Math.Cos(B) * Math.Cos(B))/*a0*/; double a4 = (0.25 + 0.00252 * Math.Cos(B) * Math.Cos(B)) * Math.Cos(B) * Math.Cos(B) - 0.04166; double a6 = (0.166*Math.Cos(B) * Math.Cos(B) - 0.084) * Math.Cos(B) * Math.Cos(B); double a3 = (0.3333333 + 0.001123 * Math.Cos(B) * Math.Cos(B)) * Math.Cos(B) * Math.Cos(B) - 0.1666667; double a5 = 0.0083 - (0.1667 - (0.1968 + 0.004 * Math.Cos(B) * Math.Cos(B)) * Math.Cos(B) * Math.Cos(B)) * Math.Cos(B) * Math.Cos(B); double L0 = 6 * 19 - 3;//已知帶號為19 double l = L + L0;//角度 l=l*3600/ρ; double x=6367558.4969*B2/ρ-(a0-(0.5+(a4+a6*l*l)*l*l)*l*l*N)*Math.Sin (B)*Math .Cos (B); double y = (1 + (a3 + a5 * l * l) * l * l) * N * l * Math.Cos(B); Console.WriteLine("{0} {1} ", x, y); textBox1.Text = Convert.ToString(x); textBox2.Text = Convert.ToString(y); } private void button2_Click(object sender, EventArgs e) { textBox1.Text = ""; textBox2.Text = ""; } private void groupBox1_Enter(object sender, EventArgs e) { } private void button3_Click(object sender, EventArgs e) { double ρ = 206264.806247096355, a = 6378245, b = 6356863.0187730473, c = 6399698.9017827110, α = 1 / 298.3, e2 = 0.006693421622966, e3 = 0.006738525414683;//先BLA是角度,BmLmAm是弧度直接進行運算,可以直接進行三角函式運算,感覺在最上邊宣告MATH後就不用再一個一個的代用 //計運算元午面的弧長 double Round = 2 * Math.PI * b + 4 * (a - b); double X=5724004.82211011; double Y=502559.918404555; //以度°為單位的bf double Bf = 27.11115372595 + 9.02468257083 * (X / 1000000 - 3) - 0.00579740442 * (X / 1000000 - 3) * (X / 1000000 - 3) - 0.00043532572 * (X / 1000000 - 3) * (X / 1000000 - 3) * (X / 1000000 - 3) + 0.00004857285 * (X / 1000000 - 3) * (X / 1000000 - 3) * (X / 1000000 - 3) * (X / 1000000 - 3); double W = Math.Sqrt(1 - e2 * Math.Sin(Bf) * Math.Sin(Bf)); double V = Math.Sqrt(1 + e3 * Math.Cos(Bf) * Math.Cos(Bf)); double N = a / W; double M = a * (1 - e2) / (W * W * W); double t = Math.Tan(Bf); double η = e3 * Math.Cos(Bf); double B = Bf - t / 2 / M / N * Y * Y + t / 24 / M / N / N / N * (5 + 3 * t * t + η * η - 9 * η * η * t * t) * Y * Y * Y * Y; Bf=Bf/180*Math .PI; double l = 1 / N / Math.Cos(Bf) * Y - 1 / 6 / N / N / N / Math.Cos(Bf) * (1 + 2 * t * t + η * η) * Y * Y * Y + 1 / 120 / N / N / N / N / N / Math.Cos(Bf) * (5 + 28 * t * t + 24 * t * t * t * t) * Y * Y * Y * Y * Y; double L0 = 6 * 19 - 3;//已知帶號為19 double y = Y - N - 500000; double L = L0 + l;//角度 double D=(int)(L); double F=(int)((L-(int)(L))*60); double MM=(L-D-F/60)*3600; textBox3.Text = Convert.ToString(B); textBox4.Text = Convert.ToString(L); Console.WriteLine("{0} {1} {2}", D, F, MM);//111° 7分 35.6402392184452秒 } private void button4_Click(object sender, EventArgs e) { textBox3.Text = ""; textBox4.Text = ""; } } }

這裡寫圖片描述
實驗截圖如上,在這裡為了簡單隻用了一個固定資料,介面好不人性化。。。需要的讀者可以自己改換為相應的資料,需要完整程式的可以聯絡我。
感謝閱讀