webPDF / RustPDF: PDF 結果に影響する物理・数式上の相違

1. 対象と結論

この文書は、現在確認できる webPDF/localPDF の calc_pdf.py と RustPDF の common.rs, pipeline.rs, zdst.rs, zoom_dst.rs を比較し、 最終的な reduced PDF (G(r)) に影響する実装上の相違だけを整理したものです。

両者は同じ物理量を目標にしています。

[ S(Q)=i(Q)+1 ]

[ G(r)=\frac{2}{\pi}\int Q\,i(Q)\sin(Qr)\,dQ ]

[ g(r)=1+\frac{G(r)}{4\pi r\rho} ]

したがって、G(r) と Rust 側の reduced_gr は物理的には同じ reduced PDF です。 差が生じる主因は G(r) の定義そのものではなく、その前段の散乱因子・規格化・Q-grid と、 short-r correction/back transform の処理です。

---

2. 原子散乱因子の二つの平均量が異なる

2.1 本来必要な量

複素原子散乱因子を

[ f_j(Q)=a_j(Q)+i b_j(Q) ]

とすると、PDF の強度規格化では二つの異なる量が現れます。

[ \langle |f|^2\rangle =\sum_j c_j |f_j|^2 =\sum_j c_j(a_j^2+b_j^2) ]

および

[ |\langle f\rangle|^2 =\left|\sum_j c_j f_j\right|^2 ]

です。ここで (c_j) は組成分率です。

2.2 webPDF

CalcFseries() は複素数に対して

np.power(f, 2.0).real

を使用しています。

したがって実際に計算される第一の量は

[ \operatorname{Re}(f_j^2)=a_j^2-b_j^2 ]

であり、

[ |f_j|^2=a_j^2+b_j^2 ]

とは異なります。

同様に <f>2 側も

(fsum2_each_q * fsum2_each_q).real

なので、

[ \operatorname{Re}(\langle f\rangle^2) ]

を使っています。

2.3 RustPDF

calc_fseries() は Complex::norm_sqr() を使用しており、

[ \langle |f|^2\rangle ]

と

[ |\langle f\rangle|^2 ]

を明示的に計算しています。

2.4 PDF への影響

(f’‘\neq0) のとき両者は一致しません。差はそのまま強度規格化と (i(Q)) の分母・ self-scattering subtraction に入り、最終的に (Q\,i(Q)) を通じて (G(r)) 全体へ伝播します。

これは単なる FFT 実装差ではなく、散乱強度として使用している複素散乱因子の物理量そのものの差です。

---

3. 強度の規格化係数の定義が異なる

(D(Q)=\langle

f

^2\rangle+C(Q)) と書きます。

3.1 webPDF

Calc_iQ_SQ() の規格化係数は

[ A_{\mathrm{web}} = \frac{1}{N} \sum_{Q\in[Q_{\min},Q_{\max}]} \frac{I_{\mathrm{corr}}(Q)}{D(Q)} ]

です。すなわち 各 Q 点で比を取ってから平均しています。

その後、

[ i_{\mathrm{web}}(Q) = \frac{ I_{\mathrm{corr}}(Q)/A_{\mathrm{web}} -\langle |f|^2\rangle -C(Q) }{ |\langle f\rangle|^2 } ]

に相当する処理を行います。ただし前節のとおり、webPDF の散乱因子二乗は 実際には modulus square ではありません。

3.2 RustPDF

calc_iq_sq() の規格化係数は

[ A_{\mathrm{Rust}} = \frac{ \sum I_{\mathrm{corr}}(Q) }{ \sum D(Q) } ]

です。すなわち 強度の和と散乱項の和の比です。

その後、

[ i_{\mathrm{Rust}}(Q) = \frac{ I_{\mathrm{corr}}(Q)/A_{\mathrm{Rust}} -\langle |f|^2\rangle -C(Q) }{ |\langle f\rangle|^2 } ]

を計算します。

3.3 両者が一致する条件

一般に

[ \frac{1}{N}\sum\frac{I_k}{D_k} \neq \frac{\sum I_k}{\sum D_k} ]

