本篇延續 PEC 球體的米氏散射解析解與 Sim4Life 數值模擬驗證(上) 一文

3.3.2 三維散射分布與數據匯出
 

令 $j=\sqrt{-1}$，由實部與虛部重建
 

$$E_{\theta }\left(\varphi \right)=Re{E}_{\theta }\left(\varphi \right)+jIm{E}_{\theta }\left(\varphi \right)$$

$$E_{\varphi }\left(\varphi \right)=Re{E}_{\varphi }\left(\varphi \right)+jIm{E}_{\varphi }\left(\varphi \right)$$
兩個分量對應互相正交的球座標方向，因此總切向電場幅度為
 

$$A\left(\varphi \right)=\sqrt{{\left|E_{\theta }\left(\varphi \right)\right|}^2+{\left|E_{\varphi }\left(\varphi \right)\right|}^2}$$

本範例在整個 $0^\circ$ 至 $360^\circ$ 範圍內，$\max|E_\theta|$ 約為 $1.72\times10^{-7}\ \mathrm{V/m}$，而 $\max|E_\phi|$ 約為 $1.94013\ \mathrm{V/m}$，顯示此切面由 $E_\phi$ 主導。在 $90^\circ$ 至 $270^\circ$ 的比較區間內，總幅度與單獨 $E_\phi$ 所得到的正規化曲線差異小於 $10^{-12}\ \mathrm{dB}$，在本文繪圖精度下可以忽略。
 

在後續解析解與模擬數值結果的比較中，Python 腳本預設設定 COMPONENT = total，代表使用總切向幅度進行比較；若只使用 $E_\phi$，可改為 COMPONENT = phi。這兩種方式在本資料中幾乎重合，但更換場源極化或觀測切面時，不能假設此關係仍成立。
 

3.3.3 米氏解析解與角度取樣對齊
 

米氏散射角以入射波傳播方向為基準，前向方向為 $\theta_s=0^\circ$。依據本範例的座標設定，在比較區間內使用
 

$${\theta }_s=|\varphi -{180}^{\circ }|$$
解析解採用第 2.5 節說明的 PEC 球體 $S_2$ 散射振幅，此選擇對應本案例的場源極化與方位角切面。
 

由模擬結果檔案 ExportedData.xlsx 的完整角度資料中，擷取 $90^\circ\le\phi\le270^\circ$ 的範圍，共取得 181 個樣本點。腳本程式將直接在這 181 個角度位置上計算 $S_2$，不以插值方式取代解析解的計算。預設繪圖中的兩條比較曲線亦採用相同的 $1^\circ$ 角度間隔進行取樣連線，以保留與原始模擬資料一致的角度解析度。
 

此外，預設設定 ENSE_MIE_STEP_DEG = None，表示解析解曲線沿用模擬數據的原始角度取樣。若希望觀察更細緻的解析零點位置，可將 DENSE_MIE_STEP_DEG 設為 0.1，此時僅將米氏解析解的繪圖網格加密至 $0.1^\circ$。此設定不會增加 Sim4Life 的實際模擬資料量，也不會改變以原始 181 個樣本點計算的誤差指標。
 

需要注意的是，當解析解的角度網格加密後，可能會觀察到更深、更接近理論零點的數值。這是由角度取樣解析度提高所造成的結果，並不代表模擬結果的準確度變差。
 

3.3.4 前向散射值正規化
 

ExportedData.xlsx 數據檔第 182 列對應 $\phi=180^\circ$，其中
 

$$E_\phi(180^\circ) = 1.491144896+j\,1.241198659\ \mathrm{V/m}.$$

因此參考幅度約為
 

$$A(180^\circ)=1.940125566\ \mathrm{V/m}.$$

將 Sim4Life 模擬數據的正規化曲線定義為
 

$$L_{\mathrm{S}\mathrm{4}\mathrm{L}}\left(\varphi \right)=20{log}_{10}\frac{A\left(\varphi \right)}{A\left({180}^{\circ }\right)}$$
也可先定義
 

$$Q={\left|E_{\theta }\right|}^2+{\left|E_{\varphi }\right|}^2$$

再計算
 

$$L_{\mathrm{S}\mathrm{4}\mathrm{L}}\left(\varphi \right)=10{log}_{10}\frac{Q\left(\varphi \right)}{Q\left({180}^{\circ }\right)}$$
其中 $Q$ 是場幅平方，不是已校正的絕對 RCS。米氏解析解則計算
 

