Swerlingターゲットモデルとレーダー検出確率の導出と実装

夜空を横切る航空機をレーダーで追尾していると、受信信号の強さがちらちらと揺らぐのに気づきます。機体までの距離も向きもほとんど変わっていないのに、あるスキャンでは強い反射が返ってきて、次のスキャンでは急に弱くなる。これは機器の故障でも測定ミスでもありません。航空機という物体が、無数のリベットや翼端やエンジンナセルといった「小さな反射体の集まり」でできているために、観測する角度がわずかに変わるだけで、それらからの反射波が強め合ったり弱め合ったりして、見かけの反射断面積(RCS, Radar Cross Section)が刻々と変動するのです。

この「RCSの揺らぎ」を確率分布としてモデル化したものが、1954年に Peter Swerling が提案した Swerlingターゲットモデル です。レーダーが目標を「検出できた」と判定する確率(検出確率 $P_d$)は、信号対雑音比 $\mathrm{SNR}$ だけでなく、このRCS変動の統計的性質に強く依存します。変動するターゲットは、変動しないターゲットに比べて、同じ平均SNRでも検出が難しくなる場合があり、その差は「変動損失(fluctuation loss)」と呼ばれ、しばしば数 dB から十数 dB にも達します。

Swerlingモデルを理解すると、次のような実務に直結する見通しが得られます。

  • レーダーの探知距離設計: 「この航空機を90%の確率で検出するには送信電力をどれだけ確保すべきか」という所要SNRの見積もりに、ターゲット変動を織り込めるようになります
  • パルス積分の設計: 複数パルスを積分すると検出性能が上がりますが、その効き方はSwerlingモデルごとに大きく違います。何パルス積分すれば十分かの判断材料になります
  • 対ステルス・対UAV評価: 形状の異なる目標(大型機・小型ドローン・ミサイル)はRCS変動の統計も異なります。どのSwerlingモデルが当てはまるかで探知性能を評価し直せます

本記事の内容

  • レーダー検出を「雑音の中の信号判定」という統計問題として定式化する
  • 非変動ターゲット(Swerling 0 / Marcum)の検出確率を導出する
  • 多数散乱体モデルからレイリー振幅分布と指数分布RCSを導く
  • Swerling I〜IVのRCS確率分布(指数分布・自由度4のカイ二乗分布)を導出する
  • 一定誤警報率 $P_{fa}$ のもとで各モデルの $P_d$ を $\mathrm{SNR}$ の関数として求める
  • 変動損失とパルス積分効果をPythonの検出確率曲線で定量的に読み取る

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

レーダー検出を統計問題として捉える

まず「レーダーが目標を検出する」とは何をしているのか、その本質を直感で押さえておきましょう。レーダーはアンテナから電波を放ち、目標で反射して戻ってきた微弱な信号を受信機で受け取ります。ところがこの受信機の出力には、目標からの信号だけでなく、いつでも避けられない 熱雑音 が乗っています。熱雑音は受信機内部の抵抗体で電子がランダムに揺らぐことから生まれる、完全にランダムな電圧です。

したがってレーダーが見ているのは「信号 + 雑音」か「雑音だけ」かのどちらかであり、両者を区別するために、出力電圧の大きさがある しきい値 を超えたかどうかで「目標あり/なし」を判定します。ここで2種類の確率が登場します。

  • 検出確率 $P_d$: 本当に目標があるとき、出力がしきい値を超えて「目標あり」と正しく判定する確率
  • 誤警報確率 $P_{fa}$: 目標がない(雑音だけの)とき、雑音がたまたま大きくなってしきい値を超え、「目標あり」と誤判定する確率

この2つはトレードオフの関係にあります。しきい値を下げれば $P_d$ は上がりますが、雑音だけでもしきい値を超えやすくなるので $P_{fa}$ も上がってしまいます。レーダー設計では、まず許容できる $P_{fa}$(典型的には $10^{-6}$ など)を決めてしきい値を固定し、その条件下で $P_d$ を最大化(あるいは所要SNRを最小化)するという考え方をとります。この「$P_{fa}$ を一定に保つ」という発想を実際の信号処理で実現したものが、前提記事で扱った CFAR検出 です。

ここで本記事の主役が登場します。目標からの信号の強さ、すなわちSNRが固定値なら話は単純ですが、冒頭で見たように実際のRCSは揺らぎます。RCSが揺らぐと信号電力も揺らぎ、ある瞬間にはしきい値を超え、別の瞬間には超えない、という確率的な振る舞いが生じます。この揺らぎ方を統計分布で表したのがSwerlingモデルであり、それによって $P_d$ の計算式が変わってくるのです。

