高 倩,章孝燦,胡 祺

(浙江大學(xué) 地球科學(xué)學(xué)院，浙江 杭州 310027)

柵格數(shù)據(jù)和矢量數(shù)據(jù)是地理信息系統(tǒng)中常用的兩種數(shù)據(jù)。柵格數(shù)據(jù)是一種離散的點(diǎn)陣數(shù)據(jù)，每一個單元的屬性值代表實(shí)體屬性，具有數(shù)據(jù)結(jié)構(gòu)簡單，易于進(jìn)行空間分析等優(yōu)點(diǎn)，但是數(shù)據(jù)存儲量太大，縮放尺度有限制。矢量數(shù)據(jù)由點(diǎn)的精確坐標(biāo)來表達(dá)點(diǎn)、線、面等地理要素，具有數(shù)據(jù)存儲量少，放大縮小不失真，精確定位等優(yōu)點(diǎn)。由于矢量數(shù)據(jù)有存儲空間小，便于觀察，具有空間對象，可以精確定位等較于柵格數(shù)據(jù)的優(yōu)點(diǎn)，在實(shí)際生產(chǎn)中常常需要將柵格數(shù)據(jù)矢量化。柵格矢量化通常分為兩種，一種是線狀數(shù)據(jù)的矢量化，比如等值線圖，另一種是圖斑的矢量化，比如遙感分類圖，專題圖等都需要柵格矢量化。本文主要討論圖斑矢量化后的問題。柵格矢量化后，由于柵格數(shù)據(jù)具有點(diǎn)陣特點(diǎn)，是離散性數(shù)據(jù)，矢量化后的數(shù)據(jù)呈現(xiàn)鋸齒狀，十分不美觀，與實(shí)際地物的連續(xù)邊界不符，且造成大量的數(shù)據(jù)冗余。為了矢量化后圖斑的美觀性、真實(shí)性，以及減少數(shù)據(jù)的冗余，節(jié)省數(shù)據(jù)的存儲空間，需要對矢量化后的數(shù)據(jù)進(jìn)行平滑。

目前多邊形平滑算法可以分為3類。一類是矢量抽稀算法[1-8]，將數(shù)據(jù)抽稀，以達(dá)到減少數(shù)據(jù)冗余且平滑的結(jié)果。主要方法有垂距法，光欄法，道格拉斯普克(Douglas-Peucker , D-P)算法，小波算法，遺傳算法等。其中D-P算法是最常用的一種，自提出以來不斷被改進(jìn)，如黃萬里等提出的基于面積保持的D-P算法[9]，顧騰等提出的D-P算法與Li-Openshaw結(jié)合改進(jìn)的方法[7]。雖然矢量抽稀算法能有效平滑鋸齒，消除數(shù)據(jù)冗余，但是往往需要一個適合的閾值，閾值的大小難以把握。第二類是曲線擬合算法[10-12]，如多項(xiàng)式擬合，B樣條擬合，貝塞爾曲線等，曲線擬合是平滑多邊形最常用方法。然而曲線擬合難以去除圖斑鋸齒，反而會出現(xiàn)“波浪”。第三類是綜合法，將多類方法結(jié)合使用，比如前面提到的D-P算法和Li-Openshaw結(jié)合的方法，D-P算法和B樣條結(jié)合的方法[13]等，此類方法因結(jié)合多種方法計算緩慢。

由于圖斑中鋸齒的存在，上述3種多邊形平滑方法無法較好地平滑圖斑。加之遙感圖斑通常具有實(shí)際意義，特別是土地利用類型圖斑的面積屬性，在實(shí)際應(yīng)用中有重要意義，平滑后圖斑與原始圖斑面積和位置的保持十分重要。此外,若相鄰圖斑公共邊處理不一致，容易出現(xiàn)“裂縫”等拓?fù)溴e誤的結(jié)果。因此，需要提出具有針對性的圖斑平滑方法。石軍南等人提出的矢量線條概率中值平滑方法[14]，宋正祥等人提出的顧及拓?fù)渑c尖角的矢量數(shù)據(jù)分組壓縮算法都是針對圖斑的平滑[15]，均達(dá)到了良好的平滑效果，但無法保持原始圖斑面積和位置。綜上圖斑的平滑存在3個問題: ①圖斑呈現(xiàn)鋸齒狀；②平滑時難以保持圖斑面積及位置信息；③平滑時圖斑公共邊易出現(xiàn)拓?fù)溴e誤。

