曹旭寧 陳思澤 俞 杰 張連鑫 張早娣 李桃生

1（中國科學院合肥物質科學研究院核能安全技術研究所 合肥 230031）

2（中國科學技術大學 合肥 230027）

中子照相技術是一種重要的無損檢測技術［1］，由于其穿透能力強、抗γ射線干擾以及對輕核素敏感等特點，彌補了傳統(tǒng)X射線成像的不足，受到廣泛關注［2］。然而由于中子相比X射線成像效率低，且中子源使用成本高昂，中子照相經(jīng)常面臨成像對比度不足的問題，限制了其在小型化、動態(tài)拍攝及三維CT成像等方向的發(fā)展和應用［3］。對比度圖像增強算法在軟件環(huán)境下直接對已成像圖像進行處理，改善原圖灰度分布，是一種解決中子照相對比度不足的有效技術方案，降低了圖像質量對中子源和成像設備的依賴，成為中子照相技術的熱點研究領域［4-5］。

目前，國內(nèi)外學者在對比度增強算法領域已進行了大量的研究，其中直方圖均衡及其改進算法和Retinex算法及其改進算法應用較為廣泛［6］。典型代表有直方圖均衡化（Histogram Equalization，HE）、自適應直方圖均衡化（Adaptive Histogram Equalization，AHE）［7］、Retinex算 法［8］、單尺度Retinex算法（Single-Scale Retinex，SSR）、多尺度Retinex算法（Multi-Scale Retinex，MSR）［9］等。然而直方圖均衡及其改進算法往往需要根據(jù)不同場景人為調(diào)整參數(shù)，并且參數(shù)的選取對圖像處理效果影響較大，導致對比度增強效果穩(wěn)定性較差，增強結果不易控制［10］。而Retinex及其改進算法大多為全局增強處理。中子照相成像系統(tǒng)由于存在中子輻射干擾，成像結果往往嵌入大量低頻背景和噪聲信號中［11］，全局增強易將背景噪聲同時放大，導致圖像增強效果不佳，甚至會引發(fā)圖像質量倒退的現(xiàn)象。因此局部自適應對比度增強算法更適合中子圖像對比度增強，可以避免大量低頻背景噪聲對提升過程的干擾，此外為了提升處理效果的穩(wěn)定性需要采用優(yōu)化算法求解參數(shù)值來代替人工設定的參數(shù)值。

針對中子圖像特點和對比度優(yōu)化需求，本文選取自適應對比度增強算法（Adaptive Contrast Enhancement，ACE）對中子圖像進行對比度優(yōu)化，這種算法可以有效緩解低頻背景的干擾。另外為了解決現(xiàn)有對比度增強算法在提升過程中的不穩(wěn)定性問題，本文將群智能技術與自適應對比度增強技術相結合，采用基于粒子群優(yōu)化的自適應對比度增強算法，找到全局最優(yōu)增益因子系數(shù)，利用全局最優(yōu)增益因子系數(shù)對中子圖像進行局部自適應對比度增強。本文算法可以還原圖像細節(jié)，提高圖像對比度，有效提高了中子照相圖像質量，且通過粒子群尋優(yōu)算法根據(jù)圖像本身特點設置最優(yōu)參數(shù)值，避免了經(jīng)驗值參數(shù)給處理效果帶來的影響，增強了算法處理效果的穩(wěn)定性。

1 相關算法

1.1 自適應對比度增強算法

自適應對比度增強（Adaptive Contrast Enhancement，ACE）算法［12］通過結合圖像自身特點，對本身對比度很強的圖像區(qū)域進行較小比例增強，對本身對比度弱的圖像區(qū)域較大幅度增強。算法基本原理是將圖像分成兩部分：低頻部分和高頻部分，低頻部分可通過均值濾波求解得到，高頻部分可通過原圖與低頻部分做差得到。利用增益因子實現(xiàn)對圖像高頻部分的增強和放大，從而達到自適應調(diào)節(jié)的目的，整體過程如下：

