吳 未，陳 華，戴 斌，張偉斌，駱吉洲，姜守旭，宗和厚

（1.中國工程物理研究院化工材料研究所，四川 綿陽621999；2.哈爾濱工業(yè)大學(xué)計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，黑龍江 哈爾濱150001）

1 引言

三氨基三硝基苯（TATB）基高聚物黏結(jié)炸藥（Polymer Bonded Explosive，PBX）是一種能量表現(xiàn)優(yōu)良的鈍感耐熱炸藥，主要由炸藥顆粒、黏結(jié)劑等組份構(gòu)成，在壓制成型后組份以及這些組份構(gòu)成的結(jié)構(gòu)是影響炸藥整體力學(xué)性能的重要因素［1］。通過對(duì)壓制成型后的TATB 基PBX 內(nèi)部顆粒邊界的精確識(shí)別，可進(jìn)一步開展裝藥顆粒間界面形變等問題的研究，建立基于真實(shí)PBX 材料的計(jì)算模型；在構(gòu)建PBX 宏觀本構(gòu)模型時(shí)，PBX 內(nèi)部顆粒的精確識(shí)別可以補(bǔ)充對(duì)顆粒和界面的描述，是研究炸藥非線性損傷、斷裂行為及本構(gòu)關(guān)系理論的重要基礎(chǔ)。當(dāng)前對(duì)于TATB 基PBX 顆粒邊界識(shí)別和提取主要采用計(jì)算機(jī)層析成像（Computed Tomography，CT）技術(shù)，但由于該技術(shù)對(duì)PBX 炸藥的檢測能力存在空間與密度分辨水平低、尺寸定量精度較差等問題，客觀上造成了CT 圖像視野較暗、對(duì)比度低、噪聲較多等情況［1］，進(jìn)而使現(xiàn)有的圖像分割技術(shù)無法對(duì)顆粒大小、位置、分布等關(guān)鍵數(shù)據(jù)進(jìn)行準(zhǔn)確的定量表征。

PBX 顆粒的邊界勾勒實(shí)質(zhì)是對(duì)PBX 炸藥的CT 圖像進(jìn)行分割。圖像分割一直是計(jì)算機(jī)視覺領(lǐng)域的研究熱點(diǎn)，但是針對(duì)PBX 等特殊材料CT 圖像的研究相對(duì)匱乏。基于深度學(xué)習(xí)的圖像分割算法利用各類卷積神經(jīng)網(wǎng)絡(luò)（如FCN［2］，U-Net［3］，Mask R-CNN［4］），通過大量經(jīng)過標(biāo)注的圖像來訓(xùn)練模型完成圖像分割的任務(wù)。而目前PBX 材料CT 圖像數(shù)量稀少且未經(jīng)標(biāo)注，難以滿足模型訓(xùn)練的需要。傳統(tǒng)的圖像分割算法主要包括基于閾值、邊緣和區(qū)域三類。OTSU［5］為代表的基于閾值的分割算法將各個(gè)像素的灰度值與閾值作比較來實(shí)現(xiàn)圖像分割。這類算法具有簡單、高效、應(yīng)用廣泛等特點(diǎn)，但它們僅考慮了灰度值而沒有考慮像素的空間分布，有對(duì)噪聲敏感的缺點(diǎn)。因此，在高噪聲的PBX炸藥CT圖像上直接應(yīng)用閾值分割算法難以準(zhǔn)確提取炸藥顆粒邊界。

