search
尋找貓咪~QQ 地點 桃園市桃園區 Taoyuan , Taoyuan

專輯論文| 宋樹華:幾種地球剖分網格分析與初探

《測繪學報》

構建與學術的橋樑 拉近與權威的距離

1

22345617

1. 蒼穹數碼技術股份有限公司, 北京 100081;

2. 61541部隊, 北京 100094;

3. 石油新疆油田公司數據公司, 克拉瑪依 834000;

4. 天繪衛星中心, 北京 102102;

5. 運載火箭技術研究院研究發展中心, 北京 100076;

6. 海參某局, 北京 100841;

7. 國家信息中心, 北京 100045

收稿日期:2016-08-10; 修回日期:2016-10-20

第一作者簡介:宋樹華 (1980—),男, 博士, 高級工程師, 研究方向為GIS理論與應用。

E-mail: [email protected]

摘要:針對全球遙感衛星影像數據缺少統一組織框架的問題,從遙感影像數據特點出發,分析歸納了遙感影像統一組織框架的設計原則和設計難點,並對當前已提出的用於遙感影像統一組織的幾種剖分格網方案進行概述。通過分析與比較,認為GeoSOT比較適合作為遙感影像的統一組織框架。

關鍵詞:遙感影像 剖分格網 GeoSOT 數據組織

Analysis and Research on Several Global Subdivision Grids

, DONG Fang22345617

Abstract: In order to solve the problem that lacking of an unified organization frame about global remote sensing satellite image data, this paper introduces serval global subdivision grids as the unified organization frame for remote sensing image. Based on the characteristics of remote sensing image data, this paper analyzes and summarizes the design principles and difficulties of the organization frame. Based on analysis and comparison with these grids, GeoSOT is more suitable as the unified organization frame for remote sensing image. To provide a reference for the global remote sensing image organization.

Key words: remote sensing image GeoSOT subdivision grid data organization

當前的空間數據組織方面,尤其是地圖數據,在統一組織方面已經形成一套完整的體系,如1992年頒布的《國家基本比例尺地形圖分幅和編號GB/T130 89—92》標準等。針對遙感影像數據的組織與管理方面,在某類或者某系列衛星的影像數據統一組織上也有很完備的體系[

1

],但缺少適合所有遙感衛星影像數據的統一組織上的網格框架,由此給基於遙感影像數據共享、檢索、全球無縫表達、分發、編目、計算以及數據更新等應用帶來諸多問題。本文從遙感影像數據特點出發,歸納了統一組織框架的設計原則和設計難點,分析了幾種剖分格網作為遙感影像統一組織框架的優勢和缺陷,為遙感影像統一組織提供參考。

1 遙感影像統一組織框架設計原則與難點

遙感數據作為空間數據的一種,通常與傳統的地圖數據結合使用。因此,全球遙感影像的統一組織框架需考慮國家標準地圖分幅體系與衛星數據分景體系之間兼容問題。

遙感衛星影像數據是關於地球表層或者地球表面附近區域信息的快照,具有與地球形狀相似的性質。因此,全球遙感影像統一組織框架既要滿足全球尺度無縫覆蓋的要求,又要滿足局部平面表達、精準量算的要求。如何在兩者之間達到平衡,成為組織框架設計的難點之一。

現階段對全球的空間數據採用經緯度進行組織,而對局部區域或者平面空間數據採用高斯坐標系下的地圖分幅框架進行組織。鑒於空間數據組織現狀,依據全球離散網格模型標準[

2

],國內外研究比較多的地球剖分模型,如基於三角形或菱形的剖分模型[

3

-

11

]、基於正六邊形剖分模型[

12

-

16

]、以及基於Voronoi圖的剖分模型[

17

-

20

],並不能很好地適合空間數據組織與管理。因此,在設計全球遙感影像數據統一組織框架時,在符合「Goodchild標準」

(1) 可與地圖分幅體系兼容。框架中各級不同大小的網格單元,又稱剖分單元,能很方便地兼容或者成為地圖分幅的有機部分,便於用於傳統空間數據整合。

(2) 劃分單元的形狀盡量為矩形或正方形。框架中的網格單元形狀最好為矩形或正方形,使其既可與影像中的像元點陣很好地對應,又便於影像的表達。

(3) 具有完整的編碼體系。框架有一個完整的編碼體系,使得大到地球、小到厘米級網格都有編碼,且層與層之間的編碼有很好的嵌套關係,以便於計算。

