高斯正反演算法程序

㈠ 高斯投影正反算 編寫c++程序 急用，按圖編寫只要正反算。真急用，求幫助。

/*
方程組維數:3
控制精度:0.0000006
系數矩陣
21-1
3-11
13-1
系數矩陣輸入完畢!
常數項矩陣
-162
常數項矩陣輸入完畢!
整理完畢，此時方程組為
1-0.3333330.3333332
01-0.40
0015

㈡ 高斯反算用c語言表示，可是找不到錯誤

整形除法結果是整形，如果不能被整除的話精度會丟失，結果出現問題，後面那個測試用例恰好沒有出現，解決方法：
把t定義和數組a定義成double不會出現該問題，演算法也可以優化下，我的實現
#include "stdafx.h"
#include#include
double a[3][4] ={123,32,12,2, 3,1,2,3, 4,5,3,3};void p(double a[3][4])
{ for(int i=0;i<3;i++)
{
for(int j=0;j<4;j++)
printf("%f\t",a[i][j]);
printf("\n");
}
}
int main()
{
int i,j,k;
double t = 1;

p(a);
for (i = 0;i < 2; i++)
{
for (j = i + 1; j < 3; j++)
{
t = a[j][i]/a[i][i];
cout<<t<<endl;
for (k = i; k < 4; k++)
{
a[j][k] = a[j][k] - a[i][k]*t;
}

}
}

p(a);
return 0;}

㈢ 高斯投影正『反算公式間接進行換帶計算的實質是什麼

高斯投影換帶計算利用了高斯投影的正、反算公式，其實質是：將橢球面上的大地坐標（B，L）作為過渡坐標。
將某投影帶（如第1帶，中央子午線為L1）內有關點的平面坐標（x1，y1）1，利用高斯投影反算公式換算成橢球面上的大地坐標（B，l1），進而得到L＝L1+l1；根據鄰帶（如第2帶，中央子午線為L2）的中央子午線來計算經差l2,即，l2＝L-L2；然後再由大地坐標（B，l2）利用高斯投影正算公式得到鄰帶的平面坐標（x2，y2）。

㈣ 求高斯公式計算的詳細步驟

過程如圖，先構建一個封閉曲面，使用高斯公式算出來等於4/3，然後算加上去的曲面積分和是10/3。最後答案是4/3－10/3=－2

這個高斯公式計算的話，你可以參照課本上的一些公式和例題進行模仿。

㈤ 80坐標高斯投影正反算，求公式，並註明其中各個字母的含義，務必用最通俗易的方式，望大神解答！

問題有點復雜，你參考下：http://wenku..com/link?url=-Mw6DgC1jqb1F9YrzFGsSgoM-

㈥ 高斯正算公式

弧長參數計算:
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

㈦ 高斯投影正反算

//高斯投影正、反算
//////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表示

㈧ 求高斯投影正反算公式以及 公式中每個字母所包含的內容 謝謝了！

這個

㈨ 誰有高斯投影正反算程序

我有。 要的話GALL我。

㈩ 我要用VB編一個能進行高斯投影正反算和鄰帶換算的程序 哪個高手幫幫忙 萬分感謝

我有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
另外,站長團上有產品團購,便宜有保證