基于邊緣檢測的分割算法其核心思想是利用微分算子檢測灰度值或者結(jié)構(gòu)突變的區(qū)域以識(shí)別目標(biāo)邊緣進(jìn)而實(shí)現(xiàn)圖像分割，有時(shí)也利用邊緣特征在同一相位具有較高頻率來識(shí)別目標(biāo)邊緣以完成分割任務(wù)［6-7］。這類算法或無法檢出連續(xù)完整的邊緣，或?qū)υ肼暜惓Ｃ舾校?］，在PBX 高噪聲高紋理的CT 圖像上難以取得良好的效果。分水嶺算法［8］為代表的基于區(qū)域的分割算法以尋找和定位分割區(qū)域?yàn)榛A(chǔ)。在低對(duì)比度高噪聲的PBX 炸藥CT 圖像上，直接利用經(jīng)典邊緣算子計(jì)算輪廓來得到種子點(diǎn)容易導(dǎo)致炸藥顆粒邊緣分界混亂的問題，難以準(zhǔn)確提取PBX 顆粒邊界。近年提出的基于超像素的圖像分割方法［9］將紋理、顏色和亮度等特征相似的像素視為具有視覺意義的不規(guī)則超像素。這種思想分別與閾值分割算法、邊緣檢測算法和分水嶺算法相結(jié) 合，產(chǎn)生了SLIC［10］、SEEDS［11］和Waterpixels［12］等圖像分割算法。這些算法運(yùn)用于PBX 的CT 圖像邊界勾勒時(shí)難以確保每個(gè)超像素恰好表示一個(gè)炸藥顆粒內(nèi)所有像素，使得所勾勒的顆粒邊界嚴(yán)重偏離真實(shí)情況。

基于上述分析，本研究提出的邊界勾勒算法，結(jié)合統(tǒng)計(jì)學(xué)數(shù)字特征和二值化方法的彈性膠囊算法以識(shí)別圖像中每個(gè)炸藥顆粒對(duì)應(yīng)的唯一種子點(diǎn)，利用選定的種子點(diǎn)調(diào)用分水嶺算法獲得TATB 基PBX 炸藥的單像素顆粒邊界，最后對(duì)單像素顆粒邊界進(jìn)行進(jìn)一步調(diào)整和優(yōu)化來獲得準(zhǔn)確清晰的邊界勾勒結(jié)果。

2 邊界勾勒算法

2.1 算法基本思想

TATB 基PBX 炸藥CT 圖像呈現(xiàn)顆粒內(nèi)部灰度值低，邊界部分灰度值高的特性，如果將灰度圖像轉(zhuǎn)變?yōu)槿S數(shù)字高程模型，顆粒內(nèi)部對(duì)應(yīng)“山谷”，邊界部分則對(duì)應(yīng)為“山脊”，與經(jīng)典分水嶺算法［8］的設(shè)計(jì)思路非常契合。但是經(jīng)典分水嶺算法流程中對(duì)于種子點(diǎn)的選取是直接利用Canny 算法［13］檢測邊緣，依據(jù)邊緣圖尋找和繪制輪廓，生成種子點(diǎn)輸入到分水嶺算法中。但由于TATB 基PBX 圖像質(zhì)量較差，用該方法尋找的種子點(diǎn)與實(shí)際炸藥顆粒對(duì)應(yīng)關(guān)系混亂，當(dāng)某個(gè)顆粒位置對(duì)應(yīng)多個(gè)種子點(diǎn)，則顆粒內(nèi)部將被分成多個(gè)小區(qū)域；如果未對(duì)應(yīng)任何種子點(diǎn)，則該顆粒區(qū)域?qū)?huì)與其他顆粒粘連。為了解決種子點(diǎn)選取不合理的問題，本文提出基于彈性膠囊的分水嶺算法［8］來實(shí)現(xiàn)TATB 基PBX 炸藥CT 圖像上顆粒邊界的自動(dòng)識(shí)別和勾勒，該算法的核心思想是將彈性膠囊作為工具，讓可以形變的“膠囊”在二值化邊界草圖中PBX 顆粒間游走，通過分析“膠囊”覆蓋下像素的類別確定該位置是否位于炸藥顆粒邊界以確定與實(shí)際炸藥顆粒一一對(duì)應(yīng)的子區(qū)域，再將子區(qū)域作為種子點(diǎn)輸入到分水嶺算法中完成PBX 顆粒的分割。最后依據(jù)彈性膠囊在圖像游走時(shí)留下的軌跡特征，適當(dāng)調(diào)整分水嶺算法輸出圖像邊界的寬度與灰度，得到最終較為清晰準(zhǔn)確的炸藥顆粒邊界圖，該算法基本流程如圖1 所示，主要包含興趣區(qū)域識(shí)別與統(tǒng)計(jì)二值化、彈性膠囊算法、邊界調(diào)整和種子點(diǎn)選取、邊界精細(xì)化4 個(gè)子步驟。

