李率杰,李 鵬,馮兆永，姚正安

(中山大學數(shù)學與計算科學學院，廣東 廣州 510275)

圖像修復技術(shù)是指針對圖像中遺失或者損壞的部分,利用周圍未被損壞的圖像信息,按照一定的規(guī)則進行填補,使修復后的圖像接近或達到原圖的視覺效果，在本質(zhì)上是一種典型的插值技術(shù)。圖像修復具有廣泛的應用前景，如填補美術(shù)作品上出現(xiàn)的裂痕，以使其清晰并恢復它的完整性；恢復舊照片中的劃痕；去除圖像中的文字或物體等特效制作；視頻通信中錯誤隱匿；圖像遮擋物體的去除等。

近來的數(shù)字圖像修復技術(shù)大多是利用偏微分方程(partial differential equation，PDE)和變分法。而這種基于PDE的數(shù)字圖像修復技術(shù)最早是由Bertalmio等[1]首先引入到圖像處理領(lǐng)域中的。在[1] 中，他們利用待修復區(qū)域的周圍信息，采用一種由粗到精的方法來估計等照度線的方向，并采用傳播機制將信息傳播到待修復區(qū)域內(nèi)。之后Balleater等[2]通過提出一個能量泛函并利用變分法極小化該能量得到非線性的偏微分方程計算出修復后的圖像。在文獻[3-4]中作者考慮圖像的等照度線方向和不可壓Navier-Stokes方程之間的聯(lián)系將圖像灰度作為二維不可壓流體的流函數(shù)，將Laplacian算子作用后的待修復圖像作為流的渦流函數(shù)?；谟闪骱瘮?shù)定義的向量場將渦流函數(shù)傳播到待修復區(qū)域，從而得到修復后的圖像。由此Tai等[5]提出了基于TV-Stokes 方程的圖像修復兩步算法,首先通過求解一個非線性TV-Stokes方程將等照度線延伸到待修復區(qū)域,然后圖像就以此等照度線方向進行修復。

不同于以上方法，Chan等[6]引入了Rudin等[7]提出的全變差(total variation，TV)圖像去噪模型的思想，得到基于TV模型的圖像修復模型，該方法通過求解最小化能量泛函而得到非線性方程，無論待修復圖像是否帶有噪聲都能夠有效地對圖像進行修復，得到了較好的效果。但是TV模型在圖像平滑區(qū)域容易產(chǎn)生階梯效應。在后面的工作如文獻[8-10]中，Chan和Shen提出了基于圖像等照度線的曲率的能量泛函，試圖在保證光滑條件下連接修復區(qū)域周圍的等照度線，由此得到的方程是非線性的四階方程，該方法可以較好的修復面積較大的非紋理區(qū)域，但是非線性的高階方程給求解帶來了困難?；谝陨螩han和Shen的工作，Esedoglu等[11]基于Mumford-Shah圖像分割模型提出了Mumford-Shah-Euler圖像修復模型[12]，利用Ambrosio和Tortorelli提出的Mumford-Shah模型橢圓逼近方法得到相應的拋物方程。方法的優(yōu)點是方程的最高階偏導是線性的可以快速的計算。但是卻有難以修復較大的區(qū)域的缺點。Grossauer等[13]利用復Ginzburg-Landau方程解決圖像修復問題，近來Bornemann等[14]提出基于一階傳輸方程的圖像修復的非迭代方法。但是這些方法計算量大且非常復雜。

另一方面紋理圖像的修復也受到廣泛的研究，在文[15]中，作者首先將待修復區(qū)域分解為紋理和結(jié)構(gòu)，然后對紋理采用紋理合成的方法，對結(jié)構(gòu)采用一般的非紋理圖像修復模型進行分別修復。另外最近紋理合成的方法也同樣被廣泛的運用在紋理圖像修復中[16-20]。

本文基于文獻[3-5]中的思想，我們在等照度線方向上向待修復區(qū)域傳輸圖像信息。在文獻[3]通過求解不可壓流體的Navier-Stokes方程，得到修復后的圖像函數(shù)。假設D為待修復區(qū)域，?D為待修復區(qū)域的邊界。一旦給定待修復區(qū)域D此模型就可以將信息沿等照度方向連續(xù)的從邊界?D傳輸?shù)酱迯蛥^(qū)域D。如此圖像修復依據(jù)連接和待修復區(qū)域邊界相交的等照度線得以完成。但是該修復算法只是基于待修復區(qū)域的邊界?D上的信息，因此難以得到與待修復區(qū)域相交的等照度線的準確方向。另外如果待修復區(qū)域D是規(guī)則的矩形區(qū)域，算法容易實現(xiàn)，但是實際上待修復區(qū)域往往是形狀不規(guī)則的區(qū)域，實現(xiàn)此算法將面臨執(zhí)行有限差分格式和邊界條件的困難，最后該修復算法對噪聲圖像的修復無能為力。

