張祥志 許子健 甄香君 王勇 郭智 嚴(yán)睿 邰仁忠

(中國科學(xué)院上海應(yīng)用物理研究所 上海 201204)

X射線顯微成像技術(shù)的分辨率高、穿透性好，可觀察厚物質(zhì)的內(nèi)部三維結(jié)構(gòu)，在生物醫(yī)學(xué)和納米材料領(lǐng)域有廣泛應(yīng)用[1]。上海光源軟X射線譜學(xué)顯微光束線站(BL08U)是我國第一條基于第三代同步輻射光源的軟X射線波段的光束線站。它用掃描透射X射線顯微術(shù)(STXM)將30 nm的高空間分辨能力和高能量分辨的近邊吸收精細(xì)結(jié)構(gòu)譜學(xué)(NEXAFS)結(jié)合在一起，有元素識別和化學(xué)態(tài)分析的能力。實驗中，掃描樣品時對樣品位置的移動能達(dá)到相應(yīng)的精度，才能保證高的空間分辨率和后期實驗數(shù)據(jù)處理結(jié)果的可靠性。目前，BL08U采用激光干涉儀確定樣品和各部件的移動位置，可確保波帶片相對樣品橫向位置的長期穩(wěn)定性，還能測量掃描樣品過程中圖像在x和y方向上的輕度失真[2]。本文根據(jù)激光干涉儀提供的掃描圖像在x和y方向的輕度失真數(shù)據(jù)，用slope函數(shù)求出樣品在移動過程中位置漂移的規(guī)律，從而修正樣品位置。盡管如此，圖像在x和y方向依然有一定漂移，給分析實驗結(jié)果造成很大誤差，需對數(shù)字圖像進(jìn)行后期配準(zhǔn)。圖像匹配技術(shù)是實現(xiàn)圖像融合、圖像校正、圖像鑲嵌及目標(biāo)識別與跟蹤的關(guān)鍵步驟之一，已廣泛應(yīng)用在圖像識別和圖像重建等領(lǐng)域中。如何自動提取穩(wěn)定可靠的特征、提高匹配準(zhǔn)確度是實現(xiàn)圖像匹配的重要環(huán)節(jié)，國內(nèi)外學(xué)者已做了大量研究工作。圖像匹配技術(shù)分為基于像素和基于特征的方法[3–6]。本文在軟件上采用基于統(tǒng)計的圖像相關(guān)性最大化方法來配準(zhǔn)圖像，并發(fā)展了一種基于圖像分窗口相關(guān)性的數(shù)字圖像自動配準(zhǔn)方法，著重針對雙能襯度圖像進(jìn)行配準(zhǔn)測試。結(jié)果顯示配準(zhǔn)誤差小于1個像素尺度。

1 用硬件對圖像漂移校正

激光干涉儀能提供不同能量掃描時每幅圖像在x和y方向的漂移。儀器本身有固有誤差，所以用激光干涉儀提供的漂移數(shù)據(jù)，通過 slope函數(shù)計算這些數(shù)據(jù)點的線性回歸擬合線方程的斜率，即樣品臺在移動過程中的漂移規(guī)律。用該規(guī)律校正樣品在移動過程中的位置漂移。

1.1 樣品

實驗中樣品是用銅質(zhì)的多級不同間隙圓環(huán)組成的圓形標(biāo)準(zhǔn)樣品，如圖 1(a)，其系經(jīng)掃描透射得到的透射圖，中心小圓區(qū)域為空，第一層圓環(huán)由寬度30 nm的空隙和銅片組成，依次類推，圓環(huán)逐層從內(nèi)到外銅片間距分別為30、60、120、240和480 nm。

1.2 實驗方案

標(biāo)準(zhǔn)樣品置于樣品架上，調(diào)整樣品至波帶片焦點處。用不同能量(N個)的X射線掃描樣品，得到N幅吸收圖，在每幅圖上確定一個有共同標(biāo)志性的點(圖1b)。比較各圖上該點相對于某幅圖像的偏移，可得到樣品在 x和 y方向上的一組漂移數(shù)據(jù)。用slope函數(shù)計算這些數(shù)據(jù)點的線性回歸擬合方程的斜率dx/dz和dy/dz，用其對掃描過程中樣品臺位置在x和y方向上的漂移作動態(tài)修正，則可確保吸收圖數(shù)據(jù)的可靠性。