それでは、まず揺らぎのない「基準ケース」である非変動ターゲットの検出確率から、丁寧に組み立てていきましょう。

受信信号の包絡線と雑音の分布

検出確率を導くには、受信機出力の電圧がどんな確率分布に従うかを知る必要があります。レーダー受信機は、受信した高周波信号を直交検波して同相成分 $I$ と直交成分 $Q$ の2つに分け、その大きさ(包絡線)$r = \sqrt{I^2 + Q^2}$ を取り出します。検出判定はこの包絡線 $r$ をしきい値 $V_T$ と比較して行います。

雑音だけのとき — レイリー分布

目標がなく雑音だけのとき、$I$ と $Q$ はそれぞれ平均0・分散 $\psi$ の独立なガウス分布に従います。これは熱雑音が多数の独立な微小寄与の和(中心極限定理)であることの帰結です。2つの独立な正規分布変数の二乗和の平方根は レイリー分布 に従うことが知られており、包絡線 $r$ の確率密度は次のようになります。

$$ p_n(r) = \frac{r}{\psi} \exp\left(-\frac{r^2}{2\psi}\right), \quad r \geq 0 $$

ここで $\psi$ は1直交成分あたりの雑音電力(分散)です。この分布から、包絡線がしきい値 $V_T$ を超える確率、すなわち誤警報確率が計算できます。$V_T$ から $\infty$ まで積分すると、

$$ P_{fa} = \int_{V_T}^{\infty} \frac{r}{\psi} \exp\left(-\frac{r^2}{2\psi}\right) dr = \exp\left(-\frac{V_T^2}{2\psi}\right) $$

この積分は被積分関数が $\exp(-r^2/2\psi)$ の微分の形になっているため、すぐに実行できます。$P_{fa}$ について逆に解くと、所望の誤警報確率を達成するしきい値が決まります。

$$ \frac{V_T^2}{2\psi} = \ln\frac{1}{P_{fa}} $$

ここで重要なのは、しきい値を雑音電力 $\psi$ で規格化して扱うと見通しがよくなる、という点です。そこで規格化しきい値を $\eta_T = V_T^2 / (2\psi)$ と定義すれば、誤警報確率は実にシンプルに $P_{fa} = e^{-\eta_T}$ となり、逆に $\eta_T = \ln(1/P_{fa})$ で必要なしきい値が決まります。たとえば $P_{fa} = 10^{-6}$ なら $\eta_T = \ln(10^6) \approx 13.8$ です。

信号があるとき — ライス分布

目標があるとき、包絡線は「信号成分 + 雑音」になります。振幅 $A$ の確定的な信号に雑音が重畳すると、包絡線 $r$ は ライス分布(Rician distribution) に従います。

$$ p_s(r) = \frac{r}{\psi} \exp\left(-\frac{r^2 + A^2}{2\psi}\right) I_0\left(\frac{rA}{\psi}\right) $$

ここで $I_0(\cdot)$ は0次の第1種変形ベッセル関数です。$A = 0$(信号なし)とすると $I_0(0) = 1$ なので、これは先ほどのレイリー分布に一致します。信号の存在によって分布のピークが $r = 0$ から正の方向へ押し出される、というのがライス分布の直感的なイメージです。

ここで信号対雑音比を $\chi = A^2 / (2\psi)$ と定義します。これは「信号電力 $A^2/2$ を1直交成分あたりの雑音電力 $\psi$ で割ったもの」に相当し、本記事を通じてSNRを表す中心的な量です。これらの分布が揃ったので、次は非変動ターゲットの検出確率を導きましょう。

非変動ターゲット(Swerling 0)の検出確率

最初に扱うのは、RCSが全く揺らがない理想的なターゲットです。これは Swerling 0(または Swerling V、あるいは Marcum ケース)と呼ばれ、すべての変動モデルを評価する際の基準になります。物理的には、完全な球体のように向きを変えても反射が一定な物体や、変動が無視できるほど短い時間スケールでの観測に相当します。

検出確率は、信号がある場合の包絡線分布 $p_s(r)$ をしきい値 $V_T$ から $\infty$ まで積分したものです。

$$ P_d = \int_{V_T}^{\infty} \frac{r}{\psi} \exp\left(-\frac{r^2 + A^2}{2\psi}\right) I_0\left(\frac{rA}{\psi}\right) dr $$

この積分を計算しやすくするため、変数を規格化します。$v = r / \sqrt{\psi}$ と置換し、SNRを $\chi = A^2/(2\psi)$、規格化しきい値を $\eta_T = V_T^2/(2\psi)$ とすると、上の積分は次の形に書き換えられます。途中、$dr = \sqrt{\psi}\, dv$ と $rA/\psi = v\sqrt{2\chi}$ を代入すると、

