王雪峰,張自豪,陳興穌,王元慶
(1.伊犁師范大學(xué)網(wǎng)絡(luò)安全與信息技術(shù)學(xué)院,新疆 伊寧 835000;2.河南工業(yè)大學(xué)人工智能與大數(shù)據(jù)學(xué)院,河南 鄭州 450001;3.南京大學(xué)電子科學(xué)與工程學(xué)院,江蘇 南京 210093)
非視域成像技術(shù)是近年來的研究熱點(diǎn)技術(shù),是剛發(fā)展起來的一種間接成像技術(shù),在醫(yī)療成像、災(zāi)難救援、軍事反恐及城市街道安全偵查等都具有重要的研究意義。非視域成像技術(shù)通過激光照射中介面,主要是通過光在中介面進(jìn)行多次反射和散射來獲取隱藏物體的信息,再使用一定的重建算法對采集的數(shù)據(jù)進(jìn)行處理,才能完成非視域物體的三維重建。美國麻省理工學(xué)院的Velten等人[1-2]利用飛秒激光和條紋相機(jī)進(jìn)行拐角成像,首先使用反投影重建法,實(shí)現(xiàn)了非視域物體的三維重建。為了提高重建的精度,對反投影法的優(yōu)化算法被相繼提出,如使用時(shí)間選通門技術(shù)及單光子相機(jī)來對隱藏物體進(jìn)行逐層探測[3-5],使用快速反投影算法來提高重建的速度[6],利用橢圓投影特點(diǎn)進(jìn)行橢圓模型消除重建的偽影[7],文獻(xiàn)[8]提出了兩個(gè)權(quán)重因子應(yīng)用在濾波反投影算法中,減少了條紋效應(yīng),并能有效應(yīng)用到噪聲數(shù)據(jù)中。采用f-k偏移方法[9]為鏡面反射對象產(chǎn)生更高質(zhì)量的重建;O′Toole等人[10-11]使用光錐變換及共聚焦光路對隱藏物體進(jìn)行探測;Xin等人[12]使用費(fèi)馬流算法得到了較高的重建精度。
由于反投影法需要采集足夠的數(shù)據(jù)(全部角度范圍),才能較好的完成三維重建,所以需要光源進(jìn)行掃描式數(shù)據(jù)采集,需要大量時(shí)間,但是在實(shí)際的應(yīng)用中,特別是對于較危險(xiǎn)的物體,如反恐及災(zāi)害救援,需要快速、有效的進(jìn)行數(shù)據(jù)采集,并能夠快速完成三維重建,這就對非視域成像技術(shù)及重建算法提出了新的挑戰(zhàn)。本文針對稀疏角度數(shù)據(jù)問題,提出了反投影修正MLEM(最大似然期望值最大化)算法,能夠有效提高非視域物體三維重建的精度。
現(xiàn)有大多數(shù)成像方式為反射式成像,也稱為拐角成像,本文采用透射式成像方式,成像設(shè)置的示意圖如圖1所示,激光照射中介面一點(diǎn)(如圖中D點(diǎn)),光透射中介面后,中介面選擇擴(kuò)散膜(不能透過擴(kuò)散膜直接觀測到后面的物體),繼續(xù)散射到擴(kuò)散膜后面的隱藏物體,光會再次被隱藏物體進(jìn)行散射和反射,并返回到擴(kuò)散膜后被探測器APD接收,其中圖1中參考光是檢測光源的出光時(shí)間,作為探測器APD接收信號的時(shí)間基準(zhǔn)。本文提出的基于APD陣列的透射式非視域激光成像系統(tǒng)具有成像結(jié)構(gòu)簡單、操作方便且無需場景掃描、能夠快速完成數(shù)據(jù)采集的特點(diǎn)。探測器布置在擴(kuò)散膜附近,直接探測隱藏目標(biāo)散射回來的光信息。一次激光照射后,APD陣列接收來自隱藏目標(biāo)多個(gè)方向的信息,這相當(dāng)于進(jìn)行了多點(diǎn)信息接收,因此不需要進(jìn)行場景掃描,使數(shù)據(jù)采集一次完成。
圖1 本文提出的透射式非視域激光成像示意圖
探測器陣列的排布方式采用矩形陣列方式,如圖2所示,以3×3的陣列為例,不同的探測器陣列(如圓環(huán)形狀和十字形狀等)采集到的數(shù)據(jù)會有不同,獲取的隱藏目標(biāo)的信息也會有所不同,對于重建算法有一定的影響。本文中使用的探測器單元由4個(gè)APD單元組成,詳細(xì)的APD單元排布方式如圖1所示。4個(gè)APD排布在擴(kuò)散膜的四個(gè)對角位置,排布成2行2列的矩陣陣列的方式,四個(gè)APD分別標(biāo)記為APD11、APD12、APD21和APD22。矩形陣列探測的數(shù)據(jù)表示為:
本文提出的透射式非視域成像的光傳輸過程如圖3所示。
圖2 矩形陣列的探測器排布方式
圖3(a)中,首先激光照射擴(kuò)散膜D點(diǎn),光經(jīng)過擴(kuò)散膜后,繼續(xù)傳輸?shù)诫[藏物體上(如點(diǎn)S1、S2和S3),傳輸路徑為D→S1、D→S2和D→S3,傳輸距離為rd1、rd2和rd3;其次光經(jīng)過隱藏物體后再次傳輸?shù)綌U(kuò)散膜,并被探測器APD接收,如其中一個(gè)APD探測單元在A點(diǎn),則傳輸路徑為S1→A、S2→A和S3→A,傳輸距離為rs1、rs2和rs3;最后在A點(diǎn)的APD探測器單元接收信號的路徑為D→S1→A、D→S2→A和D→S3→A,傳輸距離為rd1+rs1、rd2+rs2和rd3+rs3,則接收信號的三個(gè)時(shí)間點(diǎn)(ts1、ts2和ts3)如圖3(b)上面所示,對應(yīng)圖3(a)中隱藏物體的三個(gè)點(diǎn)(S1、S2和S3),則ts1=(rd1+rs1)/c、ts2=(rd2+rs2)/c和ts3=(rd3+rs3)/c,其中c為光速。
圖4 非視域成像的橢圓投影原理圖
圖4中兩個(gè)定點(diǎn)是光源D點(diǎn)和探測器A點(diǎn),隱藏物體上的點(diǎn)即是某個(gè)橢圓上的點(diǎn),如橢圓f1上的動點(diǎn)p3;隱藏物體上的兩點(diǎn)p1和p2在同一個(gè)橢圓f2上,則p1和p2會在同一個(gè)時(shí)間點(diǎn)返回信息到探測器上,即光經(jīng)過D點(diǎn)后傳輸?shù)诫[藏物體上兩點(diǎn)p1和p2,并投影到探測器上,這就是非視域成像中的橢圓投影方式。后期的圖像重建也是根據(jù)反投影方法的原理,首先將隱藏物體劃分為三維像素塊(稱為體素,如圖4中p1、p2和p3),然后根據(jù)橢圓反投影方法對隱藏物體的每個(gè)體素進(jìn)行重建,最終完成整個(gè)物體的三維重建。
MLEM算法是最大似然期望值最大化算法(maximum likelihood expectation maximization),它使用期望最大化EM(expectation maximization)算法[13-14]來完成最大似然的估計(jì),是一種統(tǒng)計(jì)迭代算法,以“參數(shù)估計(jì)理論”為基礎(chǔ),基于觀測數(shù)據(jù)統(tǒng)計(jì)模型的重建算法,廣泛應(yīng)用于CT及PETCT圖像重建[15-16],CT圖像重建的主要思想也是反投影[17],所以也需要完整的數(shù)據(jù)才能完成較好的重建,但實(shí)際中無法完成180°范圍的數(shù)據(jù)采集,且需要采集的時(shí)間間隔要小、均勻才能更好的重建圖像,針對CT圖像重建中稀疏角度數(shù)據(jù)的問題,可以使用迭代重建算法。
MLEM算法是根據(jù)最大似然估計(jì)理論,假設(shè)認(rèn)為每次橢圓投影符合泊松分布,根據(jù)似然函數(shù),求投影數(shù)據(jù)的最大概率,推導(dǎo)得出計(jì)算公式,而此公式正好為非線性方程,很難給出解析解,只能通過迭代的方法來求近似解,通常采用期望最大化(EM)來求解。由此,推導(dǎo)得出了MLEM算法,迭代公式[18]如下:
(1)
其中,重建體素設(shè)為xj;重建體素個(gè)數(shù)為M個(gè);投影矩陣中橢圓個(gè)數(shù)為i個(gè);投影矩陣元素aij表示第j個(gè)體素xj對第i個(gè)橢圓的貢獻(xiàn)值;pi表示第i個(gè)橢圓的投影值。
基于APD陣列的非視域成像MLEM算法的具體步驟如下:
(4)計(jì)算第j個(gè)重建體素xj的修正值:
(2)
(5)對重建體素xj的值進(jìn)行修正:
(3)
這里用穿過該體素的所有橢圓對它進(jìn)行修正,即完成一輪迭代。
(6)k=k+1;將以上迭代的結(jié)果作為初值,重復(fù)步驟(2)到(5)進(jìn)行新一輪迭代,直到取得符合收斂條件為止。
針對非視域成像的特點(diǎn),僅能在稀疏角度下獲取數(shù)據(jù),提出了BP-MLEM(Back Projection -MLEM)迭代重建算法,更加適用于非視域成像,改進(jìn)的方法對偽影和噪聲的抑制較好,且具有更快的收斂速度。
MLEM算法中比較重要的迭代過程是修正方法,它決定了每次迭代結(jié)果后的修正方向,為了更快更準(zhǔn)確的得到收斂結(jié)果,需要迭代方向每次都往正確的結(jié)果修正,這就需要更好的修正方法,MLEM算法修正方法的修正值為:
采用上次迭代結(jié)果和修正值進(jìn)行相乘的方法進(jìn)行迭代,即下一次迭代值是上一次迭代值與反投影乘積的結(jié)果。這樣迭代后圖像中較大的值增長越大,較小的值增長較慢,圖像的對比度越來越大,即迭代速度快,但有時(shí)速度過快,容易造成過收斂,不能得到理想的收斂結(jié)果。
本文提出的基于反投影迭代BP-MLEM算法,迭代公式如下:
(4)
其中,重建體素設(shè)為xj;重建體素個(gè)數(shù)為M個(gè);投影矩陣中橢圓個(gè)數(shù)為i個(gè);投影矩陣元素aij表示第j個(gè)體素xj對第i個(gè)橢圓的貢獻(xiàn)值;pi表示第i個(gè)橢圓的投影值。
基于APD陣列的非視域成像BP-MLEM算法的具體步驟如下:
(4)計(jì)算第j個(gè)重建體素xj的修正值為:
(5)
(5)對重建體素xj值進(jìn)行修正,方法為:
(6)
這里用穿過該體素的所有橢圓對重建體素xj進(jìn)行修正,并將反投影BPΔi加入到修正方法中,這樣就完成一輪迭代。
(6)k=k+1;將以上迭代的結(jié)果作為初值,重復(fù)步驟(2)到(5)進(jìn)行新一輪迭代,直到取得符合收斂條件為止。
實(shí)驗(yàn)設(shè)置如圖1所示,隱藏目標(biāo)選擇“十”字形狀,距離探測器220 cm,4個(gè)APD探測器的位置如圖1所示。分別使用反投影算法、MLEM算法和BP-MLEM算法進(jìn)行實(shí)驗(yàn),為了更好的進(jìn)行對比實(shí)驗(yàn),本文進(jìn)行了多次迭代重建,分別給出不同迭代次數(shù)的對比結(jié)果。
圖5顯示了使用傳統(tǒng)的反投影算法得到重建結(jié)果,反投影算法中使用了R-S-L濾波算子,從圖中可以看出,反投影的重建結(jié)果存在較大的偽影,重建結(jié)果精度不高。
圖5 濾波反投影算法重建結(jié)果圖
圖6顯示的是使用MLEM算法對非視域物體“十”字的重建結(jié)果,分別顯示了4次不同迭代次數(shù)(分別是第3、4、9、20次迭代),從圖中可以看出,當(dāng)?shù)?次迭代(圖6(a))時(shí),重建結(jié)果圖的偽影在減少,但是重建結(jié)果圖并不理想,圖像中心點(diǎn)強(qiáng)度很高,而邊緣強(qiáng)度很低(接近背景色),不能較好的重建“十”字形狀的結(jié)果;再次進(jìn)行6次迭代后,即第9次迭代(圖6(c)),重建結(jié)果的偽影比第3次迭代少,但是重建的結(jié)果仍然是中心點(diǎn)強(qiáng)度過大,而邊緣強(qiáng)度太小;當(dāng)達(dá)到第20次迭代時(shí),重建結(jié)果反而沒有最初的迭代結(jié)果好,幾乎無法重建“十”形狀結(jié)果,這就是MLEM算法迭代速度過快,從而產(chǎn)生過收斂的結(jié)果。
本文提出的BP-MLEM算法的重建結(jié)果如圖7所示,同時(shí)也給出了4個(gè)不同迭代次數(shù)(分別是第3、4、9、20次迭代),總體來看,4次迭代的重建結(jié)果都能夠重建出“十”字形狀圖,并且隨著迭代次數(shù)的增加,重建偽影在不斷減少,重建的精度也在不斷增加;當(dāng)?shù)螖?shù)達(dá)到4次時(shí)(圖7(b)),物體“十”字形狀就得到了較好的重建結(jié)果圖;當(dāng)?shù)螖?shù)不斷增加時(shí),重建結(jié)果沒有出現(xiàn)MLEM算法的過收斂結(jié)果,而是達(dá)到了更好的重建效果,重建結(jié)果圖像趨于穩(wěn)定。實(shí)驗(yàn)表明,BP-MLEM算法能夠在較少的迭代次數(shù)下,得到較好的重建結(jié)果;并且當(dāng)?shù)螖?shù)增加時(shí),重建算法能夠得到較穩(wěn)定的重建結(jié)果圖,不會出現(xiàn)過收斂等問題。
為了更加準(zhǔn)確的對比實(shí)驗(yàn)結(jié)果,本文使用了SSIM(Structural SIMilarity)結(jié)構(gòu)相似度來評價(jià)重建結(jié)果圖的圖像質(zhì)量,對比實(shí)驗(yàn)結(jié)果顯示在圖8中,分別對使用BP-MLEM算法的重建結(jié)果圖(圖7)和MLEM算法的重建結(jié)果圖(圖6)計(jì)算SSIM值,圖8中分別計(jì)算了第3、4、9、20次迭代的SSIM值,其中第一個(gè)值(橫坐標(biāo)標(biāo)識為BP)是使用反投影算法重建結(jié)果圖(圖5)。從圖8中可以看出,利用BP-MLEM算法重建結(jié)果圖的SSIM值都高于MLEM算法重結(jié)果圖的SSIM值,而且隨著迭代次數(shù)的增加,BP-MLEM算法的重建結(jié)果圖的SSIM也在不斷增加,表明重建結(jié)果圖的質(zhì)量在不斷增加;而MLEM算法的重建結(jié)果圖的SSIM值雖然剛開始也在不斷增加,但增加的趨勢較慢,且當(dāng)?shù)螖?shù)為第9次和20次時(shí),SSIM值基本沒有增加,表明重建結(jié)果圖的質(zhì)量并沒有隨著迭代次數(shù)的增加而不斷增加。
圖8 BP-MLEM算法與MLEM算法重建結(jié)果圖在 不同迭代次數(shù)下的SSIM值對比圖
圖像重建算法在非視域成像系統(tǒng)中具有重要的作用,它關(guān)乎著整個(gè)系統(tǒng)的實(shí)現(xiàn)結(jié)果,對實(shí)驗(yàn)系統(tǒng)的驗(yàn)證起到了關(guān)鍵作用。本文構(gòu)建了基于APD陣列的透射式成像系統(tǒng),快速進(jìn)行隱藏物體探測,一次完成數(shù)據(jù)采集,不需要使用光源進(jìn)行掃描成像;并針對非視域成像系統(tǒng)僅能在稀疏角度下獲取不完整數(shù)據(jù)的特點(diǎn),提出了BP-MLEM迭代算法,實(shí)驗(yàn)結(jié)果表明,相比反投影算法,能夠更好的去除偽影;并且能在較少的迭代次數(shù)下,得到較好的重建結(jié)果圖;當(dāng)重建次數(shù)不斷增加時(shí),重建結(jié)果趨于穩(wěn)定,表明BP-MLEM算法較穩(wěn)定,不會出現(xiàn)過收斂等問題。