本發明屬于數字圖像處理技術領域,具體涉及一種基于總曲率的SAR圖像變分去噪方法。
背景技術:
:
相干斑噪聲是合成孔徑雷達(Synthetic Aperture Radar,簡稱SAR)圖像的重要特征,嚴重影響SAR圖像的可解譯性。相干斑噪聲通常作為乘性噪聲來建模f=uη,f為觀察到的退化圖像,u為原始清晰圖像,η為噪聲。理想的SAR圖像去噪方法是能在去除斑點噪聲的同時保留圖像的邊緣和細節信息。去噪方法為定義一個濾波器窗口估計相干斑的局部噪聲方差,利用估計值進行濾波處理。現有技術中,常用濾波算法有均值濾波、中值濾波、局部濾波、Lee濾波、Lee-Sigma濾波、Frost濾波和Gamma-MAP濾波。研究表明,在均勻圖像區域,這些濾波方法能夠較好的削減噪聲,而在非均勻區域,圖像過于平衡或模糊,不能很好的保持邊緣細節信息;現有技術中,全變分方法將圖像去噪構建為能量函數的最小化問題,引入各項異性擴散方程,在平滑噪聲的同時保持邊緣。變分模型包括數據保真項和規則項,基于全變分TV規則項,AA模型是最早的SAR圖像Gamma分布的乘性噪聲去除模型,SST模型為Poisson分布的乘性噪聲去除模型,DTDS模型為Rayleigh分布的乘性噪聲去除模型,SO模型為綜合的乘性噪聲去除模型。常用的梯度下降法在求解乘性噪聲去噪模型時存在計算速度慢的問題,通常在求解過程中引入一些快速算法包括Split Bregman算法,對偶算法等。TV規則項能夠較好的保持邊緣,但階梯效應是其主要缺點,通常引入高階項來克服這一缺點,主要包括Hessian矩陣規則項、拉普拉斯Δu規則項和總曲率規則項。由于計算的復雜性和非線性,目前乘性噪聲變分模型還未引入高階規則項,因此設計一種基于總曲率的SAR圖像變分去噪方法,能夠將待處理的圖像既能平滑乘性噪聲又能保圖像邊緣細節信息。
技術實現要素:
:
本發明的目的在于克服現有方法存在的缺陷,尋求設計提供一種基于總曲率的SAR圖像變分去噪方法,該方法涉及的變分能量方程包括數據保真項和總曲率規則項,并且基于交替方向乘子法(Alternating Direction Method of Multipliers,縮寫為ADMM)巧妙設計輔助變量,通過L2范數約束,實現能量方程最小化極值問題的求解,求解的圖像既能平滑乘性噪聲又能保圖像邊緣細節信息。
為了實現上述目的,本發明涉及的基于總曲率的SAR圖像變分去噪方法的具體操作方法按照如下步驟進行:
a.選擇一幅待處理的原始SAR圖像f并根據該圖像f建立基于總曲率的SAR圖像變分去噪能量方程,對于輸入的原始超聲圖像f,期望得到的去噪后的圖像為u,基于總曲率的能量方程為:
其中,Ω為SAR圖像區域,α為權重系數,a、b和c為成型噪聲一階項、平方項和對數項的參數;曲率規則項的使用能夠在SAR圖像去噪過程中更好的保持邊緣細節信息;
b.對步驟a中所述的總曲率的能量方程進行轉換,步驟a建立的能量方程的數據項和規則項均為非凸非線性,因此引入u=ez進行變量替換,能量方程轉換如下:
c.步驟b建立的能量方程具有高階、非凸性,難以進行求解,引入分裂算子w、和q簡化總曲率規則項,步驟b的能量方程形式化為帶約束的極小值問題:
d.步驟c中所述的約束能夠轉換為和兩個等價約束,設計約束因此,又被轉化為和具有約束的變量是松弛的,至此,能量方程轉化為可使用增廣拉格朗日方法求解的方程:
e.步驟d中所述的約束w=z、和采用L2懲罰項,由能夠推導出因此使用L1懲罰約束這樣極小化問題轉換為以下子問題的交替優化問題:
其中,β1、β2、β3、β4和β5是正的懲罰參數,λ1、λ2、λ4和是拉格朗日乘子,能夠根據相應規則更新;
f.利用變量交替迭代優化求解分別計算步驟e中的變量z,w,q,將步驟e的極小化問題轉換為以下6個子問題:
g.分別求解步驟f中的ε1(z)、ε2(w)、ε4(q)、和的歐拉方程;ε2(w)的歐拉方程采用梯度降方法直接求解,和ε4(q)的歐拉方程采用廣義軟閾值公式求解,的歐拉方程能夠直接采用投影方法,而ε1(z)和的歐拉方程為非線性,采用快速傅里葉變換方法進行求解;
h.對步驟g中的z,w,q,進行迭代求解,當相鄰兩次迭代的能量差小于設定的閾值時停止;
i.采用u=ez得到的u即為去噪后的SAR圖像。
本發明與現有技術相比,利用總曲率規則項進行SAR圖像去噪,對于利用總曲率項建立的能量方程為了避免在求解時所產生的復雜運算,同時巧妙設計約束,引入輔助變量進行求解,不但提高了效率,而且減少了計算的復雜度,同時本發明提出的基于總曲率的SAR圖像去噪方法具有非常好的實際應用價值,對于提高SAR圖像的清晰度,提高圖像的解譯度起到了非常重要的作用,應用價值極高,市場前景廣闊。
附圖說明:
圖1為本發明涉及的基于總曲率的SAR圖像去噪方法流程圖。
圖2為本發明涉及的在圖像SAR-1得到的結果與AA模型的比較,其中圖2(a)為原始SAR-1圖像,圖2(b)為基于本發明α=0.5得到的去噪結果圖,圖2(c)為基于本發明α=1得到的去噪結果,圖2(d)為基于AA模型α=0.5得到的去噪結果圖。
圖3為本發明涉及的在圖像SAR-2得到的結果與SST模型的比較,其中圖3(a)為原始SAR-2圖像;圖3(b)為基于本發明α=0.5得到的去噪結果圖,圖3(c)為基于本發明α=1得到的去噪結果圖,圖3(d)為基于SST模型α=0.5得到的去噪結果圖。
圖4為本發明涉及的在圖像SAR-3得到的結果與DTDS模型的比較圖,其中圖4(a)為原始SAR-3圖像,圖4(b)為基于本發明α=0.5得到的去噪結果圖,圖4(c)為基于本發明α=1得到的去噪結果圖,圖4(d)為基于DTDS模型α=0.5得到的去噪結果圖。
圖5為基于本發明開發的圖像去噪應用程序,程序運行包括圖像灰度值動態模擬和圖像|能量值動態模擬兩種方式,其中圖5(a)為主界面圖,圖5(b)為圖像灰度值動態模擬圖,圖5(c)圖像|能量動態模擬圖,圖5(d)結果輸出圖。
圖6為本發明涉及的圖像SAR-1的灰度值時空變化結果圖,其中圖6(a)為原始SAR-1圖像,圖6(b)為基于本發明α=0.5得到的去噪結果的灰度值三維圖,圖6(c)為基于本發明α=1得到的去噪結果灰度值三維圖,圖6(d)為基于AA模型α=0.5得到的去噪結果的灰度值三維圖。
具體實施方式:
下面結合附圖和具體實施方式對本發明做進一步說明:
實施例1:
本實施對SAR圖像變分去噪時,具體操作方法按照如下步驟進行:
a.選擇一幅待處理的原始超聲圖像f并根據該圖像f建立基于總曲率的SAR圖像變分去噪能量方程,對于輸入的原始超聲圖像f,期望得到的去噪后的圖像為u,基于總曲率的能量方程為:
其中,Ω為SAR圖像區域,α為權重系數,a、b和c為成型噪聲一階項、平方項和對數項的參數;
b.對步驟a中所述的總曲率的能量方程進行轉換,步驟a建立的能量方程的數據項和規則項均為非凸非線性,因此引入u=ez進行變量替換,能量方程轉換如下:
c.將步驟b建立的能量方程進行求解,引入分裂算子w、和q簡化總曲率規則項,步驟b的能量方程形式化為帶約束的極小值問題:
d.步驟c中所述的約束能夠轉換為和兩個等價約束,設計約束因此,又被轉化為和具有約束的變量是松弛的,至此,能量方程轉化為可使用增廣拉格朗日方法求解的方程:
e.步驟d中所述的約束w=z、和采用L2懲罰項,由能夠推導出因此使用L1懲罰約束這樣極小化問題轉換為以下子問題的交替優化問題:
其中,β1、β2、β3、β4和β5是正的懲罰參數,λ1、λ2、λ4和是拉格朗日乘子,能夠根據相應規則更新;
f.利用變量交替迭代優化求解分別計算步驟e中的變量z,w,q,將步驟e的極小化問題轉換為以下6個子問題:
g.分別求解步驟f中的ε1(z)、ε2(w)、ε4(q)、和的歐拉方程;ε2(w)的歐拉方程采用梯度降方法直接求解,和ε4(q)的歐拉方程采用廣義軟閾值公式求解,的歐拉方程能夠直接采用投影方法,而ε1(z)和的歐拉方程為非線性,采用快速傅里葉變換方法進行求解;
h.對步驟g中的z,w,q,進行迭代求解,當相鄰兩次迭代的能量差小于設定的閾值時停止,具體迭代步驟如下:
(1)初始化參數z=logf,w=z,(β1,β2,β3,β4,β5,Δt,iteration)>0,根據噪聲分布函數確定參數a,b和c;
(2)固定wk,和求解ε1(z)的歐拉方程,采用快速傅里葉變換(Fast Fourier Transform,縮寫為FFT)方法求z;歐拉方程為
上述方程能夠形式化為如下方程:
步驟(2)中的歐拉方程的離散形式為:
(β5-β2(S1++S1--2I+S2++S2--2I))z(i,j)=g(i,j),
將離散后的歐拉方程采用離散傅里葉變換(Discrete Fourier Transform,縮寫為DFT)變換,得到如下方程:
由于dz=(β1-2β3(coszi+coszj-2))>0,采用反傅里葉變換求得z,其中
(3)固定zk+1和λ1k,求解ε2(w)的歐拉方程,先采用梯度降方法求w,其中歐拉方程為α(afe-w+bf2e-2w-c)+β1(zk+1-w)+λ1k=0;在采用梯度降求解,具體方式如下:
(4)固定zk+1、λ2k和求解的歐拉方程,采用廣義軟閾值公式求其中歐拉方程的計算方式如下:
廣義軟閾值公式求解
(5)固定和λ4k.,求解ε4(q)的歐拉方程,采用廣義軟閾值公式求q;歐拉方程為廣義軟閾值公式求解
(6)固定qk+1,和λ4k,求解的歐拉方程,采用FFT求歐拉方程為
上述歐拉方程依據移位算子表達如下:
其中
再將上述移位算子表達后的方程采用FFT變換,得到如下方程:
系數為:
確保β4β5>0,行列式D=β5-2β4β5(coszi+coszj-2)即大于0,采用反傅里葉變換求得
(7)固定λ2k和求解的歐拉方程,采用投影法求
歐拉方程為
投影法求解
i.采用u=ez得到的u即為去噪后的SAR圖像。