陳碧云 張業(yè)榮 潘鑫
(1.南京郵電大學(xué)電子科學(xué)與工程學(xué)院,南京 210003;2.空軍航空大學(xué),阜新 123000)
?
ADI-FDTD方法在乳腺癌檢測正向計算中的應(yīng)用
陳碧云1張業(yè)榮1潘鑫2
(1.南京郵電大學(xué)電子科學(xué)與工程學(xué)院,南京 210003;2.空軍航空大學(xué),阜新 123000)
傳統(tǒng)的時域有限差分(Finite-Difference Time-Domain, FDTD)算法受到穩(wěn)定性條件的制約,時間步長受限于空間網(wǎng)格的尺寸.醫(yī)學(xué)應(yīng)用講究即時性,為提高成像的速度, 文中采用無條件穩(wěn)定的交替隱式時域有限差分(Alternating-Direction Implicit Finite-Difference Time-Domain, ADI-FDTD)算法替代傳統(tǒng)的FDTD算法進(jìn)行正向計算,通過實(shí)驗(yàn)得出采用 ADI-FDTD算法在保證精度的前提下,計算時間可縮短為FDTD算法的四分之一,為乳腺癌微波即時成像提供了可能.
交替隱式時域有限差分法;微波斷層成像;乳腺癌
DOI 10.13443/j.cjors.2016013101
微波成像技術(shù)作為一種新型的乳腺癌檢測手段,因低風(fēng)險、高靈敏度和高對比度等一系列優(yōu)點(diǎn)而日益受到重視[1-2].相比于乳腺鉬靶X線檢查,微波的能量僅為幾個電子伏特,消除了對人體的健康隱患,同時非侵入式的檢查手段不會給人帶來不適;相比于磁共振成像(Magnetic Resonance Imaging, MRI)和正電子發(fā)射計算機(jī)斷層顯像(Positron Emission Computed Tomography,PECT),微波成像的低成本為普查提供了可能;相比于多普勒超聲成像,微波對腫瘤尤其是惡性腫瘤具有更高的敏感度,檢出率更高[3].
國內(nèi)外研究乳腺癌微波成像主要集中在兩種方法:一種是共焦成像(Confocal Microwave Imaging, CMI)[4-7];另一種是斷層成像(Microwave Tomography, MWT)[8-13].MWT類似于現(xiàn)代醫(yī)用中的CT,相比于CMI,其物理解釋更清晰,更具醫(yī)學(xué)診斷價值.MWT屬于電磁逆散射的研究范疇,處理逆散射問題通常采用最優(yōu)化迭代求解,大量的時間耗費(fèi)在重復(fù)的正向計算,因此在保證精度的前提下提高正向計算的效率應(yīng)該成為我們關(guān)注的焦點(diǎn).
目前MWT的研究主要集中在頻域范圍內(nèi),最早的臨床系統(tǒng)是由Meaney 等人[8]于2000年開發(fā)出的一套基于牛頓高斯迭代法微波斷層成像系統(tǒng),采用有限元法(Finite Element Method,FEM)算法和邊界元(Boundary Element Method,BEM)算法相結(jié)合進(jìn)行正向計算,然而由于頻帶的限制,獲得的分辨率不高,后來Rubk等[9]在Meaney研究的基礎(chǔ)上,將二維成像擴(kuò)展到了三維,采用矩量法(Method of Moments,MOM)進(jìn)行前向計算,然而三維成像時間過長,一次成像需要幾十個小時,為此Liu等[10]采用快速傅里葉變換(Fast Fourier Transform,FFT)算法還原了二維和三維參數(shù)分布.
針對乳腺癌成像,選擇時域方法更有優(yōu)勢:1)早期的乳腺腫瘤是毫米級的,利用寬帶高頻脈沖信號能夠提供較好的成像分辨率;2)時域MWT利用目標(biāo)的瞬態(tài)響應(yīng)波形進(jìn)行成像,包含的信息量更多,成像結(jié)果更準(zhǔn)確[14].Takenaka等人采用時域逆散射算法對二維乳腺癌模型進(jìn)行了圖像重建[11],劉廣東等人以LM算法為基礎(chǔ),分別對二維半圓形乳房[12]和三維半球乳房[13]模型進(jìn)行了時域斷層成像分析,這些算法的成像分辨率以及魯棒性都很好,然而臨床應(yīng)用普遍存在一大缺陷——成像時間過長,尤其是到了高維情形.
早期在時域上研究乳腺癌微波成像普遍采用了FDTD算法,作為經(jīng)典的電磁三大算法之一,FDTD被廣泛地運(yùn)用于電磁學(xué)的各個領(lǐng)域[15-16].然而傳統(tǒng)的FDTD算法要滿足柯朗-弗里德里希斯-列維(Courant-Friedrichs-Lewy,CFL)條件,時間步長受限于空間網(wǎng)格的大小,成像時想要得到比較好的分辨率,CPU的運(yùn)算時間將要大大增加.臨床應(yīng)用時,檢測的效率和速度不佳,因此我們需要尋找更加高效的算法來提高計算效率.近些年,ADI-FDTD算法引起了越來越多的關(guān)注[17-18],相對于傳統(tǒng)的FDTD方法,ADI-FDTD采用了交替隱式差分方法,時間步長不再受到CFL條件的制約,在保證精確度的前提下可以大大降低運(yùn)算成本.
本文首次將ADI-FDTD算法應(yīng)用到二維乳腺癌模型中,仿真實(shí)驗(yàn)采用二維半圓無限長柱狀乳房模型,分別對單個圓形,單個“十”字型和多個“十”字型的腫瘤模型進(jìn)行仿真測試,實(shí)驗(yàn)結(jié)果與傳統(tǒng)的FDTD方法進(jìn)行了比較,從而驗(yàn)證了ADI-FDTD算法的可行性與有效性.
乳房組織是非磁性介質(zhì),磁導(dǎo)率σm=0,二維TM波在乳房組織中傳播滿足如下麥克斯韋方程:
(1)
(2)
(3)
式中:r=(x,y)T表示位置矢量;ε(r)表示介電常數(shù);σ(r)為電導(dǎo)率;μ0為真空中的磁導(dǎo)率;Ez(r,t)表示電場在z軸上的分量;Hx(r,t)和Hy(r,t)分別為磁場強(qiáng)度在x軸和y軸上的分量.M根發(fā)射天線和N根接收天線分別位于r=rm(m=1,2,…,M)和r=rn(n=1,2,…,N)如圖1黑色實(shí)心點(diǎn)處.
微波斷層成像技術(shù)是利用電磁脈沖照射下產(chǎn)生的散射場,反演出目標(biāo)物體的電參數(shù)分布.通常采取的方案是將逆散射問題轉(zhuǎn)化為最優(yōu)化問題進(jìn)行迭代求解,基于最小二乘法構(gòu)建如下代價函數(shù):
(4)
Hy,m(r,t))T;
(5)
Hy,m(r,t))T.
(6)
1) 第一步n→n+1/2對方程中的x方向?qū)?shù)取隱式差分格式,對y方向?qū)?shù)取顯示,整理得:
(7)
(8)
(9)
第二步n+1/2→n+1,對方程中的y方向?qū)?shù)取隱式差分格式,對x方向?qū)?shù)取顯示,整理得:
(10)
(11)
(12)
求解類似于第一步,這里不再贅述.ADI-FDTD算法提供無條件穩(wěn)定,時間步長的選擇不再受限于空間網(wǎng)格的大小,計算時間就可以被大大的縮短.
建立二維半無限長圓柱體乳房模型如圖1,假設(shè)該模型浸沒在各向同性的耦合溶液中,該耦合溶液的電參數(shù)為εr=36,σ=0 S/m,與皮膚層的波阻抗匹配,乳房結(jié)構(gòu)簡化為一個半圓形的無限長圓柱體,最外層是厚度為2 mm的皮膚層,其電磁參數(shù)為εr=36,σ=1.0 S/m.皮膚層下面是乳房組織層,是一個半徑為50 mm的半圓柱體,電參數(shù)為εr=9.0,σ=0.15 S/m,連接乳房層的是一個20 mm厚的胸壁層,電參數(shù)為εr=50,σ=1.2 S/m.成像系統(tǒng)由距離皮膚表面8 mm放置的9根天線組成,9發(fā)8收(當(dāng)其中一根作為發(fā)射天線時,其余8根作為接收天線, 收發(fā)天線位置見圖 1中的黑色圓點(diǎn),作為發(fā)射天線時從左往右依次標(biāo)記為天線T1,天線T2,…,天線T9,同樣地接收天線從左往右依次標(biāo)記為天線R1,天線R2,…,天線R8).
圖1 乳房結(jié)構(gòu)模型圖
激勵源選取調(diào)制高斯脈沖,其時域表達(dá)式為
(13)
式中:fc為中心頻率,取2.0 GHz;τ為高斯脈沖的寬度,取1.436×10-9s;t0=1.288×10-9s.整個計算區(qū)域采用均勻網(wǎng)格剖分,由于Mur吸收條件推導(dǎo)比較簡單,在ADI-FDTD方法的運(yùn)用中較為普遍,因此此處截斷邊界采用一階Mur吸收邊界,傳統(tǒng)FDTD計算步長要滿足CFL條件:
(14)
這里取ΔtFDTD=Δx/(2c),ADI-FDTD 時間步長ΔtADI-FDTD=CFL*ΔtFDTD,這里取兩者相差的系數(shù)為CFL是為了突出ADI-FDTD與FDTD算法的區(qū)別, 這里CFL取整數(shù).
3.1 情形1
元胞尺寸取Δx=Δy=1 mm,位于點(diǎn)(52,20)有一個半徑大小為3 mm圓形腫瘤,接下來分別采用FDTD方法和ADI-FDTD方法進(jìn)行仿真實(shí)驗(yàn),FDTD時間步數(shù)Nmax取3 000,CFL取4,6,10,12,天線T1上放置激勵源,天線R1上接收的電磁場如圖2(a),圖2(b)和圖2(c)所示.
(a) Ez
(b) Hx
(c) Hy圖2 情形1(天線R1上接收到的場)
本文采用ADI-FDTD方法進(jìn)行前向計算,所涉及的信號是時域波形,相比頻域方法,時域波形由于包含了更多的信息,包括時間、強(qiáng)度以及相位,可能獲得更準(zhǔn)確的重建結(jié)果.從圖2(a),圖2(b)和圖2(c)仿真結(jié)果我們可以看出,ADI-FDTD的計算結(jié)果與FDTD基本一致,只是在后期出現(xiàn)了一些振蕩,產(chǎn)生這種振蕩的主要原因有二:一是采用的吸收邊界條件所引起的反射誤差;二是離散網(wǎng)格過程中階梯誤差,傳統(tǒng)的FDTD方法中也有相應(yīng)缺陷.同時也看到隨著CFL數(shù)值的增加,出現(xiàn)的振蕩越明顯,因?yàn)棣越大所帶來的數(shù)值色散誤差也越大,所以說CFL的選取還是有約束的.
表1中,T表示進(jìn)行一次正向計算的時間,TADI-FDTD/TFDTD表示采用ADI-FDTD方法與FDTD方法的時間比,從表1我們可以發(fā)現(xiàn)從CFL=4開始,采用ADI-FDTD的時間比FDTD短,當(dāng)CFL=12,采用ADI-FDTD計算時間僅為FDTD方法的24.9%,計算時間大大的縮短.
表1 情形1ADI-FDTD與FDTD計算時間比較
3.2 情形2
腫瘤大多都是不光滑的,為了更加接近真實(shí),采用“+”字形腫瘤模型,元胞的尺寸選取為Δx=Δy=2 mm,位于(90,30)存在一個由5個邊長為2 mm的方柱組成的“+”字型腫瘤.分別采用了FDTD方法和ADI-FDTD方法進(jìn)行實(shí)驗(yàn), FDTD時間步數(shù)Nmax取1 500,天線T1上放置激勵源,天線R1上接收的電磁場如圖3(a),圖3(b)和圖3(c)所示.
(a) Ez
(b) Hx
(c) Hy圖3 情形2(R1上接收到的場)
從圖3(a),圖3(b)和圖3(c)仿真結(jié)果我們可以看出,CFL取4,6,10,12時,FDTD和ADI-FDTD的計算結(jié)果還是比較接近的,整個起伏的走勢表現(xiàn)一致,除了后期出現(xiàn)少許振蕩.同時比較情形1和情形2,我們可以看出情形1取1 mm的誤差要小于情形2取2 mm時產(chǎn)生的誤差,網(wǎng)格剖分越細(xì),誤差越小,至于這種誤差對于成像精度的影響,在后期要結(jié)合成像結(jié)果進(jìn)行分析,本文暫不做討論.
從表2我們發(fā)現(xiàn),同樣的從CFL=6開始,ADI-FDTD的計算優(yōu)勢開始體現(xiàn),當(dāng)CFL=12時,采用ADI-FDTD方法的計算時間為FDTD的54.52%.比較表1和表2還可以發(fā)現(xiàn),采用1 mm剖分網(wǎng)格的計算效率要高于2 mm的情形,同樣當(dāng)CFL=12,采用1 mm網(wǎng)格的計算時間僅為傳統(tǒng)FDTD的24.9%,采用2 mm網(wǎng)格的計算時間為傳統(tǒng)FDTD方法的54.4%,從而可以得出采用ADI-FDTD方法另一個好處就是,分辨率越高,計算效率越高,這于成像是非常有利的.
表2 情形2ADI-FDTD與FDTD計算時間比較
3.3 情形3
存在多個腫瘤的情形,位于(90,30)和(52,20)處分別有一個5個邊長為2 mm的方柱組成的“十”字形柱體腫瘤,元胞尺寸取Δx=Δy=2 mm.仿真分別采用FDTD方法和ADI-FDTD方法進(jìn)行實(shí)驗(yàn),接下來分別采用FDTD方法和ADI-FDTD方法進(jìn)行仿真實(shí)驗(yàn),FDTD時間步數(shù)取1 500,天線T1上放置激勵源,天線R1上接收的電磁場如圖4(a),圖4(b)和圖4(c)所示.
(a) Ez
(b) Hx
(c) Hy圖4 情形3(R1上接收到的場)
從圖4(a),圖4(b)和圖4(c)仿真結(jié)果我們可以看出,FDTD和ADI-FDTD的計算結(jié)果同樣吻合得很好,多個腫瘤時采用ADI-FDTD的計算結(jié)果同樣可行,說明ADI-FDTD對于乳腺癌成像的適用性很好,為我們進(jìn)一步的成像研究奠定了基礎(chǔ).
對乳腺癌這種高對比問題進(jìn)行微波斷層成像,采用傳統(tǒng)的FDTD算法進(jìn)行分析時,由于受到CFL穩(wěn)定性條件的制約,時間步長受限于空間間隔的大小,計算時間過長,尤其是在高維情形下更為明顯.
為了提高檢測的效率,增強(qiáng)成像的臨床實(shí)用性,本文采用ADI-FDTD方法代替?zhèn)鹘y(tǒng)的FDTD方法進(jìn)行前向計算,適當(dāng)?shù)卦黾訒r間步長,縮短計算時間,采用二維半圓無限長柱狀乳房模型,分別對單個圓形、單個“十”字型和多個“十”字型的腫瘤模型進(jìn)行了測試,通過實(shí)驗(yàn)得出:1) 采用ADI-FDTD算法所得的結(jié)果和傳統(tǒng)FDTD吻合地很好,在保證現(xiàn)有成像精度的同時,大大縮短了成像時間,提供了即時成像的可能性;2) 當(dāng)網(wǎng)格剖分越細(xì)時,ADI-FDTD所體現(xiàn)出的效率越高,這于高精度成像是非常有利的,大大地提高了微波斷層成像技術(shù)的臨床可應(yīng)用性.
本文所做的研究內(nèi)容是前期的正向計算,為后續(xù)的成像問題奠定了基礎(chǔ),下一步的工作就是將ADI-FDTD方法運(yùn)用于成像分析中.當(dāng)然前期的工作還有一些問題需要考慮,比如說邊界條件的選取、天線的具體設(shè)計、激勵源的選擇等等,這些也是我們在今后的研究過程中所需考慮的實(shí)際問題.
[1]FEAR E C.Microwave imaging of the breast[J].Technology in cancer research &treatment, 2005, 4(1):69-82.
[2]LIU G D, ZHANG Y R.An overview of active microwave imaging for early breast cancer detection [J].
Journal of Nanjing University of Posts and Telecommunications, 2010, 30(1):64-70.
[3]BIN GUO.Microwave techniques for breast cancer detection and treatment [D].Gainesville:University of Florida, 2007.
[4]HAGNESS S C, TAFLOVE A, BRIDGES J E.Two-dimensional FDTD analysis of a pulsed microwave confocal system for breast cancer detection:fixed-focus and antenna-array sensors[J].IEEE transactions on biomedical engineering, 1999, 47(5):783-791.
[5]FEAR E C, L XU, HAGNESS S C.Confocal microwave imaging for breast cancer detection:localization of tumors in three dimensions[J].IEEE transactions on biomedical engineering, 2002, 49(8):812-822.
[6]LI X, HAGNESS S C.A confocal microwave imaging algorithm for breast cancer detection[J].IEEE Microwave and wireless components letters, 2001, 11(3):130-132.
[7]BOND E J, LI X, HAGNESS S C, et al.Microwave imaging via space-time beam-forming for early detection of breast cancer[J].IEEE transactions on antennas and propagation, 2003, 51(8):1690-1705.
[8]MEANEY P M, FANNING M W, LI D.A clinical prototype for active microwave imaging of the breast[J].IEEE transactions on microwave theory and techniques, 2000, 48(11):1841-1853.
[9]RUBAK T, MEANEY P M, MEINCKE P.Nonlinear microwave imaging for breast-cancer screening using Gauss-Newton's method and the CGLS inversion algorithm[J].IEEE transactions on antennas and propagation, 2007, 55(8):2320-2331.
[10]LIU Q H, YU C, STANG J.Experimental and numerical investigations of a high-resolution 3D microwave imaging system for breast cancer detection[C]//Antennas and Propagation Society International Symposium, 2007:2192.DOI:10.1109/APS.2007.4395963
[11]TAKENAKA T, JIA H, TANAKA T.Microwave imaging of electrical property distributions by a forward-backward time-stepping method[J].Journal of electromagnetic waves and applications, 2000, 14(12):1609-1626.
[12]劉廣東, 張業(yè)榮.二維有耗色散介質(zhì)的時域逆散射方法[J].物理學(xué)報, 2010, 59(10):6969-6979.
LIU G D, ZHANG Y R.Time-domain inverse scattering problem for two-dimensional frequency-dispersive lossy media[J].Acta physica sinica, 2010, 59(10):6969-6979.(in Chinese)
[13]劉廣東, 張業(yè)榮.乳腺癌檢測的三維時域微波斷層成像方法[J].電波科學(xué)學(xué)報, 2010, 25(6):1175-1181.
LIU G D, ZHANG Y R.Three-dimensional microwave tomography imaging method in the time-domain for detection of breast cancer[J].Chinese journal of radio science, 2010, 25(6):1175-1181.(in Chinese)
[14]WINTERS D W, BOND E J, BARRY D V, et al.Estimation of the frequency-dependent average dielectric properties of breast tissue using a time-domain inverse scattering technique[J].IEEE transactions on antennas and propagation, 2006, 54(11):3517-3528.
[15]唐濤, 廖成, 楊丹.FDTD求解高功率微波大氣傳播問題的可行性研究[J].電波科學(xué)學(xué)報, 2010, 25(1):122-126.
TANG T, LIAO C, YANG D.Feasibility of solving high-power microwave propagation in the atmosphere using FDTD method[J].Chinese journal of radio science, 2010, 25(1):122-126.(in Chinese)
[16]高本慶, 薛正輝, 任武.FDTD計算中關(guān)于低頻激勵源問題的研討[J].電波科學(xué)學(xué)報, 2009, 24(2):213-217.
GAO B Q, XUE Z H, REN W.Study on low frequency exciting source in FDTD simulation[J].Chinese journal of radio science, 2009, 24(2):213-217.(in Chinese)
[17]WANG S, CHEN J.Multigrid ADI method for two-dimensional electromagnetic simulations[J].IEEE transactions on antennas propagation, 2006, 54(2):715-720.
[18]張玉廷, 蔡智, 劉勝.交替隱式時域有限差分分析有耗介質(zhì)電磁波傳播[J].電波科學(xué)學(xué)報, 2011, 26(6):1088-1094.
ZHANG Y T, CAI Z, LIU S.Electromagnetic wave propagation in loss media based on ADI-FDTD[J].Chinese journal of radio science, 2011, 26(6):1088-1094.(in Chinese)
[19]湯煒, 李清亮, 焦培南, 等.ADI-FDTD在二維散射問題中的應(yīng)用[J].電波科學(xué)學(xué)報, 2003, 18(6):620-624.
TANG W, LI Q L, JIAO P N, et al.Two dimension scattering analysis using ADI-FDTD method[J].Chinese journal of radio science, 2003, 18(6):620-624.(in Chinese)
陳碧云 (1986-),女,江蘇人,南京郵電大學(xué)電子科學(xué)與工程學(xué)院博士生,主要研究方向?yàn)殡姴▊鞑?、電磁逆散射理論及其?yīng)用等.
張業(yè)榮 (1963-),男,安徽人,南京郵電大學(xué)電子科學(xué)與工程學(xué)院教授、博士生導(dǎo)師,主要研究方向?yàn)橐苿油ㄐ畔到y(tǒng)與設(shè)計、電磁逆散射、電磁場的數(shù)值計算和UWB信道等.
Application of ADI-FDTD method in the forward solver of early breast cancer detection
CHEN Biyun1ZHANG Yerong1PAN Xin2
(1.CollegeofElectronicScienceandEngineering,NanjingUniversityofPostsandTelecommunications,Nanjing210003,China;2.AviationUniversityofAirForce,Fuxin123000,China)
In the conventional finite-difference time-domain(FDTD) method, fine cells reduce the time-step size due to the Courant-Friedrich-Levy(CFL) stability condition, which results in an increase in computational effort, such as the central processing unit (CPU) time.In the alternating-direction implicit finite-difference time-domain(ADI-FDTD) method, a larger time-step size than allowed by the CFL stability condition limitation can be set because the algorithm of this method is unconditionally stable.Consequently, an increase in computational efforts caused by fine cells can be prevented.The simulation experiments data by the ADI-FDTD method were compared with that by the conventional FDTD method.The results show that simulation field by ADI-FDTD method agree quite well with that by the FDTD method, but the required CPU time can be much shorter than that for the FDTD method.
alternating-direction implicit finite-difference time-domain(ADI-FDTD);microwave tomography (MT);breast cancer
陳碧云, 張業(yè)榮, 潘鑫.ADI-FDTD方法在乳腺癌檢測正向計算中的應(yīng)用[J].電波科學(xué)學(xué)報,2016,31(5):1016-1022.
10.13443/j.cjors.2016013101
CHEN B Y, ZHANG Y R, PAN X.Application of ADI-FDTD method in the forward solver of early breast cancer detection [J].Chinese journal of radio science,2016,31(5):1016-1022.(in Chinese).DOI:10.13443/j.cjors.2016013101
2016-01-31
TN011
A
1005-0388(2016)05-1016-07
聯(lián)系人:陳碧云 E-mail:cby8612@126.com.