郭 霄 薛 睿 武慧琳 荀張媛
(1. 61363部隊, 陜西 西安 710000; 2. 自然資源部第一大地測量隊, 陜西 西安 710054;3. 自然資源部測繪標準化研究所, 陜西 西安 710054; 4. 中煤航測遙感集團公司, 陜西 西安 710199)
標準經緯度分幅地形圖接圖表作為一種項目管理不可缺少的重要數據,在測繪生產領域有著十分重要的作用。當前地方數據格式管理不一,實際工作中可能遇見數據分別采用新舊編號規則管理的圖幅,導致在數據匯交管理時,需重新修改接圖表數據屬性或重新生成接圖表,現有的制作接圖表方法操作復雜效率低。因此本文提出一種快速生成接圖表的算法就具有十分重要的意義。
目前主流制作接圖表方法依據目的分為兩類:一是通過經緯度范圍生成接圖表,主要通過ArcGIS中ArcToolbox提供的“創建漁網”工具,生成接圖表的骨架,使用者需要提供漁網的起始樞軸點,要生成對應比例尺的經差與緯差等參數共11項[1]。該方法處理不夠自動化,圖號等屬性需要后續通過字段計算器添加,費時費力,對操作者專業水平要求較高,容易出現操作失誤導致的漁網位置錯位等問題。二是通過圖號創建接圖表,主要通過Arcpy實現,分為3步:①根據文件名獲取圖幅號列表;②按照圖幅名計算坐標;③依照計算出的坐標生成接圖表多邊形,主要用到的函數有arcpy.Point()、arcpy. Array()、arcpy. Polygon()[2]。此方法一定程度的提高了生成接圖表的效率,但后續依舊需要通過字段計算器添加圖號信息,不能一步完成。
為更好地解決以上問題,本文提出了一種新的算法模型,該模型在部分以往的計算基礎上,借鑒了圖像處理深度優先搜索(depth first search,DFS)算法邏輯,減少了操作者干預接圖表骨架的生成,提高了準確度,將圖號信息在算法中直接添加,避免了后續操作,并將該算法與傳統方法進行對比實驗進行驗證。
1.1.1地形圖經緯度分幅
地形圖經緯度分幅指地圖的圖廓由經緯線構成,我國的基本比例尺就是以經緯線分幅制作的[3]。它們是以1∶1 000 000地圖為基礎,按規定的經差和緯差劃分圖幅,行列數和圖幅數成簡單的倍數關系[4]。
1.1.2新舊圖幅編號
新舊圖幅編號規則都是在標準經緯度分幅的基礎上設計,新圖幅號由10位碼組成,以1∶50 000圖號為例,J50E001002中J代表該圖幅所屬1∶1 000 000圖幅行號數字碼、50代表其所屬1∶1 000 000圖幅列號數字碼、E代表所屬比例尺代碼、001與002分表代表圖幅行號與列號數字碼[5]。舊圖幅編號系統按圖1編碼。
所有比例尺都可實現新舊圖號轉換,以1∶50 000比例尺為例,其對應公式如表1所示,其中H代表新圖幅中的行號,L代表新圖幅中的列號,X代表舊圖號中的地圖代碼值,[]代表取整,()代表取余。此為靈活生成不同編碼規則的圖號接圖表奠定了算法基礎。

圖1 我國舊式地形圖編號系統

表1 新舊圖號轉換
本文提出算法將根據經緯度范圍生成接圖表與根據圖號生成接圖表結合,整體流程如圖2所示。圖號文件可根據需求輸出為新編碼規則圖號與舊編碼規則圖號txt文件。
1.3.1深度優先搜索(DFS)[6]
DFS是計算機術語,是一種在開發爬蟲早期使用較多的方法,是搜索算法的一種。它的目的是要達到被搜索結構的葉結點。
深度優先遍歷圖的算法是,假定給定圖G的初始狀態是所有頂點均未被訪問過,在G中任選一個頂點i作為遍歷的初始點,則深度優先搜索遞歸調用步驟[7-8]。①訪問搜索到的未被訪問的鄰接點;②將此頂點標記為已訪問節點;③搜索該頂點的未被訪問的鄰接點,若該鄰接點存在,則從此鄰接點開始進行同樣的訪問和搜索,反復進行直到所有節點都被訪問為止[9]。

