基于對(duì)數(shù)域加性信號(hào)分解的時(shí)序SAR圖像相干斑抑制方法

康 健 童風(fēng)雨 白雨松 丁 翔 冀騰宇 張 柘

①(蘇州大學(xué)電子信息學(xué)院 蘇州 215006)

②(西北工業(yè)大學(xué)數(shù)學(xué)與統(tǒng)計(jì)學(xué)院 西安 710072)

③(蘇州市微波成像處理與應(yīng)用技術(shù)重點(diǎn)實(shí)驗(yàn)室 蘇州 215123)

④(蘇州空天信息研究院 蘇州 215123)

⑤(中國科學(xué)院空天信息創(chuàng)新研究院 北京 100190)

1 引言

合成孔徑雷達(dá)(Synthetic Aperture Radar,SAR)通過主動(dòng)發(fā)射、接收并處理電磁波信號(hào)，成為具有全天時(shí)、全天候?qū)Φ赜^測(cè)能力的遙感技術(shù)，在國防、自然資源等領(lǐng)域發(fā)揮著重要作用。近年來，新的衛(wèi)星及信號(hào)處理技術(shù)的飛速發(fā)展極大地促進(jìn)了SAR對(duì)地觀測(cè)技術(shù)的進(jìn)步，使得同時(shí)滿足高分辨率、大幅寬以及高重訪頻率的SAR圖像獲取成為可能，從而提升了其在高精度、大范圍地表動(dòng)態(tài)監(jiān)測(cè)等方面的能力[1–3]。通過將不同時(shí)間獲取的SAR圖像進(jìn)行配準(zhǔn)，得到的SAR圖像時(shí)間序列能提供被觀測(cè)地區(qū)在時(shí)間維度的變化信息，進(jìn)而能更好地服務(wù)于農(nóng)業(yè)監(jiān)測(cè)、災(zāi)后評(píng)估以及城市規(guī)劃等方面[4]。

SAR成像過程中的相干處理特性使得圖像中像素幅度值呈現(xiàn)出隨機(jī)擾動(dòng)的現(xiàn)象，這主要是雷達(dá)波與像素點(diǎn)內(nèi)各個(gè)散射點(diǎn)相互作用產(chǎn)生的回波相參疊加所導(dǎo)致的。與常見的加性高斯噪聲不同，這種相干斑噪聲是一種乘性噪聲，其非平穩(wěn)特性對(duì)SAR圖像的后續(xù)解譯構(gòu)成了嚴(yán)重挑戰(zhàn)[5,6]。一直以來，SAR圖像的相干斑抑制是SAR圖像處理的重要任務(wù)之一，早期的去噪方法采用簡(jiǎn)單的低通濾波器對(duì)單視SAR圖像進(jìn)行濾波，在噪聲被抑制的同時(shí)，圖像細(xì)節(jié)信息也存在一定程度的丟失。為了解決此問題，文獻(xiàn)[7,8]均利用空間域自適應(yīng)濾波的方法，在降噪的同時(shí)，可以充分保留圖像的細(xì)節(jié)特征。文獻(xiàn)[9]引入了全變分(Total Variation,TV)正則化項(xiàng)，且利用對(duì)數(shù)變換將噪聲的乘性模型轉(zhuǎn)化為滿足加性模型的數(shù)據(jù)擬合項(xiàng)，進(jìn)而提出了相干斑抑制的優(yōu)化模型并利用最優(yōu)化方法得到去噪結(jié)果。近幾年，隨著深度學(xué)習(xí)方法在圖像處理領(lǐng)域的廣泛應(yīng)用，一些工作也將其用來抑制SAR圖像相干斑，并取得了良好的效果[10–12]。其他單通道SAR圖像抑制方法詳見綜述類文獻(xiàn)[6,13,14]。