圖1 標(biāo)準(zhǔn)樣品的吸收襯度圖像(a)和一部分的吸收襯度圖像(b)(b)中箭頭所指地方是參照點，各圖間比較都以該點為基準(zhǔn)，得到位置漂移Fig.1 STXM images (a), and a part of the images (b), of a standard sample.The arrowed position is the reference point from which the image position drifts are calculated.

1.3 實驗結(jié)果和分析

圖2為一組832–907 eV(對應(yīng)于不同z軸位置)X射線掃描標(biāo)準(zhǔn)樣品的吸收圖，步長5 eV。圖中明顯可見，樣品沿z軸移動過程中在x,y平面上確實存在漂移。激光干涉儀測得樣品掃描過程中的具體位置如表1。表中SZ是樣品沿z軸移動過程中z軸上的坐標(biāo)位置，XF是二維掃描時樣品在x方向的漂移量，YF是在y方向上的漂移量。用slope函數(shù)計算這些數(shù)據(jù)點的線性回歸擬合方程的斜率dx/dz和d y/dz。slope函數(shù)返回的是由這些數(shù)據(jù)點擬合得到的直線方程的斜率(即變化率)，具體表達(dá)式為：

式中，n是參數(shù)個數(shù)，x和 y是樣本平均值A(chǔ)verage(known_x's)和 Average(known_y's)。具體形式 是 SLOPE(known_y's，known_x's)[7]， 其 中known_x's是樣品在 z方向位置數(shù)據(jù)的集合，known_y's是樣品在x或y方向上漂移量數(shù)據(jù)的集合。

把表1數(shù)據(jù)代入公式(2)和(3)求出樣品沿z方向移動過程中在x,y方向的漂移率：

圖2 832–907 eV X射線(z軸不同位置處)掃描標(biāo)準(zhǔn)樣品的吸收襯度圖(無漂移修正)Fig.2 STXM images by X-rays of 832–907 eV, before calibration with the laser interferometer.

slope函數(shù)中 known_y's是實驗中實際測得樣品在z軸運(yùn)動過程中不同位置對應(yīng)的漂移量，其本身是條曲線，由無數(shù)點依次相連，相對曲線每一點都對應(yīng)一個斜率，隨曲線延伸變化，每一點的斜率都不同，slope的作用就是把每一點的斜率作為新曲線的值，構(gòu)成一條新曲線，即樣品在z軸移動過程中y方向的漂移規(guī)律。然后把上述得到的樣品在x,y方向上的漂移率用到樣品移動中，對樣品漂移動態(tài)修正。在能量832—907 eV內(nèi)重復(fù)前面實驗，能量步長取5 eV，結(jié)果如圖3。

表1 樣品位置在x,y方向上相對參照點的漂移數(shù)據(jù)Table 1 The drift data of STXM images relative to the reference point in x and y direction.

圖3 加入動態(tài)修正后重復(fù)實驗步驟得到的吸收圖Fig.3 The STXM images calibrated dynamically in the scanning process using the calculated drift rates.

由圖 3，通過激光干涉儀定位及掃描過程動態(tài)校正后，樣品圖像在 x,y方向上漂移已很大程度上被消除，但具體的實驗掃描出來的兩幅或多幅圖像間還存在一定程度的漂移，堆棧圖像間還不十分匹配，單獨掃描兩幅圖像有時會有幾十個像素點的漂移，如不校準(zhǔn)仍會嚴(yán)重影響后續(xù)數(shù)據(jù)分析的正確性，因此還需用軟件對各幅吸收圖像進(jìn)一步配準(zhǔn)。

2 用軟件對圖像漂移校正

2.1 圖像配準(zhǔn)

