基于天地圖的局部影像渲染技術研究及應用

余萬一，張 達，冀 虎

(1.礦冶科技集團有限公司，北京 100160；2.中國-南非礦產資源綜合利用聯合研究中心，北京 102628；3.國家金屬礦綠色開采國際聯合研究中心，北京 102628)

目前，礦山安全監測系統廣泛采用WebGIS來展示監測點的狀態及空間位置信息[1-2]，并基于天地圖提供的衛星影像，為用戶提供更加直觀的可視化分析應用[3-5]。但天地圖提供的影像底圖更新周期長，一般都是1～2年，甚至更長時間。而礦山的實際情況卻已經發生了改變，無法及時呈現在衛星影像上，導致不能很好地滿足礦山企業的需求。本文通過研究天地圖的圖文數據及其開放的API，通過無人機航拍礦區現場圖像，拼接制圖，替換天地圖衛星影像局部圖層，實現了礦區影像及時更新[6-8]。

1 天地圖

天地圖是由國家基礎地理信息中心研制的國家地理信息公共服務平臺，所提供的中國區域內地理信息數據資源是最全的。影像底圖可以達到街道級的顯示分辨率。天地圖服務采用OGC WMTS標準，支持HTTP和HTTPS協議，影像底圖的投影類型有經緯度投影和球面墨卡托投影。為了適應不同的業務場景，天地圖提供了多種服務調用方式API，包括：地圖API、網頁API、WEB服務API、數據API等[9]。

本文所做的分析都是基于地圖API和網頁API來進行的，通過分析地圖API加載衛星影像數據的原理以及網頁API提供的各項接口服務說明，利用無人機航拍的礦區影像替換原有的衛星影像底圖，來達到礦區現場影像實時更新的目的[10]。

2 技術研究

2.1 天地圖底圖

天地圖底圖可以選擇顯示地圖、影像、地形，本文采用的是影像底圖，而影像底圖的來源主要包括：天地圖官方提供的衛星影像底圖、第三方供應商提供的衛星影像底圖、自建衛星影像底圖。三種方式的優缺點詳見表1。

表1 天地圖衛星影像底圖來源方式比較

結合目前礦山安全監測實際應用需求及成本控制，最終采用天地圖官方提供的衛星影像底圖，然后通過無人機航拍來解決官方衛星影像更新慢的問題。

2.2 天地圖影像底圖加載原理分析

通過分析初始化地圖及渲染地圖的方式，可知影像底圖是不同的瓦片拼接而成的，每個瓦片都是256×256像素png圖片。縮放級別是從1到20，隨著縮放級別改變分辨率也有相應變化，而每一級加載的瓦片數也是不同的，推算出：瓦片數為22n-1，列數為2n，行數為2n-1，n代表的是縮放級別。坐標原點是西經180度，北緯90度。

2.3 天地圖影像底圖加載方式

天地圖默認是不加載衛星影像底圖的，需要自定義衛星影像底圖請求地址，拼接局部瓦片完成地圖圖層的渲染，具體加載方式可以分為以下幾步：

1)地圖服務選擇影像底圖圖層服務，投影類型采用球面墨卡托投影，服務采用HTTP的方式請求，影像底圖圖層請求示例：

http：//t0.tianditu.gov.cn/img_w/wmts?SERVICE=WMTS&REQUEST=GetTile&VERSION=1.0.0&LAYER=img&STYLE=default&TILEMATRIXSET=w&FORMAT=tiles&TILEMATRIX={z}&TILEROW={x}& TILECOL={y}&tk=密鑰

其中，參數z為縮放級別，參數x為切片的行號，參數y為切片的列號，密鑰在天地圖開發者中心申請的，擁有密鑰才能使用天地圖提供的服務。

2)配置化參數，配置的內容包括：初始化縮放級別、初始化定位的中心點經緯度、是否可拖拽、是否可縮放、操作控件的加載、地圖操作事件的監聽等。

3)通過初始化、縮放、拖拽可觸發可視區域的渲染，采用HTTP的方式請求天地圖服務器獲取圖層瓦片，如果存在瓦片數據則進行拼接，同時用衛星影像圖層渲染地圖，否則用異常圖片渲染地圖，提示該區域無影像。加載流程見圖1。

圖1 地圖加載流程圖Fig.1 Map loading flow chart

2.4 行列號計算

要完成替換，需要找到經緯度與影像底圖切片行列號的對應關系，通過分析天地圖影像底圖加載的原理及開放的API資源，推算出經緯度轉換為行列號的計算公式。

1)行號計算公式

式中，Row為行號，Lat為緯度，Res為相應縮放級別的分辨率，縮放級別與分辨率對應關系見表2。

2)列號計算公式

式中，Col為列號，Lng為經度，Res為相應縮放級別的分辨率，縮放級別與分辨率對應關系見表2。

表2 縮放級別與分辨率對應關系

2.5 航拍圖校準