近年來，世界各國均著眼于SAR衛(wèi)星星座的建設(shè)，如歐空局的哨兵1號(hào)(Sentinel-1)衛(wèi)星星座、中國的高分3號(hào)(Gaofen-3)衛(wèi)星星座，使得具有全球尺度、高重訪頻率的SAR圖像處理成為領(lǐng)域內(nèi)研究的重點(diǎn)?？臻g域配準(zhǔn)得到的具有高時(shí)間分辨率的SAR圖像時(shí)間序列能充分地反映出觀測(cè)時(shí)間段內(nèi)地物變化信息。針對(duì)時(shí)序SAR圖像的相干斑抑制，國內(nèi)外學(xué)者也做了很多相關(guān)工作。文獻(xiàn)[15–17]將單通道SAR圖像搜索同質(zhì)像素點(diǎn)或者圖像塊進(jìn)行去噪的思想拓展到時(shí)序SAR圖像中，即在時(shí)空方向上均進(jìn)行同質(zhì)像素點(diǎn)或圖像塊的搜索，進(jìn)而再進(jìn)行濾波。文獻(xiàn)[18]提出的多基線干涉SAR技術(shù)(SqueeSAR)能自適應(yīng)地將同質(zhì)分布式散射點(diǎn)進(jìn)行選取并聯(lián)合濾波，進(jìn)而同時(shí)實(shí)現(xiàn)對(duì)多時(shí)相SAR幅度及相位進(jìn)行去噪。文獻(xiàn)[19]提出的多時(shí)相SAR圖像三維塊匹配(Multitemporal SAR-BM3D,MSAR-BM3D)方法利用時(shí)間維度上的冗余信息以及變換域協(xié)同濾波將單通道SAR圖像三維塊匹配(SAR Block Matching 3D algorithm,SAR-BM3D)方法[20]進(jìn)行拓展，使其能實(shí)現(xiàn)對(duì)時(shí)序SAR圖像的相干斑抑制效果。文獻(xiàn)[21]提出的基于比值的多時(shí)相SAR圖像去噪(Ratio-based Multitemporal SAR,RABASAR)方法通過利用SAR時(shí)間序列圖像與其在時(shí)間維度平均得到的均值圖像計(jì)算對(duì)數(shù)比，使SAR圖像相干斑的空間非平穩(wěn)性能被很好地抑制，從而能進(jìn)一步提升時(shí)序SAR圖像的去噪效果。

然而，除了空間域SAR圖像中原有的相干斑噪聲帶來的挑戰(zhàn)之外，時(shí)間序列中存在的變化信息，如出現(xiàn)人造目標(biāo)、建筑物的結(jié)構(gòu)變化等，均能使得觀測(cè)區(qū)域在時(shí)間維度上出現(xiàn)異常的散射強(qiáng)度，如圖1所示。由于這些沿時(shí)間維度突變信號(hào)來源于觀測(cè)區(qū)域地物變化，特別是人造目標(biāo)，現(xiàn)有的大部分針對(duì)SAR圖像時(shí)間序列進(jìn)行相干斑抑制的方法并沒有對(duì)此信號(hào)成分進(jìn)行建模，從而導(dǎo)致去噪結(jié)果中含有人為引入的噪點(diǎn)或者出現(xiàn)“虛假散射點(diǎn)”現(xiàn)象。

圖1 含有沿時(shí)間維度突變信號(hào)的SAR幅度時(shí)間序列Fig.1 SAR amplitude time series with outliers along the temporal dimension

另外，在一些應(yīng)用領(lǐng)域，如目標(biāo)檢測(cè)等，這些沿時(shí)間維度突變信號(hào)通常反映了一些人造目標(biāo)在觀測(cè)時(shí)間內(nèi)的變化情況，也能為后續(xù)的解譯工作提供有價(jià)值的信息。為此，Baier等人[22]提出了一種非局部低秩信號(hào)分解的方法對(duì)時(shí)序SAR圖像進(jìn)行相干斑抑制(Robust Nonlocal Low-Rank SAR Time Series Despeckling Considering Speckle Correlation by Total Variation Regularization,DespecKS-NLLRTV)，分別將觀測(cè)數(shù)據(jù)分解為低秩信號(hào)成分、時(shí)間相關(guān)的相干成分以及稀疏的沿時(shí)間維度突變信號(hào)成分，在相干斑抑制的同時(shí)也能降低突變信號(hào)對(duì)去噪結(jié)果產(chǎn)生的影響。雖然此方法在實(shí)際數(shù)據(jù)中的效果明顯，但是由于其加性信號(hào)分解是在原始數(shù)據(jù)域中進(jìn)行的，并沒有將信號(hào)成分與相干斑成分進(jìn)行充分解耦，而且相干斑成分與突變信號(hào)成分均被描述為稀疏項(xiàng)并由 L1范數(shù)進(jìn)行建模，也會(huì)導(dǎo)致相干斑成分與突變信號(hào)成分存在耦合現(xiàn)象，從而導(dǎo)致信號(hào)分解結(jié)果不準(zhǔn)確，并不能真實(shí)地反映觀測(cè)時(shí)間內(nèi)的地物變化所導(dǎo)致的散射強(qiáng)度突變的情況。此外，DespecKS-NLLRTV在對(duì)TV鄰近算子進(jìn)行求解時(shí)采用了Chambolle-Pock算法[23]，其運(yùn)算效率不高，進(jìn)一步限制了該方法在實(shí)際數(shù)據(jù)中的應(yīng)用。