圖1 顆粒邊界勾勒算法流程圖Fig.1 Flowchart of the granule boundary extraction algorithm

2.2 興趣區(qū)域識(shí)別與統(tǒng)計(jì)二值化

興趣區(qū)域識(shí)別是整個(gè)邊界勾勒算法流程的第一步，旨在從原始圖像（如圖2a）中提取興趣區(qū)域，排除非主體部分對(duì)后續(xù)研究的干擾。本研究樣品的興趣區(qū)域?yàn)閳D2a 中央圓形內(nèi)部，邊緣退化［14］情況較輕，故直接采用Canny 算法［13］，通過高斯平滑和非極大值抑制消除邊緣噪聲和偽影干擾，識(shí)別到準(zhǔn)確的興趣區(qū)域邊界，再通過Hough 圓變換［15］計(jì)算區(qū)域圓心和半徑，得到最終興趣區(qū)域的準(zhǔn)確空間坐標(biāo)范圍，并消除坐標(biāo)范圍外非主體部分的干擾（如圖2b）。

統(tǒng)計(jì)二值化過程是在興趣區(qū)域圖像（如圖2b）的基礎(chǔ)上，利用圖像中炸藥顆粒內(nèi)部與邊界的灰度差異，盡可能將興趣區(qū)域中粘結(jié)劑邊界分離出來，形成一張二值化邊界草圖。該子過程首先計(jì)算興趣區(qū)域內(nèi)每個(gè)像素（i，j）其20×20 鄰域內(nèi)像素的灰度均值m和標(biāo)準(zhǔn)差σ；然后對(duì)每個(gè)像素點(diǎn)進(jìn)行判斷，若像素（i，j）的灰度值g滿足g＜m-0.1σ，則判定像素（i，j）位于炸藥顆粒內(nèi)部，設(shè)為白色，反之判定為邊界，設(shè)為黑色；經(jīng)過上述操作，得到二值化邊界草圖（如圖2c）。

二值化邊界草圖具有如下特征：（1）由于密度不均和噪聲的問題，草圖中代表顆粒內(nèi)部的白色像素并非連續(xù)分布而是呈離散團(tuán)狀點(diǎn)云，內(nèi)部有孔洞；（2）草圖存在邊界不夠平滑、有孤立白色碎點(diǎn)的情況，相鄰顆粒的團(tuán)狀點(diǎn)云存在粘連。為了讓每個(gè)顆粒有且僅有一個(gè)閉合子區(qū)域與其對(duì)應(yīng)，需要在后續(xù)的操作中填充孔洞，消除碎點(diǎn)與粘連問題。

圖2 興趣區(qū)域識(shí)別與統(tǒng)計(jì)二值化過程Fig.2 Procedures of regions of interest identification and statistical binarization

2.3 彈性膠囊算法

為了讓每個(gè)顆粒都有唯一子區(qū)域與之對(duì)應(yīng)，消除團(tuán)狀點(diǎn)云之間的粘連以及碎點(diǎn)是核心的問題，而彈性膠囊算法則是為了解決這一問題，保證分水嶺算法種子點(diǎn)選取的合理性而提出。現(xiàn)有的傳統(tǒng)圖像分割算法都只是從圖像本身出發(fā)，依據(jù)像素的灰度、分布進(jìn)行圖像分割，當(dāng)圖像質(zhì)量較低時(shí)，由于局部的噪聲，其他算法無法有效將粘連的PBX 炸藥小顆粒劃分開。彈性膠囊算法除了圖像本身的信息，還考慮了PBX 炸藥顆粒壓制后類球狀的形態(tài)，將這個(gè)信息作為常識(shí)融入到算法中，通過條狀“膠囊”的旋轉(zhuǎn)，捕捉顆粒邊界語義信息，較好地消除顆粒粘連和碎點(diǎn)的問題。其次，現(xiàn)有的傳統(tǒng)圖像處理算法對(duì)該類別PBX 炸藥CT 圖像的異常像素點(diǎn)（如噪聲點(diǎn)、材料雜質(zhì)部分）應(yīng)對(duì)的能力不足，灰度或邊緣特征在受到異常像素點(diǎn)影響后都會(huì)發(fā)生變化，從而讓結(jié)果偏離真實(shí)的顆粒分布。而彈性膠囊算法利用膠囊的移動(dòng)，在分析膠囊覆蓋下的顆粒情況時(shí)，特別考慮了異常像素點(diǎn)的問題，對(duì)異常點(diǎn)有足夠的寬容度，所以在低質(zhì)量PBX 炸藥CT 圖像中，該算法能有效選取合理種子點(diǎn)，保證了最終輸出結(jié)果的準(zhǔn)確度。

