(遼寧工程技術(shù)大學(xué) 測(cè)繪與地理科學(xué)學(xué)院，遼寧 阜新 123000)

0 引言

主動(dòng)遙感技術(shù)中的合成孔徑雷達(dá)(synthetic aperture radar，SAR)具有全天時(shí)、全天候等光學(xué)傳感器所不具有的圖像采集能力，因此作為其信息載體的SAR圖像已在環(huán)境監(jiān)測(cè)，地球資源探測(cè)和災(zāi)害評(píng)估中得到了廣泛的應(yīng)用[1]。圖像分割將圖像細(xì)分為非重疊的均勻區(qū)域，是SAR圖像處理和應(yīng)用的基本步驟[2]。閾值法由于其簡(jiǎn)單直觀的特點(diǎn)，在圖像分割中占有重要地位[3]。在實(shí)際使用中，閾值的選取有多種方法，如直方圖法[4]、最大熵法[5]、Otsu法[6]、模糊聚類(lèi)法[7]等。其中，Otsu法具有原理簡(jiǎn)單、穩(wěn)定有效的優(yōu)點(diǎn)，是使用最廣泛的閾值方法之一。Otsu法由OTSU[8]提出，利用圖像的光譜測(cè)度直方圖，使用窮舉法選取使圖像類(lèi)間方差最大的閾值將圖像劃分為前景和背景2個(gè)部分。然而，由于圖像中的對(duì)象具有各自不同的光譜分布，因此需要較多的閾值才能從其背景中分割出感興趣的對(duì)象。與單閾值Otsu法相比，多閾值Otsu法可以更好地提取圖像目標(biāo)，然而，由于光譜測(cè)度直方圖的零階和一階累積矩的大量重復(fù)計(jì)算，計(jì)算量顯著增加[9]，特別是在需較多閾值的圖像分割中更為突出。

Otsu多閾值圖像分割可視為一個(gè)組合優(yōu)化問(wèn)題，其基本原理為在圖像的光譜測(cè)度范圍內(nèi)，根據(jù)特定的準(zhǔn)則和類(lèi)屬數(shù)搜索系列單閾值的最佳組合。目前，研究人員開(kāi)始嘗試將優(yōu)化方法和Otsu多閾值分割方法相結(jié)合，以減少計(jì)算量，縮短運(yùn)算時(shí)間，提高算法的魯棒性?，F(xiàn)代優(yōu)化算法主要包括蟻群算法[10]、粒子群算法[11]、遺傳算法(genetic algorithm，GA)[12]和禁忌搜索算法[13]等。其中，GA具有過(guò)程簡(jiǎn)單、尋優(yōu)能力強(qiáng)、易與其他算法相結(jié)合等優(yōu)點(diǎn)，是一種常用的優(yōu)化算法。GA是模擬自然界的進(jìn)化機(jī)制，對(duì)個(gè)體由確定的適應(yīng)度函數(shù)逐代演化，以最終演化的個(gè)體作為問(wèn)題的最優(yōu)解[14]。近幾十年來(lái)，研究人員已提出了幾種基于GA的圖像分割方法，如Lee等[15]提出的一種基于GA的圖像分割自動(dòng)閾值選取，以類(lèi)間方差為進(jìn)化方向，利用GA自動(dòng)演化能量函數(shù)參數(shù)，使輸出分割結(jié)果滿(mǎn)足定義的感知能量。趙夫群等[16]提出的一種基于GA-Otsu法的圖像閾值分割，其以圖像類(lèi)間方差為適應(yīng)度函數(shù)；通過(guò)適應(yīng)度值確定進(jìn)一步的搜索方向和范圍，最終確定最優(yōu)閾值實(shí)現(xiàn)圖像分割。但由于GA易陷入局部最優(yōu)，出現(xiàn)早熟現(xiàn)象，且迭代過(guò)程中計(jì)算量大，搜索速度較慢，降低了該類(lèi)方法的實(shí)際應(yīng)用價(jià)值。

