国务院第二次全国土地调查领导小组办公室关于统一图幅理论面积与图斑椭球面积计算要求的通知-法搜-中国法律信息搜索网

　　B1 ＝ B1 / RHO
　　L1 ＝ L1 / RHO
　　计算梯形面积
　　Double AreaVal；//梯形面积值
　　Double lDiference ；//经差
　　Double bDiference； //纬差
　　Double bSum；//纬度和
　　Double ItemValue（5）；//计算变量
　　bDiference ＝（B1 - B0）；
　　bSum ＝（B1 ＋ B0） / 2；
　　lDiference ＝（L1　＋ L） / 2；
　　ItemValue（0）＝ ParamA * Sin（bDiference / 2） * Cos（bSum）；
　　ItemValue（1）＝ ParamB * Sin（3 * bDiference / 2） * Cos（3 * bSum）；
　　ItemValue（2）＝ ParamC * Sin（5 * bDiference / 2） * Cos（5 * bSum）；
　　ItemValue（3）＝ ParamD * Sin（7 * bDiference / 2） * Cos（7 * bSum）；
　　ItemValue（4）＝ ParamE * Sin（9 * bDiference / 2） * Cos（9 * bSum）；
　　AreaVal ＝ 2 * bRadius * lDiference * bRadius * （ItemValue（0） - ItemValue（1）＋ ItemValue（2） - ItemValue（3）＋ ItemValue（4））；
　　areaSum ＝ areaSum ＋ AreaVal；
　　Next
　　End Sub
　　3、高斯坐标反解算法
　　Public Sub ComputeXYGeo（x As Double， y As Double， B As Double， L As Double， center As Double）
　　Dim y1 As Double
　　Dim bf As Double
　　y1 ＝ y - 500000
　　Dim e As Double
　　e ＝ Parak0 * x
　　Dim se As Double
　　se ＝ Sin（e）
　　bf ＝ e ＋ Cos（e） * （Parak1 * se - Parak2 * Power（se， 3）＋ Parak3 * Power（se， 5） - Parak4 * Power（se， 7））
　　Dim v As Double
　　Dim t As Double
　　Dim N As Double
　　Dim nl As Double
　　Dim vt As Double
　　Dim yn As Double
　　Dim t2 As Double
　　Dim g As Double
　　g ＝ 1
　　t ＝ Tan（bf）
　　nl ＝ ParaE1 * Power（Cos（bf）， 2）
　　v ＝ Sqr（1 ＋ nl）
　　N ＝ ParaC / v
　　yn ＝ y1 / N
　　vt ＝ Power（v， 2） * t
　　t2 ＝ Power（t， 2）
　　B ＝ bf - vt * Power（yn， 2） / 2 ＋（5 ＋ 3 * t2 ＋ nl - 9 * nl * t2） * vt * Power（yn， 4） / 24 - （61 ＋ 90 * t2 ＋ 45 * Power（t2， 2）） * vt * Power（yn， 6） / 720