張 偉,于 靜,陳儒敏,張鴻博
(北京科技大學天津學院信息工程系,天津 301830)
關鍵字:便攜定位裝置、逆地理編碼算法、傳統搜索算法、k-d樹算法
當前,隨著移動互聯網基建設施的逐步完善和智能手機、平板等終端設備的普及,基于智能終端的APP 應用程序獲得爆發式增長,其中絕大部分業務都需要使用設備的地理位置信息。同時,由于國內經濟的發展和人民生活水平的提升,私家車駕車出行的比例大幅上升。駕駛員普遍應用高德地圖、騰訊地圖等移動導航APP 軟件。通過使用此類軟件,司機往往能夠快速知曉車輛所處具體位置,獲得其經緯度坐標,并能夠進行路徑規劃等高級操作。
然而,這類應用有一個弊端,即必須在聯網的情況下才能使用。目前,在我國通訊基站實施覆蓋的平原城市、鄉村等地區,實時使用此類軟件的問題不大;但在突發或者特定情況下,通信網絡的獲取比較困難:比如,在地質勘探工作者進入深山和原始森林進行科學考察時,或者突發地質氣象災害時,由于不具備通訊基站等通訊設施,這類導航軟件往往不能正常使用。
針對這類無網絡情況下的定位需求,本文提出了一種基于k-d 樹的逆地理編碼算法,此算法應用場景定位于便攜式離線定位裝置上。同時,編寫基于歐氏距離和半正矢公式的兩種傳統搜索算法代碼,作為對照組。通過測試,不僅驗證了k-d 樹算法的正確性,而且證明在較大測試集的情況下,此算法運行效率比傳統算法高出近300倍。
逆地理編碼算法,亦稱反地理編碼算法,或地址解析算法,是指從已知的經緯度坐標推算出對應的地址描述的轉換算法。一般用于根據定位點經緯度來獲取該處的位置詳細信息,在地質考察、野生動物保護、軍事等領域有著大量的應用場景。目前,大部分系統和項目開發中,逆地理編碼算法都是通過調用互聯網廠商提供的逆地理編碼服務接口來實現的。國內市場占有率較高的廠商,例如百度地圖和高德地圖,均推出了基于開放地圖服務的地圖API接口。文獻[1]中詳細闡述了如何通過編寫JavaScript 腳本程序調用高德開放地圖接口,輸入經緯度得出地理信息的過程。國外的諸多網站,如PickPoint、GeoNames、MapQuest等服務商也提供類似的服務。
然而,在上述提到的某些應用場景中,離線逆地理編碼算法的實現又是不可回避的要求。目前,國內外對逆地理編碼具體算法實現及流程的研究很少,大多數只是調用現成的程序接口或者封裝好的模塊。高德集團的一份專利中,提及了離線的逆地理編碼的方法及其裝置,不過其側重于設備端的功能設計及實現,對具體的算法實現流程并未進行深入的闡釋,可參考價值有限。文獻[3]提出了一種基于GeoHash算法的分層調度方法和區間調度模型,可對共享單車的調度方案進行有效地規劃,通過對比已有的分層調度方案,得出了算法收斂速度快、模型時間復雜度低的結論。但是,該文章直接調用封裝好的第三方GeoHash計算模塊,未對其實現機理做詳細研究。針對上述不足,本文將逐步通過問題分析、測試集收集、測試組/對比組代碼測試等流程來研究這一問題。
離線逆地理編碼算法的流程,簡單來說,即是在未聯網的情況下,從便攜式離線定位裝置的GPS 接收器獲取一個或多個經緯度數據,再由算法與內置數據集進行搜索和比對,得出距離所測數據點最近的內置地址數據,然后將結果數據以地理名稱的字符串形式返回。該計算流程涉及兩方面關鍵問題,即:①內置數據集的來源、數據格式和參考點位的數量級;②搜索算法的計算效率和正確性。下文將逐一進行簡述。
本文算法的應用范圍先限定為國內區域,數據集采用國家標準行政地理區劃中的國內行政區劃不同粒度范圍地名劃分機制,具體信息如表1 所示。

表1 國內行政區劃分級劃分機制
根據表1定義本文所用的地址數據結構,如表2所示。

