现代数字信号处理仿真作业.docx
- 文档编号:512964
- 上传时间:2022-10-10
- 格式:DOCX
- 页数:29
- 大小:243.74KB
现代数字信号处理仿真作业.docx
《现代数字信号处理仿真作业.docx》由会员分享,可在线阅读,更多相关《现代数字信号处理仿真作业.docx(29页珍藏版)》请在冰豆网上搜索。
现代数字信号处理仿真作业
第三章仿真作业
3.17
(1)代码
clear;
N=32;
m=[-N+1:
N-1];
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:
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+noise;
uk=fft(un,2*N);
sk=(1/N)*abs(uk).^2;
r0=ifft(sk);
r1=[r0(N+2:
2*N),r0(1:
N)];
r=xcorr(un,N-1,'biased');
figure
subplot(2,2,1)
stem(m,real(r1));
xlabel('m');
ylabel('FFT估计r1实部');
subplot(2,2,2)
stem(m,imag(r1));
xlabel('m');
ylabel('FFT估计r1虚部');
subplot(2,2,3)
stem(m,real(r));
xlabel('m');
ylabel('平均估计r实部');
subplot(2,2,4)
stem(m,imag(r));
xlabel('m');
ylabel('平均估计r虚部');
仿真结果
(2)代码
clear;
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:
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+noise;
NF=1024;
spr=fftshift((1/NF)*abs(fft(un,NF)).^2);
f1=(0:
length(spr)-1)*(1/(length(spr)-1))-0.5;
M=64;
r=xcorr(un,M,'biased');
bt=fftshift(abs(fft(r,NF)));
f2=(0:
length(bt)-1)*(1/(length(bt)-1))-0.5;
figure
subplot(1,2,1)
plot(f1,10*log10(spr/max(spr)));
xlabel('w/2pi');
仿真结果
(3)
代码
clear;
N=1000;
fai1=rand(1,1)*2*pi;
fai2=rand(1,1)*2*pi;
noise=(randn(1,N)+j*randn(1,N))/sqrt
(2);
un=exp(j*0.5*pi*(0:
N-1)+j*fai1)+exp(-j*0.3*pi*(0:
N-1)+j*fai2)+noise;
p=8;
cx=xcorr(un,p,'biased');
rxx=cx(p+1:
2*p)';
R=toeplitz(rxx);
[u,s]=eig(R);
nw=128;
ww=[-128:
128]/128*pi;
e=exp(-j*ww'*[0:
p-1])%k行m列
ev=e*u(:
1:
p-2);
pw=1./real(diag(ev*ev'));
plot(ww/(2*pi),10*log10(pw)/max(pw));
仿真结果
3.20
(1)代码
clear;
N=1000;
fai1=rand(1,1)*2*pi;
fai2=rand(1,1)*2*pi;
noise=(randn(1,N)+j*randn(1,N))/sqrt
(2);
un=exp(j*0.5*pi*(0:
N-1)+j*fai1)+exp(-j*0.3*pi*(0:
N-1)+j*fai2)+noise;
p=8;
cx=xcorr(un,p,'biased');
rxx=cx(p+1:
2*p)';
R=toeplitz(rxx);
[u,s]=eig(R);
nw=128;
ww=[-128:
128]/128*pi;
e=exp(-j*ww'*[0:
p-1])%k行m列
ev=e*u(:
1:
p-2);
pw=1./real(diag(ev*ev'));
plot(ww/(2*pi),10*log10(pw)/max(pw));仿真结果
距离单位圆最近的两个解为-0.2363-0.9717i和0.3747+0.9271i,对应的归一化频率为0.1889和-0.2880
(2)代码
clear;
N=1000;
fai1=rand(1,1)*2*pi;
fai2=rand(1,1)*2*pi;
noise=(randn(1,N)+j*randn(1,N))/sqrt
(2);
un=exp(j*0.5*pi*(0:
N-1)+j*fai1)+exp(-j*0.3*pi*(0:
N-1)+j*fai2)+noise;
p=8;
cx=xcorr(un,p,'biased');
rxx=cx(p+1:
2*p)';
R=toeplitz(rxx);
[u,s]=eig(R);
nw=128;
ww=[-128:
128]/128*pi;
e=exp(-j*ww'*[0:
p-1])%k行m列
ev=e*u(:
1:
p-2);
pw=1./real(diag(ev*ev'));
plot(ww/(2*pi),10*log10(pw)/max(pw));
仿真结果
3.21
代码
clear;
N=1000;
fai1=rand(1,1)*2*pi;
fai2=rand(1,1)*2*pi;
noise=(randn(1,N)+j*randn(1,N))/sqrt
(2);
un=exp(j*0.5*pi*(0:
N-1)+j*fai1)+exp(-j*0.3*pi*(0:
N-1)+j*fai2)+noise;
p=8;
fork=1:
N-p
xs(:
k)=un(k+p-1:
-1:
k)';
end
rxx=xs(:
1:
end-1)*xs(:
1:
end-1)'/(N-p-1);
rxy=xs(:
1:
end-1)*xs(:
2:
end)'/(N-p-1);
[u,e]=svd(rxx);
ev=diag(e);
emin=ev(end);
z=[zeros(p-1,1),eye(p-1);0,zeros(1,p-1)];
cxx=rxx-emin*eye(p);
cxy=rxy-emin*z;
[U,E]=eig(cxx,cxy);
Z=diag(E);
ph=angle(Z)/(2*pi);
err=abs(abs(Z)-1);
仿真结果
最接近单位圆的两个解分别为0.5867+0.8097i和0.0349-0.9984i,对应的归一化频率为0.1502和-0.2444。
第四章仿真题作业
4-18
(1)产生随机序列
代码
clear;
data_len=512;
trials=100;
n=1:
data_len;
a1=-0.975;
a2=0.95;
sigma_v_2=0.0731;
vn=sqrt(sigma_v_2)*randn(data_len,1,trials);
u0=[00];
num=1;
den=[1a1a2];
zi=filtic(num,den,u0);
u=filter(num,den,vn,zi);
(2)单次实验估计最优权向量
代码
clear;
data_len=512;
trials=100;
n=1:
data_len;
a1=-0.975;
a2=0.95;
sigma_v_2=0.0731;
vn=sqrt(sigma_v_2)*randn(data_len,1,trials);
u0=[00];
num=1;
m=1;
den=[1a1a2];
zi=filtic(num,den,u0);
u=filter(num,den,vn,zi);
mu1=0.05;
mu2=0.005;
w1=zeros(2,data_len);
w2=zeros(2,data_len);
e1=zeros(data_len,1);
e2=zeros(data_len,1);
d1=zeros(data_len,1);
d2=zeros(data_len,1);
%form=1:
trials
forn=3:
data_len-1
w1(:
n+1)=w1(:
n)+mu1*u(n-1:
-1:
n-2,:
m)*conj(e1(n));
w2(:
n+1)=w2(:
n)+mu2*u(n-1:
-1:
n-2,:
m)*conj(e2(n));
d1(n+1)=w1(:
n+1)'*u(n:
-1:
n-1,:
m);
d2(n+1)=w2(:
n+1)'*u(n:
-1:
n-1,:
m);
e1(n+1)=u(n+1,:
m)-d1(n+1);
e2(n+1)=u(n+1,:
m)-d2(n+1);
end
figure
(1)
holdon;
plot(w1(1,:
));
plot(w1(2,:
));
plot(w2(1,:
),’g’);
plot(w2(2,:
),’g’)
仿真结果
(3)100次实验计算剩余均方误差
代码
clear;
data_len=512;
trials=100;
n=1:
data_len;
a1=-0.975;
a2=0.95;
sigma_v_2=0.0731;
vn=sqrt(sigma_v_2)*randn(data_len,1,trials);
u0=[00];
num=1;
%m=1;
den=[1a1a2];
zi=filtic(num,den,u0);
u=filter(num,den,vn,zi);
mu1=0.05;
mu2=0.005;
w1=zeros(2,data_len);
w2=zeros(2,data_len);
e1=zeros(data_len,1);
e2=zeros(data_len,1);
d1=zeros(data_len,1);
d2=zeros(data_len,1);
jn1=zeros(1,trials);
jn2=zeros(1,trials);
mse1=zeros(1,trials);
mse2=zeros(1,trials);
form=1:
trials
forn=3:
data_len-1
w1(:
n+1)=w1(:
n)+mu1*u(n-1:
-1:
n-2,:
m)*conj(e1(n));
w2(:
n+1)=w2(:
n)+mu2*u(n-1:
-1:
n-2,:
m)*
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 现代 数字信号 处理 仿真 作业