问题的引出：雷达目标检测
雷达在接收到回波信号后，需要区分目标与噪声。目标检测方法的核心是阈值法。如果雷达回波大于阈值，则显示检测到目标，否则视为噪声。
假设将当前单元的功率为 \(Y\)，噪声功率为 \(\mu\)，使用的阈值因子为 \(\alpha\)，则：
\[\begin{cases} \text{target}& Y\ge \alpha\mu\\ \text{no target}& Y< \alpha\mu \end{cases} \]眼下要解决的问题是：如何确定噪声功率 \(\mu\)。
采用恒定阈值
采用固定的门限进行目标检测时可能会产生一定的虚警。
- 虚警，没有目标时判断为有目标
- 漏警，有目标时判断为没有目标
因此，最好根据雷达杂波数据动态调整检测门限。在虚警概率保持不变的情况下实现目标检测概率最大化，这种方法称为恒虚警率（Constant False Alarm Rate，CFAR）检测技术。
采用恒虚警 CFAR
有多种 CFAR 检测算法。基本思路都是动态确定背景噪声功率。
训练单元（Cells for Mean Level Computation），用于推导噪声功率的单元；
守护单元（"Guard" Cells），放在测试单元的两侧，用于消除训练单元邻接检测单元开始和结尾的多余能量；
检测单元（Cell Under Test），一个被指定的、作为检测中心的单元。
训练单位为 5、守护单位为 1 的检测窗口就会大概长这样：
训 训 训 训 训 守 检 守 训 训 训 训 训
- CA(Cell Averaging)-CFAR，使用临近 \(2n\) 个参考单元功率的均值作为背景噪声功率水平
- GO(Greatest Of)-CFAR、SO(Smallest Of)-CFAR，前 \(n\) 个参考单元功率之和，与后 \(n\) 个参考单元功率之和，两者中的 较大 / 较小 者作为背景噪声功率水平
- OS(Order Statistic)-CFAR，对参考单元从小到大排序，取第 \(k\) 个样本作为噪声功率
其中，CA-CFAR、OS-CFAR 算法最为普遍。
虚警概率 PFA
虚警概率（PFA），是噪声或杂波峰值被 CFAR 算法误认为目标的概率。\(P_{fa}\) 值通常很小，在 \(10^{-5}\) 或 \(10^{-6}\) 量级。
阈值因子 \(\alpha\) 是 \(P_{fa}\) 的函数。
CA-CFAR
CA-CAFR 使用待检测数据周围一定范围的单元格的平均值来确定该单元格的噪声功率。
\[\mu = \frac{1}{N_{\text{train}}}\sum^{N_{\text{train}}}_{i}X_i \]门限因子 \(\alpha\) 有固定的计算公式：
\[\alpha=N_{\text{train}}(P_{fa}^{-1/N_{\text{train}}}-1) \]CA-CFAR 的主要缺点是，临近的多目标会产生严重的掩蔽效应，导致漏检。
OS-CFAR
OS-CFAR 则是将临近范围的单元格功率大小进行升序排序，取第 \(k\) 个值作为噪声功率 \(\mu\)。
门限因子 \(\alpha\) 通过以下等式获得：
\[P_{fa}=k\binom{N}{k}\frac{(k-1)!(\alpha+N-k)!}{(\alpha+N)!} \]- \(N\) 是 train 单元数量
为了便于计算，将问题转化为：求使得函数 \(\frac{N!}{(N-k)!}\frac{(\alpha+N-k)!}{(\alpha+N)!}\frac{1}{P_{fa}}\rightarrow 1\) 的 \(\alpha\) 值。并使用割线法求得 \(\alpha\)。参考 radarsimpy，Python 代码如下：
def log_factorial(n):
"""用于近似计算阶乘的对数"""
n = n + 9.0
product = 1
for i in range(1, 9):
product *= n - i
return (
1 / 2 * (np.log(2 * np.pi) - np.log(n))
+ n * (np.log(n + 1 / (12 * n - (1 / 10 / n))) - 1)
- np.log(product)
)
def os_cfar_threshold(k, n, pfa):
"""计算 alpha"""
def fun(k, n, t_os, pfa):
return (
log_factorial(n)
- log_factorial(n - k)
- np.sum(np.log(np.arange(n, n - k, -1) + t_os))
- np.log(pfa)
)
max_iter = 10000
t_max = 1e32
t_min = 1
for _ in range(0, max_iter):
m_n = t_max - fun(k, n, t_max, pfa) * (t_min - t_max) / (
fun(k, n, t_min, pfa) - fun(k, n, t_max, pfa)
)
f_m_n = fun(k, n, m_n, pfa)
if f_m_n == 0:
return m_n
if np.abs(f_m_n) < 0.0001:
return m_n
if fun(k, n, t_max, pfa) * f_m_n < 0:
t_min = m_n
elif fun(k, n, t_min, pfa) * f_m_n < 0:
t_max = m_n
else:
break
return None
参考来源
- 【雷达目标检测】恒定阈值法和恒虚警(CFAR)法及代码实现 https://blog.csdn.net/qq_43485394/article/details/123417527
- 雷达系统必备知识——恒虚警检测器CFAR https://www.zhihu.com/tardis/zm/art/269840008
- 【目标检测】雷达目标CFAR检测算法 https://blog.csdn.net/weixin_45317919/article/details/125899133
- 牛斯托洛夫斯坦，“雷达信号处理之恒虚警（CFAR）检测基础知识总结” https://zhuanlan.zhihu.com/p/652220176
- Khaled Zebiri;Amar Mezache，“Radar CFAR detection for multiple-targets situations for Weibull and log-normal distributed clutter”
- https://github.com/radarsimx/radarsimpy