楊利紅,趙變紅,張星祥,任建岳*
(1.中國科學(xué)院 長春光學(xué)精密機(jī)械與物理研究所,吉林 長春130033;2.中國科學(xué)院 研究生院,北京100039)
圖像恢復(fù)的目的是從被退化函數(shù)模糊和被噪聲污染的圖像中盡可能地恢復(fù)出真實(shí)的場景。經(jīng)典的圖像恢復(fù)方法有維納濾波、約束最小二乘法濾波、Lucy-Richardson 算法等,這些方法均需要知道退化函數(shù)的類型和參數(shù)才能取得一定的恢復(fù)效果[1-3]。由于航天遙感相機(jī)運(yùn)行的空間環(huán)境特殊,圖像獲取的大氣條件多變,因此,往往無法獲得準(zhǔn)確的退化函數(shù),只能利用成像系統(tǒng)的參數(shù)或獲取的圖像對系統(tǒng)退化函數(shù)進(jìn)行辨識,然后,根據(jù)估計(jì)出的退化函數(shù),選擇合適的算法對圖像進(jìn)行恢復(fù)。
隨著航天科技的發(fā)展,遙感圖像恢復(fù)取得了很多成果,退化函數(shù)估計(jì)是其中的一個(gè)重要研究方向。航天遙感相機(jī)成像過程中的退化函數(shù)為點(diǎn)擴(kuò)散函數(shù)( PSF) ,PSF 傅里葉變換的歸一化模值為調(diào)制傳遞函數(shù)( MTF) 。顧行發(fā)等人提出了模擬理想靶標(biāo)檢測衛(wèi)星在軌MTF,根據(jù)45°方向的MTF 數(shù)據(jù)建立插值二維MTF 矩陣來對圖像進(jìn)行恢復(fù)[4],此方法需要人工鋪設(shè)模擬點(diǎn)光源的靶標(biāo),適用性不強(qiáng);陳強(qiáng)等人利用實(shí)驗(yàn)室測定的相機(jī)系統(tǒng)MTF 值來近似遙感圖像的MTF 值,利用指數(shù)調(diào)節(jié)MTF 曲線來進(jìn)行圖像恢復(fù)[5],忽略了衛(wèi)星實(shí)際運(yùn)行環(huán)境和成像條件對MTF 的影響;張朋等人利用系統(tǒng)模型法估計(jì)MTF,此種方法沒有考慮大氣的影響,適用于獲取的圖像沒有特征地物的情況。上述方法通過估算MTF 在頻域內(nèi)對遙感圖像進(jìn)行恢復(fù),需要將MTF 對應(yīng)的成像系統(tǒng)的空間頻率轉(zhuǎn)換為傅里葉頻率,目前還沒有定論,使用受到一定的限制。
本文直接利用圖像中的特征地物估計(jì)成像系統(tǒng)的PSF,通過高斯擬合對PSF 進(jìn)行調(diào)整,使之更符合實(shí)際成像系統(tǒng)的情況。以擬合后的PSF 作為參數(shù),使用自適應(yīng)維納濾波進(jìn)行圖像恢復(fù)。在圖像恢復(fù)前,先進(jìn)行去噪預(yù)處理,去除遙感圖像上的周期條帶噪聲,增強(qiáng)圖像恢復(fù)的效果。
圖像的退化過程可以被模型化為一個(gè)退化函數(shù)和一個(gè)加性噪聲項(xiàng)。原始地物的圖像經(jīng)過成像系統(tǒng)和大氣發(fā)生了退化,在圖像的獲取和傳輸過程中又引入了大量噪聲,從而使最終獲得的圖像是降質(zhì)的圖像,存在著不同程度的模糊。
若退化函數(shù)是線性空間不變的,則圖像退化過程在空間域可表示為原圖像f(x,y) 和退化函數(shù)h(x,y) 的卷積并疊加噪聲n(x,y) ,形成退化圖像g(x,y) ,具體形式為:
由卷積定理可知,空間域兩個(gè)函數(shù)的卷積在頻率域中表示為它們傅里葉變換的乘積,因此,退化過程可在頻域表示為:
式中:F(u,v) ,H(u,v) ,N(u,v) ,G(u,v) 分別為f(x,y) ,h(x,y) ,n(x,y) ,g(x,y) 的二維傅里葉變換。
圖像恢復(fù)用于減輕或消除圖像獲取過程中發(fā)生的退化。如果知道退化函數(shù)和噪聲的相關(guān)知識,即可以采用與退化相反的過程恢復(fù)退化圖像,包括去噪和去卷積兩個(gè)過程。在頻率域中可表示為:
噪聲是一個(gè)隨機(jī)分量,很難確切估計(jì),圖像恢復(fù)會放大圖像中的噪聲,因此,圖像恢復(fù)要去除噪聲,消除噪聲對它的副作用。去噪既可以在空間域也可以在頻率域中進(jìn)行。常用的方法是根據(jù)噪聲的特點(diǎn),選擇合適的濾波器濾除噪聲。卷積運(yùn)算在頻率域可轉(zhuǎn)換為乘法運(yùn)算,其逆運(yùn)算在頻域中為除法運(yùn)算,因此,去卷積通常在頻域中進(jìn)行。最后,對得到的結(jié)果F(u,v) 進(jìn)行傅里葉反變換即可得到恢復(fù)圖像。
圖像恢復(fù)主要是提升圖像的高頻部分,增強(qiáng)圖像的邊緣細(xì)節(jié),而噪聲在圖像中也表現(xiàn)為高頻成分,恢復(fù)圖像的同時(shí)會使噪聲增大。因此,在圖像恢復(fù)之前進(jìn)行去噪處理,可以避免噪聲放大淹沒圖像細(xì)節(jié)。時(shí)間延遲積分( TDI) CCD 相機(jī)圖像附加的噪聲主要是條帶噪聲,表現(xiàn)為圖像中灰度變化緩慢的明暗相間豎條紋分布。本文的圖像恢復(fù)在頻域中進(jìn)行,因此,選擇在頻域進(jìn)行圖像去噪處理。對圖像進(jìn)行傅里葉變換時(shí),條帶噪聲在頻譜圖上表現(xiàn)為頻域橫軸上按一定規(guī)律分布的亮點(diǎn)。采用頻率域陷波濾波器來消除圖像中的條帶噪聲[7-8],過濾掉傅里葉頻譜中對應(yīng)條帶噪聲的高頻成分,然后,進(jìn)行傅里葉反變換即可得到去噪后的圖像。圖1( a) 是截取的部分遙感圖像,分布有明顯的垂直條帶噪聲,其對應(yīng)的頻譜圖為圖2( a) ;圖1( b) 是通過陷波濾波器進(jìn)行去噪處理后的圖像,其對應(yīng)的頻譜圖為圖2( b) 。通過兩圖對比可見,本文選取的去噪方法基本除去了圖像的條帶噪聲,有效地改善了圖像質(zhì)量。
圖1 圖像條帶噪聲去除Fig.1 Removing strip noise of images
圖2 圖像頻譜圖Fig.2 Fourier spectrum of the images
目前主要用3 種方法估計(jì)系統(tǒng)的PSF:
(1) 使用帶有特征地物的圖像( 例如飛機(jī)跑道、橋梁或人工靶標(biāo)圖像) 估計(jì)PSF;
(2) 高分辨率圖像對比法,需要衛(wèi)星上載有分辨率不同的多種成像儀;
(3) 根據(jù)系統(tǒng)設(shè)計(jì)標(biāo)準(zhǔn)和系統(tǒng)分析模型[9]推算PSF。
本文采用第一種方法,利用刀刃法估計(jì)成像系統(tǒng)的PSF。
由于獲取遙感圖像的條件復(fù)雜,直接估計(jì)得到的PSF 誤差很大,本文采用高斯擬合對直接估計(jì)的PSF 進(jìn)行校正,使其更符合實(shí)際的成像系統(tǒng)。
刀刃法是選取圖像中具有一定灰度差異的相鄰兩塊地物邊界作為刀刃邊緣圖像,兩塊邊界的灰度變化要均勻,并與軌道方向具有一定的傾角[10],原理如圖3 所示。詳細(xì)步驟如下:
(1) 從圖像中摳取具有刀刃邊緣特征的子圖像,對子圖像每一行提取邊緣分界點(diǎn);
(2) 采用最小均方誤差線性擬合邊緣分界點(diǎn),得到邊緣分界線,并以分界線與每一行的交點(diǎn)作為新的分界點(diǎn);
(3) 計(jì)算子圖像每一行上的其他像素相對于分界點(diǎn)的偏移量,并以0.05 的分辨率( 每兩個(gè)數(shù)據(jù)中插入20 個(gè)數(shù)據(jù)) 對偏移量數(shù)據(jù)進(jìn)行3 次樣條插值;
(4) 插值后所有行的偏移量數(shù)據(jù)取平均,得到子圖像的邊緣擴(kuò)散函數(shù)( ESF) ;
(5) 對邊緣擴(kuò)散函數(shù)微分得到線擴(kuò)散函數(shù)( LSF) 。
對于推掃式TDI-CCD 相機(jī),沿軌方向和跨軌方向的不對稱性很小,成像系統(tǒng)PSF 可以分解為沿軌和跨軌方向的兩個(gè)一維分量,分別對應(yīng)于沿軌和跨軌方向的LSF,用沿軌方向的LSF 列向量乘以跨軌方向的LSF 行向量,得到二維PSF。
圖3 刀刃法原理圖Fig.3 Schematic diagram of knife-edge method
成像系統(tǒng)的PSF 趨向于高斯型,可以表示為:
式中,σx,σy為標(biāo)準(zhǔn)差。PSF 在沿軌和跨軌兩個(gè)方向是可分的,即:
從圖3 中LSF 的形狀可見它近似為高斯型[11],但并不完全符合高斯表達(dá)式。為了減少PSF 估計(jì)的誤差,對得到的LSF 曲線進(jìn)行高斯擬合,使用擬合后的數(shù)據(jù)形成PSF 離散矩陣。
遙感圖像具有數(shù)據(jù)量龐大的特點(diǎn),對其進(jìn)行恢復(fù)的時(shí)間代價(jià)相對較高。因此,估計(jì)出成像系統(tǒng)的PSF 后,選擇維納濾波進(jìn)行圖像恢復(fù)。維納濾波性能好,沒有迭代過程,比Lucy-Richardson算法和正則濾波運(yùn)算量少,能夠以很低的計(jì)算代價(jià)獲得較好的復(fù)原效果,其頻域表達(dá)式可簡化為:
式中:H(u,v) 為高斯擬合后PSF 的二維傅里葉變換,K為恢復(fù)過程中的規(guī)整化值。
維納濾波采用固定的規(guī)整化值對整個(gè)圖像進(jìn)行恢復(fù)。較小的規(guī)整化值能突出圖像的邊緣特征,但在圖像的灰度連續(xù)區(qū)域中會產(chǎn)生寄生紋波,影響圖像恢復(fù)的效果。因此,這里對維納濾波進(jìn)行了改進(jìn),采用自適應(yīng)維納濾波算法進(jìn)行圖像恢復(fù)。自適應(yīng)維納濾波算法采用Sobel 算子對圖像進(jìn)行邊緣檢測,把圖像劃分為邊緣特征區(qū)域和灰度連續(xù)區(qū)域,在邊緣特征區(qū)域使用較小的K值,突出圖像細(xì)節(jié)的恢復(fù); 在灰度連續(xù)區(qū)域用較大的K值,避免產(chǎn)生寄生紋波。
遙感圖像中包含條帶噪聲,在圖像恢復(fù)之前進(jìn)行了有效的去噪處理。通過選取圖像中具有刀刃邊緣特征的地物估計(jì)擬合PSF,然后,分別采用維納濾波和自適應(yīng)維納濾波對整幅圖像進(jìn)行了恢復(fù)處理。圖4( a) 是截取了遙感圖像中的一個(gè)鐵塔頂部,其上分布有垂直的條帶噪聲。圖4( b) 是應(yīng)用維納濾波對去噪處理后的圖像進(jìn)行恢復(fù)處理的結(jié)果,由于進(jìn)行了去噪預(yù)處理,條帶噪聲的影響非常弱,頻域恢復(fù)中沒有受到增強(qiáng)噪聲的影響??梢娙ピ胩幚硐藞D像中的高頻條帶噪聲,有效地增強(qiáng)了圖像恢復(fù)的效果。由于維納濾波固有的不足,圖4( b) 相對于圖4( a) 雖然有一定效果,但由于存在寄生紋波,整個(gè)圖像的清晰度還不夠。圖4( c) 是采用自適應(yīng)維納濾波算法進(jìn)行恢復(fù)的結(jié)果,由于對不同的圖像區(qū)域采取不同的規(guī)整化值進(jìn)行恢復(fù),因此,恢復(fù)后圖像的細(xì)節(jié)更加突出,圖像清晰度得到增強(qiáng)。圖5( a) 是具有復(fù)雜地形的遙感圖像,圖5( b) 是先去除條帶噪聲,然后采用自適應(yīng)維納濾波算法進(jìn)行恢復(fù)的結(jié)果。原圖中明顯的條帶噪聲經(jīng)過去噪后基本消除,沒有影響最終恢復(fù)圖像的效果。
圖4 遙感圖像去噪與恢復(fù)Fig.4 Denoising and restoration of remote sensing images
圖5 復(fù)雜圖像自適應(yīng)維納濾波恢復(fù)Fig.5 Adaptive Wiener filtering restoration for complex images
除了在主觀上通過人眼觀察評價(jià)上述遙感圖像恢復(fù)效果外,還應(yīng)引入客觀評價(jià)標(biāo)準(zhǔn),對圖像的恢復(fù)效果進(jìn)行定量分析。圖像恢復(fù)消除了圖像的模糊現(xiàn)象,與原圖相比恢復(fù)圖像細(xì)節(jié)豐富,表現(xiàn)為相鄰像素的灰度變化較大,具有更大的梯度值,相關(guān)的圖像梯度函數(shù)可作為評價(jià)圖像恢復(fù)效果的標(biāo)準(zhǔn),具體如下:
(1) 方差:圖像的方差大,說明圖像灰度層次較為豐富。M、N為圖像的高度和寬度,f(i,j) 為圖像某點(diǎn)的灰度值,u為圖像的灰度均值,u =
(2) 灰度平均梯度: 反映圖像在水平和垂直方向的灰度變化。灰度平均梯度越大,表明圖像越清晰。其中:
灰度平均梯度計(jì)算公式為:
(3) Laplacian 梯度: 反映圖像每一像素附近的灰度變化。Laplacian 梯度值越大表明圖像輪廓越鮮明,相應(yīng)地圖像越清晰。其中:
由圖4 可見,待處理的遙感圖像細(xì)節(jié)不清,鐵塔的支撐結(jié)構(gòu)邊緣非常模糊。表1 的數(shù)據(jù)也驗(yàn)證了這一點(diǎn),未恢復(fù)之前的圖像各項(xiàng)評價(jià)標(biāo)準(zhǔn)值均比恢復(fù)后的圖像低很多。由于圖像恢復(fù)處理增強(qiáng)了圖像細(xì)節(jié)和邊緣特征,由表1 可見,恢復(fù)后圖像的方差和灰度平均梯度明顯增大,表明恢復(fù)后圖像包含的信息更多,圖像更清晰。結(jié)合表1,對比圖4( b) 與圖4( c) 的恢復(fù)效果可見,自適應(yīng)維納濾波比直接維納濾波的恢復(fù)圖像清晰程度高,各項(xiàng)指標(biāo)均有很大提高,采用維納濾波恢復(fù)后方差比原圖4( a) 方差增大1.866,灰度平均梯度增大1.289,Laplacian 梯度增大8.21;采用自適應(yīng)維納濾波恢復(fù)后,圖4( c) 方差比原圖4( a) 方差增大4.395,灰度平均梯度增大1.799,Laplacian梯度增大10.014,說明采用自適應(yīng)維納濾波算法可以獲得更好的圖像恢復(fù)效果。對圖5 所示的復(fù)雜遙感圖像,自適應(yīng)維納濾波同樣取得了很好的恢復(fù)效果。由表2 可知,恢復(fù)后圖像方差比原圖方差增大20.753,灰度平均梯度增大1.766,Laplacian 梯度增大5.018。
表1 遙感圖像恢復(fù)評價(jià)Tab.1 Evaluation of remote sensing image
表2 復(fù)雜遙感圖像恢復(fù)評價(jià)Tab.2 Evaluation of complex remote sensing image
針對航天遙感相機(jī)成像過程中圖像退化作用引起的圖像模糊,提出了高斯擬合的PSF 估計(jì)方法,并采用自適應(yīng)維納濾波進(jìn)行圖像恢復(fù)。考慮圖像恢復(fù)過程會放大噪聲,因此,在恢復(fù)之前先對圖像進(jìn)行去噪預(yù)處理,采用頻率域陷波濾波器消除了遙感圖像中附加的條帶噪聲。去噪后,根據(jù)所處理遙感圖像的特點(diǎn),利用圖像中的特征地物,采用刀刃法估計(jì)成像系統(tǒng)的點(diǎn)擴(kuò)散函數(shù)PSF,為了保證PSF 估計(jì)的精度,對刀刃法得到的PSF 進(jìn)行高斯擬合,將擬合后的PSF 應(yīng)用到自適應(yīng)維納濾波中進(jìn)行圖像恢復(fù)。實(shí)驗(yàn)結(jié)果表明: 去噪處理有效地消除了圖像中的條帶噪聲,自適應(yīng)維納濾波比直接維納濾波取得了更好的圖像恢復(fù)效果。細(xì)節(jié)圖像經(jīng)自適應(yīng)維納濾波恢復(fù)后方差比原圖方差增大4.395,灰度平均梯度增大1.799,Laplacian 梯度增大10.014; 復(fù)雜圖像經(jīng)自適應(yīng)維納濾波恢復(fù)后,方差比原圖方差增大20.753,灰度平均梯度增大1.766,Laplacian 梯度增大5.018。恢復(fù)后圖像中的細(xì)節(jié)更清晰,特征更易識別,利于圖像的判讀和分析。
[1] 許元男,趙元,劉麗萍,等.含噪聲模糊圖像的點(diǎn)擴(kuò)展函數(shù)參數(shù)辨識[J].光學(xué) 精密工程,2009,17(11) :2849-2856.XU Y N,ZHAO Y,LIU L P,et al.. Parameter identification of point spread function in noisy and blur images[J].Opt.Precision Eng.,2009,17(11) :2849-2856.( in Chinese)
[2] 郭永彩,王婀娜,高潮.空間自適應(yīng)和正則化技術(shù)的盲圖像復(fù)原[J].光學(xué) 精密工程,2008,16(11) :2263-2267.GUO Y C,WANG E N,GAO CH. Blind image restoration algorithm based on space-adaptive and regularization[J].Opt.Precision Eng.,2008,16(11) :2263-2267.( in Chinese)
[3] 江潔,鄧瓊,張廣軍.基于小波變換的正則化盲圖像復(fù)原算法[J].光學(xué) 精密工程,2007,15(4) :582-586.JIANG J,DENG Q,ZHANG G J. Regularization algorithm for blind image restoration based on wavelet transform[J].Opt.Precision Eng.,2007,15(4) :582-586.( in Chinese)
[4] 顧行發(fā),李小英,閔祥軍,等.CBERS-02 衛(wèi)星CCD 相機(jī)MTF 在軌測量及圖像MTF 補(bǔ)償[J].中國科學(xué)E 輯 信息科學(xué),2005,35:26-40.GU H F,LI X Y,MIN X J,et al.. On-orbit MTF Estimation and MTF compensation of CCD camera in CBERS-02 satellite[J].Science in China,Series E Information Sciences,2005,35:26-40.( in Chinese)
[5] 陳強(qiáng),戴奇燕,夏德深.基于MTF 理論的遙感圖像復(fù)原[J].中國圖像圖形學(xué)報(bào),2006,11(9) :1299-1305.CHEN Q,DAI Q Y,XIA D SH. Restoration of remote sensing images based on MTF theory[J].J. Image and Graphics,2006,11(9) :1299-1305.( in Chinese)
[6] 張朋,劉團(tuán)結(jié),王宏琦.線陣CCD 相機(jī)MTF 的系統(tǒng)模型估計(jì)法與圖像復(fù)原[J].光學(xué)技術(shù),2009,35(3) :394-398.ZHANG P,LIU T J,WANG H Q. MTF estimation based on system modal for linear CCD camera and image recovery[J].Opt. Technique,2009,35(3) :394-398.( in Chinese)
[7] GONZALEZ R C,WOODS R E. 數(shù)字圖像處理[M].2 版.北京:電子工業(yè)出版社,2007.GONZALEZ R C,WOODS R E.Digital Image Processing[M]. 2nd ed. Beijing:Publishing House of Electronics Industry,2007.( in Chinese)
[8] GONZALEZ R C,WOODS R E,EDDINS S L. 數(shù)字圖像處理:MATLAB 版[M].北京:電子工業(yè)出版社,2005.GONZALEZ R C,WOODS R E,EDDINS S L.Digital Image Processing Using MATLAB[M]. Beijing:Publishing House of Electronics Industry,2005.( in Chinese)
[9] BENSEBAA K,BANON G J F,F(xiàn)ONSECA L M G. On-orbit spatial resolution estimation of CBERS-1 CCD system from bridge images[C]. Proceedings of the XXth International Society for Photogrammetry and Remote Sensing( ISPRS) ,Istanbul,Turkey,July 2004:36-41.
[10] CHOI T. IKONOS satellite on orbit Modulation Transfer Function( MTF) measurement using edge and pulse method[D]. Brookings:South Dakota State University,2002.
[11] 周松濤,宣家斌.基于景物灰度分布特征的影像恢復(fù)技術(shù)[J].武漢測繪科技大學(xué)學(xué)報(bào),1999,24(3) :230-234,239.ZHOU S T,XUAN J B. Restoration technique based on image features[J].J. Wuhan Technical University of Surveying and Mapping,1999,24(3) :230-234,239.( in Chinese)