偏微分程數(shù)值求解以及運(yùn)算結(jié)果展示

1、第四章 偏微分程數(shù)值求解以及運(yùn)算結(jié)果展示§ 4.1 偏微分方程數(shù)值求解和水平集方法用于處理圖像問題的偏微分方程比較復(fù)雜， 往往無法求出一個(gè)解析解， 通常是求數(shù)值 逼近解。因此，理論上說，偏微分方程數(shù)值求解的一般方法都可以拿到圖像偏微分方程的 求解中來使用。從原則上說，由于圖像偏微分方程是利用微分幾何的方法來描述一幅圖像，以及一個(gè) 圖像處理問題。而圖像信息主要蘊(yùn)含在像素灰度級的跳變(圖像邊緣信息) ，以及表變得 幾何分布上。因此， 圖像偏微分方程往往具有其復(fù)雜性。而圖像偏微分方程的求解，往往 需要復(fù)雜的技巧性。常用的偏微分方程數(shù)值求解方法有：有限差分法，有限單元法，邊界元法等等，其中，

2、對 圖像偏微分方程求解主要使用有限差分法， 迎風(fēng)格式的有限差分法在圖像偏微分方程求解 中廣泛地使用。目前在偏微分圖像處理中，通常輔以使用Osher和Sethian提出的水平集 方法( Level Set Method )進(jìn)行數(shù)值求解。§ 4.1.1 圖像偏微分方程數(shù)值求解水平集方法( Level Set Method )最初由 Osher 和 Sethian 提出 ，其主要思想是將 曲線、曲面和圖像演化表示為更高維數(shù)的超平面水平集， 其演化速度是該曲線或曲面的局 部曲率的函數(shù)。 因此，水平集方法可以看作是使用歐拉方法求解隱性偏微分方程的一種具 體實(shí)現(xiàn)方式， 具體地說， 以計(jì)算的復(fù)雜度

3、換取對拓?fù)渥兓倪m應(yīng)性， 能處理比較困難的曲 線或曲面演化過程拓?fù)浣Y(jié)構(gòu)變化問題。水平集法一般需要有初始的輪廓線或輪廓面和速度圖像。 其中初始輪廓線或輪廓面一 般由用戶手工勾劃出或輸入種子點(diǎn)代替。 速度圖像可以根據(jù)圖像的特征如曲率、 平均曲率 或梯度而定。 水平集演化產(chǎn)生的圖像是一個(gè)時(shí)間的函數(shù)圖像,稱為時(shí)間距離圖 (Time-Distance Plane )。即對每個(gè)象素 , 輪廓前沿到達(dá)此象素所經(jīng)歷的最短時(shí)間。通過 指定一個(gè)停止閾值 ,當(dāng)前沿進(jìn)化迭代到一定數(shù)目時(shí)停止 ,這就產(chǎn)生了以演化時(shí)間為參數(shù)的 等值面序列 , 直至完成圖像處理任務(wù)。初始化水平集方法時(shí)， 需要區(qū)分兩維閉合曲線的內(nèi)外部， 以構(gòu)

4、造符號距離函數(shù) (SignedDistanee Function ),而如何區(qū)分任意形狀閉合曲線(閉合曲面)的內(nèi)外部是一個(gè)比較煩 瑣的問題。另外，水平集方法迭代過程中，為了計(jì)算穩(wěn)定，往往要間隔地對水平集函數(shù)重新初始 化(Re-initialization)，而常規(guī)的直接計(jì)算點(diǎn)與曲線距離的初始化方法5計(jì)算時(shí)間較長，大大增加了計(jì)算量。因此，如何有效而快速地進(jìn)行水平集函數(shù)的重新初始化，對提高 水平集方法的計(jì)算效率很有意義。水平集方法進(jìn)行數(shù)值求解，在偏微分方程圖像處理領(lǐng)域已經(jīng)得到了廣泛的應(yīng)用。比如，Alvarez9 以及Malladi , Sethian5 等用水平集去除圖像噪聲；Parogios以水

