数字信号处理编程作业机电.docx
- 文档编号:10034141
- 上传时间:2023-02-08
- 格式:DOCX
- 页数:33
- 大小:629.80KB
数字信号处理编程作业机电.docx
《数字信号处理编程作业机电.docx》由会员分享,可在线阅读,更多相关《数字信号处理编程作业机电.docx(33页珍藏版)》请在冰豆网上搜索。
数字信号处理编程作业机电
数字信号处理编程作业
姓名:
学号:
单位:
专业:
导师:
实验一离散时间系统及离散卷积
1、离散时间系统的单位脉冲响应
(1)选择一个离散时间系统;
y(n)=x(n)+2y(n-1)
y(n)=0,n<0;h(n)=0,n<0;
(2)用笔进行差分方程的递推计算;
当x(n)=
,y(n)=h(n),得
h(0)=1+2y(-1)=1;
h
(1)=0+2h(0)=2;
h
(2)=0+2h
(1)=4;
h(n)=0+2h(n-1)=2n.
从而,系统函数的单位冲击响应为
h(n)=2nu(n);
(3)编制差分方程的递推计算程序;
n=[0:
9];
y=zeros(1,length(n));
a=2;
form=1:
10
if(m==1)
y(1,1)=1;
else
y(1,m)=a*y(1,m-1);
end
end
subplot(1,1,1),stem(n,y);title('³å¼¤ÏìÓ¦');xlabel('n');ylabel('h(n)');
(4)在计算机上实现递推运算;
得到的结果,如图所示
(5)将程序计算结果与笔算的计算结果进行比较,验证程序运行的正确性;
>>y
y=
1248163264128256512
与笔算结果完全符合。
1-2参考实例验证:
离散时间系统的单位脉冲响应
%单位取样序列函数
function[x,n]=impseq(n0,n1,n2)
n=[n1:
n2];
x=[(n-n0)==0];
%单位阶跃序列函数
function[x,n]=stepseq(n0,n1,n2)
n=[n1:
n2];
x=[(n-n0)>=0];
(1)参考实例验证
>>[x,n]=impseq(0,-20,120);
>>a=[1,-1,0.9];b=1;
>>h=filter(b,a,x);
>>subplot(2,1,1),stem(n,h);title('冲激响应');xlabel('n');ylabel('h(n)');
subplot(2,1,2),zplane(b,a);title('系统零极点图');
2、离散系统的幅频、相频的分析方法
(1)参考实例验证
b=[0.0181,0.0543,0.0543,0.0181];
a=[1.000,-1.76,1.1829,-0.2781];
functionwork1(h,N)
sum=0;
w=0:
0.1:
pi;
forn=1:
N
sum=sum+h(n).*exp(-j*n*w);
end
Hm=abs(sum);
Ha=angle(sum);
subplot(2,1,1);plot(w/pi,Hm);title('幅度响应');grid;
ylabel('幅度');xlabel('以\pi为单位的频率');
subplot(2,1,2);plot(w/pi,Ha);title('相位响应');grid;
ylabel('相位/pi');xlabel('以\pi为单位的频率');
>>b=[0.0181,0.0543,0.0543,0.0181];a=[1.000,-1.76,1.1829,-0.2781];
[x,n]=impseq(0,0,100);h=filter(b,a,x);
>>work1(h);
3、离散卷积的计算
(1)参考实例验证
%矩形序列函数
function[x,n]=rectseq(n0,n1,n2,n3)
n=[n2:
n3];
x=[(n-n0>=0)&(n-n1<0)];
>>[x,n0]=rectseq(0,10,-5,50);[y,n1]=stepseq(0,-5,50);h=y.*(0.9.^n1);
c=conv(x,h);
>>subplot(3,1,1);stem(n0,x);title('输入序列');
xlabel(')');ylabel('x(n)');axis([-5,50,0,2]);
subplot(3,1,2);stem(n1,h);title('冲激响应');
xlabel('n');ylabel('h(n)');axis([-5,50,0,2]);
subplot(3,1,3);stem(-10:
100,c);title('输出');
xlabel('n');ylabel('y(n)');axis([-10,100,0,8]);
实验二离散傅立叶变换与快速傅立叶变换
一、实验内容
1、用离散傅立叶变换程序处理时间抽样信号,并根据实序列离散傅立叶变换的对称性,初步判定程序的正确性;
function[Xk]=dft(xn,N)
n=[0:
1:
N-1];
k=[0:
1:
N-1];
WN=exp(-j*2*pi/N);
nk=n'*k;
WNnk=WN.^nk;
Xk=xn*WNnk;
>>F=50;n=1:
64;T=0.000625;x=cos(2*pi*F*n*T);X=dft(x,64);
>>subplot(2,1,1);stem(n,x);xlabel('n');ylabel('x(n)');
subplot(2,1,2);stem(n,X);xlabel('k');ylabel('X(k)');
2、观察三角波和反三角波序列的时域和幅频特性,用N=8点FFT分析信号序列
和
的幅频特性,观察两者的序列形状和频谱曲线有什么异同?
绘出两序列及其幅频特性曲线。
三角波序列
反三角波序列
functionxc=xc()
%三角波序列
forn=0:
7
i=n+1;
if(n>=0&&n<=3)
xc(i)=n+1;
elseif(n>=4)
xc(i)=8-n;
else
xc(i)=0;
end
end
end
N=8;
XC=fft(xc,N);
n=0:
7;
k=0:
7;
subplot(2,2,1);stem(n,xc);title('序列xc');
xlabel('n');ylabel('xc(n)');
subplot(2,2,2);stem(k,XC);title('FFT变换');
xlabel('k');ylabel('XC(k)');
subplot(2,2,3);stem(k,abs(XC));title('幅度');
xlabel('k');ylabel('|XC(k)|');
subplot(2,2,4);stem(k,angle(XC));title('相位');
xlabel('k');
>>xc;
functionxd=xd()
%反三角波序列
forn=0:
7
i=n+1;
if(n>=0&&n<=3)
xd(i)=4-n;
elseif(n>=4)
xd(i)=n-3;
else
xd(i)=0;
end
end
end
N=8;
XD=fft(xd,N);
n=0:
7;
k=0:
7;
subplot(2,2,1);stem(n,xd);title('序列xd');
xlabel('n');ylabel('xd(n)');
subplot(2,2,2);stem(k,XD);title('FFT变换');
xlabel('k');ylabel('XD(k)');
subplot(2,2,3);stem(k,abs(XD));title('幅度');
xlabel('k');ylabel('|XD(k)|');
subplot(2,2,4);stem(k,angle(XD));title('相位');
xlabel('k');
>>xd;
从图中的幅频相位图可以看出,它们的幅度特性是偶对称的,相位特性是奇对称。
实验中的信号序列
和
,在单位圆上的Z变换频谱
和
会相同吗?
如果不同,你能说出哪一个低频分量更多一些吗?
为什么?
functionwork2(xc,xd,N)
XC=0;
XD=0;
w=0:
0.1:
pi;
forn=1:
N
XC=XC+xc(n).*exp(-j*(n-1)*w);
XD=XD+xd(n).*exp(-j*(n-1)*w);
end
Hmc=abs(XC);Hmd=abs(XD);
subplot(2,1,1);plot(w/pi,Hmc);title('Xc(w)幅度响应');grid;
subplot(2,1,2);plot(w/pi,Hmd);title('Xd(w)幅度响应');grid;
>>xc0=[1,2,3,4,4,3,2,1];xd0=[4,3,2,1,1,2,3,4];
>>work2(xc0,xd0,8);
3、已知余弦信号如下
当信号频率
,采样间隔
,采样长度
时,对该信号进行傅立叶变换。
用FFT程序分析正弦信号,分别在以下情况下进行,并且分析比较结果
(1)F=50,N=32,T=0.000625;
(2)F=50,N=32,T=0.005;
(3)F=50,N=32,T=0.0046875;
(4)F=50,N=32,T=0.004;
(5)F=50,N=64,T=0.000625
function[x,X]=work3(F,N,T)
n=1:
N;
x=cos(2*pi*F*n*T);
X=fft(x,N);
functionwork3_1(x1,x2,x3,x4,x5,X1,X2,X3,X4,X5)
n1=1:
32;
n2=1:
64;
subplot(5,2,1);stem(n1,x1);xlabel('n');ylabel('x1(n)');
subplot(5,2,3);stem(n1,x2);xlabel('n');ylabel('x2(n)');
subplot(5,2,5);stem(n1,x3);xlabel('n');ylabel('x3(n)');
subplot(5,2,7);stem(n1,x4);xlabel('n');ylabel('x4(n)');
subplot(5,2,9);stem(n2,x5);xlabel('n');ylabel('x5(n)');
subplot(5,2,2);stem(n1,X1);xlabel('k');ylabel('X1(k)');
subplot(5,2,4);stem(n1,X2);xlabel('k');ylabel('X2(k)');
subplot(5,2,6);stem(n1,X3);xlabel('k');ylabel('X3(k)');
subplot(5,2,8);stem(n1,X4);xlabel('k');ylabel('X4(k)');
subplot(5,2,10);stem(n2,X5);xlabel('k');ylabel('X5(k)');
>>[x1,X1]=work3(50,32,0.000625);
>>[x2,X2]=work3(50,32,0.005);
>>[x3,X3]=work3(50,32,0.0046875);
>>[x4,X4]=work3(50,32,0.004);
>>[x5,X5]=work3(50,64,0.000625);
>>work3_1(x1,x2,x3,x4,x5,X1,X2,X3,X4,X5);
4、选定某一时间信号进行N=64点离散傅立叶变换,并且,对同一信号进行快速傅立叶变换,并比较它们的速度。
>>formatlong;tic;n=1:
64;x=cos(2*pi/32*n);X1=dft(x,64);toc;
Elapsedtimeis0.000000seconds.
>>formatlong;tic;n=1:
64;x=cos(2*pi/32*n);X2=fft(x,64);toc;
Elapsedtimeis0.000000seconds.
程序太小,程序运算时间约为0.
一、实验要求
1、调试实验程序,并且,给参考程序加注释;
2、完成实验内容2,并对结果进行分析。
实验中的信号序列
和
,在单位圆上的Z变换频谱
和
会相同吗?
如果不同,你能说出哪一个低频分量更多一些吗?
为什么?
3、完成实验内容3,并对结果进行分析;
4、利用编制的计算卷积的计算程序,分别给出一下三组函数的卷积结果。
(1)
(2)
(3)
、
functiony1=work4()
n1=0:
14;x1=1.^n1;h1=(0.8).^n1;%
(1)
n21=0:
9;n22=0:
19;x2=1.^n21;h2=0.5*sin(0.5*n22);%
(2)
n3=0:
9;x3=1-(0.1)*n3;h3=(0.1)*n3;%(3)
%
(1)
l1=length(x1)+length(h1)-1;
X1=fft(x1,l1);H1=fft(h1,l1);
y1=ifft(X1.*H1);
%
(2)
l2=length(x2)+length(h2)-1;
X2=fft(x2,l2);H2=fft(h2,l2);
y2=ifft(X2.*H2);
%(3)
l3=length(x3)+length(h3)-1;
X3=fft(x3,l3);H3=fft(h3,l3);
y3=ifft(X3.*H3);
%showresult
figure
(1);
subplot(3,1,1);stem(n1,x1);xlabel('n');ylabel('x1(n)');
subplot(3,1,2);stem(n1,h1);xlabel('n');ylabel('h1(n)');
subplot(3,1,3);stem(1:
l1,real(y1));xlabel('n');ylabel('y1(n)');
figure
(2);
subplot(3,1,1);stem(n21,x2);xlabel('n');ylabel('x2(n)');
subplot(3,1,2);stem(n22,h2);xlabel('n');ylabel('h2(n)');
subplot(3,1,3);stem(1:
l2,real(y2));xlabel('n');ylabel('y2(n)');
figure(3);
subplot(3,1,1);stem(n3,x3);xlabel('n');ylabel('x3(n)');
subplot(3,1,2);stem(n3,h3);xlabel('n');ylabel('h3(n)');
subplot(3,1,3);stem(1:
l3,real(y3));xlabel('n');ylabel('y3(n)');
(1)
(2)
(3)
、
实验三IIR数字滤波器设计
一、实验内容如下:
1、设计一个巴特沃思数字低通滤波器,设计指标如下:
通带内
幅度衰减不大于1dB;
阻带
幅度衰减不小于15dB;
%非归一化巴特沃斯模拟滤波器原型函数
function[b,a]=buttap_o(N,omega)
[z,p,k]=buttap(N);
p=p*omega;
k=k*omega^N;
B=real(poly(z));
b0=k;b=k*B;a=real(poly(p));
%非归一化巴特沃斯模拟滤波器设计函数
function[b,a]=afd_buttap(wp,ws,Rp,As)
ifwp<=0
error('通带边缘必须大于0');
end
ifws<=0
error('阻带边缘必须大于通带边缘');
end
if(Rp<=0)|(As<0)
error('通带波动或阻带波动衰减必须大于0');
end
N=ceil((log10((10^(Rp/10)-1)/(10^(As/10)-1)))/(2*log10(wp/ws)));
fprintf('\n***ButterworthFilterOrder=%2.0f\n',N);
omega=wp/((10^(Rp/10)-1)^(1/(2*N)));
[b,a]=buttap_o(N,omega);
%计算模拟滤波器频率响应的函数
function[db,mag,pha,w]=freqs_m(b,a,wmax)
w=[0:
1:
500]*wmax/500;
H=freqs(b,a,w);
mag=abs(H);
pha=angle(H);
db=20*log10((mag+eps)/max(mag));
%冲激响应不变法设计IIR数字滤波器
function[b,a]=imp_invr(c,d,T)
[R,p,k]=residue(c,d);
p=exp(p*T);
[b,a]=residuez(R,p,k);
b=real(b');a=real(a');
%计算离散系统频率响应
function[db,mag,pha,grd,w]=freqz_m(b,a)
[H,w]=freqz(b,a,1000,'whole');
H=(H(1:
1:
501))';w=(w(1:
1:
501))';
mag=abs(H);
db=20*log10((mag+eps)/max(mag));
pha=angle(H);
grd=grpdelay(b,a,w);
2、编制计算设计的数字滤波器幅度特性和相位特性的程序,并进行实验验证。
>>wp=0.2*pi;ws=0.35*pi;Rp=1;As=15;T=1;
op=wp/T;os=ws/T;
[cs,ds]=afd_buttap(op,os,Rp,As);
[b,a]=imp_invr(cs,ds,T);
[db,mag,pha,grd,w]=freqz_m(b,a);
subplot(2,2,1);plot(w/pi,mag);title('幅度响应');grid;ylabel('幅度');xlabel('以\pi为单位的频率');
subplot(2,2,2);plot(w/pi,pha);title('相位响应');grid;ylabel('相位');xlabel('以\pi为单位的频率');
subplot(2,2,3);plot(w/pi,db);title('幅度响应(dB)');grid;ylabel('对数幅度/dB');xlabel('以\pi为单位的频率');
subplot(2,2,4);plot(w/pi,grd);title('群延迟');grid;ylabel('样本');xlabel('以\pi为单位的频率');
***ButterworthFilterOrder=5
二、实验要求
1、调试实验程序,并且,给实验程序加注释;
2、根据实验结果,给出自己设计的数字滤波器的幅度特性和相位特性;
3、用所设计的滤波器对不同频率的正弦波信号进行滤波,以说明其特性;
4、fp=0.2KHz,Rp=1dB,fs=0.3KHz,As=25dB,T=1ms;分别用脉冲响应不变法及双线性变换法设计一Butterworth数字低通滤波器,观察所设计数字滤波器的幅频特性曲线,记录带宽和衰减量,检查是否满足要求。
比较这两种方法的优缺点。
冲激响应不变法:
>>fp=200;fs=300;Rp=1;As=25;T=0.001;
wp=2*pi*fp*T;ws=2*pi*fs*T;
op=wp/T;os=ws/T;
[cs,ds]=afd_buttap(op,os,Rp,As);
[b,a]=imp_invr(cs,ds,T);
[db,mag,pha,grd,w]=freqz_m(b,a);
subplot(2,2,1);plot(w/pi,mag);title('幅度响应');grid;ylabel('幅度');xlabel('以\pi为单位的频率');
subplot(2,2,2);plot(w/pi,pha);title('相位响应');grid;ylabel('相位');xlabel('以\pi为单位的频率');
subplot(2,2,3);plot(w/pi,db);title('幅度响应(dB)');grid;ylabel('对数幅度/dB');xlabel('以\pi为单位的频率');
subplot(2,2,4);plot(w/pi,grd);title('群延迟');grid;ylabel('样本');xlabel('以\pi为单位的频率');
***ButterworthFilterOrder=9
双线性变换法:
>>fp=200;fs=300;Rp=1;As=25;T=0.001;
wp=2*pi*fp*T;ws=2*pi*fs*T;
op=(2/T)*tan(wp/2);os=(2/T)*tan(ws/2);
[cs,ds]=afd_buttap(op,os,Rp,As);
[b,a]=bilinear(cs,ds,1/T);
[db,mag,pha,grd,w]=freqz_m(b,a);
subplot(2,2,1);plot(w/pi,mag);title('幅度响应');grid;ylabel('幅度');xlabel('以\pi为单位的频率');
subplot(2,2,2);plot(w/pi,pha);title('相位响应');grid;ylabel('相位');xlabel('以\pi为单位的频率');
subplot(2,2,3);plot(w/pi,db);title('幅度响应(dB)');grid;ylabel('对数幅度/dB');xlabel('以\pi为单位的频率');
subplot(2,2,4);plot(w/pi,grd);title('群延迟');grid;ylabel('样本');xlabel('以\pi为单位的频率');
***ButterworthFilterOrder=6
实验四FIR数字滤波器设计
一、实验内容
1、设计一个FIR数字滤波器,设计指标如下:
通带内
幅度衰减不大于1dB;
阻带
幅度衰减不小于15dB;
%理想低通滤波器的计算函数
functionhd=ideal_lp(wc,M);
alpha=(M-1)/2;
n=[0:
1:
(M-1)];
m=n-alpha+eps;
hd=sin(wc*m)./(pi*m);
>>wp=0.2*pi;ws=0.35*pi;N=20;n=[0:
1:
N-1];
wc=(ws+wp)/2;%理想低通截止频率
hd=ideal_lp(wc,N);%理想低通冲激响应
w_ham=(hamming(N))';%Hamming窗
h
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数字信号 处理 编程 作业 机电
![提示](https://static.bdocx.com/images/bang_tan.gif)