航拍圖校準包括兩個方面，一方面需要確定航拍圖起始和結束經緯度，另一方面需要確定航拍圖的切片數量以及切片與行列號的對應關系。

通過分析網頁API，可以采用在影像底圖加載整張航拍圖的方式來確定起始和結束經緯度。方法如下：將無人機航拍的圖片處理成矩形，不規則部分用透明圖填充。圖片處理完成后，通過天地圖提供的添加圖片圖層功能，校準圖片在所在影像地圖區域的位置，確定航拍圖左下角和右上角的經緯度坐標，由經緯度轉換行列號算法計算出左下角和右上角所在的行列號，即可得到起始和結束的行和列進而推算出切片數量，并推算出航拍圖需要切割的切片和行列的對應關系。

1)切片數據量計算公式

T=[(Re-Rs)+1]×[(Ce-Cs)+1]

式中，T為切片數量，Re為末行值，Rs為首行值，Ce為末列值，Cs為首列值。

2)圖片編號與行列對應關系計算公式

Col=(Idxmod (Ce-Cs+1))+Cs

式中，Row為行號，Col為列號，Idx為切片編號，Re為末行值，Rs為首行值，Ce為末列值，Cs為首列值。

2.6 重新渲染影像底圖

采用專用切圖軟件的切片工具對航拍圖進行切圖，切圖大小為256×256像素的PNG圖片，不足256×256的部分用透明圖像填充。切圖排序采用原圖從左至右，從上至下，從0開始編號。通過圖片編號與行列轉換公式計算出行、列號，以此行列號命名圖片并存儲到本地，命名規則：行號_列號。

通過參考網頁API提供的天地圖圖層加載監聽事件，監聽加載的切片圖層行列號是否在航拍圖所在行列號區間內，對于在此區間內的切片圖層，利用本地存儲的影像進行替換，達到重新渲染圖層的目的。重新渲染的流程見圖2。

圖2 渲染流程圖Fig.2 Rendering flow chart

3 天地圖局部渲染技術應用

本文所述技術已經應用于某尾礦庫安全監測一張圖功能上，首先通過無人機航拍整個礦區，然后根據方案設計的算法計算出要替換的行列號和航拍圖切圖數量及與行列號的對應關系，最后用航拍圖庫區的衛星影像，以達到及時呈現尾礦庫現場情況的目的。下面以天地圖縮放級別16級為例來說明在尾礦庫安全監測系統上的應用方式及效果。

3.1 區域切片制圖

航拍圖使用的是該尾礦庫排水斜槽區域的影像圖，制圖步驟如下：

1)將該影像圖以圖片圖層的方式添加到天地圖，通過人工調整影像圖在天地圖上的空間位置，利用天地圖的坐標拾取功能，確定航拍圖左下角和右上角的經緯度坐標。得出左下角經度為117.783 690，緯度為29.028 330；右上角經度為117.791 520，緯度為29.033 130。

2)通過經緯度坐標計算得出左下角所在的行號為27241，列號為54209；右上角所在的行號為27240，列號為54210。可知首行值為27240，末行值27241，首列值54209，末列值54210。

3)通過首末行列值計算得出切片數量為4個，從左至右，從上至下，從0開始編號，得到編號值0、1、2、3。

4)通過首末行列值和圖片編號計算得出圖片編號與行列號的對應關系，見表3。

表3 圖片編號與行列號對應關系

5)利用專用切圖軟件提供的切片工具將該影像圖切成4個256×256像素的png圖片，分別命名為27240_54209.png、27240_54210.png、 27241_54209.png、27241_54210.png，并保存到本地資源文件夾下。最后影像圖的切片效果見圖3。

圖3 尾礦庫排水斜槽區域的影像圖切圖Fig.3 Image cutaway of tailings pond drainage chute area

3.2 區域局部渲染

天地圖提供圖層加載事件的監聽，通過監聽圖層加載，判斷該影像圖所在行列號是否在圖層加載切片的范圍內，如果在此范圍內則用本地資源下的影像切片替換原有衛星影像切片，不在此范圍內則用原有衛星影像渲染。相關示意代碼如下：

// 監聽圖層加載事件

map.getLayers()[0].addEventListener('tileloadstart'，function(params){

// 判斷行列號是否需要用本地圖片替換

if(coords.contains(params.coords)){

params.tile.src = '/resources/'+params.coords.x+'_'+params.coords.y+'.png' }})

在天地圖縮放至16級，該尾礦庫排水斜槽區域的拼接效果見圖4。

圖4 航拍影像拼接效果圖Fig.4 Mosaicking effect of aerial image

4 結論

基于天地圖的開放資源，通過分析地圖API和網頁API及其他相關資源，設計了相應的影像圖層局部渲染的技術，解決了尾礦庫區域衛星影像圖更新慢、無法及時呈現尾礦庫區域現有狀況的問題，并應用于礦山安全監測一張圖監測點的空間位置信息展示。通過現場應用的結果表明本技術切實可行，同時也為局部區域衛星影像實時性要求較高的業務領域提供了參考。