與上面算法不同的是，本文修復算法的思想是能夠在待修復區(qū)域進行信息填充，而在待修復區(qū)域外異性光滑圖像去除噪聲(如果存在)，且保持邊緣信息。在本文中，基于不可壓流的Navier-Stokes方程，我們主要提出一種新的圖像修復算法，目的是任給一個帶有噪聲的且信息丟失或被破壞的待修復圖像，我們可以同時的進行圖像的修復和去噪，得到清晰、完整的圖像。

1 基于Navier-Stokes方程的圖像修復

不可壓的流體一般遵循如下所謂的Navier-Stokes方程

υt+υ·▽υ=-▽p+νΔυ,▽·υ=0

(1)

其中υ為流體的速度場，p為流體的壓強，ν為黏性系數(shù)。對于二維空間，我們引進流函數(shù)Ψ，滿足

▽⊥Ψ=υ

(2)

同時設渦量ω=▽×υ，且由(1)可得滿足渦流方程

ωt+υ·▽ω=νΔω

(3)

在二維空間里，渦量ω是一個標量，并且滿足ω=ΔΨ。當黏性系數(shù)ν=0時，可以得到無黏性流的歐拉方程。在二維時黏性流或無黏性流，方程都是適定的，對于任意光滑初始條件下解都是存在的，并且連續(xù)的依賴初始條件和邊界條件[21]。

根據(jù)流函數(shù)的定義，黏性流方程(3)的穩(wěn)定態(tài)一定滿足下列條件

▽⊥Ψ·▽ΔΨ≈0

(4)

此式說明Laplacian算子作用后流函數(shù)ΔΨ，其等照度線的方向和流體的流速方向是平行的。即流體是沿著ΔΨ等照度線方向流動，如此其和在圖像修復中將周圍信息沿等照度線方向傳播到待修復區(qū)域的原理相同。

將二維不可壓流體力學中的流函數(shù)Ψ和圖像修復中待修復圖像的灰度函數(shù)I0對應起來，Bertalmi等[3]提出了基于Navier-Stokes方程的圖像修復模型。設Ω為整個圖像區(qū)域，D為待修復區(qū)域，?D為待修復區(qū)域的邊界，待修復圖像I0在ΩD是足夠光滑，并且I0和ΔI0在ΩD是已知的。用圖像灰度函數(shù)I對應流函數(shù)Ψ，將二維不可壓流體理論和圖像修復理論對應關(guān)系總結(jié)如下：

Bertalmio等[3]目的是在待修復區(qū)域上解決Navier-Stokes方程的一種形式，即渦流方程(3)。他們解決如下方程：

表1 流體力學和圖像修復物理量對應關(guān)系

ωt+υ·▽ω=ν▽·(g(|▽ω|)▽ω)

(5)

其中渦量ω=ΔI，函數(shù)g為單調(diào)遞減函數(shù)，可以使函數(shù)ω異性擴散，速度場υ=▽⊥I。由求解方程(5)所得的渦量ω,通過求解Poisson方程：

ΔI=ω，I?D=I0

(6)

從而得到修復后的圖像I。

2 圖像修復算法

本文中，我們的目的是在整個區(qū)域上求解Navier-Stokes方程的一種不同于(5)式的一種形式，可以使得在修復區(qū)域D內(nèi)進行圖像修復，而在修復區(qū)域外ΩD進行圖像去噪。我們求解下列渦流方程：

(7)

(8)

其中χD為修復區(qū)域D的特征函數(shù)，定義為χD=1,x∈D，χD=0,x∈ΩD。由方程求得的ω，我們同樣可以在區(qū)域Ω求與(6)式類似的Poisson方程，得到修復后的圖像I。

由提出的模型,可以看出它能夠在修復區(qū)域的內(nèi)部和外部起著不同的作用。其修復過程包括：修復區(qū)域內(nèi)，基于周圍信息對破壞或丟失的信息進行填充；而在修復區(qū)域外部，對ω進行異性光滑，從而可以達到去噪目的。

3 數(shù)值計算

