角分辨光電子能譜數據的快速拼接算法設計

賈小文

(陸軍軍事交通學院 基礎部，天津 300161)

1887年，赫茲發現當入射光照射到金屬表面時，會有電子從金屬表面逸出. 1905年，愛因斯坦提出了“光電子”的概念，從量子力學的角度完美地解釋了光電效應. 由于光電子攜帶了材料內部電子結構的信息，因此可以通過探測光電子來研究材料的電子結構，角分辨光電子能譜就是該原理的具體應用[1-3]. 實驗上，通過電子能量分析器探測出射的光電子，得到數據為某一方向立體角內電子的統計數. 從實驗得到的光電子統計數是按照角度分布的. 為了方便研究，需要將角度轉換為動量(實際上是動量在kx-ky面上的投影)，經過轉換后的數據為光電子統計數隨動量(kx,ky)的分布. 在角度空間分布于一條線上的光電子(目前主流探測器都是這種模式)，轉換到動量空間的分布是一條彎曲的曲線. 因此要獲得某個動量平面內所有光電子的分布，就要求對各個不同角度下(通常這個角度就是測量樣品的取向)的測量數據進行拼接[4]. 在拼接數據時，傳統的方法先由角度換算為動量，然后在動量空間進行均勻插值，這種方法存在耗時較大的缺點. 由于角分辨光電子能譜實驗要求極度清潔的樣品表面，如果耗時過長，樣品清潔度將會下降，從而導致實驗失敗. 因此無論是從數據處理效率還是從開展實驗的角度來看，都需要對傳統的方法進行改進.

本文充分考慮了角分辨光電子能譜實驗數據的相鄰傾角數據基本類似的特點，提出了快速數據拼接算法，該算法首先生成均勻分布的動量點，然后逆用角度動量轉換公式獲取對應的角度值，最后在角度空間進行插值. 由于避開了普通插值的三角剖分和幾何排序等復雜過程，可將運算速度提高2個數量級.

1 角分辨光電子能譜實驗原理及數據的拼接

角分辨光電子能譜實驗的基本原理是光電效應，如圖1所示. 探測器收集電子的窗口是1條狹縫，1次收集1個狹長立體角空間的所有電子. 原始數據為網格化的矩陣，橫坐標是動能，縱坐標是角度，數據表示光電子的譜重，如圖2所示. 處理數據時，需要將角度轉換成動量，轉換的公式取決于樣品與電子探測器相對方位的布局和定義. 以文獻[4]為例，分析器相對于樣品的角度定義取圖3所示定義：將與狹縫平行的方向定義為偏轉角度φ，與狹縫垂直的方向定義為傾角θ，樣品轉過的角度定義為旋轉角ω[4-5]. 動量(kx和ky面的投影)和實驗測量角度的轉換公式為

(1)

式(1)中，k0表示出射電子的動量，由式(2)決定：

(2)

從式(1)中可以看到，探測器探測到的動量實際上是水平方向的動量，因此角分辨光電子能譜經常用于研究電子結構二維性較好的材料.

圖1 角分辨光電子能譜實驗基本原理示意圖

圖2 典型角分辨光電子能譜原始數據圖譜

(a)繞樣品面法向旋轉角度為0時，偏轉角度φ、 傾角θ與動量的換算關系

(b)繞樣品面法向旋轉ω角度后動量的換算關系圖3 實驗中樣品、探測器狹縫取向與動量換算幾何關系

為了得到不同動量區域的光電子信息，需要不斷地改變狹縫相對于樣品的傾角θ，并將這些數據拼接在一起. 在拼接時，圖2中的原始數據在角度空間是均勻分布的，但是利用式(1)轉換到動量空間(其實是kz=C的某個平面)后，(kx,ky)就不再是均勻的網格分布了，有些動量點有數據，有些動量點則沒有數據，如圖4所示.

圖4(a)所示為鐵基超導體FeSe0.3Te0.7在20 K下費米能級處譜重隨動量的分布(光子能量22.218 eV，光電子能譜儀為Scienta公司的R4000),可以看到圖4(a)的數據分布并不均勻，有些點有數據，有些點沒有數據(相當于譜重等于0)，是非均勻分布的數據圖. 出現這種情況的原因有：1)角度到動量的轉換算法[式(1)]是非線性的，連續均勻的角度分布換算到動量空間就不再是連續均勻分布了；2)探測器在傾角方向分辨率太大(傾角人為設定，一般是1°～2°，遠大于探測器狹縫方向的角度分辨率)，數據呈現明顯的條狀分布. 這種非均勻分布的數據給后續的分析帶來不便，實踐中需要將圖4(a)的結果通過插值拼接成分布連續、光滑的數據，如圖4(b)所示. 拼接的過程中需要解決數據不均勻分布的問題，即非均勻數據網格化.

