王燕紅，程博，尤淑撐，武盟盟

（1.中國科學(xué)院遙感與數(shù)字地球研究所，北京 100094；2.中國科學(xué)院大學(xué)，北京 100049；3.中國土地勘測規(guī)劃院，北京 100035）

1 引 言

近些年，城市擴張和城市建設(shè)的迅猛發(fā)展，運用遙感技術(shù)進行城市環(huán)境研究具有重要應(yīng)用價值。星載合成孔徑雷達（Synthetic Aperture Radar，SAR）具有全天時、全天候、穿云透霧的優(yōu)勢，對我國幅員遼闊、氣象氣候和下墊面條件復(fù)雜、且急需快速獲取監(jiān)測數(shù)據(jù)的前提下，具有明顯的優(yōu)勢。2007年以來，加拿大 RADARSAT－2、德國 TerraSAR－X及意大利COSMO－SkyMed等分辨率達1m雷達衛(wèi)星的陸續(xù)成功發(fā)射，高分辨率城區(qū)SAR圖像的獲取能力大大提升，一方面為準(zhǔn)確提取建筑區(qū)提供了前提，但是另一方面卻導(dǎo)致高分辨率條件下城區(qū)環(huán)境更為復(fù)雜，給高分辨率SAR影像應(yīng)用帶來了新的挑戰(zhàn)。

建筑區(qū)提取是城區(qū)遙感圖像解譯的重要研究課題之一。建筑區(qū)是城市環(huán)境最基本的一部分，一般由住宅區(qū)、公共建筑、道路與綠地和場地組成。在高分辨率SAR圖像上，灰度分布比較復(fù)雜，用傳統(tǒng)的基于圖像灰度的提取方法很容易受到居民地內(nèi)部要素以及噪聲的影響而得不到滿意的效果，針對高分辨率SAR影像中含有豐富的紋理信息特點，紋理分析方法是對像素間空間分布關(guān)系進行研究，因而可利用紋理分析的方法來解決。國內(nèi)外學(xué)者也做了大量的研究，如Dell’Acqua等［1］利用灰度共生紋理度量來刻畫城區(qū)SAR圖像上不同區(qū)域的建筑覆蓋密度，從而區(qū)分城市中心、居住區(qū)和郊區(qū)。Tison等［2］利用馬爾科夫隨機場紋理分類方法對城區(qū)進行精細的分類。強永剛等［3］利用小波變換和數(shù)學(xué)形態(tài)學(xué)提取人工建筑區(qū)。趙凌君等［4］利用變差函數(shù)的方法，實現(xiàn)了建筑區(qū)的快速提取。變差函數(shù)方法是一種較新的紋理提取方法，在建筑區(qū)和非建筑區(qū)上的特征表現(xiàn)差異大，算法簡單、快速，具有重要應(yīng)用價值，然而目前常用的變差函數(shù)模型易受噪聲干擾，穩(wěn)健性差［5］。本文提出一種基于中值濾波的變差函數(shù)紋理分析方法提取建筑區(qū)，通過改進變差函數(shù)的計算方法，使窗口的大小不再受步長的約束，不僅保留了變差函數(shù)區(qū)分建筑區(qū)與非建筑區(qū)的優(yōu)勢，而且改善了高分辨率SAR提取建筑區(qū)時受強反射點、噪聲干擾大的缺點，提取效果更完整。

2 變差函數(shù)紋理分析

2.1 變差函數(shù)理論

變差函數(shù)理論由數(shù)學(xué)專家G.Maberon教授1962年創(chuàng)立，作為地質(zhì)統(tǒng)計學(xué)的重要工具，變差函數(shù)被應(yīng)用于空間隨機場的統(tǒng)計特性研究。隨著變差函數(shù)理論的不斷發(fā)展，其應(yīng)用領(lǐng)域從早期的空間數(shù)據(jù)分析、地理數(shù)據(jù)位置分析等逐漸擴展到遙感數(shù)據(jù)分析，尤其是在20世紀90年代，Miranda等人成功將變差函數(shù)紋理分析方法應(yīng)用于SAR圖像植被識別之后［6-7］，變差函數(shù)成為了一種重要的紋理特征分析方法，被廣泛應(yīng)用于遙感圖像信息提取領(lǐng)域。

變差函數(shù)又稱半變差函數(shù)（Semivariogram Function），定義為區(qū)域化變量Z（x）和Z（x＋h）（同時包含兩點距離和方向信息）兩點之差的方差之半，即

對于離散的柵格數(shù)據(jù)，實驗中的變差函數(shù)定義為：