5、平集方 法進(jìn)行紋理分割以及運(yùn)動目標(biāo)分割7,8 ； Bertalmio 等將水平集方法應(yīng)用于圖像變形和 破損圖像修復(fù)(Inpainting )中10,11 ； Mansouri將水平集運(yùn)用于運(yùn)動目標(biāo)跟蹤?quán)徲?12 ； Fuaeras和Keriven 13以水平集方法解決立體匹配問題。在計(jì)算過程中，所有網(wǎng)格點(diǎn)上的函數(shù)值都有更新的計(jì)算。因此計(jì)算量往往是很龐大的， 尤其是對于計(jì)算用于圖像處理的偏微分方程時(shí)，更是如此。復(fù)雜的計(jì)算方法所換來的優(yōu)勢是，曲線的拓?fù)浣Y(jié)構(gòu)變化處于整個(gè)算法的控制之下。因?yàn)椴徽撉€如何變化，它始終處于平面或空間網(wǎng)格中，而演化后的輪廓跟蹤步驟可以得到平面上任何位置的曲線。因此，水平集方

6、法求解也大大促進(jìn)了曲率方程曲線演化模型的迅速發(fā)展66。§ 水平集方法圖像微分方程作的求解，通常使用的必要工具是圖像的水平線 (Level Set Con tours)， 即，圖像的水平集各自對應(yīng)的邊界線。在這個(gè)框架下，一幅圖像就由一系列的零水平集函 數(shù)表示。圖像處理的過程就表現(xiàn)為相應(yīng)的圖像的水平線的移動/演化過程。因此，在流體力學(xué)中表達(dá)各種曲率流的方程在圖像處理中得到了廣泛的應(yīng)用。水平集方法中，圖像偏微分方程可以用一個(gè)普遍的表達(dá)式來表示。一幅圖像可以由一 系列的零水平集函數(shù)表示。圖像處理的過程表現(xiàn)為水平集函數(shù)組所約束的曲線(稱為水平 線)的演化過程。我們可以把以上所述的圖像偏微分方程

7、用一個(gè)普遍的表達(dá)式來表示。設(shè)叮I :山2表示一幅灰度圖像，而叮°(x,y)是灰度值。引入迭代次數(shù)參量t (“偽 時(shí)間”參量)，圖像的變化可以用偏微分演化方程 "(x,y,t):山2 0, )是隨時(shí)間變化的圖像，而:：-0是初始條件，>是表示所使用的算子。偏微分方程的解叮（x,y,t）給出了偽時(shí)刻t時(shí)（即迭代次數(shù)為t）圖像的狀態(tài)。對于矢量值圖像，可以建立相應(yīng)的偏微 分方程組。在這個(gè)框架下表示曲線（輪廓線，水平線）的演化方程。對于曲線 G是山2:的函 數(shù)。曲線的演化視為在法線方向上的變形， 而變形速度與相應(yīng)曲線或曲面的曲率有關(guān)， 方 程如下 C i）N其中： i是曲率，N

8、是曲線G的法矢量，J是曲線演化算子。有.t時(shí)，方程中可能也同時(shí)出現(xiàn)切向速度，但這一速度分量不影響曲線幾何形狀的變化。這一方法可以方便地推廣到曲面演化模型的情況下。水平集方法（Level Set Method）最初由Osher和Sethian提出28，其主要思想是將曲線、曲面和圖像演化表示為更高維數(shù)的超平面水平集，其演化速度是該曲線或曲面的 局部曲率的函數(shù)。因此，水平集方法可以看作是使用歐拉方法求解隱性偏微分方程的一種 具體實(shí)現(xiàn)方式，具體地說，以計(jì)算的復(fù)雜度換取對拓?fù)渥兓倪m應(yīng)性，能處理比較困難的曲線或曲面演化過程拓?fù)浣Y(jié)構(gòu)變化問題。 在計(jì)算過程中，所有網(wǎng)格點(diǎn)上的函數(shù)值都有更新 的計(jì)算。復(fù)雜的計(jì)算

9、方法所換來的優(yōu)勢是，曲線的拓?fù)浣Y(jié)構(gòu)變化處于整個(gè)算法的控制之下。 因?yàn)椴徽撉€如何變化，它始終處于平面或空間網(wǎng)格中，而演化后的輪廓跟蹤步驟可以得 到平面上任何位置的曲線。Osher-Sethain提出的水平集方法有效地解決了在活動圍線演化過程中曲線拓?fù)浣Y(jié)構(gòu) 變化的問題。因此，水平集的方法在包括圖像處理在內(nèi)的諸多領(lǐng)域內(nèi)，迅速地得到了廣泛應(yīng)用。但是傳統(tǒng)的水平集方法對圖像只能作兩個(gè)分片區(qū)域的分割，而且在數(shù)值計(jì)算中常常需要引入對零水平集函數(shù)的重新初始化過程，造成計(jì)算的繁瑣。針對這些弊端，適合圖像 分割的新的水平集方法不斷被提出。典型的技術(shù)是Vese-Chan多相水平集方法和P-D測地 活動區(qū)域組，其他的