針對上述問題，以矢量化后的圖斑為研究對象，考慮到相鄰圖斑的拓?fù)潢P(guān)系，以及圖斑平滑前后面積和位置保持的問題，提出基于弧段的遙感信息圖斑分段平滑方法。通過對矢量化后圖斑的弧段進(jìn)行分段后再平滑，達(dá)到保持其拓?fù)潢P(guān)系，也保持原始圖斑面積和位置的結(jié)果。

1 基于弧段的分段平滑方法

基于弧段的遙感專題信息圖斑的分段平滑方法利用深度搜索的方式提取弧段，防止后續(xù)平滑出現(xiàn)拓?fù)洳灰恢碌慕Y(jié)果，再根據(jù)鋸齒梯度對弧段分段，然后進(jìn)行平滑，以降低矢量化后數(shù)據(jù)平滑造成的面積和位置的偏差。

1.1 弧段提取

柵格矢量化后的數(shù)據(jù)沒有拓?fù)浣Y(jié)構(gòu)，平滑時相鄰圖斑的公共邊平滑不一致，容易產(chǎn)生如圖1所示的“裂縫”，“重疊”等拓?fù)洳灰恢碌慕Y(jié)果。但是在拓?fù)浣Y(jié)構(gòu)下，相鄰圖斑的公共邊共用一條弧段，僅平滑一次，則避免了平滑不一致產(chǎn)生的拓?fù)溴e誤。為消除平滑時容易出現(xiàn)的拓?fù)溴e誤，采取深度搜索的方式提取弧段，為后續(xù)平滑做準(zhǔn)備。

圖1 有公共邊的平滑

實(shí)驗(yàn)數(shù)據(jù)為常用ArcGIS的shp矢量文件，shp文件沒有拓?fù)浣Y(jié)構(gòu)，按數(shù)組順序存儲多邊形，多邊形外環(huán)的點(diǎn)順時針存儲，多邊形內(nèi)環(huán)的點(diǎn)逆時針存儲。針對shp文件的特點(diǎn)，以及深度搜索方法[2-3, 16]的特點(diǎn)，提出以下方面修正，使其更適用于shp文件，流程如圖2所示。

圖2 弧段提取流程

1)深度搜索采取“順序比較”和“逆序比較”的方法搜索公共邊。在shp文件中，公共邊在相鄰圖斑中必然是順序相反，因此只需要逆序比較。

2)搜索時起始點(diǎn)為第一個存儲點(diǎn)，不是弧段的端點(diǎn)。如圖3所示，圖斑AB在公共邊均有圖斑存儲的起始點(diǎn)與終止點(diǎn)，于是除了“逆序比較”找到end點(diǎn)，還要將start點(diǎn)向前移動。即A多邊形的點(diǎn)1，2，3，…，比較B多邊形i,i-1,i-2，…，找到end點(diǎn)，還需要A多邊形的點(diǎn)n，n-1，n-2,…,比較B多邊形i,i+1,i+2,…,找到真正的start點(diǎn)。

圖3 圖斑中點(diǎn)的存儲

3)搜索到一段弧段時立即平滑替換掉原始弧段，并標(biāo)記，再次搜索到該弧段時，不需要二次對比，直接跳過。

4)弧段僅需平滑一次，若是公共邊的弧段，在另一圖斑中只需將平滑后弧段逆序存儲替換掉原始弧段。

1.2 弧段分段平滑

柵格矢量化后，由于遙感圖像的點(diǎn)陣特點(diǎn)，柵格數(shù)據(jù)對邊界的表達(dá)有限，矢量化后的圖斑存在大量鋸齒。圖斑不美觀和數(shù)據(jù)冗余都由鋸齒造成，于是平滑圖斑只需對鋸齒平滑，即可達(dá)到平滑效果。采取將弧段分段為平滑鋸齒單元的方式找到需要平滑的鋸齒，采用鋸齒中線代替鋸齒的方式進(jìn)行平滑，達(dá)到既平滑圖斑，又保持圖斑面積和位置的結(jié)果。

1.2.1 平滑方式