其中，N（h）表示觀測數(shù)據(jù)中間距為h的點對數(shù)目，估計值γ＊（h）通常稱為實驗變差函數(shù)。變差函數(shù)用于度量區(qū)域化變量的空間相關(guān)性，能充分反映圖像數(shù)據(jù)的隨機性和結(jié)構(gòu)性。圖1是一個理想的變差函數(shù)模型，當(dāng)h＝0時的變差函數(shù)值C0為塊金值（Nugget Value），代表不確定的變化，隨著h的增大，空間樣本之間的相關(guān)性越來越?。ǚ窍嚓P(guān)性越來越大），當(dāng)樣本變?yōu)橄嗷オ毩r的距離a稱為變程，相應(yīng)的變差函數(shù)值逐漸增大并在間距為a時趨于平穩(wěn)值C0＋C1，這一平穩(wěn)值稱為基臺值（Sill）。

圖1 理想的變差函數(shù)模型

2.2 改進的變差函數(shù)計算方法

變差函數(shù)用于紋理分析的常用計算方式是首先確定步長h、窗口w、計算方向，然后計算窗口w內(nèi)所有間距為h的點對的半方差值，然后取平均（如公式（2））作為窗口中心點的變差函數(shù)值，遍歷全圖即得到影像的變差函數(shù)特征圖。對于高分辨率的SAR圖像而言，窗口w的選取取決于h，要保證窗口w內(nèi)間距h的點對數(shù)目足夠，窗口w至少應(yīng)為3h～5h［4］。然而，窗口w過大會造成圖像整體模糊、邊緣虛警率高，且窗口內(nèi)取平均的方法易受噪聲、孤立強反射點干擾，算法穩(wěn)健性差。

根據(jù)變差函數(shù)計算方法所存在的弊端，本文提出一種穩(wěn)健的變差函數(shù)計算方法，該算法繼承了變差函數(shù)和標(biāo)準(zhǔn)中值濾波方法的優(yōu)點，窗口w取值不再受h約束，主要思想是：在計算像素點的變差函數(shù)值時，不再取窗口w內(nèi)間距為h的點對來求取，而是計算像素點（x，y）與固定方向（如0°方向）上間距h的點（x＋h，y）的半方差值，將此值作為點（x，y）的半方差值；然后取窗口內(nèi)所有像素的半方差值的中值代替均值作為變差函數(shù)值。

計算公式如下：

其中，w為二維模板。

下面以窗口w＝3，h＝1，全方向為例，如圖2所示。

圖2 變差函數(shù)紋理特征計算示意圖

3 基于改進的變差函數(shù)紋理分析的建筑區(qū)提取方法

3.1 高分辨率SAR影像的建筑區(qū)紋理特征

遙感影像上可以分為建筑區(qū)與非建筑區(qū)兩大類。建筑區(qū)可分為城鎮(zhèn)建筑區(qū)和農(nóng)村建筑區(qū)。在高分辨率SAR中，建筑物由于角反射器及屋頂強反射而產(chǎn)生的亮斑，通?；叶龋炼龋┮哂谥苓叺貐^(qū)，然而城鎮(zhèn)建筑區(qū)和農(nóng)村建筑區(qū)由于建筑結(jié)構(gòu)、材質(zhì)、分布的差異，紋理特征也不盡相同。城鎮(zhèn)建筑物分布整齊，樓房之間間隔大，大多數(shù)為平頂整齊的高層樓房，使用材料大多具有良好反射率，當(dāng)屋頂在垂直于雷達觀測角的情況下產(chǎn)生的鏡面發(fā)射是強反射，因此多表現(xiàn)為亮斑。而高層樓房之間道路、草坪等散射回波較弱，在影像上表現(xiàn)為亮斑之間的暗斑塊，加之建筑物通常排列較為整齊，形成有一定周期規(guī)律的明暗相間的紋理特征。農(nóng)村建筑區(qū)房屋間隔小，分布相對密集，通常為房屋聚集的村莊，在影像上表現(xiàn)為棋盤式分布的亮斑，紋理特征不如城鎮(zhèn)建筑區(qū)有規(guī)律。非建筑區(qū)主要包括植被、水體等，整體亮度較低。

3.2 建筑區(qū)的變差函數(shù)曲線