$$L_{\mathrm{M}\mathrm{i}\mathrm{e}}\left(\varphi \right)=20{log}_{10}\frac{\left|S_2\left({\theta }_s\left(\varphi \right)\right)\right|}{\left|S_2\left(0\right)\right|}$$
兩條曲線分別使用自身的前向值正規化。此方式比較的是相對角度響應；即使兩條曲線在 $180^\circ$ 均為 0 dB，也不代表絕對散射幅度已通過驗證。若需比較 dBsm，還必須確認 Sim4Life 輸出的入射場參考、距離因子及 RCS 校正定義，不能僅由這份 V/m 資料直接推算。
 

3.3.5 執行 Python 腳本與圖形輸出
 

本節使用 mie_pec_sphere_scripter.py 腳本來比較 PEC 球體的米氏散射解析解與 Sim4Life EM-FDTD 模擬結果。腳本從 ExportedData.xlsx 讀取指定觀測切面的複數遠場電場，依教材設定重建模擬曲線，再計算 PEC 球體的米氏 S2 散射振幅。兩條曲線分別以前向散射方向作為 0 dB 參考，最後輸出圖形、逐點資料及誤差指標。

腳本按照下列順序執行:

• 載入 Python 套件
• 讀取使用者設定與物理常數

• 以 read_export() 讀取並驗證 Excel

• 以 mie_amplitudes() 計算 PEC 球體米氏散射振幅
• 以 power_to_db() 將線性功率正規化為 dB
• 以 analyze() 建立模擬與解析曲線並計算誤差
• 以 write_outputs() 輸出 CSV、JSON、PNG 與 SVG
• 由 main() 串接上述流程並完成執行

(1) 執行環境與與套件

腳本開頭匯入標準函式庫與第三方套件:

各套件的用途如下:

• SimpleNamespace: 將分散的設定集中成 args 物件，供各函式傳遞使用
• csv: 寫出逐角度比較資料
• hashlib: 計算原始 Excel 的 SHA-256，建立結果追溯資訊
• json: 保存誤差指標與執行設定
• Path: 管理 Windows 輸入與輸出路徑
• NumPy: 執行複數陣列、遮罩、排序及誤差統計
• openpyxl: 以唯讀方式讀取 Excel
• spherical_jn、spherical_yn: 計算球面 Bessel 函數及其導數，用於建立 PEC 米氏係數

在 Sim4Life Scripter 執行前，可先測試套件是否可正常匯入:

import numpy
import scipy
import matplotlib
import openpyx

若出現 ModuleNotFoundError，則必須將缺少的套件安裝至 Sim4Life Scripter 實際使用的 Python 環境，而不是其他獨立的 Python 環境。(可參考【Python API】如何在 Sim4Life 環境中安裝其他 Python 套件?)

(2) 使用者設定區

腳本以固定變數取代命令列參數，使其可直接在 Sim4Life Scripter 視窗中執行:

a.. 輸入與輸出路徑

INPUT_XLSX 指向 Sim4Life 匯出的 Excel，OUTPUT_DIRECTORY 則指定所有輸出檔案的位置。兩者使用絕對路徑。

b. 工作表與場分量

SHEET_NAME = None 表示由程式搜尋包含五個必要欄位的唯一工作表。若活頁簿內存在多張相同格式工作表，應明確指定: SHEET_NAME = EM RCS E(phi)
 

COMPONENT 可設為:

• total：使用總切向場幅平方 $|E_\phi|^2+|E_\theta|^2$
• phi：只使用 $|E_\phi|^2$

本案例中 $E_\phi$ 明顯遠大於 $E_\theta$，因此兩種設定所得的正規化曲線幾乎重合。但更換入射極化或觀測切面後，不能假設 $E_\theta$ 仍可忽略。

c. 物理參數與角度範圍

SPHERE_RADIUS_M 與 FREQUENCY_HZ 必須與 Sim4Life 模型一致，本例為半徑 0.127 m、頻率 10 GHz 的 PEC 球體。
PHI_MIN_DEG 至 PHI_MAX_DEG 定義比較範圍，FORWARD_PHI_DEG = 180.0 則指定前向散射方向。

DENSE_MIE_STEP_DEG 只控制解析曲線的繪圖密度:

• None: 解析曲線沿用 Excel 原始角度取樣
• 0.1: 以 0.1° 間隔重算更密集的解析曲線

加密解析曲線不會增加 Sim4Life 模擬資料，也不會改變誤差統計，因為誤差仍在 Excel 的原始角度樣本上計算。

PLOT_Y_MIN_DB 只設定圖形 y 軸下限，不會裁切一般 RMSE 或 MAE 的資料。
$$E_\phi=\operatorname{Re}(E_\phi)+j\operatorname{Im}(E_\phi).$$

