羅孝楠 周廣勇 謝全遠(yuǎn)

摘要：因各個(gè)國(guó)家衛(wèi)星導(dǎo)航使用的坐標(biāo)系不一致，GPS工程測(cè)量中坐標(biāo)轉(zhuǎn)換是實(shí)際應(yīng)用的重要步驟，針對(duì)目前坐標(biāo)轉(zhuǎn)換軟件數(shù)據(jù)錄入不便、過(guò)程不可控等問(wèn)題，文中探析四參數(shù)法坐標(biāo)轉(zhuǎn)換參數(shù)求解方法，并在Excel下利用少量?jī)?nèi)部函數(shù)，設(shè)計(jì)實(shí)現(xiàn)了四參數(shù)求解和坐標(biāo)轉(zhuǎn)換電子表格。通過(guò)工程實(shí)例坐標(biāo)進(jìn)行驗(yàn)證，求解結(jié)果準(zhǔn)確，計(jì)算過(guò)程數(shù)據(jù)可視，編輯錄入數(shù)據(jù)方便，方便在測(cè)量工程、教學(xué)實(shí)驗(yàn)中使用。

關(guān)鍵詞：衛(wèi)星導(dǎo)航；四參數(shù)；坐標(biāo)轉(zhuǎn)換；Excel

中圖分類(lèi)號(hào)：TP311.1? ? 文獻(xiàn)標(biāo)識(shí)碼：A

文章編號(hào)：1009-3044（2023）35-0118-03

開(kāi)放科學(xué)（資源服務(wù)）標(biāo)識(shí)碼（OSID）

0 引言

伴隨著衛(wèi)星技術(shù)深入發(fā)展，在實(shí)踐中得到了廣泛應(yīng)用，特別是美國(guó)全球定位系統(tǒng)（GPS） 的成功建立和應(yīng)用，我國(guó)北斗衛(wèi)星導(dǎo)航測(cè)量技術(shù)也取得了前所未有的進(jìn)步和發(fā)展。衛(wèi)星定位測(cè)量技術(shù)具有精準(zhǔn)、快速、高效和全天候等顯著的特點(diǎn)，使得測(cè)繪行業(yè)經(jīng)歷了一場(chǎng)深刻的變革。豐富的應(yīng)用場(chǎng)景已經(jīng)被應(yīng)用到了現(xiàn)代社會(huì)經(jīng)濟(jì)中的很多領(lǐng)域，尤其是在軍事、交通、測(cè)量等方面和領(lǐng)域取得了顯著的成效，凸顯了重要的應(yīng)用價(jià)值。在因此全球?qū)Ш蕉ㄎ恍l(wèi)星系統(tǒng)，在全球范圍內(nèi)的眾多行業(yè)中，尤其是在測(cè)繪行業(yè)有著越來(lái)越廣泛的應(yīng)用。

全球定位系統(tǒng)（GPS） 所采用的坐標(biāo)系是WGS-84坐標(biāo)系，因此直接采集的點(diǎn)位坐標(biāo)屬于WGS-84坐標(biāo)系，而目前我們國(guó)家工程測(cè)量成果普遍使用的是以CGCS2000坐標(biāo)系或是地方（任意）獨(dú)立坐標(biāo)系為基礎(chǔ)的坐標(biāo)數(shù)據(jù)。因此必須將WGS-84坐標(biāo)轉(zhuǎn)換到CSCS2000坐標(biāo)系或地方（任意）獨(dú)立坐標(biāo)系[1]。坐標(biāo)轉(zhuǎn)換是從一種坐標(biāo)系統(tǒng)映射到另一種坐標(biāo)系統(tǒng)的過(guò)程。通過(guò)建立兩個(gè)坐標(biāo)系統(tǒng)之間一一對(duì)應(yīng)關(guān)系來(lái)實(shí)現(xiàn)。一般四參數(shù)應(yīng)用于平面轉(zhuǎn)換和七參數(shù)應(yīng)用于三維轉(zhuǎn)換。