です。

(D(Q)) が規格化区間内でほぼ一定、または (I(Q)/D(Q)) が特殊な条件を満たす場合にのみ 差が小さくなります。

3.4 PDF への影響

これは (i(Q)) の scale と baseline に直接作用するため、(G(r)) の振幅および termination/background 成分に直接影響します。

---

4. background subtraction 後の負値処理が異なる

webPDF

subtract_bg() は

[ I_{\mathrm{sub}}=I_{\mathrm{raw}}-aI_{\mathrm{BG}} ]

を計算した後、

[ I_{\mathrm{sub}}<0 \quad\Rightarrow\quad I_{\mathrm{sub}}=0 ]

として負値を clipping します。

RustPDF

subtract_bg() は

[ I_{\mathrm{sub}}=I_{\mathrm{raw}}-aI_{\mathrm{BG}} ]

をそのまま返し、負値を 0 に clipping しません。

PDF への影響

background が測定強度を局所的に上回る領域では、webPDF はその情報を 0 に置換し、 RustPDF は負の残差を保持します。

その差は polarization correction、規格化、(i(Q)) を経由して (G(r)) に入ります。

---

5. Q-grid と補間を行う位置が異なる

5.1 webPDF

webPDF は概略として

raw 2theta
 -> raw Q
 -> polarization / form factor / Compton
 -> normalization
 -> i(Q)
 -> i(Q) を等間隔 Q-grid に補間
 -> transform

です。

Interpolate() では

[ Q=0,\Delta Q,2\Delta Q,\ldots,\pi/dr ]

に相当する長い grid を作り、norm_max より高 Q の (i(Q)) を 0 にしてから変換します。 測定 Q 範囲より低い補間外領域には fill_value=-1 が指定されています。

5.2 RustPDF

RustPDF は

raw 2theta
 -> raw Q
 -> q_cut ... norm_max の等間隔 Q-grid を構築
 -> intensity / background をその grid へ補間
 -> form factor / Compton / polarization
 -> normalization
 -> i(Q)
 -> Zoom-DST

です。

つまり 強度を先に共通 Q-grid へ載せ、その grid 上で物理補正と規格化を行う構成です。

5.3 PDF への影響

補間と非線形な除算・規格化は一般に交換できません。

[ \mathcal{I}\left[\frac{I(Q)}{D(Q)}\right] \neq \frac{\mathcal{I}[I(Q)]}{D(Q)} ]

ここで (\mathcal{I}) は補間演算です。

また low-Q の扱いも異なるため、長距離側 baseline や termination 成分にも差が入り得ます。

---

6. Q→r 変換の離散化・規格化が異なる

両者が目標としている連続式は同じです。

[ G(r)=\frac{2}{\pi}\int Q\,i(Q)\sin(Qr)\,dQ ]

ただし離散化が異なります。

6.1 webPDF

CalcGrFFT() は SciPy の DST-I を使用し、

fft.dst(for_FT, type=1) / sqrt(nn_max + 1) / pi * 2

で G(r) を作ります。

SciPy DST-I の離散和と nn_max に基づく規格化を使用し、r-grid もその FFT/DST grid から決まります。

6.2 RustPDF

RustPDF の zoom_dst() は sine transform の離散和に入力 grid 間隔

[ \Delta Q ]

を明示的に掛けています。calc_gr() はその結果へ

[ \frac{2}{\pi} ]

を掛けます。

したがって実装上は積分

[ \frac{2}{\pi}\sum_j Q_j i(Q_j)\sin(Q_jr_k)\Delta Q ]

を直接近似する形です。

また Zoom-DST により、r-grid は FFT の natural grid に固定せず dr, r_max から指定します。

6.3 PDF への影響

連続極限で目標とする物理式は同じですが、有限 grid での quadrature、 normalization、Q/r sampling が異なるため、振幅・peak shape・termination ripple に 数値的な差が生じます。

これは「G(r) の物理的定義が違う」という意味ではなく、 同じ Fourier sine integral の離散表現が違うという差です。