每個(gè)彈性膠囊是寬為2W+1 長為L的矩形塊，矩形塊中每個(gè)像素記錄其與膠囊最左側(cè)中點(diǎn)的偏移，可以用二維點(diǎn)集C=｛（x，y）||y|≤W，0≤x＜L且x，y∈Z｝來保存。彈性膠囊的旋轉(zhuǎn)是圍繞其最左側(cè)中點(diǎn)進(jìn)行的，故稱該點(diǎn)為“軸點(diǎn)”（如圖3 中紅色像素，對(duì)應(yīng)點(diǎn)集C 中（0，0）點(diǎn)）。將膠囊繞軸點(diǎn)逆時(shí)針旋轉(zhuǎn)角度α（0°≤α＜180°）得到的點(diǎn)集C'=｛（［xcosα-ysinα］，［xsinα+ycosα］）|（x，y）∈C｝，圖3 中綠色膠囊就是橙色膠囊逆時(shí)針旋轉(zhuǎn)45°得到的。本文根據(jù)樣品特征，選取彈性膠囊寬度范圍｛5，7，9｝，長度最大為47，最小為6；通過在0°與180°之間均勻選取16 個(gè)角度，對(duì)初始膠囊進(jìn)行旋轉(zhuǎn)，每個(gè)初始膠囊將對(duì)應(yīng)16 種不同方向的新膠囊。

圖3 膠囊及其旋轉(zhuǎn)示意圖Fig.3 Schematic diagram of a elastic capsule and its rotation

在膠囊生成后，按面積從大到小的順序選取膠囊，二值化草圖中的每個(gè)像素都要與膠囊經(jīng)過一輪比較，當(dāng)像素（i，j）與膠囊進(jìn)行比較時(shí)，就將膠囊“軸點(diǎn)”與像素（i，j）重合，再對(duì)膠囊覆蓋下的像素進(jìn)一步分析以確定像素（i，j）的類別。具體像素分析過程如下：（1）統(tǒng)計(jì)膠囊覆蓋下邊界點(diǎn)像素?cái)?shù)nResult；（2）統(tǒng)計(jì)膠囊覆蓋下顆粒內(nèi)部點(diǎn)像素?cái)?shù)nPFrag；（3）統(tǒng)計(jì)膠囊覆蓋下顆粒內(nèi)部點(diǎn)構(gòu)成的連通分支數(shù)nFrag；（4）若nResult大于等于20，nPFrag小于等于膠囊面積的1/20，nFrag小于等于5，三個(gè)條件若同時(shí)滿足則當(dāng)前像素為邊界的概率很大，置像素（i，j）的灰度值為128+nResult，同時(shí)這三個(gè)條件的設(shè)定也能接受一定量異常像素點(diǎn)的存在，使算法對(duì)噪聲不敏感；若三個(gè)條件不同時(shí)滿足，則換用更小的膠囊進(jìn)行比較；若所有膠囊比較后仍然無法滿足條件，則認(rèn)為該點(diǎn)無法通過，置灰度值為0。經(jīng)過上述彈性膠囊算法，將會(huì)得到一張彈性膠囊邊界圖（如圖4a）。

2.4 邊界調(diào)整和種子點(diǎn)選取

彈性膠囊算法得到的炸藥顆粒粗略邊界仍含有一些噪聲，為了消除這些噪聲對(duì)后續(xù)操作的影響，需要進(jìn)行邊界調(diào)整操作。先計(jì)算膠囊軌跡像素（如圖4a 中白色區(qū)域）的灰度均值m和標(biāo)準(zhǔn)差σ，然后將灰度值小于m-σ/2 的膠囊軌跡像素置為黑色以縮減邊界，所得結(jié)果如圖4b。