圖2 整體算法流程
1.3.2計算圖號
本文將經緯度范圍圍成的矩形四角頂點類比作圖G的頂點i,這里稱為基點,以基點為根結點搜索與其相鄰的圖幅,相鄰圖幅作為其子節點[10-11]。
當輸入數據為經緯度范圍時,根據下式計算四角其對應的新式圖號編號,a為1∶1 000 000圖幅所在經度帶的字符所對應的數字碼;b為1∶1 000 000圖幅所在緯度帶的字符所對應的數字碼[12];c為所求比例尺地形圖在1∶1 000 000地形圖編號后的行號;d為所求比例尺地形圖在1∶1 000 000地形圖編號后的列號;λ為某點的經度;φ為某點的緯度;Δλ為所求比例尺的經差;Δφ為所求比例尺的緯差[13-14]。
以矩形4個交點圖號為控制搜索生成邊框圖號,當搜索結束后,計算四角點圖號內對角圖號,再以內對角圖號為基準再次進行搜索。當搜索進行到最后一個循環時會有6種情況,分別為(a)田字分布(b)豎向一字分布(c)橫向一字分布(d)單幅分布(e)豎向兩行分布(f)橫向兩行分布。
針對(a)情況,當四角圖號的對角圖號為其本身,則跳出搜索;針對(b)、(c)情況以兩邊圖幅為基準向下或向右搜索圖號至完畢;(d)中當四角圖號的對角圖號為同一圖號時跳出循環,結束搜索;(e)(f)判斷向下或向右搜索是否右圖幅如果沒有,結束搜索。
在已知新圖號的前提下,通過(5)(6)式可計算出該圖幅在西南圖廓點的經緯度坐標[4]。
通過該比例尺所屬經差與緯差可計算出其他3個角點的經緯度坐標,獲取的XY坐標,寫出文件時將四角點坐標與其對應圖號一起寫出,為后續直接在接圖表面數據中增加圖號信息做準備。根據實際需求轉換圖號編碼實現需要圖號的直接輸入。
將XY坐標轉為點要素并通過spatRef = arcpy.SpatialReference(4326)設置其空間參考,再將點要素生成線,線轉面時將之前寫出的圖號寫入面要素屬性中。完成接圖表制作,具體流程如圖3所示。

圖3 通過圖號計算接圖表算法
2.1.1計算對角圖號實現
搜索計算四角圖號由4個主要函數組成,分別是ldown(string number)計算基準圖幅左下角圖號;lup(string number)計算基準圖幅左上角圖號;rdown(string number)計算基準圖幅右下角圖號;rup(string number)計算基準圖幅右上角圖號。每一個計算方法都有4種情況:(1)a變b不變;(2)a不變b變;(3)a、b同時變化;(4)a、b都不變,最后返回計算圖號。
2.1.2循環獲取邊框圖號
循環獲取圖號主要借助toright(string num, string destination);todown(string num,string destination),采用while(true)死循環,當搜索圖號num與目標圖號destination相等,或num圖號排在destination下方和左方時循環跳出,搜索結束。
最后一個循環借助hashset結構,在預先存儲了搜索基準圖號的前提下,hashset的長度為4,hashset主要特點為無序的不重復,因此借助hashset可以了解最后一個循環屬于上文6種情況中的哪一種,利用while(true)死循環,當為(a)情況時添加其對角圖號后hashset長度為仍為4,跳出循環;(b)(c)情況時hashset長度為6,借助todown(string num,string destination)、toright(string num, string destination)完成搜索后跳出;(d)情況時長度為5,break跳出循環。
所有圖號最終都存儲在hashset中通過foreach循環遍歷借助stream類寫出。通過python實現角點圖號與邊框圖號函數編寫,并將腳本插入Arctoolbox中,方便加入modelbuilder使用。
2.1.3圖號生成接圖表實現
借助modelbuilder,實現圖3流程如圖4所示。其中“計算xy坐標”與“xy轉為點要素”工具為插入的腳本工具。

圖4 Modelbuilder實現算法
2.2.1參數對比
本文計算經緯度范圍內的圖號,只需輸入上下緯度與左右經度值,與保存路徑。也可以將該腳本拖入的Model中,一步實現生成。
依據圖號生成接圖表界面參數包括存儲圖號文件的路徑,判斷是否為一條線上的點的字段name,以及輸出路徑。
而傳統方法所需參數繁雜,后續還需要在屬性表中創建字段,使用字段計算器,計算幾何還有投影設置等設置接圖表圖號、投影屬性[15]。
2.2.2實驗結果對比
以生成1∶50 000數據接圖表為例,假定使用傳統方法熟練的情況下,與本文方法對比如圖5所示,其中粗線為本文方法,細線為傳統方法。
從圖5可以看出傳統方法人工操作執行時間受操作者的影響較大,且基本運行時間是本文方法的1~2倍。

圖5 整體對比分析
本文借鑒深度優先搜索算法邏輯,類比接圖表格網為柵格數據的像素格網,以接圖表四頂點開始搜索經緯度范圍內的所有圖號,并在通過圖號生成接圖表的過程中將圖號屬性寫入接圖表矢量數據中,與現有制作接圖表的方法相比減少了處理步驟與人工操作,克服了因人工誤操作導致的誤差問題。通過多次對比實驗分析證明,該方法大大提高了接圖表生成的效率,有效地降低了制作的復雜度。但該方法還存在識別圖號格式不夠靈活的限制,例如在圖號中如果出現其他連接符需更改為標準圖號格式算法才可使用,未來將對算法進行改進使其更加靈活適應圖號的不同格式,擴大其適用領域。