高昭良

(1.武漢大學(xué)測(cè)繪遙感信息工程國(guó)家重點(diǎn)實(shí)驗(yàn)室，湖北武漢 430079;2.福州市勘測(cè)院，福建福州 350003)

1 引言

從2008年7月1日起正式啟用中國(guó)大地坐標(biāo)系統(tǒng)2000(CGCS 2000)作為國(guó)家法定的坐標(biāo)系，作為我國(guó)新一代的平面基準(zhǔn)［1］。該系統(tǒng)有利于GPS快速、精確地獲取高精度城市坐標(biāo)和高程成果，有利于GIS與GPS應(yīng)用的結(jié)合，進(jìn)一步提升數(shù)字城市的綜合服務(wù)能力。

目前全國(guó)各地城市基本還是沿用原有的1954北京坐標(biāo)系、1980西安坐標(biāo)系統(tǒng)及地方坐標(biāo)系統(tǒng)，因此，各城市都需要進(jìn)行采用CGCS2000建立城市平面坐標(biāo)系統(tǒng)。陳俊勇院士等探討國(guó)家小比例尺地圖的“西安80坐標(biāo)系”更換為“大地參考系1980”(GRS80)，給出這些變化量的大小區(qū)間并討論它對(duì)地圖表示可能產(chǎn)生的影響［2］。李江衛(wèi)、解斌等人以武漢市為例，提出了以投影改正平方和最小和投影變形大于 2.5 cm/km的點(diǎn)數(shù)最少為原則進(jìn)行最佳中央子午線和抵償高程面選取的方法［3］。武豐雷、劉曦燦等人提出建立斜軸墨卡托投影濟(jì)南市坐標(biāo)系［4］。

在討論如何確定城市平面坐標(biāo)系統(tǒng)選擇時(shí)，相關(guān)文獻(xiàn)基本采用兩個(gè)方面來驗(yàn)證:一是采用有限的控制點(diǎn)作為計(jì)算結(jié)果驗(yàn)證，二是采用最小二乘法。由于城市區(qū)域范圍較大，地形起伏分布不均勻，因此，有限的控制點(diǎn)無法全面準(zhǔn)確概括。并且少量的地形突變會(huì)導(dǎo)致最小二乘法的最優(yōu)評(píng)價(jià)結(jié)果畸變。本文提出的基于DEM的CGCS2000城市平面坐標(biāo)系統(tǒng)確定方法可以兼顧解決這兩方面問題。

2 城市平面坐標(biāo)系統(tǒng)定義模型

一個(gè)城市坐標(biāo)系統(tǒng)確定應(yīng)根據(jù)以上基本原則，還得考慮城市具體地理范圍和形狀，最后來確定城市平面坐標(biāo)系統(tǒng)的定義模型。根據(jù)《城市測(cè)量規(guī)范》，坐標(biāo)系統(tǒng)應(yīng)遵循以下原則:

(1)投影變形盡量小，投影長(zhǎng)度變形值不大于2.5 cm/km(1/4 萬)［5］。

(2)坐標(biāo)系與國(guó)家(省)級(jí)控制點(diǎn)應(yīng)聯(lián)測(cè)，并建立嚴(yán)密的轉(zhuǎn)換關(guān)系［6］。

2.1 高斯投影長(zhǎng)度變形分析

由于定義國(guó)家大地坐標(biāo)系的橢球面是一個(gè)凸起的、不可展平的曲面，當(dāng)采用高斯正形投影將這個(gè)曲面上的元素投影到平面上時(shí)，投影后就會(huì)發(fā)生長(zhǎng)度變形問題。高斯投影時(shí)，先把地面觀測(cè)值歸化至參考橢球面上，再把參考橢球面上的觀測(cè)值歸化至高斯平面上。兩次歸化產(chǎn)生的長(zhǎng)度變形值稱為投影變形。投影變形△D計(jì)算過程如下:

上式中各變量說明如下:

Hm:測(cè)量?jī)啥藶楦叱鰠⒖紮E球面的平均高程，單位:m;

H0:抵償面高程，單位:m;

D:測(cè)量長(zhǎng)度兩端點(diǎn)平均高程面的水平距離，單位:m;

RA:歸算邊方向參考橢球法截弧的曲率半徑，可取概值 6 371 000 m，單位:m;

ym:測(cè)距邊兩端點(diǎn)投影后近似橫坐標(biāo)平均值，單位:m;

Rm:測(cè)距邊中點(diǎn)的平均曲率半徑，可取概值6 371 000 m，單位:m。

2.2 采用格網(wǎng)占比法確定最佳關(guān)系

對(duì)于面積較大的城市或地區(qū)，在確定CGCS2000框架下的平面坐標(biāo)系時(shí)，還應(yīng)根據(jù)區(qū)域地理特點(diǎn)以及建設(shè)發(fā)展的實(shí)際情況，合理選取最佳中央子午線和抵償高程面。不斷調(diào)整高程參考面的大地高為H，使其投影變形滿足要求，條件式可近似表示為:

若確定的中央子午線與高程抵償面在城市坐標(biāo)系統(tǒng)覆蓋的范圍內(nèi)，能夠較好地滿足式(2)或式(3)的關(guān)系，則該中央子午線與高程抵償面關(guān)系最佳。

3 基于DEM模型數(shù)據(jù)的計(jì)算模型

大部分學(xué)者計(jì)算最優(yōu)值采用最小二乘法∑△D2=min［9，10］模型，由于選取的控制點(diǎn)有限，不能全面均勻代表該區(qū)域投影變形，因此，條件通過遍歷搜索獲得∑△D2=min的方案是否可以確定為最優(yōu)，值得商榷。本文基于福州地區(qū) 1∶1萬DEM數(shù)據(jù)，提取分布均勻的高程網(wǎng)格，利用以上的投影計(jì)算模型，計(jì)算最佳高程抵償面。最佳高程抵償面與中央子午線的計(jì)算采取遍歷搜索尋優(yōu)法:搜索范圍在中央經(jīng)度［Lmin，Lmax］，抵償面高程［Hmin，Hmax］內(nèi)，按適合精度的步長(zhǎng)，以投影改正平方和最小為原則，以變形值小于 2.5 cm的格網(wǎng)占比(符合變形要求的格網(wǎng)個(gè)數(shù)/坐標(biāo)定義覆蓋區(qū)域的所有格網(wǎng)個(gè)數(shù):m/n，if g(i)＜2.5)then m=m+1，n為所有格網(wǎng)數(shù))作為參考觀察指標(biāo)，按逐步縮小搜索區(qū)間，直到理論最優(yōu)值所在范圍內(nèi)達(dá)到預(yù)期設(shè)定要求，最后進(jìn)行幾種最優(yōu)方案比選。本文將這種方法簡(jiǎn)稱為格網(wǎng)占比法。

3.1 基于福州地區(qū)DEM提取500 m格網(wǎng)間距高程點(diǎn)

為了避免采用有限的控制點(diǎn)不夠均勻，代表性不夠，從而導(dǎo)致最優(yōu)關(guān)系計(jì)算不夠嚴(yán)密。因此，利用福州地區(qū) 1∶1萬DEM數(shù)據(jù)，福州地區(qū) 1∶1萬DEM數(shù)據(jù)格網(wǎng)尺寸為 12.5 m，高程精度 1 m。坐標(biāo)系統(tǒng)為1980西安坐標(biāo)系，1985黃海高程基準(zhǔn)。根據(jù)福州山地地形的變化情況，選擇邊長(zhǎng) 250 m格網(wǎng)來提取高程，基本可以滿足計(jì)算需要。按一定的格網(wǎng)間距提取三維坐標(biāo)點(diǎn)，進(jìn)行最優(yōu)關(guān)系計(jì)算。

3.2 最優(yōu)中央子午線和抵償高程面確定算法流程

(1)基礎(chǔ)數(shù)據(jù)準(zhǔn)備:按投影經(jīng)度步長(zhǎng)2分(約3 km)、抵償面高程步長(zhǎng) 50 m設(shè)置，生成所有計(jì)算方案表1、基于CGCS2000坐標(biāo)系統(tǒng)的格網(wǎng)中心點(diǎn)數(shù)據(jù)表2。

(2)按遍歷計(jì)算每個(gè)方案每個(gè)格網(wǎng)點(diǎn)的變形值△D，存入表3。

(3)完成一個(gè)方案的變形值計(jì)算后，將∑△D2和△D≤2.5 cm的格網(wǎng)占比結(jié)果計(jì)算存入表3的對(duì)應(yīng)方案符合。

(4)最后遍歷搜索法對(duì)比最優(yōu)方案。

具體算法流程圖如圖1所示:

圖1 計(jì)算流程圖

3.3 最優(yōu)中央子午線和抵償高程面計(jì)算實(shí)驗(yàn)結(jié)果

以福州地區(qū)數(shù)據(jù)為例進(jìn)行實(shí)驗(yàn)計(jì)算，選擇精度118°40'～119°40'范圍內(nèi)，按 2 分步長(zhǎng)遞增，高程抵償面從 0 m～600 m，按 20 m步長(zhǎng)遞增，對(duì)福州地區(qū)范圍內(nèi)的格網(wǎng)數(shù)據(jù)進(jìn)行遍歷計(jì)算。計(jì)算結(jié)果得到最優(yōu)方案不同，如下:

(1)格網(wǎng)占比法:按△D≤2.5 cm的網(wǎng)格占比最大計(jì)算結(jié)果及對(duì)比數(shù)據(jù)。

變形值符合要求的格網(wǎng)占比計(jì)算結(jié)果及對(duì)比數(shù)據(jù)表 表1

變形值符合要求的格網(wǎng)占比模型計(jì)算結(jié)果表1列出最大的5個(gè)方案，最優(yōu)的方案是中央子午線為119°34'，高程抵償面為 50 m，符合每千米變形△D≤2.5 cm的網(wǎng)格占福州區(qū)域內(nèi)格網(wǎng)總數(shù)的55.62%。但是每千米變形值平方和為39.78并非最小。空間分布專題圖如圖2(b)所示。