量子計(jì)算[17]可以解決傳統(tǒng)計(jì)算領(lǐng)域的許多難題。Narayanan和Moore[18]首次將量子力學(xué)的概念和原理引入到GA中，提出了量子啟發(fā)遺傳算法，并成功應(yīng)用于求解TSP(traveling salesman problem)問(wèn)題。由于這種類(lèi)型的算法不使用量子比特來(lái)編碼染色體，因此它仍然屬于傳統(tǒng)意義上的GA。Han和Kim[19]將量子比特和量子門(mén)的概念引入到GA，并提出了量子遺傳算法(quantum genetic algorithm，QGA)，其采用量子比特編碼染色體增加種群多樣性和全局尋優(yōu)能力；采用量子門(mén)演化量子染色體加快算法收斂，有效地避免早熟現(xiàn)象。目前，QGA已成功用于許多領(lǐng)域，如信號(hào)處理[20]、機(jī)組組合問(wèn)題[21]等。

然而，由于SAR的特殊成像機(jī)制，使SAR圖像具有強(qiáng)相干斑噪聲，降低了圖像質(zhì)量[1]。使得直接對(duì)SAR圖像進(jìn)行Otsu多閾值分割的結(jié)果中含有大量的雜散點(diǎn)[22]，這無(wú)疑會(huì)對(duì)后續(xù)的信息提取造成不可逆的錯(cuò)誤。而提前對(duì)圖像進(jìn)行去噪則是一種消除這種現(xiàn)象的有效操作，但傳統(tǒng)的去噪方法都存在不同程度的局限性，如Lee濾波[23]和曲波濾波法[24]在有限的迭代次數(shù)內(nèi)很難去除原圖中無(wú)用的紋理信息，無(wú)法起到去除大量雜散點(diǎn)的效果，而增加迭代次數(shù)則會(huì)導(dǎo)致對(duì)象的精細(xì)邊緣全部被平滑掉，降低分割精度。SRAD濾波則可在精確地保持邊緣的前提下以很少的迭代次數(shù)平滑掉紋理信息，最大限度地去除雜散點(diǎn)，特別適用于具有復(fù)雜精細(xì)結(jié)構(gòu)的目標(biāo)。但傳統(tǒng)的SRAD由于無(wú)法自動(dòng)控制迭代次數(shù)使得其會(huì)過(guò)度平滑掉部分精細(xì)結(jié)構(gòu)。

綜上，為了解決Otsu多閾值法對(duì)SAR圖像分割具有的噪聲敏感且計(jì)算量大的問(wèn)題，本文采用圖像間平均結(jié)構(gòu)相似性指標(biāo)來(lái)自適應(yīng)地控制SRAD的迭代次數(shù)，以最小迭代次數(shù)達(dá)到既保持精細(xì)邊緣又不過(guò)度平滑掉紋理結(jié)構(gòu)的目的。后經(jīng)自適應(yīng)QGA優(yōu)化閾值組合，獲得最佳的圖像分割效果。

1 算法描述

1.1 改進(jìn)的SRAD濾波

設(shè)I0={I0(x，y)，(x，y)∈Ω}為定義在圖像域Ω、大小為M×N的原始SAR圖像，其中，(x，y)表示像素位置，I0(x，y)∈{0，…，L}為像素(x，y)強(qiáng)度的光譜測(cè)度，L為圖像的光譜測(cè)度量化級(jí)。將I0(x，y)作為圖像初始輸入I(x，y；0)，通過(guò)偏微分方程建立SRAD模型，可表示為：

(1)

(2)

式中：q(x，y；t)為濾波的瞬時(shí)擴(kuò)散系數(shù)，定義為：

(3)

(4)

q0(t)≈q0exp(-ρt)

(5)

式中：ρ為指數(shù)衰減常數(shù)；q0是為初始擴(kuò)散閾值。

對(duì)SRAD模型，使用雅可比迭代法進(jìn)行求解。假設(shè)較小的時(shí)間步長(zhǎng)Δt和在坐標(biāo)x和y方向上較小的空間步長(zhǎng)h，則時(shí)間和空間坐標(biāo)可被離散為：t=nΔt，n=0，1，2，…；x=ih，i=0，1，2，…，M′-1；y=jh，j=0，1，2，…，N′-1；M′h×N′h=M×N。

(6)

(7)

(8)

(9)

(10)

其中，j∈{0，…，N′-1}，i∈{0，…，M′-1}。

(11)

(12)

其滿(mǎn)足對(duì)稱(chēng)邊界條件，

(13)

其中，j∈{0，1，…，N′-1}，i∈{0，1，…，M′-1}。

圖1 擴(kuò)散示意圖

