レーダーにおける目標検出は、受信信号から目標の有無を判定する統計的な意思決定問題です。最も単純なアプローチは固定閾値と比較する方法ですが、クラッタレベルが時間的・空間的に変動する実環境では、誤警報率（False Alarm Rate）が大きく変動してしまいます。この問題を解決するのがCFAR（Constant False Alarm Rate: 定誤警報率）検出器です。
CFARは、検出セル周辺の参照セルからクラッタ+雑音レベルを推定し、閾値を適応的に設定することで、環境が変化しても誤警報率を一定に保ちます。本記事では、ネイマン=ピアソン検定に基づく統計的定式化から、CA-CFAR（Cell Averaging CFAR）の閾値係数 $\alpha$ の導出、クラッタエッジ対策としてのGO-CFAR/SO-CFAR、マルチターゲット環境に強いOS-CFAR（Order Statistic CFAR）まで解説し、Pythonで実装・可視化します。
本記事の内容
- レーダー目標検出の統計的定式化（ネイマン=ピアソン検定）
- 固定閾値の問題とCFARの必要性
- CA-CFARのアルゴリズムと閾値係数 $\alpha$ の導出
- GO-CFAR/SO-CFARのクラッタエッジ対策
- OS-CFAR（順序統計量ベース）のアルゴリズム
- マルチターゲット環境での問題
- Pythonでの各種CFAR実装とシミュレーション
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
レーダー目標検出の統計的定式化
仮説検定としての定式化
レーダーの目標検出は、各レンジセルにおいて以下の2つの仮説のどちらが正しいかを判定する二項仮説検定です。
$$ \begin{cases} H_0: & r = n \quad (\text{目標なし: 雑音のみ}) \\ H_1: & r = s + n \quad (\text{目標あり: 信号 + 雑音}) \end{cases} $$
ここで $r$ は受信信号、$s$ は目標からの反射信号、$n$ は雑音（クラッタを含む）です。
検出性能の指標
仮説検定の結果は4通りに分類されます。
| $H_0$ 真 (目標なし) | $H_1$ 真 (目標あり) | |
|---|---|---|
| $H_0$ 採択 | 正棄却 (TN) | 見落とし (Miss) |
| $H_1$ 採択 | 誤警報 (FA) | 検出 (Detection) |
誤警報確率 $P_{fa}$: 目標がないのに「あり」と判定する確率
$$ P_{fa} = P(\text{判定} = H_1 \mid H_0) $$
検出確率 $P_d$: 目標があるとき正しく「あり」と判定する確率
$$ P_d = P(\text{判定} = H_1 \mid H_1) $$
ネイマン=ピアソンの基準
ネイマン=ピアソン検定は、誤警報確率 $P_{fa}$ を指定値以下に制約しつつ、検出確率 $P_d$ を最大化する最適な判定方式を与えます。
尤度比検定（Likelihood Ratio Test）として定式化すると、次の判定規則が最適です。
$$ \Lambda(r) = \frac{p(r \mid H_1)}{p(r \mid H_0)} \underset{H_0}{\overset{H_1}{\gtrless}} \gamma $$
ここで $p(r \mid H_i)$ は仮説 $H_i$ のもとでの $r$ の確率密度関数、$\gamma$ は閾値です。
包絡線検波後の検出
レーダー受信機の包絡線検波後の出力電力 $z = |r|^2$ を考えます。
$H_0$（雑音のみ）: 複素ガウス雑音の振幅の2乗 $z$ は指数分布に従います。
$$ p(z \mid H_0) = \frac{1}{\mu} \exp\!\left(-\frac{z}{\mu}\right), \quad z \geq 0 $$
ここで $\mu = E[z \mid H_0]$ は雑音の平均電力です。
固定閾値 $T$ に対する誤警報確率は次の通りです。
$$ \begin{align} P_{fa} &= P(z > T \mid H_0) \\ &= \int_T^{\infty} \frac{1}{\mu} \exp\!\left(-\frac{z}{\mu}\right) dz \\ &= \exp\!\left(-\frac{T}{\mu}\right) \end{align} $$
固定閾値の問題
誤警報率の変動
上式を変形すると、所望の $P_{fa}$ を達成するための閾値は次のようになります。
$$ T = -\mu \ln P_{fa} = \mu \cdot \alpha_{\mathrm{fixed}} $$
ここで $\alpha_{\mathrm{fixed}} = -\ln P_{fa}$ は閾値係数です。
問題は、この閾値が雑音電力 $\mu$ に依存することです。$\mu$ が既知で一定であれば、固定閾値で所望の $P_{fa}$ を維持できます。しかし、現実のレーダーではクラッタレベルが空間的・時間的に変動するため、$\mu$ は一定ではありません。
$\mu$ が変化したときの誤警報確率を計算してみましょう。閾値を $T = \mu_0 \cdot \alpha_{\mathrm{fixed}}$（設計時の $\mu_0$ に基づく）に固定し、実際の雑音電力が $\mu_1 = k\mu_0$ に変化した場合の $P_{fa}$ は次の通りです。
$$ P_{fa} = \exp\!\left(-\frac{T}{\mu_1}\right) = \exp\!\left(-\frac{\mu_0 \alpha_{\mathrm{fixed}}}{k\mu_0}\right) = \exp\!\left(-\frac{\alpha_{\mathrm{fixed}}}{k}\right) $$
例えば、$P_{fa} = 10^{-6}$ を設計値として $\alpha_{\mathrm{fixed}} = 6\ln 10 \approx 13.8$ とし、クラッタレベルが3 dB上昇（$k = 2$）した場合を計算します。
$$ P_{fa} = \exp\!\left(-\frac{13.8}{2}\right) = \exp(-6.9) \approx 10^{-3} $$
誤警報確率が $10^{-6}$ から $10^{-3}$ へと1000倍に悪化しました。このように、固定閾値はクラッタ変動に対して非常に脆弱です。
CA-CFAR（Cell Averaging CFAR）
アルゴリズム
CA-CFARは、検出セル（CUT: Cell Under Test）の周囲にある参照セル（Reference Cells）の電力値を平均して雑音レベルを推定し、適応的に閾値を設定します。
構成は以下の通りです。
[ガード] [参照セル N/2 個] [ガード] [CUT] [ガード] [参照セル N/2 個] [ガード]
ガードセル: CUT近傍のセルには目標信号が漏れ込む可能性があるため、参照セルから除外します。
参照セル: CUTの両側に $N/2$ 個ずつ、計 $N$ 個の参照セルを配置します。
雑音レベル推定
$N$ 個の参照セルの電力値を $z_1, z_2, \ldots, z_N$ とすると、雑音レベルの推定値は次の通りです。
$$ \hat{\mu} = \frac{1}{N} \sum_{i=1}^{N} z_i $$
閾値は次のように設定されます。
$$ T = \alpha \cdot \hat{\mu} = \alpha \cdot \frac{1}{N} \sum_{i=1}^{N} z_i $$
CUTの電力 $z_{\mathrm{CUT}}$ が閾値を超えれば目標と判定します。
$$ z_{\mathrm{CUT}} \underset{H_0}{\overset{H_1}{\gtrless}} T = \alpha \cdot \hat{\mu} $$
閾値係数 $\alpha$ の導出
CA-CFARの閾値係数 $\alpha$ を、所望の誤警報確率 $P_{fa}$ から導出しましょう。
$H_0$ のもとで各セルの電力 $z_i$ は独立同一な指数分布 $\mathrm{Exp}(\mu)$ に従います。$N$ 個の参照セルの合計を $Z = \sum_{i=1}^{N} z_i$ とすると、$Z$ はガンマ分布 $\Gamma(N, \mu)$ に従います。
$$ p_Z(z) = \frac{z^{N-1}}{\mu^N \Gamma(N)} \exp\!\left(-\frac{z}{\mu}\right) $$
誤警報確率は、CUTの電力 $z_{\mathrm{CUT}}$（$\mathrm{Exp}(\mu)$）が閾値 $\alpha Z / N$ を超える確率です。
$$ P_{fa} = P\!\left(z_{\mathrm{CUT}} > \frac{\alpha}{N} Z \bigg| H_0\right) $$
$z_{\mathrm{CUT}}$ と $Z$ は独立なので、条件付き確率を $Z$ に関して積分します。$Z = z$ を与えたときの条件付き誤警報確率は次の通りです。
$$ P(z_{\mathrm{CUT}} > \alpha z / N \mid Z = z) = \exp\!\left(-\frac{\alpha z}{N\mu}\right) $$
$Z$ に関して期待値をとります。
$$ \begin{align} P_{fa} &= \int_0^{\infty} \exp\!\left(-\frac{\alpha z}{N\mu}\right) \cdot \frac{z^{N-1}}{\mu^N \Gamma(N)} \exp\!\left(-\frac{z}{\mu}\right) dz \\ &= \frac{1}{\mu^N \Gamma(N)} \int_0^{\infty} z^{N-1} \exp\!\left(-\frac{z}{\mu}\left(1 + \frac{\alpha}{N}\right)\right) dz \end{align} $$
変数置換 $u = z(1 + \alpha/N)/\mu$ により次のようになります。
$$ \begin{align} P_{fa} &= \frac{1}{\mu^N \Gamma(N)} \cdot \frac{\mu^N}{(1 + \alpha/N)^N} \int_0^{\infty} u^{N-1} e^{-u} du \\ &= \frac{\Gamma(N)}{(1 + \alpha/N)^N \cdot \Gamma(N)} \\ &= \frac{1}{(1 + \alpha/N)^N} \end{align} $$
$$ \boxed{P_{fa} = \left(1 + \frac{\alpha}{N}\right)^{-N}} $$
この式から閾値係数 $\alpha$ を求めると次のようになります。
$$ \boxed{\alpha = N\left(P_{fa}^{-1/N} – 1\right)} $$
数値例
$N = 24$（参照セル数）、$P_{fa} = 10^{-6}$ の場合の閾値係数を計算します。
$$ \alpha = 24 \times \left((10^{-6})^{-1/24} – 1\right) = 24 \times \left(10^{0.25} – 1\right) = 24 \times (1.778 – 1) = 24 \times 0.778 = 18.67 $$
固定閾値の場合は $\alpha_{\mathrm{fixed}} = -\ln(10^{-6}) = 13.82$ でしたので、CA-CFARでは推定誤差を補うために閾値が高くなっています。この差はCFAR損失と呼ばれ、参照セル数 $N$ が多いほど小さくなります。
CFAR損失
CFAR損失（CFAR Loss）は、固定閾値に対するCA-CFAR閾値の増加分であり、目標検出に必要なSNRの増加量に対応します。
$$ L_{\mathrm{CFAR}} \, [\mathrm{dB}] = 10\log_{10}\!\left(\frac{\alpha}{\alpha_{\mathrm{fixed}}}\right) = 10\log_{10}\!\left(\frac{N(P_{fa}^{-1/N}-1)}{-\ln P_{fa}}\right) $$
$N$ が十分大きいとき、$P_{fa}^{-1/N} \approx 1 + (-\ln P_{fa})/N$ となり、$\alpha \approx -\ln P_{fa} = \alpha_{\mathrm{fixed}}$ に近づきます。つまり、参照セル数を増やすほどCFAR損失は減少します。
GO-CFARとSO-CFAR
クラッタエッジの問題
CA-CFARは均一なクラッタ環境では良好に機能しますが、CUTの片側にクラッタ境界（クラッタエッジ）がある場合に問題が生じます。
ケース1: CUTがクラッタ内、片側が低雑音領域
低雑音側の参照セルが全体の平均を引き下げるため、閾値が低くなりすぎて誤警報が増大します。
ケース2: CUTが低雑音領域、片側がクラッタ
クラッタ側の参照セルが全体の平均を引き上げるため、閾値が高くなりすぎて検出確率が低下します。
GO-CFAR（Greatest Of CFAR）
GO-CFARは、CUTの両側の参照セルそれぞれの平均を計算し、大きい方を雑音推定値として使用します。
$$ \hat{\mu}_{\mathrm{lead}} = \frac{1}{N/2} \sum_{i=1}^{N/2} z_i, \quad \hat{\mu}_{\mathrm{lag}} = \frac{1}{N/2} \sum_{i=N/2+1}^{N} z_i $$
$$ \hat{\mu}_{\mathrm{GO}} = \max(\hat{\mu}_{\mathrm{lead}}, \hat{\mu}_{\mathrm{lag}}) $$
閾値は $T_{\mathrm{GO}} = \alpha_{\mathrm{GO}} \cdot \hat{\mu}_{\mathrm{GO}}$ です。
長所: クラッタエッジで大きい方を選ぶため、低雑音側に引きずられることなく誤警報を抑制できます。
短所: 均一環境では両側の大きい方を選ぶことで閾値が高くなりがちで、CFAR損失が増大します。
SO-CFAR（Smallest Of CFAR）
SO-CFARは、GO-CFARとは逆に小さい方を使用します。
$$ \hat{\mu}_{\mathrm{SO}} = \min(\hat{\mu}_{\mathrm{lead}}, \hat{\mu}_{\mathrm{lag}}) $$
長所: 均一環境での損失がCA-CFARより小さい場合があります。
短所: クラッタエッジで低雑音側を選んでしまうため、誤警報が増大するリスクがあります。
OS-CFAR（Order Statistic CFAR）
アルゴリズム
OS-CFARは、$N$ 個の参照セルの電力値を昇順にソートし、$k$ 番目の順序統計量を雑音推定値として使用します。
$$ z_{(1)} \leq z_{(2)} \leq \cdots \leq z_{(N)} $$
$$ \hat{\mu}_{\mathrm{OS}} = z_{(k)} $$
閾値は $T_{\mathrm{OS}} = \alpha_{\mathrm{OS}} \cdot z_{(k)}$ です。典型的には $k \approx 3N/4$ が使用されます。
マルチターゲット環境での利点
OS-CFARの最大の利点は、マルチターゲット環境での耐性です。CA-CFARでは参照セル内に目標が含まれると平均値が上昇し、閾値が不当に高くなってCUT内の目標を見落とす問題（マスキング効果）があります。
OS-CFARでは $k$ 番目の順序統計量を使用するため、参照セル内に $N – k$ 個以下の目標（外れ値）があっても雑音推定値に影響しません。
閾値係数の計算
$H_0$ のもとで各セルが独立な $\mathrm{Exp}(\mu)$ に従う場合、$k$ 番目の順序統計量 $z_{(k)}$ の期待値は次の通りです。
$$ E[z_{(k)}] = \mu \sum_{i=1}^{k} \frac{1}{N – i + 1} $$
閾値係数 $\alpha_{\mathrm{OS}}$ は、所望の $P_{fa}$ を達成するように数値的に決定されます。$k$ 番目の順序統計量の分布からの正確な導出は複雑ですが、モンテカルロシミュレーションで容易に求められます。
Pythonでの実装
テスト用レーダーデータの生成
import numpy as np
import matplotlib.pyplot as plt
def generate_radar_data(N_cells=500, seed=42):
"""クラッタ + 目標 + 雑音のテストデータを生成"""
np.random.seed(seed)
# 雑音（指数分布: |複素ガウス|^2）
noise_power = 1.0
data = np.random.exponential(noise_power, N_cells)
# クラッタエッジ: セル200-350で雑音レベルが10倍
clutter_region = slice(200, 350)
data[clutter_region] = np.random.exponential(noise_power * 10, 150)
# 目標の追加 (セル番号, SNR [線形])
target_info = [
(80, 15), # 孤立目標
(160, 20), # 孤立目標（強い）
(205, 12), # クラッタエッジ付近の目標
(280, 25), # クラッタ内の目標
(400, 10), # 孤立目標（弱い）
(402, 12), # 近接目標
]
for (cell, snr) in target_info:
data[cell] += noise_power * snr
return data, target_info
# データ生成と表示
data, target_info = generate_radar_data()
fig, ax = plt.subplots(figsize=(14, 5))
ax.plot(10 * np.log10(data + 1e-10), 'b-', linewidth=0.8, label='Received Power')
for (cell, snr) in target_info:
ax.axvline(x=cell, color='red', linestyle='--', alpha=0.5, linewidth=1)
ax.set_xlabel('Range Cell', fontsize=12)
ax.set_ylabel('Power [dB]', fontsize=12)
ax.set_title('Radar Data: Targets + Clutter Edge + Noise', fontsize=14)
ax.legend(fontsize=11)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
CA-CFARの実装
import numpy as np
import matplotlib.pyplot as plt
def ca_cfar(data, N_ref, N_guard, Pfa):
"""
CA-CFAR（Cell Averaging CFAR）の実装
Parameters:
data: 入力電力データ（線形スケール）
N_ref: 片側の参照セル数
N_guard: 片側のガードセル数
Pfa: 所望の誤警報確率
Returns:
detections: 検出フラグ配列
threshold: 閾値配列
"""
N = len(data)
N_total_ref = 2 * N_ref # 両側合計の参照セル数
# 閾値係数の計算
alpha = N_total_ref * (Pfa**(-1.0/N_total_ref) - 1)
detections = np.zeros(N, dtype=bool)
threshold = np.zeros(N)
for i in range(N_guard + N_ref, N - N_guard - N_ref):
# 参照セルのインデックス
lead_start = i - N_guard - N_ref
lead_end = i - N_guard
lag_start = i + N_guard + 1
lag_end = i + N_guard + N_ref + 1
# 参照セルの電力平均
ref_cells = np.concatenate([data[lead_start:lead_end],
data[lag_start:lag_end]])
noise_estimate = np.mean(ref_cells)
# 閾値設定
threshold[i] = alpha * noise_estimate
# 検出判定
if data[i] > threshold[i]:
detections[i] = True
return detections, threshold
# データ生成
def generate_radar_data(N_cells=500, seed=42):
np.random.seed(seed)
noise_power = 1.0
data = np.random.exponential(noise_power, N_cells)
data[200:350] = np.random.exponential(noise_power * 10, 150)
target_info = [
(80, 15), (160, 20), (205, 12), (280, 25), (400, 10), (402, 12),
]
for (cell, snr) in target_info:
data[cell] += noise_power * snr
return data, target_info
data, target_info = generate_radar_data()
# CA-CFAR適用
N_ref = 12 # 片側の参照セル数
N_guard = 3 # 片側のガードセル数
Pfa = 1e-4 # 所望の誤警報確率
det, thresh = ca_cfar(data, N_ref, N_guard, Pfa)
# 可視化
fig, ax = plt.subplots(figsize=(14, 6))
ax.plot(10*np.log10(data + 1e-10), 'b-', linewidth=0.8, label='Received Power')
ax.plot(10*np.log10(thresh + 1e-10), 'r-', linewidth=1.5, label='CA-CFAR Threshold')
det_indices = np.where(det)[0]
ax.plot(det_indices, 10*np.log10(data[det_indices] + 1e-10),
'go', markersize=8, label='Detections')
ax.set_xlabel('Range Cell', fontsize=12)
ax.set_ylabel('Power [dB]', fontsize=12)
ax.set_title(f'CA-CFAR Detection (N_ref={2*N_ref}, N_guard={2*N_guard}, Pfa={Pfa})', fontsize=14)
ax.legend(fontsize=11)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
print(f"閾値係数 α = {2*N_ref*(Pfa**(-1.0/(2*N_ref))-1):.2f}")
print(f"検出されたセル: {det_indices}")
GO-CFAR / SO-CFAR の実装
import numpy as np
import matplotlib.pyplot as plt
def go_cfar(data, N_ref, N_guard, Pfa):
"""GO-CFAR: 両側の平均の大きい方を使用"""
N = len(data)
N_half = N_ref # 片側の参照セル数
alpha = N_half * (Pfa**(-1.0/N_half) - 1) # 片側のNで計算
detections = np.zeros(N, dtype=bool)
threshold = np.zeros(N)
for i in range(N_guard + N_ref, N - N_guard - N_ref):
lead = data[i - N_guard - N_ref : i - N_guard]
lag = data[i + N_guard + 1 : i + N_guard + N_ref + 1]
mu_lead = np.mean(lead)
mu_lag = np.mean(lag)
noise_estimate = max(mu_lead, mu_lag)
threshold[i] = alpha * noise_estimate
if data[i] > threshold[i]:
detections[i] = True
return detections, threshold
def so_cfar(data, N_ref, N_guard, Pfa):
"""SO-CFAR: 両側の平均の小さい方を使用"""
N = len(data)
N_half = N_ref
alpha = N_half * (Pfa**(-1.0/N_half) - 1)
detections = np.zeros(N, dtype=bool)
threshold = np.zeros(N)
for i in range(N_guard + N_ref, N - N_guard - N_ref):
lead = data[i - N_guard - N_ref : i - N_guard]
lag = data[i + N_guard + 1 : i + N_guard + N_ref + 1]
mu_lead = np.mean(lead)
mu_lag = np.mean(lag)
noise_estimate = min(mu_lead, mu_lag)
threshold[i] = alpha * noise_estimate
if data[i] > threshold[i]:
detections[i] = True
return detections, threshold
# データ生成
def generate_radar_data(N_cells=500, seed=42):
np.random.seed(seed)
noise_power = 1.0
data = np.random.exponential(noise_power, N_cells)
data[200:350] = np.random.exponential(noise_power * 10, 150)
target_info = [
(80, 15), (160, 20), (205, 12), (280, 25), (400, 10), (402, 12),
]
for (cell, snr) in target_info:
data[cell] += noise_power * snr
return data, target_info
data, target_info = generate_radar_data()
N_ref = 12
N_guard = 3
Pfa = 1e-4
det_go, thresh_go = go_cfar(data, N_ref, N_guard, Pfa)
det_so, thresh_so = so_cfar(data, N_ref, N_guard, Pfa)
# 可視化
fig, axes = plt.subplots(2, 1, figsize=(14, 10))
for ax, det, thresh, name, color in [
(axes[0], det_go, thresh_go, 'GO-CFAR', 'orange'),
(axes[1], det_so, thresh_so, 'SO-CFAR', 'purple'),
]:
ax.plot(10*np.log10(data + 1e-10), 'b-', linewidth=0.8, label='Received Power')
ax.plot(10*np.log10(thresh + 1e-10), color=color, linewidth=1.5,
label=f'{name} Threshold')
det_idx = np.where(det)[0]
ax.plot(det_idx, 10*np.log10(data[det_idx] + 1e-10), 'go', markersize=8,
label='Detections')
ax.set_xlabel('Range Cell', fontsize=11)
ax.set_ylabel('Power [dB]', fontsize=11)
ax.set_title(f'{name} Detection', fontsize=13)
ax.legend(fontsize=11)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
OS-CFARの実装
import numpy as np
import matplotlib.pyplot as plt
def os_cfar(data, N_ref, N_guard, k, Pfa):
"""
OS-CFAR（Order Statistic CFAR）の実装
Parameters:
data: 入力電力データ
N_ref: 片側の参照セル数
N_guard: 片側のガードセル数
k: 使用する順序統計量のランク（1-indexed）
Pfa: 所望の誤警報確率
"""
N = len(data)
N_total = 2 * N_ref
# 閾値係数（モンテカルロで推定するか、近似式を使用）
# ここでは近似的に計算
# P_fa = k * C(N,k) * B(alpha/(1+alpha), k, N-k+1) の数値逆解法
# 簡易的に固定値を使用（実用上は事前計算テーブル）
from scipy import special
# 正確な閾値係数をモンテカルロで推定
N_mc = 100000
ref_samples = np.random.exponential(1.0, (N_mc, N_total))
sorted_samples = np.sort(ref_samples, axis=1)
os_values = sorted_samples[:, k-1]
test_samples = np.random.exponential(1.0, N_mc)
# αを二分探索で求める
alpha_low, alpha_high = 0.1, 100.0
for _ in range(50):
alpha_mid = (alpha_low + alpha_high) / 2
pfa_est = np.mean(test_samples > alpha_mid * os_values)
if pfa_est > Pfa:
alpha_low = alpha_mid
else:
alpha_high = alpha_mid
alpha = (alpha_low + alpha_high) / 2
detections = np.zeros(N, dtype=bool)
threshold = np.zeros(N)
for i in range(N_guard + N_ref, N - N_guard - N_ref):
lead = data[i - N_guard - N_ref : i - N_guard]
lag = data[i + N_guard + 1 : i + N_guard + N_ref + 1]
ref_cells = np.concatenate([lead, lag])
sorted_ref = np.sort(ref_cells)
noise_estimate = sorted_ref[k - 1] # k番目の順序統計量
threshold[i] = alpha * noise_estimate
if data[i] > threshold[i]:
detections[i] = True
return detections, threshold, alpha
# データ生成
def generate_radar_data(N_cells=500, seed=42):
np.random.seed(seed)
noise_power = 1.0
data = np.random.exponential(noise_power, N_cells)
data[200:350] = np.random.exponential(noise_power * 10, 150)
target_info = [
(80, 15), (160, 20), (205, 12), (280, 25), (400, 10), (402, 12),
]
for (cell, snr) in target_info:
data[cell] += noise_power * snr
return data, target_info
data, target_info = generate_radar_data()
N_ref = 12
N_guard = 3
Pfa = 1e-4
k = 18 # 3/4 * 24 = 18
det_os, thresh_os, alpha_os = os_cfar(data, N_ref, N_guard, k, Pfa)
fig, ax = plt.subplots(figsize=(14, 6))
ax.plot(10*np.log10(data + 1e-10), 'b-', linewidth=0.8, label='Received Power')
ax.plot(10*np.log10(thresh_os + 1e-10), 'm-', linewidth=1.5,
label=f'OS-CFAR Threshold (k={k})')
det_idx = np.where(det_os)[0]
ax.plot(det_idx, 10*np.log10(data[det_idx] + 1e-10), 'go', markersize=8,
label='Detections')
ax.set_xlabel('Range Cell', fontsize=12)
ax.set_ylabel('Power [dB]', fontsize=12)
ax.set_title(f'OS-CFAR Detection (k={k}, α={alpha_os:.2f})', fontsize=14)
ax.legend(fontsize=11)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
print(f"OS-CFAR 閾値係数 α = {alpha_os:.2f}")
print(f"検出されたセル: {det_idx}")
全手法の比較
import numpy as np
import matplotlib.pyplot as plt
def generate_radar_data(N_cells=500, seed=42):
np.random.seed(seed)
noise_power = 1.0
data = np.random.exponential(noise_power, N_cells)
data[200:350] = np.random.exponential(noise_power * 10, 150)
target_info = [
(80, 15), (160, 20), (205, 12), (280, 25), (400, 10), (402, 12),
]
for (cell, snr) in target_info:
data[cell] += noise_power * snr
return data, target_info
data, target_info = generate_radar_data()
# 各CFARの閾値を計算して比較プロット
fig, ax = plt.subplots(figsize=(14, 7))
data_dB = 10 * np.log10(data + 1e-10)
ax.plot(data_dB, 'b-', linewidth=0.8, label='Received Power', zorder=1)
# 真の目標位置をマーク
for (cell, snr) in target_info:
ax.axvline(x=cell, color='gray', linestyle=':', alpha=0.3)
# 各CFARの結果をプロット
N_ref = 12
N_guard = 3
Pfa = 1e-4
# CA-CFAR
N_total = 2 * N_ref
alpha_ca = N_total * (Pfa**(-1.0/N_total) - 1)
thresh_ca = np.zeros(len(data))
for i in range(N_guard+N_ref, len(data)-N_guard-N_ref):
lead = data[i-N_guard-N_ref:i-N_guard]
lag = data[i+N_guard+1:i+N_guard+N_ref+1]
thresh_ca[i] = alpha_ca * np.mean(np.concatenate([lead, lag]))
# GO-CFAR
alpha_go = N_ref * (Pfa**(-1.0/N_ref) - 1)
thresh_go = np.zeros(len(data))
for i in range(N_guard+N_ref, len(data)-N_guard-N_ref):
lead = data[i-N_guard-N_ref:i-N_guard]
lag = data[i+N_guard+1:i+N_guard+N_ref+1]
thresh_go[i] = alpha_go * max(np.mean(lead), np.mean(lag))
ax.plot(10*np.log10(thresh_ca + 1e-10), 'r-', linewidth=1.5, label='CA-CFAR', alpha=0.8)
ax.plot(10*np.log10(thresh_go + 1e-10), 'orange', linewidth=1.5, label='GO-CFAR', alpha=0.8)
ax.set_xlabel('Range Cell', fontsize=12)
ax.set_ylabel('Power [dB]', fontsize=12)
ax.set_title('Comparison of CFAR Thresholds', fontsize=14)
ax.legend(fontsize=11)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
上記のコードを実行すると、クラッタエッジと複数目標を含むレーダーデータに対して、CA-CFAR、GO-CFAR、SO-CFAR、OS-CFARそれぞれの検出結果が可視化されます。クラッタエッジ付近では閾値の振る舞いが大きく異なり、各手法の特性の違いが明確に確認できます。
まとめ
本記事では、CFAR（定誤警報率）検出の理論を数式で丁寧に導出しました。
- レーダーの目標検出はネイマン=ピアソン検定として定式化され、誤警報確率 $P_{fa}$ を制約して検出確率 $P_d$ を最大化する
- 固定閾値ではクラッタレベルの変動に対して誤警報率が大きく変動するため、適応的閾値が必要
- CA-CFARは参照セルの平均から雑音推定を行い、閾値係数は $\alpha = N(P_{fa}^{-1/N} – 1)$ で与えられ、誤警報確率は $P_{fa} = (1 + \alpha/N)^{-N}$ となる
- GO-CFARは両側の大きい方を選ぶことでクラッタエッジでの誤警報を抑制する
- SO-CFARは小さい方を選ぶことで均一環境での損失を抑えるが、クラッタエッジに弱い
- OS-CFARは順序統計量ベースであり、マルチターゲット環境でマスキング効果に耐性がある
- 環境に応じた適切なCFAR手法の選択が実用上重要である
次のステップとして、以下の記事も参考にしてください。