$$ P_d = \int_{\sqrt{2\eta_T}}^{\infty} v \exp\left(-\frac{v^2}{2} – \chi\right) I_0\left(v\sqrt{2\chi}\right) dv $$

この積分は初等関数では書けませんが、MarcumのQ関数 $Q_M(\alpha, \beta)$ という特殊関数を使うと閉じた形で表せます。MarcumのQ関数は次のように定義されます。

$$ Q_M(\alpha, \beta) = \int_{\beta}^{\infty} v \exp\left(-\frac{v^2 + \alpha^2}{2}\right) I_0(\alpha v)\, dv $$

この定義と先ほどの積分を見比べると、$\alpha = \sqrt{2\chi}$、$\beta = \sqrt{2\eta_T}$ と対応づけられます。したがって非変動ターゲットの検出確率は、

$$ \boxed{P_d = Q_M\left(\sqrt{2\chi},\ \sqrt{2\eta_T}\right)} $$

と簡潔に表せます。ここで $\eta_T = \ln(1/P_{fa})$ で誤警報確率と結びついているので、$P_{fa}$ と $\chi$(SNR)を与えれば $P_d$ が一意に決まります。

直感的に言うと、$Q_M(\alpha, \beta)$ は「中心が $\alpha$ だけずれた2次元ガウス分布が、半径 $\beta$ の円の外側にある確率」です。信号が強い($\chi$ が大きい)ほど分布全体が外側へ押し出され、円の外に出る確率すなわち $P_d$ が増えます。これがレーダー検出の最も基本的な定量関係です。MarcumのQ関数は SciPy の scipy.stats.ncx2(非心カイ二乗分布)や scipy.special.chndtr を通じて数値計算できるので、後ほどPythonで実装します。

この非変動ケースを基準として、いよいよRCSが揺らぐ場合に進みます。揺らぎを表すには、まず「なぜ揺らぐのか」を物理から導く必要があります。

RCS変動はなぜ生じるか — 多数散乱体モデル

航空機のような複雑な形状の目標は、機首・主翼・尾翼・エンジン・無数の突起など、たくさんの「散乱体(scatterer)」の集まりと見なせます。レーダーから見た合成RCSは、これらすべての散乱体からの反射波をベクトル的に足し合わせたものです。各散乱体までの距離は波長 $\lambda$(数センチ)のオーダーでわずかに異なるため、反射波の位相はバラバラで、互いに強め合ったり打ち消し合ったりします。

目標の姿勢が少し変われば、各散乱体までの距離が波長の何分の一か変化し、位相関係が一変します。その結果、合成RCSが大きく揺らぐのです。これを数式で表しましょう。$N$ 個の散乱体からの反射波の複素振幅の和を $\bm{E}$ とすると、

$$ \bm{E} = \sum_{k=1}^{N} a_k e^{j\phi_k} $$

ここで $a_k$ は各散乱体の反射振幅、$\phi_k$ は位相です。多数の散乱体があり、位相 $\phi_k$ が $[0, 2\pi)$ に一様にランダムに分布していると仮定します。この和の実部と虚部はそれぞれ、多数の独立なランダム量の和なので、中心極限定理により平均0のガウス分布に近づきます。

実部 $E_I = \sum_k a_k \cos\phi_k$ と虚部 $E_Q = \sum_k a_k \sin\phi_k$ が独立な平均0・等分散のガウス変数になると、その包絡線(振幅)$|\bm{E}| = \sqrt{E_I^2 + E_Q^2}$ は レイリー分布 に従います。これは雑音の包絡線がレイリー分布になったのと全く同じ数学的構造です。

そして RCS は振幅の二乗 $\sigma \propto |\bm{E}|^2$ に比例します。レイリー分布する変数の二乗は 指数分布 に従うことが知られています。すなわち、多数の同程度の散乱体からなる目標のRCSは、平均 $\bar\sigma$ の指数分布

$$ p(\sigma) = \frac{1}{\bar\sigma} \exp\left(-\frac{\sigma}{\bar\sigma}\right), \quad \sigma \geq 0 $$

に従うのです。これがSwerling I・II の基礎になる分布です。

一方、目標の中に1つだけ飛び抜けて大きな支配的散乱体(dominant scatterer)があり、その周りに多数の小さな散乱体がある場合(例えばミサイルの先端や大きな平面)には、分布の形が変わります。この場合のRCS分布は 自由度4のカイ二乗分布 で近似され、

$$ p(\sigma) = \frac{4\sigma}{\bar\sigma^2} \exp\left(-\frac{2\sigma}{\bar\sigma}\right), \quad \sigma \geq 0 $$