使用前向差分來(lái)近似時(shí)間導(dǎo)數(shù)，則SRAD模型的數(shù)值近似為：

(14)

為了在平滑圖像的同時(shí)保持圖像的結(jié)構(gòu)特征，在迭代過(guò)程中引入平均結(jié)構(gòu)相似性指標(biāo)(mean structural similarity index measurement，MSSIM)[25]作為定量的標(biāo)準(zhǔn)來(lái)控制迭代過(guò)程。表達(dá)式為：

(15)

(16)

|MSSIM(n)-MSSIM(1)|/MSSIM(1)≤ε

(17)

其中，ε為迭代終止閾值，一般取ε=0.5。

當(dāng)滿(mǎn)足式(17)時(shí)，得到平滑圖像z={z(x，y)，(x，y)∈Ω }，其中，z(x，y)為像素(x，y)的光譜測(cè)度值，z(x，y)∈{1，…，L}。

1.2 基于自適應(yīng)QGA的Otsu多閾值圖像分割

利用QGA在函數(shù)尋優(yōu)方面具有計(jì)算量小，收斂速度快等特點(diǎn)，將其在大小為L(zhǎng)的搜索空間中，以圖像類(lèi)間方差最大為準(zhǔn)則，優(yōu)化得到C個(gè)最優(yōu)閾值，并據(jù)此閾值將圖像域Ω分割為C+1個(gè)同質(zhì)區(qū)域，即Ω={Ωk，k=1，…，C+1}，其中，k為同質(zhì)區(qū)域索引，Ωk表示第k個(gè)同質(zhì)區(qū)域，滿(mǎn)足Ωk∩ Ωk′=?，?k，k′∈{1，…， C+1}，k≠ k′， Ω1∪… ∪ ΩC+1=Ω。在尋優(yōu)過(guò)程中，包含以下步驟：量子染色體編碼、量子染色體測(cè)量和解碼、適應(yīng)度評(píng)價(jià)、量子染色體進(jìn)化。

(1)量子染色體編碼。QGA中的染色體用量子比特表示，一個(gè)量子比特的狀態(tài)可以表示為：

|ψ〉=α|0〉+β|1〉

(18)

式中：α和β為概率幅，滿(mǎn)足以下歸一化條件：

|α|2+|β|2=1

(19)

式中：|α|2和|β|2分別表示量子比特坍縮到態(tài)|0〉和|1〉的概率。

(20)

(21)

對(duì)二進(jìn)制串進(jìn)行解碼，得到閾值的十進(jìn)制值，其原理為在第K個(gè)二進(jìn)制串中，第d個(gè)二進(jìn)制值1表示十進(jìn)制值2d-1。則閾值與二進(jìn)制串的關(guān)系可表示為：

(22)

(23)

式中：Tk是第k個(gè)閾值，其滿(mǎn)足T1

(24)

整個(gè)圖像光譜測(cè)度級(jí)的平均值μ為：

(25)

則基于自適應(yīng)QGA的多閾值Otsu的適應(yīng)度函數(shù)F(·)為：

(26)

圖2 量子旋轉(zhuǎn)門(mén)調(diào)整量子比特

(27)

式中：φK，v為第K個(gè)染色體中第v個(gè)量子比特的旋轉(zhuǎn)角，其轉(zhuǎn)換過(guò)程如下：

(28)

因此，第t代種群的更新過(guò)程可描述為：

Q(t+1)=G(t)×Q(t)

(29)

式中：G(t)為第t代的量子門(mén)；Q(t)和Q(t+1)分別為第t代及其進(jìn)化后第t+1代的種群。

表1 旋轉(zhuǎn)角設(shè)置方案

為了使量子染色體更快地向適應(yīng)度高的個(gè)體方向進(jìn)化和避免局部最優(yōu)，當(dāng)量子染色體與歷史最優(yōu)染色體的適應(yīng)度差異較大時(shí)，旋轉(zhuǎn)角將適當(dāng)?shù)卦龃?，反之，將適當(dāng)減小。同時(shí)，根據(jù)相鄰兩代量子染色體的差異，逐代地調(diào)整量子旋轉(zhuǎn)角的大小。具體的調(diào)整策略為：

(30)

1.3 算法流程

算法流程圖如圖3所示。

圖3 流程圖

綜上，算法的基本步驟可總結(jié)為：