圖像配準(zhǔn)是將不同時間、不同傳感器(成像設(shè)備)或不同條件下(天候、照度、攝像位置和角度等)獲取的兩幅或多幅圖像進(jìn)行匹配、對準(zhǔn)和疊加，已廣泛應(yīng)用于遙感數(shù)據(jù)分析、計算機(jī)視覺、圖像處理等領(lǐng)域[8]?，F(xiàn)在常用的圖像配準(zhǔn)方法，有基于特征和基于統(tǒng)計的圖像配準(zhǔn)[9]，前者提取圖像信息的特征并以這些特征為模型配準(zhǔn)，圖像的特征點遠(yuǎn)少于圖像的像素點，配準(zhǔn)過程的計算量就大為減少，但特征提取的計算代價較大，也不便實時應(yīng)用；后者是指最大互信息的圖像配準(zhǔn)方法，與前者相比，其突出優(yōu)點是魯棒性(robustness)好、配準(zhǔn)精度高、人工干預(yù)少[10]?；诨バ畔⒌膱D像配準(zhǔn)，是用兩幅圖像的聯(lián)合概率分布與完全獨立時的概率分布的廣義距離來估計互信息。理論上，上述兩法均適用于STXM圖像配準(zhǔn)。實際應(yīng)用中，STXM實驗通常在感興趣元素吸收邊附近能量上對樣品掃描，各圖像襯度差別有可能較大，在能量掃描范圍不大情況下，基于統(tǒng)計的兩幅圖相關(guān)性最大化的圖像配準(zhǔn)方法可滿足要求；但若能量掃描范圍大，或兩幅圖襯度差別甚大，則基于統(tǒng)計相關(guān)性的配準(zhǔn)方法有可能失效。鑒此，我們發(fā)展了一種新的基于圖像分窗口相關(guān)的數(shù)字圖像自動配準(zhǔn)方法，可對吸收襯度差別非常大的兩幅或多幅圖準(zhǔn)確配準(zhǔn)。

2.2 整體相關(guān)最大化圖像配準(zhǔn)方法

根據(jù)STXM工作模式和數(shù)據(jù)特征的具體情況，用基于統(tǒng)計的兩幅圖相關(guān)性最大化圖像配準(zhǔn)方法。由于掃描兩幅圖像的時間間隔很短，實驗站的光強(qiáng)變化很小，即使光強(qiáng)變化很大，也可通過計算像素點的光密度對數(shù)據(jù)歸一化處理，認(rèn)為環(huán)境無變化，只有能量變化，導(dǎo)致樣品不同部位的吸收發(fā)生變化。這些都會反映在數(shù)字圖像上，因此用兩幅圖像的相關(guān)性來進(jìn)行圖像配準(zhǔn)。兩幅圖的相關(guān)性用互相關(guān)(CC)表示[11]：

將此法應(yīng)用到測量苔蘚植物頸部細(xì)胞切片中的錳元素分布實驗中結(jié)果見圖4。其中，圖4(a)和(b)是在錳元素吸收邊能量E1=640.3 eV和吸收邊前能量E2=638.6 eV測得的吸收圖像，都未經(jīng)過圖像漂移校正，圖像尺寸20 μm×20 μm，步長50 nm，像素數(shù)目401×401。把圖4(b)作參照圖，用基于統(tǒng)計的兩幅圖相關(guān)性最大化圖像配準(zhǔn)算法計算得出圖4(a)相對(b)的漂移量為(–3,1)，即把圖4(a)向 x 負(fù)方向和y正方向平移3和1個像素點，得到圖4(c)。把圖4(c)平移的部分和圖4(b)與之對應(yīng)的區(qū)域都裁掉，得到對齊的公共區(qū)域。圖4(d)和(e)是與(a)和(b)分別對應(yīng)的、經(jīng)漂移校正及裁剪的兩幅圖。裁剪后兩幅圖大小都是19.85 μm×19.95 μm，步長50 nm，像素數(shù)目為398×400。

圖4 用兩幅圖相關(guān)性最大化方法的圖像漂移校正過程(a) E1處校正前圖像；(b) E2處校正前圖像；(c) 根據(jù)計算出漂移參數(shù)對圖像(a)平移的結(jié)果；(d) 圖像(a)配準(zhǔn)后公共區(qū)域；(e) 圖像(b)配準(zhǔn)后公共區(qū)域Fig.4 Image registration process using global correlation maximizing technique.(a) before calibration at E1=640.3 eV; (b) before calibration at E2=638.6 eV; (c) pixel-moved Image (a) based on the drift parameters;(d) the registered area of Image (a); (e) the registered area of Image (b)

以上結(jié)果表明，基于統(tǒng)計的兩幅圖相關(guān)性最大化的圖像配準(zhǔn)方法在掃描能量范圍不大情況下可滿足圖像配準(zhǔn)的要求。

2.3 分窗口相關(guān)最大化圖像配準(zhǔn)方法