(a) 直接將角度轉換成動量

(b)插值拼接以后的結果 圖4 費米能級處譜重隨動量的分布

2 非均勻數據網格化

非均勻數據網格化是常見的數據處理方法，有非常成熟的數值算法，其過程是先對自變量空間進行三角剖分，再利用三角剖分結果重新抽樣. 很多軟件(如Igor，Matlab等)都有相應的函數完成上述運算. 這里介紹Igor下的做法[6-7].

首先創建N×3維的矩陣M，對每個傾角下的譜應用式(1)計算(kx,ky)，分別記錄在M的第1列和第2列，M的第3列存放對應的測量值. 然后利用ImageInterpolate命令進行插值運算：

Imageinterpolate/S={x0,dx,x1,y0,dy,y1}

Voronoi,M

以上做法是目前角分辨光電子能譜實驗中數據拼接的傳統做法，圖4(b)的結果就是利用該方法計算得到的，其缺點是運算時間過長，通常需要1 min左右(取決于譜的數量)，甚至更長. 耗時來源于三角剖分過程及以三角剖分為條件的插值幾何排序過程.

3 快速數據拼接算法

三角剖分法適用的范圍廣，但某些應用場合效率較低. 實際上，角分辨光電子能譜實驗數據有自己的特點：數據的自變量并非完全隨機分布，相鄰狹縫之間除了傾角不同，其他的實驗條件完全相同. 根據這個特點構建了快速拼接算法：逆用式(1)從動量值反過來推算角度，利用角度值在原始數據中插值，將插值作為該動量處的結果. 快速拼接算法的優點是：

1)在動量空間，可以預先定義均勻分布的網格數據，省去了重新抽樣的過程；

2)在角度空間，由于數據的分布是完全均勻的，可直接通過線性插值來獲取插值結果，省去了三角剖分的計算過程.

快速拼接算法的復雜度為O(n)，因此運算速度非常快. 實驗結果表明:快速拼接算法要比傳統算法運算速度提高2個數量級.

快速拼接算法的詳細描述如下：

1)創建矩陣M(kx,ky).

2)取出(kx,ky)，逆用式(1)計算(θ,φ).

3)對原始測量數據通過偏轉角、傾角進行排序，找出相鄰傾角θ1和θ2，使得θ1<θ<θ2，φmin<φ<φmax.φmin和φmax為原始譜最小角和最大角加上偏轉角φ.

4)設θ1對應的原始數據是data1，θ2對應的原始數據是data2，利用雙線性插值分別獲取2個值z1=data1(0)(φ)，z2=data2(0)(φ). 這里0表示費米能級.

5)通過線性插值得到(θ,φ)處數據：

6)將M中(kx,ky)處相應的值設置為z.

下面是算法具體實現(Igor代碼)：

variable N

Make/O/N=(N,N) M_FS

〈setscale〉

〈get wtheta, wphi, womega, wdatas〉

Sort {wtheta, wphi} wtheta, wphi, womega, wdatas

M_FS=render_FS(x, y, wtheta, wphi, wdatas)

Function render_FS(kx, ky, wtheta, wphi, womega, wdatas)

variable kx, ky

wave wtheta, wphi, womega wave/wave wdatas

variable omega=womega[0]

variable kx1, ky1, theta, phi

variable k0=sqrt(2*me*e0)/hb*sqrtEk/1E10*LA/pi

kx1=kx*cos(omega)-ky*sin(omega)

ky1=kx*sin(omega)+ky*cos(omega)

phi=asin(kx1/k0)

theta=asin(ky1/k0/cos(phi))

theta=theta*180/pi

phi=phi*180/pi

〈get ipos using (theta, phi) from wtheta, wphi〉

wave w1=wdatas[ipos]

wave w2=wdatas[ipos+1]

z1=w1(0)(phi)

z2=w2(0)(phi)

return (z2-z1)/(theta2-theta1)*(theta-theta1)+z1

end

4 實驗結果

以圖4所示數據拼接為例，對上述算法進行測試. 測試數據一共有77張譜，每張譜在角度方向的數據長度為360.

圖5是快速拼接算法給出的插值結果. 可以看到圖5(a)和圖4(b)的結果幾乎沒有區別. 為了進一步確認快速拼接算法的可靠程度，圖5(b)給出了沿不同方向切面[ProfileLine，方向如圖5(a)中紅色虛線箭頭所指]的動量分布曲線(Momentum distribution curves，MDC)，為了方便比較，曲線都進行了歸一化. 可以看出快速拼接算法在曲線數據特征如峰位、半高全寬等方面與傳統算法是一致的. 盡管如此，如果仔細觀察，傳統算法和快速拼接算法還是存在一些微小的差異. 出現這種現象的原因主要是插值機制的區別. Igor下傳統算法使用的泰森多邊形插值(Voronoi interpolation)是平滑插值方法，而本文中的快速拼接算法則完全采用線性插值. 另外，由于計算過程不同，計算機的舍入誤差也會給插值結果帶來一些差別.

