A. 高斯投影正反算 编写c++程序 急用,按图编写只要正反算。真急用,求帮助。
/*
方程组维数:3
控制精度:0.0000006
系数矩阵
21-1
3-11
13-1
系数矩阵输入完毕!
常数项矩阵
-162
常数项矩阵输入完毕!
整理完毕,此时方程组为
1-0.3333330.3333332
01-0.40
0015
B. 高斯算法是什么
一次数学课上,老师让学生练习算数。于是让他们一个小时内算出1+2+3+4+5+6+……+100的得数。全班只有高斯用了不到20分钟给出了答案,因为他想到了用(1+100)+(2+99)+(3+98)……+(50+51)…………一共有50个101,所以50×101就是1加到一百的得数。后来人们把这种简便算法称作高斯算法。
具体的方法是:
首项加末项乘以项数除以2
项数的计算方法是末项减去首项除以项差(每两项之间的差)加1.
1+2+3+4+5+······+n
字母表示:n(1+n)/2
等差数列求和公式Sn=(a1+an)n/2Sn=n(2a1+(n-1)d)/2; d=公差Sn=An2+Bn; A=d/2,B=a1-(d/2)
C. 高斯正算公式
弧长参数计算:
i_a = x1
i_a_ = x2
//if i_a = 6378245.000 and i_a_ = 298.3 then
// i_a0 = 111134.8611
// i_a1 = 32005.7799
// i_a2 = 133.9238
// i_a3 = 0.6973
// i_a4 = 0.0039
//elseif i_a = 6378140.000 and i_a_ = 298.257 then //xi80
// i_a0 = 111133.0047
// i_a1 = 32009.8575
// i_a2 = 133.9602
// i_a3 = 0.6976
// i_a4 = 0.0039
//else
// return -1
//end if
i_b = i_a - i_a/i_a_
i_c = i_a*i_a/i_b
i_e1 = (i_a + i_b)*(i_a - i_b)/i_a/i_a
i_e2 = (i_a + i_b)*(i_a - i_b)/i_b/i_b
i_a0 = 1.0 + (3.0/4.0 + (45.0/64.0 + ( 525.0/768.0 + (33075.0/49152.0 + ( 654885.0/983040.0)*i_e1)*i_e1)*i_e1)*i_e1)*i_e1
i_a1 = (3.0/4.0 + (15.0/16.0 + (1575.0/1536.0 + (6615.0/6144.0 + (1091475.0/983040.0)*i_e1)*i_e1)*i_e1)*i_e1)*i_e1
i_a2 = ((15.0/64.0 + ( 315.0/768.0 + (6615.0/12288.0 + ( 155925.0/245760.0)*i_e1)*i_e1)*i_e1)*i_e1)*i_e1
i_a3 = (((105.0/1536.0 + (945.0/6144.0 + (467775.0/1966080.0)*i_e1)*i_e1)*i_e1)*i_e1)*i_e1
i_a4 = ((((945.0/49152.0 + ( 51975.0/983040.0)*i_e1)*i_e1)*i_e1)*i_e1)*i_e1
i_a5 = ((((( 10395.0/1966080.0)*i_e1)*i_e1)*i_e1)*i_e1)*i_e1
double d
d = i_a*(1.0 - i_e1)
i_a0 = i_a0*d
i_a1 = -i_a1*d/2.0
i_a2 = i_a2*d/4.0
i_a3 = -i_a3*d/6.0
i_a4 = i_a4*d/8.0
i_a5 = -i_a5*d/10.0
---------------------------------------------------------------------------------------------
计算弧长:
i_a0* x + &
i_a1* sin( 2*x) + &
i_a2* sin( 4*x) + &
i_a3* sin( 6*x) + &
i_a4* sin( 8*x) + &
i_a5* sin(10*x)
----------------------------------------------------------------------------
正算公式:
GX_L0 = L0
GX_L = L - L0
GX_B = B
i_tanB = tan(gx_B)
i_tanBB = i_tanB*i_tanB
i_cosB = cos(gx_B)
i_cosBBLL = i_cosB*i_cosB*GX_L*GX_L
i_nn = i_e2*i_cosB*i_cosB
i_N = i_c/sqrt(1.0 + i_nn)
GX_X =getarc(GX_B) + i_N*i_tanB*i_cosBBLL*(0.5 + &
i_cosBBLL*((5.0 - i_tanBB +9.0*i_nn + 4.0*i_nn*i_nn)/24.0 + &
i_cosBBLL*((61.0 -58.0*i_tanBB + i_tanBB*i_tanBB)/720.0)))
GX_Y = i_N*i_cosB*GX_L*(1.0 + i_cosBBLL*((i_nn - i_tanBB + 1.0)/6.0+ &
i_cosBBLL*((5.0 + (i_tanBB - 18.0)*i_tanBB+ &
(14.0 -58.0*i_tanBB)*i_nn)/120.0)))
http://www.gisforum.net/bbs/dispbbs.asp?BoardID=44&ID=55396
D. 80坐标高斯投影正反算,求公式,并注明其中各个字母的含义,务必用最通俗易的方式,望大神解答!
问题有点复杂,你参考下:http://wenku..com/link?url=-Mw6DgC1jqb1F9YrzFGsSgoM-
E. 高斯投影正反算
//高斯投影正、反算
//////6度带宽 54年北京坐标系
//高斯投影由经纬度(Unit:DD)反算大地坐标(含带号,Unit:Metres)
void GaussProjCal(double longitude, double latitude, double *X, double *Y)
{
int ProjNo=0; int ZoneWide; ////带宽
double longitude1,latitude1, longitude0,latitude0, X0,Y0, xval,yval;
double a,f, e2,ee, NN, T,C,A, M, iPI;
iPI = 0.0174532925199433; ////3.1415926535898/180.0;
ZoneWide = 6; ////6度带宽
a=6378245.0; f=1.0/298.3; //54年北京坐标系参数
////a=6378140.0; f=1/298.257; //80年西安坐标系参数
ProjNo = (int)(longitude / ZoneWide) ;
longitude0 = ProjNo * ZoneWide + ZoneWide / 2;
longitude0 = longitude0 * iPI ;
latitude0=0;
longitude1 = longitude * iPI ; //经度转换为弧度
latitude1 = latitude * iPI ; //纬度转换为弧度
e2=2*f-f*f;
ee=e2*(1.0-e2);
NN=a/sqrt(1.0-e2*sin(latitude1)*sin(latitude1));
T=tan(latitude1)*tan(latitude1);
C=ee*cos(latitude1)*cos(latitude1);
A=(longitude1-longitude0)*cos(latitude1);
M=a*((1-e2/4-3*e2*e2/64-5*e2*e2*e2/256)*latitude1-(3*e2/8+3*e2*e2/32+45*e2*e2
*e2/1024)*sin(2*latitude1)
+(15*e2*e2/256+45*e2*e2*e2/1024)*sin(4*latitude1)-(35*e2*e2*e2/3072)*sin(6*l
atitude1));
xval = NN*(A+(1-T+C)*A*A*A/6+(5-18*T+T*T+72*C-58*ee)*A*A*A*A*A/120);
yval = M+NN*tan(latitude1)*(A*A/2+(5-T+9*C+4*C*C)*A*A*A*A/24
+(61-58*T+T*T+600*C-330*ee)*A*A*A*A*A*A/720);
X0 = 1000000L*(ProjNo+1)+500000L;
Y0 = 0;
xval = xval+X0; yval = yval+Y0;
*X = xval;
*Y = yval;
}
//高斯投影由大地坐标(Unit:Metres)反算经纬度(Unit:DD)
void GaussProjInvCal(double X, double Y, double *longitude, double *latitude) 字串9
{
int ProjNo; int ZoneWide; ////带宽
double longitude1,latitude1, longitude0,latitude0, X0,Y0, xval,yval;
double e1,e2,f,a, ee, NN, T,C, M, D,R,u,fai, iPI;
iPI = 0.0174532925199433; ////3.1415926535898/180.0;
a = 6378245.0; f = 1.0/298.3; //54年北京坐标系参数
////a=6378140.0; f=1/298.257; //80年西安坐标系参数
ZoneWide = 6; ////6度带宽
ProjNo = (int)(X/1000000L) ; //查找带号
longitude0 = (ProjNo-1) * ZoneWide + ZoneWide / 2;
longitude0 = longitude0 * iPI ; //中央经线
X0 = ProjNo*1000000L+500000L;
Y0 = 0;
xval = X-X0; yval = Y-Y0; //带内大地坐标
e2 = 2*f-f*f;
e1 = (1.0-sqrt(1-e2))/(1.0+sqrt(1-e2));
ee = e2/(1-e2);
M = yval;
u = M/(a*(1-e2/4-3*e2*e2/64-5*e2*e2*e2/256));
fai = u+(3*e1/2-27*e1*e1*e1/32)*sin(2*u)+(21*e1*e1/16-55*e1*e1*e1*e1/32)*sin(
4*u)
+(151*e1*e1*e1/96)*sin(6*u)+(1097*e1*e1*e1*e1/512)*sin(8*u);
C = ee*cos(fai)*cos(fai);
T = tan(fai)*tan(fai);
NN = a/sqrt(1.0-e2*sin(fai)*sin(fai)); 字串1
R = a*(1-e2)/sqrt((1-e2*sin(fai)*sin(fai))*(1-e2*sin(fai)*sin(fai))*(1-e2*sin
(fai)*sin(fai)));
D = xval/NN;
//计算经度(Longitude) 纬度(Latitude)
longitude1 = longitude0+(D-(1+2*T+C)*D*D*D/6+(5-2*C+28*T-3*C*C+8*ee+24*T*T)*D
*D*D*D*D/120)/cos(fai);
latitude1 = fai -(NN*tan(fai)/R)*(D*D/2-(5+3*T+10*C-4*C*C-9*ee)*D*D*D*D/24
+(61+90*T+298*C+45*T*T-256*ee-3*C*C)*D*D*D*D*D*D/720);
//转换为度 DD
*longitude = longitude1 / iPI;
*latitude = latitude1 / iPI;
}
NN卯酉圈曲率半径,测量学里面用N表示
M为子午线弧长,测量学里用大X表示 字串2
fai为底点纬度,由子午弧长反算公式得到,测量学里用Bf表示 字串4
R为底点所对的曲率半径,测量学里用Nf表示
F. 高斯算法怎么算
高斯算法怎么算
以首项加末项乘以项数除以2用来计算“1+2+3+4+5+···+(n-1)+n”的结果。 这样的算法被称为高斯算法。
G. 高斯投影正‘反算公式间接进行换带计算的实质是什么
高斯投影换带计算利用了高斯投影的正、反算公式,其实质是:将椭球面上的大地坐标(B,L)作为过渡坐标。
将某投影带(如第1带,中央子午线为L1)内有关点的平面坐标(x1,y1)1,利用高斯投影反算公式换算成椭球面上的大地坐标(B,l1),进而得到L=L1+l1;根据邻带(如第2带,中央子午线为L2)的中央子午线来计算经差l2,即,l2=L-L2;然后再由大地坐标(B,l2)利用高斯投影正算公式得到邻带的平面坐标(x2,y2)。
H. 求高斯投影正反算公式以及 公式中每个字母所包含的内容 谢谢了!
这个
I. 《大地测量学基础》里面的高斯投影正反算公式及换带计算的VB或者C语言编程
我有VB的,自己很多年前写的,一直用,但是正算->反算->正算后,Y坐标与原来的差了0.5-0.7mm,不知道怎么回事,这两年工作忙也没有时间再深究,但是这样的计算精度做控制足够了,如果楼主或是者是哪位同仁见此贴能顺便把这个问题解决了,咱们就一起进步了!代码如下:
'高斯坐标正算
Private Sub DadiZs()
Dim t As Double, Itp As Double, X0 As Double, N As Double, L0 As Double
Dim V As Double, ll As Double, W As Double, M As Double
Lat = Radian(Lat)
Lon = Radian(Lon)
L0 = Radian(Lo)
If Tq = 0 Then
a = 6378245 '54椭球参数
b = 6356863.01877305
ep = 0.006693421622966
ep1 = 0.006738525414683
f = (a - b) / a
c = a ^ 2 / b
d = b ^ 2 / a
X0 = 111134.8611 * (Lat * 180# / Pi) - (32005.7799 * Sin(Lat) + 133.9238 * (Sin(Lat)) ^ 3 + 0.6973 * (Sin(Lat)) ^ 5 + 0.0039 * (Sin(Lat)) ^ 7) * Cos(Lat)
'X0 = 111134.8611 * (Lat * 180# / Pi) - (32005.7798 * Sin(Lat) + 133.9238 * (Sin(Lat)) ^ 3 + 0.6972 * (Sin(Lat)) ^ 5 + 0.0039 * (Sin(Lat)) ^ 7) * Cos(Lat)
Else
a = 6378140 '75椭球参数
b = 6356755.28815753
ep = 0.006694384999588
ep1 = 0.006739501819473
f = (a - b) / a
c = a ^ 2 / b
d = b ^ 2 / a
X0 = 111133.0047 * (Lat * 180 / Pi) - (32009.8575 * Sin(Lat) + 133.9602 * (Sin(Lat)) ^ 3 + 0.6976 * (Sin(Lat)) ^ 5 + 0.0039 * (Sin(Lat)) ^ 7) * Cos(Lat)
End If
ll = Lon - L0
t = Tan(Lat)
Itp = ep1 * Cos(Lat) ^ 2
W = Sqr(1 - ep * Sin(Lat) ^ 2)
V = Sqr(1 + ep1 * Cos(Lat) ^ 2)
M = c / V ^ 3
N = a / W
'x = X0 + N * t * (Cos(Lat)) ^ 2 * ll ^ 2 / 2 + N * t * (5 - t * t + 9 * Itp + 4 * Itp * Itp) * (Cos(Lat)) ^ 4 * ll ^ 4 / 24 + N * t * (61 - 58 * t ^ 2 + t ^ 4 + 270 * Itp - 330 * t ^ 2 * Itp) * (Cos(Lat)) ^ 6 * ll ^ 6 / 720 + N * t * (1385 - 3111 * t ^ 2 + 543 * t ^ 4 - t ^ 6) * Cos(Lat) ^ 8 * ll ^ 8 / 40320
x = X0 + N * t * (Cos(Lat)) ^ 2 * ll ^ 2 / 2 + N * t * (5 - t * t + 9 * Itp ^ 2 + 4 * Itp ^ 4) * (Cos(Lat)) ^ 4 * ll ^ 4 / 24 + N * t * (61 - 58 * t ^ 2 + t ^ 4 + 270 * Itp ^ 2 - 330 * t ^ 2 * Itp ^ 2) * (Cos(Lat)) ^ 6 * ll ^ 6 / 720 + N * t * (1385 - 3111 * t ^ 2 + 543 * t ^ 4 - t ^ 6) * Cos(Lat) ^ 8 * ll ^ 8 / 40320
y = N * Cos(Lat) * ll + N * (1 - t * t + Itp) * (Cos(Lat)) ^ 3 * ll ^ 3 / 6 + N * (5 - 18 * t * t + t ^ 4 + 14 * Itp - 58 * Itp * t * t) * (Cos(Lat)) ^ 5 * ll ^ 5 / 120 + N * (61 - 479 * t ^ 2 + 179 * t ^ 4 - t ^ 6) * Cos(Lat) ^ 7 * ll ^ 7 / 5040
r = Sin(Lat) * ll + Sin(Lat) * (Cos(Lat)) ^ 2 * ll ^ 3 * (1 + 3 * Itp + 2 * Itp ^ 2) / 3 + Sin(Lat) * (Cos(Lat)) ^ 4 * ll ^ 5 * (2 - t * t) / 15
r = Degree(r)
y = y + 500000#
End Sub
'高斯反算
Private Sub DadiFs()
Dim t As Double, Itp As Double, X0 As Double, Bf As Double, N As Double
Dim v As Double, ll As Double, W As Double, M As Double, L0 As Double
L0 = Radian(Lo)
X0 = x * 0.000001
y = y - 500000#
If Tq = 0 Then
a = 6378245 '54椭球参数
b = 6356863.01877305
ep = 0.006693421622966
ep1 = 0.006738525414683
f = (a - b) / a
c = a ^ 2 / b
d = b ^ 2 / a
If X0 < 3 Then
Bf = 9.04353301294 * X0 - 0.00000049604 * X0 ^ 2 - 0.00075310733 * X0 ^ 3 - 0.00000084307 * X0 ^ 4 - 0.00000426055 * X0 ^ 5 - 0.00000010148 * X0 ^ 6
ElseIf X0 < 6 Then
Bf = 27.11115372595 + 9.02468257083 * (X0 - 3) - 0.00579740442 * (X0 - 3) ^ 2 - 0.00043532572 * (X0 - 3) ^ 3 + 0.00004857285 * (X0 - 3) ^ 4 + 0.00000215727 * (X0 - 3) ^ 5 - 0.00000019399 * (X0 - 3) ^ 6
End If
Else
a = 6378140 '75椭球参数
b = 6356755.28815753
ep = 0.006694384999588
ep1 = 0.006739501819473
f = (a - b) / a
c = a ^ 2 / b
d = b ^ 2 / a
If X0 < 3 Then
Bf = 9.04369066313 * X0 - 0.00000049618 * X0 ^ 2 - 0.00075325505 * X0 ^ 3 - 0.0000008433 * X0 ^ 4 - 0.00000426157 * X0 ^ 5 - 0.0000001015 * X0 ^ 6
ElseIf X0 < 6 Then
Bf = 27.11162289465 + 9.02483657729 * (X0 - 3) - 0.00579850656 * (X0 - 3) ^ 2 - 0.00043540029 * (X0 - 3) ^ 3 + 0.00004858357 * (X0 - 3) ^ 4 + 0.00000215769 * (X0 - 3) ^ 5 - 0.00000019404 * (X0 - 3) ^ 6
End If
End If
Bf = Bf * Pi / 180#
t = Tan(Bf)
Itp = ep1 * Cos(Bf) ^ 2
W = Sqr(1 - ep * Sin(Bf) ^ 2)
v = Sqr(1 + ep1 * Cos(Bf) ^ 2)
M = c / v ^ 3
N = a / W
Lat = Bf - 0.5 * v ^ 2 * t * ((y / N) ^ 2 - (5 + 3 * t * t + Itp - 9 * Itp * t * t) * (y / N) ^ 4 / 12 + (61 + 90 * t * t + 45 * t ^ 4) * (y / N) ^ 6 / 360)
ll = ((y / N) - (1 + 2 * t * t + Itp) * (y / N) ^ 3 / 6 + (5 + 28 * t * t + 24 * t ^ 4 + 6 * Itp + 8 * Itp * t * t) * (y / N) ^ 5 / 120) / Cos(Bf)
r = y * t / N - y ^ 3 * t * (1 + t * t - Itp) / (3 * N ^ 3) + y ^ 5 * t * (2 + 5 * t * t + 3 * t ^ 4) / (15 * N ^ 5)
Lat = Degree(Lat)
Lon = Degree(L0 + ll)
r = Degree(r)
End Sub
有了正反算,换带也就完成了!
用到的子程序:
Public Const Pi = 3.14159265358979, p = 206264.806
Public Cktq As String
'角度化弧度
Public Function Radian(a As Double) As Double
Dim Ro As Double
Dim c As Double
Dim Fs As Double
Dim Ib As Integer
Dim Ic As Integer
If a < 0 Then a = -a: t = 1
Ro = Pi / 180#
Ib = Int(a)
c = (a - Ib) * 100#
Ic = Int(c + 0.000000000001)
Fs = (c - Ic) * 100#
If t = 1 Then Radian = -(Ib + Ic / 60# + Fs / 3600#) * Ro Else Radian = (Ib + Ic / 60# + Fs / 3600#) * Ro
End Function
'弧度化角度
Public Function Degree(a As Double) As Double
Dim Bo As Double
Dim Fs As Double
Dim Im As Integer
Dim Id As Integer
If a < 0 Then a = -a: t = 1
Bo = a
Call DMS(Bo, Id, Im, Fs)
If t = 1 Then Degree = -(Id + Im / 100# + Fs / 10000#) Else Degree = Id + Im / 100# + Fs / 10000#
End Function
Public Sub DMS(a As Double, Id As Integer, Im As Integer, Fs As Double)
Dim Bo As Double
Dim c As Double
c = a
c = 180# / Pi * c
Id = Int(c)
Bo = (c - Id) * 60
Im = Int(Bo)
Fs = (Bo - Im) * 60
End Sub
'取位计算
Public Function Qw(a As Double, Ws As Integer) As Double
Qw = Int(a * 10 ^ Ws + 0.5) / 10 ^ Ws
End Function
J. 哪里能够找到高斯投影反算公式
这是一个高斯反算的算法,你可以参考以下,这个程序是针对克拉索夫斯基椭球的,转成其他的椭球体的坐标,只需要换一下参数就可以了。
#define PI 3.1415926535
main()
{
float P;
float x,y,z;
float a,b,bf,Bf,B,l1,l2,L,L0;
float b2,b3,b4,b5;
float N,i,j,k;
int B1,B2,L1,L2;float B3,L3;
P=(180.0/PI)*3600.0;
printf("please enter x,y:\n");
scanf("%f,%f",&x,&y);
y-=500000;
a=x*P/6367558.4969;
b=(a*PI/180.0)/3600.0;
i=cos(b)*cos(b);
bf=a+(50221746+(293622+(2350+22*i)*i)*i)*sin(b)*cos(b)*P*(1E-10);
Bf=(bf*PI/180.0)/3600.0;
j=cos(Bf)*cos(Bf);
N=6399698.902-(21562.267-(108.973-0.612*j)*j)*j;
b2=(0.5+0.003369*j)*sin(Bf)*cos(Bf);
b3=0.333333-(0.166667-0.001123*j)*j;
b4=0.25+(0.16161+0.00562*j)*j;
b5=0.2-(0.1667-0.0088*j)*j;
z=(y/N)/cos(Bf);
k=z*z;
printf("Please enter L0:\n");
scanf("%f",&L0);
B=bf-(1-(b4-0.12*k)*k)*k*b2*P;
B=B/3600.0;
B1=(int)B;
B2=(int)((B-B1)*60);
B3=((B-B1)*60-(int)((B-B1)*60))*60.;
l1=(1-(b3-b5*k)*k)*z*P;
L=(l1/3600)+L0;
L1=(int)L;
L2=(int)((L-L1)*60);
L3=((L-L1)*60-(int)((L-L1)*60))*60.;
}
}