大地主题解算 正反算matlab程序

此资源内容涉及的是大地测量领域的大地主题解算，是基于白塞尔模型进行的正算与反算。matlab简单易学代码。

%白塞尔大地主题反算
clcclear;%清理命令窗口和工作区间
format long g;%规定输出数据为长整型
a=6378245;b=6356863.0187730473;%克氏椭球
%a=6378140;b=6356755.2881575287;%1975国际椭球
c=a^2/b;e1=sqrt（a^2-b^2）/b;
raw=1;%计算相关椭球参数
B1=（47+46/60+52.6470/3600）*pi/180;
L1=（35+49/60+36.3300/3600）*pi/180;
B2=（48+04/60+09.6384/3600）*pi/180;
L2=（36+14/60+45.0505/3600）*pi/180;%输入已知值（将角度值变成弧度值）
dL=L2-L1;
dB=B2-B1;
Bm=1/2*（B1+B2）;
Vm=sqrt（1+（e1*cos（Bm））^2）;
Nm=c/Vm;
tm=tan（Bm）;
gm=sqrt（e1^2*（cos（Bm））^2）;%计算大地线中点M处相关参数
r01=（Nm/raw）*cos（Bm）;
r21=（Nm*cos（Bm）/（24*raw^3*Vm^4））*（1+gm^2-9*（gm^2）*（tm^2）*（gm^4））;
r03=（-Nm/（24*raw^3））*（cos（Bm）^3）*tm^2;
S10=Nm/（raw*Vm^2）;
S12=（Nm/（24*（raw^3）*（Vm^2）））*（cos（Bm））^2*（+2+3*tm^2+2*gm^2）;
S30=Nm/（8*（raw^3）*（Vm^6））*（gm^2-（tm^2）*（gm^2））;%公式（4-227）
U=r01*dL+r21*dB^2*dL+r03*dL^3;
V=S10*dB+S12*dB*dL^2+S30*dB^3;%公式（4-226）
t01=tm*cos（Bm）;
t21=（1/（24*raw^2*Vm^4））*cos（Bm）*tm*（2+7*gm^2+9*tm^2*gm^2+5*gm^4）;
t03=（1/（24*raw^2））*（（cos（Bm））^3）*tm*（2+tm^2+2*gm^2）;%公式（4-229）
dA=t01*dL+t21*（dB^2）*dL+t03*dL^3;%公式（4-228）
Am=atan（U/V）;%公式（4-230）
S=U/sin（Am）;%计算大地线长度
A1=Am-1/2*dA;
A2=Am+1/2*dA;
if （A2>pi）
    A2=A2-pi;
else
    A2=A2+pi;
end
%计算方位角
SA1A2%输出的A1，A2为弧度值

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件       1293  2016-11-26 15:40  Bessel fansuan.m

     文件       1104  2016-11-26 15:34  Bessel zhengsuan.m

     文件          0  2018-05-03 08:44  CUST.txt

----------- ---------  ---------- -----  ----

                 2397                    3