式中：M(i，j)為以像素點(i，j)為中心的局部均值；I(i，j)為像素點(i，j)對應的初始像素值；G為增益系數(shù)；I0(i，j)為結果圖像中像素點(i，j)對應的像素值。

由于標準差表示的是圖像像素均勻性，局部區(qū)域標準差越大，則表明其像素值越不均勻，對比度越強，反之標準差越小，其像素值越均勻，對比度越弱。所以將標準差作為增益因子G的動態(tài)取值計算參數(shù)之一，如式（2）：

式中：M為全局均值；σ(i，j)為以像素點(i，j)為中心的局部標準差；α為增益因子系數(shù)，選?。?，1）任意固定值；ε為極小正常數(shù)項，用于避免因出現(xiàn)局部標準差為零而導致的公式計算錯誤（本文ε=10-7）。

由式（2）可以看出，自適應對比度增強算法中增益因子的計算需要一個增益因子系數(shù)α調(diào)節(jié)高頻部分的增益程度，通常增益因子系數(shù)α往往選取一個固定經(jīng)驗值，但是這樣就會給增強過程帶來不穩(wěn)定性，因此對于自適應對比度增強算法有必要選取一個最優(yōu)增益因子系數(shù)α，以便能夠保證最佳的對比度提升效果。

1.2 粒子群優(yōu)化算法

粒子群優(yōu)化算法（Particle Swarm Optimization，PSO）［13］由Kennedy和Eberhart提出，模擬鳥群在覓食過程中的集群行為而設計。相比于遺傳算法，它省去了復雜的“交叉”、“變異”等過程，PSO算法從隨機解出發(fā)，通過不斷迭代的過程在解空間中搜索最優(yōu)解。其核心思想延續(xù)了智能優(yōu)化算法，即通過群體中個體之間的協(xié)作和信息共享來尋找最優(yōu)解。

PSO算法［14］通過設計一種無質量粒子模擬覓食鳥群中的鳥，這種無質量粒子具有兩個屬性：速度和位置。速度代表粒子移動的快慢，位置代表粒子移動的方向。每個無質量粒子在搜索空間中單獨尋找最優(yōu)解，記為當前個體最優(yōu)解（pbest），和整個粒子群中的其他粒子共享這個個體最優(yōu)解，整個粒子群中最優(yōu)的個體值視為當前全局最優(yōu)解（gbest）。粒子群中的所有粒子根據(jù)當前找到的個體最優(yōu)解和整個粒子群共享的當前全局最優(yōu)解不斷迭代更新自己的速度和位置，最終所有粒子都收斂于最優(yōu)解周圍，從而實現(xiàn)最優(yōu)化。

由于PSO算法簡潔，易于實現(xiàn)，且在非線性連續(xù)、組合和混合整數(shù)非線性等優(yōu)化問題上具有較好性能，目前已經(jīng)在函數(shù)優(yōu)化、神經(jīng)網(wǎng)絡訓練、模糊系統(tǒng)控制以及其他遺傳算法等領域有著非常廣泛的應用［15］?？紤]到PSO算法在函數(shù)優(yōu)化方面的優(yōu)勢，本文選取該算法用于自適應對比度增強算法的增益因子系數(shù)α的尋優(yōu)。

2 改進的自適應對比度增強算法

本文算法主要為按照自適應對比度增強算法（ACE）思想將圖像分成高頻部分和低頻部分，對于代表細節(jié)的高頻部分進行增益計算。以圖像信息熵和標準差作為目標函數(shù)，通過粒子群優(yōu)化算法（PSO）迭代優(yōu)化找到圖像的最優(yōu)增益因子系數(shù)，將所求解的全局最優(yōu)增益因子系數(shù)代入ACE算法中進行局部自適應對比度增強，克服了ACE算法采用固定系數(shù)值增益所導致的處理效果不穩(wěn)定，對中子圖像進行了自適應的對比度增強，減小了大量暗區(qū)對中子檢測的干擾，還原圖像細節(jié)，提高圖像質量。本文算法基本流程如圖1所示。

