%%   對數(shù)正態(tài)分布的雜波序列
xxi=conv(b,xi);
clear all;
azi_num=2000;   %雷達(dá)回波幀數(shù)
fr=1e3;        %脈沖重復(fù)頻率
lamda0=1/30;    %波長
sigmav=1.0;
sigmaf=2*sigmav/lamda0;
rand('state',sum(100*clock));%產(chǎn)生服從u（0，1）的隨機(jī)分布
d1=rand(1,azi_num);
rand('state',7*sum(100*clock)+3);
d2=rand(1,azi_num);
xi=2*sqrt(-2*log(d1)).*cos(2*pi*d2);%正交獨(dú)立的高斯序列 服從N（0，1）
xq=2*sqrt(-2*log(d1)).*sin(2*pi*d2);
coe_num=12;    %求濾波器系數(shù)，用傅立葉級(jí)數(shù)展開法
for n=0:coe_num
    coeff(n+1)=2*sigmaf*sqrt(pi)*exp(-4*sigmaf^2*pi^2*n^2/fr^2)/fr;
end
for n=1:2*coe_num+1
    if n<=coe_num+1
        b(n)=1/2*coeff(coe_num+2-n);
    else
        b(n)=1/2*coeff(n-coe_num);
    end
end
xxi=xxi(coe_num*2+1:azi_num+coe_num*2);%目的是去掉暫態(tài)響應(yīng)
xisigmac=std(xxi);
ximuc=mean(xxi);
yyi=(xxi-ximuc)/xisigmac;
muc=10;              %中位值
sigmac=0.6;          %雜波的標(biāo)準(zhǔn)差
yyi=sigmac*yyi+log(muc);      %使雜波具有指定的標(biāo)準(zhǔn)差
xdata=exp(yyi); %對數(shù)正態(tài)分布的雜波序列
figure
plot(xdata);title('對數(shù)正態(tài)分布雜波時(shí)域波形');
%求概率密度函數(shù)的參數(shù)
num=100;
maxdat=max(abs(xdata));
mindat=min(abs(xdata));
NN=hist(abs(xdata),num);
xpdf1=num*NN/((sum(NN))*(maxdat-mindat));%用直方圖估計(jì)概率密度函數(shù)
xaxis1=mindat:(maxdat-mindat)/num:maxdat-(maxdat-mindat)/num;
th_val=lognpdf(xaxis1,log(muc),sigmac);
% th_val=(xaxis1./sigmac.^2).*exp(-xaxis1.^2./(2*sigmac.^2));%概率密度函數(shù)理論值
figure
plot(xaxis1,xpdf1);hold;%做出仿真結(jié)果的概率密度函數(shù)曲線
plot(xaxis1,th_val,'r');grid on;%做出理論噶率密度函數(shù)曲線
title('雜波幅度分布');xlabel('幅度');ylabel('概率密度');
signal=xdata;
signal=signal-mean(signal); %求功率譜密度，先去掉直流分量
%用Burg法估計(jì)功率譜密度
M=128;
psd_dat=pburg(real(signal),16,M,fr);
psd_dat=psd_dat/(max(psd_dat));      %歸一化處理
freqx=0:0.5*M;
freqx=freqx*fr/M;
powerf=exp(-freqx.^2/(2*sigmaf.^2));%理想高斯曲線
figure
plot(freqx,psd_dat);title('雜波頻譜');xlabel('頻率');ylabel('功率譜密度');hold;
plot(freqx,powerf,'r');grid on;