高分辨率SAR圖像中建筑區(qū)、植被區(qū)和水域的實驗變差函數(shù)曲線如圖3所示（A為建筑區(qū)，B為植被區(qū)，C為水域）。從圖中可以看出，首先，建筑區(qū)的變差函數(shù)曲線整體上都遠遠高于植被區(qū)和水域，且隨著間距h的增加周期性地出現(xiàn)峰值和谷點，但是在數(shù)值上，峰值與峰值、谷點與谷點不一定嚴格相等，隨著h的增大會表現(xiàn)出衰減的現(xiàn)象［8-9］，這是由于建筑區(qū)是富含大量強散射點，且建筑區(qū)與街道等間隔分布；植被區(qū)的變差函數(shù)曲線隨著間距h的增大，變差函數(shù)曲線先緩慢上升，隨后趨于平穩(wěn)。水體由于鏡面反射回波很弱，像素灰度值很低且非常均勻，變差函數(shù)曲線近似為一條直線。建筑區(qū)的變差函數(shù)曲線具有特殊性，與非建筑區(qū)特征有明顯差異，因此可以作為有效區(qū)分建筑區(qū)和非建筑區(qū)的依據(jù)。

圖3 實驗區(qū)各地物的變差函數(shù)曲線

3.3 變差函數(shù)參數(shù)選取

（1）步長h

如圖3（b）所示，變差函數(shù)曲線在h1＝6的位置達到第一個峰值，h2＝10的位置達到波谷。在0～h1區(qū)間內(nèi)，曲線上升非?？?，不具有穩(wěn)定性，因此步長h應(yīng)大于h1；另外考慮到隨著步長增大，曲線逐漸衰減，相關(guān)性減小，因此本文選取稍大于h1，小于h2的位置作為步長。這樣既具有穩(wěn)健性，又能保證相關(guān)性比較強。

（2）窗口w

一般來說，窗口選擇太大可能包含多余的無用信息；窗口太小又不能有效、全面地描述地物的特征，而且當(dāng)窗口選擇的較小時，噪聲對紋理的計算影響較大。所以，合適的窗口大小對于紋理的計算十分重要。在實際應(yīng)用中，一般根據(jù)影像地物的實際特征來確定窗口的大小。

（3）方向

要考慮的因素是變差函數(shù)的各向同性和各向異性，即變差函數(shù)計算方向的選擇。一般采用全方向。

3.4 提取流程

建筑區(qū)提取流程如圖4所示，首先進行變差函數(shù)實驗，確定變差函數(shù)參數(shù)h、w，然后計算目標(biāo)影像的變差函數(shù)紋理圖，利用EM算法對得到的紋理圖進行閾值分割，得到初步的提取結(jié)果，利用形態(tài)學(xué)方法進行后處理，包括填洞，去除小面積區(qū)域等，最后得到建筑區(qū)提取結(jié)果。

圖4 算法提取流程

4 實 驗

4.1 實驗數(shù)據(jù)

為驗證算法的有效性，選取了高分辨率圖像中農(nóng)村建筑區(qū)（實驗區(qū)1）和城鎮(zhèn)建筑區(qū)（實驗區(qū)2）進行了實驗。實驗區(qū)1和實驗區(qū)2均采用2007年8月19日獲取北京市升軌右視的TerraSAR－X數(shù)據(jù)，極化方式為HH，空間分辨率為1.25m，大小均為1600×1200像素。

4.2 變差函數(shù)參數(shù)的確定

根據(jù)圖3的變差曲線特征（本文所有實驗影像均來自同一景影像），本文設(shè)定步長h＝8。一般來說，窗口w的選定根據(jù)實驗區(qū)特征選取，目前多數(shù)依靠經(jīng)驗選取。為了定量的分析窗口大小對提取結(jié)果的影響，本文通過固定步長h，改變窗口w的大小，計算實驗區(qū)的檢測率、虛警率、漏警率，以此選擇最優(yōu)的窗口w。

圖5為固定h＝8，檢測率、虛警率、漏警率隨窗口變化的曲線。兩實驗區(qū)的曲線變化大致相同，檢測率曲線先上升后下降，漏警率呈下降趨勢，虛警率則是先下降后上升。當(dāng)窗口在0～20時，檢測率一直上升，大約達到23～28的時候，檢測率最高，此時虛警率和漏警率較小。因此本文選取w＝23作為最優(yōu)窗口。

4.3 實驗結(jié)果及評價

圖5 檢測率、虛警率、漏警率隨窗口變化曲線圖

根據(jù)上節(jié)提出的方法對原始影像進行建筑區(qū)提?。╤＝8，窗口w＝23，全方向），本文將結(jié)合光學(xué)影像人工提取的建筑區(qū)作為標(biāo)準(zhǔn)結(jié)果。實驗區(qū)的初始提取結(jié)果如圖6（b）、圖7（b）所示，然后將初始結(jié)果進行了填洞、去除小區(qū)域（小于300像素），與標(biāo)準(zhǔn)影像進行疊加如圖6（c）、圖7（c）所示，并且與未改進變差函數(shù)的方法的后處理結(jié)果進行對比（白色部分為正確區(qū)域，藍色部分為漏警區(qū)，黃色部分為虛警區(qū)）。

