多元回歸格網改正量模型的大區域坐標轉換方法

楊國創，尹 彤，王金月，高 嚴，孫澤宇

多元回歸格網改正量模型的大區域坐標轉換方法

楊國創1，尹 彤2，王金月3，高 嚴2，孫澤宇2

（1. 廣東省地質測繪院，廣州 510800；2. 自然資源部 測繪標準化研究所，西安 710054；3. 自然資源部 地圖技術審查中心，北京 100036）

為了解決我國存量地理信息成果向2000中國大地坐標系（CGCS2000）的整體無縫轉換問題，提出一種多元回歸格網改正量模型的大區域坐標轉換方法：采用格網改正量法建立統一的轉換模型，格網節點改正量采用多元回歸算法實現，待轉點改正量采用雙線性內插法獲取，研究整套轉換模型及其使用方法。實驗結果表明，多元回歸模型計算精度較高，且可保證轉換結果的唯一性；格網改正量法是大區域轉換的統一模型，能夠實現大區域的無縫轉換；格網轉換方法效率高、可實施性強、準確度高。

2000中國大地坐標系（CGCS2000）；多元回歸；格網改正量；雙線性內插；坐標轉換

0 引言

中國當前最新的國家大地坐標系2000國家大地坐標系（China geodetic coordinate system 2000，CGCS2000）是采用空間技術手段建立的高精度地心動態大地坐標系[1]。2008年7月啟用之后，很多大地測量和測繪工程工作均基于CGCS2000框架展開。

由于之前我國先后采用1954北京坐標系和1980西安坐標系，所以在我國基礎測繪及工程應用領域存在大量往期坐標系統成果資料，這些資料在我國的建設領域大量存在并產生著關鍵作用。CGCS2000的啟用必然存在一個關鍵問題，即往期地理信息成果如何轉換至CGCS2000坐標系下，如何采用一種通用方法完成全國海量的地理信息數據轉換工作，可以考慮將格網改正量法應用到大區域成果坐標系統轉換中。

目前，我國通過與空間大地網聯合平差已獲取了將近5萬多個不同級別天文大地網點的CGCS2000坐標[2]，這些點在我國中東部地區分布較為均勻，西部地區依據地理環境和經濟情況也整體覆蓋；所以基于這些公共點可以實現1954北京坐標系、1980西安坐標系與CGCS2000之間的相互轉換。

在文獻[3-4]中，將反距離加權法和德勞奈（Delaunay）三角網法應用到我國速度場應用中，采用鄰近點綜合特征表述任一點特征并構建三角網的方法，效果很好。本文將重點研究如何利用公共點，采用多元回歸模型，實現基于方形格網改正量的坐標轉換方法，并結合某市的情況進行實驗驗證和精度分析。

1 基于多元回歸模型的格網改正量計算方法

1.1 多元回歸方程

多元回歸方程是含有多個回歸變量的多項式回歸模型，在大地測量領域可用來轉換不同參考框架之間的大地坐標或直角坐標[5]。在格網模型的格網節點改正量計算中可采用多元回歸模型，通過搜索一定范圍內的公共點及其坐標改正值來計算節點處的坐標改正值。其輸入數據可以是搜索半徑內符合條件的點，也可以是全部數據。因輸入數據的數量影響模型計算速度，故可視情況而定。在每個格網節點處均須單獨構建多項式曲面模型，根據每個節點半徑范圍內的所有公共點計算回歸方程參數。

如圖1所示，輸入數據為節點周邊半徑范圍內搜索的公共點，為搜索到的公共點個數，特征值為每個公共點的原坐標系和目標坐標系的經緯度差異（?B,?L）。

圖1 格網節點處半徑搜索（圓點為格網節點，三角為搜索到的公共點）

首先確定多元方程模擬量，計算節點半徑范圍內的每個公共點處2套坐標系統在緯度和經度位置的差異（?B,?L）（0<≤）；然后計算搜索范圍內所有公共點的中心點位置（0,0）；最后根據以上數據在每個節點處構建多元回歸多項式。節點處緯度和經度坐標改正量表達式為：

式中：A（0≤≤）為緯度特征項多項式系數；C為經度特征項多項式系數；U=B-0；V=L-0（B、L為節點在原坐標系的緯度和經度值）；為回歸方程階數。

將每個公共點處的經緯度差（?B,?L）以及U和V帶入多項式中，形成組多項式方程，并采用最小二乘擬合計算A和C。最終將節點的經緯度帶入多項式中，即可獲取節點處的改正量值。0和0即為格網節點的坐標改正量估值[6]。

1.2 格網節點改正量數據存儲