となります。これがSwerling III・IV の基礎です。指数分布(自由度2のカイ二乗分布)に比べてピークがはっきりしており、極端に小さなRCSになる確率が下がる、つまり「揺らぎが穏やか」な分布です。

これで2つの基礎分布が揃いました。次に、揺らぎの「速さ」の違いによって、これら同じ分布が4つのモデルに枝分かれする様子を見ていきます。

Swerling I〜IVの分類

Swerlingモデルは、上で導いた2種類のRCS分布(指数分布か自由度4カイ二乗分布か)と、揺らぎの 時間スケール(相関時間)の2つの軸で4つに分類されます。揺らぎの時間スケールには2種類あります。

  • スキャン間変動(slow fluctuation): 1回のスキャン(目標に向けたビーム照射)の間ではRCSは一定で、スキャンごとに独立に変わる。レーダーがアンテナを回して再びその目標を照射するまでに姿勢が変わるイメージ
  • パルス間変動(fast fluctuation): 1回のスキャン内で送る各パルスごとにRCSが独立に変わる。目標の姿勢変化や回転が速い、あるいはパルス繰り返し周期が長い場合

この時間スケールは、複数のパルスを積分(コヒーレント/ノンコヒーレント積分)して検出する際に決定的な違いを生みます。スキャン間変動では $n$ パルスがすべて同じRCS値を共有するのに対し、パルス間変動では $n$ パルスがそれぞれ独立なRCS値を持つため、平均化(ダイバーシティ)の効果が効いて検出性能が改善します。

この2軸を組み合わせると、次の表のように4つのモデルが定義されます。

モデル RCS分布 物理的描像 変動の速さ
Swerling I 指数分布(自由度2 $\chi^2$) 多数の同程度の散乱体 スキャン間(遅い)
Swerling II 指数分布(自由度2 $\chi^2$) 多数の同程度の散乱体 パルス間(速い)
Swerling III 自由度4の $\chi^2$ 1つの支配的散乱体 + 多数の小散乱体 スキャン間(遅い)
Swerling IV 自由度4の $\chi^2$ 1つの支配的散乱体 + 多数の小散乱体 パルス間(速い)

ここで「自由度」という言葉が出てきました。指数分布は自由度2のカイ二乗分布、もう一方は自由度4のカイ二乗分布です。一般に、RCS分布を自由度パラメータ $2K$(あるいは形状パラメータ $K$)を持つガンマ分布族として書くと統一的に扱えます。

$$ p(\sigma) = \frac{K^K}{(K-1)!\,\bar\sigma}\left(\frac{\sigma}{\bar\sigma}\right)^{K-1} \exp\left(-\frac{K\sigma}{\bar\sigma}\right) $$

Swerling I/II は $K = 1$(指数分布)、Swerling III/IV は $K = 2$(自由度4カイ二乗分布)に対応します。$K \to \infty$ の極限では分布がデルタ関数に近づき、非変動ケース(Swerling 0)に一致するという美しい統一が成り立ちます。$K$ が大きいほど揺らぎが小さく、$K=1$ が最も激しく揺らぐケースです。

分類が整理できたので、次は実際に各モデルの検出確率を導出します。鍵になるのは、SNRがRCSに比例して揺らぐとき、$P_d$ をRCS分布で平均化するという操作です。

単一パルスでの検出確率の導出

平均化の考え方

非変動ケースでは、与えられたSNR $\chi$ に対して検出確率 $P_d(\chi) = Q_M(\sqrt{2\chi}, \sqrt{2\eta_T})$ でした。ターゲットが変動するということは、$\chi$ 自体が確率的に揺らぐということです。瞬時SNR $\chi$ はRCS $\sigma$ に比例するので、平均RCS $\bar\sigma$ に対応する平均SNRを $\bar\chi$ とすれば、$\chi$ も平均 $\bar\chi$ の同じ形の分布に従います。

したがって変動ターゲットの検出確率は、各瞬時SNRでの検出確率を、SNRの確率分布で重み付け平均(期待値)したものになります。

$$ \bar{P_d} = \int_0^{\infty} P_d(\chi)\, p(\chi)\, d\chi = \int_0^{\infty} Q_M\left(\sqrt{2\chi}, \sqrt{2\eta_T}\right) p(\chi)\, d\chi $$

これが変動ターゲットの検出確率を求める一般的な枠組みです。あとは各モデルの $p(\chi)$ を代入して積分するだけです。

Swerling I・II(指数分布, $K=1$)

Swerling I と II は同じ指数分布のRCSを持ちます。単一パルス($n=1$)の場合、両者の検出確率は同一の式になります(パルス積分のときに初めて差が出ます)。瞬時SNRが平均 $\bar\chi$ の指数分布

