獨立坐標系高斯投影實用算法程序實現及應用

姜韶,張坤先,匡志威

(1.四川省冶金地質勘查局測繪工程大隊,四川 成都 610212; 2.長沙市規劃勘測設計研究院,湖南 長沙 410007)

0 引 言

目前GNSS技術已經基本取代了傳統大地測量技術,在城市測量、工程測量中已經成為主要的控制網測量技術。GNSS技術直接確定的是大地坐標系坐標(B、L、H),在實際應用中需要通過高斯投影正算計算出平面直角坐標系坐標。根據高斯投影邊長改化公式,離中央子午線的距離、高程均會引起長度變形,因此在海拔較高的城市、地區以及線性工程的控制網測量中,均應考慮選擇合適的中央子午線、相對參考橢球面來限制投影變形,精度要求應滿足長度變形值不應大于 25 mm/km的要求。

隨著無人機航空攝影測量技術的廣泛應用,利用無人機機載RTK實時獲取航攝相機曝光瞬間的高精度坐標信息,加上無人機IMU慣性測量單元的高度集成化、高性能化、組合化發展,現階段無人機在小于1∶1 000比例尺的地形圖測量中可以完全實現免像控。但是POS數據中的坐標一般都是地理坐標和橢球高,在海拔較高城市、地區以及線性工程為了限制變形,仍需要在空三前將坐標系統轉換到獨立坐標系和1985國家高程基準。少數的商業軟件雖然具備獨立坐標系投影計算功能,但是對數據格式有固定要求,往往需要單獨將機載POS數據中的大地坐標數據提取,轉換到獨立坐標系后再替換到原POS數據中進行空三處理,數據格式的往返和人工編輯替換比較費事,不方便且容易出錯。因此,經過程序開發實現高精度無人機機載POS數據的坐標批量自動轉換可有效減少工作環節,有助于提高工作效率。

1 橢球膨脹法

橢球膨脹法是保持參考橢球的定位、定向和扁率不變,對橢球進行縮放,使得縮放后的參考橢球面與獨立坐標系所選定的平面相切。因此,采用橢球膨脹法投影后得到的平面坐標在數值上與國家參考橢球的橢球面上的平面坐標接近,只是存在縮放關系,不會引起經度(東坐標)的變化,只要選擇合適的投影面,就可以把高程所產生的變形降到最小,繼而忽略不計。因此在建立獨立坐標系時,通常采用橢球膨脹法。

1.1 算法流程

橢球膨脹法關鍵是計算膨脹后的新橢球的長半軸改正值dα,以及緯度改正值dB,具體的算法流程如圖1所示。

圖1 橢球膨脹法算法流程

1.2新橢球的長半軸變化量計算公式

要滿足橢球膨脹法的縮放條件,要求項目區的平均曲率半徑R變化dR=H,其中H為區域的大地高。根據平均曲率半徑計算公式:

(1)

(2)

(3)

將式(2)、式(3)代入(1)求微分得到:

(4)

式中:M為子午圈曲率半徑,N為卯酉圈曲率半徑,B為大地緯度,e為橢球的第一偏心率,α為原橢球長半軸,dα為縮放后橢球相對于原橢球的長半軸改正數,H為大地高。

1.3大地坐標改正公式

從國家參考橢球到縮放后的參考橢球的坐標由莫洛金斯基公式可得:

L=L0

(5)

(6)

上式中,L、B為縮放后的參考橢球對應的大地坐標,L0、B0為原國家參考橢球對應的大地坐標,H為大地高。

2 高斯投影通用算法

高斯投影正反算算法在很多文獻中均給出了詳細公式,在本文中將采用適用于不同橢球的高斯平面坐標正算的實用公式[3],本算法可以在編程語言或者Excel表格中簡單實現。

(7)

(8)

X=c0B-cosB(c1sinB+c2sin3B+c3sin5B)

(9)

式中,c0=Aa(1-e2)/ρ;c1=(C-B-D)a(1-e2)

(10)

3 算法編程實現及精度驗證

3.1 橢球膨脹法任意投影面高斯投影算法實現

