微分スペクトルを計算する

オージェ電子分光装置で測定を行うと、ダイレクトスペクトルの他に微分スペクトルを使うと便利な場合があります。しかし、たとえば JAMP-9500F（日本電子）の付属ソフトでは、微分スペクトルデータのエクスポートはできないように見えます。自分で計算しようと思ったら、どうすればいいでしょうか？

表示される微分スペクトルの仕様
計算する方法
次回予告
参考文献

表示される微分スペクトルの仕様

JAMP-9500F（の付属ソフト、以下省略）で表示される微分スペクトルは、以下の仕様になっているようです。

Savitzky-Golay法[1]によって微分スペクトルを計算している。
Savitzky-Golay法と一口にいってもパラメータがいろいろあるが、2次関数近似を用いて、ウィンドウサイズを7で計算している。これを簡単に2次7点のSavitzky-Golay法、と呼んだりする？

計算する方法

計算式

いま、入力 $x_i$ に対して出力 $y_i$ を得たデータ列があるとします。オージェ（分析のダイレクト）スペクトルであれば、電子のエネルギーが $x_i$ で、そのカウント数（もしくはカウント率）が $y_i$ になります。
Savitzky-Golay法は入力 $x_i$ が一定間隔であることを前提としており、その間隔を $\epsilon$ とします。 $x_i$ における微分値は、ターゲット近傍の測定値列に対して線形結合をするだけで求められます。上述の2次、7点の場合に具体的な係数を含めて示すと、以下のようになります。

$y'_{x=x_i} = \frac{-3y_{i-3}-2y_{i-2}-y_{i-1}+y_{i+1}+2y_{i+2}+3y_{i+3}}{28 \epsilon}$

ちなみに分母の28というのは分子の各項の係数の2乗和です。ウィンドウサイズを9にすると（容易に推測できる通りに）分子に $-4y_{i-4}$ と $+4y_{i+4}$ の項が増え、分母の係数が28ではなくて60になります。

両端の処理

データ列の端の方では、「近傍」にデータがないということが起こります。その場合の処理方法は、データの性質によっていろいろな方法がありますが、JAMP-9500Fでは最寄りの測定点（すなわち端点）の出力の値をそのまま使っています。データのインデックスが $i=0$ から始まっているときに $x_2$ における微分値を求めようとすると $y_{-1}$ の値が必要になりますが、それは測定していないので代わりに $y_0$ を使って、次のように計算します。

$y'_{x=x_2} = \frac{-3y_0-2y_0-y_1+y_3+2y_4+3y_5}{28 \epsilon}$

Pythonで計算する

Pythonを使っている場合は、科学技術計算用のライブラリ Scipyにそれ用の関数(scipy.signal.)savgol_filterが用意されています。yが出力データだけの1次元配列であるときに今の仕様で微分スペクトルを求めるには、（もちろんimport scipyした上で）以下のように使います。
y_d = scipy.signal.savgol_filter(y, 7, 2, 1, epsilon, mode = 'nearest')

次回予告

上記の計算式を導出してみたところ、具体的なケースだけを考えるのであれば意外と簡単でした。ということで、その雰囲気を味わう記事を投稿するかもしれませんし、しないかもしれません。

参考文献

Savitzky, A.; Golay, M.J.E. Smoothing and Differentiation of Data by Simplified Least Squares Procedures, Analytical Chemistry 1964 36(8) 1627-1639 DOI: 10.1021/ac60214a047 現時点では、弘前大学附属図書館の第2書庫に所蔵されているようです。