(3)  read_export(): 讀取與驗證 Excel
 

read_export(path, sheet=None) 是資料入口。函式先以 read_only=True 開啟活頁簿，避免修改原始檔案，並以 data_only=True 讀取公式的快取數值。

函式依序執行下列檢查:

• 若使用者指定工作表，確認該名稱存在
• 檢查第一列是否完整包含五個必要欄名
• 確認符合格式的工作表只有一張
• 從第 2 列開始讀取數值，略過五欄全空的列
• 將每個儲存格轉換為浮點數
• 排除缺值、文字、NAN與無限大
• 確認至少存在兩筆資料
• 依方位角由小到大排序
• 拒絕重複或幾乎重複的角度樣本。

函式最後回傳一個字典:

{
    "sheet": 工作表名稱,
    "row": 原始 Excel 列號,
    "phi": 方位角陣列,
    "etheta": Etheta 複數陣列,
    "ephi": Ephi 複數陣列,
}

保留 Excel 原始列號的目的，是讓輸出的 CSV 可追溯回原始資料位置。即使讀取或驗證中途失敗，Finally 區塊仍會關閉活頁簿。

(4) mie_amplitudes(): 計算 PEC 球體米氏散射振幅

mie_amplitudes() 是解析計算的核心。函式輸入方位角、半徑、頻率、前向角與可選的最高階數，輸出 S1、S2、自由空間波數 $k$ 與實際階數 $N$。

a. 波數、尺寸參數與級數階數

自由空間波數與尺寸參數分別為: 
 

$$k=\frac{2\pi f}{c_0}$$
$$x=ka$$

若未指定 n_max，程式採用 Wiscombe 型估算:
 

$$N=\lceil x+4x^{1/3}+2\rceil$$
對本案例 $x\approx26.62$，所需最高階數約為 41。階數過低時，旁瓣數量、零點位置與前向主瓣寬度可能尚未收斂。
 

b. 球面 Bessel、Hankel 函數與 PEC 米氏係數
 

SciPy 提供第一、第二類球面 Bessel 函數 $j_n(x)$ 與 $y_n(x)$。程式將兩者組成第一類球面 Hankel 函數:
 

$$h_n^{\left(1\right)}\left(x\right)=j_n\left(x\right)+i{y}_n\left(x\right)$$

對 PEC 球體，程式使用:
 

$$a_n=-\frac{j_n+x{j}_n\prime }{h_n+x{h}_n\prime }$$
$$b_n=-\frac{j_n}{h_n}$$
第一式等價於 $-[xj_n(x)]'/[xh_n^{(1)}(x)]'$。

a_n 與 b_n 都是複數，包含各階部分波的振幅與相位資訊。

c. Sim4Life 方位角與物理散射角
 

赤道切面的角度關係為:
 

$$\mu=\cos\theta_s=\cos(\phi-\phi_{\mathrm{forward}}).$$

在本例中 $\phi_{\mathrm{forward}}=180^\circ$。使用餘弦可同時處理前向角左右兩側，並保持球體的對稱性。此轉換依賴 theta = 90° 切面；若改變觀測平面，必須重新確認角度定義。
 

d. 角向函數遞迴與 S1、S2
 

程式從
 

$$\pi_0=0,\qquad \pi_1=1$$

開始遞迴。第 $n$ 階的 $\tau_n$ 為:
 

$$\tau_n=n\mu\pi_n-(n+1)\pi_{n-1}.$$

接著以權重
 

$$\frac{2n+1}{n(n+1)}$$

累加兩個正交極化散射振幅:
 

$$S_1=\sum_{n=1}^{N}\frac{2n+1}{n(n+1)}(a_n\pi_n+b_n\tau_n),$$

$$S_2=\sum_{n=1}^{N}\frac{2n+1}{n(n+1)}(a_n\tau_n+b_n\pi_n).$$

最後利用
 

$$\pi_{n+1}=\frac{2n+1}{n}\mu\pi_n-\frac{n+1}{n}\pi_{n-1}$$

更新下一階角向函數。遞迴方式可避免直接計算伴隨 Legendre 函數時，在前向角附近因除以很小的$\sin\theta_s$ 而產生數值不穩定。
 

本範例的入射極化與觀測切面對應標準定義下的 S2，因此後續比較使用 $|S_2|^2$; S1 雖同時計算，但不納入本案例誤差指標。
 

(5) power_to_db(): 線性功率正規化


power_to_db(power, reference)使用:
 