格網改正數據存儲可以采用國際標準的格網文件格式，如：美國或加拿大大地測量局建立的標準格網格式[7-8]；一般研究中可以根據編程需要來設定。本文采用txt文本存儲，存儲格式為行號、列號、緯度改正量?B、經度改正量?L，同時配合元數據文件來說明格網起始經緯度和經緯度方向的格網間隔。

2 雙線性內插算法

雙線性插值是先后在2個方向或者2個變量上進行線性插值的算法。雙線性內插法在矢量和影像轉換中使用普遍，具有計算簡單、轉換結果平滑、過渡自然等特點[9-10]。本文涉及到如何通過格網4個節點經緯度坐標及對應的特征值，采用差值方法來獲取格網內任意位置的特征值，即在經度和緯度方向上依次實施線性插值。如圖2所示。

圖2 格網內任意位置的雙線性內插（圓圈為格網點，五星為內插點）

圖2中，已知格網4個節點經緯度坐標和節點經緯度坐標轉換改正量分別為1（1,1,?1,?1）、2（1,2,?2,?2）、3（2,1,?3,?3）、4（2,2,?4,?4）。內插點為待定點（B,L,?b,?l），其中位置參數（B,L）已知，需要計算特征參數坐標轉換改正量（?b,?l）。

雙線性內插算法首先在方向（東方向）上內插，然后在方向（北方向）內插。即先在緯度方向內插點1,2和3,4，然后由1,2和3,4內插點。具體流程如下（其中所有點的位置坐標已知）：

1）在1方向上進行線性內插，得到1,2（1,L,?1,2,?1,2），即：

2）在2方向上進行線性內插，得到3,4（2,L, ?3,4, ?3,4），即：

3）在B方向上進行線性內插，得到（B,p,?b,?l），即：

在這個模型中，只要輸入區域內的任一點（格網內部、格網線上）的經緯度就可以插值該點的坐標改正量。坐標轉換過程只需要在原坐標上加入改正量既可獲取在目標坐標系的新坐標，從而可實現基于格網的坐標轉換工作。轉換后新的坐標為：

同時，在區域上插值獲取的改正量具有連續性和唯一性特點[11]，這一點十分重要。

3 基于格網的坐標轉換流程

改正量格網的坐標轉換主要包括5個步驟：1）按照合適的經緯度間隔制作格網；2）基于公共點采用半徑搜索法和多元回歸方程計算格網節點處的坐標改正量；3）格網改正量文件制作；4）格網改正量文件精度驗證；5）格網制作完畢后，基于格網文件采用雙線性內插方法實現逐點坐標轉換。

這個轉換模型即顧及了公共點的空間相關性，又保證了區域轉換的連續性，避免了圖形間的接邊誤差，從而實現高精度、無縫隙的大區域坐標轉換，其優勢是小區域四參數、七參數等轉換模型無法相提并論的。

3.1 區域改正量格網制作

格網制作的主要目的是確定格網上各個節點的經緯度值，首先確定格網起始經度和緯度值(0,0)，然后確定格網經緯度間隔（G,G），以及確定各個節點行列號（h,l）。

1）起始經緯度(0,0)的確定依據。依照轉換區域的范圍，一般采用長方形覆蓋，須在所轉換區域的邊界向外延伸5 km。

2）格網間隔確定依據。所述格網點間隔是指建立格網的縱橫相鄰格網的經緯度等差異量，主要依據控制點的平均間隔。一般格網間隔要小于控制點的平均間隔，可以采用1″、1.5″或3″。

區域坐標改正量格網網形如圖3所示。

圖3 區域格網網絡覆蓋（三角代表公共點）

3.2 基于多元回歸計算格網節點改正值

采用多元回歸方程，利用格網點周邊公共點來計算格網節點處坐標改正量，如圖4所示。

首先采用半徑加四象限搜索法找到節點周邊的公共點；然后將搜索到的公共點經緯度和原坐標系和目標坐標系的經緯度差異作為輸入值，采用多元回歸模型計算節點處的改正量。

遍歷所有的格網節點之后，將節點及其改正量保存在文件中，以備坐標轉換使用。可以繪制格網轉換改正量分量圖，來評估改正量模型的連續性和平緩性——從數值圖來觀察其連續性，一般情況下是連續的、平滑的、平緩的[12]。

圖4 基于多元回歸的格網節點改正數計算（圓點為格網點，三角為公共點）

3.3 格網轉換精度驗證

為了保證格網改正量模型的正確性，須對格網文件進行連續性與符合性評估[13]。可以采用內符合精度和外符合精度驗證。

1）內符合精度驗證。為了保證格網改正文件的精度，需要采用內插方法對參與格網計算的公共點進行轉換，并統計其轉換后的精度。

2）外符合精度驗證。采用外業實測的方法采集原坐標和目標系下的全球衛星導航系統（global navigation satellite system，GNSS）接收機測得的坐標，然后采用格網改正模型將原坐標下的GNSS成果轉換至目標系下并與實測目標系成果進行比較。