$$ p(\chi) = \frac{1}{\bar\chi} e^{-\chi/\bar\chi} $$

に従うとして、上の積分を実行します。指数分布RCSの場合、MarcumのQ関数を含む積分が解析的に実行でき、驚くほど簡単な閉じた形が得られます。

$$ \boxed{\bar{P_d} = \exp\left(-\frac{\eta_T}{1 + \bar\chi}\right)} $$

この式を導くには、$Q_M$ の積分表現と指数分布を掛けて積分順序を入れ替え、ガウス積分を実行します。結果は非常にきれいで、$P_{fa} = e^{-\eta_T}$ を使うと次のようにも書けます。

$$ \bar{P_d} = (P_{fa})^{1/(1+\bar\chi)} $$

この形は実用上とても便利です。$P_{fa}$ と平均SNR $\bar\chi$ さえ与えれば、特殊関数なしで $P_d$ が一発で計算できます。$\bar\chi \to \infty$ で $\bar{P_d} \to 1$、$\bar\chi \to 0$ で $\bar{P_d} \to P_{fa}$ となり、極限の振る舞いも妥当です。

Swerling III・IV(自由度4カイ二乗分布, $K=2$)

Swerling III と IV は自由度4のカイ二乗分布、すなわち $K=2$ のガンマ分布をRCSが持ちます。瞬時SNRの分布は

$$ p(\chi) = \frac{4\chi}{\bar\chi^2} \exp\left(-\frac{2\chi}{\bar\chi}\right) $$

です。これを同様に積分すると、単一パルスでの検出確率は次の閉じた形で得られます。

$$ \boxed{\bar{P_d} = \left(1 + \frac{2}{\bar\chi} \cdot \frac{\eta_T}{(1 + \bar\chi/2)^2}\cdot \frac{1}{\cdots}\right)\cdots} $$

実際の導出は項が増えますが、結果として Swerling III/IV の単一パルス検出確率は

$$ \bar{P_d} = \left(1 + \frac{2 \eta_T}{2 + \bar\chi} – \cdots\right) \exp\left(-\frac{\eta_T}{1 + \bar\chi/2}\right) $$

の形になります。複雑なので、本記事では実装上は数値積分または非心カイ二乗分布の関係を用いて計算します。重要なのは、$K=2$ の分布は $K=1$ の指数分布よりもピークが鋭く、極端に小さなSNRになる確率が低いため、同じ平均SNRでも検出が容易になるという定性的事実です。

ここまでで単一パルスの $P_d$ が出揃いました。しかし実際のレーダーは1回の観測で複数パルスを送って積分します。次にその効果を見ましょう。

パルス積分の効果

レーダーがアンテナビームを目標に向けている間、通常は1発ではなく $n$ 発のパルスを送り、その受信信号を積み重ねます。これを 積分(integration) と呼びます。雑音はパルスごとに独立なので積分すると平均化されて相対的に小さくなり、信号成分は積み上がるため、実効的なSNRが向上します。

ノンコヒーレント積分(包絡線検波後に足し合わせる)では、$n$ パルスを積分したときの検出確率は、$n$ 個の独立な信号+雑音サンプルの和がしきい値を超える確率として計算されます。雑音だけの和は自由度 $2n$ のカイ二乗分布に従い、しきい値も $n$ パルス用に設定し直します。

ここで Swerlingモデルの違いが本質的に効いてきます。

  • スキャン間変動(Swerling I, III): $n$ パルス全体で1つのRCS値を共有します。つまり「不運にもRCSが小さいスキャン」では $n$ パルスすべてが弱く、積分してもダメなものはダメ。積分による平均化はSNRに対してのみ効き、RCS変動に対しては効きません。
  • パルス間変動(Swerling II, IV): $n$ パルスがそれぞれ独立なRCS値を持ちます。すると「全パルスが同時に弱くなる」確率が低く、いくつかは強いパルスが混ざるため、変動による不利が積分で平均化されます。これが ダイバーシティ利得 です。

この違いは検出確率曲線の形にはっきり現れます。少数パルスでは Swerling I(スキャン間)の方が II(パルス間)より検出が困難で、所要SNRが高くなります。逆に高い $P_d$ を目指すほど変動損失は大きくなり、パルス間変動モデルではダイバーシティでそれが緩和されます。

一般化すると、ガンマ分布パラメータ $K$ と積分パルス数 $n$ を使って、各モデルの「実効的な自由度」が決まります。スキャン間変動では揺らぎの自由度が $2K$ のまま、パルス間変動では $2Kn$ に増えると考えると見通しがよくなります。自由度が大きいほど分布が集中し、変動損失が小さくなる、という理解です。