(4) 每個單元格具有明確的地理位置。框架中每個網格都有固定的空間範圍,即每個網格中心、四個角點都有明確的坐標位置,可方便地與影像中坐標對應;每個網格編碼在地球上具有準確的地理含義。

(5) 便於多源數據融合。在採用統一組織框架網格來切割影像時,可根據不同的解析度,便於多源數據融合。

2 幾種剖分格網分析

鑒於國內外學者對柏拉圖立體、正多邊形、Voronoi圖等剖分模型進行了深入研究,並根據統一組織框架設計原則,本文只對等距、等緯差、正六面體、GeoSOT等幾種剖分格網用於遙感影像統一組織的方案進行分析,具體如下。

2.1 等距格網

等距剖分網格即基於CGCS2000大地坐標系統的地球橢球體,沿經線方向採用等距、沿緯線方向採用與緯圈相等的正方形建立全球等距格網框架。因此,等距剖分格網又分如下兩種情況。

2.1.1 基於基點的等距剖分格網框架

基於基點的等距剖分格網是指剖分單元以地球表面某點作為起始點形成全球無縫覆蓋的等距剖分格網。比較典型的基點是地球兩極點、本初子午線與赤道的交點、任意其他經線與赤道的交點。為了便於描述與計算,本文將其基點設為北極點,並假設正方形剖分單元的邊長為d。基於北極的等距剖分格網示意圖如圖 1,圖中剖分單元的邊長與兩個緯圈之間的距離相等,並採用不同的填充對不同距離的剖分面片進行區別。

圖 1 基於北極的等距剖分格網示意圖Fig. 1 The equidistant subdivision grid based on the arctic pole

,若假設以正方形邊長

d

為從北極點到某緯度之間沿大圓經線的弧長,則通過式 (1) 可算出該緯線圈所處的緯度,並由式 (2) 計算出緯圈的半徑,由此得計算緯圈的周長。

(1)

(2)

式中,

a

為長半軸;

e

為第一偏心率;

ϕ

為緯度,

E