圖1 算法基本流程圖Fig.1 The flow chart of algorithm

2.1 增益因子系數(shù)尋優(yōu)過程

為解決固定系數(shù)帶來的不穩(wěn)定性，本文采用PSO算法尋找增益因子系數(shù)α的最優(yōu)解，將增益因子系數(shù)看作粒子，首先初始化一個固定大小的粒子群（本算法中構造了一個含有24個粒子的粒子群），并對粒子群中的每一個粒子進行速度和位置的隨機初始化，粒子在每一次迭代中都會根據(jù)目標函數(shù)評估出每個粒子的適應度值（fitness），若fitness是粒子的當前個體最優(yōu)解，則將它存儲為pbest，若fitness是當前粒子群中的全局最優(yōu)解，則將它存儲為gbest，然后根據(jù)pbest和gbest對各粒子進行位置和速度的更新，更新公式如式（3）和式（4）所示：

式中：i表示當前粒子；t表示迭代次數(shù)；表示第i個粒子在第t次迭代時的速度；xti表示第i個粒子在第t次迭代時的位置。C1、C2為加速常數(shù)，前者為每個粒子的個體學習因子，后者為每個粒子的社會學習因子，通常令C1=C2=2，可以在不易陷入局部最優(yōu)的前提下加快收斂。表示第i個粒子在第t次迭代時的個體最優(yōu)解，表示第i個粒子在第t次迭代時的全局最優(yōu)解。

式（3）中ω為慣性因子，用來衡量全局尋優(yōu)性能和局部尋優(yōu)性能。當ω較大時，全局尋優(yōu)能力較強，局部尋優(yōu)能力較弱。考慮到在迭代早期較大的ω有助于粒子群的全局快速收斂，而迭代后期粒子群已經(jīng)較集中于最優(yōu)范圍內(nèi)，此時應通過較小的ω增強局部搜索能力，在最優(yōu)范圍內(nèi)更加精細地搜索最優(yōu)解。因此采用線性遞減法自適應地選擇ω，其計算公式如式（5）所示：

式中：ωmax和ωmin分別為慣性權重最大值和最小值；tmax為最大迭代次數(shù)。

通過不斷地迭代和更新，當滿足迭代終止條件（迭代次數(shù)達到最大迭代次數(shù)tmax）時各粒子收斂于最優(yōu)值附近，即增益因子系數(shù)獲得了最優(yōu)解。

2.2 選擇評價函數(shù)

在對增益因子系數(shù)進行粒子群優(yōu)化的過程中，評價函數(shù)的選擇是至關重要的，每個粒子根據(jù)評價函數(shù)評估所得到的適應值fitness是評價中子圖像增強效果的重要標準。中子圖像由于其中子源和成像系統(tǒng)的限制，導致圖像分辨率較差，存在大量暗區(qū)，細節(jié)丟失嚴重。因此要求對比度增強后的中子圖像應具有信息量大、對比度高等特點［16］。在本文算法中，將灰度標準方差和圖像信息熵融入評價函數(shù)中，作為每個粒子的尋優(yōu)目標函數(shù)，其評價函數(shù)如式（6）所示：

式中：α1和α2分別表示信息熵和標準差在評估中所占比例，此處α1和α2均等于1/2，表示H和Std占比相同。

H代表圖像信息熵，在圖像中，圖像的熵代表整張圖像的不確定性，可以用圖像熵來衡量圖像所包含信息量的大小，也可以表示圖像灰度值分布的聚集特征，熵越大則圖像中的灰度值在一定范圍內(nèi)越均衡，表現(xiàn)在視覺上就是圖像對比度越大，圖像質量越高，計算過程如式（7）所示：

其中：pn表示圖像中灰度值n所占比例。

Std代表圖像標準差，標準差表示的是圖像的像素值的均勻性，可以認為標準差越大的局部區(qū)域，其像素值越不均勻，對比度越強；反之，標準差越小的局部區(qū)域，其像素值越均勻，對比度越弱。