パルス積分の定性的な効き方が分かったので、いよいよPythonで検出確率曲線を描いて、これらの議論を定量的に確かめましょう。

Pythonによる検出確率曲線の実装

単一パルスでの各モデル比較

まず単一パルス($n=1$)について、Swerling 0・I/II・III/IV の検出確率を平均SNRの関数として描きます。Swerling I/II は閉じた式 $\bar{P_d} = (P_{fa})^{1/(1+\bar\chi)}$ を使い、非変動(Swerling 0)は MarcumのQ関数(非心カイ二乗分布で計算)、Swerling III/IV は瞬時SNRのガンマ分布での数値積分で計算します。

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import ncx2, gamma
from scipy.integrate import quad

# 誤警報確率としきい値(単一パルス)
Pfa = 1e-6
eta_T = -np.log(Pfa)   # 規格化しきい値 eta_T = ln(1/Pfa)

# SNR 範囲 [dB]
snr_db = np.linspace(-5, 30, 200)
snr_lin = 10 ** (snr_db / 10)   # 平均SNR chi_bar (真値)

# --- Swerling 0 (非変動): Marcum Q関数 = 非心カイ二乗の上側確率 ---
def pd_swerling0(chi, eta):
    # 包絡線二乗 2r^2/(2psi) は自由度2・非心パラメータ 2*chi の非心カイ二乗
    # Pd = P(X > 2*eta), X ~ ncx2(df=2, nc=2*chi)
    return ncx2.sf(2 * eta, df=2, nc=2 * chi)

pd0 = np.array([pd_swerling0(c, eta_T) for c in snr_lin])

# --- Swerling I/II (指数分布, 単一パルスは同一): 閉じた式 ---
pd12 = np.exp(-eta_T / (1 + snr_lin))

# --- Swerling III/IV (自由度4カイ二乗, 単一パルスは同一): 数値積分 ---
def pd_swerling34(chi_bar, eta):
    # 瞬時SNR chi はガンマ分布 (形状K=2, スケール chi_bar/2)
    # Pd = ∫ Pd0(chi) p(chi) dchi
    def integrand(chi):
        p_chi = gamma.pdf(chi, a=2, scale=chi_bar / 2)
        return pd_swerling0(chi, eta) * p_chi
    val, _ = quad(integrand, 0, chi_bar * 30 + 50, limit=200)
    return val

pd34 = np.array([pd_swerling34(c, eta_T) for c in snr_lin])

ここまでで各モデルの検出確率を計算しました。Swerling 0 は非心カイ二乗分布の上側確率としてMarcumのQ関数を評価し、Swerling I/II は導いた閉じた式、Swerling III/IV は瞬時SNRをガンマ分布(形状 $K=2$)で平均化する数値積分で求めています。次にこれらを1枚のグラフにまとめて描画します。

# 検出確率曲線の描画
plt.figure(figsize=(9, 6))
plt.plot(snr_db, pd0,  'k-',  lw=2, label='Swerling 0 (non-fluctuating)')
plt.plot(snr_db, pd12, 'b--', lw=2, label='Swerling I/II (exponential RCS)')
plt.plot(snr_db, pd34, 'r-.', lw=2, label='Swerling III/IV (chi-sq 4 dof)')

plt.axhline(0.9, color='gray', ls=':', lw=1)
plt.text(-4, 0.91, 'Pd = 0.9', fontsize=10, color='gray')
plt.xlabel('Average SNR [dB]')
plt.ylabel('Detection Probability  $P_d$')
plt.title(f'Single-Pulse Detection Probability (Pfa = {Pfa:.0e})')
plt.legend(loc='lower right')
plt.grid(alpha=0.3)
plt.ylim(0, 1.02)
plt.tight_layout()
plt.show()

このグラフから、変動損失の本質が読み取れます。第一に、$P_d = 0.9$ という高い検出確率を達成するのに必要なSNRは、非変動(Swerling 0)が最も低く、Swerling III/IV がその次、Swerling I/II(指数分布)が最も高くなります。$P_{fa}=10^{-6}$ では、Swerling 0 が約13 dB なのに対し、Swerling I/II は約21 dB と、8 dB 以上の変動損失 が生じます。第二に、揺らぎの激しい指数分布($K=1$)の方が自由度4カイ二乗($K=2$)より損失が大きく、$K$ が大きいほど非変動に近づくという理論通りの順序になっています。第三に、低い $P_d$ 領域(左側)では逆転が起き、変動ターゲットの方が検出されやすい場合がある点も見て取れます。これは、揺らぎによってたまたま大きなRCSになる「幸運な」確率があるためです。

