基于matlab实现随机性信号检测（能量检测器）

% 已知参数的随机性信号检测 （能量检测器）% 经典的假设检验问题，通常假设实信号：% H0：只有高斯白噪声，均值为0，方差σw^2% H1：高斯随机信号+白噪声，均值为0，方差σs^2+σw^2% 这时候用似然比准则可以得到检测器为能量检测器% 但是阈值不知道。本程序中给出了检测阈值。% 另外，对于复信号理论上并没有给出推导，虽然也是能量% 检测器，但是实际上和实信号还是有一定的区别的，% 尤其在阈值上。本程序也给出了复信号的检测阈值，% 并给出了较为简洁的分析方法。% 在理论分析的基础上，利用蒙特卡洛仿真验证了结果。% 给出的PDF文件中，有从书上截图，可以对比结果。clear all; close all; clc;N=25;   % 信号长度Pfa=1e-5;   % 虚警概率p = 1-Pfa;   % 卡方分布的累积概率密度为1-Pfa% 注意：实数据和复数据的门限是不一样的% 原因在于，实数据是N个x^2相加，每个x~N(0,σw^2)% Tx=sum(x^2)，那么Tx/σw^2就服从于N自由度的卡方分布% 而复数据，虽然也是N个数据的平方和相加，但如果把实部和虚部拆开来看% 就变成了2N个数据的平方和相加，每个数据x~N(0,σw^2/2)，% Tx=sum(x^2)，那么Tx/（σw^2/2）就服从于2N自由度的卡方分布% 所以应该按照下面的两个公式分别求解实数和复数的检测门限sigma_w=1;   % 高斯白噪声的标准差 σwthres_shi_theory = chi2inv(p,N)*sigma_w^2;   %根据逆卡方分布的累积函数，求解实序列的检测阈值thres_com_theory = chi2inv(p,2*N)*sigma_w^2/2;   %根据逆卡方分布的累积函数，求解复序列的检测阈值SNR_dB=-20:10;  % 信噪比 σs^2/σw^2，单位dBSNR_len=length(SNR_dB);  SNR=10.^(SNR_dB/10);     % σs^2sigma_s=SNR.^0.5;   % 标准差σssigma_all=sqrt(sigma_s.^2+sigma_w.^2);   % H1假设下的标准差，因为总方差是σs^2+σw^2Mont_num=10000;  % 蒙特卡洛仿真次数for ii=1:SNR_len    sigma_s=SNR(ii)^0.5;   % 信号分量的标准差σs    num_detect_shi=0;   % 实信号检测时的计数器，初始化    num_detect_complex=0;    for jj=1:Mont_num        %% 复信号        x_real=sigma_all(ii)/(2^0.5)*randn(1,N);    % 复信号的实部        x_imag=sigma_all(ii)/(2^0.5)*randn(1,N);   % 复信号的虚部        x_complex=x_real+1i*x_imag;  % 复信号        Tx_complex=x_complex*x_complex';   % 计算出来的Tx        if Tx_complex>thres_com_theory  % 复信号检测            num_detect_shi=num_detect_shi+1;   % 计数器加1        end        %% 实信号        x_shi=sigma_all(ii)*randn(1,N);   % 实信号        Tx_shi=x_shi*x_shi';    % 计算出来的Tx        if Tx_shi>thres_shi_theory            num_detect_complex=num_detect_complex+1;        end    end    detect_rate_complex(ii)=num_detect_shi/Mont_num;  % 复信号仿真的检测概率    detect_rate_shi(ii)=num_detect_complex/Mont_num;    % 实信号仿真的检测概率endPD_complex= 1-chi2cdf(2*thres_com_theory/sigma_w.^2./(SNR+1),2*N);   %复信号的理论检测概率PD_shi= 1-chi2cdf(thres_shi_theory/sigma_w.^2./(SNR+1),N);   %实信号的理论检测概率figureplot(SNR_dB,PD_shi,'r-','Linewidth',2); hold on;plot(SNR_dB,detect_rate_shi,'ro');hold on;plot(SNR_dB,PD_complex,'b-','Linewidth',2); hold on;plot(SNR_dB,detect_rate_complex,'bo');grid on;legend('理论值(实)','仿真值(实)','理论值(复)','仿真值(复)')axis([SNR_dB(1) SNR_dB(end) 0 1]);