4 實驗與結果分析

4.1 實驗數據概述

本次實驗數據采用某市C級和D級點作為已知公共點，這些點均具有1980西安坐標系和CGCS2000坐標系成果，共計107個點（包括46個C級點、61個D級點）。公共點的經緯度坐標范圍為34°20'N—36°07' N、108°23'E—109°55'E。

C級點和D級點的點位精度統計如表1所示。

表1 C級和D級點精度統計 mm

由于公共點分布比較均勻，確定本市改正量格網范圍起始位置為108°20' E—109°58' E、34°17' N—36°10' N，格網經緯方向間隔均為1″，格網采用搜索半徑為3″。公共點分布如圖5所示。

圖5 公共點分布（藍色點將參與格網計算、紅色點作為驗證點）

4.2 格網改正量模型計算

根據前面確定的格網參數，計算并存儲起始點位置及各個節點的行列號，保存在txt文件中。然后采用編寫的多元回歸格網算法軟件逐個節點計算。其中選用97個點參與格網制作，其他10個點（C級點3個、D級點7個）作為外符合精度驗證點。

輸入文件為格網行列號文件、公共點文件、計算參數（搜索半徑、起始點坐標等參數設置文件），經過約1.5 h的計算，獲得格網改正量文件coordinate.dll。

4.3 采用雙線性內插方法計算待轉換點的改正量

讀取轉換圖形或者文本的各個要素，并逐個根據待轉換點的原始坐標確定其所在格網及4個節點的特征參數；采用雙線性內插法以及格網中對應4個格網點的格網改正量確定待轉換點的改正量；然后將待轉點原坐標與改正量相加即可得到目標坐標系下的新坐標，并逐點逐要素實現轉換。

4.4 精度分析

分別采用內符合和外符合檢驗的方法進行精度分析。

采用格網改正文件將公共點從原始坐標系轉換至目標坐標系，然后與其目標坐標系的已知坐標進行比較。

內符合驗證采用97個公共點，將其轉換誤差統計結果為：方向最小誤差為1.1 mm，方向最大誤差為9.2 mm，方向中誤差為4.5 mm；方向最小誤差為2.9 mm，方向最大誤差為11.3 mm，方向中誤差為6.2 mm。

同時，對轉換精度誤差區間分布進行統計：方向轉換誤差在0～5 mm區間的有65個點，方向轉換誤差在5～10 mm區間的有32個點，方向轉換誤差大于10 mm的有0個點；方向轉換誤差在0～5 mm區間的有33個點，方向轉換誤差在5～10 mm區間的有62個點，方向轉換誤差大于10 mm的有2個點。

外符合精度驗證基于沒參與格網改正模型計算的3個C級點和7個D級點，采用格網改正文件將公共點從原始坐標系轉換至目標坐標系，然后跟其目標坐標系的已知坐標進行比較。外符合驗證點轉換誤差統計為：

3個C級點外符合轉換誤差統計結果：方向最小誤差為1.9 mm，方向最大誤差為7.4 mm，方向中誤差為4.8 mm；方向最小誤差為4.1 mm，方向最大誤差為5.8 mm，方向中誤差為5.3 mm。

7個D級點外符合轉換誤差統計結果：方向最小誤差為2.4 mm，方向最大誤差為9.4 mm，方向中誤差為6.9 mm；方向最小誤差為5.3 mm，方向最大誤差為11.3 mm，方向中誤差為7.3 mm。

綜合統計10個外符合點的坐標轉換誤差：方向坐標轉換中誤差為4.3mm，方向轉換最小差異為2.0 mm，最大差異為9 mm；方向坐標轉換中誤差為7.1 mm，方向轉換最小差異為4 mm，最大差異為11 mm。符合轉換精度要求。

5 結束語

本文采用多元回歸格網改正量模型進行地理信息數據的逐點無縫轉換，經過理論分析和實驗驗證，得到以下幾點結論：

1）國家坐標系統轉換涉及到全國或各省海量數據的整體統一轉換，格網改正量模型避免了傳統的四參數或七參數轉換的不足，能夠實現轉換前后圖形元素的無縫接邊；在格網模型所涉及的格網特征量計算方法、差值方法、格網間隔等參數固定的情況下，任意點的坐標轉換結果是唯一的，確保了圖件中點、線、面等元素的坐標轉換結果的唯一性和統一性。

2）格網轉換的本質是通過待轉換點周邊的控制點來實現坐標轉換，距離轉換點越近的公共點的特征表達權重越大，這種基于控制點的直接轉換方法轉換準確度和現勢性更高；采用多元回歸方程計算格網節點改正量的方法精度及其運算速度比較適合生產需求，且內外符合精度較好。