為了對(duì)時(shí)序SAR圖像進(jìn)行相干斑抑制的同時(shí)，準(zhǔn)確地提取出觀測(cè)時(shí)間段內(nèi)產(chǎn)生的沿時(shí)間維度突變散射點(diǎn)，本文提出了一種基于對(duì)數(shù)域加性信號(hào)分解的方法，首先，通過對(duì)對(duì)數(shù)域SAR強(qiáng)度值進(jìn)行統(tǒng)計(jì)建模，得到時(shí)序SAR數(shù)據(jù)的對(duì)數(shù)似然，再引入針對(duì)不同成分的先驗(yàn)知識(shí)，提出相應(yīng)的正則項(xiàng)，進(jìn)而得到各個(gè)成分的后驗(yàn)概率，通過最大后驗(yàn)概率(Maximum A Posterior,MAP)的估計(jì)方法，再得出問題對(duì)應(yīng)的優(yōu)化模型，最后利用交替方向乘子法(Alternating Direction Method of Multipliers,ADMM)[24]進(jìn)行求解優(yōu)化。所提方法在仿真以及實(shí)際數(shù)據(jù)中均驗(yàn)證了效果，通過對(duì)比試驗(yàn)分析了所提方法在去噪性能以及突變信號(hào)成分分離效果上的提升。

2 問題建模與求解方法

本文提出的方法總體框架如圖2所示，為了保持時(shí)序SAR圖像中永久散射點(diǎn)(Persistent Scatterer,PS)的散射特性，首先計(jì)算振幅離差指數(shù)(Amplitude Dispersion Index,ADI)[25]剔除可能的永久散射點(diǎn)，對(duì)于每個(gè)非永久散射點(diǎn)，再通過Kolmogorov-Smirnov (K-S)檢測(cè)[18]將非局部搜索窗口中的同質(zhì)像素點(diǎn)進(jìn)行篩選，進(jìn)而將所有同質(zhì)點(diǎn)與其排列并構(gòu)成待分解矩陣，應(yīng)用所提出的對(duì)數(shù)域信號(hào)分解方法得到相干斑抑制成分以及沿時(shí)間維度突變信號(hào)成分，最后采用聚合求均值策略得到最后的去噪以及突變信號(hào)結(jié)果。

圖2 本文所提方法流程圖Fig.2 Flow chart of the proposed method

2.1 問題建模

SAR圖像的強(qiáng)度值服從以下Gamma分布[9,26]：

其中，I ∈R+表示像素強(qiáng)度值，R ∈R+為待恢復(fù)的信號(hào)強(qiáng)度值，L>0為視數(shù)個(gè)數(shù)，Γ(·)為Gamma函數(shù)。在乘性噪聲模型下，強(qiáng)度I表示為

其中，n為均值，為E(n)=1的獨(dú)立同分布的隨機(jī)變量，進(jìn)而得到強(qiáng)度的均值與方差為

從式(4)可以看出，噪聲在空間方向上是非平穩(wěn)的，并與信號(hào)真實(shí)強(qiáng)度值有關(guān)。為了方便處理，通常采用對(duì)數(shù)函數(shù)將乘性噪聲模型進(jìn)行轉(zhuǎn)化，得到滿足加性噪聲模型的Fisher-Tippett分布[26]：

其中，g=lg(I),z=lg(R)。根據(jù)上述對(duì)數(shù)域中的強(qiáng)度似然模型，在MAP框架下，待恢復(fù)的信號(hào)可以通過求解以下模型來得到：

其中，lgpg(g|z)為對(duì)數(shù)域中的信號(hào)對(duì)數(shù)似然，R(z)為待恢復(fù)信號(hào)的先驗(yàn)知識(shí)對(duì)應(yīng)的正則化項(xiàng)，λ為正則化項(xiàng)參數(shù)。對(duì)于時(shí)序SAR圖像，空間像素點(diǎn)中的時(shí)間序列對(duì)數(shù)域強(qiáng)度值服從以下對(duì)數(shù)似然：

其中，t為時(shí)間索引，C 為常量，在對(duì)變量zt進(jìn)行優(yōu)化時(shí)，可以將其忽略，并不影響優(yōu)化結(jié)果。

上述MAP問題需要對(duì)對(duì)數(shù)域時(shí)間序列強(qiáng)度值的先驗(yàn)知識(shí)進(jìn)行建模，本文主要考慮以下先驗(yàn)：