對方程數(shù)值計算，我們假設m×n的圖像I是在連續(xù)區(qū)域[0,m]×[0,n]上，離散為m×n網(wǎng)格，定義時間步長Δt。我們用有限差分半隱格式逼近此問題。

Δ

Δtν▽·(g(|▽ωk|)▽ωk+1)i,j+

Δtλ(1-χD)(ω0-ωk+1)i,j

(9)

▽·(g(|▽ωk|)▽ωk+1)i,j=

(10)

(11)

其中K是取定的常數(shù)。

對于Poisson方程ΔIk+1=ωk+1，我們采用五點中心離散格式，將Ik+1和ωk+1重組成mn的向量，得到線性方程組AIk+1=ωk+1，這里A為一個稀疏的mn×mn的矩陣，求解這個線性方程組由ωk+1我們快速的得到修復后的圖像Ik+1。詳細的計算流程如下：

3)求解Poisson方程得到Ik；

其中k=1,2,…，直到穩(wěn)定狀態(tài)，得到的Ik既是修復后的圖像。

4 仿真實驗

在本節(jié)將用數(shù)值實驗結(jié)果展示本模型的有效性和實用性。我們選了不同類型的圖像，他們有若干區(qū)域丟失信息，劃痕，或者同時帶有噪聲，我們對其進行修復。

圖1 圖像修復結(jié)果比較

實驗2 本文算法針對兩幅合成圖像修復結(jié)果。如圖2所示， (a)和(d)分別為原始圖像，(b)和(e)為丟失信息的圖像，為了清晰顯示修復區(qū)域，我們同樣進行了隨機像素填充。通過修復，我們可以分別得到它們修復后的圖像(c)和(f)，從結(jié)果可以看出修復效果非常明顯。由于二值合成圖像灰度對比較強，在邊緣處和原始圖像有所差別，但是圖像信息基本得到了恢復。本實驗中參數(shù)選擇如表2所示。

表2 實驗中參數(shù)選擇

圖2 圖像修復結(jié)果

實驗3 本文算法針對兩幅真實圖像的修復結(jié)果，所用參數(shù)和迭代次數(shù)見表2。在圖3中，(a)為一幅原始圖像，(b)為其帶有劃痕的圖像，(c)為修復后的圖像。從效果上來看修復圖像和原始圖像非常相近，可以證明本文算法對這種劃痕是有效的。(d)為選擇的另一幅原始圖像，(e)為在其不同位置有信息缺失的圖像，其缺失區(qū)域大多分布在圖像邊緣上，(f)為修復后的圖像。我們可以看出修復結(jié)果仍然十分理想。

圖3 圖像修復結(jié)果

實驗4 本文算法同時去除噪聲和修復圖像信息的仿真實驗。本實驗利用本章所提的修復算法同時去除噪聲和修補圖像中缺失的信息，結(jié)果如圖4所示。我們同樣選擇了兩幅圖像，其中(a)和(d)為原始圖像，(b)和(e)為破損且加噪圖像，(c)和(f)為本文算法修復的圖像。從圖中可以看出，本文算法可以有效地修復破損的噪聲圖像。

圖4 圖像修復結(jié)果

實驗5 本文算法移除圖中不需要的景物，修復圖畫的仿真實驗。結(jié)果如圖5所示，(a)和(d)為原始圖像，(b)和(e)是在圖中標出修復區(qū)域的圖像，(c)和(f)為修復結(jié)果。從修復后的圖像可以看出，雖然移除景物面積較大但是本文算法仍然得到了很好的修復效果，使觀察者難以察覺圖像曾經(jīng)被修改過。但是修復后的圖像亮度發(fā)生了輕微的變化，相關(guān)計算過程值得我們研究改進。

圖5 圖像修復結(jié)果

5 結(jié) 論

本文基于Navier-Stokes方程，提出了一個新的圖像修復的數(shù)學模型。該模型能夠在待修復區(qū)域?qū)D像丟失或破損的信息進行修補，而在待修復區(qū)域外去除噪聲。

該模型的思想來自于不可壓流體力學的Navier-Stokes方程的物理背景，以及流體速度方向和圖像的等照度線方向的之間的一致性。對于圖像信息丟失或損壞，并且?guī)в性肼暤膱D像，修復將被自適應方程完成，在修復區(qū)域內(nèi)表現(xiàn)為信息由邊界向內(nèi)部傳播，在修復區(qū)域外表現(xiàn)為異性光滑。實驗結(jié)果顯示了模型在處理破損圖像(圖2，圖3和圖5)或噪聲圖像(圖4)時的有效性和實用性。但是模型也存在某些不足，如在保持修復區(qū)域圖像邊緣方面如圖2，圖4，如何避免這些缺點將是我們以后繼續(xù)研究的工作。