在許多測(cè)量數(shù)據(jù)處理軟件中，如徠卡LGO、天寶TBC、南方SGO、Pinncle等數(shù)據(jù)處理軟件，都有坐標(biāo)系轉(zhuǎn)換模塊，有些功能比較齊全，如在TBC軟件中包含了七參數(shù)法、格網(wǎng)法、多元回歸法；LGO軟件中有格網(wǎng)法、七參數(shù)法、三參數(shù)法、格網(wǎng)與參法結(jié)合法，有三維轉(zhuǎn)換也有二維轉(zhuǎn)換。在實(shí)際應(yīng)用中，可以結(jié)合測(cè)區(qū)大小、范圍內(nèi)已知點(diǎn)的數(shù)量與分布情況決定采用哪一種方法。這些軟件在學(xué)習(xí)、教學(xué)過(guò)程中還是存在一些問(wèn)題，軟件比較復(fù)雜龐大、學(xué)習(xí)成本比較高，計(jì)算過(guò)程不可見(jiàn)，不利于對(duì)原理的理解。

Excel是微軟公司推出的一個(gè)廣泛應(yīng)用于商業(yè)和科學(xué)領(lǐng)域的電子表格程序。它具有一些強(qiáng)大的功能和特點(diǎn)，使其成為各種應(yīng)用程序中不可缺少的工具。簡(jiǎn)潔的界面、優(yōu)秀的計(jì)算功能和圖表工具，再加上出色的市場(chǎng)營(yíng)銷(xiāo)，使Excel成為流行的個(gè)人計(jì)算機(jī)電子表格軟件。它的主要功能有制作表格、數(shù)字格式的轉(zhuǎn)換、制作圖表、函數(shù)功能、數(shù)據(jù)處理功能等[2]。強(qiáng)大的數(shù)據(jù)編輯功能，使得利用Excel編制坐標(biāo)轉(zhuǎn)換程序表格方便可行。

1 四參數(shù)坐標(biāo)轉(zhuǎn)換原理

四參數(shù)是用于兩個(gè)平面直角坐標(biāo)系之間的互相轉(zhuǎn)換，而七參數(shù)是用于兩個(gè)三維空間直角坐標(biāo)系之間的轉(zhuǎn)換。四參數(shù)可以由兩個(gè)及以上具有平面二維坐標(biāo)的已知等級(jí)控制點(diǎn)求出，高程可以用擬合方法求出，求解較為簡(jiǎn)單，也較容易理解；七參數(shù)轉(zhuǎn)換是空間坐標(biāo)系的轉(zhuǎn)換需要在測(cè)區(qū)布設(shè)一定密度的控制網(wǎng)點(diǎn)，需要三個(gè)以上已知三維坐標(biāo)的控制點(diǎn)，利用整個(gè)網(wǎng)的WGS-84坐標(biāo)系下的三維約束平差結(jié)果和當(dāng)?shù)刈鴺?biāo)系統(tǒng)的二維約束平差結(jié)果及各點(diǎn)的高程解算，參數(shù)求解比較煩瑣，理解起來(lái)相對(duì)困難。

1.1 平面坐標(biāo)轉(zhuǎn)換

點(diǎn)位的平面坐標(biāo)可以用高斯投影取得。當(dāng)?shù)孛鎯牲c(diǎn)的距離小于十千米時(shí)，地球曲率的影響可以忽略，此時(shí)可以采用四參數(shù)來(lái)完成兩種坐標(biāo)系的轉(zhuǎn)換，精度也能達(dá)到預(yù)定的需求。四參數(shù)轉(zhuǎn)換計(jì)算如公式（1） 。

[x2y2=?x?y+mcosα-sinαsinαcosαx1y1]? ? （1）