(1) 非局部自相似性：在一定觀測(cè)范圍內(nèi)，存在地物目標(biāo)尺寸遠(yuǎn)大于SAR圖像分辨率或者不同目標(biāo)屬于同種地物類別等情況，使得區(qū)域內(nèi)的目標(biāo)散射特性有相似性，進(jìn)而可以利用非局部搜索方法將一些同質(zhì)像素點(diǎn)對(duì)應(yīng)的時(shí)間序列找出并進(jìn)行聯(lián)合處理。為此，本文利用K-S檢測(cè)方法將一定搜索范圍內(nèi)(例如，21×21窗口范圍)的同質(zhì)點(diǎn)進(jìn)行搜索選取，再將候選時(shí)間序列構(gòu)成觀測(cè)矩陣。

(2) 平滑性：由于觀測(cè)矩陣是由同質(zhì)像素點(diǎn)對(duì)應(yīng)的時(shí)間序列構(gòu)成的，其在時(shí)間以及空間維度上(表示為矩陣的行和列)具有平滑特性，即散射強(qiáng)度變化不大，本文采用TV范數(shù)來對(duì)此進(jìn)行建模。

(3) 稀疏性：對(duì)于時(shí)序信號(hào)中存在的沿時(shí)間維度突變信號(hào)，其本身具有稀疏特性，本文利用 L1范數(shù)對(duì)其進(jìn)行刻畫。

根據(jù)上述似然以及正則化項(xiàng)，本文提出了如下模型對(duì)時(shí)序SAR圖像進(jìn)行相干斑抑制并對(duì)觀測(cè)數(shù)據(jù)中的穩(wěn)定信號(hào)及突變信號(hào)成分進(jìn)行分離：

其中，元素gti所構(gòu)成的矩陣是將K-S檢測(cè)方法選取到的同質(zhì)像素點(diǎn)對(duì)應(yīng)的時(shí)序列向量進(jìn)行按行排列得到的，其時(shí)間維度及空間維度的索引分別為t,i。式(8)中X為相干斑抑制后的信號(hào)成分，S為沿時(shí)間維度突變信號(hào)成分，α,β分別為不同正則化項(xiàng)對(duì)應(yīng)的參數(shù)。給出矩陣A，其TV范數(shù)∥A∥TV定義為

其中，i1,i2分別為矩陣列和行的索引。

2.2 求解方法

為了有效求解所提模型式(8)，本文首先利用變量替換及分解方法[9]將上述模型進(jìn)行轉(zhuǎn)化，得到：

其中，D(·):=[Dt(·);Di(·)]表示二維1階導(dǎo)數(shù)算子。再將有約束的優(yōu)化模型轉(zhuǎn)化為無約束優(yōu)化模型，得到增廣拉格朗日函數(shù)：

其中，〈·〉表示為內(nèi)積算子，∥·∥F為矩陣F-范數(shù)，T1,T2為優(yōu)化引入的輔助變量，μ為懲罰項(xiàng)對(duì)應(yīng)的參數(shù)。進(jìn)而利用ADMM方法對(duì)其進(jìn)行優(yōu)化求解，求解過程如下。

2.2.1Z子問題

類似于文獻(xiàn)[26]，上述問題可以通過牛頓法(Newton’s method)進(jìn)行逐像素多次迭代求解：

2.2.2X子問題

可以通過計(jì)算式(14)對(duì)于變量X的梯度，并且令梯度為0，得到線性方程：

其中，D?(·)表示D(·)的伴隨算子，利用頻域快速求解方法[27]得到式(14)的解。

2.2.3Y,S 子問題

式(16)與式(17)的解可以由軟閾值(Soft Thresholding)[28]算子得出。

2.2.4 乘子T1,2更新

按照式(18)，式(19)子問題求解規(guī)則，在迭代一定次數(shù)并滿足收斂條件后，得到的與為當(dāng)前像素點(diǎn)及其同質(zhì)點(diǎn)時(shí)間序列的相干斑抑制矩陣和突變信號(hào)矩陣。

由于上述分解過程是在對(duì)數(shù)域中進(jìn)行的，需要將這兩種成分轉(zhuǎn)換到原始信號(hào)域，令為去噪后的SAR圖像像素強(qiáng)度值，根據(jù)所提方法的分解模型可得到：

則在原始信號(hào)域中，exp()為去噪后的穩(wěn)定信號(hào)成分，exp()(exp()-1)為沿時(shí)間維度突變信號(hào)成分。

2.3 永久散射點(diǎn)剔除