パルス積分による検出性能の改善

次に、ノンコヒーレント積分のパルス数 $n$ を変えたとき、検出性能がどう改善するかを見ます。ここでは比較を明確にするため、$n$ パルス積分後の実効しきい値を自由度 $2n$ のカイ二乗分布で設定し、スキャン間変動(全パルス同一RCS)とパルス間変動(パルス独立RCS)の違いをモンテカルロ・シミュレーションで確かめます。

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import chi2

rng = np.random.default_rng(0)

def simulate_pd(snr_db, n_pulses, model, Pfa=1e-6, trials=200000):
    """ノンコヒーレント積分の検出確率をモンテカルロで推定"""
    snr = 10 ** (snr_db / 10)            # 平均SNR
    # n パルス積分後のしきい値: 雑音和は自由度2nのカイ二乗(各I/Qで)
    VT = chi2.isf(Pfa, df=2 * n_pulses)  # 上側Pfa点
    detections = 0
    for _ in range(trials):
        if model == 'swerling1':   # スキャン間: 全パルス同一RCS(指数)
            a = rng.exponential(snr)          # 1つのRCS倍率
            amps = np.full(n_pulses, a)
        elif model == 'swerling2': # パルス間: 各パルス独立RCS(指数)
            amps = rng.exponential(snr, n_pulses)
        else:                       # 非変動
            amps = np.full(n_pulses, snr)
        # 各パルスの信号+雑音(包絡線二乗) ~ 非心カイ二乗近似を簡易生成
        # I,Q に信号平均 sqrt(amp), 雑音分散1/2 を与える
        sig = np.sqrt(amps)
        I = rng.normal(sig, np.sqrt(0.5), n_pulses)
        Q = rng.normal(0,   np.sqrt(0.5), n_pulses)
        total = np.sum(2 * (I**2 + Q**2))     # 自由度2nに規格化
        if total > VT:
            detections += 1
    return detections / trials

このシミュレーション関数は、各試行でモデルに応じたRCS倍率を生成し(スキャン間は全パルス同一値、パルス間はパルスごとに独立)、信号+雑音の包絡線二乗を $n$ パルス分合計してしきい値と比較します。雑音だけのときに合計が自由度 $2n$ のカイ二乗分布に従うよう規格化しているため、しきい値は chi2.isf(Pfa, df=2n) で設定できます。次に複数のSNRと積分数で実行して曲線を描きます。

# パルス数を変えて Swerling I (スキャン間) と II (パルス間) を比較
snr_range = np.arange(0, 26, 2)
n_list = [1, 10]
plt.figure(figsize=(9, 6))
styles = {'swerling1': ('o-', 'b'), 'swerling2': ('s--', 'r')}
labels = {'swerling1': 'Swerling I (scan-to-scan)',
          'swerling2': 'Swerling II (pulse-to-pulse)'}

for n in n_list:
    for model in ['swerling1', 'swerling2']:
        pd = [simulate_pd(s, n, model, trials=40000) for s in snr_range]
        st, col = styles[model]
        plt.plot(snr_range, pd, st, color=col, alpha=0.5 + 0.5*(n == 10),
                 label=f'{labels[model]}, n={n}')

plt.xlabel('Average SNR per pulse [dB]')
plt.ylabel('Detection Probability  $P_d$')
plt.title('Pulse Integration: Scan-to-scan vs Pulse-to-pulse Fluctuation')
plt.legend(loc='lower right', fontsize=9)
plt.grid(alpha=0.3)
plt.ylim(0, 1.02)
plt.tight_layout()
plt.show()

このグラフから、パルス積分とSwerlingモデルの相互作用が明確に読み取れます。第一に、$n=10$ パルス積分の曲線(濃い線)は $n=1$(薄い線)より左にシフトしており、積分が検出性能を向上させることが確認できます。第二に、$n=1$ では Swerling I と II の曲線がほぼ重なります(単一パルスでは両者は同一でした)。第三に、$n=10$ では両者が明確に分離し、Swerling II(パルス間変動)の方が高い $P_d$ 領域で有利 になります。これは独立なRCS値を10個平均することでダイバーシティ利得が得られ、「全パルスが同時に弱くなる不運」が起こりにくくなるためです。スキャン間変動(Swerling I)では全パルスが同じRCSを共有するため、この恩恵を受けられません。

変動損失の定量化

最後に、目標とする検出確率 $P_d = 0.9$ を達成するために各モデルで必要な所要SNRを求め、非変動ケースに対する「変動損失」を dB で定量化します。

import numpy as np
from scipy.stats import ncx2, gamma
from scipy.integrate import quad
from scipy.optimize import brentq

