peakfit.py 技術ドキュメント
プログラムの動作
peakfit.py は、X線回折(XRD)パターンやその他のスペクトルデータに含まれるピークを自動的に検出し、その特性(中心位置、強度、半値幅)を推定するためのPythonプログラムです。
主な機能:
Excelデータ読み込み: 指定されたExcelファイル(
.xlsx)から、X-Y形式のスペクトルデータを読み込みます。データフィルタリング: X軸の指定された範囲に基づいてデータをフィルタリングし、解析対象の領域を限定します。
データ平滑化:
peaksearchモジュールを使用して、入力データにSavitzky-Golayフィルターなどの平滑化処理を適用し、ノイズを低減します。ピーク検出: 平滑化されたデータに対して、ピーク検出アルゴリズムを適用し、スペクトル内のピークの位置、強度、半値幅を自動的に特定します。
パラメータ初期値の生成: 検出されたピークの情報に基づいて、後続のピークフィッティング処理で使用するためのパラメータ(背景、ピークの強度、中心位置、半値幅)の初期値を構築します。
結果の可視化: 入力データ、平滑化されたデータ、および検出されたピークの形状(現在の実装ではローレンツ関数を使用)をグラフとしてプロットし、視覚的に確認できるようにします。
データ保存: フィルタリングされた入力データを新しいExcelファイルに保存します。
解決する課題:
複数のピークが重なり合う複雑なスペクトルデータにおいて、手動でのピーク特定は時間と手間がかかり、主観的な判断が入りがちです。peakfit.py は、このピーク特定プロセスを自動化し、客観的かつ効率的なピーク解析の初期段階を提供することで、ユーザーの負担を軽減します。これにより、より高度なフィッティングや解析に進むための堅牢な基盤を築くことができます。
原理
本プログラムは、主にデータ処理、ピーク検出、およびピーク形状の数学的表現に基づいています。
データ処理
入力データの読み込みとフィルタリング:
pandasライブラリを使用してExcelファイルからデータを読み込みます。その後、指定されたX軸の最小値 (xmin) と最大値 (xmax) に基づいてデータをフィルタリングします。これにより、解析対象の領域に絞り込むことが可能です。
ピーク検出アルゴリズム
ピーク検出は、外部の peaksearch モジュールによって行われます。このモジュールは、以下のステップでピークを特定すると推測されます。
データ平滑化:
nsmooth(平滑化ウィンドウサイズ) とnorder(多項式の次数) を用いたSavitzky-Golayフィルターなどの手法でデータを平滑化し、高周波ノイズを除去します。これにより、ピークの識別が容易になります。微分解析: 平滑化されたデータに対して1次微分および2次微分を計算します。
1次微分がゼロで、2次微分が負である点がピークの中心位置の候補となります。
ピークの裾野や形状を判断するために、1次微分の変化率 (
ydiff1_threshold) も利用されます。
閾値処理:
threshold(ピークの最小強度) とydiff1_threshold(1次微分の変化率の閾値) を用いて、ノイズによる誤検出を除外し、有意なピークのみを識別します。ピークパラメータ推定: 検出された各ピークについて、その中心位置 (
xc)、ピーク高さ(平滑化後のデータ値I0)、および半値幅 (w) を概算します。半値幅はピークの形状(ガウスまたはローレンツ)を仮定して推定されます。
ピーク形状関数
本プログラムは、スペクトルピークの形状を記述するためにガウス関数とローレンツ関数を実装していますが、現在の実装ではローレンツ関数が使用されています。
ガウス関数 (Gaussian): $\(f(x, x_0, \text{whalf}, A) = A \cdot \exp\left(-\left(\frac{x - x_0}{a}\right)^2\right)\)\( ただし、\)a = \frac{\text{whalf}}{\sqrt{\ln 2}}\( であり、\)\text{whalf}\( は半値幅、\)x_0\( はピークの中心位置、\)A\( はピークの高さまたは面積に関連する振幅です。コード内の \)A\( (デフォルト値) は、\)A = \frac{\sqrt{\ln 2 / \pi}}{\text{whalf}}$ と設定されており、これは特定の正規化条件を満たすための係数です。
ローレンツ関数 (Lorentzian): $\(f(x, x_0, \text{whalf}, A) = A \cdot \frac{1}{1 + \left(\frac{x - x_0}{\text{whalf}}\right)^2}\)\( ここで、\)\text{whalf}\( は半値幅、\)x_0\( はピークの中心位置、\)A\( はピークの高さまたは面積に関連する振幅です。コード内の \)A\( (デフォルト値) は、\)A = \frac{1}{\text{whalf} \cdot \pi}$ と設定されており、これも特定の正規化条件を満たすための係数です。
peak_func 関数では、Lorentzian 関数が使用され、I0 (強度) が直接乗じられます。したがって、解析結果のピーク関数は、以下の形式となります。
ここで \(I_0\) はピークの高さ、 \(x_c\) はピークの中心位置、\(w\) は半値幅を示します。
必要な非標準ライブラリとインストール方法
peakfit.py を実行するには、以下の非標準Pythonライブラリが必要です。
numpy: 数値計算を効率的に行うための基盤ライブラリ。pandas: データ構造とデータ解析ツールを提供するライブラリ。特にExcelファイルの読み書きに使用されます。matplotlib: グラフ描画のためのライブラリ。openpyxl:pandasがExcelファイル (.xlsx形式) を読み書きする際に内部的に利用するライブラリ。peaksearch: ピーク検出アルゴリズムを提供するカスタムモジュール。このモジュールは別途提供されるか、同じディレクトリに配置されている必要があります。
以下の pip コマンドを使用して、必要なライブラリをインストールできます。
pip install numpy pandas matplotlib openpyxl
peaksearch モジュールについては、Pythonの標準的なパッケージインストール方法ではインストールできません。peakfit.py と同じディレクトリに peaksearch.py ファイルが存在することを確認してください。
必要な入力ファイル
プログラムは、デフォルトで test/xrd.xlsx というパスにあるExcelファイル(.xlsx 形式)を入力として期待します。
ファイル形式:
ファイルの種類: Microsoft Excelワークブック (
.xlsx)データ構造: 少なくとも2つの列を持つデータで構成されている必要があります。
1列目: X軸データ(例: 2\(\theta\)、回折角度、エネルギー、波長など)。
2列目: Y軸データ(例: 強度、カウント数、吸光度など)。
ヘッダー: 1行目に各列のラベル(ヘッダー)が含まれていることを想定しています。
例 (xrd.xlsx の内容):
2theta Intensity
20.0 100
20.1 150
20.2 200
... ...
30.0 5000
30.1 6000
30.2 5500
... ...
パスの変更:
プログラムのソースコード内で、infile = 'test/xrd.xlsx' の行を編集することで、入力ファイルのパスを変更できます。
生成される出力ファイル
peakfit.py は、解析結果の一部をExcelファイルとして出力します。
ファイル名: デフォルトでは
simulated.xlsx内容: 入力データ (
infileで指定されたファイルから読み込まれ、xminおよびxmaxでフィルタリングされたデータ) と全く同じ内容が保存されます。具体的には、フィルタリング後のX軸データとY軸データが、元の列ラベルとともにExcelワークシートに書き込まれます。 現在の実装では、ピークサーチで検出されたピークのパラメータ(中心位置、強度、半値幅など)は、この出力ファイルには保存されません。これらの情報はプログラムの実行中にコンソールに表示され、グラフに可視化されます。
パスの変更:
プログラムのソースコード内で、outfile = 'simulated.xlsx' の行を編集することで、出力ファイルのパスと名前を変更できます。
コマンドラインでの使用例 (Usage)
peakfit.py プログラムは、コマンドライン引数を直接受け取りません。プログラムの動作を制御するための設定(入力ファイル名、出力ファイル名、解析範囲、ピーク検出パラメータなど)は、ソースコード内のグローバル変数や cfg オブジェクトに直接記述されています。
基本的な実行コマンドは以下の通りです。
python peakfit.py
このコマンドを実行すると、プログラムは内部で定義された設定値に基づいて動作します。
主な設定値(main 関数内の cfg オブジェクトで定義):
cfg.infile: 入力Excelファイルのパス (デフォルト:'test/xrd.xlsx')cfg.outfile: 出力Excelファイルのパス (デフォルト:'simulated.xlsx')cfg.xmin,cfg.xmax: X軸データの解析範囲 (デフォルト:20.0,40.0)cfg.norder: Savitzky-Golayフィルターの多項式の次数 (デフォルト:5)cfg.nsmooth: Savitzky-Golayフィルターのウィンドウサイズ(奇数である必要あり) (デフォルト:7)cfg.threshold: ピーク検出の最小強度閾値 (デフォルト:3000)cfg.ydiff1_threshold: 1次微分の変化率に基づくピーク検出閾値 (デフォルト:1.0)cfg.figsize,cfg.fontsize,cfg.legend_fontsize: グラフ描画に関する設定
これらの設定値を変更するには、peakfit.py のソースコードを直接編集する必要があります。
コマンドラインでの具体的な使用例
peakfit.py はコマンドライン引数を受け取らないため、実行方法は一つです。
実行コマンド:
python peakfit.py
実行結果の説明:
コンソール出力: プログラムが起動すると、以下の情報がコンソールに順次表示されます。
プログラムのタイトルと入出力ファイル名。
入力データの読み込み状況。
ピークサーチの進行状況 (
norder,nsmooth,threshold,ydiff1_threshold,ndataなど、設定されたパラメータが表示されます)。peaksearchモジュールによる詳細なピーク検出過程(例:Peak search start,Peak search endedなど)。検出されたピークのリスト: 各ピークについて、そのインデックス、中心位置 (
xc)、ピーク高さ (I0、平滑化後のY値)、および半値幅 (w) が整形されて表示されます。
Peak fit test program / odata module infile=test/xrd.xlsx outfile=simulated.xlsx Peak search mode : plot infile : test/xrd.xlsx xlabel : 0 ylabel : 1 x range : 20.0 40.0 norder : 5 nsmooth : 7 Warning: nsmooth must be odd. Changed to 7 (※もしnsmoothが偶数だった場合) threshold : 3000 dy/dx threshold : 1.0 ndata = [データポイント数] Peak search start ... Peak search ended. Npeak = [検出されたピーク数] peak list: 000: xc= 28.0000 I0= 3500.0 w= 0.5000 001: xc= 32.5000 I0= 6200.0 w= 0.7000 ...
ファイル出力: プログラムが正常に終了すると、指定されたパス(デフォルトは
./simulated.xlsx)に新しいExcelファイルが生成されます。このファイルには、test/xrd.xlsxから読み込まれ、X軸範囲 (xminからxmaxまで) でフィルタリングされたデータがそのまま書き込まれています。グラフ表示:
matplotlibによってグラフウィンドウがポップアップ表示されます。このグラフには以下の情報が含まれます。入力データ: 曲線または点でプロットされます(デフォルトでは青い丸)。
検出されたピーク: 各ピークの中心位置に垂直線、半値幅を示す水平線、および推定されたローレンツ関数曲線(デフォルトでは赤線)が重ねて描画されます。これにより、検出結果を視覚的に確認できます。
X軸とY軸にはそれぞれ入力データの列ラベルが表示されます。
グラフウィンドウは、コンソールで「Press ENTER to terminate>>」と表示され、ユーザーがEnterキーを押すまで開いたままになります。Enterキーを押すと、プログラムが終了し、グラフウィンドウも閉じます。
注意点:
プログラムを実行する前に、test/xrd.xlsx ファイルが実行ディレクトリまたは指定されたパスに存在すること、および peaksearch.py が同じディレクトリに存在することを確認してください。