厚积薄发：毫米波雷达开发手册之大话恒虚警率检测
严肃声明
严格禁止未经本人允许转载本文工作成果（包括引言、证明、代码分析等），一经发现被用于学术或商业用途，将予以法律警告。欢迎和笔者深度合作，探讨学术话题。
本人发现部分网站（如https://it.cha138.com/等野鸡网站）未经本人允许私自将本人博客内容转发公开，本人严禁杜绝此类行为，一经发现将投诉举报相关行为并追究违法责任！！！
写在前面
深知新手在接触毫米波雷达板硬件时需要花费的沉没成本，因此在行将告别毫米波雷达之际，总结这两年以来在毫米波雷达上的一些经验和教训。
本文档用于为实现基于AWR1243BOOST等单板毫米波雷达开发提供参考指南与解决方案，主要包括硬件配置、基础参数、信号模型、应用DEMO开发以及可深入研究方向思考等；为更好地匹配后续级联雷达应用的学习路线，在本手册中会尽可能同化单板雷达和级联雷达中的相关表述。
本指南作者信息：Xl，联系方式：xxxxx@zju.edu.cn。未经本人允许，请勿用于商业和学术用途。
希望后者在使用本指南时可以考虑引用作者在毫米波雷达旅途中的相关工作，如本文参考文献[1].
本章节为可深入研究方向思考章节之恒虚警率算法，主要讨论均值类、有序统计类、自动筛选类、逻辑类等恒虚警率检测算法，并给出了多目标和单目标场景下的不同算法关于SNR和虚警率的性能实例程序。
欢迎各位读者通过邮件形式与笔者交流讨论，本章节完整程序请私信笔者，希望使用本代码时能够提供一份引用和Star，以表示对笔者工作的尊重，谢谢！在后续将定时维护更新。
https://github.com/DingdongD/TDMA-MIMO
https://github.com/DingdongD/CFAR_Lib
往期内容：
一统江湖：毫米波雷达开发手册之大话线谱估计
炉火纯青：毫米波雷达开发手册之大话空间谱估计
登堂入室：毫米波雷达开发手册之信号模型
初出茅庐：毫米波雷达开发手册之基础参数
扬帆起航：毫米波雷达开发手册之硬件配置
眼观四海：自动驾驶&4D成像毫米波雷达 如今几何？
简介
恒虚警率检测算法是现代雷达信号处理和目标检测系统中较为关键的部分。恒虚警率检测算法是在系统中给检测策略提供检测阈值并且保持虚警概率稳定、不随外界杂波和干扰而变化的信号处理算法。(本部分内容摘自笔者的本科毕业设计)
目前，为满足不同特定环境下的检测需求，恒虚警率检测算法研究已经出现了多个研究方向，但恒虚警率检测算法的分类至今仍没有一个规范的国际国内标准。针对大多数恒虚警率检测算法面临在海杂波环境下多目标背景检测性能显著下降的问题，本节主要从经典恒虚警率检测算法、海杂波环境下恒虚警率检测算法、多目标背景恒虚警率检测算法这三者展开论述。
经典恒虚警率检测算法主要有均值类恒虚警率（ML CFAR）检测算法、有序统计类恒虚警率（OS CFAR）检测算法。均值类恒虚警率检测算法中较为经典的是单元平均恒虚警率（CA CFAR）检测算法，该算法基于均匀背景且各参考单元间独立同分布的假设，将检测单元两侧的参考单元样本平均值作为背景功率估计值[2]，能够在均匀背景下实现较好的目标检测性能，但在多目标背景中存在目标掩蔽效应以及在杂波边缘背景中存在虚警和漏警。为改进CA CFAR杂波边缘背景下的检测缺陷，Hansen提出的单元平均选大（GO CFAR）检测算法[3]，通过将较大参考窗中参考单元采样值的平均作为背景功率估计值，减少了杂波边缘强杂波区域发生漏警。Trunk提出的单元平均选小恒虚警率（SO CFAR）检测算法[4]，通过将较小参考窗中参考单元采样值的平均作为背景功率估计值，减少了杂波边缘弱杂波区域发生虚警，同时能够克服一侧参考窗中存在多目标干扰时带来的目标遮蔽效应，但当其两侧参考窗中存在多目标干扰时其检测性能略显逊色。为更好抑制两侧参考窗中存在多目标时的目标遮蔽效应，Rohling提出了有序统计量恒虚警率（OS CFAR）检测算法[5]，将排序后的参考窗中第k个单元功率值作为背景功率估计值，相比SO CFAR，该算法在两侧参考窗中存在多目标时具有更好的包容性，但当参考窗中干扰目标数目超过算法所能容纳最大干扰阈值N-k时其检测性能会开始衰减。此外，有序统计类检测算法的检测性能十分依赖选择序号k的值，k值大小可能会影响目标的命中率。
现有的关于海杂波环境下恒虚警率检测算法研究较少，本节主要从海杂波的非均匀性入手来论述在非均匀环境下恒虚警率检测算法的研究与发展现状。由于CA CFAR检测算法在非瑞利环境下无法克服海杂波的拖尾效应，因而Goldstein等提出了适合于对数正态分布杂波模型的对数-T恒虚警率（Log-t CFAR）检测算法[6]，通过引入辅助单元实现了检测门限的调整，该算法在对数正态分布和威布尔分布杂波环境下能够控制虚警概率来维持良好的目标检测性能。后来，Guida等人又提出了适合于威布尔分布杂波模型的最优线性无偏估计恒虚警率（BLUE CFAR）检测算法[7]，通过对数变换将Weibull概率密度函数退化为Gumbel概率密度函数，并利用位置尺度参数的最佳线性无偏估计调整检测门限，该算法具有一定的抵抗局部不均匀性能力。近年，刘怡等人依据信号的中心极限定理和对数压缩准则提出了综合恒虚警率（COMP CFAR）检测算法[8]，其中主要包含对信号的对数压缩处理和信号幅度平均值累积这两个步骤，该算法在非均匀环境中具有较好的鲁棒性，并且在多目标与杂波边缘背景中均具有很好的检测性能。
多目标背景恒虚警率检测算法，顾名思义，是研究当参考窗中存在一侧或两侧多目标干扰时如何有效提升检测算法虚警调节性能和目标分辨性能的问题。在有序统计类算法的基础上，Rickard和Dillard提出了删除平均（CMLD CFAR）[9]恒虚警率检测算法，在对参考窗中样本的功率值进行排序后来删除前r个较大的样本值，再对剩余的样本采样加权作为检测单元的杂波背景功率估计值，与OS CFAR相比，该算法在多目标背景中实现相同的检测性能时损失的信噪比更小。但OS CFAR与CMLD CFAR分别受参数k、r的限制，这极大程度上限制了雷达回波的实时运算，因此Gálvez等人提出了神经网络单元平均有序统计类恒虚警率（NNCAOS CFAR）检测算法[10][11]，通过神经网络对雷达回波进行先验分析，设置两个神经网络块NN1和NN2来搜索并判断背景属于多目标背景、杂波边缘背景、均匀背景的哪一种，根据背景类型从经典CFAR检测算法中选择最优策略完成门限估计，该算法既能提高回波的实时运算能力，又能保持检测性能不随杂波的均匀性变化而严重衰减。为降低经典CFAR检测算法的计算复杂度，ChiaHung等人提出了深度学习恒虚警率（DL CFAR）检测算法[12]，通过训练卷积神经网络来学习距离多普勒图中目标的结构，在剔除目标结构后得到具有纯噪声的距离多普勒图并估计门限水平，该算法在多目标背景下不同信杂比的情况能够保持鲁棒的检测性能与虚警调节能力，缓解了目标掩蔽效应。
均值类
均值类算法包括单元平均、单元平均选大、单元平均选小等经典算法。
%% CA-CFAR标称因子计算
function alpha = ca_threhold(Pfa,N)
% Pfa:虚警概率
% N：参考单元数
alpha = N .* (Pfa .^ (-1 / N) - 1 );
end
%% CA-CFAR恒虚警率检测算法的函数实现
function result = func_cfar_ca(x,alpha,NSlide,Pro_cell)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
if isempty(left)
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
end
T = zeros(1,len); %检测阈值
% target = java.util.LinkedList; %利用Java链表来实现目标的存储
target = [];
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
T(1,i) = mean(cell_right) * alpha;
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Slide - Half_Pro_cell:i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
T(1,i) = (mean(cell_left) + mean(cell_right)) / 2 * alpha; %求解门限
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
for i = right + 1:len %右边界
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i-Half_Pro_cell-1);
T(1,i) = mean(cell_left) * alpha;
if T(1,i) < x(i)
target = [target, i];
end
end
result = {'CFAR_CA',T,target};
end
%% SO-CFAR的虚警概率计算
function Pfa = SO_Pfa(alpha,N)
Pfa = 0;
n = ceil(N / 2);
for i = 0:n-1
Pfa = Pfa + 2 * gamma(n+i) ./ gamma(i+1) ./ gamma(n) .* (2 + alpha ./ n) .^ (-(n+i));
end
end
%% SO-CFAR恒虚警率检测算法的函数实现
function result = func_cfar_so(x,alpha,NSlide,Pro_cell)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
end
T = zeros(1,len); %检测阈值
% target = java.util.LinkedList; %利用Java链表来实现目标的存储
target = [];
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
T(1,i) = mean(cell_right) * alpha;
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Slide - Half_Pro_cell:i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
T(1,i) = min(mean(cell_left),mean(cell_right)) * alpha; %求解门限
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
for i = right + 1:len %右边界
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i-Half_Pro_cell-1);
T(1,i) = mean(cell_left) * alpha;
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
result = {'CFAR_SO',T,target}; %采用字典类型
end
%% GO-CFAR标称因子计算
function [alpha] = go_threhold(Pfa,N)
% Pfa:虚警概率
% N：参考单元数
scope = [0,100]; %区间
precision = 0.01; %精度
func = @GO_Pfa;
parameter = [N,Pfa];
[alpha] = binary_solution(scope,precision,func,parameter);
end
%% GO-CFAR恒虚警率检测算法的函数实现
function result = func_cfar_go(x,alpha,NSlide,Pro_cell)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
end
T = zeros(1,len); %检测阈值
target = [];
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
T(1,i) = mean(cell_right) * alpha;
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Slide - Half_Pro_cell:i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
T(1,i) = max(mean(cell_left),mean(cell_right)) * alpha; %求解门限
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
for i = right + 1:len %右边界
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i-Half_Pro_cell-1);
T(1,i) = mean(cell_left) * alpha;
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
result = {'CFAR_GO',T,target}; %采用字典类型
end
有序统计类
有序统计类算法包括有序统计、有序统计选大、有序统计选小等经典算法。
%% OS-CFAR的虚警概率计算
function Pfa = OS_Pfa(alpha,N,rate)
k = ceil(N .* rate);
Pfa = gamma(N+1) .* gamma(N-k+alpha+1) ./ gamma(N-k+1) ./ gamma(N+alpha+1);
end
%% OS-CFAR恒虚警率检测算法的函数实现
function result = func_cfar_os(x,alpha,NSlide,Pro_cell,rate)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
%rate:有序比例，OS中选择序号如k=3/4N时rate=0.75
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
persistent K;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
K = round(rate * NSlide); %有序序号
end
T = zeros(1,len); %检测阈值
target = [];
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
cell_right = sort(cell_right);
T(1,i) = cell_right(K) * alpha;
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Slide - Half_Pro_cell:i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
cell = sort([cell_left,cell_right]); %合并左右侧滑动窗
T(1,i) = cell(K) * alpha; %求解门限
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
for i = right + 1:len %右边界
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i-Half_Pro_cell-1);
cell_left = sort(cell_left);
T(1,i) = cell_left(K) * alpha;
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
result = {'CFAR_OS',T,target}; %采用字典类型
end
%% OSGO-CFAR的虚警概率计算
function Pfa = OSGO_Pfa(alpha,N,rate)
k = ceil(N * rate / 2);
pfa = 0;
for j = 0:N/2-k
for i = 0:N/2-k
pfa = pfa + nchoosek(N/2-k,j)*nchoosek(N/2-k,i)*((-1)^(N-2*k-j-i)*gamma(N-j-i)*gamma(alpha+1))/((N/2-i)*gamma(N-j-i+alpha+1));
end
end
Pfa = pfa * 2 * k * k * nchoosek(N/2,k) * nchoosek(N/2,k);
end
%% 本程序主要实现osgo-CFAR恒虚警率检测算法的函数形式
function result = func_cfar_osgo(x,alpha,NSlide,Pro_cell,rate)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
%rate:有序比例，OS中选择序号如k=3/4N时rate=0.75
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
persistent K;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
K = round(rate * Half_Slide); %有序序号
end
T = zeros(1,len); %CMLD检测阈值
target = java.util.LinkedList; %利用Java链表来实现目标的存储
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
cell_right = sort(cell_right);
T(1,i) = cell_right(K) * alpha; %OSGO CFAR的门限
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Pro_cell - Half_Slide : i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
cell_left = sort(cell_left);
cell_right = sort(cell_right);
T(1,i) = max(cell_left(K),cell_right(K)) * alpha; %OSGO CFAR的门限
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = right + 1 : len
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i - Half_Pro_cell - 1);
cell_left = sort(cell_left);
T(1,i) = cell_left(K) * alpha; %OSGO CFAR的门限
if T(1,i) < x(i)
target.add(i);
end
end
result = {'CFAR_OSGO',T,target}; %采用字典类型
end
%% OSSO-CFAR的虚警概率计算
function Pfa = OSSO_Pfa(alpha,N,rate)
k = ceil(N * rate / 2);
pfa = 0;
for j = 0:N/2-k
for i = 0:N/2-k
pfa = pfa + nchoosek(N/2-k,j)*nchoosek(N/2-k,i)*((-1)^(N-2*k-j-i)*gamma(N-j-i)*gamma(alpha+1))/((N/2-i)*gamma(N-j-i+alpha+1));
end
end
Pfa = 2 * k * nchoosek(N/2,k) * ((gamma(k) * gamma(alpha+N/2-k+1)/gamma(alpha+N/2+1)) - k * nchoosek(N/2,k) * pfa);
end
%% 本程序主要实现osso-CFAR恒虚警率检测算法的函数形式
function result = func_cfar_osso(x,alpha,NSlide,Pro_cell,rate)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
%rate:有序比例，OS中选择序号如k=3/4N时rate=0.75
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
persistent K;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
K = round(rate * Half_Slide); %有序序号
end
T = zeros(1,len); %CMLD检测阈值
target = java.util.LinkedList; %利用Java链表来实现目标的存储
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
cell_right = sort(cell_right);
T(1,i) = cell_right(K) * alpha; %OSSO CFAR的门限
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Pro_cell - Half_Slide : i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
cell_left = sort(cell_left);
cell_right = sort(cell_right);
T(1,i) = min(cell_left(K),cell_right(K)) * alpha; %OSSO CFAR的门限
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = right + 1 : len
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i - Half_Pro_cell - 1);
cell_left = sort(cell_left);
T(1,i) = cell_left(K) * alpha; %OSSO CFAR的门限
if T(1,i) < x(i)
target.add(i);
end
end
result = {'CFAR_OSSO',T,target}; %采用字典类型
end
自动筛选类
自动筛选类算法包括CMLD、TM、ICCA/ICOS（自适应）、Switch-OS等经典算法。
%% CMLD-CFAR的虚警概率计算
function Pfa = CMLD_Pfa(alpha,N,rate)
k = ceil(N * rate);
Pfa = 0;
for i = 1:N-k
ai = nchoosek(N,k) * nchoosek(N-k,i-1) * (-1)^(i-1) * ((N-i+1-k)/k)^(N-k-1);
ci = (N-i+1) / (N-k-i+1);
Pfa = Pfa + ai / (ci + alpha);
end
end
%% 本程序主要实现CMLD-CFAR恒虚警率检测算法的函数形式
function result = func_cfar_cmld(x,alpha,NSlide,Pro_cell,rate)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
%rate:有序比例，OS中选择序号如k=3/4N时rate=0.75
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
persistent K;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
K = round(rate * NSlide); %有序序号
end
T = zeros(1,len); %CMLD检测阈值
target = [];
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
cell_right = sort(cell_right);
cell_right(end-K+1:end) = []; %删除第K个单元
T(1,i) = mean(cell_right) * alpha; %CMLD CFAR的门限
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Pro_cell - Half_Slide : i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
cell = [cell_left,cell_right];
cell = sort(cell); %对参考单元排序
cell(end-K+1:end) = []; %删除两侧参考窗排序后的第K个单元
T(1,i) = (mean(cell)) * alpha; %CMLD CFAR的门限
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
for i = right + 1 : len
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i - Half_Pro_cell - 1);
cell_left = sort(cell_left);
cell_left(end-K+1:end) = [];
T(1,i) = mean(cell_left) * alpha; %CMLD CFAR的门限
if T(1,i) < x(i)
target = [target, i]; %加入目标
end
end
result = {'CFAR_CMLD',T,target}; %采用字典类型
end
%% TM-CFAR的虚警概率计算
function Pfa = TM_Pfa(alpha,N,rate1,rate2)
k1 = ceil(N * rate1);
k2 = ceil(N * rate2);
init = 0;
for j = 0:k1
init = init + (N - k1 - k2) * nchoosek(k1,j) * (-1)^(k1-j) /(N - j + alpha);
end
vi1 = factorial(N) / (factorial(k1) * factorial(N - k1 - 1) * (N - k1 - k2)) * init;
vi2 = 1;
for i = 2:N - k1 - k2
ai = (N - k1 - i + 1) / (N - k1 - k2 - i + 1);
vi2 = vi2 * ai / (ai + alpha / ( N - k1 - k2));
end
Pfa = vi1 * vi2;
end
%% 本程序主要实现TM-CFAR恒虚警率检测算法的函数形式
function result = func_cfar_tm(x,alpha,NSlide,Pro_cell,rate1,rate2)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
%rate:有序比例，OS中选择序号如k=3/4N时rate=0.75
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
persistent K1;
persistent K2;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
K1 = round(rate1 * NSlide); %有序序号
K2 = round(rate2 * NSlide); %有序序号
end
T = zeros(1,len); %CMLD检测阈值
target = java.util.LinkedList; %利用Java链表来实现目标的存储
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
cell_right = sort(cell_right);
cell_right(1:K1) = [];
cell_right(end-K2+1:end) = [];
T(1,i) = mean(cell_right) * alpha; %TM CFAR的门限
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Pro_cell - Half_Slide : i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
cell = [cell_left,cell_right];
cell = sort(cell);
cell(1:K1) = [];
cell(end+1-K2:end) = [];
T(1,i) = mean(cell) * alpha; %TM CFAR的门限
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = right + 1 : len
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i - Half_Pro_cell - 1);
cell_left = sort(cell_left);
cell_left(1:K1) = [];
cell_left(end-K2+1:end) = [];
T(1,i) = mean(cell_left) * alpha; %TM CFAR的门限
if T(1,i) < x(i)
target.add(i);
end
end
result = {'CFAR_TM',T,target}; %采用字典类型
end
%% ICCA-CFAR恒虚警率检测算法的函数实现
function result = func_cfar_icca(x,Pfa,NSlide,Pro_cell)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
end
T = zeros(1,len); %检测阈值
target = java.util.LinkedList; %利用Java链表来实现目标的存储
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
alpha = 1 - Pfa^(1 / (NSlide-1)); %初始alpha0
for time = 1:30 %IC迭代更新次数
t = sum(cell_right) * alpha;%阈值系数提取
cell_right = cell_right(cell_right <= t);
if 1 - Pfa^(1 / (length(cell_right)-1)) == alpha
break; %如果alpha 更新仍相等,说明迭代收敛
else
alpha = 1 - Pfa^(1 / (length(cell_right)-1)); %否则更新alpha
end
end
T(1,i) = sum(cell_right) * alpha;
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Slide - Half_Pro_cell:i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
cell = [cell_left,cell_right]; %合并两侧窗
alpha = 1 - Pfa^(1 / (NSlide-1)); %初始alpha0
for time = 1:30 %IC迭代更新次数
t = sum(cell) * alpha;%阈值系数提取
cell = cell(cell <= t);
if 1 - Pfa^(1 / (length(cell)-1)) == alpha
break; %如果alpha 更新仍相等,说明迭代收敛
else
alpha = 1 - Pfa^(1 / (length(cell)-1)); %否则更新alpha
end
end
T(1,i) = sum(cell) * alpha; %求解门限
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = right + 1:len %右边界
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i-Half_Pro_cell-1);
alpha = 1 - Pfa^(1 / (NSlide-1)); %初始alpha0
for time = 1:30 %IC迭代更新次数
t = sum(cell_left) * alpha;%阈值系数提取
cell_left = cell_left(cell_left <= t);
if 1 - Pfa^(1 / (length(cell_left)-1)) == alpha
break; %如果alpha 更新仍相等,说明迭代收敛
else
alpha = 1 - Pfa^(1 / (length(cell_left)-1)); %否则更新alpha
end
end
T(1,i) = sum(cell_left) * alpha;
if T(1,i) < x(i)
target.add(i);
end
end
result = {'CFAR_ICCA',T,target}; %采用字典类型
end
%% ICOS-CFAR恒虚警率检测算法的函数实现
function result = func_cfar_icos(x,Pfa,NSlide,Pro_cell,rate)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
%rate:有序比例，OS中选择序号如k=3/4N时rate=0.75
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
end
T = zeros(1,len); %检测阈值
target = java.util.LinkedList; %利用Java链表来实现目标的存储
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
alpha = os_threhold(Pfa,length(cell_right),rate); %初始alpha
for time = 1:30 %IC迭代更新次数
cell_right = sort(cell_right);
t = cell_right(ceil(rate*length(cell_right))) * alpha; %阈值系数提取
cell_right = cell_right(cell_right <= t);
if os_threhold(Pfa,length(cell_right),rate) == alpha
break; %如果alpha 更新仍相等,说明迭代收敛
else
alpha = os_threhold(Pfa,length(cell_right),rate); %否则更新alpha
end
end
T(1,i) = cell_right(ceil(rate*length(cell_right))) * alpha;
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Slide - Half_Pro_cell:i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
cell = [cell_left,cell_right]; %合并两侧窗
alpha = os_threhold(Pfa,length(cell),rate); %初始alpha
for time = 1:30 %IC迭代更新次数
cell = sort(cell);
t = cell(ceil(rate*length(cell))) * alpha; %阈值系数提取
cell = cell(cell <= t);
if os_threhold(Pfa,length(cell),rate) == alpha
break; %如果alpha 更新仍相等,说明迭代收敛
else
alpha = os_threhold(Pfa,length(cell),rate); %否则更新alpha
end
end
T(1,i) = cell(ceil(rate*length(cell))) * alpha; %求解门限
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = right + 1:len %右边界
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i-Half_Pro_cell-1);
alpha = os_threhold(Pfa,length(cell_left),rate); %初始alpha
for time = 1:30 %IC迭代更新次数
cell_left = sort(cell_left);
t = cell_left(ceil(rate*length(cell_left))) * alpha; %阈值系数提取
cell_left = cell_left(cell_left <= t);
if os_threhold(Pfa,length(cell_left),rate) == alpha
break; %如果alpha 更新仍相等,说明迭代收敛
else
alpha = os_threhold(Pfa,length(cell_left),rate); %否则更新alpha
end
end
T(1,i) = cell_left(ceil(rate*length(cell_left))) * alpha;
if T(1,i) < x(i)
target.add(i);
end
end
result = {'CFAR_ICOS',T,target}; %采用字典类型
end
详细原理可以见笔者的深入分析：恒虚警率检测算法之Switch-CFAR
%% 本程序主要实现S-CFAR恒虚警率检测算法的函数形式
function result = func_cfar_sos(x,beta,NSlide,Nt,Pro_cell,rate)
%x:原始杂波数据
%alpha：标称因子
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
%Nt:整数门限
%rate:选择序号
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
K = rate * NSlide; %选择序号
end
T = zeros(1,len); %CMLD检测阈值
target = java.util.LinkedList; %利用Java链表来实现目标的存储
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
S0 = cell_right(cell_right <= mean(cell_right));
S1 = cell_right(cell_right > mean(cell_right));
if length(S0) > Nt
S0 = sort(S0);
T(1,i) = S0(K) * beta;
else
T(1,i) = mean(cell_right) * beta;
end
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Pro_cell - Half_Slide : i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
cell = [cell_left,cell_right];
S0 = cell(cell <= mean(cell));
S1 = cell(cell > mean(cell));
if length(S0) > Nt
S0 = sort(S0);
T(1,i) = S0(K) * beta;
else
T(1,i) = mean(cell) * beta;
end
if T(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = right + 1 : len
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i - Half_Pro_cell - 1);
S0 = cell_left(cell_left <= mean(cell_left));
S1 = cell_left(cell_left > mean(cell_left));
if length(S0) > Nt
S0 = sort(S0);
T(1,i) = S0(K) * beta;
else
T(1,i) = mean(cell_left) * beta;
end
if T(1,i) < x(i)
target.add(i);
end
end
result = {'CFAR_SOS',T,target}; %采用字典类型
end
逻辑类CFAR
逻辑类算法包括AND/OR等经典算法。
%% AND-CFAR标称因子计算
function [alpha] = and_threhold(Pfa,N,rate)
% Pfa:虚警概率
% N：参考单元数
% alpha：门限系数，需要考虑CA和OS
alpha = zeros(1,2);
alpha(1,1) = ca_threhold(Pfa,N);
alpha(1,2) = os_threhold(Pfa,N,rate);
end
%% AND-CFAR恒虚警率检测算法的函数实现
function result = func_cfar_and(x,alpha,NSlide,Pro_cell,rate)
%x:原始杂波数据
%alpha：标称因子,alpha(1,1)用于CA，alpha(1,,2)用于OS
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
%rate:有序比例，OS中选择序号如k=3/4N时rate=0.75
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
persistent K;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
K = round(rate * NSlide); %有序序号
end
T_CA = zeros(1,len); %CA检测阈值
T_OS = zeros(1,len); %OS检测阈值
T = zeros(1,len); %AND检测阈值
target = java.util.LinkedList; %利用Java链表来实现目标的存储
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
cell_right = sort(cell_right);
T_CA(1,i) = mean(cell_right) * alpha(1,1);
T_GO(1,i) = cell_right(K) * alpha(1,2);
T(1,i) = max(T_CA(1,i),T_GO(1,i)); %AND CFAR的门限
if T_CA(1,i) < x(i) && T_GO(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Slide - Half_Pro_cell:i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
cell = sort([cell_left,cell_right]); %合并左右侧滑动窗
T_CA(1,i) = (mean(cell_left) + mean(cell_right)) / 2 * alpha(1,1);
T_GO(1,i) = cell(K) * alpha(1,2);
T(1,i) = max(T_CA(1,i),T_GO(1,i)); %AND CFAR的门限
if T_CA(1,i) < x(i) && T_GO(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = right + 1:len %右边界
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i-Half_Pro_cell-1);
cell_left = sort(cell_left);
T_CA(1,i) = mean(cell_left) * alpha(1,1);
T_GO(1,i) = cell_left(K) * alpha(1,2);
T(1,i) = max(T_CA(1,i),T_GO(1,i)); %AND CFAR的门限
if T_CA(1,i) < x(i) && T_GO(1,i) < x(i)
target.add(i);
end
end
result = {'CFAR_AND',T,target}; %采用字典类型
end
%% OR-CFAR标称因子计算
function [alpha] = or_threhold(Pfa,N,rate)
% Pfa:虚警概率
% N：参考单元数
% alpha：门限系数，需要考虑CA和OS两种
alpha = zeros(1,2);
alpha(1,1) = ca_threhold(Pfa,N);
alpha(1,2) = os_threhold(Pfa,N,rate);
end
%% OR-CFAR恒虚警率检测算法的函数实现
function result = func_cfar_or(x,alpha,NSlide,Pro_cell,rate)
%x:原始杂波数据
%alpha：标称因子,alpha(1,1)用于CA，alpha(1,,2)用于OS
%Nslide：滑动窗大小
%Pro_cell:保护单元大小
%rate:有序比例，OS中选择序号如k=3/4N时rate=0.75
persistent left;
persistent right;
persistent Half_Slide;
persistent Half_Pro_cell;
persistent len;
persistent K;
if isempty(left)
left = 1 + Half_Pro_cell + Half_Slide; %设置左边界
right = length(x) - Half_Pro_cell - Half_Slide; %设置右边界
Half_Slide = NSlide / 2; %半滑动窗
Half_Pro_cell = Pro_cell / 2; %一侧保护单元长度
len = length(x); %杂波单元
K = round(rate * NSlide); %有序序号
end
T_CA = zeros(1,len); %CA检测阈值
T_OS = zeros(1,len); %OS检测阈值
T = zeros(1,len); %OR检测阈值
target = java.util.LinkedList; %利用Java链表来实现目标的存储
for i = 1:left - 1 %左边界
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Slide * 2 + Half_Pro_cell);
cell_right = sort(cell_right);
T_CA(1,i) = mean(cell_right) * alpha(1,1);
T_GO(1,i) = cell_right(K) * alpha(1,2);
T(1,i) = min(T_CA(1,i),T_GO(1,i)); %OR CFAR的门限
if T_CA(1,i) < x(i) || T_GO(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = left:right %中间区域
cell_left = x(1,i - Half_Slide - Half_Pro_cell:i - Half_Pro_cell - 1);
cell_right = x(1,i + Half_Pro_cell + 1 : i + Half_Pro_cell + Half_Slide);
cell = sort([cell_left,cell_right]); %合并左右侧滑动窗
T_CA(1,i) = (mean(cell_left) + mean(cell_right)) / 2 * alpha(1,1);
T_GO(1,i) = cell(K) * alpha(1,2);
T(1,i) = min(T_CA(1,i),T_GO(1,i)); %OR CFAR的门限
if T_CA(1,i) < x(i) || T_GO(1,i) < x(i)
target.add(i); %加入目标
end
end
for i = right + 1:len %右边界
cell_left = x(1,i - Half_Pro_cell - Half_Slide * 2 : i-Half_Pro_cell-1);
cell_left = sort(cell_left);
T_CA(1,i) = mean(cell_left) * alpha(1,1);
T_GO(1,i) = cell_left(K) * alpha(1,2);
T(1,i) = min(T_CA(1,i),T_GO(1,i)); %OR CFAR的门限
if T_CA(1,i) < x(i) || T_GO(1,i) < x(i)
target.add(i);
end
end
result = {'CFAR_OR',T,target}; %采用字典类型
end
算法信噪比和ROC分析
%% 本程序主要研究多目标情况下不同CFAR算法关于信噪比的检测性能
clc;clear;close all;
%% 参数初始化
N = 64; %滑动窗长度
n = N / 2; %前后沿参考单元长度
Pfa = 1e-6; %虚警概率
SNR_dB = 1:1:35; %信噪比
SNR = 10 .^ (SNR_dB / 10);
Lth = length(SNR); %仿真信噪比区间长度
monte_num = 1e4; % 蒙特卡洛仿真次数，值越大越平滑
rate = 3 / 4;
k = 3 * N / 4; % 有序选择序号
tm_rate1 = 1 / 32;
tm_rate2 = 1 / 16;
cmld_rate = 1 / 4;
%% 门限系数估计
T_CA = ca_threhold(Pfa,N); %CA-CFAR门限因子
T_GO = go_threhold(Pfa,N); %GO-CFAR门限因子
T_SO = so_threhold(Pfa,N); %SO-CFAR门限因子
T_OS = os_threhold(Pfa,N,rate); % OS-CFAR 的门限因子
T_TM = tm_threhold(Pfa,N,tm_rate1,tm_rate2); % TM-CFAR的门限因子
T_CMLD = cmld_threhold(Pfa,N,cmld_rate); % CMLD-CFAR的门限因子
% 检测概率矩阵
Pd_CA = zeros(1, Lth);
Pd_SO = zeros(1, Lth);
Pd_GO = zeros(1, Lth);
Pd_OS = zeros(1, Lth);
Pd_TM = zeros(1, Lth);
Pd_CMLD = zeros(1, Lth);
%% 多目标程序（此处设置了3个目标，可根据自己的需求改）
for i = 1:Lth
detect_ca = 0; %设置CA检测信号数目变量
detect_so = 0; %设置SO检测信号数目变量
detect_go = 0; %设置GO检测信号数目变量
detect_os = 0; %设置OS检测信号数目变量
detect_tm = 0; %设置TM检测信号数目变量
detect_cmld = 0; %设置CMLD检测信号数目变量
for j = 1:monte_num
lambda = 1;
u = rand(1,N);
exp_noise = log(u) * (-lambda);
lambda = SNR(i) + 1;
u = rand(1,3); % 产生目标
exp_target = log(u(1)) * (-lambda); % 目标
% 注意 此处可以根据需求选择干扰的位置放在哪里 和 控制干扰的数目 使用时注释掉其他部分
%=================================================================================
InterPos = randi([1+n,N],1,2); %后沿参考窗产生两个干扰
exp_noise(InterPos(1)) = log(u(2)) * (-lambda);
exp_noise(InterPos(2)) = log(u(3)) * (-lambda);
%=================================================================================
% InterPos = randi([1,n],1,2); %前沿参考窗产生两个干扰
% exp_noise(InterPos(1)) = log(u(2)) * (-lambda);
% exp_noise(InterPos(2)) = log(u(3)) * (-lambda);
%=================================================================================
% InterPos1 = randi([1,n],1,1); %前沿参考窗产生1个干扰
% InterPos2 = randi([1+n,N],1,1); %后沿参考窗产生1个干扰
% exp_noise(InterPos(1)) = log(u(2)) * (-lambda);
% exp_noise(InterPos(2)) = log(u(3)) * (-lambda);
%=================================================================================
% 不同检测器的估计
% CA-CFAR
cfar_ca = exp_target / mean(exp_noise);
if cfar_ca > T_CA
detect_ca = detect_ca + 1;
end
% GO-CFAR
cfar_go = exp_target / max(mean(exp_noise(1:N/2)),mean(exp_noise(N/2+1:N)));
if cfar_go > T_GO
detect_go = detect_go + 1;
end
% SO-CFAR
cfar_so = exp_target / min(mean(exp_noise(1:N/2)),mean(exp_noise(N/2+1:N)));
if cfar_so > T_SO
detect_so = detect_so + 1;
end
% OS-CFAR
os_noise = sort(exp_noise);
cfar_os = exp_target / os_noise(k);
if cfar_os > T_OS
detect_os = detect_os + 1;
end
% TM-CFAR
tm_noise = sort(exp_noise);
tm_noise(1:N*tm_rate1) = [];
tm_noise(end+1-N*tm_rate2:end) = [];
cfar_tm = exp_target / mean(tm_noise);
if cfar_tm > T_TM
detect_tm = detect_tm + 1;
end
% CMLD-CFAR
cmld_noise = sort(exp_noise);
cmld_noise(end-cmld_rate*N+1:end) = [];
cfar_cmld = exp_target / mean(cmld_noise);
if cfar_cmld > T_CMLD
detect_cmld = detect_cmld + 1;
end
end
Pd_CA(1,i) = detect_ca / monte_num;
Pd_SO(1,i) = detect_so / monte_num;
Pd_GO(1,i) = detect_go / monte_num;
Pd_OS(1,i) = detect_os / monte_num;
Pd_TM(1,i) = detect_tm / monte_num;
Pd_CMLD(1,i) = detect_cmld / monte_num;
end
plot(SNR_dB,Pd_CA,'Color','#FFD700','LineWidth',2,'Marker','<');hold on;
plot(SNR_dB,Pd_GO,'Color','#DC143C','LineWidth',2,'Marker','+');hold on;
plot(SNR_dB,Pd_SO,'Color','#FF1493','LineWidth',2,'Marker','*');hold on;
plot(SNR_dB,Pd_OS,'Color','#8B008B','LineWidth',2,'Marker','v');hold on;
plot(SNR_dB,Pd_TM,'Color','#483D8B','LineWidth',2,'Marker','^');hold on;
plot(SNR_dB,Pd_CMLD,'Color','#0000FF','LineWidth',2,'Marker','o');hold on;
grid minor;
xlabel('\fontname{Times New Roman}SNR/dB');
ylabel('\fontname{宋体}检测概率\fontname{Times New Roman}Pd');
h = legend('CA-CFAR','GO-CFAR','SO-CFAR','OS-CFAR','TM-CFAR','CMLD-CFAR','Location','SouthEast','NumColumns',1);
set(h,'edgecolor','none');
clc;clear;close all;
% 本程序主要研究多目标情况下不同CFAR算法关于信噪比的检测性能
%% 参数初始化
N = 64; %滑动窗长度
n = N / 2; %前后沿参考单元长度
Pfa = [1e-6:2e-6:9e-6, 1e-5:2e-5:9e-5,1e-4:2e-4:9e-4,1e-3:2e-3:9e-3,1e-2:2e-2:9e-2,1e-1:2e-2:9e-2]; %虚警概率
SNR_dB = 15; %信噪比
SNR = 10 .^ (SNR_dB / 10);
Lth = length(Pfa); %仿真信噪比区间长度
monte_num = 1e5;
k = 3 * N / 4; %有序样本选择序号
T_CA = zeros(1, Lth);
T_GO = zeros(1, Lth);
T_SO = zeros(1, Lth);
T_OS = zeros(1, Lth);
for i = 1 : Lth
T_CA(1,i) = Pfa(i) ^ (-1 / N) - 1; %门限因子
T_GO(1,i) = Pfa(i) ^ (-1 / N) - 1; %门限因子
T_SO(1,i) = Pfa(i) ^ (-1 / N) - 1; %门限因子
syms T_os
g = Pfa(i) - k * nchoosek(N,k) * gamma(k) * gamma(N-k+1+T_os) / gamma(N+T_os +1);
x = solve(g);
T_os = double(x);
T_OS(1,i) = T_os(T_os == abs(T_os));
end
% 蒙特卡洛仿真
Pd_CA = zeros(1, Lth);
Pd_SO = zeros(1, Lth);
Pd_GO = zeros(1, Lth);
Pd_OS = zeros(1, Lth);
%% 多目标程序
for i = 1:Lth
detect_ca = 0; % 设置CA检测信号数目变量
detect_so = 0; % 设置SO检测信号数目变量
detect_go = 0; % 设置GO检测信号数目变量
detect_os = 0; % 设置OS检测信号数目变量
for j = 1:monte_num
lambda = 1;
u = rand(1,N-5);
exp_noise = log(u) * (-lambda);
lambda = SNR + 1;
u = rand(1,6); %
exp_target = log(u(1)) * (-lambda); % 目标
exp_noise(N-4) = log(u(2)) * (-lambda); % 干扰目标1
exp_noise(N-3) = log(u(3)) * (-lambda); % 干扰目标2
exp_noise(N-2) = log(u(4)) * (-lambda); % 干扰目标3
exp_noise(N-1) = log(u(5)) * (-lambda); % 干扰目标4
exp_noise(N) = log(u(6)) * (-lambda); % 干扰目标5
cfar_ca = exp_target / sum(exp_noise);
if cfar_ca > T_CA(i)
detect_ca = detect_ca + 1;
end
cfar_go = exp_target / max(sum(exp_noise(1:N/2)),sum(exp_noise(N/2+1:N)));
if cfar_go > T_GO(i)
detect_go = detect_go + 1;
end
cfar_so = exp_target / min(sum(exp_noise(1:N/2)),sum(exp_noise(N/2+1:N)));
if cfar_so > T_SO(i)
detect_so = detect_so + 1;
end
exp_noise = sort(exp_noise);
cfar_os = exp_target / exp_noise(k);
if cfar_os > T_OS(i)
detect_os = detect_os + 1;
end
end
Pd_CA(1,i) = detect_ca / monte_num;
Pd_SO(1,i) = detect_so / monte_num;
Pd_GO(1,i) = detect_go / monte_num;
Pd_OS(1,i) = detect_os / monte_num;
end
% %% 单目标程序
% for i = 1:Lth
% detect_ca = 0; %设置CA检测信号数目变量
% detect_so = 0; %设置SO检测信号数目变量
% detect_go = 0; %设置GO检测信号数目变量
% detect_os = 0; %设置OS检测信号数目变量
% for j = 1:monte_num
% lambda = 1;
% u = rand(1,N);
% exp_noise = log(u) * (-lambda);
% lambda = SNR + 1;
% u = rand(1,1); %
% exp_target = log(u(1)) * (-lambda); % 目标
%
% cfar_ca = exp_target / sum(exp_noise);
% if cfar_ca > T_CA(i)
% detect_ca = detect_ca + 1;
% end
%
% cfar_go = exp_target / max(sum(exp_noise(1:N/2)),sum(exp_noise(N/2+1:N)));
% if cfar_go > T_GO(i)
% detect_go = detect_go + 1;
% end
%
% cfar_so = exp_target / min(sum(exp_noise(1:N/2)),sum(exp_noise(N/2+1:N)));
% if cfar_so > T_SO(i)
% detect_so = detect_so + 1;
% end
%
% exp_noise = sort(exp_noise);
% cfar_os = exp_target / exp_noise(k);
% if cfar_os > T_OS(i)
% detect_os = detect_os + 1;
% end
%
% end
% Pd_CA(1,i) = detect_ca / monte_num;
% Pd_SO(1,i) = detect_so / monte_num;
% Pd_GO(1,i) = detect_go / monte_num;
% Pd_OS(1,i) = detect_os / monte_num;
% end
semilogx(Pfa,Pd_CA,'rd-','LineWidth',1.5);hold on;
semilogx(Pfa,Pd_GO,'g-*','LineWidth',1.5);hold on;
semilogx(Pfa,Pd_SO,'b-.','LineWidth',1.5);hold on;
semilogx(Pfa,Pd_OS,'k-x','LineWidth',1.5);hold on;
grid minor;
xlabel('\fontname{宋体}虚警概率');
ylabel('\fontname{宋体}检测概率\fontname{Times New Roman}Pd');
h = legend('CA-CFAR','GO-CFAR','SO-CFAR','OS-CFAR','Location','SouthEast','NumColumns',1);
set(h,'edgecolor','none');
参考文献
[1] X. Yu, Z. Cao, Z. Wu, C. Song, J. Zhu and Z. Xu, "A Novel Potential Drowning Detection System Based on Millimeter-Wave Radar," 2022 17th International Conference on Control, Automation, Robotics and Vision (ICARCV), Singapore, Singapore, 2022, pp. 659-664, doi: 10.1109/ICARCV57592.2022.10004245.
[2] 付康. 复杂背景下雷达恒虚警检测方法研究 [硕士]: 西安电子科技大学; 2019.
[3] Hansen V. Constant false alarm rate processing in search radars. In; 1973.
[4] Trunk GV. Range Resolution of Targets Using Automatic Detectors. IEEE Transactions on Aerospace and Electronic Systems 1978,AES-14:750-755.
[5] Rohling H. Radar CFAR Thresholding in Clutter and Multiple Target Situations. IEEE Transactions on Aerospace and Electronic Systems 1983,AES-19:608-621.
[6] Goldstein GB. False-Alarm Regulation in Log-Normal and Weibull Clutter. IEEE Transactions on Aerospace and Electronic Systems 1973,AES-9:84-92.
[7] Guida M, Longo M, Lops M. Biparametric linear estimation for CFAR against Weibull clutter. IEEE Transactions on Aerospace and Electronic Systems 1992,28:138-151.
[8] Liu Y, Zhang S, Suo J, Zhang J, Yao T. Research on a New Comprehensive CFAR (Comp-CFAR) Processing Method. IEEE Access 2019,7:19401-19413.
[9] Rickard JT, Dillard GM. Adaptive Detection Algorithms for Multiple-Target Situations. IEEE Transactions on Aerospace and Electronic Systems 1977,AES-13:338-343.
[10] Gálvez N, Cousseau J, Pasciaroni Jl, Agamennoni O. Improved Neural Network Based CFAR Detection for non Homogeneous Background and Multiple Target Situations. Latin American applied research Pesquisa aplicada latino americana = Investigación aplicada latinoamericana 2012,42.
[11] Akhtar J, Olsen KE. A Neural Network Target Detector with Partial CA-CFAR Supervised Training. In: 2018 International Conference on Radar (RADAR); 2018. pp. 1-6.
[12] Lin C, Lin Y, Bai Y, Chung W, Lee T, Huttunen H. DL-CFAR: A Novel CFAR Target Detection Method Based on Deep Learning. In: 2019 IEEE 90th Vehicular Technology Conference (VTC2019-Fall); 2019. pp. 1-6.
