本ページは、本書第10章〜第12章の 1/3 オクターブ受音点レベル L_p(Z)・L_p(A)・L_p(C)・L_p(G) および 重み付け総合値を算出した Python ソースコード一式を整理したものです。各ファイルは Python 標準ライブラリのみを使用しており、NumPy・SciPy 等の外部依存はありません。
計算は次の 3 系統からなります。
iso_compute.pycompute_weather.pycompute_cylindrical.pyPython 3 のみで動作します。python/ ディレクトリに移動して各メインスクリプトを直接実行してください。
cd python
python3 iso_compute.py # ISO 9613-2 要約、9 ケース
python3 iso_compute.py --per-band # 1/3 オクターブ詳細
python3 compute_weather.py # ケース1〜4 気象補正、9 × 4 = 36 総合値
python3 compute_weather.py --per-band # 1/3 オクターブ詳細
python3 compute_cylindrical.py # 円筒拡散による上限値、9 ケース
python3 iso_emit_tables.py # 付録の標準参照テーブル(α、A/C/G/Z 重み)を LaTeX 形式で出力
各ファイルの目的と本書中の位置づけを以下に示す。ソースコードのリンクから生のスクリプトを参照できる。
iso9613_1.pyISO 9613-1:1993 の大気吸収係数 α(f, T, RH, p_a) を解析公式で計算する。第6.2節の式 (3)〜(5) をそのまま実装。
iso9613_2.pyISO 9613-2 の点音源屋外伝搬モデル。幾何拡散 A_div、大気吸収 A_atm、地表効果 A_gr を算出する。地表項は海域・陸域の2 セグメントに分けて処理する。
iso_weightings.py1/3 オクターブバンド中心周波数における周波数重み付けカーブ:IEC 61672-1:2013 の A 特性・C 特性、ISO 7196:1995 の G 特性。
iso_vestas_spectrum.pyVestas V236-15 MW の 1/3 オクターブ音響パワースペクトル(Z・A)。エストニア Hiiumaa EIA(Ramboll、2022)で採用された安全マージン +2 dB を含む(L_WA = 120 dB(A), L_WZ = 134 dB)。
iso_compute.pyISO 9613-2 ベースラインのメインスクリプト。9 ケース(3 配置 × 3 受音点)の要約と、オプション --per-band による 1/3 オクターブバンド詳細を出力する。表 12.2 を再現。
iso_emit_tables.py本書付録の標準参照テーブル(10 °C / 70 % および 20 °C / 80 % における α、ならびに A・C・G・Z 重み付け値)を LaTeX 形式で出力する。
weather_corrections.pyケース1〜4 の帯域別気象補正を実装する。屈折集束 ΔL_refrac(表 11.2.4)、ダクト内拡散 ΔL_duct(実効拡散指数 n_i、表 11.2.3・表 11.2.6)、および海面反射補正 ΔL_sea(表 11.2.5)の 3 補正からなる。
cylindrical_spreading.py第11.3節の円筒拡散モデル。ISO 9613-2 の球面拡散 20·log10(r) を、ダクト高 h = 150 m の円筒拡散 10·log10(2π r h) に置き換える。大気吸収(ISO 9613-1)と海面反射 +3 dB はそのまま用いる。
英文の総合 README(インストール条件・モデル定式化・検証結果):python/README.md
各風車 k および各 1/3 オクターブバンド i について、受音点における音圧レベル L_p,i,k を次式で算出する。
L_p,i,k = L_W,i − A_div(r_k) − A_atm,i(r_k) − A_gr,i(r_k)
A_div(r) = 20·log10(r) + 11 — 幾何拡散(ISO 9613-2 §7.1, eq. 7)。A_atm,i(r) = α_i · r — 大気吸収(§7.2, eq. 8)。α_i は ISO 9613-1 解析公式から計算。A_gr,i — 地表効果(§7.3)。海域セグメント(G = 0)では全帯域一律 −3 dB、陸域セグメントは帯域別の値を割り当てる:
f ≤ 200 Hz → A_gr,land = −2 dB(ほぼ反射的)200 Hz < f ≤ 1000 Hz → A_gr,land = 0 dBf > 1000 Hz → A_gr,land = +1 dBr_sea / r_total、r_land / r_total で線形混合する。
全風車の寄与を帯域ごとにエネルギー和で集約し、最後に重み付けを加える。
L_p,i = 10·log10 Σ_k 10^(L_p,i,k / 10)
L_p(W) = 10·log10 Σ_i 10^((L_p,i + W_i) / 10)
重み W_i は標準値:A・C は IEC 61672-1:2013 表 3(0.1 dB 単位)、G は ISO 7196:1995 表 2、Z は全帯域 0 dB。
ISO 9613-2 §7.3.1 注 10 では、海域上の気温逆転条件は適用外と明記されている。本書第11.2章では、ISO ベースライン値の上に 3 つの帯域別補正を重ねることで、この問題に対応する。
ΔL_refrac(表 11.2.4)— 風速勾配および温度勾配により下方に屈折する音響エネルギー。距離非依存・帯域別・ケース依存。ΔL_duct — 表面ダクトに捕捉された音は r^n(n < 2)で拡散する。実効拡散指数 n_i を表 11.2.3(海域)・表 11.2.6(陸域)から取り、海域・陸域セグメントを別々に計算したのち合算する。距離依存。ΔL_sea(表 11.2.5)— 強風時(ケース 1, 2)は海面散乱で +3 dB 反射が中高域でやや減少、無風時(ケース 3, 4)は鏡面反射に近づく。4 ケースは独立した 2 つの気象軸の組み合わせとして定義される。
3 補正はすべて風車ごとのエネルギー和の内側で適用されるため、ダクト項は風車ごとの実際の経路長を反映する。
L_p,i,k_corr = L_p,i,k_ISO + ΔL_refrac,i + ΔL_duct,i,k + ΔL_sea,i
L_p,i_corr = 10·log10 Σ_k 10^(L_p,i,k_corr / 10)
ISO 9613-2 の球面拡散 20·log10(r) を、表面ダクト高 h = 150 m の円筒拡散 10·log10(2π r h) に置き換えた物理的上限値を算出する。大気吸収(ISO 9613-1)と海面反射 +3 dB は維持される。
本書本文および付録に掲載されている計算結果は、すべて以下のスクリプト群の出力を 0.1 dB に丸めたものである。本ページから取得したコードを再実行することで、本書記載値と次の精度範囲で同一の数値が得られる。
iso_compute.py の出力は表 12.2 の L_p(Z・A・C・G)全 9 ケースに対応する(0.1 dB 丸めによる差 ±0.05 dB 以内)。L_p(Z) の出力は付録 ISO ベースラインのテーブル群(tab:app_*_winter, tab:app_*_3km, tab:app_*_5km, tab:app_C_coast)の数値と 0.0〜0.1 dB(付録の丸め)で対応する。compute_weather.py の出力は表 12.3 および付録 C の値、全 9 × 4 = 36 通りに対応する(同じく ±0.05 dB の丸め範囲)。compute_cylindrical.py の出力は表 12.4 の円筒拡散行(C 配置・海岸線:L_p(Z) = 83.0、L_p(A) = 62.2、L_p(C) = 79.0、L_p(G) = 86.6)に対応する。第11.3節で指定された 10 °C / 70 % 条件下の ISO 9613-1 α を用いる。すべての計算で次の条件を用いる(ISO 9613-2 参照条件)。
L_WA = 120 dB(A)、L_WZ = 134 dB。3 配置案(A・B・C)で同一スペクトルを用い、C 案 18 MW 級の出力差は本書第10章記載の +2 dB 補正で吸収する。