①對(duì)原始SAR圖像進(jìn)行改進(jìn)的SRAD濾波；

②定義Otsu多閾值分割時(shí)，圖像的類(lèi)間方差；

③利用量子比特將閾值組合編碼為種群Q(t)；

④對(duì)種群Q(t)進(jìn)行測(cè)量和解碼，得到閾值的十進(jìn)制值；

⑤將②的類(lèi)間方差作為適應(yīng)度函數(shù)，根據(jù)適應(yīng)度大小評(píng)價(jià)各個(gè)個(gè)體的優(yōu)劣，保留適應(yīng)度最高的個(gè)體作為進(jìn)化目標(biāo)；

⑥利用量子旋轉(zhuǎn)門(mén)G(t)作用于Q(t)，得到新種群Q(t+1)；

⑦對(duì)新種群Q(t+1)執(zhí)行④、⑤和⑥；

⑧將算法執(zhí)行到最大進(jìn)化代數(shù)或相鄰兩代的最大類(lèi)間方差的差值小于0.05作為終止條件，如滿(mǎn)足終止條件，則輸出分割結(jié)果，否則繼續(xù)計(jì)算。

2 實(shí)驗(yàn)與分析

為了驗(yàn)證該方法的可行性和有效性，對(duì)模擬和真實(shí)SAR強(qiáng)度圖像進(jìn)行了實(shí)驗(yàn)。

2.1 模擬SAR圖像

圖4(a)示出了尺度為128像素×128像素的模擬圖像模板，其中，數(shù)字Ⅰ-Ⅴ分別表示不同的均勻區(qū)域。表2示出了對(duì)應(yīng)于模擬SAR圖像的均勻區(qū)域的Gamma分布形狀參數(shù)γ和尺度參數(shù)λ，得到的模擬SAR圖像如圖4(b)所示。

圖4 模擬圖像

表2 模擬圖像各同質(zhì)區(qū)域的Gamma分布參數(shù)

利用該方法對(duì)圖4(b)進(jìn)行實(shí)驗(yàn)，分割結(jié)果如圖5所示，其中，圖5(a)和圖5(b)分別為SRAD濾波結(jié)果和分割結(jié)果。為了驗(yàn)證提出方法的優(yōu)越性，利用文獻(xiàn)[16]和文獻(xiàn)[22]的方法對(duì)圖4(b)進(jìn)行實(shí)驗(yàn)，分割結(jié)果如圖5(c)和圖5(d)所示。從圖5(b)～圖5(d)的比較可以看出，文獻(xiàn)[16]的分割結(jié)果含有較多的分割噪聲，而文獻(xiàn)[22]的方法在區(qū)域Ⅱ-V內(nèi)的分割精度有所提高，但在區(qū)域Ⅰ內(nèi)仍較低。另外，3種方法對(duì)圖4(b)的運(yùn)算時(shí)間分別為6.35 s、22.52 s、35.81 s，因此，本文方法能較快得到結(jié)果。

圖5 分割結(jié)果

對(duì)于定性評(píng)估，提取結(jié)果的輪廓線(圖6(a))，將其疊加在圖4(b)的模擬圖像上，如圖6(b)所示。由疊加結(jié)果可以看出，所提輪廓線與實(shí)際輪廓線非常一致，表明本文方法能較好地處理各均勻區(qū)域的邊界。

圖6 視覺(jué)評(píng)價(jià)

對(duì)于定量評(píng)估，根據(jù)模板圖像分別計(jì)算分割結(jié)果圖5(b)～圖5(d)的混淆矩陣[26]，如表3所示。在此基礎(chǔ)上，計(jì)算產(chǎn)品精度、用戶(hù)精度、總體精度和Kappa值(表4)，其中各項(xiàng)指標(biāo)越高，表明其分割精度越高。

從表4可以看出，本文方法得到的產(chǎn)品精度，用戶(hù)精度和總精度均最高，都達(dá)到91.1%以上，Kappa值為0.966，而一般分類(lèi)器Kappa值達(dá)到0.8以上就被視為高質(zhì)量分類(lèi)器[26]。上述定量分析結(jié)果表明了該算法的有效性和準(zhǔn)確性。

表3 模擬圖像分割結(jié)果的混淆矩陣

表4 產(chǎn)品精度、用戶(hù)精度、總精度和Kappa值

2.2 真實(shí)SAR圖像

