高斯正反算程序代碼

1、 . #include<math.h>#include<stdio.h>#define PI 3.97932#define P 206264.5void main()long double AngleToRadian(long double alfa); long double RadianToAngle(long double alfa); long double Bf(long double a,long double b,long double x); void GSZS(long double a,long double b,long double l,long

2、 double B,long double *x,long double *y); void GSFS(long double a,long double b,long double Bf,long double y,long double *B,long double *l); char XZ,XZZF,XZSL; int N,num; long double a,b,B,L,x,y; char YN='O' long double L1,B1,L0,l; long double *pointer_x,*pointer_y; long double FSB,FSL;long

3、double *pointer_B,*pointer_L;long double DH; pointer_B=&FSB;pointer_L=&FSL;pointer_x=&x;pointer_y=&y;printf("n=歡迎使用高斯投影計(jì)算程序=n");printf(" 說(shuō)明n");printf("1.程序可實(shí)現(xiàn)在任意橢球上進(jìn)行高斯三度帶和六度帶正算和反算;n");printf("2.經(jīng)緯度的輸入輸出格式為：度分秒,例如：113°05'13.68466"輸入時(shí)為:

4、113.051368466;n");printf("3.y坐標(biāo)的輸入輸出格式為：帶號(hào)*10E6+500公里+自然坐標(biāo),例如：y=18637682.388,18為帶號(hào),自然坐標(biāo)為137682.388;n");printf("4.其余可按照提示完成,如有疑問(wèn)可發(fā)送Email至gys_long126.,我們將盡快為您解答。");printf("n=n");printf("請(qǐng)選擇計(jì)算方式：A.三度帶正算;tB.三度帶反算;tC.六度帶正算;tD.六度帶反算n");printf("請(qǐng)選擇字母后回車(chē)：&qu

5、ot;);scanf("%c",&XZ);if(XZ='A'|XZ='a')XZZF='Z'XZSL='S'printf("注意：您將要進(jìn)行的是三度帶正算!");else if(XZ='B'|XZ='b')XZZF='F'XZSL='S'printf("注意：您將要進(jìn)行的是三度帶反算!");else if(XZ='C'|XZ='c')XZZF='Z'XZ

6、SL='L'printf("注意：您將要進(jìn)行的是六度帶正算!");else if(XZ='D'|XZ='d')XZZF='F'XZSL='L'printf("注意：您將要進(jìn)行的是六度帶反算!");else printf("警告：您操作有誤!"); printf("是否繼續(xù)Y/N:");scanf("%s",&YN);printf("n");doif(YN='N'|YN=&#

7、39;n')printf("歡迎您繼續(xù)使用此程序,如有疑問(wèn)可發(fā)送Email至gys_long126.,我們將盡快為您解答!n");break; if(XZZF='Z'|XZZF='z')printf("請(qǐng)選擇橢球：n1.克拉索夫斯基橢球; 2.1975年國(guó)際橢球體; 3.WGS-84橢球體; 4.自定義橢球n");printf("請(qǐng)選擇數(shù)字后回車(chē)：");scanf("%d",&num);switch(num)case 1: a=6378245.0000000000;b

8、=6356863.0187730473;printf("注意：您選擇的是克拉索夫斯基橢球!n請(qǐng)輸入經(jīng)度L和緯度B:(輸入時(shí)經(jīng)度和緯度之間用逗號(hào)隔開(kāi))n");scanf("%lf,%lf",&L,&B);L0=(XZSL='L'|XZSL='l')?(N=(int)(L+3)/6.+0.5)*6-3):(N=(int)(L/3.+0.5)*3);printf("t=中間計(jì)算過(guò)程結(jié)果=nn");printf("ttN=%dnttL0=%fn",N,L0);B1=Angle

