基于改進(jìn)Runge-Kutta型landweber的電容層析成像圖像重建算法

陳宇， 陳德運(yùn)

(1.東北林業(yè)大學(xué)信息與計(jì)算機(jī)工程學(xué)院，黑龍江哈爾濱 150040;2.哈爾濱理工大學(xué)計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，黑龍江哈爾濱 150080)

0 引言

隨著20世紀(jì)80年代中期過程層析成像技術(shù)(process tomography，PT)的興起和發(fā)展，該技術(shù)之一的電容層析成像技術(shù)(electrical capacitance tomography，ECT)成為眾多科研工作者研究的對(duì)象［1－2］。ECT技術(shù)以兩相流或者多相流為主要研究對(duì)象，可實(shí)現(xiàn)過程參數(shù)在線實(shí)時(shí)監(jiān)測的功能，已廣泛用于電力、能源、醫(yī)藥、石油等眾多多相流檢測領(lǐng)域［3］。因ECT技術(shù)可靠性高、成本低、結(jié)構(gòu)簡單、非入侵式、使用范圍廣、安全性強(qiáng)等優(yōu)點(diǎn)，國內(nèi)外研究者對(duì)ECT圖像重建進(jìn)行了大量研究并取得了較好的成果，成為當(dāng)今圖像重建領(lǐng)域研究發(fā)展的主流［4－5］。由于ECT技術(shù)是通過有限個(gè)數(shù)的電容值來重建圖像，其“軟場”效應(yīng)和非線性的特點(diǎn)，使ECT系統(tǒng)的解極易隨多項(xiàng)流體的變化而改變，穩(wěn)定性差，很難用具體的數(shù)學(xué)解析式描述，致使圖像重建難度增大［6－7］。因此深入研究探索更好的圖像重建算法是當(dāng)務(wù)之急。

國內(nèi)外學(xué)者經(jīng)過長期的研究，提出了很多ECT圖像重建的算法，目前較常使用的方法有:線性反投影算法LBP［8］、landweber迭代算法［9］、基于模型的MOR算法、人工神經(jīng)網(wǎng)絡(luò)算法、共軛梯度算法［10］等。LBP算法是最早使用的簡單成像法，運(yùn)算速度快，但重建圖像失真嚴(yán)重，效果較差。landweber迭代算法是在landweber迭代法的基礎(chǔ)上求解欠定方程的圖像重建算法，是目前應(yīng)用最為廣泛的算法，極大地提高了成像速度，但介質(zhì)分界面仍存在明顯的平滑效應(yīng)?；谀Ｐ偷腗OR算法比LBP算法精確率高，但過程中需設(shè)定參數(shù)來描述介電常數(shù)分布情況，參數(shù)越多圖像重建時(shí)間越長，因此消耗時(shí)間大是該方法的最大缺點(diǎn)［11］。神經(jīng)網(wǎng)絡(luò)法是通過建立電容值和圖像像素灰度值之間的映射關(guān)系來實(shí)現(xiàn)圖像重建，質(zhì)量較高、響應(yīng)速度快，但對(duì)于訓(xùn)練樣本的完整性要求較高，且訓(xùn)練過程復(fù)雜、耗時(shí)大，在實(shí)際應(yīng)用中難度較大。共軛梯度法(CG)收斂速度較快、穩(wěn)定性高，但只適合系數(shù)矩陣為對(duì)稱正定的情況，對(duì)復(fù)雜流型的圖像重建效果不理想。各種方法都有自身的優(yōu)點(diǎn)但也存在各自的局限性，還需不斷研究加以完善。

本文在電容層析成像算法的基礎(chǔ)上，提出一種改進(jìn)R－K型landweber算法。該算法可以增強(qiáng)圖像重建的穩(wěn)定性，提高重建圖像的質(zhì)量。實(shí)驗(yàn)中采用多種算法對(duì)本文考察對(duì)象進(jìn)行圖像重建，結(jié)果表明本文提出的算法獲得的重建圖像的質(zhì)量高于landweber、LBP、共軛梯度算法和最速下降法，為ECT技術(shù)圖像

重建領(lǐng)域提供了一種有效可行的新方法。

1 電容層析成像系統(tǒng)基本原理

電容層析成像ECT系統(tǒng)由3部分組成:電容傳感器、數(shù)據(jù)采集系統(tǒng)、成像計(jì)算機(jī)，如圖1所示。其基本原理為由于不同物質(zhì)的介電常數(shù)不同，介質(zhì)分布濃度的變化會(huì)引起混合物的等價(jià)介電常數(shù)變化，致使測量的電容值隨之改變。由此可依據(jù)實(shí)際測量的電容值來重建管道內(nèi)介質(zhì)的分布情況。

