陳宇， 陳德運

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

0 引言

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

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

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

重建領域提供了一種有效可行的新方法。

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

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

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

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

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

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

2 算法原理

非線性算子方程為

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

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

式中，δ表示誤差水平。

當誤差δ≠0時，通過求解初值問題，即

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

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

為了闡述方便，記

因為F(x*)=y，φ(x*)=x*，因此式(7)可以表示為

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

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

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

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

式(11)可以改寫為

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

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

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

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

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

下對本文所提出算法的收斂性進行分析:設

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

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

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

又因為

由式(26)有

由式(28)有

按廣義偏差原則，有

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

進一步有

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

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

由Morozov偏差原則

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

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

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

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

3 仿真與實驗結果

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

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

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

成像完成后需要分析重建圖像的質量。實驗選用空間圖像誤差作為評判圖像質量的標準，定義為

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

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

圖2 成像結果的比較Fig.2 Comparison of reconstructed tomograms

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

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

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

4 結語

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

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

［1］郭志恒，律德才，邵富群.基于差分靈敏度模型的電容層析成像圖像重建方法［J］.中國電機工程學報，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電容層析成像結構參數(shù)分析與優(yōu)化［J］.哈爾濱理工大學學報，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］.哈爾濱理工大學學報，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è)自動化，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.