第1页 / 共11页
第2页 / 共11页
第3页 / 共11页
第4页 / 共11页
第5页 / 共11页
第6页 / 共11页
第7页 / 共11页
第8页 / 共11页
数字信号实验代码2.pdf
1、用 MATLAB 实现下列序列，并画出图形：
① 单位采样序列移位，
；
提示：实现单位采样序列：
，可通过以下语句实现：x=zeros(1,N);x(1)=1;
n=0:10;
x=[zeros(1,3),1,zeros(1,7)];
figure(1);
stem(n,x);
② 单位阶跃序列移位，
提示：实现单位阶跃序列：
，可通过以下语句实现：x=ones(1,N);
，其中 A=2；f=10； =0.005；
n=0:10;
x=[zeros(1,3),1,ones(1,7)];
figure(2);
stem(n,x);
③ 正弦序列，
n=0:10;
A=2;f=10;Ts=0.005;
xn=2*sin(A*pi*f*n*Ts);
figure(3);
stem(n,xn);
④ 指数序列，
n=0:10;
xn=0.9.^n;
figure(4);
stem(n,xn);
⑤ 复指数序列，
，画出该序列的实部、虚部，幅值和相位。
提示：可通过下列语句实现：
实部 real(x)，虚部 imag(x)，幅值 abs(x)，相位 angle(x)
n=-20:20;
xn=exp(0.05+j*pi/4*n);
figure(5);
subplot(411);
stem(n,real(xn));
subplot(412);
stem(n,imag(xn));
subplot(413);
stem(n,abs(xn));
subplot(414);
stem(n,angle(xn));
2、用 MATLAB 实现两个序列相加：
序列 1：x1=[1 0.5 0.3 0]，n1=1:4；
序列 2：x2=[0.2 0.3 0.4 0.5 0.8 1]，n2=1:6；实现 x=x1+x2，n=1:8，并画出 x 的图
形。
提示：MATLAB 中可用算术运算符“+”实现序列相加，但两个序列的长度必须相等。如果序
列长度不等，或者长度虽然相等但采样的位置不同，就不能运用“+”了。当两序列的长度不
等或位置不对应时，首先应使两者位置对齐，然后通过 zeros 函数左右补零使其长度相等后
再进行相加。
x1=[1 0.5 0.3 0];x1(8)=0;
x2=[0.2 0.3 0.4 0.5 0.8 1];x2(8)=0;n=1:8;
100),3()(−=nnnx0001{)(==nnn100),3()(−=nnunx0001{)(==nnnu100),****2sin(*)(=nTnfAnxssT100,9.0)(=nnxn0.05+j*pi/4*()e,2020nxnn=−
x=x1+x2;
figure(6);
stem(n,x);
3、用 MATLAB 实现序列的反转：
实现
，序列 x(n)采用
，并画出 y(n)的图形。
提示：可利用 fliplr(x) 函数，例如：x=[1 2 3 4];y=fliplr(x);结果为 y=[4 3 2 1]，要实现
，需要对 fliplr(x) 函数进行合理运用。
n=0:10;
xn=[zeros(1,3),1,ones(1,7)];
yn=fliplr(xn);%·×ª
m=-fliplr(n);
figure(7);
stem(m,yn);
4、序列的尺度变换，实现插值和抽取：
已知序列
，用 MATLAB 分别实现下列尺度变换。
提示：可对序列 x(n)的下标进行取余计算，余数为零即为插值和抹去的点，函数如下：
mod(nx,m)，nx 为序列 x(n)的下标，m 为插值或抽取的倍数。
xn=[1 2 3 4 5 6 7 8];
n=-4:3;
figure(8);
subplot(211);
y1=xn(1:2:length(n));n1=-2:1;
stem(n1,y1);
subplot(212);
y2=zeros(1,16);
for m=1:8
y2(2*m)=xn(m);n2(2*m)=n(m);
end
stem(-8:6,y2(2:end));
1.求两个序列的卷积（离散线性卷积）
对于一个 LSI 系统，设其输入序列为
，其中
。系统为单
位冲激响应
，编程求出输出序列，并画出图形。
提示：输出序列应是输入序列和系统单位冲激相应序列的卷积，在 MATLAB 中可用 conv 函
数来实现两个序列的卷积。但是使用 conv 函数求卷积时，两序列卷积后的输出序列的下标
不是任意的（具体是多少，查函数），所以需要自己构造下标。
n=-5:20;
u1=stepfun(n,0);u2=stepfun(n,10);
u=u1-u2;
h=0.9.^n.*u1;
y=conv(u,h);
figure(1);
subplot(311);
stem(n,u1);
subplot(312);
)()(nxny−=100),3()(−=nnunx)()(nxny−=(){1,2,3,4,5,6,7,8},-43xnn=12()(2){()(/2)ynxnynxn==)10()()(−−=nununx205−n)(*9.0)(nunnh=
stem(n,u);
subplot(313);
stem(-10:40,y);
2.求两个序列的卷积（离散周期卷积）
已知序列
，
，计算它们的周期卷积。
提示：将
和
的主值序列 x(n)和 y(n)作线性卷积，然后再将线性卷积序列以 N 为周
期进行周期延拓。
x=[1 0 1 2 1];
y=[1 1 0 1 2];
z=conv(x,y);
h1=[z(6:9),0]; h2=h1+z(1:5); h3=[h2 h2 h2];
figure(2);
stem(h3);
3.已知 LSI 系统的差分方程为：
，其中
求单位冲激响应
，并画出图形。
提示：可以使用 filter 函数来完成，其调用格式为
，注意 b，a，x 参数都
代表什么，如何设置。
n=-5:30;
b=[1]; a=[1 -1 0.9];
x=[zeros(1,5),1,zeros(1,30)];
H=filter(b,a,x);
figure(3);
stem(n,H);
4.已知 LSI 系统的差分方程为：
求单位阶跃响应
，并画出图形。
，其中
提示：可以使用 filter 函数来完成，其调用格式为
，注意 b，a，x 参数都
n=-5:30;
b=[1]; a=[1 -1 0.9];
x=stepfun(n,0);
g=filter(b,a,x);
figure(4);
stem(n,g);
1. 求下列序列的 z 变换：
①
②
提示：将序列 x(n)符号化，然后利用 ztrans()函数实现对 x(n)进行 z 变换，其结果为 X(z)；
simple()函数用于获取符号变量的最简形式；pretty()函数用于将符号表达式化简成与高等数
学课本上显示符号表达式类似的形式。
syms n;
%un=stepfun(n,0);
xn1=cos(n.*pi/2);
~(){1,0,1,2,1}xn=~y(){1,1,0,1,2}n=~()xn~y()n)()2(9.0)1()(nxnynyny=−+−−305−n)(nh),,(xabfilterH=)()2(9.0)1()(nxnynyny=−+−−305−n)(ng),,(xabfilterH=)()2cos()(nunnx=)(2)1()(nunnnx−=
xn2=(n.*(n-1))/2;
Xz1=ztrans(xn1);
Xz2=simplify(ztrans(xn2));
pretty(Xz1);
pretty(Xz2);
2. 求 z 反变换：
①
②
提示：将序列 x(n)符号化，然后利用 iztrans()函数实现求函数 X(z)的 z 反变换 x(n)；然后利用
pretty()函数实现与书本上显示一致的符号表达。
syms z;
Xz1=z/(z-2)^2;
Xz2=z*(z-3)/(z^2+2*z+1);
xn1=iztrans(Xz1);
xn2=iztrans(Xz2);
pretty(xn1);
pretty(xn2);
3. 离 散 时 间 线 性 时 不 变 系 统 的 时 域 实 现 ： 已 知 因 果 LSI 系 统 的 差 分 方 程 为 ：
，写出系统传递函数 H(z)及其收敛域，编程实现以下
内容：
①求系统的零极点，并绘图；
②求系统冲激响应 h(n)，并绘图；
③求系统的单位阶跃响应 g(n)，并绘图；
④求系统函数的幅频响应 H(ejw)，并绘出幅频和相频特性。
提示：b 是差分方程中输入的系数，a 是差分方程中输出的系数。注意缺项的系数前要添加
系数 0。①zplane(b,a)用来实现画出系统的零极点图；②impz(b,a,n)函数，如果 n 是一个整
形变量，impz()计算的是在这些整数点位置上的冲激响应，且从 0 开始。单位阶跃响应可用
stepz(b,a,n)函数，与 impz()函数用法相似。③幅频响应可利用 freqz(b,a,w)函数，其中
w=[0:1:500]*pi/500 或 freqz(b,a,n)函数，可查看 MATLAB Help 了解它的相关用法。③幅频响
应可用 freqz(b,a,w)，其中 w=[0:1:500]*pi/500 和 freqz(b,a,n)分别求，可查看 MATLAB Help
了解它的相关用法。
m=-5:30;
b=[1 0 -1]; a=[1 0 -0.81];
%(1)
figure(1);
zplane(b,a);
%(2)
figure(2);
impz(b,a,100);
%(3)
figure(3);
stepz(b,a,100);
%(4)
w=[0:1:100]*pi/25;
figure(4);
freqz(b,a,100,'whole');
figure(5);
2)2()(−=zzzX12)3()(2++−=zzzzzX)2()()2(81.0)(−−+−=nxnxnyny
freqz(b,a,w);
（1）若
是一个 N=12 的有限长序列，利用 MATLAB 计算它的 DFT 并画出
图形。
提示：在计算机上实现信号的频谱分析及其他方面的处理工作时，对信号的要求是：在时域
和频域都应该是离散的，且都应是有限长的，只有 DFS（离散周期序列的傅立叶级数）在时
域和频域都是离散的。因此，DFT 实际上来是自于 DFS 的，只不过仅在时域、频域各取一个
周期而已。本题要求自己编写一个
函数，其中，xn 代表序列 x(n)，N 代表有限长
序列的点数。
函数的算法采用 DFT 变换的基本公式，即：
，其中，
。由于 X(k)为复数，所以画
图时，应画 X(k)的幅度特性 abs()和相位特性 angle()。
N=12; n=0:N-1; k=0:N-1;
x=cos(n*pi/6);
WN=exp(-j*(2*pi/N)*n'*k);
X=x*WN;figure(1);
X1=abs(X); X2=angle(X);
subplot(121);stem(n,X1);
subplot(122);stem(n,X2);
（2）对于实序列
① 分解成偶部
和奇部
；
n=0:19; N=length(n); m=mod(-n,N)+1;
x1=0.9.^n; x2=x1(m);
xe=1/2*(x1+x2); xo=1/2*(x1-x2);figure(2);
subplot(321);stem(n,xe);title('x1(n)ż');
subplot(322);stem(n,xo);title('x1(n)Æ æ');
② 验证实序列的性质：
提 示 ： 任 意 一 个 实 序 列 都 可 以 表 示 成 偶 对 称 序 列 与 奇 对 称 序 列 之 和 ， 即 ：
。令
，
，
验证实序列的上述性质。注：因为 MATLAB 的计算精度问题，理论上，
只包含
实部，
只包含虚部，但实际计算中仅是
的虚部很小，
的实部很小。因此，比较时：
本题要求自己编写一个 circevod(x)函数，它的功能是将一个实序列分解成
与
之
和，并分别求出
与
。其中需要有判断语句，判断序列 x(n)是否是实序列并有相应
文字输出；
的求法可参考 mod(-n,N)函数。
)6cos()(nnx=),(Nxndft),(Nxndft−=−==−=1010)()2exp()()(NnNnnkNWnxnkNjnxkX1,,1,0−=Nk200,)9.0()(=nnxn)(nxe)(nxo)](Re[)]([kXnxDFTe=)](Im[)]([kXjnxDFTo=)()()(nxnxnxoe+=]))(()([21)(NenNxnxnx−+=]))(()([21)(NonNxnxnx−−=)]([nxDFTe)]([nxDFTo)]([nxDFTe)]([nxDFTo)])([()])([(nxDFTrealnxDFTreale=)])([()])([(nxDFTimagnxDFTimago=)(nxe)(nxo)(nxe)(nxoNnNx))((−
Xe=dft(xe,N); Xo=dft(xo,N);
X=dft(x1,N); ReX=real(X); ImX=imag(X);
subplot(323);stem(n,Xe);title('xe(n)µÄDFT');
subplot(324);stem(n,ReX);title('Xʵ²¿');
subplot(325);stem(n,imag(Xo));title('xoµÄDFT');
subplot(326);stem(n,ImX);title('XÐ鲿');
（3）给定模拟信号
，以
秒进行采样，
用 DFT 进行信号频谱分析：
①若要能区别该信号中的两个频率分量，试问取样信号的长度 N 至少为多少？
②DFT 的点数 N 分别等于 64、256、1024，画出信号的 N 点 DFT 的幅度谱和相位谱。
N1=64; n1=0:N1-1;
N2=256;n2=0:N2-1;
N3=1024;n3=0:N3-1;
t=0.01*n1;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
X=dft(x,N1);figure(3)
subplot(321);
plot(n1*100/N1,abs(X));axis([0 100 0 150]);
title('x(t)--N=64');
subplot(322);plot(n1*100/N1,angle(X));
title('abs(X)--N=64');
t=0.01*n2;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
X2=dft(x,N2);
subplot(323);plot(n2*100/N2,abs(X2));
title('x(t)--N=256');
subplot(324);plot(n2*100/N2,angle(X2));
title('abs(X)--N=256');
t=0.01*n3;
x=2*sin(4*pi*t)+5*cos(8*pi*t);
X3=dft(x,N3);
subplot(325);plot(n3*100/N3,abs(X3));
title('x(t)--N=1024');
subplot(326);plot(n3*100/N3,angle(X3));
title('abs(X)--N=1024');
DFT
function Xk= dft(xn,N)
n=0:N-1; k=0:N-1;
WN=exp(-j*(2*pi/N)*n'*k);
Xk=xn*WN;
end
（1）求序列 x=[5 2 7 4 1 1 3 9]的快速傅立叶变换，并绘出幅频特性曲线；
提示：采用 fft(x,N)函数，N 为 FFT 的点数，当序列长度小于 N 时，系统自动在序列末尾补
零；当序列长度大于 N 时，系统自动截断序列多余的部分。分别画出幅频特性曲线 abs()和
相频特性曲线 angle()。
x=[5 2 7 4 1 1 3 9];
N=length(x);
)8cos(5)4sin(2)(tttx+=)1:0(01.0−==Nnnt
X=fft(x,N);
figure(1);
subplot(211);
stem(0:N-1,abs(X));
subplot(212);
stem(0:N-1,angle(X));
（2）设一序列中含有两种频率成分，f =2 Hz，f =2.05Hz，采样频率取为 fs=10Hz，即：
x(n)=sin(2 f n/f
)+sin(2 f n/f
)，根据公式
<
，要区分出这两种频率成分，N
必须满足多少？
①取 x(n)(0 n<128)时，计算 128 点 FFT，并绘出幅频特性曲线；
②取 x(n)(0 n<128)时，补 384 个 0，计算 512 点 FFT，并绘出幅频特性曲线；
③取 x(n)(0 n<512)时，计算 512 点 FFT，并绘出幅频特性曲线。
f1=2; f2=2.05; fs=10;
figure(2);
%1
N=128;n=0:N-1;
xn=sin(2*pi*f1*n/fs)+sin(2*pi*f2*n/fs);
X=fft(xn,N);
subplot(321);
plot(n*fs/N,abs(X));
subplot(322);
plot(n*fs/N,angle(X));
%2
N=512;n=0:127;
xn=sin(2*pi*f1*n/fs)+sin(2*pi*f2*n/fs);
xn(512)=0;
X=fft(xn,N);
subplot(323);
n=0:N-1;
plot(n*fs/N,abs(X));
subplot(324);
plot(n*fs/N,angle(X));
%3
N=512;n=0:N-1;
xn=sin(2*pi*f1*n/fs)+sin(2*pi*f2*n/fs);
X=fft(xn,N);
subplot(325);
plot(n*fs/N,abs(X));
subplot(326);
plot(n*fs/N,angle(X));
（3）计算序列 y(n)的自相关函数，其中：
x(n)=sin(2 n/8)，y(n)=e(n)+x(n)，式中 e(n)为一随机噪声，N=100。
提示：随机噪声 e(n)可通过 rand(1,N)函数构造，rand 函数产生[0,1]分布的随机数，均值为
0.5，(1,N)表示的是 1 行 N 列；本程序要求编写产生均值为 0，分布在[-0.5,0.5] 的随机数，
则 e= rand(1,N)-0.5；求 y(n)的自相关函数，可采用 ry=xcorr(yn,yn,'unbiased')函数实现。本
实验要求取
，画出 x(n)，y(n)和 y(n)的自相关函数 ry 的波形。
n=0:1:99; N=100;
en=rand(1,N)-0.5;
xn=sin(2*pi*n/8);
yn=en+xn;
121s2sNfs221ff−990n
ry=xcorr(yn,yn,'unbiased');%unbiased-ÎÞÆ «¹À¼Æ
figure(3);
subplot(311);
plot(n,xn);
subplot(312);
plot(n,yn);
subplot(313);
plot(0:2*length(n)-2,ry);
（1）归一化模拟低通滤波器的设计
①巴特沃斯低通滤波器的设计：3dB 通带截止频率和阶数 N 可以完全描述一个巴特沃
斯归一化模拟低通滤波器。编程实现：滤波阶次分别取 N=2、5、10、20，在同一图中绘出
巴特沃斯归一化模拟低通原型滤波器的平方幅频响应曲线。
提示：采用[z0,p0,k0]=buttap(N)：N 为巴特沃斯归一化模拟低通滤波器的阶数，返回值为滤
波器组的零点数 z0、极点数组 p0 和增益 k0。[b,a]=zp2tf(z0,p0,k0)该函数从系统的零极点及
增益得到滤波器系统函数的分析多项式系数向量 B 和分母多项式系数向量 A。w=0:2π，
H=freqs(b,a,w)函数求模拟系统的频率响应，用法类似于 freqz。画出平方幅频特性曲线采用
abs(H).^2。在同一图中
N=[2 5 10 20];
w=0:2*pi*100:2*pi*30000;
figure(1);
for i=1:4
[z,p,k]=buttap(N(i));
[b,a]=zp2tf(z,p,k);
[h,w]=freqs(b,a,length(w));
plot(w,abs(h).^2);axis([0 2 0 1]);
hold on;
end
②切比雪夫 I 型低通滤波器的设计：编程实现：滤波阶次分别取 N=2、4、6、8，在同
一图中绘出切比雪夫 I 型模拟低通原型滤波器的平方幅频响应曲线。
提示：采用[z0,p0,k0]=cheb1ap(N,Rp)，N 为切比雪夫 I 型归一化模拟低通滤波器的阶数，Rp
为通带最大衰减。[B,A]=zp2tf(z0,p0,k0) ；w=0:2π，H=freqs(b,a,w)函数；其中 Rp=1dB。
N=[2 4 6 8];
w=0:2*pi*100:2*pi*30000;Rp=1;
figure(2);
for i=1:4
[z,p,k]=cheb1ap(N(i),Rp);
[b,a]=zp2tf(z,p,k);
[h,w]=freqs(b,a,length(w));
subplot(2,2,i);
plot(w,abs(h).^2);
end
（2）IIR 数字滤波器设计
IIR 数字滤波器的设计过程
①冲激响应不变法：编程实现：设计巴特沃斯低通滤波器，通带衰减为 Rp=1dB，通带
上限角频率 =0.2π，阻带下限角频率 =0.3π，阻带最小衰减 =15dB，画出设计后的数
实际数字滤波器指标归一化低通滤波器设计指标Ha(p)Ha(s)H(z)MATLAB或查表频率变换脉冲响应不变或双线性变换pssa
相关推荐
2023年江西萍乡中考道德与法治真题及答案.doc
2012年重庆南川中考生物真题及答案.doc
2013年江西师范大学地理学综合及文艺理论基础考研真题.doc
2020年四川甘孜小升初语文真题及答案I卷.doc
2020年注册岩土工程师专业基础考试真题及答案.doc
2023-2024学年福建省厦门市九年级上学期数学月考试题及答案.doc
2021-2022学年辽宁省沈阳市大东区九年级上学期语文期末试题及答案.doc
2022-2023学年北京东城区初三第一学期物理期末试卷及答案.doc
2018上半年江西教师资格初中地理学科知识与教学能力真题及答案.doc
2012年河北国家公务员申论考试真题及答案-省级.doc
2020-2021学年江苏省扬州市江都区邵樊片九年级上学期数学第一次质量检测试题及答案.doc
2022下半年黑龙江教师资格证中学综合素质真题及答案.doc