10、水平集方法還有1,4-6,23,40。§ 4.2運(yùn)算結(jié)果展示以及對于水平集方法求解技巧的討論§ 運(yùn)算結(jié)果展示本論文工作主要跟蹤了偏微分方程圖像處理方面的一些感興趣的問題，尤其對偏微分方程圖像分割和圖像修補(bǔ)問題的研究方法和研究動態(tài)有了較為清晰的認(rèn)識。本工作主要完成了兩個(gè)方面的圖像修補(bǔ)工作。本工作的第一個(gè)方面是利用 R-O-F圖像復(fù)原模型和BSBC圖像修補(bǔ)模型，在圖像子空 間中進(jìn)行修補(bǔ)的算法5, 48。工作包括如下三個(gè)步驟：改進(jìn)R-O-F圖像復(fù)原模型，從而使得原模型中被摒除的部分v成為輸入圖像I的紋理子圖像?；诟倪M(jìn)了的R-O-F圖像復(fù)原模型進(jìn)行圖像分解，得到輸入圖像I的結(jié)構(gòu)子

11、圖像u 和紋理子圖像v。在圖像分解子空間中進(jìn)行圖像修補(bǔ)。其中對結(jié)構(gòu)子圖像u做和紋理子圖像v。本實(shí)驗(yàn)中所用的測試圖像有兩個(gè)。一個(gè)是使用 MATLAB成的簡單的幾何測試圖像， 圖像中也包含了簡單的紋理信息。另一幅是自然圖像白塔。所使用的彩色圖像首先被分解 為R,G,B,子圖像，然后分別通過本修補(bǔ)程序處理；最后重新合成為彩色圖像。幾何測試圖像的測試圖像結(jié)果如下。圖4-1.1原始圖像，缺損部分正在被填補(bǔ)inpain ted紅色圓形有缺損fig.4-1.2 in termediateimagefig.4-1.1 in put imageresult圖4-1.3輸出圖像，紅色圖4-1.2中間結(jié)果，圖像圓形

12、缺損部分基本被填補(bǔ)自然圖像白塔的測試圖像結(jié)果如下。圖4-2.1原始圖像，白塔將被剔除fig.4-2.1 original image, with tower to be removedfig.4-2. 2b inpainting maskfig.4-2. 2a tailored image to reduce calculati on圖4-2. 2b標(biāo)記圖像填補(bǔ)區(qū)域圖4-2.4 a,b,c R-,G-,B-子空間圖像的中間結(jié)果fig.4-2.4a,b,c in termediate result of R-,G-,B-subspace of tailored imagefig.4-2. 5b

13、result image of 7圖4-2. 5a輸出圖像 圖4-2. 5b輸出圖像7 fig.4-2. 5a result image本工作的第二個(gè)方面是再現(xiàn)了 EBI最佳相似復(fù)制法24。這是一個(gè)非基于偏微分方程的圖像修補(bǔ)方法，其主要算法包括如下三個(gè)步驟:對待填補(bǔ)區(qū)域邊界上的像素點(diǎn) P進(jìn)行優(yōu)先級排序。其中，鄰近的已知像素信息越豐富的點(diǎn)更具有優(yōu)先權(quán)在圖像的已知區(qū)域內(nèi)尋找與待填補(bǔ)區(qū)域邊界上的像素點(diǎn) P相對應(yīng)的最佳相似點(diǎn)Q, 并將相應(yīng)的像素信息向待填補(bǔ)區(qū)域內(nèi)部推進(jìn)。依次反復(fù)推進(jìn)，直至修補(bǔ)完畢。過程中，有必要進(jìn)行人工干預(yù)以保證最佳相似點(diǎn) 的選取和像素點(diǎn)優(yōu)先級排序沒有出現(xiàn)明顯的差錯(cuò)。實(shí)驗(yàn)中的測試圖像結(jié)