表2 地址數據結構
其中,scale_int 取值范圍為1~6,表征該條數據的地址級別;lon_f、lat_f,使用浮點數記錄地址的經緯度信息;add_1~add_6,記錄地址要素的字符串,比如add_3取值街道名稱,上限為256個字符。
本文的數據集取自國家統計局2020 年公布的數據,分別獲取3/4/5 級行政區劃數據用于后續測試。經過清洗后的數據個數如下:3級數據3552 個,4 級 數 據42358 個,5 級 數 據758050個。從3 級數據中隨機抽取296 個作為測試數據,以便驗證算法的正確性。因此,3級行政區劃數據用于內置數據集的個數變為3256 個。下文中,主要應用3 級行政數據驗證算法的正確性,4級、5級行政數據測試算法的執行效率。
在計算給定坐標點具體屬于哪個地理區劃時,最直接的辦法是計算此點與內置數據集全部參考數據點的距離,并排序選出最小的點,以此作為所屬地理區塊劃分。當然,考慮到行政區劃邊界不是嚴格的與經緯線平行的直線,還需要綜合附近多點的計算結果進行判斷。本文將問題簡化,以距離計算為基本算法進行逆地理編碼研究。
圖1 中,點是已知點,點是所求點,點是地球球心,、兩點間的弧線表示球面距離。

圖1 球面兩點間距離計算示意圖
計算距離的算法有多種,一為歐拉公式(Euclidean Distance),計算與點的直線距離,如公式(1)所示。

考慮到地球是個球面,第二種采用半正矢公式(Haversine Formula)計算球面上兩點間的距離,推導過程如下所示。
對于球面上任意兩點,圓心角的半正矢值()可用公式(2)進行計算。

公式(2)中,代表弧度制圓心角;代表地球半徑,取值6378 千米;(lon,lat) 和(lon,lat)分別代表球面上任意兩點和的弧度制經緯度坐標,其值采用WGS 84 標準GPS 坐標,取自便攜式裝置中的GPS 定位模塊。半正矢公式hav()的具體計算過程見公式(3)。

將公式(3)帶入公式(2)第二個等號右邊式子,得到公式(4)。

反半正矢函數的計算公式如公式(5)所示。

將公式(4)和公式(5)聯立,得出距離的計算公式如下。

本文采用的傳統搜索算法主要通過計算測試點與內置數據集各點的距離,選擇數值最小的點作為輸出結果。為了加快搜索速度,將經過Euclidean 或Haversine 公式計算后的距離集合進行排序,比較大小時采用傳統的二分法進行搜索。詳細步驟如下:
(1)算法初始化,將測試數據與內置數據集讀入內存;
(2)利用Euclidean 或Haversine 公式計算內置數據集各點與測試數據間的距離,得到距離集合并排序;
(3)采用二分法搜索距離最小的點;
(4)結果返回對應內置數據集,輸出地址字符串;
(5)算法結束。
應用5級行政區劃數據進行測試。為排除計算過程中的隨機影響,每次重復運行100次,取平均計算結果作為記錄值。表3 中,傳統1 和傳統2 分別對應Euclidean 和Haversine 計算公式。傳統1 代表采用Euclidean 公式計算單次運行時間,傳統2 代表采用Haversine 公式計算方法單次運行時間。

表3 傳統算法在不同行政級別數據集情況下的計算情況
從表3可以看出,整體來看,隨著內置數據集的增大,兩種傳統算法的單次計算運行時間逐漸增長。分別來看,采用Euclidean 公式的執行時間稍短,但是正確率低于采用Haversine 公式的傳統2 算法。經手動排查后,發現傳統1 算法判斷失誤的測試點與相鄰最近的內置數據集點間距離較遠,歐拉直線距離與球體表面的圓弧距離相對誤差較大。總的來說,采用傳統算法執行效率比較低,單個測試點的查詢時間在4秒鐘以上,因此在多數場景應用中是無法滿足時效要求的。鑒于采用Haversine 公式的正確率較高,故選其作為下文的對照組算法。
k-d 樹,即k-dimensional tree,常用作空間劃分及近鄰搜索算法,是二叉空間劃分樹的一個特例。通常,對于維度為、數據點數為N 的數據集,k-d 樹適用于?2的情形。以采用3級行政區劃做內置測試集為例,維度=2,數據點數=3256,滿足3256?4 的應用條件。k-d 樹算法可以分為兩部分,一是為構建k-d 樹本身這種數據結構而建立的算法,即如何把包含大量數據的內置數據集構建成一棵k-d 樹;二是如何在建立的k-d 樹上進行最鄰近查找的算法,即查找k-d樹上距離測試數據最近的節點。
構建內測數據集經緯度的k-d 樹算法簡述如下:循環依序取數據點的經緯度數據作為切分維度,分別取數據點在該維度的經緯度中值作為切分超平面,將中值左、右側的數據點分別掛在其左子樹和右子樹。遞歸處理其子樹,直至所有數據點掛載完畢。下面以3級行政區劃數據集中的甲~庚共7個數據點為例,只保留其1級行政區劃名稱及經緯度數據;目標數據點是地處河北省石家莊市橋西區的新百廣場商業綜合體。具體數據見表4。