圖112 電極電容層析成像系統(tǒng)的組成Fig.1 Composition chart of the 12-electrode ECT system

測量過程如下:以12個(gè)極板中的任意極板為起點(diǎn)，逆時(shí)針標(biāo)記12個(gè)極板編號(hào)。選取電極板1為源極板(即公共電極)，剩余的2，3，…，12電極板為檢測電極板。對(duì)源極板施加大小為U的固定電壓，測量電極板對(duì)1－2，1－3，…，1－12之間的電容值，每次測量時(shí)閑置電極均需接地。上述操作后選取電極板2作為公共電極，其余為檢測電極，依照之前的方法測量電極對(duì)2－3，2－4，…，2－12的電容值。依此類推，直到測量得到電極對(duì)11－12的電容值，完成整個(gè)測量過程。最終總共得到66個(gè)獨(dú)立的電極對(duì)電容值。

目前，多數(shù)ECT成像算法是在介電常數(shù)到電容映射的線性模型基礎(chǔ)上，經(jīng)離散化、線性化和歸一化建立的模型為

式中:C∈Rm×1是歸一化電容向量;S∈Rm×n為系數(shù)矩陣(又稱靈敏度矩陣);G∈Rn為歸一化介質(zhì)分布圖像向量。ECT圖像重建的關(guān)鍵在于根據(jù)給定的電容值C來求解介電常數(shù)的分布情況G。

2 算法原理

非線性算子方程為

其中:F是Hilbert空間上的一個(gè)非線性算子;H是Hilbert空間，具有相應(yīng)的內(nèi)積和范數(shù)［12］。

對(duì)于式(3)右邊的觀測數(shù)據(jù)不能得到準(zhǔn)確值，只能得到具有一定誤差擾動(dòng)的yδ，且有

式中，δ表示誤差水平。

當(dāng)誤差δ≠0時(shí)，通過求解初值問題，即

式中，0＜t＜T，xδ(0)=x0，可以獲得解x*的正則化估計(jì)值xδ(T)，T代表正則化參數(shù)。用R級(jí)龍格－庫塔(Runge－Kutta)方法求解式(5)，采用的離散化方程為

稱式(6)為Runge－Kutta型landweber方法。當(dāng)誤差水平δ=0時(shí)，有

為了闡述方便，記

因?yàn)镕(x*)=y，φ(x*)=x*，因此式(7)可以表示為

當(dāng)誤差水平δ≠0時(shí)，式(8)和式(9)可分別表述為

由上述理論可知，ECT重建算法中的電容向量C和圖像分布向量G之間存在非線性關(guān)系為

將計(jì)算值和電容測量值的誤差范數(shù)平方作為目標(biāo)函數(shù)，則有

用靈敏度矩陣S取代F'(G)，以SGk取代F(G)，所以式(10)可以表示為

式(11)可以改寫為

對(duì)于ECT系統(tǒng)，式(17)的迭代公式具有很高的收斂性，但需要的步長較多，屬于小步長搜索。

考慮到ECT系統(tǒng)的傳統(tǒng)landweber迭代公式為

式(18)的右邊的梯度用hk表示為

由式(18)和式(19)，則知圖像向量G的第k+1次迭代格式為

可知landweber算法的每次迭代的搜索方向就是負(fù)梯度方向，由于該步長較大，可以很快到達(dá)收斂區(qū)域。但到達(dá)收斂區(qū)域后，在向真正的解逼近時(shí)，會(huì)由于算法自身的缺陷，而只能得到一個(gè)解的近似。且當(dāng)跳過收斂域后，迭代解和真解之間的距離會(huì)變大。故考慮每次迭代時(shí)先用landweber算法進(jìn)行搜索，然后再用Runge－Kutta型landweber算法的迭代步進(jìn)行精細(xì)搜索，以得到更準(zhǔn)確的解，同時(shí)提高了收斂速度。算法的最終迭代步可以寫成

下對(duì)本文所提出算法的收斂性進(jìn)行分析:設(shè)

對(duì)于Runge－Kutta迭代公式(17)，有

假設(shè)C=SG的解G*的鄰域內(nèi)，算子滿足非線性條件

式中，0＜η＜0.5，并設(shè)δ代表誤差水平，L代表范數(shù)上界，則有

又因?yàn)?/p>

由式(26)有

由式(28)有

按廣義偏差原則，有

式中，τ是依賴于η的正數(shù)，且τ滿足

進(jìn)一步有

