吳 凱 馬藝文

(1. 濟(jì)南市勘察測繪研究院, 山東 濟(jì)南 250013; 2. 武漢市國土資源和規(guī)劃信息中心, 湖北 武漢 430014)

0 引言

在現(xiàn)實(shí)世界中,兩個(gè)目標(biāo)物之間可以直接通行的情況較少,一般都需要繞過若干障礙物才可到達(dá)。地理空間問題中有很多典型的障礙空間問題,例如:軍事問題中如何選擇合適的進(jìn)軍路線以避開難以逾越的地理障礙及人為障礙[1];日常生活中交通工具繞開擁堵路段或施工路段的最短路徑規(guī)劃[2];災(zāi)害救援中“安全島”問題[3];無人駕駛飛機(jī)的防撞避障和導(dǎo)航設(shè)計(jì)問題[4];科學(xué)研究中機(jī)器人避障路徑規(guī)劃問題[5];迷宮求解問題[6-7];障礙空間最近鄰查詢問題[8]。障礙空間,指空間中具有各種障礙物,此處的障礙物并非單指點(diǎn)或線或多邊形,而是指全形態(tài)圖形,它是點(diǎn)、線、面以及它們的組合,障礙本身不計(jì)算距離,并且不傳遞距離。

距離變換是計(jì)算并標(biāo)識(shí)空間點(diǎn)(對參照體)距離的變換或過程。障礙距離變換,即在空間中存在障礙物時(shí),距離的傳播需要繞過障礙物,它屬于條件距離變換[9-10]。這類問題比較抽象卻來源于實(shí)際地理世界的需要,障礙距離變換是解決此類地理空間問題的有效途徑。本文主要對比了兩種障礙距離變換柵格生成方法,并進(jìn)行了精度分析和對比。

1 障礙距離變換(DTO)

距離變換是將包含實(shí)體特征和空間背景兩種像元的二值圖像轉(zhuǎn)變?yōu)榫嚯x圖像的變換。在距離圖像中,每一個(gè)像素值表示該像素到其最近的一個(gè)實(shí)體像素的距離,具體體現(xiàn)為每個(gè)實(shí)體的距離波不斷地往外空間進(jìn)行擴(kuò)張,直到與鄰近實(shí)體的距離波相遇。距離變換的結(jié)果往往是空間等距離線,任一點(diǎn)的趨源梯度線(即該線上任一點(diǎn)的切線方向都垂直于過該點(diǎn)的等距離線),即該點(diǎn)趨源的最短路徑。

當(dāng)考慮障礙空間問題時(shí),上述距離變換則應(yīng)擴(kuò)展為DTO(Distance Transformation with Obstacles),即在障礙空間中進(jìn)行距離變換?？臻g中有生成元還有若干障礙,生成元傳播的距離波需要繞過障礙進(jìn)行傳播。由以上定義可知,一般意義上的距離變換可以看作是障礙物個(gè)數(shù)為零的障礙距離變換。

2 地圖代數(shù)障礙距離變換算法(MA-DTO)

地圖代數(shù)的障礙距離變換(MA-DTO)在二維障礙空間中是通過柵格路徑距離變換實(shí)現(xiàn)的?？筛鶕?jù)精度要求確定開窗大小,采用3×3,5×5等奇數(shù)(2n+1)×(2n+1)模板。模板柵格路徑定義為模板的中心至其他各個(gè)柵格中心間的有限個(gè)方向線段的組合,模板柵格路徑所對應(yīng)的長度稱為模板柵格路徑長度。

3 ArcGIS障礙距離變換的實(shí)現(xiàn)方法(ArcGIS-DTO)

本文主要利用ArcGIS的成本距離(Cost Distance)工具,來計(jì)算每個(gè)像元到成本面上最近源的最小累計(jì)成本距離,即ArcGIS-DTO。在使用Cost Distance工具進(jìn)行障礙距離變換時(shí),生成元(Source Data)可以為矢量數(shù)據(jù)也可以為柵格數(shù)據(jù),如果是柵格數(shù)據(jù)則具有 NoData 值的像元不包括在源集內(nèi),成本柵格數(shù)據(jù)(Cost Raster)中的NoData部分充當(dāng)成本面中的障礙。