其中，x1、y1為原始坐標(biāo)，x2、y2為目標(biāo)坐標(biāo)，單位米，Δx、Δy為兩個(gè)坐標(biāo)平移參數(shù)，即兩個(gè)平面坐標(biāo)系的坐標(biāo)原點(diǎn)之間的坐標(biāo)差值。α為兩個(gè)坐標(biāo)系的坐標(biāo)軸間的旋轉(zhuǎn)角度，通過(guò)旋轉(zhuǎn)α，可以使兩個(gè)坐標(biāo)系的X和Y軸重合在一起。M為尺度縮放因子，即兩個(gè)坐標(biāo)系內(nèi)的同一段直線的長(zhǎng)度比值，實(shí)現(xiàn)尺度的比例轉(zhuǎn)換。通常m值非常接近于1。因?yàn)橛兴膫€(gè)參數(shù)所以稱之為四參數(shù)坐標(biāo)轉(zhuǎn)換。這樣通過(guò)坐標(biāo)的平移、坐標(biāo)軸的旋轉(zhuǎn)和尺度的縮放，可以實(shí)現(xiàn)平面坐標(biāo)的轉(zhuǎn)換。而這四個(gè)參數(shù)如何求解是坐標(biāo)轉(zhuǎn)換的關(guān)鍵步驟。

1.2 四參數(shù)的求解

四參數(shù)求解是坐標(biāo)轉(zhuǎn)換的一個(gè)逆過(guò)程，求解需要至少兩對(duì)已知點(diǎn)，分別知道其在目標(biāo)平面的坐標(biāo)和原始平面的坐標(biāo)，根據(jù)間接平差方法，令[a=m*cosα-1] ， [b=m*sinα]， 按最小二乘原理推導(dǎo)[3]，限于篇幅直接給出公式如下：

[L=x21-x11y21-y11…x2i-x1iy2i-y1i…x2n-x1ny2n-y1n，B=1001…10…10…01…01x11-y11y11x11…x1iy1i…x1ny1n…-y1ix1i…-y1nx1n，]

[X=?x?yab，P=1001000000001001]? ?（2）

因?yàn)橐粋€(gè)點(diǎn)只能建立兩個(gè)誤差方程，要解四個(gè)未知數(shù)，至少需要兩個(gè)點(diǎn)建立四個(gè)誤差方程，如果有多的點(diǎn)，就使用最小二乘法。

求得結(jié)果X=（BTPB）-1（BTPL），可以得到Δx、Δy、a、b的值，進(jìn)而可以解得α=acrtan[b/（a+1）]、m=sqrt[（a+b）2+b2]，得到Δx、Δy、α、m四個(gè)參數(shù)值[4-5]。

1.3 四參數(shù)轉(zhuǎn)換的可靠精度

轉(zhuǎn)換參數(shù)的求取時(shí)，所采用的已知點(diǎn)的位置分布需要注意如下幾個(gè)方面：已知點(diǎn)的選取位置應(yīng)分布合理均勻，分布于測(cè)區(qū)四周和中心。為達(dá)到既定轉(zhuǎn)換精度，盡量采用多個(gè)已知點(diǎn)，讓這些點(diǎn)位能完全并均勻覆蓋整個(gè)轉(zhuǎn)換區(qū)域。留取幾個(gè)檢查點(diǎn)，作為檢核。為了保證參數(shù)求解精度，需要注意以下幾點(diǎn)：

1） 在參數(shù)求解前，對(duì)已知坐標(biāo)進(jìn)行預(yù)處理，保證數(shù)據(jù)準(zhǔn)確。

2） 控制坐標(biāo)分布、密度。

3） 殘差分析。參數(shù)求解后對(duì)殘差進(jìn)行綜合分析，剔除誤差較大的點(diǎn)，提高可靠性。

2 Excel的電子表格編輯設(shè)計(jì)

