付五洲,錢傳俊,陸 彬,李 濤
(長江水利委員會長江口水文水資源勘測局,上海 200136)
合成孔徑雷達干涉測量(Synthetic Aperture Radar Interferometry,InSAR)是基于雷達遙感發(fā)展起來的一種新型主動地表探測技術(shù),能夠全天時、全天候、大范圍、低成本的獲取數(shù)字高程模型、地表變形等。InSAR數(shù)據(jù)處理過程主要包括影像配準(zhǔn)、干涉圖生成、去平地效應(yīng)、相位解纏、絕對相位獲取等,其中關(guān)健的相位解纏步驟是在相位干涉圖基礎(chǔ)上處理的,相位干涉圖質(zhì)量直接影響數(shù)據(jù)處理精度和準(zhǔn)確度[1]。因SAR影像質(zhì)量、時間、空間失相關(guān)影響,干涉圖包含大量噪聲,增加了相位解纏的難度,因此在相位解纏前對干涉圖予以濾波去噪,是InSAR數(shù)據(jù)處理的重要一環(huán)[2]。
常用圖像濾波方法有均值濾波、中值濾波、自適應(yīng)濾波、LEE濾波、Frost濾波等,由于這些方法并不是為干涉圖去噪而設(shè)計,且存在噪聲抑制與分辨率保持之間的矛盾,并不能很好地去除干涉圖中的噪聲[3-7]。本文將小波變換與圓周期中值濾波進行組合研究InSAR干涉圖去噪方法,并與圓周期中值濾波、梯度自適應(yīng)濾波對比,分析小波-圓周期中值組合濾波的去噪效果。
小波變換是一種時間窗和頻率窗都可改變的時頻局域化分析方法,根據(jù)圖像中信號和噪聲在小波域中的不同形態(tài),構(gòu)造相應(yīng)規(guī)則,將高頻系數(shù)中含有噪聲的那一部分濾除,并對小波進行重構(gòu),得到去噪后的干涉圖[6-9]。干涉圖小波變換去噪方法如下:
(1)干涉圖小波分解
利用小波函數(shù)將干涉圖進行信號分解,分解成低頻部分和高頻部分,其中噪聲主要分布在高頻部分。
(2)小波閾值處理
確定閾值大小,對小波分解的高頻部分在水平、垂直和對角3個方向進行閾值處理。
(3)圖像重構(gòu)
由小波分解的低頻部分的系數(shù)和經(jīng)過閾值處理后的高頻部分系數(shù)進行干涉圖重構(gòu),得到小波去噪后的干涉圖。
圓周期中值濾波是一種在干涉相位圖中采用局部統(tǒng)計方法來實現(xiàn)濾波的技術(shù),其先決條件是:局部區(qū)域內(nèi)的地形起伏可以被認(rèn)為是緩變的,具有一定的相關(guān)性,而其中的噪聲則是統(tǒng)計獨立的。若干涉相位圖表示為φ(m,n)(1≤m≤M,1≤n≤N)(m代表干涉圖行方向指標(biāo),n代表干涉圖列方向指標(biāo)),濾波窗口w大小為(2Mw+1)×(2Nw+1),目標(biāo)相位在窗口中心。φ′(m,n)為濾波后的干涉相位,其圓周期中值濾波算法表示為:
φ′(m,n)=
median{arg[exp(jφ(km,kn))/dm,n]}+arg(dm,n)
(1)
式中,median(m,n)是以(m,n)為中心的窗口進行中值運算;(km,kn)為窗口中元素坐標(biāo),jφ(km,kn)為相位點(km,kn)處相位矢量;dm,n為樣點(m,n)主矢量。
圓周期中值濾波較好的保持了相位連續(xù)性和邊緣信息[3,5]。
小波-圓周期中值濾波利用小波變換去除干涉圖中的尖峰或突變,結(jié)合圓周期中值濾波保持相位連續(xù)和邊緣信息,發(fā)揮兩種濾波的優(yōu)點。算法如下:
(1)采用小波函數(shù)sym4干涉圖進行2層分解。
(2)確定閾值大小,對高頻系數(shù)中水平、垂直、對角方向的高頻分量分別進行閾值處理。
(3)將低頻系數(shù)和水平、垂直、對角方向的閾值處理后的高頻系數(shù)進行相位重構(gòu),得到小波變換后的干涉圖。
(4)將經(jīng)小波變換后的干涉圖進行圓周期中值濾波,得到最終的去噪后的干涉圖。
定性評價主要是通過對干涉圖的目視效果來判斷濾波的去噪能力、保持圖像紋理和邊緣信息的能力,定性評價較為直觀。
定量評價是通過數(shù)值指標(biāo)來客觀評價濾波去噪效果[2]。常見評價指標(biāo)有:
(1)均值
均值表示圖像灰度的變化,濾波前后的均值應(yīng)基本保持不變。
(2)均方根誤差
均方根誤差(RMS)是用來衡量濾波后的干涉相位圖和原干涉相位圖之間的偏差,其值越小表示濾波后的干涉圖保真性越好。計算公式如下:
(2)
式中,φs(m,n)是濾波后的相位值,φ0(m,n)是原干涉圖相位值,N是滑動窗口像元個數(shù)。
(3)等效視數(shù)
等效視數(shù)(ENL)是衡量干涉圖斑點噪聲相對強度,是評定濾波器濾波去噪的一項指標(biāo)。ENL越大,濾波器對圖像的平滑效果越好。計算公式如下:
(3)
式中,μ和σ表示干涉相位圖的均值和標(biāo)準(zhǔn)差。
(4)邊緣保持指數(shù)
邊緣保持指數(shù)(EPI)是用來描述因相干抑制處理使圖像中邊緣發(fā)生變化的程度,即濾波器對圖像邊緣信息的保持能力[10]。EPI越接近1表示濾波器的邊緣保持能力越好。計算公式如下:
(4)
(5)殘差點
殘差點的多少是反映濾波器去噪效果的好壞,殘差點數(shù)量越少,濾波器去噪效果越好。計算干涉相位圖中兩兩相鄰的4個像素點中按順時針(或逆時針)對相鄰兩像素點之間求相位差值Δ,然后對Δ求模和,即:
(5)
若S=0,則存在殘差點,否則不存在殘差點。
實驗采用歐空局Sentinel-1A雷達衛(wèi)星于2015年3月25日和2015年4月18號獲取西藏鐵路當(dāng)雄至羊八井段兩景Sentinel-1數(shù)據(jù)。兩景影像時間基線為24 d,空間基線為33.174 m。對兩景雷達影像進行配準(zhǔn)、干涉、去平地效應(yīng)得到干涉相位圖(圖1a),分別采用圓周期中值濾波、梯度自適應(yīng)濾波以及小波-圓周期中值濾波進行去噪。濾波后圖像如圖1中的(b)、(c)、(d)所示。
圖1 原始干涉圖和濾波后的相位干涉圖
對圖1中圓周期中值濾波、梯度自適應(yīng)濾波和小波-圓周期中值濾波的目視判讀,3種濾波均能起到一定的去噪效果,濾波后的干涉圖斑點減少、條紋明顯可見。梯度自適應(yīng)濾波去噪后的干涉相位圖條紋顏色不清晰,依然含有較多的斑點噪聲。圓周期中值濾波和小波-圓周期中值濾波去噪后的干涉條紋清晰,斑點噪聲較少,說明這兩種濾波去噪效果更優(yōu)。對比圖1(b)和圖1(d),圓周期中值濾波去噪后的干涉相位圖存在少量斑點噪聲,條紋不是很連續(xù),而小波-圓周期中值濾波條紋變化較為連續(xù)、平滑,圖像邊緣細(xì)節(jié)信息保持較好。分別對原始干涉相位圖以及濾波后干涉圖第500列的相位數(shù)據(jù)做剖面,繪制剖面散點如圖2所示。
圖2 原始干涉圖和濾波去噪后干涉圖第500列相位數(shù)據(jù)剖面圖
梯度自適應(yīng)濾波去噪后的干涉相位圖第500列的剖面曲線變化平穩(wěn)且幅度較小,但毛刺較多,去噪效果不是很好,邊緣保持能力不強,邊緣、紋理信息流失嚴(yán)重;圓周期中值濾波和小波-圓周期中值濾波去噪后曲線變化幅度較大,但毛刺較少,邊緣保持能力較好。為了更好地比較圓周期中值濾波和小波-圓周期中值濾波的去噪效果,選取剖面圖中A、B、C處去噪后相位值,小波-圓周期中值濾波去除毛刺效果更強,說明小波-圓周期中值濾波在保持邊緣紋理信息的同時可以達到更好的去噪效果。
為進一步評價濾波去噪效果,采用均值(u)、均方根誤差(RMS)、等效視數(shù)(ENL)、邊緣保持指數(shù)(EPI)以及殘差點來定量評價濾波去噪能力(表1)。
表1 干涉圖濾波質(zhì)量評價
從表1可發(fā)現(xiàn),濾波后的相位均值與原始干涉圖相位均值的差異不明顯;圓周期濾波的均方根誤差(RMS)較大,等效視數(shù)(ENL)、邊緣保持指數(shù)(EPI)和去殘差點能力表現(xiàn)適中,說明其保真性相比其他兩種濾波較差、去噪能力一般;梯度自適應(yīng)濾波的均方根誤差最小,但等效視數(shù)、邊緣保持指數(shù)較小,去殘差點能力較弱,這說明梯度自適應(yīng)濾波保證性較好,但去噪和邊緣信息保持能力較差;小波-圓周期中值濾波均方根誤差適中,等效視數(shù)、邊緣保持指數(shù)以及殘差點去除能力表現(xiàn)最優(yōu),說明了小波-圓周期中值濾波的平滑效果最好,邊緣保持能力最強,去噪效果最佳。
本文采用圓周期中值濾波、梯度自適應(yīng)濾波、小波-圓周期中值濾波對歐空局Sentinel-1A衛(wèi)星干涉相位圖進行濾波去噪實驗,采用定性和定量評價。實驗結(jié)果表明:
(1)3種濾波算法均能抑制噪聲干擾,改善干涉圖質(zhì)量。
(2)3種濾波算法中,去噪效果梯度自適應(yīng)濾波較弱,圓周期中值濾波次之,小波-圓周期中值濾波最優(yōu)。
(3)3種濾波算法中,相位保持能力和邊緣細(xì)節(jié)完整性方面,小波-圓周期中值濾波最優(yōu)。
小波-圓周期中值濾波,既有效去除了噪聲,又保持相位條紋的連續(xù)性和邊緣細(xì)節(jié)信息的完整性。