---

7. short-r cutoff / back-Fourier correction の物理処理が異なる

ここは最終 PDF、とくに補正後 (G(r)) と Lorch (G(r)) に直接影響します。

7.1 webPDF

Rcut_and_bFT() はまず

[ p_0=\mathrm{int}(r_{\mathrm{cut}}/dr) ]

から低 r 側へ探索し、g(r) の符号変化

[ g(r_i)g(r_{i-1})\le0 ]

を見つけます。

その位置より低 r の g(r) を hard zero にします。

[ g_{\mathrm{bFT}}(r)=0 ]

その後、

[ G_{\mathrm{bFT}}(r) = 4\pi r\rho\,[g_{\mathrm{bFT}}(r)-1] ]

を DST で Q-space に戻します。

さらに raw (i(Q)) と back-transform (i_{\mathrm{bFT}}(Q)) の差

[ \Delta i(Q)=i(Q)-i_{\mathrm{bFT}}(Q) ]

へ moving average を掛け、

[ i_{\mathrm{smooth}}(Q) = i(Q)-\mathrm{MA}[\Delta i(Q)] ]

を最終的な再計算用 (i(Q)) とします。

現在の norm_for_bg は 10 です。

7.2 RustPDF

RustPDF は zero-crossing search のコードを残していますが、production path では

let cut_idx = 0;

として無効化されています。

実際には

[ r_t=\max(r_{\mathrm{cut}}-w,\ r_{\min}) ]

として cosine taper

[ W(r)= \begin{cases} 1 & r\le r_t
\frac12\left[1+\cos\left( \pi\frac{r-r_t}{r_{\mathrm{cut}}-r_t} \right)\right] & r_t<r<r_{\mathrm{cut}}
0 & r\ge r_{\mathrm{cut}} \end{cases} ]

を使用します。

除去する low-r 成分は

[ g_{\mathrm{loss}}(r)=W(r)g(r) ]

で、残す側は

[ g_{\mathrm{after}}(r)=g(r)-W(r)g(r) ]

です。

除去成分を逆変換して

[ i_{\mathrm{diff}}(Q) ]

を求め、

[ i_{\mathrm{bFT}}(Q) = i(Q)-i_{\mathrm{diff}}(Q) ]

とします。

webPDF の moving-average smoothing に相当する処理は入っていません。

7.3 PDF への影響

両者は「非物理的な short-r 成分を Q-space correction として戻す」という目的は同じですが、

• cutoff 位置の決め方
• hard cutoff / cosine taper
• Q-space correction の smoothing 有無

が異なります。

したがって補正後 (i(Q))、補正後 (G(r))、Lorch (G(r)) は同一にはなりません。

---

8. Lorch 処理へ入力する i(Q) が異なる

Lorch window 自体は両者とも基本的に

[ L(Q)=\frac{\sin(Q\Delta)}{Q\Delta}, \qquad \Delta=\frac{\pi}{Q_{\max}} ]

です。

そして

[ G_{\mathrm{Lorch}}(r) = \frac{2}{\pi} \int Q\,i(Q)L(Q)\sin(Qr)\,dQ ]

を計算します。

しかし入力 (i(Q)) が異なります。

webPDF

Rcut_and_bFT() → moving average correction → UseSmoothBG() の後の

[ i_{\mathrm{smooth}}(Q) ]

へ Lorch window を掛けます。

RustPDF

cosine-taper で得た loss を直接差し引いた

[ i_{\mathrm{bFT}}(Q)=i(Q)-i_{\mathrm{diff}}(Q) ]

へ Lorch window を掛けます。

したがって Lorch window の式がほぼ同じでも、窓を掛ける前の物理的な Q-space signal が違う ため、G(r)_Lorch も異なります。

---

9. 同じ、または本質的差ではない部分

Q の定義

両者とも

[ Q=\frac{4\pi}{\lambda}\sin\theta ]

です。webPDF は 12.3984 を用いた波長換算、RustPDF は物理定数 (h,c,e) から換算します。 これは概念上の物理式の差ではなく、定数表現の微小な数値差です。