在 Excel中，函數(shù)實(shí)際上是一個(gè)預(yù)先編寫(xiě)的特定計(jì)算公式。按照這個(gè)特定的計(jì)算公式對(duì)一個(gè)或多個(gè)值執(zhí)行計(jì)算，并得出一個(gè)或多個(gè)運(yùn)算結(jié)果，叫作函數(shù)值。使用這些函數(shù)不僅可以完成許多復(fù)雜的計(jì)算，而且還可以簡(jiǎn)化公式的繁雜程度。使用函數(shù)可以簡(jiǎn)化或縮短工作表中的公式，使數(shù)據(jù)處理簡(jiǎn)單方便。

Excel每一個(gè)函數(shù)都由一系列的參數(shù)來(lái)控制其輸入輸出結(jié)果。函數(shù)的參數(shù)可以是數(shù)值、文本、單元格引用等數(shù)據(jù)類(lèi)型，將函數(shù)的參數(shù)傳入函數(shù)中進(jìn)行處理得到結(jié)果。其中，一些函數(shù)會(huì)有循環(huán)嵌套，即通過(guò)一個(gè)函數(shù)調(diào)用另一個(gè)函數(shù)的方式進(jìn)行計(jì)算，從而實(shí)現(xiàn)更為復(fù)雜的計(jì)算需求。

2.1 所用函數(shù)

實(shí)現(xiàn)坐標(biāo)轉(zhuǎn)換需要用到矩陣運(yùn)算，Excel中有矩陣運(yùn)算相關(guān)函數(shù)可以方便得到所需結(jié)果，常見(jiàn)的矩陣函數(shù)有：1） 矩陣轉(zhuǎn)置：TRANSPOSE（原矩陣）；2） 矩陣乘法：A、B矩陣相乘MMULT（矩陣A，矩陣B）；3） 求逆矩陣：MINVERSE（原矩陣）；4） 求矩陣行列式的值 MDETERM（矩陣）；5） INDIRECT函數(shù)。在使用矩陣函數(shù)時(shí)候在公式編輯欄中輸入公式：如“=MINVERSE（矩陣A）”，按下“Shift+Ctrl+Enter”組合鍵，即可得到對(duì)應(yīng)的逆矩陣。

2.2 具體編輯設(shè)計(jì)

為了防止數(shù)據(jù)編輯出錯(cuò)、界面簡(jiǎn)潔，在Excel中新建兩個(gè)工作表并命名為“輸入”和“過(guò)程”。在“輸入”工作表中編輯輸入已知點(diǎn)和坐標(biāo)轉(zhuǎn)換如圖1：

界面數(shù)據(jù)錄入分為目標(biāo)坐標(biāo)和源坐標(biāo)，各自對(duì)應(yīng)X、Y坐標(biāo)。按行從上到下排列分為P1、P2、P3…Pn點(diǎn)。右側(cè)F列編輯殘差列，對(duì)應(yīng)每個(gè)點(diǎn)轉(zhuǎn)換后的殘差，來(lái)監(jiān)測(cè)已知點(diǎn)的數(shù)據(jù)質(zhì)量。殘差的計(jì)算方式為：1） 由已知點(diǎn)求得轉(zhuǎn)換參數(shù)；2） 源坐標(biāo)按第一步的轉(zhuǎn)換參數(shù)求出轉(zhuǎn)換后的坐標(biāo)；3） 轉(zhuǎn)換后的坐標(biāo)和已知目標(biāo)坐標(biāo)X、Y進(jìn)行求差，再求平方和、開(kāi)方即為殘差。

在“過(guò)程”工作表中編輯生成過(guò)程矩陣，其中L為2n行1列的矩陣，B為2n行4列的矩陣（n為已知點(diǎn)的數(shù)量）。實(shí)際應(yīng)用中已知點(diǎn)數(shù)量為2個(gè)及以上，數(shù)量并不是固定的，因此L、B兩個(gè)矩陣宜豎向排列，以便在點(diǎn)數(shù)增加時(shí)進(jìn)行矩陣擴(kuò)展如圖2所示。