圖4 彈性膠囊算法和分水嶺算法步驟Fig.4 Procedures of the elastic capsule algorithm and watershed algorithm

調(diào)整后邊界變窄，連通性變?nèi)跚疫吔鐑?nèi)部可能存在孤立噪點(diǎn)。利用邊界形態(tài)學(xué)變換在圖4b 所示的圖像上進(jìn)行兩次開操作和一次膨脹操作，在圖中形成多個(gè)較為可靠的黑色子區(qū)域（如圖4c 黑色部分）。其中黑色像素表示炸藥顆粒內(nèi)部，每個(gè)炸藥顆粒近似與圖4c 中黑色子區(qū)域一一對(duì)應(yīng)。本研究采用Canny 算法［13］提取黑色子區(qū)域的邊緣，再依據(jù)提取到的邊緣繪制每個(gè)子區(qū)域的輪廓（所得結(jié)果如圖4d），將每個(gè)輪廓內(nèi)部所包含的獨(dú)立區(qū)域作為種子點(diǎn)，調(diào)用分水嶺算法［8］產(chǎn)生單像素寬度的顆粒邊界勾勒圖像（如圖4e 所示）。

2.5 邊界精細(xì)化

邊界精細(xì)化旨在根據(jù)彈性膠囊算法的輸出（如圖4a），計(jì)算單像素顆粒邊界（如圖4e）中的每個(gè)邊界像素的粗細(xì)權(quán)值，并根據(jù)權(quán)值向外擴(kuò)張不同的寬度以獲得更準(zhǔn)確、真實(shí)的邊界勾勒結(jié)果。首先在圖4a 中標(biāo)記單像素邊界位置；然后循環(huán)地用顆粒內(nèi)部黑色像素腐蝕圖4a 中彈性膠囊的邊界，當(dāng)腐蝕操作到達(dá)標(biāo)記的單像素顆粒邊界點(diǎn)（i，j）時(shí)將所用循環(huán)的輪數(shù)roundij作為圖4e 中白色邊界點(diǎn)（i，j）的權(quán)值，如果像素（i，j）兩側(cè)腐蝕操作到達(dá)時(shí)循環(huán)輪數(shù)不同則保留較大的值；接下來根據(jù)所有單像素邊界點(diǎn)權(quán)值roundij將所有單像素邊界點(diǎn)的權(quán)值線性映射到［1，3］的集合中得到最終權(quán)值weightij，再在圖4e中從每個(gè)白色邊界點(diǎn)依據(jù)對(duì)應(yīng)的權(quán)值向外擴(kuò)張1+2weightij個(gè)像素，最后再將圖像興趣區(qū)域中白色像素的灰度設(shè)置成與原始圖像（圖2a）對(duì)應(yīng)位置一致，產(chǎn)生最終的邊界勾勒結(jié)果（如圖4f）。

3 實(shí)驗(yàn)部分

3.1 實(shí)驗(yàn)樣品與CT 檢測

實(shí)驗(yàn)樣品：TATB 基PBX，樣品為等靜壓成型的PBX 試件經(jīng)過機(jī)械加工后制備成Φ10 mm×10 mm 的圓柱形試件。試驗(yàn)樣品的CT 掃描成像實(shí)驗(yàn)在nano-Voxel 設(shè)備（Sanying Precision Engineering Research Center，Tianjin，China）上展開，如圖5 所示。實(shí)驗(yàn)所采用的參數(shù)為電壓100 kV，電流100 μA，單幅投影圖采樣曝光時(shí)間0.3 ms，樣品旋轉(zhuǎn)360 度，投影圖采樣數(shù)為1080，檢測空間分辨率為15 μm/pixel。

圖5 CT 原理示意圖Fig.5 Principle of X-ray CT

3.2 不同算法的對(duì)比實(shí)驗(yàn)