polarization correction

両者とも

[ P=p+(1-p)\cos^2(2\theta) ]

に対応する係数で強度を除算します。RustPDF はその逆数を correction factor として掛けています。

Compton recoil

webPDF は組成平均した Compton 項を

[ C(Q)^{\mathrm{recoil}} ]

とします。RustPDF も calc_compton() の後に recoil_correction() を pipeline で適用しています。 したがって、ここは現在の標準 pipeline では本質的な相違ではありません。

G(r) と g(r) の関係

両者とも

[ g(r)=1+\frac{G(r)}{4\pi r\rho} ]

です。Rust の reduced_gr は webPDF の _Gr に対応する (G(r)) です。

---

10. 最終 PDF に差を生む経路のまとめ

                         webPDF                         RustPDF
                         ------                         -------

complex f(Q)        Re[f^2]                       |f|^2
                    Re[<f>^2]                     |<f>|^2
                         |                              |
                         +----------> normalization <---+
                                      |
                       mean(I/D)     vs     sum(I)/sum(D)
                                      |
                                     i(Q)
                                      |
                 interpolation order / Q-grid differs
                                      |
                    DST-I             vs          Zoom-DST
                                      |
                                     G(r)
                                      |
       zero-crossing + hard cut       vs      fixed cosine taper
                                      |
       back transform + moving avg    vs      direct back correction
                                      |
                              corrected i(Q)
                                      |
                         Lorch transform
                                      |
                               G_Lorch(r)

実際に LocalPDF と RustPDF の .dat を比較したときに差を追跡する場合は、 まず G(r) の Fourier engine だけを見るのではなく、

1. (\langle

   f

   ^2\rangle) と (

   \langle f\rangle

   ^2)

2. normalization factor
3. raw (i(Q))
4. short-r correction 前後の (i(Q))
5. Lorch 入力 (i(Q))

の順に比較する必要があります。

---

11. 参照した実装

webPDF / localPDF

• calc_pdf.py
  • CalcComp()
  • CalcForm()
  • CorrPol()
  • CalcFseries()
  • Calc_iQ_SQ()
  • Interpolate()
  • CalcGrFFT()
  • Rcut_and_bFT()
  • CalcGrLorch()
  • subtract_bg()
• localPDF は webPDF の Calc_PDF_from_xy を bridge 経由で使用するため、 PDF の物理計算についてはこの calc_pdf.py が基準。

RustPDF

• src/compute/common.rs
  • subtract_bg()
  • make_q()
  • make_pol_from_q()
  • calc_fseries()
  • calc_compton()
  • recoil_correction()
  • raw_q_process()
  • interpolate()
  • calc_iq_sq()
• src/pipeline.rs
  • standard pipeline の処理順
• src/compute/zdst.rs
  • calc_gr()
  • cut_gr()
  • calc_back_iq()
  • lorch_window()
  • calc_ggr_lorch()
• src/zoom_dst.rs
  • zoom_dst() 内で sine-transform 離散和へ delta_q を乗算

---

12. 要点

PDF の物理結果に関して、現在のコードから明確に確認できる重要な差は次の通りです。

1. 複素原子散乱因子の二乗: webPDF は Re(f^2)、RustPDF は物理的 modulus square |f|^2。
2. 強度規格化: webPDF は ratio の平均、RustPDF は sum の ratio。
3. background の負値: webPDF は 0 clipping、RustPDF は保持。
4. Q-grid/補間順序: webPDF は (i(Q)) 算出後、RustPDF は強度段階で先に補間。
5. Fourier sine transform の離散化: 目標式は同じだが、DST-I と積分重み付き Zoom-DST で有限 grid の数値表現が異なる。
6. short-r correction: webPDF は zero-crossing/hard cutoff + moving average、RustPDF は fixed cosine taper + direct correction。
7. Lorch 入力: window 自体は同系統だが、その前の corrected (i(Q)) が異なる。

これらはすべて、最終的な (G(r)) または (G_{\mathrm{Lorch}}(r)) に伝播する差です。