從對比圖可以明顯看出，本文提出的方法更能準(zhǔn)確、完整地提取出建筑區(qū)的邊界，邊緣虛警區(qū)小，建筑區(qū)內(nèi)部的空洞很少，一定程度改善了受孤立強反射點、噪聲干擾影響的缺點，提取效果優(yōu)于未改進的方法。為了使得對比結(jié)果清晰，本文在實驗區(qū)1截取建筑區(qū)邊緣部分進行了放大，如圖8所示，首先觀察建筑區(qū)邊界，本文的方法邊界的虛警很小，而未改進方法明顯地擴大的邊界范圍，其次對于影像中的道路旁邊的強散射點，未改進方法錯誤地提取為建筑區(qū)，造成了虛警，而本文方法一定程度改善了這一缺陷。

為了定量地評價提取結(jié)果，本文采用檢測率（Detection Rate），虛警率（False Alarm Rate）［10］和漏警率（Missing Alarm Rate）3個度量值進行結(jié)果評價，評價結(jié)果如表1所示。從結(jié)果可以看出，改進的方法檢測率較高，虛警率和漏警率較為改進方法均有明顯降低，總體明顯優(yōu)于未改進方法提取的結(jié)果。

圖6 實驗區(qū)1提取結(jié)果

圖7 實驗區(qū)2提取結(jié)果

圖8 提取結(jié)果細節(jié)圖

表1 實驗區(qū)精度評價表

5 結(jié)束語

本文提出了一種改進變差函數(shù)的建筑區(qū)提取方法，該方法以建筑區(qū)及非建筑區(qū)的變差函數(shù)特征為依據(jù)，通過建筑區(qū)變差函數(shù)實驗提取計算間距h，進而確定窗口w，計算此參數(shù)下的變差函數(shù)值作為分類特征，采用EM算法實現(xiàn)建筑區(qū)和非建筑區(qū)的兩分類問題。相比于傳統(tǒng)的變差函數(shù)計算方法，本文在計算方法上提出改進：計算像素點（x，y）與固定方向（如0°方向）上間距h的點（x＋h，y）的半方差值，將此值作為點（x，y）的半方差值，然后計算窗口w內(nèi)的所有像素點的半方差值的中值作為變差函數(shù)值。結(jié)果表明，改進的方法能夠更完整、更高效地提取建筑區(qū)。但是對于紋理特征不明顯，灰度值較低的建筑區(qū)，提取的效果不是很理想，下一步將嘗試引入其他特征對紋理特征不明顯的建筑區(qū)進行提取算法研究。

［1]DELL’ACQUA F，GAMBA P.Texture－based characterization of urban environments on satellite SAR images［J］.IEEE Transactions on Geoscience and Remote Sensing，2003，41（1）：153－159.

［2]TISON C，NICOLAS J M，TUPIN F，et al.A new statistical model for Markovian classification of urban areas in highresolution SAR images［J］.IEEE Transactions on Geosicence and Remote Sensing，2004，42（10）：2046－2057.

［3]強永剛，殷建平，祝恩，等.基于小波變換和數(shù)學(xué)形態(tài)學(xué)的遙感圖像人工建筑區(qū)提?。跩］.中國圖象圖形學(xué)報，2008，13（8）：1459－1464.

［4]趙凌君，高貴，匡綱要.基于變差函數(shù)紋理特征的高分辨率SAR圖像建筑區(qū)提?。跩］.信號處理，2009，（25）9：1433－1442.

［5]金毅，潘懋，姚凌青，等.一種穩(wěn)健變差函數(shù)計算方法［J］.北京大學(xué)學(xué)報（自然科學(xué)版），2009（6）：1033－1038.

［6]MIRANDA F P，F(xiàn)ONSECA L E N，CARR J R，et al.Analysis of JERS－1（Fuyo－1）SAR data for vegetation discrimination in northwestern Brazil using the semivariogram textural classifier［J］.International Journal of Remote Sensing，1996，17（17）：3523－3529.

［7]MIRANDA F P，F(xiàn)ONSECA L E N，CARR J R.Semivariogram textural classification ofJERS－1（Fuyo－1）SAR data obtained over a flooded area of the Amazonrainforest［J］.International Journal of Remote Sensing，1998，19（3）：549－556.

［8]CURRAN P J.The semivariogram in remote sensing：An introduction［J］.Remote Sensing of Environment，1988（24）：493－507.

［9]OLIVER M，WEBSTER R，GERRARD J.Geostatistics in physical geography.Part I：Theory［J］.Transactions of the Institute of British Geographer，1989，41（3）：259－269.

［10]SHUFELT A.Performance evaluation and analysis of monocular building extraction from aerial imagery ［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1999，21（4）：311－326.