參考文獻：

[1]BERTALMIO M,SAPIRO G,BALLESTER C,et al.Image inpainting[C]∥Siggraph 2000,Computer Graphics proceedings,K.Akeley,Ed ACM Press/ ACM SIGGRAPH/ Addison Wesley Longman,2000：417-424.

[2]BALLESTER C,BERTALMIO M,CASELLES V.Filling in by joint interpolation of vector fields and gray levels[J].IEEE Trans Image Processing,2001,10(8): 1200-1211.

[3]BERTALMIO M,BERTOZZI A L,SAPIRO G.Navier-Stokes,fluid dynamics and image and video inpainting[J].IEEE Computer Vision and Pattern Recognition (CVPR),2001,1: 355-362.

[4]WILSON A,RYO T.Inpainting with the Navier-Stokes equations [R].http://www.math.ucla.edu/～rrtakei/gradProj/930project.pdf.

[5]TAI X C,OSHER S,HOLM R.Image inpainting using TV-Stokes equation[C]∥Image Processing Based on Partial Differential Equations,Springer,Heidelberg,2006:3-22.

[6]CHAN T F,SHEN J.Mathematical models for local non-texture inpaintings[J].SIAM J Appl Math,2002,62(3): 1019-1043.

[7]RUDIN L,OSHER S，F(xiàn)ATEMI E.Nonlinear total variation based noise removal algorithms[J].Physica D,1992,60:259-268.

[8]CHAN T F,KANG S H,SHEN J.Euler’s elastica and curvature-based inpainting[J].SIAM J Appl Math,2002,63(2): 564-592.

[9]CHAN T F,SHEN J,VESE L.Variational PDE models in image processing[J].Notices Am Math Soc,2003,50(1): 14-26.

[10]CHAN T,SHEN J.Non-texture inpainting by curvature-driven diffusions (CDD)[J].J Visual Comm.Image Rep,2001,12(4):436-449.

[11]ESEDOGLU S,SHEN J.Digital inpainting based on the Mumford-Shah-Euler image model[J].European Journal of Applied Mathematics,2002,13:353-370.

[12]MUMFORD D,SHAH J.Optimal approximations by piecewise smooth functions and associated variational problems[J].Communications on Pure and Applied Mathematics,1989,42: 577-685.

[13]GROSSAUER H,SCHERZER O.Using the complex Ginzburg-Landau equation for digital inpainting in 2d and 3d,scale space methods in computer vision[J].Lecture Notes in Computer Science 2695,2003: 225-236.

[14]BORNEMANN F,MRZ T.Fast image inpainting based on coherence transport[J].Journal of Mathematical Imaging and Vision,2007,28(3): 259-278.

[15]BERTALMIO M,VESE L,SAPIRO G,et al.Simultaneous texture and structure image inpainting[J].IEEE Tran Image Process.2003,12(8): 882-889.

[16]CRIMINISI A,PEREZ P,TOYAMA K.Region filling and object removal by exemplar-based image inpainting[J].IEEE Trans Image Processing,2004,13(9): 1200-1212.

[17]WEXLER Y,SHECHTMAN E,IRANI M.Space-time video completion[J].IEEE Trans Pattern Anal Mach Intell,2007,29(3): 1463-1476.

[18]FADILI M,STARCK J,MURTAGH F.Inpainting and zooming using sparse representations[J].Comput J,2009,52(1):64-79.

[19]KOMODAKIS N,TZIRITAS G.Image completion using efficient belief propagation via priority scheduling and dynamic pruning[J].IEEE Transactions on Image Processing,2007,16(11): 2649-2661.

[20]AUJOL J,LADJAL S,MASNOU S.Exemplar-based inpainting from a variational point of view[J].SIAM Journal on Mathematical Analysis,2010,42(3): 1246-1285.

[21]MAJDA A,BERTOZZI A.Vorticity and incompressible flow [M].Cambridge Univ Press,2001.

[22]TAI X.Global extrapolation with a parallel splitting method[J].Numerical Algorithm,1991,3:527-440.

[23]WEICKERT J,ROMENY B,VIERGEVER M.Efficient and reliable schemes for nonlinear diffusion filtering[J].IEEE Trans on Image Processing,1998,7(3):398-410.

[24]KUHNE G,WEICKERT J,VIERGEVER M.Fast implicit active contours models lecture[J].Notes on Computer Science,2002,2449: 133-140.