通過從信息熵和標準差兩個維度對粒子進行評估，將多目標函數(shù)線性排列來解決多目標優(yōu)化的問題，從而使選取的增益因子系數(shù)能夠最大化保證圖像信息量和像素均勻性，滿足增強中子圖像細節(jié)，提高中子圖像對比度的要求。

2.3 改進算法過程

改進的自適應對比度增強優(yōu)化算法具體過程為：

1）輸入原始中子圖像。

2）將增益因子系數(shù)α看作粒子，初始化增益因子系數(shù)的種群和參數(shù)。包括種群個數(shù)N、粒子位置x和速度v；慣性權重最大值最小值ωmax和ωmin；學習因子C1和C2、最大迭代次數(shù)tmax等。

3）評價每個增益因子系數(shù)α。分別計算熵值H、計算灰度標準方差Std，得到每個增益因子系數(shù)的適應度值（fitness）。

4）更新增益因子系數(shù)最優(yōu)值。如果更新后的適應度值大于原來個體最優(yōu)值pbest，則將該增益因子系數(shù)位置作為當前個體最優(yōu)值pbest，如果更新后增益因子系數(shù)適應度值大于原來全局最優(yōu)值gbest，則將該增益因子系數(shù)位置作為當前全局最值優(yōu)gbest。

5）判斷是否達到迭代停止條件（即達到迭代次數(shù)tmax），若未達到，則更新每個增益因子系數(shù)的速度與位置，跳轉到3）。若滿足迭代停止條件，退出循環(huán)。

6）使用迭代中存儲下來的最優(yōu)增益因子系數(shù)α（gbest）對圖像進行ACE算法增強處理。

7）輸出經(jīng)算法處理之后的中子圖像。

3 結果與分析

為了驗證新算法的實際效果，本文對多幅中子照相成像圖片進行優(yōu)化測試。其中圖2和圖3分別為來自德國慕尼黑工業(yè)大學FRM-II實驗堆的鐘表樣品和渦輪葉片樣品中子照相圖像，圖4為基于HINEG加速器中子源獲取的線對樣品和階梯樣品的快中子照相圖像，圖5為中國原子能研究院先進堆中子照相裝置獲取的汽車火花塞中子照相圖像。將本文算法優(yōu)化結果分別與原圖、MSR算法和傳統(tǒng)ACE算法優(yōu)化結果進行效果對比，并計算給出了信息熵、灰度平均梯度和圖像對比度三項客觀評價參數(shù)結果。

3.1 效果對比

不同算法對于鐘模型中子圖像的增強效果如圖2所示?？梢钥闯鲈紙D像存在大量暗區(qū)，樣品細節(jié)模糊不清。圖2（b）存在由于圖像整體同尺度提升而導致的對比度提升不均勻的現(xiàn)象。圖2（c）雖然細節(jié)處較為清晰但存在嚴重的噪聲放大現(xiàn)象。采用本文算法對原始圖像數(shù)據(jù)處理之后的結果，在提高對比度的同時齒輪細節(jié)處也得到了很好的保持，具有局部自適應性，對比度均勻。迭代優(yōu)化過程中，最優(yōu)增益因子系數(shù)α對應的目標函數(shù)值隨迭代次數(shù)更新的曲線圖（見下文圖6（a）），函數(shù)值在大約第25次迭代時收斂于最優(yōu)值附近（α=0.000 102 419，fitness=2.255 52）。

圖2 鐘表樣品中子成像原始圖(a)、MSR算法處理圖(b)、ACE算法處理圖(c)和本文算法處理圖(d)Fig.2 Neutron radiography of clock sample's original image(a),image processed by MSR algorithm(b),image processed by ACE algorithm(c)and image processed by improved adaptive contrast enhancement algorithm(d)

渦輪葉片中子圖像用不同算法增強后效果如圖3所示，可以看出原始圖像中冷卻通道和邊緣小孔處較暗，圖3（b）對比度略有提高但圖像對比度不均勻，圖3（c）對比度提升效果不明顯。圖3（d）表明本文算法對渦輪葉片中子圖像進行處理后，更容易看清葉片的邊緣小孔和冷卻通道，對比度得到了較好的增強。迭代曲線（見下文圖6（b））表明其函數(shù)值在大約第70次迭代時函數(shù)值收斂于最優(yōu)值附近（α=0.000 509 162，fitness=2.252 34）。