基于兩幅圖相關(guān)性最大化的圖像配準(zhǔn)方法在STXM圖像配準(zhǔn)實際應(yīng)用中有時會失效，特別是成像特性差異較大的圖像間配準(zhǔn)。為此，綜合基于灰度和基于特征的兩類配準(zhǔn)方法，我們發(fā)展了一種基于圖像分窗口相關(guān)的數(shù)字圖像自動配準(zhǔn)方法[12]，它結(jié)合了STXM具體特點，可對STXM系列圖像配準(zhǔn)。

此方法是將數(shù)字圖像看作一個二維灰度分布函數(shù)，選取兩幅或多幅圖中襯度最高的一幅圖像作參照圖，其它的以它為參照圖配準(zhǔn)。在STXM中一般選取感興趣元素吸收邊能量處的掃描圖像作參照圖，然后在參照圖上選擇一個局部區(qū)域作臨時窗口，此窗口也可看作是基于特征配準(zhǔn)方法中的一個特征。在另一目標(biāo)圖上尋找與其對應(yīng)最相似的窗口位置(移動窗口)。STXM 圖像漂移只發(fā)生在垂直光路方向的x,y平面內(nèi)，沒有旋轉(zhuǎn)，空間分辨率也固定，所以窗口配準(zhǔn)只需窗口(象元)在x和y方向上的變化(n, m)，兩參數(shù)的基本單位是一個像素點。具體算法處理過程如下：

(1) 首先定義窗口大小W×W，其中W單位為像素點，W大小能滿足將參照圖劃分為整數(shù)個窗口的條件。然后將參照圖及目標(biāo)圖按該尺寸劃分成若干個臨時窗口，同時給出目標(biāo)圖像上各窗口的起始坐標(biāo)值(x0，y0)和坐標(biāo)移動范圍±r。r為像素數(shù)目[13]，實際應(yīng)用中r一般取30。

(2) 移動目標(biāo)圖像上一個分窗口的位置并計算它與參照圖像上對應(yīng)窗口的相關(guān)值，再用循環(huán)比較計算方法搜尋一定范圍內(nèi)最大相關(guān)值的位置，并記錄該位置相對參照圖像上對應(yīng)窗口位置的移動參數(shù)(漂移量 xshift和 yshift)。

(3) 移動并計算下一窗口對應(yīng)相關(guān)值最大時的窗口位置，直至完成所有分窗口匹配，并產(chǎn)生包括所有分窗口位置移動參數(shù)的配準(zhǔn)文件。此處結(jié)合STXM圖像的特點把各分離窗口作為獨立窗口分別匹配。計算相關(guān)值采用公式類似求兩幅圖像的互相關(guān)(CC)公式。

(4) 計算所有參考窗和移動窗間的相對位移參數(shù)，找出相同漂移量最多的漂移值，此即圖像的漂移。

將該方法應(yīng)用于STXM吸收襯度圖像，兩幅圖相關(guān)性很小時也能很好地配準(zhǔn)，如圖5。其中圖5(a)和(b)是未配準(zhǔn)的STXM掃描得到的原始圖像，兩幅圖的位置發(fā)生了變化，相關(guān)性很差。用基于分窗口相關(guān)的數(shù)字圖像自動配準(zhǔn)方法算出圖5(b)相對于(a)的漂移量為(3,7)，把(b)向x和y正方向平移3和7個像素點，得到圖 5(c)。將目標(biāo)圖像按照計算出來的漂移量平移后，把平移的那部分在兩幅圖的對應(yīng)區(qū)域都裁掉，保留配準(zhǔn)的那部分共同區(qū)域，得到圖5(d)和(e)。兩幅圖裁剪后的大小都是 2.85 μm×2.6 μm，空間步長是25 nm，像素數(shù)目為115×114。這樣，就可對這兩幅STXM圖像進(jìn)行實驗數(shù)據(jù)處理。

該方法能對STXM圖像進(jìn)行較精確地配準(zhǔn)，很適合STXM的圖像處理。此算法在STXM圖像處理中準(zhǔn)確、有效，能解決成像差異較大的圖像間的配準(zhǔn)問題，有效消除圖像漂移的影響，特別是消除了雙能襯度圖像中的浮雕現(xiàn)象，確保了后續(xù)實驗數(shù)據(jù)處理的準(zhǔn)確性。