如圖4所示AB兩圖斑，公共邊為鋸齒。提到柵格數(shù)據(jù)轉(zhuǎn)換為矢量數(shù)據(jù)時，線條形狀變化，其特點(diǎn)分布概率不確定，柵格數(shù)據(jù)矢量化的曲線是一條柵格線條的中線時方差最小。在圖斑中同樣，如圖4所示，鋸齒為AB兩圖斑的交接處，線1和線3包圍區(qū)域?yàn)殇忼X所屬范圍，同時為既有圖斑A也有圖斑B的混淆地帶。在線1和線3包圍區(qū)域內(nèi)的線都能將此區(qū)域等分為兩份來替代鋸齒，但是中間線線2為鋸齒邊最佳擬合的直線，其殘差最小。線2可使該交叉地帶一分為二，使斑AB位置不偏移，面積也能保持一致。不難發(fā)現(xiàn)，若需要平滑后AB兩圖斑面積保持一致，中線需平分鋸齒，為保證左右分到的面積相等，中線的起始點(diǎn)在水平線(垂直線)上，終止點(diǎn)也應(yīng)在水平線(垂直線)上。如圖4所示，可以保證AB圖斑面積不變，這樣中線左右兩邊的點(diǎn)數(shù)相同，面積也相同。

圖4 鋸齒邊界分析

1.2.2 弧段分段

1)平滑鋸齒單元。單純地使用中線進(jìn)行鋸齒平滑，保持原狀的邊也進(jìn)行處理，且未考慮到平滑后圖斑面積保持的問題。如圖5(b)中弧段的2～3處，此處數(shù)據(jù)不進(jìn)行處理，也能去除鋸齒且更能保持原始數(shù)據(jù)。為了減小平滑后數(shù)據(jù)的偏差，在平滑時此類數(shù)據(jù)不應(yīng)平滑。由此，提出將鋸齒弧段分段為平滑鋸齒單元，平滑鋸齒單元為連續(xù)的相同梯度的鋸齒，僅對平滑鋸齒單元進(jìn)行平滑處理，這樣就能在去除鋸齒的同時盡量減少平滑對數(shù)據(jù)造成的偏差。

圖5 一般去鋸齒算法

相關(guān)概念定義如下，平滑鋸齒單元：相同梯度的鋸齒邊分割為一個平滑鋸齒單元，如圖5(b)中弧段的1～2處為一個平滑鋸齒單元。鋸齒邊：相鄰三邊，中間邊為單位像元長，前后兩邊走向一致，則為鋸齒邊，如圖6(a)所示。鋸齒梯度：相鄰三邊，中間邊為單位像元長，前一條邊與中間邊長度的比值為鋸齒梯度。

2)平滑鋸齒單元分割方法。由上述思想，提出一種平滑鋸齒單元分割方法。設(shè)弧段上點(diǎn)的坐標(biāo)為(Xi,Yi)，相鄰點(diǎn)橫坐標(biāo)差值|Xi=Xi+1|X，縱坐標(biāo)差值|Yi=Yi+1|Y。因?yàn)殇忼X弧段都是由水平或者豎直線段組成，所以相鄰兩點(diǎn)有且僅有一個為0，因此相鄰兩點(diǎn)組成的線段長度L(i,i+1)=|ΔXi|+Yi,，且L(i,i+1)的正負(fù)代表其走向。平滑鋸齒單元分割方法如下：

步驟1：對弧段進(jìn)行冗余點(diǎn)去除。遍歷弧段上所有點(diǎn)，連續(xù)3個點(diǎn)在一條直線上時，中間點(diǎn)為冗余點(diǎn)，需去除，然后進(jìn)入步驟2。

步驟2：依次遍歷弧段上的點(diǎn)，計算相鄰4個點(diǎn)的關(guān)系，當(dāng)條件：L(i-1,i)×L(i+1,i+2)>0,且|L(i,i+1)|的值為單位長時，分以下3種情況，不滿足條件時為線段標(biāo)注0：

②如圖6(b)所示，若|L(i+1,i)

圖6 3種情況

步驟3：標(biāo)注完成后，當(dāng)標(biāo)注連續(xù)數(shù)字相同且大于0的邊分為一個平滑鋸齒單元。如圖7(a)弧段可以分為三段，平滑鋸齒單元1～2，平滑鋸齒單元3～4，非平滑鋸齒單元2～3。

圖7 弧段分段平滑

步驟4：平滑。平滑鋸齒單元進(jìn)行平滑處理。非平滑鋸齒單元保持原始數(shù)據(jù)。如圖7(b)所示弧段分段平滑后。

2 實(shí)驗(yàn)與討論

將分辨率為30 m的Landsat8遙感影像經(jīng)過非監(jiān)督分類為5類，得到分類圖見圖8(a)，將其柵格矢量化后得到圖8(b)，以圖8(b)為實(shí)驗(yàn)圖。

圖8 實(shí)驗(yàn)數(shù)據(jù)