(a) 快速拼接算法對費米面附近譜重隨動量分布的拼接結果

(b)動量分布曲線圖5 快速拼接算法結果

表1給出了Igor Pro (Ver 6.37)下二者所用時間的對比. 從表1可以看出快速拼接算法的運算速度要比傳統算法提高了近200倍.

利用多線程技術，還可以將快速拼接算法并行分配到多個處理器(多線程)，能夠進一步將運算速度提高到300多倍. 多線程實現的方法如下(Igor代碼)：

1)將渲染函數render_FS聲明為線程安全函數.

Threadsafe function render_FS(…)〈function body〉End

2)利用multithread指令調用渲染函數對M_FS賦值.

Multithread M_FS=render_FS(x, y, wtheta, wphi, wdatas)

快速拼接算法(多線程)和傳統三角剖分算法運行消耗時間比較如表1所示.

表1 快速拼接算法、多線程和傳統三角剖分算法運行消耗時間比較

注：CPU I3-3320 4核處理器，3.3 GHz. RAM：2.00 G.

5 討 論

原始數據譜的數量不同，每張譜的分辨率不同，無論是傳統算法還是快速拼接算法，所用的時間都是不同的. 大量的實驗結果表明，對于常見的角分辨光電子能譜實驗，傳統算法插值需要的時間集中在30 s～3 min，而快速拼接算法時間基本不會超過1 s，因此快速拼接算法對于速度的提升和本文中的例子(～200倍)基本一致，大概2個數量級. 在實踐中，采用快速拼接算法，可以實時獲取插值結果，這在使用傳統算法時是不可能實現的.

快速拼接算法的適用范圍要比傳統算法的適用范圍窄. 當測量不同動量空間是通過旋轉樣品(以樣品面法向為旋轉軸)進行時，快速拼接算法就不再適用了.

對于常見的角分辨光電子能譜實驗，通常是通過改變樣品相對于分析器的傾角來測量不同的動量空間，相鄰傾角(通常少于2°)之間的數據差別不是很大. 快速拼接算法正是利用了這一特點對相鄰傾角之間不存在的角度用兩相鄰數據線性插值結果來代替，這和傳統插值算法從插值原理上來說并沒有本質區別(此時平滑插值可以用線性插值近似)，因此快速拼接算法和傳統插值算法的插值結果一致. 但是由于跳過了三角剖分和插值幾何排序等復雜耗時的計算過程，快速拼接算法提高了運算速度.

需要注意的是無論是快速拼接算法還是傳統算法，都是通過已測數據推斷未測量點的數據，從這點來說，本文所述傳統方法亦不能作為絕對的參考，真正的插值效果有賴于待插值的數據環境. 事實上，還有其他的不同的插值方法，如最近鄰法、克里金法(Kriging)等，這些方法同樣能將非均勻分布的數據拼接為均勻分布的數據，但是插值的結果之間必然會存在一些差別.

利用快速拼接算法，可以實時查看不同結合能處譜重隨動量的分布，還可以動態顯示譜重動量分布隨結合能的變化，這可以提高角分辨光電子能譜數據處理效率.

[1] 黃昆,韓汝琦. 固體物理學[M]. 北京:高等教育出版社,2010:213-223.

[2] Damascelli A, Shen Z X, Hussain Z. Angle-resolved photoemission of the cuprate superconductors [J]. Rev. Phys., 2002,75(2):473-541.

[3] Zhou X J, Cuk T, Devereaux T, et al. Angle-resolved photoemission spectroscopy on electronic structure and electron-phonon coupling in cuprate superconductors in handbook of high-temperature superconductivity [M]. New York: Springer, 2007:87-144.

[4] Liu G D, Wang G, Zhu Y, et al. Development of a vacuum ultraviolet laser-based angle-resolved photoemission system with a superhigh energy resolution better than 1 meV [J]. Rev. Sci. Instrum., 2008,79:023105.

[5] 賈小文. 角分辨光電子能譜對重費米子超導體CeCoIn5和鐵基超導體Fe(SexTe1-x)的研究[D]. 北京:中國科學院物理研究所,2017.

[6] 丁永祥,夏巨諶. 任意多邊形的Delaunay三角剖分[J]. 計算機學報,1994(4):270-275.

[7] Wavemetrics, Inc. The manual of IGOR Pro [Z]. 2015:1516-1518.