(2)最小二乘法:∑△D2最小值計(jì)算結(jié)果及對(duì)比數(shù)據(jù)。

每公里變形值平方和計(jì)算結(jié)果及對(duì)比數(shù)據(jù)表表2

∑△D2最小值模型計(jì)算結(jié)果表2列出最小的5個(gè)方案，最優(yōu)的方案是中央子午線為119°30'，高程抵償面為 150 m，符合每千米變形值平方和最小為38.00。但是符合每千米變形△D≤2.5 cm的網(wǎng)格占福州區(qū)域內(nèi)格網(wǎng)總數(shù)為37.90%?？臻g分布專題圖如圖2(b)所示。

圖2 格網(wǎng)占比法與最小二乘法最優(yōu)方案對(duì)比圖

從圖2很明顯看出左圖的符合每千米變形△D≤2.5 cm的網(wǎng)格覆蓋量比右圖合理，左圖福州市中心及沿海經(jīng)濟(jì)開發(fā)熱點(diǎn)平原區(qū)域基本符合要求，而右圖福州市中心城區(qū)大東邊區(qū)域與沿海平潭島及羅源區(qū)域大部分都不符合變形要求。顯然，按格網(wǎng)占比法計(jì)算模型結(jié)果更為科學(xué)合理、更符合實(shí)際應(yīng)用需要。按福州地方平面坐標(biāo)系統(tǒng)選擇的投影中央子午線為119°18'21″方 案［11］，計(jì) 算 結(jié) 果 為 符 合 每 千 米 變形△D≤2.5 cm的網(wǎng)格占比為46.04%，每千米變形值平方和值為64.15，并不是最優(yōu)的方案。福州地區(qū)西北面基本為丘陵高山，許多高程基本在 400 m以上，經(jīng)過遍歷搜索計(jì)算，驗(yàn)證了該地形變化規(guī)律不適合采用高程抵償面長(zhǎng)度歸化改正與投影變形進(jìn)行抵消的計(jì)算模型。對(duì)于該高山區(qū)域的測(cè)圖，只能采用抬高投影面的方法。

4 結(jié)論

本文采用DEM數(shù)據(jù)與地理信息系統(tǒng)計(jì)算方法，對(duì)CGCS2000城市平面坐標(biāo)系統(tǒng)定義的兩個(gè)計(jì)算模型進(jìn)行研究，得出采用格網(wǎng)占比法方法得出的最優(yōu)高程抵償面與中央子午線的選擇方案更具嚴(yán)密性，該方法對(duì)城市區(qū)域CGCS2000平面坐標(biāo)系統(tǒng)定義具有較大的參考意義。

［1］陳俊勇，張鵬.國(guó)家測(cè)繪基準(zhǔn)“十二五”重大項(xiàng)目的思考［J］.武漢大學(xué)學(xué)報(bào)·信息科學(xué)版，2009，34(10):1136～1138.

［2］陳俊勇.中國(guó)采用地心三維坐標(biāo)系統(tǒng)對(duì)現(xiàn)有地圖的影響［J］.測(cè)繪學(xué)報(bào)，2004，32(4):283 ～288.

［3］李江衛(wèi)，白潔，郭際明等.基于CGCS2000的城市平面坐標(biāo)系最佳選?。跩］.城市勘測(cè)，2011，3(4):120～123.

［4］武豐雷，劉曦燦，王方等.濟(jì)南市獨(dú)立坐標(biāo)系的改造探討［J］.城市勘測(cè)，2010，2(2):105 ～106.

［5］CJJ T8－2011.城市測(cè)量規(guī)范［S］.

［6］吳云孫，晁定波，楊堂堂.淺談廣州市平面控制測(cè)量坐標(biāo)系統(tǒng)的選?。跩］.測(cè)繪通報(bào)，2005，23(8):42～43.

［7］陸鵬程，林冬偉.斜軸墨卡托投影模型及其應(yīng)用分析［J］.鐵道勘察，2010，12(4):26～29.

［8］呂忠剛，許世寧.關(guān)于抵償高程面與移動(dòng)中央子午線最佳選取問題的研究［J］.東北測(cè)繪，2003，25(2):3～8.

［9］謝征海，張澤烈.城市獨(dú)立坐標(biāo)系平面控制網(wǎng)之?dāng)U建及改造［J］.北京測(cè)繪，2000，23(4):22～23.

［10］孔建，姚宜斌，吳寒.整體最小二乘的迭代解法［J］.武漢大學(xué)學(xué)報(bào)·信息科學(xué)版，2010，35(6):711～714.

［11］高昭良.基于相似變換的城市坐標(biāo)系統(tǒng)轉(zhuǎn)換研究［J］.福州大學(xué)學(xué)報(bào)·自然科學(xué)版，2010，38(1):54～57.