3）格網改正量模型可以根據區域進行任意切割使用，為轉換參數的分區管理、加密處理和生產使用提供了方便，可操作性強，是我國開展大區域 CGCS2000坐標系統轉換的實用方法。

[1] 郭充, 呂志平. 基于Bursa模型的格網坐標轉換[C]// 測繪通報. 測繪科學前沿技術論壇摘要集. 北京: 測繪出版社, 2008: 979-983.

[2] 郭充. 面向CGCS2000的格網坐標轉換方法及應用研究[D]. 鄭州: 信息工程大學, 2009.

[3] 任營營, 王虎, 王解先, 等. 基于K-Means++的省內子塊體劃分及中國大陸水平相對運動速度場模型的建立與分析[J]. 地球物理學報, 2020, 63(7): 2516-2533.

[4] 任營營, 王解先, 王虎, 等. 基于局部無縫Delaunay三角網反距離加權法構建中國大陸速度場[J]. 武漢大學學報(信息科學版), 2021, 46(7): 1071-1080.

[5] 許琳琳. 分布式環境下的格網坐標轉換方法與應用[D]. 鄭州: 信息工程大學, 2012.

[6] 王惠文. 多元回歸模型評價、模型分類以及模型預測理論研究及其應用[D]. 北京: 北京航空航天大學, 2013.

[7] 王凱時, 馮楊民, 劉立. 格網法在2000國家大地坐標系基準轉換中的關鍵技術[J]. 測繪地理信息, 2022, 47(2): 15-18.

[8] 蔣光偉, 焦立芬, 田曉靜, 等. 基于三角剖分與方向線拼接的高精度坐標轉換格網模型的建立[J]. 大地測量與地球動力學, 2015, 35(2): 222-224.

[9] 朱毅, 丁云鵬, 李正會, 等. 利用三角剖分內插法精化似大地水準面模型[J]. 測繪科學技術學報, 2014, 31(6): 584-586.

[10] 趙輝, 張顏秋, 樊琦, 等. 基于格網模型的地理信息數據坐標轉換[J]. 測繪與空間地理信息, 2020, 43(10): 201-203.

[11] 王謙, 趙俊三, 王濤, 等. 改進的高精度CGCS2000格網改正量自動計算模型設計與實現[J]. 地質災害與環境保護, 2019, 30(3): 92-95.

[12] 郭充, 呂志平. 基于Bursa模型的格網坐標轉換[C]// 測繪通報. 測繪科學前沿技術論壇摘要集. 北京: 測繪出版社, 2008: 979-983.

[13] 呂志平, 魏子卿, 李軍, 等. CGCS2000高精度坐標轉換格網模型的建立[J]. 測繪學報, 2013, 42(6): 791-797.

Large area coordinate transformation method based on multivariate regression grid correction model

YANG Guochang1, YIN Tong2, WANG Jinyue3, GAO Yan2, SUN Zeyu2

（1. Geology Surveying and Mapping Institute of Guangdong, Guangzhou 510800, China;2. Institute of Surveying and Mapping Standardization, Ministry of Natural Resources, Xi’an 710054, China;3. Supervision Center of Map, Ministry of Natural Resources, Beijing 100036, China）

In order to solve the overall seamless transformation of China's previous geographic information data results to the China geodetic coordinate system 2000 (CGCS2000), the paper proposed a large area coordinate transformation method based on multivariate regression grid correction model: the unified transformation model was established by using the grid correction model, by which the grid node correction was realized by multiple regression algorithm, and the correction of the points to be transformed was obtained by bilinear interpolation, to perform the calculation and application of the whole set of the transformation model. Experimental result showed that the multiple regression model could have high accuracy and ensure the uniqueness of the conversion result; and the grid correction method would be a unified conversion model of large area to help realize seamless conversion of large areas; moreover, the grid conversion method would have high efficiency, strong implementability and high accuracy.

China geodetic coordinate system 2000 (CGCS2000); multiple regression; grid correction; bilinear interpolation; coordinate transformation

楊國創, 尹彤, 王金月, 等. 多元回歸格網改正量模型的大區域坐標轉換方法[J]. 導航定位學報, 2023, 11(3): 185-190.（YANG Guochang, YIN Tong, WANG Jinyue, et al. Large area coordinate transformation method based on multivariate regression grid correction model[J]. Journal of Navigation and Positioning, 2023, 11(3): 185-190.）

10.16547/j.cnki.10-1096.20230325.

P228

A

2095-4999(2023)03-0185-06

2022-05-13

楊國創（1972—），男，廣東陽春人，大學學歷，測繪高級工程師，研究方向為大地測量、攝影測量技術應用及測繪管理。

尹彤（1989—），男，寧夏固原人，碩士，工程師，研究方向為定位與導航技術應用及測繪地理信息標準化。