由式(33)和式(34)，可以得出式(32)不等式右邊為負(fù)數(shù)，說明Gk+1比Gk更單調(diào)收斂接近于真解。

對(duì)于landweber的迭代公式(18)，有

由Morozov偏差原則

則可知式(35)的不等式右邊為負(fù)數(shù)，說明式(18)中的Gk+1比Gk更單調(diào)收斂接近于真解。綜上可知式(21)迭代公式是單調(diào)的，且逐漸逼近真解。

每次迭代過程都要對(duì)圖像從物理意義的角度進(jìn)行修正，即在每次迭代循環(huán)中引入Gk的估計(jì)值應(yīng)在0和1.0之間的先驗(yàn)信息。為此，式(21)修改成投影迭代為

式中，P+是非負(fù)凸集上的投影，且有

其中，G(k)代表列向量G第k次迭代的第i個(gè)分量。為保證每次迭代得到的解都為有限非負(fù)值，投影算子需確保每次迭代都收斂于一個(gè)凸集。實(shí)驗(yàn)表明，通過施加投影算子來引入了物理意義上的約束的方法，在一定程度上可以加快重建速度、提高重建圖像的質(zhì)量。

3 仿真與實(shí)驗(yàn)結(jié)果

本研究采用12電極系統(tǒng)應(yīng)用Matlab軟件進(jìn)行仿真實(shí)驗(yàn)，驗(yàn)證該算法性能。實(shí)驗(yàn)中，用32×32的網(wǎng)格將管道截面分為1024個(gè)像素，其中856個(gè)成像單元為該界面的有效區(qū)域。實(shí)驗(yàn)時(shí)需先為極低位層流、柱狀流、小半徑核心流和低位層流進(jìn)行預(yù)設(shè)置，再采用本文提出的(improved Runge－Kutta type landweber，簡稱IRKL)算法進(jìn)行圖像重建。同時(shí)采用線性反投影法(LBP)、共軛梯度法(CG)、最速下降法(SD)和landweber迭代法重建圖像，并與IRKL算法的重建圖像質(zhì)量進(jìn)行比較，所有成像圖像采用閾值進(jìn)行了濾波處理。

圖像重建的速度用迭代次數(shù)N表示，N越大表示重建時(shí)間越長，重建速度越慢。由于LBP算法為單步處理，故N=0。迭代次數(shù)N的選取需由數(shù)值實(shí)驗(yàn)確定。典型方法為迭代誤差滿足

就停止迭代。大量實(shí)驗(yàn)表明，本文IRKL算法的迭代次數(shù)N通常選取十幾(簡單模型)到幾十(復(fù)雜模型)之間就可得到較好的成像效果。

成像完成后需要分析重建圖像的質(zhì)量。實(shí)驗(yàn)選用空間圖像誤差作為評(píng)判圖像質(zhì)量的標(biāo)準(zhǔn)，定義為

其中:gimg為重建圖像向量;ginit為介質(zhì)分布原型圖像向量;i為成像區(qū)域剖分單元索引;n為成像區(qū)域單元總數(shù)。

成像結(jié)果的比較如圖2所示(黑色區(qū)域是水，白色區(qū)域是變壓器油)。

圖2 成像結(jié)果的比較Fig.2 Comparison of reconstructed tomograms

從成像結(jié)果可知，LBP算法的重建圖像大體接近原流型，而landweber、SD算法和本文的IRKL算法的重建圖像則十分接近原流型。IRKL算法對(duì)于極底位層流的成像效果相對(duì)要好，最接近原始流型，且明顯減緩了邊界模糊效應(yīng)。表1為重建圖像誤差比較，結(jié)合圖2結(jié)果可知，對(duì)于低位層流和柱狀流，LBP誤差最大;對(duì)于極低位層流，CG法誤差最大;對(duì)于小半徑核心流，SD算法誤差最大;而對(duì)于這四種流型，IRKL算法產(chǎn)生的誤差最小，且電容值混入了5%的高斯白噪聲后，該算法仍然能較好的進(jìn)行圖像重建，表現(xiàn)出了較強(qiáng)的抗噪聲能力。表2表明，IRKL算法同landweber法和SD算法的迭代步數(shù)相近，而同CG迭代步數(shù)相差較大。從以上分析可以看出，使用改進(jìn)Runge－Kutta型landweber的圖像重建算法，其成像的精度和質(zhì)量比LBP、landweber、CG和SD算法要好。

表1 圖像誤差Table 1 Image error/%

表2 迭代步數(shù)Table 2 Number of iteration

4 結(jié)語