2.2 節(jié)所提方法是針對(duì)每個(gè)非永久散射點(diǎn)進(jìn)行的，由于永久散射點(diǎn)幾乎不受相干斑噪聲的影響，這些散射點(diǎn)的時(shí)序散射信息需要被保留。本文采用ADI方法對(duì)其進(jìn)行剔除，通過計(jì)算SAR時(shí)間序列幅度值的均值與標(biāo)準(zhǔn)差之間的比值，即mA/σA，并將大于一定閾值的像素點(diǎn)進(jìn)行剔除。

對(duì)于每個(gè)非永久散射點(diǎn)，所提方法均會(huì)找到相對(duì)應(yīng)的同質(zhì)點(diǎn)時(shí)序信號(hào)，并對(duì)所構(gòu)成的矩陣進(jìn)行分解，在遍歷所有像素點(diǎn)后，再利用求平均的策略將每個(gè)像素點(diǎn)對(duì)應(yīng)的去噪成分以及沿時(shí)間維度突變信號(hào)成分進(jìn)行聚合，從而得到最終的相干斑抑制結(jié)果以及突變信號(hào)分離結(jié)果。

3 仿真與實(shí)測(cè)數(shù)據(jù)實(shí)驗(yàn)

3.1 仿真數(shù)據(jù)實(shí)驗(yàn)結(jié)果

本文采用的仿真數(shù)據(jù)與文獻(xiàn)[22]一致，其幅度圖像與其隨時(shí)間變化如圖3所示。產(chǎn)生均值為1的高斯隨機(jī)變量以及值為0.3的指數(shù)相參損失用來模擬SAR圖像中的相干斑效應(yīng)。共生成24幅圖像作為仿真的時(shí)間序列。此外，分別隨機(jī)選取了0.5%以及1.0%比例的SAR像素，用服從Rice分布的隨機(jī)數(shù)將其替換，用來模擬沿時(shí)間維度突變信號(hào)干擾現(xiàn)象。利用K-S檢測(cè)選擇同質(zhì)點(diǎn)的閾值設(shè)置為0.25，正則化項(xiàng)參數(shù)α,β分別選取為5.0,10.0。為了驗(yàn)證所提方法的有效性，本文選擇了以下4種時(shí)序SAR圖像去噪方法進(jìn)行對(duì)比：(1) SqueeSAR[18];(2) MSAR-BM3D[19];(3) RABASAR[21];(4) DespecKS-NLLRTV[22]，并采用峰值信噪比(Peak Signal-to-Noise Ratio,PSNR)、平均結(jié)構(gòu)相似性指標(biāo)(Mean Structure Similarity Index Measure,MSSIM)以及平均互相關(guān)(Mean Cross Correlation,MCC)來定量驗(yàn)證方法在相干斑抑制上的效果。

圖3 本文采用的仿真數(shù)據(jù)Fig.3 Simulation data used in this paper

圖4(a)—圖4(f)、圖4(g)—圖4(l)分別展示了所有對(duì)比方法在沒有沿時(shí)間維度突變信號(hào)以及在0.5%像素的突變信號(hào)干擾下得到的相干斑抑制結(jié)果。相應(yīng)的定量結(jié)果如表1所示。在沒有突變信號(hào)的情況下，現(xiàn)有的方法如MSAR-BM3D和RABASAR能取得比DespecKS-NLLRTV和本文方法更好的相干斑抑制效果。其中一個(gè)原因是MSAR-BM3D和RABASAR均利用了相似圖像塊對(duì)噪聲進(jìn)行抑制，這種方法在例子中給出的同質(zhì)區(qū)域較多的場(chǎng)景中能更好地在空間維度對(duì)信號(hào)進(jìn)行平滑處理，然而，在MSSIM指標(biāo)下，所提出的方法能取得比MSARBM3D更好的效果，主要原因是基于同質(zhì)像素點(diǎn)的濾波方法能更好地保持場(chǎng)景的空間結(jié)構(gòu)信息。當(dāng)突變信號(hào)的數(shù)量開始增加時(shí)，DespecKS-NLLRTV和本文方法均能魯棒地對(duì)相干斑噪聲進(jìn)行抑制，其恢復(fù)效果很大程度上并不受突變信號(hào)的數(shù)量影響。由于其他方法并沒有考慮時(shí)序數(shù)據(jù)中存在的沿時(shí)間維度突變信號(hào)，其圖像恢復(fù)效果隨著突變信號(hào)的增多而變差。與本文所提出的對(duì)數(shù)域信號(hào)分解方法相比，直接在原始信號(hào)域進(jìn)行信號(hào)分解的方法DespecKS-NLLRTV在PSNR等指標(biāo)上有大約3 dB的差距。由于原始信號(hào)域中的信號(hào)與相干斑噪聲成分存在耦合，直接對(duì)觀測(cè)得到的數(shù)據(jù)進(jìn)行信號(hào)分解并不能消除耦合產(chǎn)生的影響，使得信號(hào)的恢復(fù)結(jié)果并不理想。此外，對(duì)于非局部自相似性構(gòu)建的時(shí)序矩陣，DespecKS-NLLRTV引入了低秩先驗(yàn)進(jìn)行約束，這會(huì)導(dǎo)致恢復(fù)得到的時(shí)序信號(hào)能量有所丟失，也會(huì)使得PSNR評(píng)價(jià)指標(biāo)降低。與此相比，本文所提出的方法充分考慮到了相干斑噪聲的統(tǒng)計(jì)特性，并在對(duì)數(shù)域中將信號(hào)成分與其他成分進(jìn)行分解，能在噪聲抑制的同時(shí)，很好地保持恢復(fù)信號(hào)的原始結(jié)構(gòu)。由于本文方法并沒有采用低秩約束，并不會(huì)導(dǎo)致恢復(fù)出的時(shí)序信號(hào)能量丟失，并且可以大幅提升算法的計(jì)算效率。值得注意的是，在沿時(shí)間維度突變信號(hào)數(shù)量較少或不存在的情況下，所提出算法的PSNR值反而較低，而隨著突變信號(hào)數(shù)量增加，PSNR有所增加。其原因在于，沿時(shí)間維度突變信號(hào)的數(shù)量增加使得觀測(cè)信號(hào)的本質(zhì)特性更加符合所提方法的信號(hào)分解模型，從而使得分解得到的穩(wěn)定信號(hào)成分更加接近于真值。如果不存在突變信號(hào)，仍然會(huì)有一部分信號(hào)成分被分解成為稀疏信號(hào)項(xiàng)，反而會(huì)影響相干斑抑制得到的信號(hào)成分恢復(fù)精度。

