数字信号处理胡广书例题作业程序.docx
- 文档编号:12100900
- 上传时间:2023-04-17
- 格式:DOCX
- 页数:70
- 大小:621.37KB
数字信号处理胡广书例题作业程序.docx
《数字信号处理胡广书例题作业程序.docx》由会员分享,可在线阅读,更多相关《数字信号处理胡广书例题作业程序.docx(70页珍藏版)》请在冰豆网上搜索。
数字信号处理胡广书例题作业程序
1、
%---filter求卷积,B(Z)/A(Z)=H(Z),已知B(Z)和A(Z),求y(n)=x(n)*h(n)-----
clear;
x=ones(100);
t=1:
100;
b=[.001836,.007344,.011016,.007374,.001836];
a=[1,-3.0544,3.8291,-2.2925,.55075];
%
y=filter(b,a,x);
%求所给系统的输出,本例实际上是求所给系统的阶跃响应;
plot(t,x,'r.',t,y,'k-');gridon;
ylabel('x(n)andy(n)')
xlabel('n')
1、
%---filter求卷积,B(Z)/A(Z)=H(Z),已知B(Z)和A(Z),求y(n)=x(n)*h(n)-----
clear;
x=ones(100);
t=1:
100;
b=[.001836,.007344,.011016,.007374,.001836];
a=[1,-3.0544,3.8291,-2.2925,.55075];
%
y=filter(b,a,x);
%求所给系统的输出,本例实际上是求所给系统的阶跃响应;
plot(t,x,'r.',t,y,'k-');gridon;
ylabel('x(n)andy(n)')
xlabel('n')
第一章产生信号,求卷积和自相关函数
1、
%信号产生
n=0:
100;
%工频
f0=50;A=220;fs=400;
x1=A*sin(2*pi*f0*n/fs);
subplot(321);plot(n,x1);xlabel('n');ylabel('x1(n)');gridon;
%率减正弦
f0=2;A=2;alf=0.5;fs=16;
x2=A*exp(-alf*n/fs).*sin(2*pi*f0*n/fs);
subplot(323);plot(n,x2);xlabel('n');ylabel('x2(n)');gridon;
%谐波信号
f0=5;A1=1.0;A2=0.5;A3=0.2;fs=100;
x3=A1*sin(2*pi*f0*n/fs)+A2*sin(2*pi*2*f0*n/fs)+A3*sin(2*pi*3*f0*n/fs);
subplot(322);plot(n,x3);xlabel('n');ylabel('x3(n)');gridon;
%哈明窗
f0=10;fs=1000;
x4=0.54-0.46*cos(2*pi*f0*n/fs);
subplot(324);plot(n,x4);xlabel('n');ylabel('x4(n)');gridon;
%采样
n=-50:
50;
f0=10;fs=400;
w=2*pi*f0*n/fs;
x5=sinc(w);
subplot(325);plot(n,x5);xlabel('n');ylabel('x5(n)');gridon;
2、
%产生均匀分布的白噪信号,使均值为0,功率为p
%-----------------------------------------------------------------
clear;
p=0.01;
N=50000;
u=rand(1,N);
u=u-mean(u);
a=sqrt(12*p);
u1=u*a;
power_u1=dot(u1,u1)/N
subplot(211)
plot(u1(1:
200));gridon;
ylabel('u(n)')
xlabel('n')
3、
%产生高斯分布的白噪信号,使功率为p,并观察数据分布的直方图
%-----------------------------------------------------------------
clear;
p=0.1;
N=500000;
u=randn(1,N);
a=sqrt(p);
u=u*a;
power_u=var(u);
subplot(211)
plot(u(1:
200));gridon;
ylabel('u(n)');
xlabel('n')
subplot(212)
hist(u,50);gridon;
ylabel('histogramofu(n)');
4、
%产生一sinc函数;
%-----------------------------------------------------------------
clear;
n=200;
stept=4*pi/n;
t=-2*pi:
stept:
2*pi;
y=sinc(t);
plot(t,y,t,zeros(size(t)));
ylabel('sinc(t)');
xlabel('t=-2*pi~2*pi');
gridon;
5、
%产生一chirp信号;
%chirp(T0,F0,T1,F1):
%T0:
信号的开始时间;F0:
信号在T0时的瞬时频率,单位为Hz;
%T1:
信号的结束时间;F1:
信号在T1时的瞬时频率,单位为Hz;
%-----------------------------------------------------------------
clear;
t=0:
0.001:
1;
x=chirp(t,0,1,125);
plot(t,x);
ylabel('x(t)')
xlabel('t')
6、
%计算两个序列的线性卷积;
%-----------------------------------------------------------------
clear;
N=5;%第一个序列的长度
M=6;%第二个序列的长度
L=N+M-1;
x=[1,2,3,4,5];
h=[6,2,3,6,4,2];
y=conv(x,h);
nx=0:
N-1;
nh=0:
M-1;
ny=0:
L-1;
subplot(231);
%绘制x
stem(nx,x,'.k');xlabel('n');ylabel('x(n)');gridon;
subplot(232);
%绘制h
stem(nh,h,'.k');xlabel('n');ylabel('h(n)');gridon;
subplot(233);
%绘制卷积
stem(ny,y,'.k');xlabel('n');ylabel('y(n)');gridon;
7、
%求两个序列的互相关函数,或一个序列的自相关函数;
%-----------------------------------------------------------------
clear;
N=500;
p1=1;
p2=0.1;
f=1/8;
Mlag=50;%自相关的单边长度
u=randn(1,N);
n=[0:
N-1];
s=sin(2*pi*f*n);
%混有高斯白噪的正弦信号的自相关
u1=u*sqrt(p1);%高斯白噪声
x1=u1(1:
N)+s;%混合信号
rx1=xcorr(x1,Mlag,'biased');%自相关,无偏估计
subplot(221);
plot(x1(1:
Mlag));title('信号x1');
xlabel('n');
ylabel('x1(n)');gridon;
subplot(223);
plot((-Mlag:
Mlag),rx1);title('x1自相关');gridon;
xlabel('m');ylabel('rx1(m)');
%高斯白噪功率由原来的p1减少为p2,再观察混合信号的自相关
u2=u*sqrt(p2);%改变高斯白噪声
x2=u2(1:
N)+s;%新的混合信号
rx2=xcorr(x2,Mlag,'biased');
subplot(222);
plot(x2(1:
Mlag));title('信号x2');
xlabel('n');ylabel('x2(n)');gridon;
subplot(224);
plot((-Mlag:
Mlag),rx2);title('x2自相关');
gridon;xlabel('m');ylabel('rx2(m)');
8、
%求序列的自相关函数
clear
N=500;
Mlag=50;%单边长度
nx=0:
N-1;
x=exp(-nx*0.1);
rx=xcorr(x,Mlag,'biased');
nrx=-Mlag:
Mlag;%自相关序列的程度
subplot(211);plot(nx,x);xlabel('n');ylabel('x(n)');gridon;
subplot(212);plot(nrx,rx);xlabel('n');ylabel('rx(n)');gridon;
9、
%正弦加白噪声,自相关
p=0.1;N=5000;Mlag=100;
u=rand(1,N);u=u-mean(u);
a=sqrt(12*p);u=u*a;
power_u=dot(u,u)/N
nx=1:
1000;
x=1.414*sin(nx*pi/16.0);
x1=x(1:
1000)+5*u(1:
1000);
rx=xcorr(x1,Mlag,'biased');
nrx=-Mlag:
Mlag;
subplot(211);plot(x1(1:
200));xlabel('n');ylabel('x(n)');gridon;
subplot(212);plot(nrx,rx);xlabel('n');ylabel('rx(n)');gridon;
1、
%产生信号,求卷积,FFT,求平均
clearall;
N=1024;%采样点数
fs=100.0;%采样频率
alf1=-1.0;f1=5;
alf2=-1.5;f2=8;
alf3=-0.7;f3=10;
%产生x和w两个信号
%产生x
u=rand(1,N);
u=u-mean(u);%均值为0的白噪声
t=[0:
1/fs:
(N-1)/fs];
x=1.0*exp(alf1*t).*sin(2*pi*f1*t)+1.0*exp(alf2*t).*sin(2*pi*f2*t)+1.0*exp(alf3*t).*sin(2*pi*f3*t);
x=x+u;
x=x/max(x);
%产生w
alf4=-1.0;
w=1.0*exp(alf4*t);
%x=x-mean(x);
figure
(1);
subplot(211);
plot(t,x,t,w,'r');title('x(t)w(t)')
%应用FFT求频谱;
f=0:
fs/N:
fs/N*(N-1);
X=fft(x,N);
X=abs(X);
%X=20*log10(X);
subplot(212);
plot(f(1:
N/2),X(1:
N/2));title('x(t)频谱')
%stem(f,X,'.');gridon;
xlabel('Hz')
%求卷积-----------------------------------------------
y=x.*w;%时域相乘,频域卷积
figure
(2);
subplot(211);
plot(t,y);title('正弦加白噪声后与w时域相乘')
Y=fft(y,N);
Y=abs(Y);
subplot(212);
plot(f(1:
N/2),Y(1:
N/2));title('正弦加白噪声后与w时域相乘的FFT')
xlabel('Hz')
xlabel('Hz')
%平均1000次----------------------------------------------
figure(3);
Y=zeros(1,N);
u=rand(1,1000*N);u=u-mean(u);%零均值白噪声
fori=0:
999
x=1.0*exp(alf1*t).*sin(2*pi*f1*t)+1.0*exp(alf2*t).*sin(2*pi*f2*t)+1.0*exp(alf3*t).*sin(2*pi*f3*t);
x=x+10*u(1+i*N:
i*N+N);x=x/max(x);
X=fft(x,N);
X=abs(X);
Y=Y+X;
end
Y=Y/1000;
subplot(211);plot(f(1:
N/2),Y(1:
N/2));title('正弦加白噪声的FFT1000次平均')
xlabel('Hz')
Y=zeros(1,N);
fori=0:
999
x=1.0*exp(alf1*t).*sin(2*pi*f1*t)+1.0*exp(alf2*t).*sin(2*pi*f2*t)+1.0*exp(alf3*t).*sin(2*pi*f3*t);
x=x+10*u(1+i*N:
i*N+N);x=x/max(x);
x=x.*w;
X=fft(x,N);
X=abs(X);
Y=Y+X;
end
Y=Y/1000;
subplot(212);plot(f(1:
N/2),Y(1:
N/2));title('正弦加白噪声后与w时域相乘的FFT1000次平均')
xlabel('Hz')
2、
clearall;
%三个正弦信号相加,分段函数,进行频谱分析
%产生三个正弦相加的函数;
N=512;
f0=10;fs=100.0;
t=[0:
N-1];
x=1.0*sin(2*pi*f0*t/fs)+1.0*sin(2*pi*2*f0*t/fs)+1.0*sin(2*pi*3*f0*t/fs);
subplot(211);
plot(t(1:
N),x(1:
N));title('x(t)')
%加窗
w=1-1.93*cos(2*pi*t/N)+1.29*cos(4*pi*t/N)-0.388*cos(6*pi*t/N)+0.0322*cos(8*pi*t/N);
%w=1.0-cos(2*pi*t/N);
x=w.*x;%加窗等于时域点乘
%应用FFT求频谱;
f=0:
fs/N:
fs/N*(N-1);
X=fft(x,N);%先点乘再进行傅里叶变换
X=abs(X)/N;
subplot(212);
plot(f(1:
N/2),X(1:
N/2));title('x(t)加窗之后的傅里叶变换')
xlabel('Hz')
%分段函数-----------------------------------------------
M=170;L=N-2*M;
x(1:
M)=1.0*sin(2*pi*f0*t(1:
M)/fs);
x(M+1:
2*M)=1.0*sin(2*pi*2*f0*t(1:
M)/fs);
x(2*M+1:
N)=1.0*sin(2*pi*3*f0*t(1:
L)/fs);
figure;
subplot(211);
plot(t(1:
N),x(1:
N));title('x(t)为分段函数')
%w=1-1.93*cos(2*pi*t/M)+1.29*cos(4*pi*t/M)-0.388*cos(6*pi*t/M)+0.0322*cos(8*pi*t/M);
%x(1:
M)=w(1:
M).*x(1:
M);
X=fft(x,N);
X=abs(X)/N;
subplot(212);
plot(f(1:
N/2),X(1:
N/2));title('x(t)频谱分析')
xlabel('Hz')
3、
%采样长度不同对FFT的影响----------------------------------------------------
clearall;
%观察数据长度N的变化对DTFT分辨率的影响
f1=2;f2=2.02;f3=2.07;fs=10;
w=2*pi/fs;
N=256;%N=256
x=sin(w*f1*(0:
N-1))+sin(w*f2*(0:
N-1))+sin(w*f3*(0:
N-1));
f=0:
fs/N:
fs/2-1/N;
X=fft(x);
X=abs(X);
subplot(221);
plot(f(45:
60),X(45:
60));title('N=256');gridon;
subplot(223)
plot(f(1:
N/2),X(1:
N/2));title('N=256');gridon;
xlabel('Hz')
%
N=N*4;%N=1024
x=sin(w*f1*(0:
N-1))+sin(w*f2*(0:
N-1))+sin(w*f3*(0:
N-1));
f=0:
fs/N:
fs/2-1/N;
X=fft(x);
X=abs(X);
subplot(222)
plot(f(45*4:
4*60),X(4*45:
4*60));title('N=1024');gridon;
xlabel('Hz')
4、
%补零的影响-------------------------------------------------------------------
clear;
%计算长度为N的原始信号的DTFT
f1=2.67;f2=3.75;f3=6.75;fs=20;w=2*pi/fs;
N=16;
x=sin(w*f1*(0:
N-1))+sin(w*f2*(0:
N-1)+pi/2)+sin(w*f3*(0:
N-1));
f=0:
fs/N:
fs/2-1/N;
X=fft(x);
X=abs(X);
f=fs/N*(0:
N/2-1);
subplot(221)
stem(f,X(1:
N/2),'.');title('不补零');gridon;
xlabel('Hz')
%在数据末补N个零
x(N:
2*N-1)=0;
X=fft(x);X=abs(X);
f=fs*(0:
N-1)/(2*N);
subplot(222)
stem(f,X(1:
N),'.');title('补N个零');gridon;
xlabel('Hz')
%在数据末补7*N个零
x(N:
8*N-1)=0;
X=fft(x);X=abs(X);
f=fs*(0:
4*N-1)/(8*N);
subplot(223)
stem(f,X(1:
4*N),'.');title('补7N个零');gridon;
xlabel('Hz')
%在数据末补29*N个零
x(N:
30*N-1)=0;
X=fft(x);X=abs(X);
f=fs*(0:
15*N-1)/(30*N);
subplot(224)
plot(f,X(1:
15*N));title('补29N个零');gridon;
xlabel('Hz')
5、
%x(n)是两个正弦信号和一个白噪声相加,FFT和IFFT------------------
-----------
clearall;
%产生两个正弦加白噪声;
N=256;
f1=.1;f2=.2;fs=1;
a1=5;a2=3;
w=2*pi/fs;
x=a1*sin(w*f1*(0:
N-1))+a2*sin(w*f2*(0:
N-1))+randn(1,N);
%应用FFT求频谱;
subplot(3,1,1);
plot(x(1:
N/4));title('x(n)')
f=-0.5:
1/N:
0.5-1/N;
X=fft(x);
y=ifft(X);
subplot(3,1,2);
plot(f,fftshift(abs(X)));title('fft(x)')
subplot(3,1,3);
plot(real(x(1:
N/4)));title('x(n)实部')
6、
%fftfilt和conv比较
n=0:
2;
h=1./2.^n;
fori=1:
51
x(i)=(i-1)/5;
end
fori=52:
100
x(i)=20-(i-1)/5;
end
y=fftfilt(h,x);%y=x*h
z=conv(h,x);
subplot(311)
title('x(t)')
hold;
plot(x);
subplot(312)
title('fftfilt叠接相加法')
hold;
plot(y);
subplot(313)
plot(z);
axis([0,100,0,20]);
title('conv卷积')
7、
clear;
%补零后对频谱的影响
%计算长度为N的原始信号的DFT
f0=50;%信号频率
fs=200;%采样频率
N=16;%抽样点数
w=2*pi/fs;%数字角频率
x=sin(w*f0*(0:
N-1));
X=fft(x);
X=abs(X);
f=fs/N*(0:
N/2-1);%N/2的数据
subplot(211)
stem(f,X(1:
N/2),'.');title('未补零');gridon;
xlabel('Hz')
%在数据末补N个零
x(N:
2*N-1)=0;
X=fft(x);X=abs(X);
f=fs*(0:
N-1)/(2*N);
subplot(212)
stem(f,X(1:
N),'.');title('补N个零');gridon;
xlabel('Hz')
8、
%抽样频率不同的影响
%计算长度为N的原始信号的DFT
%fs=100
f0=50;fs=100;w=2*pi/fs;
N=16;
x=sin(w*f0*(0:
N-1));
X=fft(x);
X=abs(X);
f=fs/N*(0:
N/2-1);
subplot(221)
stem(f,X(1:
N/2),'.');gridon;
xlabel('Hz')
%fs=150
f0=50;fs=150;w=2*pi/fs;
N=16;
x=sin(w*f0*(0:
N-1));
X=fft(x);
X=abs(X);
f=fs/N*(0:
N/2-1);
subplot(222)
stem(f,X(1:
N/2),'.');gridon;
xlabel('Hz')
%fs=200
f0=50;fs=200;w=2*pi/fs;
N=16;
x=sin(w*
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数字信号 处理 胡广书 例题 作业 程序