為了分析本研究提出的算法對(duì)比現(xiàn)有傳統(tǒng)算法的表現(xiàn)，我們選取使用廣泛且具有代表性的5 種圖像分割算法進(jìn)行比較。圖6a 的Canny 算法和圖6b 的相位一致性算法是典型基于邊緣的圖像分割算法；圖6c 的分水嶺算法是用途最為廣泛的基于區(qū)域算法；圖6d 的SLIC 算法和圖6e 的SEEDS 算法是具有代表性的基于超像素算法；上述算法靈活性強(qiáng)、用途廣泛，且分屬傳統(tǒng)圖像分割算法的不同子類別，具有典型性，所以選取了這幾個(gè)算法作為算法表現(xiàn)的參照。本研究將以上5 種算法與本文提出的算法在PBX 炸藥樣品圖像中運(yùn)行，實(shí)驗(yàn)結(jié)果表明：（1）Canny 算法（如圖6a）邊界識(shí)別不完整且噪聲嚴(yán)重，相位一致性算法（如圖6b）幾乎無法識(shí)別炸藥顆粒邊界；（2）直接調(diào)用分水嶺算法時(shí)由于種子點(diǎn)選擇不合理，出現(xiàn)了嚴(yán)重的顆粒粘連情況（如圖6c）；（3）SLIC 算 法（如 圖6d）和SEEDS 算 法（如圖6e）等超像素圖像分割算法的目的是為了降低圖像維度以及剔除異常像素點(diǎn)，并不適合直接應(yīng)用于PBX圖像顆粒實(shí)例的分割，所以分割結(jié)果與顆粒邊界差距較大；（4）本文提出的算法由于確保了可以恰當(dāng)?shù)貜念w粒內(nèi)部選取分水嶺算法的種子點(diǎn)，故邊界勾勒的效果明顯優(yōu)于其他算法。

3.3 準(zhǔn)確性分析

為了定性分析本文提出的算法在TATB 基PBX 炸藥CT 圖像邊界勾勒的表現(xiàn)，我們在樣品圖像上運(yùn)行我們的算法，得到的結(jié)果如圖7，利用圖7a 調(diào)整得到的算法參數(shù)，在該批次其他圖像中均有良好的效果（如圖7b～7d）。實(shí)驗(yàn)結(jié)果表明：本研究提出的算法在不同圖像上得到了連續(xù)且與實(shí)際貼合程度較高的邊界勾勒結(jié)果。同時(shí)，依據(jù)某一圖像情況設(shè)定的算法參數(shù)，依然可以使得算法在同批次的其他樣品圖像中有良好的效果，同批次圖像數(shù)據(jù)執(zhí)行算法時(shí)，無需反復(fù)針對(duì)圖像調(diào)整參數(shù)，具有較高的自動(dòng)化程度。

圖6 各算法在圖2a 所示的原始圖像上的輸出結(jié)果對(duì)比Fig.6 Comparison of the outputs of different algorithms applied to original image given in Fig.2a

圖7 本文提出的算法在不同輸入圖像上的結(jié)果對(duì)比Fig.7 Comparison of the outputs of the proposed algorithm applied to different original images

4 結(jié)論

（1）通過CT 掃描得到的TATB 基PBX 炸藥CT 圖像較為模糊，難以獲取準(zhǔn)確邊界輪廓，進(jìn)而造成分水嶺算法種子點(diǎn)選取精度較差，本文提出的基于彈性膠囊的處理流程，可以輔助選取準(zhǔn)確的種子點(diǎn)，實(shí)驗(yàn)證明，利用彈性膠囊算法選取種子點(diǎn)后，分水嶺算法的邊界提取精度有明顯的提升，其邊界提取準(zhǔn)確性要明顯好于其他具有代表性的5 種傳統(tǒng)圖像分割算法。

（2）針對(duì)TATB 基PBX 炸 藥CT 掃描圖像 質(zhì)量較差，通過常用的圖像處理算法無法有效提取炸藥顆粒邊界特征的問題，提出了基于統(tǒng)計(jì)數(shù)字特征和彈性膠囊的分水嶺算法，從PBX 炸藥CT 掃描原始圖像上可以準(zhǔn)確、完整地提取炸藥顆粒邊界，算法自動(dòng)化程度較高，能縮減人力成本，有一定的實(shí)用價(jià)值。

（3）基于彈性膠囊的分水嶺算法提高了TATB 基PBX 內(nèi)部細(xì)觀結(jié)構(gòu)的表征能力，為進(jìn)一步研究TATB基PBX 宏觀力學(xué)性能提供了研究基礎(chǔ)。