表4 簡化k-d樹算法所需數據點
具體步驟如下:
(1)構建根節點時切分維度為緯度,上述7個數據點按緯度從小到大進行排序,取中值乙點天津(39.14393,117.21081)為根節點;
(2)石家莊、邯鄲、邢臺三點在根節點的左半子樹,北京、唐山、秦皇島三點在根節點的右半子樹;
(3)構建根節點的左子樹時,切分維度為經度,中值點邢臺(37.06953,114.52048)作為分割平面,邯鄲為左子葉、石家莊為右子葉;
(4)根節點的右子樹構建方法如上。
至此,k-d 樹構建完成,如圖2 所示。圖中,每個數據點上的“經-N”或“緯-N”代表劃分樹或子葉時,以本點的經緯度作為切分維度的劃分標準。

圖2 7個數據點k-d樹的構建示意圖
上述構建過程結合圖3,可以看出,構建一個k-d 樹即是將一個二維平面逐步劃分的過程。分割所用的超平面即是各個數據點的某一維度值。圖3背景中的各條實線,即為分割所用的超平面,為方便辨識,經度線用粗實線標注,緯度線用細實線標注。

圖3 k-d樹搜索算法的“回溯”計算示意圖
k-d樹的搜素算法步驟如下:
(1)首先對k-d 樹做深度遍歷,取出距離最近的點,以圖4標注的①→②→③為確定的搜索方向;
(2)為了明確丙葉子節點是否為真正的最鄰近點,還需要進行“回溯”操作:具體為算法沿搜索路徑反向查找是否有距離查詢點更近的數據點;在圖4 中,即為反向從③→④→⑤的過程;

圖4 k-d樹搜索算法的實現
(3)“回溯”算法實現由圖3 所示,圖中的左下角實心五角星為需要查詢的目標地址。圖3中的“經-1/緯-3”分別為在丙點按經度、在己點按緯度做分割用的兩個超平面。計算目標地址與丙點間的距離,以目標地址為圓心,距離為半徑做圓。可以看到,圓與丙點的上一節點(己點)所確定的“經-1”存在交點,說明己點所在的左半樹需要被“回溯”計算,完成過程④。而唐山所處的整棵樹的右半樹,由于未與所做出的圓出現交點,故不需要被計算,圖中用陰影覆蓋以示區分。同理,二級回溯完成過程⑤。
(4)再比較己點、庚點和兩點與目標點的具體距離,選出距離目標點最短的那個點,k-d樹搜索算法完成。得出丙點是所查詢地址的結論,即新百廣場位置在石家莊。
(5)以上步驟演示的是需要進行“回溯”計算的情況。如圖3中右下方,虛線圓圈中的實心三角為目標點的另外一種情況:當目標點與最近的內置數據點之間所作的圓,與其他超平面無交點時,則不必進行“回溯”操作。
應用上述的k-d 樹算法編寫測試程序,在其他軟硬件測試條件與第4節一致的情況下,得到的測試結果如圖5所示。圖中橫坐標刻度值為log 值,測試集數量在第一個測試點為296 個,隨后每次增加10 倍,最后一個測試點的數量達到2.96*10個。折線圖每個點附近的數值代表其算法執行時間。

圖5 k-d樹算法在5級行政區劃下的計算結果
在第一個測試點,同樣輸入296個測試數據和5 級行政數據的情況下,k-d 樹算法的執行時間僅為5.213 秒,換算成單點計算時間為0.017611 秒,相比于表2 中的傳統2 算法4.851秒,算法效率提升了274.45 倍。同時,算法正確率達到100%。
針對離線逆地理編碼算法的研究,本文首先采用兩種傳統搜索算法進行測試,并對由半正矢公式進行球面距離計算的過程進行詳細推導,通過計算得到了傳統算法時效性低的結論。然后,通過搭建7個數據點的簡單k-d 樹模型與算法演示,闡述了k-d 樹算法的基本原理與計算步驟。測試結果表明k-d 樹算法能夠大幅提高計算效率和準確率。