以L矩陣第1、2行為例，第一行輸入“=IF（F1>=1，INDIRECT（"輸入！D3"）-INDIRECT（"輸入！B3"），0）”，第二行輸入“=IF（F1>=1，INDIRECT（"輸入！E3"）-INDIRECT（"輸入！C3"），0）”。F1單元格為已知點(diǎn)數(shù)量可以用COUNT函數(shù)取得，即已知點(diǎn)數(shù)量大于等于1時(shí)第一行值為x21-x11否則為0。因?yàn)镋xcel的公式會(huì)隨著復(fù)制剪切同步變化，這會(huì)導(dǎo)致已知點(diǎn)編輯時(shí)引起后面公式變化，所以加入INDIRECT函數(shù)防止剪切單元格時(shí)公式跟隨變化。以此類(lèi)推根據(jù)公式2中L矩陣編輯其余單元格。

矩陣B第一行第一列輸入“=IF（F1>=1，1，0）”即已知點(diǎn)數(shù)量大于等于1，第一列值為1，否則為0;第二列固定輸入0，第三列輸入“=IF（F1>=1，INDIRECT（"輸入！B3"），0）”，即已知點(diǎn)大于等于1時(shí)，第三列單元格值為“輸入”B3單元格的值，否則為零；第四列輸入“=IF（F1>=1，INDIRECT（"輸入！C3"）*-1，0）”，即已知點(diǎn)F1數(shù)量大于1時(shí)，第四列單元格值為“輸入”C3單元格值乘以-1，否則為零。據(jù)此根據(jù)需要擴(kuò)展矩陣的行數(shù)，此例為5個(gè)點(diǎn)，當(dāng)輸入點(diǎn)數(shù)量少于5時(shí)，空白的地方利用IF函數(shù)根據(jù)F1單元格的點(diǎn)數(shù)自動(dòng)填充為0。即已知點(diǎn)數(shù)量小于等于5個(gè)時(shí)，矩陣會(huì)自動(dòng)把不足5的地方填充為零，不影響最終結(jié)果的計(jì)算?？梢赃_(dá)到已知點(diǎn)數(shù)量動(dòng)態(tài)使用的效果。如果需要更多的已知點(diǎn)數(shù)量，只需要手動(dòng)把矩陣范圍擴(kuò)展到相應(yīng)點(diǎn)數(shù)即可。

編輯求得B矩陣的轉(zhuǎn)置矩陣BT，因Excel的TRANSPOSE函數(shù)使用時(shí)需要固定單元格數(shù)量，所以矩陣BT采用手動(dòng)引用編輯生成，以便矩陣的跟隨矩陣B動(dòng)態(tài)擴(kuò)展。

矩陣BT第一行從左到右等于矩陣B的第一列從上到下，矩陣BT第二行從左到右等于矩陣B的第二列從上到下，依次編輯好矩陣BT，結(jié)果如圖3所示。

編輯構(gòu)造矩陣BTPB，權(quán)系數(shù)矩陣P默認(rèn)為1的單位陣，在此忽略，BTPB結(jié)果為4行4列的矩陣，它的大小是不變的，和已知點(diǎn)數(shù)量無(wú)關(guān)。在空白處選中4行4列單元格，在編輯欄輸入“=MMULT（矩陣BT，矩陣B）”，然后同按Ctrl+Shift+回車(chē)鍵返回得到結(jié)果圖4所示：

求逆矩陣（BTPB）-1，逆矩陣大小和源矩陣大小相同，在空白處選中4行4列單元格，在編輯欄輸入“=MINVERSE（矩陣BTPB）” 然后同按Ctrl+Shift+回車(chē)鍵返回得到結(jié)果圖5所示：

