自信号处理仿真题作业Word格式.docx
- 文档编号:5259962
- 上传时间:2023-05-04
- 格式:DOCX
- 页数:18
- 大小:123.97KB
自信号处理仿真题作业Word格式.docx
《自信号处理仿真题作业Word格式.docx》由会员分享,可在线阅读,更多相关《自信号处理仿真题作业Word格式.docx(18页珍藏版)》请在冰点文库上搜索。
f2=0.17;
f3=0.26;
SNR1=30;
>>SNR2=30;
SNR3=27;
A1=10^(SNR1/20);
A2=10^(SNR2/20);
A3=10^(SNR3/20);
>signal1=A1*exp(j*2*pi*f1*(0:
N-1));
signal2=A2*exp(j*2*pi*f2*(0:
N-1));
signal3=A3*exp(j*2*pi*f3*(0:
N-1));
un=signal1+signal2+signal3+y
基于FFT的自相关函数快速计算方法:
N=32;
Uk=fft(un,2*N);
Sk=(1/N)*abs(Uk).^2;
r0=ifft(Sk);
r1=[r0(N+2:
2*N),r0(1:
N)];
>figure
(1);
>stem(real(r1));
>figure
(2);
stem(imag(r1))
输出结果为:
图1基于FFT的自相关函数快速计算
实部:
虚部:
教材中式(3.1.2)估计自相关函数
r=xcorr(un,N-1,'
biased'
);
figure(1);
stem(real(r))
>figure
(2);
stem(imag(r))
图 2教材式(3.1.2)估计的自相关函数
(2)令信号观测样本长度N=256,试用BT法和周期图法估计
的功率谱,这里设BT法中所用自相关函数的单边长度M=64。
周期图法
N=256;
noise=(randn(1,N)+j*randn(1,N))/sqrt
(2);
f1=0.15;
f2=0.17;
f3=0.26;
SNR1=30;
SNR2=30;
SNR3=27;
A1=10^(SNR1/20);
A2=10^(SNR2/20);
A3=10^(SNR3/20);
signal1=A1*exp(j*2*pi*f1*(0:
signal2=A2*exp(j*2*pi*f2*(0:
signal3=A3*exp(j*2*pi*f3*(0:
N-1));
xn=signal1+signal2+signal3+noise
NF=1024;
Spr=fftshift((1/NF)*abs(fft(xn,NF)).^2);
Spr= Spr/max(Spr);
t=[-0.5:
1/NF:
0.5-1/NF];
plot(t,10*log10(Spr));
图
3 周期图法所得功率谱
BT法
M=64;
r=xcorr(xn,M,'
biased'
);
NF=1024;
BT=fftshift(fft(r, NF));
BT=BT/max(BT);
>>plot((-511:
512)/1024,10*log10(BT))图4 BT法所得功率谱
(3)令信号观测样本长度N=256,试用Levinson-Durbin迭代算法求解AR模型的系数并估计
的功率谱,模型的阶数取为p=16。
p=16;
r0=xcorr(xn,p,'
biased'
r=r0(p+1:
2*p+1);
a(1,1)=-r
(2)/r
(1);
sigma
(1)=r(1)-(abs(r
(2))^2)/r
(1);
form=2:
p
k(m)= -(r(m+1)+sum(a(m-1,1:
m-1).*r(m:
-1:
2)))/sigma(m-1);
a(m,m)=k(m);
fori=1:
m-1
a(m, i)=a(m-1,i) +k(m) *conj(a(m-1,m-i));
end
sigma(m)=sigma(m-1)*(1-abs(k(m))^2);
end
NF=1024;
Par =sigma(p)./fftshift(abs(fft([1,a(p,:
)], NF)).^2);
Par=Par/max(Par);
plot((-511:
512)/1024,10*log10(Par))
图5 16阶AR模型的功率谱估计
3.18设随机过程
为
其中,
是零均值、方差为1的白噪声,
、
是相互独立并在
上服从均匀分布的随机相位。
使用MVDR方法进行信号频率估计的方针实验,画出频率估计谱线,并给出正弦信号频率的估计值。
(要求:
信号样本数取1000,估计的自相关矩阵为8阶)
解:
clearall;
%产生0均值,方差为1的复高斯白噪声序列v(n)
N=1000;
FS=1;
noise=(randn(1,N)+j*randn(1,N))/sqrt(2);
%产生带噪声的信号样本u(n)
signal1=exp(j*0.5*pi*(0:
N-1)/FS+j*2*pi*rand);
signal2=exp(-j*0.3*pi*(0:
N-1)/FS+j*2*pi*rand);
un=signal1+signal2+noise;
%利用un来估计自相关函数
r=zeros(1,N);
for m=1:
N
forn=m:
r(m)=r(m)+un(n)*(un(n-(m-1)))'
;
end
r(m)=r(m)/N;
end
M=8;
%设定横向滤波器的阶数
%利用自相关矩阵和自相关函数的关系,构建自相关矩阵
R=zeros(M,M);
fori=1:
M
forj=1:
ifi<
=j
R(i,j)=r(1+j-i);
%得到了M阶的矩阵R
else
R(i,j)=(r(1+i-j))'
end
end
end
N3=2048;
%设定画图时描点的数目。
d=1/(N3-1);
%求画图用的横坐标的间隔。
h=zeros(1,N3);
fori=1:
N3
h(i)=-0.5+(i-1)*d;
y=zeros(1,N3);
forj=1:
w=h(j)*2*pi;
y(j)=10*log10(abs(getPMVDR(w,M,R)));
end
plot(h,y);
title('MVDR谱估计'
xlabel('
\omega/2\pi');
ylabel('归一化功率谱/dB'
axis([-0.5 0.5-16 1]);
图
6MVDR谱估计
3.19复正弦加白噪声随机过程
同题3.18中所给。
试用MUSIC方法进行信号频率估计的仿真实验。
要求:
信号样本数取1000,估计的自相关矩阵为8阶)
(1)采用AIC和MDL准则估计信号源个数;
noise=(randn(1,N) +j*randn(1, N))/sqrt
(2);
signal1 =exp(j *0.5*pi*(0:
N-1)+j*2*pi*rand);
signal2=exp(-j*0.3*pi *(0:
N-1)+j*2* pi *rand);
un=signal1+signal2+noise;
M=8;
fork=1:
N-M
xs(:
,k)=un(k+M-1:
-1:
k).'
end
R=xs*xs'
/(N-M);
%自相关矩阵的特征值分解
[U,E]=svd(R);
% U是特征矢量组成的矩阵,E是对角阵,对角元素是由大到小排列的特征值
ev=diag(E);
%提取对角元素上的特征值;
%根据AIC准则进行信号源个数的估计
for k=1:
M
dec=prod(ev(k:
M).^(1/(M-k+1)));
% 计算第一项中对数的自变量的分子
nec=mean(ev(k:
M));
% 计算第一项中对数的自变量的分母
lnv=(dec/nec)^((M-k+1)*N);
%计算第一项中对数的自变量
AIC(k)=-2*log(lnv)+2*(k-1)*(2*M-k+1);
%计算AIC准则
end
[Amin, K]=min(AIC);
N1=K-1;
%信号源个数估计
%根据MDL准则进行信号源个数的估计
for k=1:
dec=prod(ev(k:
M).^(1/(M-k+1)));
nec=mean(ev(k:
M));
lnv=(dec/nec)^((M-k+1)*N);
MDL(k)=-log(lnv)+(k-1)/2*(2*M-k+1)*log(N);
%计算DML准则
end
[Amin, K]=min(MDL);
%寻找使DML准则最小的索引
N2=K-1;
%信号源个数的估计
(2)根据
(1)中信源个数画出相应的MUSIC频率估计谱线。
%计算MUSIC谱
En=U(:
N1+1:
M);
%噪声子空间的向量组成的矩阵
NF=2048;
forn=1:
NF
Aq=exp(-j*2*pi*(n-1)/NF*(0:
M-1)');
Pmusic(n)=1/(Aq'
*En*En'
*Aq);
%music谱
Pmusic(n)=10*log10(Pmusic(n));
end
f=-0.5:
1/NF:
0.5-1/NF;
plot(f,Pmusic);
图7MUSIC谱估计
3.20复正弦加白噪声随机过程
同题3.18中所给。
试使用Root-MUSIC方法进行信号频率估计的仿真实验。
信号样本数取1000,估计的自相关矩阵为8阶。
)
(1)计算正弦信号的频率估计值。
N=1000;
noise=(randn(1,N)+j*randn(1,N))/sqrt(2);
signal1= exp(j*0.5* pi*(0:
N-1)+j*2*pi*rand);
signal2 =exp(-j*0.3*pi* (0:
N-1) +j*2*pi* rand);
un=signal1+signal2+noise;
M=8;
fork=1:
N-M
xs(:
,k)=un(k+M-1:
-1:
end
R=xs*xs'
/(N-M);
[U,E]=svd(R);
ev=diag(E);
G=U(:
3:
Gr=G*G'
co=zeros(2*M-1,1);
form=1:
co(m:
m+M-1)=co(m:
m+M-1)+Gr(M:
1,m);
end
z=roots(co);
ph=angle(z)/(2*pi);
err=abs(abs(z)-1);
输出结果为:
最接近单位圆的两个根是:
-0.2707680+1.169i
-0.5924671+0.997740903255518i
上述根的归一化频率为:
0.2563
0.2559
(2)与MUSIC算法的估计结果比较。
En=G;
NF=2048;
forn=1:
Aq=exp(-j*2*pi*(n-1)/NF*(0:
M-1)'
Pmusic(n)=1/(Aq'*En*En'
*Aq);
Pmusic(n)=10*log10(Pmusic(n));
f=-0.5:
1/NF:
0.5-1/NF;
plot(f,Pmusic);
图8 MUSIC谱估计
3.21复正弦加白噪声随机过程
同题3.18中所给。
试使用ESPRIT算法进行信号频率估计的仿真实验,给出正弦信号频率的估计值(要求:
信号样本数取1000,估计的自相关矩阵为8阶。
N=1000;
noise=(randn(1,N) +j*randn(1,N))/sqrt
(2);
signal1=exp(j* 0.5* pi* (0:
N-1)+ j*2*pi*rand);
signal2=exp(-j*0.3*pi*(0:
N-1) +j* 2*pi *rand);
un=signal1+signal2+noise;
M=8;
fork=1:
N-M
xs(:
k)=un(k+M-1:
Rxx=xs(:
,1:
end-1)*xs(:
1:
end-1)'
/(N-M-1);
%计算自相关矩阵Rxx
Rxy=xs(:
end-1)*xs(:
2:
end)'/(N-M-1);
%计算互相关矩阵Rxy
%—————相关矩阵的特征值分解,寻找最小特征值————————%
[U,E]=svd(Rxx);
%矩阵U是特征矢量组成的矩阵,E是对角阵,对角元素是由大到小排列的特征值
ev=diag(E);
%提取对角元素上的特征值;
emin=ev(end);
%获取最小特征值;
%——————构造矩阵对{Cxx,Cxy}————————%
Z=[zeros(M-1,1),eye(M-1);
0,zeros(1,M-1)];
%构造矩阵Z;
Cxx=Rxx-emin*eye(M);
%计算Cxx;
Cxy=Rxy-emin*Z;
%计算Cxy;
%——————矩阵对{Cxx,Cxy}的广义特征值分解————————%
[U,E]=eig(Cxx,Cxy);
%广义特征值分解
z=diag(E);
%提取对角矩阵中的特征值
ph=angle(z)/(2*pi);
%求所有特征值的相位对应的归一化频率
err=abs(abs(z)-1);
%与单位元的距离
单次ESPRIT算法中最接近单位圆的两个特征值为:
0.588516696017445-0.8638i
-0.09186 +0.999262092748743i
上述特征值的相位对应的归一化频率为:
-0.149944811600338
0.2573
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 信号 处理 仿真 作业