成本距離的算法如下:像元八鄰域從像元0移動(dòng)到四個(gè)與其邊連接(1、3、5、7四像元)的近鄰之一,取像元1時(shí),跨越連接線移動(dòng)到相鄰結(jié)點(diǎn)的成本為a1=(cost0+cost1)/2,其中:cost0代表像元0的成本,cost1代表像元1的成本,a1為從像元0到像元1連接線的總成本。如果沿對角線移動(dòng),例如從像元0到像元8,則連接線上的行程成本為a2=1.414 214(cost0+cost8)/2,其中cost0代表像元0的成本,cost8代表像元8的成本,a2為從像元0到像元8連接線的總成本，如圖1所示。

圖1 像元八鄰域

4 兩種DTO算法的精度對比

4.1 理論分析

MA-DTO算法是一個(gè)統(tǒng)一的、通用的算法,具有優(yōu)良的數(shù)據(jù)初始化、可視化性能,并且障礙物和發(fā)生元為全形態(tài)圖。柵格空間下兩定點(diǎn)之間的距離由其所屬的兩柵格中心間的距離定義。在障礙空間下由于連通概念的改變,障礙空間下的柵格路徑定義為指定兩定點(diǎn)a、b所屬柵格中心間并由相通柵格中心間線段組成的有序線段集。柵格路徑所對應(yīng)的累計(jì)長度稱為柵格路徑長度。一般情況下,兩定點(diǎn)之間的柵格路徑并不唯一,往往有多條,所有路徑中最短累計(jì)長度的柵格路徑最有意義。實(shí)際上基于地圖代數(shù)的距離變換給出了整個(gè)障礙空間所有點(diǎn)的趨源距離,并且給出了多生成元的整個(gè)障礙空間所有點(diǎn)的趨最近元距離。

ArcGIS-DTO的實(shí)現(xiàn)工具,成本加權(quán)距離通過將距離等同為成本因子來修改“歐氏”距離,該成本即為經(jīng)過任何指定像元所需的成本。成本分配工具可以基于累積行程成本來識(shí)別最近(或成本最小的)源的像元。成本回溯鏈接工具將提供一個(gè)道路地圖,用于標(biāo)識(shí)從任何像元開始,沿著最小成本路徑,返回到最近源的路徑。

4.2 可視化對比

下面以一個(gè)點(diǎn)的距離變換為例,對比其距離波的形態(tài),如圖2所示,圖2(a)為一個(gè)點(diǎn)的歐氏距離變換圖,其產(chǎn)生的距離波為同心圓;圖2(b)為ArcGIS中Cost Distance所使用的3×3模板距離變換,可以看出由于3×3模板的八鄰域結(jié)構(gòu),生成的距離波是以等邊八邊形的方式逼近歐氏距離變換的圓形;圖2(c)為地圖代數(shù)方法一個(gè)點(diǎn)的距離變換結(jié)果,使用的是17×17模板,從距離波的形態(tài)上來看,已經(jīng)十分精確地逼近歐氏距離變換結(jié)果。

圖2 一個(gè)點(diǎn)的距離變換對比

用ArcGIS的距離變換結(jié)果疊加歐氏距離變換的同心圓,來觀察ArcGIS實(shí)現(xiàn)結(jié)果的誤差分布情況。如圖3所示,為ArcGIS的3×3模板距離變換結(jié)果。

圖3 誤差分布示意圖

疊加歐氏距離變換距離波同心圓可以看出,在八方向上,如線段OP所在位置,ArcGIS的結(jié)果和歐氏距離相等;在兩個(gè)相鄰八方向的中點(diǎn),如線段OR所在位置,誤差最大;誤差隨著向八方向的靠近而減小。設(shè)線段OP、OQ、OR所處位置處的誤差分別為ΔOP、ΔOQ、ΔOR,則ΔOP<ΔOQ<ΔOR。

以線狀生成元在障礙空間中的距離變換試驗(yàn)來進(jìn)行兩種DTO實(shí)現(xiàn)方法的精度對比,如圖4所示,圖4(a)中直線段生成元,螺旋區(qū)域內(nèi)部為距離變換空間,外部為障礙空間;圖4(b)為ArcGIS實(shí)現(xiàn)的DTO圖;圖4(d)為地圖代數(shù)方法實(shí)現(xiàn)的DTO圖。圖4(c)和圖4(e)分別為兩種方法實(shí)現(xiàn)結(jié)果的局部放大圖,放大部分為螺旋彎曲單元的轉(zhuǎn)彎區(qū)域。對比圖4(c)和4(e)可以看出:地圖代數(shù)方法的距離波更平滑、更細(xì)膩,而ArcGIS的距離波在轉(zhuǎn)彎區(qū)域容易形成折線,距離波紋不夠平滑。