基于本算法,在C#語言環境中開發小程序“基于橢球膨脹法任意投影面高斯正反算工具”,小程序界面如圖2所示。選取四川西部某山區城市為例,投影面大地高為 2 600 m,投影中央子午線為102°,城市中心緯度為31°54′,現有GNSS控制點6個(CGCS2000坐標系),利用該程序轉換到獨立坐標系。

圖2 基于橢球膨脹法任意投影面高斯正反算軟件界面

3.2 POS數據坐標自動轉換功能實現

為方便將無人機POS數據中的坐標轉換到工程獨立坐標系和1985國家高程基準,并按照規范要求進行航片編號等,本文利用C#語言實現POS數據坐標轉換和航攝文件整理功能,程序界面如圖3所示。其中平面投影轉換采用上文算法,由于無人機航飛區域一般跨度較小,一般僅僅幾平方公里,因此高程采用GNSS高程平面擬合即可。POS數據高程轉換目標值(正常高)h由下式計算得到:

圖3 無人機POS數據坐標轉換及航攝文件整理工具界面

圖4 HGO內置CoordTool軟件界面

h=H-ζ

(11)

其中:

ζ=a0+a1×(B-B0)+a2×(L-L0)

(12)

式中:B、L為待轉POS坐標點地理坐標;H為待轉POS坐標點橢球高;B0、L0為測區中心地理坐標;a0、a1、a2為轉換系數,由地面控制點計算得到。

3.3 橢球膨脹法任意投影面高斯投影算法精度驗證

(1)HGO、CosaGPS軟件投影計算

為驗證本文實用算法投影后的長度精度,采用相同的參數用CoordTool軟件(HGO軟件內置系統工具)對以上六點進行坐標轉換。打開CoordTool,“參數設置”界面上選擇“源橢球”“當地橢球”均為“國家2000”;“投影”設置界面輸入中央子午線為“102:00:00.00000E”,“投影面高程”為“2600”,“平均緯度”為“031:54:00.00000N”,其他默認;“選項”界面“橢球變形方法”勾選“膨脹”;返回程序主界面,打開源坐標數據進行高斯投影。

本文算法與HGO軟件計算結果如表1所示。

表1 本文算法與HGO軟件計算結果比較

采用相同的參數用在CosaGPS軟件里面設置好坐標系統、中央子午線、測區平均緯度、投影面大地高等參數進行坐標轉換,具體設置如圖5所示。

圖5 CosaGPS軟件坐標轉換設置界面

本文算法與CosaGPS軟件計算結果如表2所示。

表2 本文算法與CosaGPS軟件計算結果比較

由于與HGO、CosaGPS軟件在計算新橢球的長半軸變化量及大地坐標改正值采用的公式、參數精度存在差異,因而計算結果存在細微差別,這種差值在隨著投影面的高程增大而變大,投影點的緯度與中心緯度差值增大則變大。根據表1、表2可知,這種差值對邊長影響較小,一般可以忽略不計。

(2)全站儀測距

為驗證本算法的精度,采用全站儀對C1-C2、C3-C4、C5-C6邊長進行觀測,經過大氣改正的邊長觀測值與本文計算結果反算的邊長值對照如表3所示,可見均明顯優于《工程測量標準》《城市測量規范》中對四等控制網最弱邊相對中誤差的要求。

表3 邊長精度驗算成果

4 結 論

(1)本文采用的獨立坐標系高斯投影算法簡單實用,適用于各種橢球,不僅可以通過計算機編程語言實現,也可以利用Excel表格進行投影計算,投影精度可以滿足一般城市地區、工程項目控制網測量要求。

(2)在進行獨立坐標系高斯投影時,應選擇合適的投影中央子午線、投影高程面、中心緯度。對于線性工程,比如高差較大的線路工程測量,應進行分帶分區投影,確保投影后的邊長變形滿足工程精度要求。

(3)高精度無人機機載RTK技術在線路工程測量中已得到廣泛應用,特別是帶狀地形圖測量基本由無人機航測技術替代,將POS數據中的高精度坐標統一到工程獨立坐標系中,有利于后期的空三加密處理工作,提高工作效率。

(4)本文算法有一定的局限性,單次投影的范圍不宜過大,特別是高差較大區域應分區進行投影計算,注意控制計算范圍。