為了表明方法的有效性，實(shí)驗(yàn)對比3種用在圖斑平滑中的算法，分別是概率中值算法[14]、D-P算法[2]( D-P算法使用的壓縮閾值為30)、去尖角算法[15]。3種算法及文中提出的方法均在Visual Studio 2010下用C++語言利用GDAL庫編程實(shí)現(xiàn)，對 shp 文件格式的矢量實(shí)驗(yàn)圖進(jìn)行處理。平滑結(jié)果如圖9所示。

圖9 實(shí)驗(yàn)結(jié)果

為了比較算法間的優(yōu)缺點(diǎn)，實(shí)驗(yàn)對比分析面積保持和位置保持效果以及平滑程度和壓縮率這4個方面的效果，統(tǒng)計結(jié)果見表1與圖10，分析方法如下：

1)平滑程度分析。平滑程度除了可以從目視效果看出，還可以從尖角(0°～90°或者270°～360°的角)數(shù)量得出，尖角數(shù)越少越平滑。

2)面積保持效果分析。對平滑前后的面積作比較，計算各個圖斑的面積變化率，用平滑后的面積減去平滑前的面積得到面積差，面積差除以原始面積得到圖斑的面積變化率。實(shí)驗(yàn)統(tǒng)計面積變化率的最大、最小值和均值，統(tǒng)計結(jié)果見圖10。

3)位置誤差分析。通過緩沖區(qū)限差[17]對位置誤差進(jìn)行分析，緩沖區(qū)限差是對原始弧段建立緩沖區(qū)，統(tǒng)計平滑后弧段完全落入對應(yīng)原始弧段緩沖區(qū)里的弧段數(shù)占總弧段數(shù)的比例，用此比例來評價位置的保持情況，比例值越大，位置保持越好。

緩沖區(qū)限差：

(1)

式中：Ln為落入相應(yīng)原始弧段緩沖區(qū)的平滑后弧段數(shù)，L為總弧段數(shù)。

在統(tǒng)計落入緩沖區(qū)的弧段時，為保證相鄰弧段的緩沖區(qū)不會相交，對原始弧段建立半個像元距離的緩沖區(qū)。

4)壓縮程度。以平滑后的總點(diǎn)數(shù)與原始總點(diǎn)數(shù)之比的壓縮率來評價平滑前后的壓縮程度。

5)時間效率。統(tǒng)計4個方法計算同一數(shù)據(jù)的計算時長，以此來比較4個方法的時間效率。計算時長越短計算效率越高。

表1 4種方法平滑程度、位置誤差、壓縮率、時間效率的評價結(jié)果

圖10 面積變化率

觀察上述實(shí)驗(yàn)結(jié)果，圖9可以清楚觀察到4個算法平滑結(jié)果，僅目視效果而言，D-P算法呈現(xiàn)大量尖角，并不是很美觀。去尖角算法雖然去除了D-P算法的大部分尖角，但是由不同公共邊共同形成的尖角卻未能去除。通過表1中的對尖角的統(tǒng)計，可以看出分段平滑方法平滑程度最高，D-P算法平滑程度最低。由圖10 和表1中面積和位置的統(tǒng)計結(jié)果可以看到，分段平滑方法平滑前后面積保持了一致，而且弧段均在緩沖區(qū)內(nèi)，分段平滑方法面積與位置保持是最優(yōu)的。D-P算法是位置和面積保持最差的，但是在壓縮程度上D-P算法有著絕對的優(yōu)勢。可見D-P算法對圖斑的壓縮平滑，在面積和位置信息上有所犧牲。另外兩個算法的壓縮率和面積位置保持都介于D-P算法和分段平滑方法之間。時間效率上去尖角算法是在D-P算法基礎(chǔ)上進(jìn)行了尖角去除，時間效率上比D-P低，本文算法較另外3種方法稍稍高一點(diǎn)。綜上所述,幾種方法雖然都能有效平滑，但是在平滑程度，面積保持，位置保持以及計算效率方面文中提出的方法最優(yōu)。

3 結(jié)束語

本文針對遙感信息圖斑邊界呈現(xiàn)鋸齒狀，提出一種基于弧段的遙感專題信息圖斑分段平滑方法，用深度搜索的方式提取圖斑弧段，根據(jù)鋸齒弧段的梯度分段為平滑鋸齒單元，對平滑鋸齒單元進(jìn)行平滑處理，以達(dá)到平滑弧段、圖斑面積與位置保持之間平衡，為遙感分類圖，專題圖等矢量化后提供一種平滑方法。該方法減輕了以往方法平滑后對原始數(shù)據(jù)造成的偏差，但平滑后數(shù)據(jù)冗余消除程度遠(yuǎn)不及矢量抽稀算法。