Pfa = 1e-6
eta_T = -np.log(Pfa)
target_pd = 0.9

def pd0(chi):
    return ncx2.sf(2 * eta_T, df=2, nc=2 * chi)

def pd_exp(chi_bar):          # Swerling I/II 単一パルス
    return np.exp(-eta_T / (1 + chi_bar))

def pd_chi4(chi_bar):         # Swerling III/IV 単一パルス
    f = lambda c: pd0(c) * gamma.pdf(c, a=2, scale=chi_bar / 2)
    return quad(f, 0, chi_bar * 30 + 50, limit=200)[0]

# 所要SNR(真値)を二分法で求める
def required_snr(pd_func):
    g = lambda snr_db: pd_func(10 ** (snr_db / 10)) - target_pd
    return brentq(g, -10, 40)

snr0   = required_snr(pd0)
snr12  = required_snr(pd_exp)
snr34  = required_snr(pd_chi4)

print(f"所要SNR (Pd=0.9, Pfa=1e-6, 単一パルス):")
print(f"  Swerling 0   : {snr0:6.2f} dB  (基準)")
print(f"  Swerling I/II: {snr12:6.2f} dB  (変動損失 {snr12-snr0:5.2f} dB)")
print(f"  Swerling III/IV: {snr34:6.2f} dB  (変動損失 {snr34-snr0:5.2f} dB)")

このコードの出力では、Swerling 0 の所要SNRが約13 dB なのに対し、Swerling I/II は約21 dB(変動損失 約8 dB)、Swerling III/IV は約16〜17 dB(変動損失 約3〜4 dB)となります。ここから読み取れる実務的な結論は明快です。RCSが激しく揺らぐ指数分布ターゲット(多数散乱体型、典型的には複雑形状の航空機)を高確率で検出するには、非変動ターゲットの想定より8 dB 近く、すなわち送信電力にして約6倍以上を上乗せして設計する必要があるということです。一方、支配的散乱体を持つターゲット($K=2$)では損失は半分以下に収まります。レーダーの探知距離はSNRの4乗根でしか伸びないため、この数 dB の差が探知距離の数十パーセントの違いに直結します。

これらの定量結果が、Swelingモデルを学ぶ最大の実用的価値です。最後に全体を整理しましょう。

まとめ

本記事では、レーダーターゲットのRCS変動を表すSwerlingモデルと、それが検出確率に与える影響を、確率分布の導出から定量評価まで一貫して解説しました。要点を振り返ります。

  • 検出は統計判定: レーダー検出は「信号+雑音」と「雑音のみ」を、誤警報確率 $P_{fa}$ を固定したしきい値で判別する問題であり、検出確率 $P_d$ はSNRとRCS変動統計の両方に依存する
  • 非変動の基準: 揺らがないターゲット(Swerling 0)の検出確率は $P_d = Q_M(\sqrt{2\chi}, \sqrt{2\eta_T})$ とMarcumのQ関数で表され、すべての評価の基準になる
  • 変動の物理的起源: 多数の同程度の散乱体からの反射のベクトル和は、中心極限定理によりレイリー振幅・指数分布RCS($K=1$)を生む。支配的散乱体があると自由度4カイ二乗分布RCS($K=2$)になる
  • 4モデルの分類: RCS分布(指数 vs 自由度4カイ二乗)と変動の速さ(スキャン間 vs パルス間)の2軸で Swerling I〜IV が定義される。Swerling I/II は指数分布、III/IV は自由度4カイ二乗分布に対応する
  • 閉じた式: Swerling I/II の単一パルス検出確率は $P_d = \exp(-\eta_T/(1+\bar\chi)) = (P_{fa})^{1/(1+\bar\chi)}$ という簡潔な形になる
  • 変動損失: $P_d=0.9$ を得る所要SNRは、非変動に比べ指数分布で約8 dB、自由度4カイ二乗で約3〜4 dB 増える。揺らぎが激しいほど($K$ が小さいほど)損失が大きい
  • パルス積分とダイバーシティ: パルス間変動(Swerling II/IV)では各パルスのRCSが独立なため、積分によるダイバーシティ利得で変動損失が緩和される。スキャン間変動(I/III)ではこの恩恵が小さい

Swerlingモデルは、レーダー方程式から得られる「平均的な探知距離」を、現実的な変動を考慮した「ある検出確率を保証する探知距離」へと格上げするための要です。次のステップとして、しきい値を雑音環境に適応させる実際の信号処理である CFAR検出 や、探知距離設計の基礎となる レーダー方程式、変動の物理的起源を担う レーダー断面積(RCS) を合わせて読むと、レーダー検出性能の全体像がつながって見えてきます。

次のステップとして、以下の記事も参考にしてください。