=(

l

000012

表 1 基於北極的等距剖分格網分析統計表Tab. 1 The statistical list about equidistant subdivision grid based on the arctic pole

距北極弧長/m緯度/(°)緯圈半徑/m緯圈周長/m面片數面片周長/m緯圈、面片周長之差/m
25 00089.77625 000.399157 082.1374200 00042 917.863
50 00089.55249 999.291314 154.8116400 00085 845.19
75 00089.32974 998.522471 229.61136600 000128 770.389
100 00089.10599 995.468628290.05760600 00028 290.057
125 00088.881124 991.975785 347.73888600 000185 347.738
150 00088.657149 986.535942 393.191132800 000142 393.191
175 00088.433174 977.6431 099 416.957172800 000299 416.957
200 00088.209199 967.1441 256 430.6212241 000 000256 430.621
225 00087.986224 953.5311 413 424.7242841 000 000413 424.724
250 00087.762249 936.4161 570 396.823441 200 000370 396.82

表 1可知,在第一圈剖分面片的周長與對應緯圈周長之差已經達到了42 917.863 m,並隨著正方形單元從北極向赤道逼近,它們之間的差值越大。當d取不同的值時,也有類似的規律。因此,基於基點的等距剖分格網不能很好滿足框架設計原則。

2.1.2 基於經線的等距剖分格網框架

基於經線的等距剖分格網與基於基點的等距格網相似,只是在剖分面片的展開方式上,基於經線的等距剖分格網以某一經線作為「起跑線」,每一等距緯圈的剖分面片覆蓋地球並由此展開 (如圖 2),剖分單元的邊長與兩個緯圈之間的距離相等,並採用不同的填充對不同緯圈的剖分面片進行區別。

圖 2 基於本初子午線的等距剖分格網示意圖Fig. 2 The equidistant subdivision grid based on the prime meridian

採用式 (1) 可計算出不同緯圈的緯度,並由式 (2) 計算出各緯圈的半徑。在此基礎上,計算出各個緯圈內正方形面片的個數,即可計算在此高緯部分的各個正方形單元的重疊度,計算結果見表 2

表 2 基於本初子午線的等距剖分格網分析統計表Tab. 2 The statistical list about equidistant subdivision grid based on the prime meridian

距北極點的弧長/m緯度/(°)緯圈半徑/m緯圈周長/m上下緯圈周長之差/m面片數高緯重疊百分比/(%)
25 00089.77625 000.399157 082.137157 082.1376.28100
50 00089.55249 999.291314 154.810157 072.67412.5750.00
75 00089.32974 998.522471 229.611157 074.80018.8533.33
100 00089.10599 995.468628 290.057157 060.44625.1325.00
125 00088.881124 991.975785 347.738157 057.68131.4120.00
9 950 0000.4706 377 923.88040 073 677.6111 598.3211 602.950.004 0
9 975 0000.2446 378 079.61340 074 656.110978.4981 602.990.002 4
10 000 0000.0186 378 136.69540 075 014.769358.6601 603.000.000 9

在每個緯圈內,正方形面片均不能完整地覆蓋整個緯度帶,且在高緯區域,越靠近極地,正方形單元之間的重疊度越大。因此,這種剖分方式也不能使得正方形無縫無疊地覆蓋整個地球。

2.2 等緯差格網

在傳統的測繪數據管理上,多是採用等經差、等緯差方式對空間數據進行組織和管理,例如標準地形圖分幅管理機制即採用這種方式。且當經差或者緯差小到一定程度時,其所代表地球表面的地理空間範圍可以近似視為平面,如在赤道附近經度相差10′×10′的區域,其範圍約為一個18.5 km×18.5 km。

若遙感影像統一組織框架採用緯差為φ′、以及緯差φ′對應的地表距離為正方形的邊長為剖分單元進行等緯差剖分,其中,φ′=1、2、3、…、30。

基於式 (2) 的緯圈半徑,可算出每一緯度所對應的緯圈周長、相鄰兩緯圈之間的距離,以及兩者的比值,即為緯差φ′對應的正方形網格個數。為了計算方便,取φ′=10,其計算結果如表 3所示。

表 3 基於等緯差剖分格網分析統計表Tab. 3 The statistical list about equal difference of longitude and latitude subdivision grid based on the prime meridian

緯度緯圈上下緯圈周長之差/m上下緯圈的距離/m面片數
半徑/m周長/m
10′6 378 110.19640 074 848.27505.24118 429.0462 174.548
20′6 378 029.78440 074 343.03842.06418 429.052 174.52
30′6 377 895.76640 073 500.971 178.88018 429.0562 174.474
87°10′316 334.6861 987 589.451116 828.20218 615.176106.773
87°20′297 740.900 11 870 761.249116 844.22818 615.23100.496
87°30′279 144.563 61 753 917.021116 859.61418 615.28194.219
87°40′260 545.778 41 637 057.407116 873.62318 615.32987.941
87°50′241 944.763 71 520 183.784116 886.62218 615.37381.663
88°223 341.681 403 297.162116 898.98318 615.41475.384

根據表 3結果分析,按此種剖分方案,某一地圖圖幅內剖分面片的分佈如圖 3。出現這種現象的原因是,由於從赤道到兩極等緯差之間的地表距離並不相等,導致不同緯度帶的正方形單元大小不相等,也不能完全部覆蓋該緯度帶對應的整個圖幅。

圖 3 等緯差標準剖分框架示意圖Fig. 3 The equal difference of longitude and latitude subdivision grid2.3 正六面體格網

2.3.1 基於正六面體的剖分框架

針對剖分單元為矩形或正方形的原則,在柏拉圖立體中,正六面體內接地球,有可能使得地球表面剖分為由多個四邊形組成,如圖 4

圖 4 基於正六面體的剖分框架示意圖Fig. 4 The subdivision grid based on cube

(1) 設地球球面內接的正六面體棱長為

d

,赤道面與六面體的一面平行。六面體對角線等於球面直徑,即,因此六面體的8個頂點經緯度分別為

即8點為 (35.264 390°,0°)、(-35.264 390°,0°)、(35.264 390°,90°)、(-35.264 390°,90°)、(35.264 390°,180°)、(-35.264 390°,180°)、(35.264 390°,-90°)、(-35.264 390°,-90°),其相應的直角坐標容易計算,此處略。

(2) 從步驟 (1) 任選一面進行討論,不妨選如下四邊形:(35.264 390°,0°)、(35.264 390°,90°)、(-35.264 390°,90°)、(-35.264 390°,0°)。由此可知0級剖分對應的邊長為赤道周長的四分之一,即約10 008 014.6 m。

(3) 將上一步得到的四邊形,逐個進行 (近似)4等分:分別取各邊中點,並用測地線 (大圓弧) 連接對邊的中點 (如圖 4(c))。

(4) 對同一解析度的每個四邊形,重複步驟 (3)。

(5) 當解析度逐級提高 (四邊形邊長逐級減半) 時,對上一級的四邊形逐個重複步驟 (3)、(4)。

由此即可得到由不同層級、不同尺度的剖分面片無縫覆蓋地球表面的剖分體系。

2.3.2 變形分析

下面對基於正六面體剖分框架中的剖分面片的面積變形、邊長變形和角度變形進行分析。

2.3.2.1 剖分面片的面積變形分析

該剖分體系中的剖分面片在地球表面為一個橢球面四邊形 (如圖 5)。設任意一個橢球面四邊形ABCD,其面積的計算方法:

圖 5 球面四邊形ABCD示意圖Fig. 5 The ellipsoidal quadrilateral ABCD

先用餘弦公式求出四邊形的4個球面角,然後代入式 (3) 即可得到面積:

(3)

式中,ABCD分別為橢球四邊形ABCD頂點所對應的內側球面角。球面角的計算方法:

AD

AB

分別為過A點沿著曲線AD、AB的切線,相應的切面法向量分別為

nAD

nAB

,則球面角∠

DAB

(即角

A

) 滿足

(4)

通過式 (3) 和式 (4),即可分析出不同剖分層級上剖分面片的面積變形情況,如圖 6表 4所示。

圖 6 部分等級剖分面片面積變形分布圖Fig. 6 The distribution curve about deformation area of some subdivision cell

表 4 部分等級剖分面片的面積變形統計表Tab. 4 The statistical list about deformation area of some subdivision cell

參數剖分等級 (將初始面的中線逐級二等分)234567891011121314151617標準邊長/m2 502 003.61 251 001.8625 500.9312 750.4156 375.278 187.639 093.819 546.99 773.44 886.72 443.31 221.7610.8305.4152.776.3正變形極值/(%)0.240 6120.240 6780.553 6850.669 1330.665 8180.643 8060.637 8770.631 110.628 1150.626 5580.625 7650.625 0490.624 6520.624 4390.624 3230.624 323負變形極值/(%)-2.351 87-7.874 66-13.677 1-17.325 6-19.334 4-20.384 4-20.920 8-20.452 5-21.328-21.302 8-21.290 2-20.935 1-21.222 4-21.366 7-21.439-21.438 7

圖 6可知,以初始剖分 (0級) 中心線為參照 (下同),沿水平方向,剖分面片的面積變形隨著經度的增加而增加;沿豎直方向,隨著緯度的增加而增加。可歸納為:距離初始剖分中心線越遠,剖分面片的面積變形越大。

2.3.2.2 剖分面片邊長變形分形

剖分面片的邊長,等於過該球面曲線兩端點大圓的劣弧之長,用球面距離公式即可得到。圖 7表 5是部分等級剖分面片的邊長變形分析結果。

圖 7 部分等級剖分面片邊長變形分布圖Fig. 7 The distribution curve about deformation side of some subdivision cell

表 5 部分等級剖分面片的邊長變形統計表Tab. 5 The statistical list about deformation side of some subdivision cell

標準邊長/m2 502 003.61 251 001.8625 500.9312 750.4156 375.278 187.639 093.819 546.99 773.44 886.72 443.31 221.7610.8305.4152.776.3正變形極值/%0.109 0460.109 0160.514 2980.109 0460.514 2960.516 5950.517 3350.517 3350.517 3350.517 3350.517 3350.517 0530.516 8620.516 7580.516 7050.516 666負變形極值/%-5.68 488-11.129 7-15.942 1-18.648 9-20.075 2-20.806 3-21.176 3-21.549 3-21.455 8-21.409 1-21.385 7-21.025 6-21.31-21.452 8-21.524 4-21.523 7

表選項

2.3.2.3 剖分面片角度變形分形

剖分面片內側夾角,即為該四邊形頂點所對應的球面角,計算方法見式 (4)。圖 8表 6是部分等級剖分面片的角度變形分析結果。

圖 8 部分等級剖分面片角度變形分布圖 (與90°對比)Fig. 8 The distribution curve about deformation angle of some subdivision cell (compare with 90°)

表 6 部分等級剖分面片角度變形統計表 (與90°對比)Tab. 6 The statistical list about deformation angle of some subdivision cell (compare with 90°)

標準邊長/m2 502 003.61 251 001.8625 500.9312 750.4156 375.278 187.639 093.819 546.99 773.44 886.72 443.31 221.7610.8305.4152.776.3正變形極值/%6.456 2124.942 1494.781 0481.557 3111.217 5921.035 4481.028 9261.025 6671.027 2961.028 1111.028 5181.024 8741.030 6361.033 5221.034 9661.034 936負變形極值/%-2.423 564-2.422 564-2.040 443-1.222 417-1.035 448-1.035 448-1.035 448-1.028 926-1.028 926-1.028 926-1.028 926-1.024 394-1.030 396-1.033 402-1.034 906-1.034 906

根據上述分析可知,採用正六面體剖分地球的方式,雖然有可能到一定尺度后的剖分面片為正方形、且不同層級之間的剖分面片具有較好的遞歸關係,但存在如下問題:

(1) 雖然正六面體的6個面為正方形,但根據其等分遞歸生成的剖分面片並不都是正方形,只是特定的層級上的剖分面片接近正方形;

(2) 剖分面片與傳統測繪中的經緯度並沒有很明確的關係,尤其是兩極剖分面片,如圖 4(a)中的⑤、⑥對應的區域。因為在0級剖分面片⑤和⑥區域,在細分到1級剖分面片時,雖然可以與地球±45°和±135°經線重合,然而在繼續按四等分剖分時,不再可能與經線、緯線重合。因此會增加與傳統測繪數據之間轉換的難度。

(3) 由於正六面體的8個頂點所在地球表面的緯度比較特殊,南北緯35.264 390°,則在與地圖分幅的對應關係上,必然會割裂該區域的地圖圖幅。

因此,基於正六面體的格網剖分框架並不能很好地滿足遙感數據統一組織框架設計原則。

2.4 GeoSOT格網

為了解決跨部門間、部門內各業務階段遙感影像數據組織基準不統一的問題,以北京大學程承旗教授為首的研究團隊提出了GeoSOT剖分網格[

24

-

26

]。GeoSOT網格首先將地球表面空間擴展為512°×512°,面片中心與赤道和本初子午線的交點重合。在此基礎上遞歸四叉剖分至1°剖分網格。然後將1°網格擴展為64′,再對其遞歸四叉剖分至1′剖分網格。最後將1′剖分網格擴展至64″,再遞歸四叉剖分至 (1/2048)″(如

圖 9

)。在此基礎上利用64 bit的反Z序對稱編碼對每個剖分面片進行唯一標識

圖 9 GeoSOT網格多級剖分示意圖[24]Fig. 9 Multi-level subdivision of GeoSOT[24]圖 10 GeoSOT網格單元的編碼[24]Fig. 10 GeoSOT cells' codes[24]

通過GeoSOT劃分框架可知,GeoSOT網格面片不僅可實現全球覆蓋而且可對每個剖分面片進行唯一標識。但由於GeoSOT網格是基於等經緯度進行劃分,GeoSOT網格剖分面片沒有角度變形,其邊長長度沿經線方向沒有變形 (圖 11(b))、沿緯向方向隨著緯度增加逐漸變小 (圖 11(a)),剖分面片的面積也是隨著緯度的增加逐漸變小 (圖 12)。其中,剖分面片面積計算公式與式 (3) 相同。

圖 11 GeoSOT網格剖分面片長度變形分布圖Fig. 11 The distribution curve about side of GeoSOT subdivision cells圖 12 GeoSOT網格剖分面片面積變形分布圖Fig. 12 The distribution curve about area of GeoSOT subdivision cells

且經過研究證明,GeoSOT網格與Google Earth、Worldwind、Bing Maps和天地圖等格比較,GeoSOT網格與測繪數據地圖圖幅具有很好的同構性[

24

],更適合對傳統空間數據和遙感數據進行組織管理。

綜上,GeoSOT網格符合全球遙感數據統一組織框架設計的各項原則,可用於遙感影像的統一組織。

3 結論

通過對等距、等緯差、正六面體等幾種剖分格網方案分析,採用矩形或正方形理論上不能無縫覆蓋地球表面或者其剖分面片編碼與傳統測繪數據之間進行快速地轉換,而GeoSOT剖分網格能夠很好地與測繪數據地圖圖幅框架同構,可作為全球遙感影像統一組織框架。

【引文格式】宋樹華,董芳,萬元嵬,等。 幾種地球剖分網格分析與初探[J]. 測繪學報,2016,45(S1):48-55. DOI: 10.11947/j.AGCS.2016.F006

圖片| 太空俯瞰地球,衛星視角下的農田竟如此壯觀!

權威 | 專業 | 學術 | 前沿

微信投稿郵箱 | [email protected]

歡迎加入《測繪學報》作者QQ群: 297834524

進群請備註:姓名+單位+稿件編號



熱門推薦

本文由 yidianzixun 提供 原文連結

寵物協尋 相信 終究能找到回家的路
寫了7763篇文章,獲得2次喜歡
留言回覆
回覆
精彩推薦