編輯構(gòu)造矩陣BTPL， 權(quán)系數(shù)矩陣P默認(rèn)為1的單位陣，在此忽略。矩陣BTPL結(jié)果為4行1列的矩陣，它的大小也是不變的，和已知點(diǎn)數(shù)量無(wú)關(guān)?？瞻滋庍x中4行1列單元格，在編輯欄輸入“=MMULT（矩陣BT，矩陣L）”，然后同按Ctrl+Shift+回車(chē)鍵返回得到結(jié)果圖6左邊BTPL：

求方程的最小二乘解X=（BTPB）-1*BTPL，其結(jié)果為4行1列固定大小矩陣。空白處選中4行1列單元格，在編輯欄輸入“=MMULT（矩陣（BTPB）-1， 矩陣BTPL）”，引用好各個(gè)參數(shù)矩陣的位置，然后同時(shí)按Ctrl+Shift+Enter鍵返回得到矩陣X結(jié)果圖6右邊所示，分別為Δx、Δy、a、b的值。按照公式2最終解得α=acrtan[b/（a+1）]、m=sqrt[（a+b）2+b2]，得到Δx、Δy、α、m四個(gè)參數(shù)值。

3 實(shí)例與結(jié)果分析

現(xiàn)有某工程已知點(diǎn)坐標(biāo)如表1。

如圖7所示，輸入P1、P2、P3三個(gè)點(diǎn)得到四個(gè)轉(zhuǎn)換參數(shù)平移Δx、Δx、α、m分別為53.193m、-50.134m、6.578E-6弧度、1.0000003。殘差均小于1mm。并有殘差列可供使用者評(píng)估各個(gè)已知點(diǎn)的精度。轉(zhuǎn)換后的參數(shù)可以用于實(shí)際測(cè)量中。

4 結(jié)束語(yǔ)

本文針對(duì)工程測(cè)量中重要的坐標(biāo)轉(zhuǎn)換實(shí)現(xiàn)問(wèn)題進(jìn)行研究，探討二維平面坐標(biāo)轉(zhuǎn)換的實(shí)現(xiàn)原理、精度控制和注意事項(xiàng)。經(jīng)過(guò)矩陣運(yùn)算、實(shí)現(xiàn)使用Excel 內(nèi)部函數(shù)進(jìn)行坐標(biāo)轉(zhuǎn)換。經(jīng)頁(yè)面設(shè)計(jì)和調(diào)試測(cè)試達(dá)到精度要求。該方法相對(duì)于專(zhuān)用軟件和Excel VBA編程方法，實(shí)現(xiàn)簡(jiǎn)單，無(wú)需復(fù)雜的編程過(guò)程、學(xué)習(xí)成本?。贿^(guò)程數(shù)據(jù)可見(jiàn)，過(guò)程中生成的矩陣可以分析數(shù)據(jù)，更能方便理解其轉(zhuǎn)換的意義、錄入保存數(shù)據(jù)方便，適用于多種測(cè)量、教學(xué)場(chǎng)景。

參考文獻(xiàn)：

[1] 李征航，黃勁松.GPS測(cè)量與數(shù)據(jù)處理[M].3版.武漢：武漢大學(xué)出版社，2016.

[2] 馮注龍，丘榮茂.Excel之光：高效工作的Excel完全手冊(cè)[M].北京：電子工業(yè)出版社，2019.

[3] 武漢大學(xué)測(cè)繪學(xué)院測(cè)量平差學(xué)科組.誤差理論與測(cè)量平差基礎(chǔ)[M].2版.武漢：武漢大學(xué)出版社，2009.

[4] 王萌，高遠(yuǎn)，范咪娜.區(qū)域坐標(biāo)轉(zhuǎn)換模型對(duì)比分析[J].礦山測(cè)量，2021，49（2）：71-76.

[5] 謝偉杰.平面四參數(shù)法坐標(biāo)轉(zhuǎn)換的實(shí)例應(yīng)用探討[J].工程建設(shè)與設(shè)計(jì)，2020（13）：46-47，50.

【通聯(lián)編輯：梁書(shū)】