9、ToRadian(B);/printf("HD=%lf",B1);L0=AngleToRadian(L0);L1=AngleToRadian(L);l=L1-L0;GSZS(a,b,l,B1,pointer_x,pointer_y);printf("nt=計(jì)算結(jié)果=n");printf("nttx=%-20.4lftty=%-20.4lfn",x,(y+500000+N*1000000);printf("nt=n");break;case 2: a=6378140.0000000000;b=6356755.2881

10、575287;printf("注意：您選擇的是1975年國(guó)際橢球體!n請(qǐng)輸入經(jīng)度L和緯度B:(輸入時(shí)經(jīng)度和緯度之間用逗號(hào)隔開(kāi))n");scanf("%lf,%lf",&L,&B);L0=(XZSL='L'|XZSL='l')?(N=(int)(L+3)/6.+0.5)*6-3):(N=(int)(L/3.+0.5)*3);printf("t=中間計(jì)算過(guò)程結(jié)果=nn");printf("ttN=%dnttL0=%fn",N,L0);B1=AngleToRadian(B)

11、;L0=AngleToRadian(L0);L1=AngleToRadian(L);l=L1-L0;GSZS(a,b,l,B1,pointer_x,pointer_y);printf("nt=計(jì)算結(jié)果為=n");printf("nttx=%-20.3lftty=%-20.3lfn",x,(y+500000+N*1000000);printf("nt=n");break;case 3: a=6378137.0000000000;b=6356752.3142;printf("注意：您選擇的是WGS-84橢球體!n請(qǐng)輸入經(jīng)度L和緯

12、度B:(輸入時(shí)經(jīng)度和緯度之間用逗號(hào)隔開(kāi))n");scanf("%lf,%lf",&L,&B);L0=(XZSL='L'|XZSL='l')?(N=(int)(L+3)/6.+0.5)*6-3):(N=(int)(L/3.+0.5)*3);printf("t=中間計(jì)算過(guò)程結(jié)果=nn");printf("ttN=%dnttL0=%fn",N,L0);B1=AngleToRadian(B);L0=AngleToRadian(L0);L1=AngleToRadian(L);l=L1-L

13、0;GSZS(a,b,l,B1,pointer_x,pointer_y);printf("nt=計(jì)算結(jié)果為=n");printf("nttx=%-20.3lftty=%-20.3lfn",x,(y+500000+N*1000000);printf("nt=n");break;case 4: printf("注意：您選擇的是自定義橢球體!n請(qǐng)輸入長(zhǎng)半軸a、短半軸b、經(jīng)度L和緯度B:(輸入時(shí)數(shù)據(jù)之間用逗號(hào)隔開(kāi))n");scanf("%lf,%lf,%lf,%lf",&a,&b,&am

14、p;L,&B);L0=(XZSL='L'|XZSL='l')?(N=(int)(L+3)/6.+0.5)*6-3):(N=(int)(L/3.+0.5)*3);printf("t=中間計(jì)算過(guò)程結(jié)果=nn");printf("ttN=%dnttL0=%fn",N,L0);B1=AngleToRadian(B);L0=AngleToRadian(L0);L1=AngleToRadian(L);l=L1-L0;GSZS(a,b,l,B1,pointer_x,pointer_y);printf("nt=計(jì)算結(jié)果為

15、=n");printf("nttx=%-20.3lftty=%-20.3lfn",x,(y+500000+N*1000000);printf("nt=n");break;default:printf("您選擇了錯(cuò)誤的數(shù)字n");printf("是否繼續(xù)Y/N:");scanf("%s",&YN);if(YN='N'|YN='n')printf("歡迎您下次使用此程序,如有疑問(wèn)可發(fā)送Email至gys_long126.,我們將盡快為您解答

16、!n");break;else if(XZZF='F'|XZZF='f')printf("請(qǐng)選擇橢球：n1.克拉索夫斯基橢球; 2.1975年國(guó)際橢球體; 3.WGS-84橢球體; 4.自定義橢球n");printf("請(qǐng)選擇數(shù)字后回車(chē)：");scanf("%d",&num);/printf("請(qǐng)輸入中央子午線(xiàn)L0：");/scanf("%d",&N);switch(num)case 1:a=6378245.0000000000;b=63

17、56863.0187730473;printf("注意：您選擇的是克拉索夫斯基橢球!n請(qǐng)輸入x坐標(biāo)和y坐標(biāo):(輸入時(shí)x坐標(biāo)和y坐標(biāo)之間用逗號(hào)隔開(kāi))n");scanf("%lf,%lf",&x,&y);printf("t=中間計(jì)算過(guò)程結(jié)果=nn");FSB=Bf(a,b,x);printf("ttBf=%-20.10lfn",FSB); /DH=floor(y/1000000); /y=y-DH*1000000-500000;/GSFS(a,b,FSB,y,pointer_B,pointer_L);F

18、SB=RadianToAngle(FSB);FSL=RadianToAngle(FSL);L0=(XZSL='L')?(DH*6-3):(DH*3);printf("ttDH=%-20.0lfnttL0=%-20.0lfnttl=%-20.10lfn",DH,L0,FSL);/printf("nt=計(jì)算結(jié)果=n");printf("nttB=%-20.10lfttL=%-20.10lfn",FSB,L0+FSL);printf("nt=n");break;case 2: a=6378140.0000

19、000000;b=6356755.2881575287;printf("注意：您選擇的是1975年國(guó)際橢球體!n請(qǐng)輸入x坐標(biāo)和y坐標(biāo):(輸入時(shí)x坐標(biāo)和y坐標(biāo)之間用逗號(hào)隔開(kāi))n");scanf("%lf,%lf",&x,&y);printf("t=中間計(jì)算過(guò)程結(jié)果=nn");FSB=Bf(a,b,x);printf("ttBf=%-20.10lfn",FSB); /DH=floor(y/1000000); /y=y-DH*1000000-500000;/GSFS(a,b,FSB,y,pointer_B

20、,pointer_L);FSB=RadianToAngle(FSB);FSL=RadianToAngle(FSL);L0=(XZSL='L')?(DH*6-3):(DH*3);/printf("ttDH=%-20.0lfnttL0=%-20.0lfnttl=%-20.10lfn",DH,L0,FSL);/printf("nt=計(jì)算結(jié)果=n");printf("nttB=%-20.10lfttL=%-20.10lfn",FSB,L0+FSL);printf("nt=n");break;case 3:a

21、=6378137.0000000000;b=6356752.3142;printf("注意：您選擇的是WGS-84橢球體!n請(qǐng)輸入x坐標(biāo)和y坐標(biāo):(輸入時(shí)x坐標(biāo)和y坐標(biāo)之間用逗號(hào)隔開(kāi))n");scanf("%lf,%lf",&x,&y);printf("t=中間計(jì)算過(guò)程結(jié)果=nn");FSB=Bf(a,b,x);printf("ttBf=%-20.10lfn",FSB); /DH=floor(y/1000000); /y=y-DH*1000000-500000;/GSFS(a,b,FSB,y,poi

22、nter_B,pointer_L);FSB=RadianToAngle(FSB);FSL=RadianToAngle(FSL);L0=(XZSL='L')?(DH*6-3):(DH*3);/printf("ttDH=%-20.0lfnttL0=%-20.0lfnttl=%-20.10lfn",DH,L0,FSL);/printf("nt=計(jì)算結(jié)果=n");printf("nttB=%-20.10lfttL=%-20.10lfn",FSB,L0+FSL);printf("nt=n");break;ca

23、se 4: printf("注意：您選擇的是自定義橢球體!n請(qǐng)輸入長(zhǎng)半軸a、短半軸b、x坐標(biāo)和y坐標(biāo):(輸入時(shí)數(shù)據(jù)之間用逗號(hào)隔開(kāi))n");scanf("%lf,%lf,%lf,%lf",&a,&b,&x,&y);printf("t=中間計(jì)算過(guò)程結(jié)果=nn");FSB=Bf(a,b,x);printf("ttBf=%-20.10lfn",FSB); DH=floor(y/1000000); y=y-DH*1000000-500000;GSFS(a,b,FSB,y,pointer_B,p

24、ointer_L);FSB=RadianToAngle(FSB);FSL=RadianToAngle(FSL);L0=(XZSL='L')?(DH*6-3):(DH*3);printf("ttDH=%-20.0lfnttL0=%-20.0lfnttl=%-20.10lfn",DH,L0,FSL);/printf("nt=計(jì)算結(jié)果=n");printf("nttB=%-20.10lfttL=%-20.10lfn",FSB,L0+FSL);printf("nt=n");break;default:pri

25、ntf("您選擇了錯(cuò)誤的數(shù)字n");printf("是否繼續(xù)Y/N:");scanf("%s",&YN);if(YN='N'|YN='n')printf("歡迎您下次使用此程序,如有疑問(wèn)可發(fā)送Email至gys_long126.,我們將盡快為您解答!n");break;else printf("警告：您選擇錯(cuò)誤!n");printf("是否繼續(xù)Y/N:");scanf("%s",&YN);while(YN=&

26、#39;Y'|YN='y'); long double AngleToRadian(long double alfa) long double alfa1,alfa2;long double HS;alfa = alfa + 0.001; alfa1=floor(alfa)+floor(alfa-floor(alfa)*100.)/60;alfa2=(alfa*100.-floor(alfa*100.0)/36;alfa1+=alfa2;HS=alfa1/180.*PI;return(HS);long double RadianToAngle(long double a

27、lfa)long double alfa1,alfa2;alfa=alfa*180./PI;alfa1=floor(alfa)+floor(alfa-floor(alfa)*60.)/100.;alfa2=(alfa*60.-floor(alfa*60.)*0.006;alfa1+=alfa2;return(alfa1);void GSZS(long double a,long double b,long double l,long double B,long double *x,long double *y) long double e,e1,t,u2; long double m0,m2,

28、m4,m6,m8,a0,a2,a4,a6,a8,X; long double n0,n2,n4,n6,n8,sin2B,sin4B,sin6B,sin8B,N; e=(sqrt(a*a-b*b)/a; e1=(sqrt(a*a-b*b)/b; t=tan(B);u2=e1*e1*cos(B)*cos(B);m0=a*(1-e*e);m2=3.0/2.0*e*e*m0;m4=5.0/4.0*e*e*m2;m6=7.0/6.0*e*e*m4;m8=9.0/8.0*e*e*m6;a0=m0+m2/2.0+3.0/8.0*m4+5.0/16.0*m6+35.0/128.0*m8;a2=m2/2.0+m

29、4/2.0+15.0/32.0*m6+7.0/16.0*m8;a4=m4/8.0+3.0/16.0*m6+7.0/32.0*m8;a6=m6/32.0+m8/16.0;a8=m8/128.0;X=a0*B-a2/2.0*sin(2*B)+a4/4.0*sin(4*B)-a6/6.0*sin(6*B)+a8/8.0*sin(8*B);printf("ttX=%20.10lfn",X);n0=a;n2=1.0/2.0*e*e*n0;n4=3.0/4.0*e*e*n2;n6=5.0/6.0*e*e*n4;n8=7.0/8.0*e*e*n6;sin2B=pow(sin(B),2);

30、sin4B=pow(sin(B),4);sin6B=pow(sin(B),6);sin8B=pow(sin(B),8);N=n0+n2*sin2B+n4*sin4B+n6*sin6B+n8*sin8B;printf("ttN=%20.10lfn",N);*x=X+N/2.0*t*cos(B)*cos(B)*l*l+N/24.0*t*(5-t*t+9*u2+4*pow(u2,2)*pow(cos(B),4)*pow(l,4)+N/720.0*t*(61-58*t*t+pow(t,4)*pow(cos(B),6)*pow(l,6); *y=N*cos(B)*l+N/6.0*(1

31、-t*t+u2)*pow(cos(B),3)*pow(l,3)+N/120.0*(5-18*t*t+pow(t,4)+14*u2-58*u2*t*t)*pow(cos(B),5)*pow(l,5);long double Bf(long double a,long double b,long double x)long double e,m0,m2,m4,m6,m8; long double a0,a2,a4,a6,a8;long double B,FB,Bn1,Bn;int n=0;e=(sqrt(a*a-b*b)/a; m0=a*(1-e*e);m2=3.0/2.0*e*e*m0;m4=5.0/4.0*e*e*m2;m6=7.0/6.0*e*e*m4;m8=9.0/8.0*e*e*m6;a0=m0+m2/2.0+3.0/8.0*m4+5.0/16.0*m6+35.0/128.0*m8;a2=m2/2.0+m4/2.0+15.0/32.0*m6+7.0/16.0*m8;a4=m