表1 仿真數(shù)據(jù)相干斑抑制結(jié)果的定量分析Tab.1 Quantitative analysis of speckle suppression results in simulation data

圖4 仿真數(shù)據(jù)相干斑抑制結(jié)果Fig.4 Speckle suppression results of the simulated data

此外，為了測(cè)試提出方法對(duì)于不同參數(shù)設(shè)置下的敏感程度，以PSNR值為例，圖5給出了在不同α以及β值下的時(shí)序圖像恢復(fù)效果?？梢钥闯?，當(dāng)α>1,β >2時(shí)，所提方法的恢復(fù)效果受參數(shù)值變化敏感度不高，而且當(dāng)突變信號(hào)增多時(shí)，給予 L1范數(shù)項(xiàng)更大的權(quán)重可以使得恢復(fù)效果有進(jìn)一步提升。為了進(jìn)一步驗(yàn)證所提方法受圖像數(shù)目的影響，分別從原始的24張圖像中選取3,5,10,15,20張圖像構(gòu)成新的時(shí)序數(shù)據(jù)，用所提方法對(duì)其進(jìn)行相干斑抑制，相應(yīng)的超參數(shù)設(shè)置保持不變，即α=5.0,β=10.0，進(jìn)而計(jì)算得到結(jié)果的PSNR值如圖6所示?？梢钥闯?，隨著圖像數(shù)目的增加，所提方法的相干斑抑制結(jié)果也隨之提升。如果圖像數(shù)目過少，會(huì)影響K-S檢測(cè)方法對(duì)于同質(zhì)點(diǎn)選取的精度，進(jìn)而使得所構(gòu)成的時(shí)序矩陣自相似性較差，導(dǎo)致所提方法的去噪效果降低。因此，所提出的方法適用于圖像數(shù)目較多(多于20張)的時(shí)序SAR數(shù)據(jù)去噪。此外，盡管實(shí)驗(yàn)中突變信號(hào)的數(shù)目增加，所提方法并不受到突變信號(hào)數(shù)目不同而產(chǎn)生效果上的影響，從而進(jìn)一步驗(yàn)證了所提方法的魯棒性。

圖5 本文方法的參數(shù)敏感性分析Fig.5 Parameter sensitivity analyses of the proposed method

3.2 仿實(shí)測(cè)數(shù)據(jù)實(shí)驗(yàn)結(jié)果

本文分別選取了上海近海以及浦東機(jī)場(chǎng)區(qū)域的垂直極化干涉寬幅哨兵1號(hào)(Sentinel-1 IW VV)數(shù)據(jù)，兩組時(shí)間序列獲取時(shí)間段均為2019-12-02到2021-02-06，每組包含36張圖像。其距離向及方位向分辨率分別約為2.33 m和13.94 m。所提方法的參數(shù)設(shè)置與3.1節(jié)仿真數(shù)據(jù)一致。