圖3 渦輪葉片中子成像原始圖(a)、MSR算法處理圖(b)、ACE算法處理圖(c)和本文算法處理圖(d)Fig.3 Neutron radiography of turbine plate's original image(a),image processed by MSR algorithm(b),image processed by ACE algorithm(c)and image processed by improved adaptive contrast enhancement algorithm(d)

圖4給出了不同算法在線對、階梯樣品中子圖像上的效果對比。原始圖像中樣品信息嵌入在大量低頻背景信號中，圖像整體較暗，對比度差，不利于對樣品的分析。圖4（b）所示MSR算法使原始圖像對比度小幅提升，但由于是全局增強，沒有突出樣品細節(jié)。圖4（c）雖在細節(jié)處提高了樣品分辨率，但整體對比度仍不足。本文算法對樣品細節(jié)和圖像整體的對比度均有明顯提升，肉眼更容易識別線對樣品中不同間隔的明暗線對信息以及階梯樣品的灰度漸變信息。迭代曲線（見下文圖6（c）表明，大約在第60次迭代時函數(shù)值收斂于最優(yōu)值附近（α=0.072 763，fitness=2.263 52）。

圖4 線對、階梯樣品中子成像原始圖(a)、MSR算法處理圖(b)、ACE算法處理圖(c)和本文算法處理圖(d)Fig.4 Neutron radiography of lines and ladder samples'original image(a),image processed by MSR algorithm(b),image processed by ACE algorithm(c)and image processed by improved adaptive contrast enhancement algorithm(d)

圖5為汽車火花塞中子圖像通過不同算法增強之后的效果對比圖，可以看出原圖中螺旋條紋細節(jié)不清晰，圖5（b）所示MSR算法雖提高了火花塞中子圖像整體對比度，但細節(jié)部分并沒有得到突出顯示。相比于圖5（b），圖5（c）細節(jié)處清晰度略有提升。本文算法提高了火花塞中子圖像的整體和局部對比度，尤其是火花塞底部細節(jié)得到了較大程度還原。整個尋優(yōu)過程的迭代曲線，大約在90次迭代時函數(shù)值收斂于最優(yōu)值附近（α=0.846 859，fitness=2.236 8）。圖2～圖5樣品迭代尋優(yōu)曲線結果如圖6。

圖5 汽車火花塞中子成像原始圖(a)、MSR算法處理圖(b)、ACE算法處理圖(c)和本文算法處理圖(d)Fig.5 Neutron radiography of sparkplug's original image(a),image processed by MSR algorithm(b),image processed by ACE algorithm(c)and image processed by improved adaptive contrast enhancement algorithm(d)

圖6 最優(yōu)目標函數(shù)值隨迭代次數(shù)更新曲線圖(a)圖2，(b)圖3，(c)圖4，(d)圖5Fig.6 The optimal objective function value updates with the number of iterations(a)Fig.2,(b)Fig.3,(c)Fig.4,(d)Fig.5

3.2 評價參數(shù)對比

3.2.1 圖像信息熵

由前文可知，圖像信息熵可以衡量圖像信息量的多少，因此本文根據(jù)式（7）計算得到圖像信息熵數(shù)值，用來評價測試圖像處理效果［17］。表1為測試圖像2、圖像3、圖像4和圖像5的熵值比較，可以看出本文算法對于4幅圖像的信息熵值均有明顯提升，表明本文算法可以還原圖像細節(jié)，提高圖像質量。

表1 測試圖2～圖5的熵值比較Table 1 Comparison of Fig.2～Fig.5's entropy value

3.2.2 灰度平均梯度

