本文汇集CA、SO、GO、OS、杂波图等恒虚警算法的门限因子求解方法及其matlab实现程序。
CA-CFAR
非常简单,可以直接求解。
%% 均值恒虚警_门限因子计算公式
%% 版本:v1
%% 时间:2019.11.08
%% 终版【不在优化】
function [ alpha ] = form_ALPHA_ca( PFA,N )
%FORM_ALPHA_CA 此处显示有关此函数的摘要
% PFA:虚警概率
% N:参考单元个数
alpha=N.*(PFA.^(-1./N)-1);
end
SO-CFAR
这个门限因子求解是一个非常复杂的高阶函数,基本上无法通过反函数求解,这里使用二分求解法解决。
PFA计算函数
%% 虚警概率_均值选小
%% 时间:2019.11.27
%% 版本:v1【不在优化】
function [ PFA ] = form_PFA_so( ALPHA,N )
%FORM_PFA_GO 此处显示有关此函数的摘要
% ALPHA:门限因子
% N:参考窗个数
PFA=0;
n=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
GO-CFAR
PFA计算函数
%% 虚警概率_均值选大
%% 时间:2019.11.27
%% 版本:v1【不在优化】
function [ PFA ] = form_PFA_go( ALPHA,N )
%FORM_PFA_GO 此处显示有关此函数的摘要
% ALPHA:门限因子
% N:参考窗个数
PFA=0;
n=N/2;
for i=0:n-1
PFA=PFA-2.*gamma(n+i)./gamma(i+1)./gamma(n).*(2+ALPHA./n).^(-(n+i));
end
PFA=PFA+2.*(1+ALPHA./n).^(-n);
end
OS-CFAR
%% 虚警概率_有序
%% 时间:2019.11.27
%% 版本:v1【不在优化】
function [ PFA ] = form_PFA_os( ALPHA,N,Rate )
%FORM_PFA_GO 此处显示有关此函数的摘要
% 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
杂波图
点参数
%% 虚警概率_杂波图_点参数
%% 时间:2019.11.27
%% 版本:v1【不在优化】
function [ PFA ] = form_PFA_cm_point( ALPHA,m,r )
%FORM_PFA_CM 此处显示有关此函数的摘要
% ALPHA:门限因子
% m:天线旋转周期
PFA=1;
for n=0:m-1
PFA=PFA.*(1+ALPHA.*r.*(1-r).^n).^-1;
end
end
面参数
%% 虚警概率_杂波图_面参数
%% 时间:2019.11.27
%% 版本:v1【不在优化】
function [ PFA ] = form_PFA_cm_surface( ALPHA,m,r,M )
%FORM_PFA_CM 此处显示有关此函数的摘要
% ALPHA:门限因子
% m:天线旋转周期
% M:参考单元数
PFA=1;
for n=0:m-1
PFA=PFA.*(1+ALPHA.*r.*(1-r).^n./M).^-1;
end
end
二分求解法核心函数
%% 二分法求解方程的解
%% 时间:2019.11.27
%% 版本:v1【不在优化】
function [ d1 ] = func_SOLUTION_binary( d1_scope,d2_precision,func,parameter)
%FUNC_SOLUTION_BINARY 此处显示有关此函数的摘要
% d1: 已知数值
% d2: 目标数值
% d2_scope: 目标初始范围【小值,大值】
% d1_precision: 目标数值精度【精度值】
% func: 公式
% parameter: 相关参数【类型自定】 1=ALPHA, end=d2
shape=size(parameter);
if shape(2)==2 %只有一个参数额外参数需要输入
while 1
d1=mean(d1_scope);
d2=func(d1,parameter(1,1)); % 本质上是PFA_pre
d2_difference=1/d2-1/parameter(1,end);
if abs(d2_difference)<d2_precision || abs(d1_scope(1,1)-d1_scope(1,2))<0.001
return;
elseif d2_difference<0
d1_scope(1,1)=d1;
else
d1_scope(1,2)=d1;
end
end
elseif shape(2)==3 % 用于OS_CFAR门限因子的计算
while 1
d1=mean(d1_scope);
d2=func(d1,parameter(1,1),parameter(1,2)); % 本质上是PFA_pre
d2_difference=1/d2-1/parameter(1,end);
if abs(d2_difference)<d2_precision || abs(d1_scope(1,1)-d1_scope(1,2))<0.001
return;
elseif d2_difference<0
d1_scope(1,1)=d1;
else
d1_scope(1,2)=d1;
end
end
end
end
测试样例及其输出结果
clear;
PFA=10^(-4);
N=36;
Rate=0.5;
[ ALPHA ] = form_ALPHA_ca( PFA,N );
[ PFA1 ] = form_PFA_ca( ALPHA,N );
d1_scope=[0,100];
d2_precision=1;
func1=@form_PFA_ca;
func2=@form_PFA_so;
func3=@form_PFA_go;
func4=@form_PFA_os;
func5=@form_PFA_cm_point;
parameter=[N,PFA];
parameter4=[N,Rate,PFA];
parameter5=[100,1./512,PFA];
[ ALPHA1 ] = func_SOLUTION_binary( d1_scope,d2_precision,func1,parameter);
[ ALPHA2 ] = func_SOLUTION_binary( d1_scope,d2_precision,func2,parameter);
[ ALPHA3 ] = func_SOLUTION_binary( d1_scope,d2_precision,func3,parameter);
ALPHA3_1=form_ALPHA_os (N,Rate,PFA);
[ ALPHA4 ] = func_SOLUTION_binary( d1_scope,d2_precision,func4,parameter4);
[ ALPHA5 ] = func_SOLUTION_binary( d1_scope,d2_precision,func5,parameter5);
输出结果展示
参考
雷达无线电系列(三)经典CFAR算法门限因子alpha计算(matlab)https://www.cnblogs.com/Mufasa/p/11950154.html
万分感谢博主万雨的无私分享。
写于关雎
新浪微博:https://weibo.com/tianzhejia
CSDN博客:https://blog.csdn.net/qq_35605018
博客网站:http://www.zhijiadeboke.xyz
GitHub:https://github.com/ZhijiaTian
QQ邮箱:2461824656@qq.com
126邮箱:tianzhejia@126.com
Outlook邮箱:tianzhejia@outlook.com
以上均可与本人取得联系,欢迎探讨。^ v ^