$$L=10\log_{10}\left(\frac{P}{P_{\mathrm{ref}}}\right)$$

將線性功率轉換成相對 dB。參考值必須大於零且為有限數。函式不套用繪圖下限，因此若某個角度的功率精確為零，結果會自然成為 $-\infty$。

保留 $-\infty$ 比任意抬高理論零點更符合物理意義。後續 analyze() 會將非有限配對排除於一般 RMSE 與 MAE 之外，並在 JSON 中記錄排除筆數。

(6) analyze() 建立曲線與誤差指標
 

analyze() 將 Excel 資料、模型參數與比較範圍整合成完整分析結果。

a. 範圍與前向樣本檢查
 

函式首先確認:

• 下限小於上限
• 前向角位於比較範圍內
• Excel 實際涵蓋完整比較範圍
• Excel 中恰好存在一筆前向角樣本

腳本不會以插值或外插建立前向參考值。若缺少 $180^\circ$ 樣本，應重新匯出資料或修正角度設定。
 

b. 模擬功率與解析功率

選擇 COMPONENT = "phi" 時:

$$P_{\mathrm{s}\mathrm{i}\mathrm{m}}={\left|E_{\varphi }\right|}^2$$

選擇 COMPONENT = "total" 時:
 

$$P_{\mathrm{s}\mathrm{i}\mathrm{m}}={\left|E_{\varphi }\right|}^2+{\left|E_{\theta }\right|}^2$$
解析端則使用: 
 

$$P_{\mathrm{M}\mathrm{i}\mathrm{e}}={\left|S_2\right|}^2$$
兩條曲線必須分別使用自己的前向值正規化: 
 

$$L_{\mathrm{S}\mathrm{4}\mathrm{L}}\left(\varphi \right)=10{log}_{10}\frac{P_{\mathrm{s}\mathrm{i}\mathrm{m}}\left(\varphi \right)}{P_{\mathrm{s}\mathrm{i}\mathrm{m}}\left({180}^{\circ }\right)}$$
$$L_{\mathrm{M}\mathrm{i}\mathrm{e}}\left(\varphi \right)=10{log}_{10}\frac{{\left|S_2\left(\varphi \right)\right|}^2}{{\left|S_2\left({180}^{\circ }\right)\right|}^2}$$

c. 有限值篩選與誤差

逐點差值定義為:
 

$$\Delta L(\phi)=L_{\mathrm{S4L}}(\phi)-L_{\mathrm{Mie}}(\phi).$$

正值代表 Sim4Life 曲線高於解析曲線，負值代表 Sim4Life 曲線較低。一般統計只使用模擬與解析值皆為有限數的角度點。

腳本計算下列主要指標:

• rmse_db: 未裁切逐點差值的均方根誤差
• mae_db: 未裁切逐點差值的平均絕對誤差
• max_abs_error_db: 最大絕對誤差
• max_error_phi_deg: 最大誤差所在方位角
• excluded_nonfinite_pairs: 因散射零點而排除的非有限配對數
• max_total_vs_phi_normalized_difference_db: 總切向場與單獨 $E_\phi$ 曲線的最大差異

此外，程式提供一個獨立的穩健指標: 先把兩條曲線都限制在不低於 −40 dB，再計算 RMSE。此指標可降低極深零點對單一統計量的支配，但不取代未裁切 RMSE，也不改變圖形資料。
 

d. 米氏級數收斂檢查
 

腳本以原階數 $N$ 與增加十階的 $N+10$ 分別計算 S2，再求:
 

$$\varepsilon_N= \frac{\max_\phi|S_2^{(N+10)}(\phi)-S_2^{(N)}(\phi)|} {|S_2^{(N)}(180^\circ)|}.$$

此值記錄於 n_plus_10_relative_amplitude_change。數值越小，表示加入更多高階模式後結果幾乎不變，米氏級數截斷越充分。
 

(7) write_outputs()：輸出資料與圖形


write_outputs() 使用 Matplotlib 的物件導向介面及離屏 Agg 畫布，不切換 Sim4Life 主程式的 pyplot 後端，因此適合在 Scripter 中重複執行。

a. figure3_comparison.csv

CSV 保存比較區間內每一筆資料，包括: 

• 原始 Excel 列號與方位角
• $E_\theta$、$E_\phi$ 的實部與虛部
• Sim4Life 正規化曲線
• 米氏 S2 正規化曲線
• Sim4Life − Mie 的逐點差值

此檔案適合用於後續繪圖、抽查特定零點或追溯原始 Excel 儲存格。
 

b. comparison_metrics.json