選取4幅真實(shí)SAR強(qiáng)度圖像如圖7(a)～圖7(d)所示。圖7(a)和圖7(b)是分辨率為25 m，HV極化、尺度為128像素×128像素的Radarsat-Ⅱ圖像，其中，圖7(a)是河口圖像，人為判讀含有2個(gè)同質(zhì)區(qū)域，圖7(b)是城市圖像，人為判讀含有3個(gè)同質(zhì)區(qū)域。圖7(c)和圖7(d)是分辨率為30 m、VV極化、尺度分別為128像素×128像素和256像素×256像素的Radarsat-I海冰圖像，人為判讀含有3個(gè)同質(zhì)區(qū)域。

從圖7可以看出，原始圖像含有大量斑點(diǎn)噪聲，為此，利用改進(jìn)的SRAD濾波原始SAR圖像，其設(shè)置參數(shù)Δt=0.1；ρ=0.3；h=1，由圖8濾波結(jié)果可發(fā)現(xiàn)，圖像斑點(diǎn)噪聲被很好的濾除且邊緣得到了較好的保留。

圖7 真實(shí)SAR圖像

圖9顯示了4幅真實(shí)SAR強(qiáng)度圖像的分割結(jié)果。其中圖9(a1)～圖9(e1)為文獻(xiàn)[16]的分割結(jié)果；圖9(a2)～圖9(e2)為文獻(xiàn)[22]的分割結(jié)果；圖9(a3)～圖9(e3)為本文方法的分割結(jié)果。從圖9可以看出，文獻(xiàn)[16]的分割結(jié)果中含有較多的分割噪聲，文獻(xiàn)[22]的分割結(jié)果的分割噪聲有所降低，但在邊界處存在誤分割現(xiàn)象，而本文方法在分割噪聲和邊界方面均能得到較好地實(shí)現(xiàn)。次外，對(duì)圖7中4幅圖像分割的3種方法的平均運(yùn)行時(shí)間分別為23.22 s、35.42 s及6.03 s，該方法的速度明顯高于其他2種方法。

對(duì)于定性評(píng)估，提取分割結(jié)果的同質(zhì)區(qū)域的邊界線并將其疊加在原始圖像上，如圖10所示。由圖10中可以看出，文獻(xiàn)[16]和文獻(xiàn)[22]中獲得的輪廓線與實(shí)際輪廓線不一致，而本文方法獲得的輪廓線可精確地匹配實(shí)際輪廓。因此，驗(yàn)證了該方法的可行性和有效性。

圖9 分割結(jié)果

圖10 視覺(jué)評(píng)價(jià)

3 結(jié)束語(yǔ)

Otsu法是一個(gè)簡(jiǎn)單而有效的閾值法，具有良好的性能，被廣泛應(yīng)用于圖像分割。然而，它對(duì)噪聲敏感，并且當(dāng)它從單閾值推廣到多閾值時(shí)，因?yàn)樾栌?jì)算每一類(lèi)的累積概率和平均值，存在計(jì)算量大，耗時(shí)過(guò)長(zhǎng)的缺點(diǎn)，這使其在實(shí)際應(yīng)用中受到限制。為了提高分割精度，同時(shí)實(shí)現(xiàn)更好的時(shí)間性能，本文提出一種結(jié)合降斑各向異性擴(kuò)散和自適應(yīng)量子遺傳算法的Otsu多閾值SAR圖像分割方法。該方法克服了傳統(tǒng)方法對(duì)噪聲敏感以及邊緣細(xì)節(jié)丟失等問(wèn)題，通過(guò)在邊緣兩側(cè)以異性擴(kuò)散來(lái)最大限度地去除同質(zhì)區(qū)域的斑點(diǎn)噪聲，并采用平均相似度有效的控制迭代過(guò)程，同時(shí)兼顧了去除斑點(diǎn)噪聲與保持精度。將改進(jìn)的自適應(yīng)量子遺傳算法作為Otsu多閾值的優(yōu)化工具，快速精確地確定最佳閾值組合，使Otsu多閾值圖像分割更加實(shí)用。但本文需人為確定類(lèi)別數(shù)，使其分割自動(dòng)化降低，因此，在未來(lái)工作中，擬對(duì)本文方法進(jìn)行改進(jìn)，快速地實(shí)現(xiàn)SAR圖像的可變類(lèi)分割。