圖7展示了所有對(duì)比方法分別在2020-09-15以及2021-02-06獲取的SAR圖像上的相干斑抑制結(jié)果?？梢钥闯觯^(qū)域的艦船產(chǎn)生的強(qiáng)散射點(diǎn)構(gòu)成了時(shí)間序列中的沿時(shí)間維度突變信號(hào)。傳統(tǒng)的方法如SqueeSAR對(duì)于艦船區(qū)域的去噪效果差，其平均結(jié)果導(dǎo)致了此區(qū)域產(chǎn)生模糊的亮斑，對(duì)于MSARBM3D，可以發(fā)現(xiàn)艦船附近區(qū)域的去噪結(jié)果存在人為引入的噪聲，使得其并不能真實(shí)反映出海面散射特性。盡管RABASAR能對(duì)大部分同質(zhì)區(qū)域進(jìn)行很好的相干斑抑制，但對(duì)于艦船目標(biāo)，其散射結(jié)果易出現(xiàn)模糊現(xiàn)象，從而損失了真實(shí)的船體結(jié)構(gòu)信息。相比于上述幾種方法，采用信號(hào)分解思想的DespecKS-NLLRTV和本文方法均能通過沿時(shí)間維度突變信號(hào)分離的策略降低其對(duì)于相干斑抑制的影響。相比于DespecKS-NLLRTV，本文方法更好地利用了SAR圖像的統(tǒng)計(jì)特性，減小了直接在原始圖像域中做分解存在的信號(hào)與相干斑噪聲的耦合效應(yīng)。另外，由于DespecKS-NLLRTV中的模型對(duì)于突變信號(hào)以及相干斑成分均采用L1范數(shù)進(jìn)行建模，會(huì)導(dǎo)致兩種成分分離結(jié)果互相干擾。相比之下，本文方法能更加準(zhǔn)確地對(duì)突變信號(hào)進(jìn)行提取。為定量評(píng)估不同方法的相干斑抑制結(jié)果，我們分別在兩次觀測(cè)場(chǎng)景中選取了同質(zhì)區(qū)域A1和A2，并根據(jù)所有方法得到的去噪結(jié)果計(jì)算相應(yīng)區(qū)域內(nèi)的等效視數(shù)(Equivalent Number of Looks,ENL)。從表2可以看出，所提方法可以實(shí)現(xiàn)良好的去噪效果，可以對(duì)于同質(zhì)區(qū)域內(nèi)的噪聲進(jìn)行有效抑制。雖然MSAR-BM3D在這些同質(zhì)區(qū)域內(nèi)取得了最好的等效視數(shù)，但其在散射特性復(fù)雜的區(qū)域并不能得到很好的去噪結(jié)果，尤其在具有沿時(shí)間維度突變散射點(diǎn)附近，其去噪結(jié)果通常會(huì)引入人為噪聲干擾。與上述分析一致，相比于原始信號(hào)域信號(hào)分解方法DespecKS-NLLRTV，本文提出的信號(hào)分解策略能在同質(zhì)區(qū)域內(nèi)獲得更高的ENL。此外，針對(duì)得到的突變信號(hào)圖像，本文進(jìn)一步計(jì)算其香農(nóng)熵，從而定量地反映了圖像的復(fù)雜程度。表3列出了兩種方法得到的突變信號(hào)圖像所計(jì)算的香農(nóng)熵，可以看出本文方法提取得到的突變信號(hào)圖像質(zhì)量較高，而DespecKS-NLLRTV并不能很好地將相干斑噪聲與突變信號(hào)成分進(jìn)行有效的分離，使得提取得到的突變信號(hào)中含有噪聲干擾，從而導(dǎo)致了其圖像熵值較高。

表2 4塊同質(zhì)區(qū)域計(jì)算得到的等效視數(shù)Tab.2 The equivalent apparent number calculated by four homogeneous regions

表3 沿時(shí)間維度突變信號(hào)的熵值分析Tab.3 Entropy analysis of the outliers along the temporal dimension

圖7 上海近海區(qū)域哨兵1號(hào)數(shù)據(jù)的相干斑抑制結(jié)果Fig.7 Speckle suppression results for Sentinel-1 data in the Shanghai offshore region