JSON 保存模型參數、資料範圍、角度間距、參考 Excel 列、前向複數場、誤差指標與級數收斂結果。程式也會計算輸入 Excel 的 SHA-256；若兩次報告的雜湊值不同，表示使用的原始檔案內容不同。
 

c. PNG 與 SVG 圖形
 

腳本輸出:

• figure3_mie_vs_sim4life.png:  200 dpi 畫素圖檔
• figure3_mie_vs_sim4life.svg: 向量格式圖檔

解析曲線使用藍色實線，Sim4Life 曲線使用橘色點畫線。圖形標題顯示球體半徑與工作頻率；y 軸標示為正規化雙站 RCS，避免與絕對 dBsm 混淆。
 

若啟用 DENSE_MIE_STEP_DEG，程式只在圖形中加密米氏解析曲線。前向參考值會另外在精確的 FORWARD_PHI_DEG 計算，不依賴密集網格是否恰好包含該角度。

(8) main(): 主程式與防呆檢查

main()在執行分析前檢查

• 輸入及輸出是否為絕對路徑
• Excel 是否存在
• 輸出位置是否為資料夾
• COMPONENT 是否為 total 或 phi
• DENSE_MIE_STEP_DEG 是否為 None 或正值

通過檢查後，設定會被包裝成 SimpleNamespace:

args = SimpleNamespace(
    ......
)

主流程依序呼叫:

data = read_export(args.xlsx, args.sheet)

result = analyze(
    data,
    args.radius,
    args.frequency,
    args.component,
    args.phi_min,
    args.phi_max,
    args.forward_phi,
)

write_outputs(data, result, args)

檔案最後使用:

if __name__ == "__main__":
    main()

因此，直接執行腳本時會啟動分析; 若日後由其他 Python 程式匯入此檔，則不會自動執行 main()。

3.3.6 比較結果與誤差解讀

figure3_mie_vs_sim4life.png 圖示:

本案例 $a=0.127\ \mathrm{m}$、 $f=10\ \mathrm{GHz}$ 及 $S_2$ 設定的計算，尺寸參數為 $ka=26.617231779$，截斷階數為 $N=41$。兩條曲線在原始 181 個角度上的比較結果如下:
 

其中一般 RMSE 定義為
 

$$\mathrm{RMSE}=\sqrt{\frac{1}{M}\sum_{i=1}^{M}\left[L_{\mathrm{S4L}}(\phi_i)-L_{\mathrm{Mie}}(\phi_i)\right]^2},\qquad M=181.$$

限制至 −40 dB 的指標屬於額外敏感度分析: 先將兩條曲線各自取 $\max(L,-40)$，再計算 RMSE，不取代未截斷指標。181 組比較值皆為有限值，沒有因無效資料而排除角度。
 

由圖可觀察到，兩者的前向主瓣與多數旁瓣結構相近，但部分局部極小值仍有可辨識差異。最大誤差出現在 171°，其米氏值約為 −24.97 dB，並非僅發生於低於 −40 dB 的深零點。因此，不能將全部誤差都歸因於極低場值的對數放大。
 

將級數由 41 階提高至 51 階後，各角度振幅變化相對於前向振幅的最大值約為 $1.8\times10^{-11}$，顯示本次差異不是由米氏級數截斷不足主導。不過，僅憑單一 Excel 匯出仍無法判定 FDTD 誤差來源或宣稱網格已收斂; 仍應分別檢查網格加密、模擬時間與遠場角度取樣。
 

因此，本文的比較不只是在同一張圖上疊加兩條曲線，而是依序確認模型條件、觀測角度、極化分量、正規化方式與取樣網格的一致性。完成這些核對後，剩餘差異才適合用來討論數值離散與模擬收斂; 也不能因為前向值同為 0 dB，就宣稱所有角度或絕對 RCS 均已驗證。

4. 總結

本文以理想電導球體為例，建立由散射物理、解析理論到數值驗證的完整流程。第一章介紹雙站雷達截面的物理意義與觀測方向，說明散射強度如何隨接收位置改變; 第二章闡述米氏散射以球面波展開電磁場，並透過理想電導體的表面邊界條件求得散射係數，呈現不同階次散射波疊加形成的角度分布；第三章將理論落實於 Sim4Life 的幾何建模、平面波激發、EM-FDTD 模擬與後處理，利用匯出的複數電場及 Python 腳本完成解析解比較。可靠的驗證須確保模型參數、極化、觀測角度、取樣與正規化方式一致，並配合網格及時間收斂檢查。以前向值正規化的曲線可用於評估相對散射分布，但不能取代絕對 RCS 與相位的驗證。