武 晶,陸 明,陳 平*,潘晉孝
(1.中北大學(xué) 信息與探測處理山西省重點實驗室,太原 030051;2.中北大學(xué) 理學(xué)院,太原 030051;3.中國航天科工集團有限公司六院601所,呼和浩特 010076)
固體火箭發(fā)動機結(jié)構(gòu)簡單、運行可靠、維護使用方便,已成為當今戰(zhàn)略、戰(zhàn)術(shù)導(dǎo)彈武器系統(tǒng)首選動力裝置[1]。但由于環(huán)境、載荷作用以及工藝制作等原因?qū)е碌乃幹毕莺芸赡軙斐苫鸺l(fā)射失敗。其中,裂紋會在推進劑燃燒過程中產(chǎn)生額外的燃燒面積,對發(fā)動機內(nèi)彈道性能造成嚴重影響[2-3]。因此,為進一步觀察裂紋對燃面影響,從而判斷發(fā)動機安全性能,開展點火試驗環(huán)境下裂紋缺陷檢測研究是至關(guān)重要的。
為了實現(xiàn)藥柱燃燒退移過程中的燃面退移檢測,王晶、王世輝等[4-5]利用微波CT技術(shù)對燃面進行重建,有效表征了燃面退移變化,但其重建圖像分辨率較低,難以實現(xiàn)缺陷的有效檢測。X射線高速實時熒屏分析技術(shù)(RTR)以其具有穿透能力強、非接觸、實時記錄數(shù)字圖像等優(yōu)點,被用于研究不透明高壓燃燒器或發(fā)動機內(nèi)復(fù)雜工作過程,近幾十年在固體火箭發(fā)動機點火試驗內(nèi)部工作的研究中得到了廣泛的應(yīng)用[6-8]。但由于其為直接投影成像,缺陷的空間位置無法確定[9]。相比RTR技術(shù),CT成像能更直觀、準確地獲取物體內(nèi)部結(jié)構(gòu),實現(xiàn)發(fā)動機內(nèi)部結(jié)構(gòu)的有效可視化[10]。但是,由于試車環(huán)境下無法像常規(guī)的CT成像模式,實現(xiàn)發(fā)動機旋轉(zhuǎn)或布置多視角射線——探測器對,只能采取超稀疏成像方式。因此,亟需提出一種能夠?qū)崟r對發(fā)動機內(nèi)結(jié)構(gòu)進行超稀疏CT三維重建的方法,從而在點火試驗環(huán)境下對藥柱缺陷進行更直觀地檢測。
考慮到傳統(tǒng)CT重建算法需要對物體進行全方位掃描來完成重建,在醫(yī)學(xué)、工業(yè)應(yīng)用上都有很大限制,研究人員進行了稀疏角CT重建算法研究。自1970年,GORDON等[11]提出了迭代重建算法(ART)之后,為解決其重建圖像噪聲大、含有偽影、細節(jié)丟失等問題。SIDKY等[12]在此基礎(chǔ)上引入了全變分(total variation,TV)正則項,在一定程度上消除了條狀偽影。BUADES等[13]提出了非局部均值(Non-local Means,NLM)去噪算法來消除噪聲。JING等[14]在文獻[12-13]基礎(chǔ)上,提出了ART-NLM算法,進一步改善了重建圖像質(zhì)量。但是,這些方法無法兼顧稀疏性與成像質(zhì)量,在超稀疏CT重建時,得到的重建圖像可能會產(chǎn)生細節(jié)丟失和邊緣平滑,無法滿足實時成像需求。近年來,深度學(xué)習因其強大的學(xué)習能力與計算能力,在圖像處理、目標檢測等領(lǐng)域被廣泛使用,研究人員便將其用來解決稀疏角重建的病態(tài)問題。HENZLER等[15]基于深度卷積神經(jīng)網(wǎng)絡(luò),從單張X射線圖像中重建獸頭骨三維結(jié)構(gòu)。YING等[16]利用生成對抗網(wǎng)絡(luò)(Generative Adversarial Network,GAN)模型,從兩張正交的胸腔X射線圖像中重建CT圖像。SHEN等[17]利用數(shù)據(jù)驅(qū)動的方法,實現(xiàn)了單視圖重建患者體積數(shù)據(jù)。但是,基于深度學(xué)習的方法是通過模型在訓(xùn)練階段學(xué)習數(shù)據(jù)中隱含的先驗知識來進行預(yù)測,而在發(fā)動機未裝備時無法獲取點火試驗環(huán)境下內(nèi)結(jié)構(gòu)缺陷的先驗信息。因此,現(xiàn)有深度學(xué)習方法無法實時重建出缺陷,使得設(shè)計人員不能準確判別缺陷的性質(zhì)和對發(fā)動機可能造成的危害。
針對以上問題,本文結(jié)合深度學(xué)習技術(shù)和視差原理,提出一種可以快速、實時重建發(fā)動機藥柱缺陷的雙視角CT重建算法。本文提出的算法模型僅利用兩個視角下的投影信息,即可實時重建帶有缺陷信息的三維體。與現(xiàn)有的超稀疏深度重建網(wǎng)絡(luò)模型相比,無需缺陷區(qū)域的先驗知識,且通過實驗表明,該算法對自然噪聲、發(fā)動機振動影響具有魯棒性,有效實現(xiàn)了點火試驗環(huán)境下對固體火箭發(fā)動機藥柱缺陷的實時檢測。
在本文提出的方法中,首先考慮到實際工程環(huán)境下采集系統(tǒng)的布置成本以及單一視角的信息遮擋問題,采用相互正交的雙視角成像系統(tǒng)來保證物體特征分布的正交性,利用一個基于深度學(xué)習的雙視角CT重建模型完成雙視角投影數(shù)據(jù)到三維結(jié)構(gòu)體的重建。同時,將實時采集到的兩張正交的投影數(shù)據(jù)通過一個基于無監(jiān)督學(xué)習的異常檢測網(wǎng)絡(luò)進行缺陷檢測與定位,得到缺陷在二維圖像中的坐標信息。最后,利用視差原理,獲取缺陷在三維中的深度信息,從而在三維結(jié)構(gòu)體上重建缺陷。該算法流程圖如圖1所示。
圖1 算法流程圖
雙視角投影數(shù)據(jù)采集系統(tǒng)如圖2所示,采用兩個射線源正交的排布方式,同時同步進行成像。
圖2 雙視角投影數(shù)據(jù)采集系統(tǒng)
為了實現(xiàn)雙視角下對發(fā)動機內(nèi)結(jié)構(gòu)的CT重建,本文采用基于數(shù)據(jù)驅(qū)動的稀疏投影快速重構(gòu)理論,在文獻[17]基礎(chǔ)上對其模型進行改進,考慮到實際工程環(huán)境下,采集系統(tǒng)數(shù)據(jù)采集速度快、數(shù)據(jù)量大,將生成網(wǎng)絡(luò)中反卷積層的三維反卷積運算用二維反卷積代替,即將三維深度信息與通道信息編碼在同一維度中,同時減少網(wǎng)絡(luò)層數(shù),實現(xiàn)網(wǎng)絡(luò)輕量化設(shè)計,以保證實時重建效率。本文采用的網(wǎng)絡(luò)架構(gòu)如圖3所示,包括表征網(wǎng)絡(luò)、轉(zhuǎn)換網(wǎng)絡(luò)、生成網(wǎng)絡(luò)三部分,圖中每一層下面的數(shù)字為通道數(shù)。
圖3 雙視角CT網(wǎng)絡(luò)架構(gòu)圖
表征網(wǎng)絡(luò)由3個下采樣塊構(gòu)成,每個下采樣塊由卷積核為4×4,步長為2的2D卷積層、2D批標準化層、ReLU層、卷積核為3×3,步長為1的2D卷積層、2D批標準化層和ReLU層按順序組成,從輸入的雙視角投影數(shù)據(jù)中提取特征信息。此外,在應(yīng)用第二個ReLU層之前,將每個下采樣塊的輸入與輸出通過快捷路徑連接起來,進行簡單的同等映射。
轉(zhuǎn)換網(wǎng)絡(luò)是由卷積核為1×1,步長為1的2D卷積層、ReLU層以及卷積核為1×1,步長為1的2D反卷積層、ReLU層構(gòu)成。對表征網(wǎng)絡(luò)的輸出進行跨維度映射,將二維特征與三維特征聯(lián)系起來,以便于生成網(wǎng)絡(luò)進行重建。
生成網(wǎng)絡(luò)由3個上采樣塊構(gòu)成,每個上采樣塊由卷積核為4×4,步長為2的2D反卷積層、2D批標準化層、ReLU層、卷積核為3×3,步長為1的2D反卷積層、2D批標準化層和ReLU層按順序組成,與下采樣塊相同,每個上采樣塊的輸入輸出之間也采用快捷路徑鏈接,生成網(wǎng)絡(luò)將學(xué)習到的空間特征生成對應(yīng)的CT圖像。
目前在工業(yè)生產(chǎn)中,異常檢測與定位算法大多用于工業(yè)產(chǎn)品的表面缺陷檢測[18-20]。但為了在點火試驗時可以實時從二維投影圖像中獲取到發(fā)動機內(nèi)部缺陷的信息,本文在文獻[18]基礎(chǔ)上,在其異常檢測與定位模塊之后,增加了缺陷角點輸出模塊,用于發(fā)動機內(nèi)部缺陷的檢測與定位。
在訓(xùn)練階段,如圖4所示,該算法可分為特征提取和特征學(xué)習兩部分。在特征提取階段,為避免繁瑣的神經(jīng)網(wǎng)絡(luò)優(yōu)化,使用預(yù)訓(xùn)練過的ResNet網(wǎng)絡(luò)來提取訓(xùn)練集中用于異常檢測的特征。在特征學(xué)習階段,為了學(xué)習圖像位置(i,j)處的正常圖像特征,將在特征提取時從N張訓(xùn)練圖像位置(i,j)處得到的特征向量集進行處理,使其符合多元高斯分布。假設(shè)向量集Xij符合多元高斯分布N(μij,εij),其中μij為Xij的樣本均值,則樣本的協(xié)方差εij如式(1)所示:
圖4 異常檢測算法原理圖
(1)
其中,正則化項φI使樣本協(xié)方差矩陣εij滿秩且可逆。最后,通過以上計算,使得每個圖像塊都與一個多元高斯分布相關(guān)聯(lián)。在測試階段,該算法利用馬氏距離M(xij)計算測試圖像中位置(i,j)處的異常分數(shù)。M(xij)是通過衡量測試圖像塊特征向量xij與訓(xùn)練階段學(xué)習正常樣本得到的分布N(μij,εij)之間的距離來與合適的閾值γ做對比,從而檢測并定位異常,其計算公式如下:
(2)
因此,通過計算得到M=(M(xij))1
(3)
完成異常定位后,假設(shè)在二維投影圖中檢測到的異常區(qū)域如圖5所示,該算法進一步輸出了兩幅不同視角投影異常區(qū)域四個角點的坐標(a1,b1)、(a2,b2)、(a3,b3)、(a4,b4)。
圖5 模擬投影數(shù)據(jù)缺陷定位圖
需要注意的是,該算法中閾值γ在沒有真實缺陷標注作標簽時,需要手動調(diào)節(jié)來尋找最優(yōu)參數(shù),手動選取參數(shù)會在一定程度上導(dǎo)致缺陷角點坐標輸出不準確。
在機器視覺三維重建領(lǐng)域中,由于計算機接受的客觀世界信息以及攝像機拍攝的圖像都是二維的,需要通過立體視覺技術(shù)模擬人眼來獲取三維真實物體的深度信息。其中,雙目視覺由于其計算簡單、有效性高,其實用性顯著高于其他方法[21]。雙目視覺是利用視差原理從攝像機拍攝的兩幅不同視角圖像中獲取對應(yīng)像素點的視差,再結(jié)合攝像機標定參數(shù)來求取深度信息,從而重建真實場景三維結(jié)構(gòu),在航空測繪、工業(yè)檢測等多個領(lǐng)域中具有廣泛應(yīng)用[22]。
雖然X射線立體成像與可見光立體視覺有一定差異,但其成像系統(tǒng)的核心都是不同角度下二維圖像的獲取。因此,雙目視覺三維重建原理也可用于本文雙視角下缺陷區(qū)域的CT重建中[23]。
由以上分析可知,利用視差原理,可從兩幅不同角度二維投影中恢復(fù)深度信息,從而重建三維立體圖像。本文提出的固體火箭發(fā)動機點火試車環(huán)境下采用雙視角投影數(shù)據(jù)采集系統(tǒng),利用方向正交的兩束射線源同時對物體進行投影來獲取雙視角投影圖,成像示意圖如圖6所示。
圖6 成像示意圖
圖6中,O為物體幾何中心,兩個射線源之間的夾角θ=90°,O1、O2為各自的投影中心,Z軸與XY坐標系垂直且為物體的旋轉(zhuǎn)軸,物體中某點(x,y)在兩幅投影中的位置坐標(x1,y1,z1)、(x2,y2,z3)之間的關(guān)系為
(4)
如果該點在兩幅投影中的橫坐標x1、x2已知,則該點的縱坐標y1、y2的計算公式為
(5)
y2=-x1sinθ+y1cosθ
(6)
由于平行投影,因此z1=z2,所以該點在三維中坐標(x1,y1,z1)已知。
在本文中,異常檢測與定位模型對實時采集的雙視角X射線斷層圖像進行檢測與定位后,獲取到了兩個相互正交視角下缺陷區(qū)域的各角點坐標,即對于每個角點來說,其三維坐標中x與z已知,通過以上公式計算得到深度信息y,即可獲得缺陷在三維結(jié)構(gòu)中的具體定位與大小。將其在上述生成的無缺陷三維結(jié)構(gòu)體中進行重建后,即可得到帶有缺陷的發(fā)動機三維圖像。
為了驗證本文所提方法的有效性,分別對仿真數(shù)據(jù)與模擬工件數(shù)據(jù)進行實驗。
仿真數(shù)據(jù)是通過對藥柱燃燒退移過程進行數(shù)學(xué)建模,并結(jié)合CT仿真平臺進行全角度投影數(shù)據(jù)采集,投影圖像大小為512×512,投影間隔1°,對采集到的360張投影數(shù)據(jù)用FDK算法重建作為真實CT圖像,重建大小為256×256×256。某一時刻藥柱仿真退移模型如圖7所示,模型大小為240×240×450。
(a)The simulation model (b)Top view of simulation model
為了構(gòu)造含有缺陷的藥柱模型,在上述模型200層到230層構(gòu)建長為22像素、寬為6像素的裂縫,同時為了保證模型對噪聲與發(fā)動機振動模糊具有一定魯棒性,加入0.01的高斯噪聲與運動位移為10像素的徑向振動,在同等條件下進行投影圖像采集。圖8中,第一、二行分別為構(gòu)造缺陷前后模型在0°、90°的投影圖像及第210層切片圖,圖8中紅色標注為缺陷區(qū)域。
(a)0°project image without defect (b)90°project image without defect (c)The 210th of slice without defect
此外,考慮到真實發(fā)動機CT影像受殼體材料、絕熱層等影響,在上述藥柱退移模型基礎(chǔ)上,粘接了環(huán)寬為5像素的殼體和環(huán)寬為3像素的絕熱層,得到的模型大小為256×256×450,缺陷構(gòu)造方式與投影數(shù)據(jù)采集條件同上。圖9中,第一、二行分別為構(gòu)造缺陷前后模型在0°、90°的投影圖像及第210層切片圖,圖9中紅色標注為缺陷區(qū)域。
(a)0°project image without defect (b)90°project image without defect (c)The 210th of slice without defect
模擬工件是采用玻璃纖維制作的發(fā)動機藥柱形態(tài)的假藥模型,使用YXLON FF20微焦點CT系統(tǒng)對其進行投影間隔為1°的全角度投影數(shù)據(jù)采集,并對其重建作為真實CT圖像,實驗參數(shù)設(shè)置如表1所示。
表1 實驗條件
圖10(a)為CT系統(tǒng)重建的假藥模型,為了模擬發(fā)動機點火試驗下產(chǎn)生的自然缺陷,在該模體中人為構(gòu)建了退移不均勻和裂縫兩種常見缺陷,在同等實驗條件下進行投影數(shù)據(jù)采集與重建。圖10(b)為CT系統(tǒng)重建的帶缺陷的假藥模型,兩種缺陷分別分布于A、B兩處位置,A處為退移不均勻?qū)е碌娜毕?B處為裂縫。
(a)CT image of non-defect (b)CT image with defect
為了保證模型具有一定魯棒性,在帶有缺陷投影圖像中加入運動位移為10像素的振動模糊及0.0001的高斯噪聲。上述的假藥模體投影數(shù)據(jù)如圖11所示,圖11中紅色標注為缺陷區(qū)域。
在圖11(c)、(d)中,由于圖像對比度較低,缺陷區(qū)域不明顯。為了更直觀的觀察缺陷,圖12分別給出了0°、90°下缺陷投影數(shù)據(jù)灰度值曲線,圖12中紅色標注為缺陷A、B處列灰度值曲線。
(a)0°project image without defects (b)90°project image without defects
(a)Column gray scale curve of 0°projection data
在對仿真數(shù)據(jù)與模擬工件數(shù)據(jù)進行實驗時,實驗環(huán)境為Windows10,處理器為Intel i5 6300HQ 4核,顯卡為NVIDIA GeForce GTX 960M,基于Pytorch框架,構(gòu)建雙視角CT模型與異常檢測模型,模型參數(shù)設(shè)置如表2所示。在MATLAB平臺上,構(gòu)建缺陷CT重建模型。實驗采用平均絕對誤差(MAE)、均方根誤差(RMSE)、結(jié)構(gòu)相似性(SSIM)作為評價指標,將本文方法重建的帶有缺陷的發(fā)動機三維體與真實缺陷CT圖像進行定量分析,MAE、RMSE值越小代表模型預(yù)測準確度越高,SSIM取值范圍為[0,1],越接近1,重建質(zhì)量越好。
表2 模型設(shè)置
本實驗選取仿真數(shù)據(jù)集中缺陷模型在0°、90°下的投影數(shù)據(jù)進行驗證。為了體現(xiàn)對文獻[12]中PatRecon模型改進的效果,統(tǒng)計了改進前后訓(xùn)練參數(shù)量、每秒浮點運算次數(shù)(Floating-point operations per second,FLOPS))以及測試時RMSE、SSIM值和重建時間,對比結(jié)果如表3所示。可看出,在對PatRecon模型改進前后,雖然成像質(zhì)量稍顯降低,但大幅度減少了模型參數(shù)計算量與時間消耗,有效提高了發(fā)動機點火試驗環(huán)境下實時采集數(shù)據(jù)的重建效率。
表3 模型效果對比
在進行異常檢測時,閾值γ取值為0.7,仿真藥柱與仿真發(fā)動機數(shù)據(jù)異常檢測可視化結(jié)果分別如圖13、圖14所示。由圖14可以看出,在有殼體與絕熱層粘接的情況下,也能有效檢測到藥柱缺陷。
(a)Defect detection results of 0°projection data
(a)Defect detection results of 0°projection data
圖15(a)~(c)與圖16(a)~(c)分別為藥柱仿真數(shù)據(jù)與發(fā)動機仿真數(shù)據(jù)經(jīng)過雙視角CT模型、本文方法、FDK算法重建缺陷CT的第120層切片圖??梢钥闯?僅通過基于深度學(xué)習的雙視角CT重建模型進行測試,缺陷區(qū)域無法被重建,而本文方法基本可以重建出裂縫。
(a)Dual-view CT model (b)The proposed method (c)FDK
以FDK算法重建結(jié)果作為真實缺陷CT圖像衡量本文結(jié)果,定量指標結(jié)果分別如表4、表5所示。由表中MAE、RMSE指標可看出,本文方法對仿真數(shù)據(jù)的預(yù)測準確度較高,且SSIM指標表明,本文方法受噪聲與發(fā)動機振動的影響較小,可通過兩張角度正交的投影圖像重建出質(zhì)量較高的發(fā)動機三維結(jié)構(gòu)。
表4 藥柱實驗結(jié)果定量指標
表5 發(fā)動機實驗結(jié)果定量指標
本實驗選取模擬工件數(shù)據(jù)集中缺陷模型在0°、90°下的投影數(shù)據(jù)進行驗證。在進行異常檢測時,閾值γ=0.6,異常檢測可視化結(jié)果如圖17所示。
(a)Defect detection results of 0°projection data
圖18、圖19分別為雙視角CT模型、本文方法、CT系統(tǒng)重建結(jié)果的第16層、第60層切片圖??梢钥闯?僅通過基于深度學(xué)習的雙視角CT重建模型進行測試,兩種缺陷均無法被重建,而本文方法基本可以重建出兩種缺陷。
(a)Dual-view CT model(b)The proposed method(c)CT system
(a)Dual-view CT model (b)The proposed method (c)CT system
以CT系統(tǒng)重建結(jié)果作為真實缺陷CT圖像衡量本文結(jié)果,定量指標結(jié)果如表6所示。
由表6結(jié)果可知,對于有多個缺陷且加噪加振動的模擬工件數(shù)據(jù),本文方法可以利用兩張角度正交的投影圖像較準確地重建出三維CT圖像,為后續(xù)實驗測試系統(tǒng)的設(shè)計提供了理論支撐。
表6 實驗結(jié)果定量指標
在以上實驗中,X射線雙視角成像時間與圖像三維重建時間如表7所示。
表7 成像與重建時間
由表7可知,三維重建時間小于成像時間,滿足實時檢測需求,可實現(xiàn)固體火箭發(fā)動機點火試驗時的動態(tài)檢測。
(1)本文選用仿真數(shù)據(jù)和模擬工件數(shù)據(jù)進行實驗驗證,實驗結(jié)果表明,該方法可以從雙視角投影數(shù)據(jù)中有效重建三維結(jié)構(gòu)。對于有多個缺陷的模擬工件,也可以較準確地進行重建,可實現(xiàn)對固體火箭發(fā)動機點火試驗時的動態(tài)監(jiān)測目的。
(2)為接近真實工業(yè)環(huán)境,在對投影數(shù)據(jù)進行缺陷檢測前,加入了高斯噪聲和振動模糊。從異常檢測模型可視化結(jié)果來看,該模型對噪聲與振動模糊具有一定的魯棒性,為進一步實時觀察藥柱缺陷對燃面燃燒的影響奠定了基礎(chǔ)。
(3)本文方法與傳統(tǒng)超稀疏重建算法相比,只需兩個角度的投影數(shù)據(jù),即可重建三維結(jié)構(gòu),大大減少了工程成本。與現(xiàn)有的深度重建網(wǎng)絡(luò)相比,其無需缺陷區(qū)域的先驗信息,能夠在實時采樣情況下完成重建,使設(shè)計人員較準確地獲取發(fā)動機藥柱缺陷信息。
(4)本文方法為固體火箭發(fā)動機點火情況下藥柱缺陷的檢測提供了一種實現(xiàn)途徑,并進行了實驗驗證。下一步,圍繞實際點火實驗開展缺陷重建結(jié)果的精細化表征,以及缺陷對發(fā)動機燃面影響的分析研究。