圖8給出了所有方法在上海浦東機(jī)場(chǎng)區(qū)域2019-12-02以及2020-09-27獲取到的SAR圖像相干斑抑制結(jié)果。場(chǎng)景中由于飛機(jī)目標(biāo)的出現(xiàn)(放大區(qū)域)，其強(qiáng)后向散射特性構(gòu)成了機(jī)場(chǎng)區(qū)域內(nèi)的SAR像素強(qiáng)度突變，使其成為時(shí)間序列中的沿時(shí)間維度突變信號(hào)。在這些區(qū)域內(nèi)，SqueeSAR方法中的簡(jiǎn)單空間平均模糊掉了強(qiáng)散射目標(biāo)，從而使得去噪結(jié)果中含有虛假亮斑。盡管MSAR-BM3D和RABASAR能對(duì)大部分區(qū)域進(jìn)行有效的相干斑抑制，但由于這些方法沒有特別考慮時(shí)序圖像中的沿時(shí)間維度突變信號(hào)，使得去噪結(jié)果中包含人為引入的噪點(diǎn)或者導(dǎo)致突變信號(hào)的目標(biāo)結(jié)構(gòu)被過度平滑。針對(duì)選取的同質(zhì)區(qū)域A3和A4，所提方法均能在取得良好去噪效果的同時(shí)，實(shí)現(xiàn)對(duì)于突變目標(biāo)的提取，并在所提取的圖像中完整地保留目標(biāo)的結(jié)構(gòu)。表3得到的突變信號(hào)圖像熵值計(jì)算結(jié)果進(jìn)一步說明了兩種分解方法對(duì)于異常目標(biāo)分離效果差異。

圖8 上海浦東機(jī)場(chǎng)區(qū)域哨兵1號(hào)數(shù)據(jù)的相干斑抑制結(jié)果Fig.8 Speckle suppression results for Sentinel-1 data in Shanghai pudong airport region

為了驗(yàn)證不同方法對(duì)于圖像細(xì)節(jié)的保持特性，分別選取了近海區(qū)域(如圖9所示)以及海上艦船區(qū)域(如圖10所示)，并利用邊緣檢測(cè)技術(shù)對(duì)所有算法的去噪結(jié)果進(jìn)行圖像結(jié)構(gòu)信息提取?？梢钥闯?，本文提出的方法能很好地保持近海區(qū)域的地物的結(jié)構(gòu)信息，并沒有因相干斑噪聲抑制導(dǎo)致圖像細(xì)節(jié)信息丟失。與之對(duì)比，RABASAR等方法并沒有很好地保持圖像的結(jié)構(gòu)信息，從而使得邊緣提取結(jié)果并不準(zhǔn)確。在艦船區(qū)域結(jié)果中，可以看出所提方法很好地將艦船與背景噪聲進(jìn)行分離，使得艦船邊緣提取結(jié)果并沒有受到背景噪聲的干擾。與之相比，RABASAR等方法邊緣提取結(jié)果模糊，不能準(zhǔn)確地反映目標(biāo)的結(jié)構(gòu)特性。

圖9 近海區(qū)域圖像邊緣提取結(jié)果對(duì)比Fig.9 Comparison of image edge extraction results in offshore areas

圖10 海上艦船區(qū)域圖像邊緣提取結(jié)果對(duì)比Fig.10 Comparison of image edge extraction results of ship region

4 結(jié)語

本文針對(duì)時(shí)序SAR圖像中的相干斑噪聲，提出了一種基于對(duì)數(shù)域加性信號(hào)分解技術(shù)的相干斑抑制方法，所提方法不僅能很好地對(duì)相干斑進(jìn)行抑制，并且能對(duì)時(shí)間序列中的穩(wěn)定信號(hào)成分以及沿時(shí)間維度突變信號(hào)成分進(jìn)行有效分離，得到的突變信號(hào)能展示出觀測(cè)時(shí)間內(nèi)由于人造目標(biāo)出現(xiàn)或地物變化產(chǎn)生的強(qiáng)散射現(xiàn)象，為后續(xù)的圖像解譯提供了有效的參考數(shù)據(jù)。為了驗(yàn)證方法的有效性，本文選取了現(xiàn)有的4種主流時(shí)序SAR圖像相干斑抑制方法進(jìn)行對(duì)比，并且用PSNR等幾種指標(biāo)定量地分析了不同方法在仿真數(shù)據(jù)上的去噪效果。此外，本文選取了上海近海以及浦東機(jī)場(chǎng)區(qū)域的哨兵1號(hào)時(shí)序SAR數(shù)據(jù)，并將上述方法所得結(jié)果進(jìn)一步進(jìn)行對(duì)比分析，均能驗(yàn)證所提方法的有效性。在后續(xù)工作中，會(huì)進(jìn)一步驗(yàn)證本文所提方法在觀測(cè)角度大幅變化情況下獲取的SAR時(shí)序數(shù)據(jù)上的相干斑抑制效果。