圖4 線狀生成元障礙距離變換精度對比

4.3 量化對比

如圖5所示,圖5(a)中點(diǎn)為生成元,共7個(gè),線和面為障礙物;圖5(b)為ArcGIS實(shí)現(xiàn)的多點(diǎn)生成元DTO圖;圖5(c)為地圖代數(shù)方法實(shí)現(xiàn)的DTO圖。兩種方法所得距離波的形態(tài)不同之處前面已經(jīng)進(jìn)行了討論。在此,布設(shè)距離變換圖上n個(gè)隨機(jī)點(diǎn)作為檢查點(diǎn),通過對比它們所在位置處兩種距離變換結(jié)果的距離值,分析距離差值的變化規(guī)律。

圖5 多點(diǎn)生成元DTO

在圖幅范圍內(nèi)隨機(jī)選取40個(gè)檢查點(diǎn),如圖6所示,三角形點(diǎn)在此用P1,P2,…,P40來表示。圖6(a)為檢查點(diǎn)與ArcGIS-DTO結(jié)果的疊加圖,圖6(b)為檢查點(diǎn)與MA-DTO結(jié)果的疊加圖。

圖6 隨機(jī)點(diǎn)與DTO結(jié)果的疊加

采集40個(gè)檢查點(diǎn)所在位置,ArcGIS-DTO和MA-DTO所得到的兩個(gè)距離變換結(jié)果圖中的距離值。如表1所示,D1為ArcGIS所得結(jié)果的距離值,D2為地圖代數(shù)方法所得結(jié)果的距離值。再求出兩種距離值的差值(ΔD),ΔD=D1-D2。對40個(gè)檢查點(diǎn)按照距離值從小到大的順序進(jìn)行排序,排序后的結(jié)果如表2所示。

表1 隨機(jī)點(diǎn)處距離值及距離差值 單位:m

表2 按距離值從小到大排序 單位:m

由于MA-DTO方法已十分接近于歐氏距離變換,所以本文將MA-DTO方法設(shè)為標(biāo)準(zhǔn),考察ArcGIS-DTO的誤差。則表2中距離差值ΔD可表示為ArcGIS-DTO的誤差值。利用誤差值(ΔD)和ArcGIS-DTO的距離值(D1)繪制成折線圖,即誤差隨距離值變大的變化趨勢圖。如圖7所示,圖7(a)中的0值點(diǎn)即為上述的八方向上的點(diǎn),在統(tǒng)計(jì)變化趨勢前我們先將0值點(diǎn)去除,統(tǒng)計(jì)其余37個(gè)點(diǎn)的趨勢圖,如圖7(b)所示。

圖7 ArcGIS-DTO結(jié)果誤差值隨距離值的變化趨勢圖

對圖7(b)所示的誤差值隨距離值的變化趨勢圖進(jìn)行線性回歸,如圖8所示,線性回歸函數(shù)為:

圖8 誤差隨距離值變化趨勢線性回歸

y=1.926 54+0.062 38x(x為距離值,y為誤差值),相關(guān)系數(shù)R=0.811 48?？梢钥闯?雖然距離差值有波動(dòng),但整體上來講,隨著距離值的增大,也就是距離生成元越遠(yuǎn),ArcGIS-DTO的誤差越大。

5 結(jié)束語

通過可視化和量化兩種方式的精度對比,本文得出以下結(jié)論:ArcGIS-DTO和MA-DTO在實(shí)現(xiàn)障礙空間的距離變換時(shí)都是模擬無障礙空間中歐氏距離變換的思想。地圖代數(shù)算法使用的模板大小可以為3×3、5×5、奇數(shù)(2n+1)×(2n+1),隨著n的增大誤差會(huì)越來越小,越接近歐氏距離,精度也越高,但算法的復(fù)雜性也會(huì)增加;ArcGIS中使用的算法是固定的3×3模板,所以靈活性比較差,精度相對也較低些。