圖5 用分窗口相關(guān)配準(zhǔn)方法的圖像配準(zhǔn)過程(a) E1處校正前圖像；(b) E2處校正前圖像；(c) 根據(jù)計算出漂移參數(shù)對(a)平移的結(jié)果；(d) (a)配準(zhǔn)后公共區(qū)域；(e) (b)配準(zhǔn)后公共區(qū)域Fig.5 Image registration process using multi-window cross correlation technique.(a) before calibration at E1=707.2 eV; (b) before calibration at E2=703.9 eV; (c) pixe-movedI Image (a) based on the drift parameters;(d) the registration area of Image (a); (e) the registration area of Image (b)

3 結(jié)論

本文從硬件和軟件兩方面探討了STXM對樣品掃描出來的兩幅或多幅圖像空間位置漂移的校正問題。硬件上，通過激光干涉儀對標(biāo)準(zhǔn)樣品的精確定位及用slope函數(shù)計算得到樣品在掃描過程中漂移的規(guī)律，從而對樣品的位置漂移進(jìn)行動態(tài)修正，實驗結(jié)果表明用硬件對樣品位置動態(tài)校正切實可行。軟件上，基于統(tǒng)計的兩幅圖相關(guān)性最大化圖像配準(zhǔn)方法對大部分圖像都能很好配準(zhǔn)，但也有失效的時候，因此用遙感領(lǐng)域圖像配準(zhǔn)中的基于圖像分窗口相關(guān)性的配準(zhǔn)方法基本理論，發(fā)展出一種適用于STXM掃描圖像基于圖像分窗口相關(guān)的數(shù)字圖像自動配準(zhǔn)方法。對兩幅或多幅STXM圖像進(jìn)行配準(zhǔn)時，該方法對每一子區(qū)的自動相關(guān)匹配來達(dá)到對整幅圖像的精確配準(zhǔn)。具體處理結(jié)果表明該方法能對STXM圖像精確地配準(zhǔn)，且能解決成像差異較大圖像間的配準(zhǔn)問題，有效消除圖像漂移的影響，特別是消除了雙能襯度圖像中的浮雕現(xiàn)象，精度可達(dá)小于1個像素尺度。實際上此方法可達(dá)亞像素級配準(zhǔn)，但目前STXM只針對像素級計算，所以算法上采用限定在亞像素級范圍內(nèi)最多的漂移量取整作為配準(zhǔn)參數(shù)。圖像實際漂移量不僅像素級，有些是亞像素級，該方法確保了后續(xù)實驗數(shù)據(jù)處理的精度。

1 Chao W I, Harteneck B D, Liddi E J A, et al.Nature,2005,435: 1210–1213

2 Kilcoyne A L D, Tyliszczak T, Steele W F, et al.J Synchrotron Radiat, 2003, 10: 125–136

3 Keller Y, Averbuch A, Israeli M.IEEE Trans Image Process, 2005, 14(1): 12–22

4 Siggelkow S.Feature histograms for content-based image retrieval.Frieiburg: Albert Ludwigs University of Frieiburg, 2002

5 Mikolajczyk K, Schmid C.An affine invariant interest point detector.Proceedings of the 7th European Conference on Computer Vision

6 楊曉敏, 吳 煒, 卿粼波, 等.光學(xué)精密工程, 2009, 17(9): 2276–2282 YANG Xiaomin, WU Wei, QING Linbo, et al.Opt Precis Eng, 2009, 17(9): 2276–282

7 Hohenwarter M, Preiner J.J Online Math Appl, 2007, 7:Article ID 1574

8 Zitov′a B, Flusser J.Image Vision Comput, 2003, 21:977–1000

9 吳 畏, 趙文杰, 劉 輝.紅外, 2009, 30(10): 37–43 WU Wei, ZHAO Wenjie, LIU Hui.Infrared, 2009, 30(10):37–43

10 Viola P, Wells W M.Int J Comput Vision, 1997, 24(2):137–154

11 Pratt W K.Digital image processing.2nd ed., New York:Wiley, 1991

12 李 震, 范湘濤, 施建成.中國圖像圖形學(xué)報, 2001,6A(2): 129–132 LI Zhen, FAN Xiangtao, SHI Jiancheng.J Image Graphics, 2001, 6A(2): 129–132

13 沈 楠, 曹計昌.計算機(jī)工程與應(yīng)用, 2004,40(3): 75–78 SHEN Nan, CAO Jicang.Comput Eng Appl, 2004,40(3):75–78