本文提出了一種改進(jìn)的Runge－Kutta型landweber的正則化迭代電容層析成像算法(IRKL)，在分析ECT反問題病態(tài)性的基礎(chǔ)上，采用迭代求解時(shí)廣搜和精搜同時(shí)進(jìn)行的策略，并利用范數(shù)殘量誤差的性質(zhì)分析了改進(jìn)的Runge－Kutta型landweber算法的收斂性。該算法程序簡單、所需存儲(chǔ)量小、成像精度高、抗噪聲能力強(qiáng)且易于滿足收斂條件。數(shù)值實(shí)驗(yàn)結(jié)果表明該算法的重建圖像更接近原流型，且質(zhì)量遠(yuǎn)遠(yuǎn)優(yōu)于傳統(tǒng)的LBP、SD、landweber和CG算法，為ECT圖像重建技術(shù)提供了一種有效可行的新方法。

未來將把研究重心放在復(fù)雜流型的自動(dòng)辯識(shí)領(lǐng)域，在迭代的過程中加入補(bǔ)償算法以此減少測量誤差，進(jìn)一步提高圖像重建質(zhì)量，使ECT圖像重建技術(shù)更廣泛地應(yīng)用于生產(chǎn)實(shí)踐。

［1］郭志恒，律德才，邵富群.基于差分靈敏度模型的電容層析成像圖像重建方法［J］.中國電機(jī)工程學(xué)報(bào)，2012，32(23):75－82.

GUO Zhiheng，LU Decai，SHAO Fuquan.Image reconstruction method for electrical capacitance tomography based on the different sensitivity model［J］.Proceedings of the CSEE，2012，32(23):75－82.

［2］RIMPILAINEN V，HEIKKINEN L M，VAUHKONEN M.Moisture distribution and hydrodynamics of wet granules during fluidized－bed drying characterized with volumetric electrical capacitance tomography［J］.Chemical Engineering Science，2012，75(18):220－234.

［3］李巖，朱艷丹，袁小花，等.基于ANSYS電容層析成像結(jié)構(gòu)參數(shù)分析與優(yōu)化［J］.哈爾濱理工大學(xué)學(xué)報(bào)，2012，17(1):54－58.

LI Yan，ZHU Yandan，Yuan Xiaohua，et al.Analysis and optimization of structural parameters based on ANSYS in ECT［J］.Journal of Harbin University of Science and Technology，2012，17(1):54－58.

［4］周云龍，衣得武，高云鵬.基于混合模型電容層析成像系統(tǒng)敏感矩陣建立［J］.哈爾濱理工大學(xué)學(xué)報(bào)，2012，17(2):40－44.

ZHOU Yunlong，YI Dewu，GAO Yunpeng.Sensitivity matrix construction in electrical capacitance tomography system based on mixed models［J］.Journal of Harbin University of Science and Technology，2012，17(2):40－44.

［5］OSTROWSKI K L，LUKE S P，WILLIAMS R A.Simulation of the performance of electrical capacitance tomography for muasurement of dense phase pneumatic conveying［J］.Chemical Engineering Journal，1997，68(2/3):197－205.

［6］HJERTAKER B T.Static characterization of a dual sensor flow imaging system［J］.Flow Measurement and Instrumentation，1998，9(3):183－191.

［7］LIU S，YANG W Q，WANG H，et al.Investigation of square fluidized beds using capacitance tomography:preliminary results［J］.Measurement Science and Technology，2001，12(8):1120－1125.

［8］XI Shaolin，ZHAO Fengzhi.Computation Methods for Optimization［M］.Shanghai:Shanghai Scientifie＆Technical Publishers.1983:113－117.

［9］YANG W Q，PENG L H.Image reconstruction algorithms for electrical capacitance tomography［J］.Measurement Science and Technology，2003，14(1):R1－R3.

［10］周云龍，衣得武，高云鵬.基于INGA的電容層析成像圖像重建算法的研究［J］.制造業(yè)自動(dòng)化，2012，34(11):81－86.

ZHOU Yunlong，YI Dewu，GAO Yunpeng.Imaging reconstruction algorithm based on improved niche genetic algorithms for electrical capacitance tomography system［J］.Manufacturing Automation，2012，34(11):81－86.

［11］LIU S，F(xiàn)U L，YANG W Q，et al.Prior-online iteration for image reconstruction with electrical capacitance tomography［J］.IEE Proceedings Science，Measurement and Technology，2004，151(3):195－200.

［12］HANKE M，NEUBAUER A，SCHERZER O.A convergence analysis of the landweber iteration for nonlinear ill-posed problems［J］.Numerische Mathematik，1995，72(1):21－37.