灰度平均梯度(Grayscale Mean Gradient,GMG)［18］指圖像的邊緣細節(jié)或影線兩側附近有明顯差異的灰度值，可以反映圖像細微處反差變化的速率，即灰度變化率。這種灰度變化率的大小可用來表征圖像的相對清晰度。平均梯度越大，圖像層次越多，圖像細節(jié)越清晰。計算公式如式（8）所示：

由表2可以看出，本文算法對于4幅圖像的灰度平均梯度均有明顯提升，表明經(jīng)過本文算法處理后圖像清晰度提高。

表2 測試圖2～圖5的GMG值比較Table 2 Comparison of Fig.2～Fig.5's GMG value

3.2.3 圖像對比度

對比度計算［19］是衡量圖像對比度性能的最直觀和最基本的重要指標，對比度數(shù)值越高，說明圖像對比度越好，展示的細節(jié)就越豐富，其計算公式如式（9）所示：

式中：δ(i，j)為濾波窗口內(nèi)相鄰像素之間的灰度差值（本文中求取對比度的濾波窗口大小取3×3），Pδ(i，j)為相鄰像素間灰度差為δ的像素分布概率。計算結果如表3所示，通過對多幅圖像、不同算法的對比度數(shù)值對比，證明了本文算法在對比度提升方面的有效性。

表3 測試圖2～圖5的對比度比較Table 3 Comparison of Fig.2～Fig.5's contrast value

3.3 自適應對比度增強算法中窗口選值分析

自適應對比度增強算法中窗口大小的取值可能會對圖像增強效果帶來影響，窗口內(nèi)每個像素值與其均值相似，式（1）中I(i，j)-M(i，j)很小，因此對圖像增益效果不是很好，而當窗口較大時，在較大圖像塊內(nèi)式（1）中M(i，j)相同，圖像塊內(nèi)得到的新像素值較為接近，因此邊緣和細節(jié)部分增強效果略差，易發(fā)生局部模糊。為證明窗口大小設置為7×7的有效性，本文對圖2～圖5在濾波窗口尺寸大小分別設置為3、5、7、9、11時運用本文算法增強后圖像的信息熵值、GMG值和對比度值進行了對比分析，結果分別如圖7的（a）、（b）、（c）所示，可知當窗口大小取7×7時，能夠更好地兼顧整體和局部的處理效果，其各項評價指標計算總體優(yōu)于其他窗口大小對應數(shù)值。因此本文選取7×7作為自適應對比度增強算法局部處理過程中的窗口尺寸。

圖7 4幅圖像取不同窗口大小時對應的信息熵值(a)、GMG值(b)和對比度值(c)Fig.7 The entropy value(a),GMG value(b)and contrast value(c)corresponding to different window sizes of the four images

4 結語

本文提出了一種基于粒子群優(yōu)化的自適應對比度增強算法，并將該算法運用到中子圖像上進行分析。通過圖像處理之后的效果和指標的對比，可以得出以下結論：

1）新算法采用粒子群優(yōu)化算法獲取最優(yōu)增益因子系數(shù)α，從圖像本身獲取算法參數(shù)，再進行局部自適應的圖像對比度優(yōu)化，圖像處理具有較好的穩(wěn)定性。

2）與傳統(tǒng)的MSR算法、ACE算法相比，新算法在提高中子照相圖像對比度方面具有一定優(yōu)勢，其處理圖像的信息熵、灰度平均梯度和對比度參數(shù)均更為優(yōu)異。

基于上述研究結果，新算法在中子照相圖像處理領域具有明確的應用前景，并已開始廣泛應用于本單位基于加速器中子源的中子照相圖像優(yōu)化。由于該算法采用粒子群優(yōu)化算法獲取最優(yōu)增益因子系數(shù)α的過程需要多次迭代耗時較長，算法目前還無法滿足中子CT或動態(tài)成像等高速多幀圖像實時處理需求，下一代的高速對比度優(yōu)化算法目前正在開發(fā)中。

致謝感謝中國原子能科學研究院物理所中子散射團隊及德國慕尼黑工業(yè)大學Burkhard Schillinger教授提供了相關中子照相圖像用于本文算法測試。