14、果如下。圖4-3.1原始圖像，圖像前景人物將被剔除fig.4-3.1in put image ， woma n to beremoved圖4-3.2輸出圖像，人物用圖像背景成功地填補(bǔ)fig.4-3.2 image inpain ted在附錄中，列出了我用MATLA語言實(shí)現(xiàn)的兩個(gè)圖像修補(bǔ)算法的源程序片段： BSBC在圖像子空間中修補(bǔ) ，和一個(gè)非基于偏微分方程的圖像修補(bǔ)方法，EBI最佳相似復(fù)制法 24 。§ 4.2.2 對于水平集方法求解技巧的討論我的工作主要是試圖理解偏微分方程求解的基本方法， 并通過跟蹤現(xiàn)有的方 法和實(shí)現(xiàn)其提出的水平集方法輔助求解的基本策略， 試圖理解水平集方法在偏微

15、 分方程圖像處理中的使用技巧。附錄中的 BSBC偏微分方程圖像處理程序，是利 用平均曲率移動 ( M ean Curvature Motion) ，或稱為平均曲率流( MeanCurvature Flow)的模型，使用的是迎風(fēng)格式(Upwind Scheme)的有限差分法。平均曲率 移動是一個(gè)數(shù)學(xué)性質(zhì)很好的偏微分方程， 因而，比較容易用簡單的有限差分法實(shí) 現(xiàn)。重新初始化( Re-initialization )和粘性延伸( Velocity Extension )， 是現(xiàn)有水平集( Level Set )方法中兩個(gè)重要的技巧。然而，這些都是為了保證 曲線演化( Evolution )的穩(wěn)定和可

16、用的結(jié)果，人們不得不采取的補(bǔ)救措施 (Remedy。當(dāng)然，絕大多數(shù)水平集的文章都不說這是缺陷，相反，這些補(bǔ)救措 施的實(shí)現(xiàn)方法都成為一個(gè)重要的課題。 比如說，單單關(guān)于如何做重新初始化這方 面的文章就有很多。所謂重新初始化( Re-initialization ，顧名思義就是反復(fù)初始化。在迭代 之前就要把水平集函數(shù)( Level Set Function )初始化為一個(gè)數(shù)學(xué)性質(zhì)很好的函 數(shù)。通常是利用符號距離函數(shù)( SDF, Signed Distance Function )，然后每迭代 幾步就要停下來， 修改這個(gè)函數(shù)， 把它改成一個(gè)新的符號距離函數(shù)； 再接著迭代， 周而復(fù)始 . 這個(gè)每迭代幾步

17、就停下來修改水平集函數(shù)的過程就是重新初始化。在修改水平集函數(shù)的同時(shí)， 還往往要修改或重建速度函數(shù) (SpeedFunction ) (就是所求解的水平集函數(shù)表示下的偏微分方程的系數(shù)) 。這個(gè)過程叫粘性延伸 ( Velocity Extension )。每做一次重新初始化， 就需要重新計(jì)算符號距離函數(shù) ( SDF, Signed Distance Function )，以確定各個(gè)像素點(diǎn)演化的優(yōu)先順序。符號距離函數(shù)，其符號是指圖 像中的一個(gè)像素點(diǎn)與零水平曲線的相對位置， 由 Heaviside 函數(shù)確定； 其大小是 指該像素點(diǎn)到一個(gè)閉合曲線零水平曲線上的最短距離，這個(gè)問題比較麻煩 ( 盡管 有些文

18、章也提到了一些算法 ) 。如果是指給定封閉曲線，算相應(yīng)的符號距離函數(shù)，還是很直接的?？墒窃趯?水平集函數(shù)進(jìn)行重新初始化的時(shí)候， 即每迭代幾步， 停止迭代， 計(jì)算所得到的水 平集函數(shù)的零水平集曲線( Zero Level Contours )的符號距離函數(shù)，這就不容 易了?，F(xiàn)在主要有三種方法：1。 確定水平集函數(shù)(以下記為u)的零水平集曲線(這可以用Matlab函數(shù) contour 得到)，然后算這些零水平集曲線的符號距離函數(shù)。 這種方法的思路 很直接，但是在實(shí)際應(yīng)用中，如果 u的零水平集函數(shù)(Zero Level Set )由 幾個(gè)零水平集曲線構(gòu)成， 而且出現(xiàn)大零水平集曲線里面還有小零水平集曲線 的情況，這時(shí)候在判斷哪些點(diǎn)是這些零水平集曲線的內(nèi)部外部(即確定 SDF 的符號)會相當(dāng)麻煩。2。用偏微分方程( PDE， Partial Differential Equation )的方法，即解 一個(gè)基于時(shí)間(time dependent)的偏微分方程，其平穩(wěn)穩(wěn)定解(