02356数字信号处理

第二章2.1 判断下列序列是否是周期序列。

若是，请确定它的最小周期。

（1）x(n)=Acos(685ππ+n ) （2）x(n)=)8(π-ne j（3）x(n)=Asin(343ππ+n ) 解 (1)对照正弦型序列的一般公式x(n)=Acos(ϕω+n )，得出=ω85π。

因此5162=ωπ是有理数，所以是周期序列。

最小周期等于N=)5(16516取k k =。

（2）对照复指数序列的一般公式x(n)=exp[ωσj +]n,得出81=ω。

因此πωπ162=是无理数，所以不是周期序列。

（3）对照正弦型序列的一般公式x(n)=Acos(ϕω+n )，又x(n)=Asin(343ππ+n )＝Acos(-2π343ππ-n )＝Acos(6143-n π)，得出=ω43π。

因此382=ωπ是有理数，所以是周期序列。

最小周期等于N=)3(838取k k =2.2在图2.2中，x(n)和h(n)分别是线性非移变系统的输入和单位取样响应。

计算并列的x(n)和h(n)的线性卷积以得到系统的输出y(n)，并画出y(n)的图形。

(a)1111(b)(c)111110 0-1-1-1-1-1-1-1-1222222 33333444………nnn nnnx(n)x(n)x(n)h(n)h(n)h(n)21u(n)u(n)u(n)a n ===22解 利用线性卷积公式y(n)=∑∞-∞=-k k n h k x )()(按照折叠、移位、相乘、相加、的作图方法，计算y(n)的每一个取样值。

(a) y(0)=x(O)h(0)=1y(l)=x(O)h(1)+x(1)h(O)=3y(n)=x(O)h(n)+x(1)h(n-1)+x(2)h(n-2)=4,n ≥2 (b) x(n)=2δ(n)-δ(n-1)h(n)=-δ(n)+2δ(n-1)+ δ(n-2)y(n)=-2δ(n)+5δ(n-1)= δ(n-3) (c) y(n)=∑∞-∞=--k kn k n u k u a)()(=∑∞-∞=-k kn a=aa n --+111u(n)2.3 计算线性线性卷积 (1) y(n)=u(n)*u(n) (2) y(n)=λnu(n)*u(n)解：(1) y(n)=∑∞-∞=-k k n u k u )()(=∑∞=-0)()(k k n u k u =(n+1),n ≥0即y(n)=(n+1)u(n) (2) y(n)=∑∞-∞=-k k k n u k u )()(λ=∑∞=-0)()(k kk n u k u λ=λλ--+111n ,n ≥0 即y(n)=λλ--+111n u(n)2.4 图P2.4所示的是单位取样响应分别为h 1(n)和h 2(n)的两个线性非移变系统的级联，已知x(n)=u(n), h 1(n)=δ(n)-δ(n-4), h 2(n)=a nu(n),|a|<1,求系统的输出y(n).解 ω(n)=x(n)*h 1(n) =∑∞-∞=k k u )([δ(n-k)-δ(n-k-4)]=u(n)-u(n-4)y(n)=ω(n)*h 2(n) =∑∞-∞=k kk u a )([u(n-k)-u(n-k-4)]=∑∞-=3n k ka,n ≥32.5 已知一个线性非移变系统的单位取样响应为h(n)=an-u(-n),0<a<1 用直接计算线性卷积的方法，求系统的单位阶跃响应。

《数字信号处理》研究生复试大纲南京理工大学研究生复试大纲学院（系）：电光学院课程名称：数字信号处理执笔人：宋耀良修（制）订日期： 2003年3月一、课程的考试目的与基本要求本课程的考试目的主要在于考查学生掌握离散时间信号与系统的分析方法、离散傅里叶变换的概念以及离散傅里叶变换的快速算法、掌握数字滤波器的各种设计方法，以及了解有限字长效应对系统性能的影响等方面内容的情况；通过考查同时检验学生分析问题，解决实际问题的能力。

要求学生较全面深刻地掌握数字信号变换的基本理论及其在数字系统分析的应用，全面掌握数字系统的基本结构和系统设计和实现的基本方法。

二、考查内容以教学大纲为依据，具体内容包括:1 绪论数字信号和数字信号处理的基本概念2离散时间信号与系统2.1 离散时间信号与系统的频域表示2.1.1 系统的频率响应；2.1.2 系统频率响应的两个性质；2.1.3 系统频率响应与单位取样响应的关系；2.1.4 序列的频域表示法；2.1.5 输出序列与输入序列的傅氏变换间的关系。

2.2 傅里叶变换的对称性质2.2.1 序列的共轭对称和共轭反对称2.2.2 偶序列与奇序列2.2.3 傅里叶变换的共轭对称和共轭反对称2.2.4 傅里叶变换的对称性质3连续时间信号的采样3.1 周期采样3.2 采样的频域表示3.2.1 奈奎斯特采样定理(△)3.3 由样本重构带限信号（★）3.4 连续时间信号的离散时间处理3.4.1 线性时不变离散时间系统3.4.2 冲激响应不变3.5 离散时间信号的连续时间处理3.5.1 非整数延时3.5.2 滑动平均3.6 利用离散时间处理改变采样率（★）3.6.1 采样率按整数因子减小3.6.2 采样率按整数因子增加3.6.3 采样率按非整数因子变化3.7 实际应用考虑的问题3.7.1 消除混叠的预滤波3.7.2 模拟到数字(A/D)转换3.7.3 量化误差分析3.7.4 D/A转换3.7.5 A/D和D/A转换中抽取和内插的应用4 z 变换4.1 z变换的定义和性质4.2 z 变换的存在条件及系统的收敛域4.3 逆z变换5线性时不变 ( LTI ) 系统的变换分析5.1 LTI系统的频率响应5.1.1 理想频率选择性滤波器5.1.2 相位失真和延迟5.2 用线性常系数差分方程所表征的系统函数5.2.1 稳定性和因果性；5.2.2 逆系统（★）5.2.3 有理系统函数的冲激响应5.3 有理系统函数的频率响应5.3.1 单个零点或极点的频率响应5.3.2 多重零极点的例子5.4 幅度和相位特性之间的关系5.5 全通系统5.6 最小相位系统5.6.1 频率响应的补偿5.6.2 最小相位系统的性质5.7 广义线性相位的线性系统（△）5.7.1 线性相位系统5.7.2 广义线性相位5.7.3 因果广义线性相位系统5.7.4线性相位系统与最小相位系统的关系6 离散时间系统结构6.1 线性常系数差分方程的方框图表示6.2 线性常系数差分方程的信号流图表示6.3 IIR系统的基本结构6.3.1 直接型(△)6.3.2 级联型(△)6.3.3 并联型(△)6.3.4 IIR系统中的反馈6.4 转置形式6.5 FIR系统的基本结构6.5.1 直接型(△)6.5.2 级联型(△)6.5.3 线性相位FIR系统的结构(△)6.6 有限精度数值效应概述6.6.1 数的表示法6.6.2 在系统中实现的量化6.8 系数量化的影响6.7.1 IIR系统中系数量化的影响(△)6.7.2 FIR系统中系数量化的影响6.7.3 其它灵敏度分析方法6.8 数字滤波器中舍入噪声的影响6.8.1 直接型IIR结构分析6.8.2 IIR系统定点实现中的幅度加权6.8.3 级联和并联IIR结构分析6.8.4 直接型FIR结构分析6.8.5 离散时间系统的浮点实现6.9 IIR数字滤波器定点实现中的零输入极限环6.9.1 由于舍入和截尾引起的极限环6.9.2 由于溢出而出现的极限环7滤波器设计方法7.1 由连续时间波器设计离散时间IIR滤波器7.1.1 滤波器设计的冲激响应不变法（△）7.1.2 双线性变换法（△）7.1.3 双线性变换法设计举例7.2 低通IIR滤波器的频率变换法7.3 用窗函数法设计FIR滤波器7.3.1 常用窗函数的性质（△）7.3.2 广义线性相位的合并7.3.3 凯泽窗滤波器设计法（△）7.3.4 凯泽窗与其它窗之间的关系7.4 凯泽窗法设计FIR滤波器举例7.4.1 高通滤波器举例7.4.2 离散时间微分器7.6 FIR滤波器的最佳逼近（★）7.6.1 最佳I型低通滤波器7.6.2 最佳II型低通滤波器7.6.3 帕克斯-麦克莱伦算法7.6.4 最佳FIR滤波器的特性7.7 FIR等波纹逼近举例7.7.1 低通滤波器7.7.2 零阶保持器的补偿7.7.3 带通滤波器8离散傅里叶变换8.1 周期序列的表示：离散傅里叶级数8.2 离散傅里叶级数的性质（△）8.2.1 线性；8.2.2 序列的移位8.2.3 对偶性8.2.4 对称性8.2.5 周期卷积8.3 周期序列DFS表示的性质汇总8.4 周期信号的傅里叶变换8.5 对傅里叶变换采样8.5 有限长序列的傅里叶表示：离散傅里叶变换（△）8.7 离散傅里叶变换的性质（△）8.7.1 线性8.7.2 序列的循环移位（★）8.7.3 对偶性8.7.4对称性8.7.5 循环卷积（★）8.8 离散傅里叶变换的性质汇总8.9 用离散傅里叶变换实现线性卷积（△）8.9.1 两个有限长序列的线性卷积8.9.2 循环卷积作为带有混叠的线性卷积8.9.3 用DFT实现线性时不变系统9 离散傅里叶变换的计算9.1 离散傅里叶变换的高效计算9.2 戈泽尔算法9.3 按时间抽取的FFT算法9.3.1 算法原理(△)9.3.2 算法流图(△)9.3.3 算法特点(△)9.3.4 变形的算法流图(★)9.4 按频率抽取的FFT算法9.4.1 算法原理(△)9.4.2 算法流图(△)9.4.3 算法特点(△)9.4.4 变形的算法流图(★)9.5 FFT算法的实现9.5.1 标号排列9.5.2 系数9.5.3 实现FFT的算法流程图9.5.4 实序列FFT算法9.5 复合数N的FFT算法（★）9.5.1 库利-图基算法9.5.2 素因子算法9.6 用卷积实现DFT线性调频变换算法（△）9.7 离散傅里叶变换计算中有限寄存器长度的影响9.7.1 在FFT的直接法计算中量化的分析9.7.2 定点FFT算法中量化效应的分析9.7.3 浮点FFT算法中量化效应小结9.7.4 在FFT中系数量化的影响Alan V. Oppenheim & Ronald W. Schafer,三. 主要参考教材Discrete-Time Signal Processing,(2and ) Prince Hall. （教材）1997.注:课程复试大纲填写说明:重点△，难点★。

实验五谱分析一，试验目的：1 研究不同类型的窗函数，研究一些不同的方法来测试窗的性能；2 专注于有关窄带信号的几个不同的情形。

二，实验原理信号是无限长的，而在进行信号处理时只能采用有限长信号，所以需要将窄带信号“截断”。

在信号处理中，“截断”被看成是用一个有限长的“窗口”看无限长的信倍号，或者从分析的角度是无限长的信号x(t)乘以有限长的窗函数w(t),由傅立叶变换性质可知：x(t)w(t)==1/2πX(jw)W(jw)如果x(t)是频带有限信号，而w(t)是频带无限函数，截断后的信号也必是频带无限信号，从而产生所谓的频谱泄露。

频谱泄露是不可避免的，但是尽量减小，因此设计了不同的窗函数满足不同的要求。

从能量的角度，频谱泄露也是能量泄露，因为加窗后，是原来的信号集中在宅频带内的能量分散到无限的频带范围。

三，实验内容1．用MATLAB编程绘制各种窗函数的形状（1）矩形窗程序N=20;w=boxcar(N);nn=0:N-1;plot(nn,w)（2）汉宁窗程序N=20;w=hanning(N);nn=0:N-1;plot(nn,w)（3）汉明窗程序N=20;w=hamming(N);nn=0:N-1;plot(nn,w)（4）巴特利特窗程序N=20;w=bartlett(N);nn=0:N-1;plot(nn,w)（5）布莱克曼窗程序N=20;w=blackman(N);nn=0:N-1;plot(nn,w)（6）Triang窗程序N=20;w=triang(N); nn=0:N-1;plot(nn,w)（7）Kaiser窗程序N=20;w=kaiser(N,10); nn=0:N-1;plot(nn,w)（8）切比雪夫窗程序N=20;w=chebwin(N,30);nn=0:N-1;plot(nn,w)2，用MATLAB编程绘制各种窗函数的形状及其幅度响应。

(1)矩形窗N=20;w=boxcar(N);[H,W]=dtft(w,512);subplot(111),plot(W/2/pi,20*log(abs(H)));grid,title('MAGNITUDE RESPONSE')xlabel('NORMALIZED FREQUENCY'),ylabel('|H(w)|/dB')（2）汉宁窗N=20;w=hanning(N);[H,W]=dtft(w,512);subplot(111),plot(W/2/pi,20*log(abs(H)));grid,title('MAGNITUDE RESPONSE')xlabel('NORMALIZED FREQUENCY'),ylabel('|H(w)|/dB')（3）汉明窗N=20;w=hamming(N);[H,W]=dtft(w,512);subplot(111),plot(W/2/pi,20*log(abs(H)));grid,title('MAGNITUDE RESPONSE')xlabel('NORMALIZED FREQUENCY'),ylabel('|H(w)|/dB')（4）巴特利特窗N=20;w=bartlett(N);[H,W]=dtft(w,512);subplot(111),plot(W/2/pi,20*log(abs(H)));grid,title('MAGNITUDE RESPONSE')xlabel('NORMALIZED FREQUENCY'),ylabel('|H(w)|/dB')（5）布莱克曼窗N=20;w=blackman(N);[H,W]=dtft(w,512);subplot(111),plot(W/2/pi,20*log(abs(H)));grid,title('MAGNITUDE RESPONSE')xlabel('NORMALIZED FREQUENCY'),ylabel('|H(w)|/dB')（6）Triang窗N=20;w=triang(N);[H,W]=dtft(w,512);subplot(111),plot(W/2/pi,20*log(abs(H)));grid,title('MAGNITUDE RESPONSE')xlabel('NORMALIZED FREQUENCY'),ylabel('|H(w)|/dB')（7）Kaiser窗N=20;w=kaiser(N,10);[H,W]=dtft(w,512);subplot(111),plot(W/2/pi,20*log(abs(H)));grid,title('MAGNITUDE RESPONSE')xlabel('NORMALIZED FREQUENCY'),ylabel('|H(w)|/dB')（8）切比雪夫窗N=20;w=chebwin(N,30);[H,W]=dtft(w,512);subplot(111),plot(W/2/pi,20*log(abs(H)));grid,title('MAGNITUDE RESPONSE')xlabel('NORMALIZED FREQUENCY'),ylabel('|H(w)|/dB')3，绘制矩形窗的幅频响应，窗长度分别为：N=10，N=20，N=50，N=100。

程序清单及波形显示： clc;close all;clear all;%======内容1：调用filter 解差分方程，由系统对u(n)的响应判断稳定性======A=[1,-0.9];B=[0.05,0.05]; %系统差分方程系数向量B 和A x1n=[1 1 1 1 1 1 1 1 zeros(1,50)]; %产生信号x1(n)=R8(n) x2n=ones(1,128); %产生信号x2(n)=u(n)y1n=filter(B,A,x1n); %求系统对x1(n)的响应y1(n) n=0:length(y1n)-1;subplot(2,2,1);stem(n,y1n,'.'); title('(a) 系统对R8(n)的响应y1(n)');xlabel('n');ylabel('y1(n)');y2n=filter(B,A,x2n); %求系统对x2(n)的响应y2(n) n=0:length(y2n)-1;subplot(2,2,2);stem(n,y2n,'.'); title('(b) 系统对u(n)的响应y2(n)');xlabel('n');ylabel('y2(n)');hn=impz(B,A,58); %求系统单位脉冲响应h(n) n=0:length(hn)-1;subplot(2,2,3);y=hn;stem(n,hn,'.'); title('(c) 系统单位脉冲响应h(n)');xlabel('n');ylabel('h(n)');(a) 系统对R8(n)的响应y1(n)ny 1(n )(b) 系统对u(n)的响应y2(n)ny 2(n )(c) 系统单位脉冲响应h(n)nh (n )%===内容2：调用conv 函数计算卷积============================x1n=[1 1 1 1 1 1 1 1 ]; %产生信号x1(n)=R8(n) h1n=ones(1,10); h2n=[1 2.5 2.5 1 ];y21n=conv(h1n,x1n); y22n=conv(h2n,x1n); figure(2)n=0:length(h1n)-1;subplot(2,2,1);stem(n,h1n); title('(d) 系统单位脉冲响应h1n');xlabel('n');ylabel('h1(n)'); n=0:length(y21n)-1;subplot(2,2,2);stem(n,y21n); title('(e) h1(n)与R8(n)的卷积y21n');xlabel('n');ylabel('y21(n)'); n=0:length(h2n)-1;subplot(2,2,3);stem(n,h2n); title('(f) 系统单位脉冲响应h2n');xlabel('n');ylabel('h2(n)'); n=0:length(y22n)-1;subplot(2,2,4);stem(n,y22n); title('(g) h2(n)与R8(n)的卷积y22n');xlabel('n');ylabel('y22(n)');nh 1(n )ny 21(n )(f) 系统单位脉冲响应h2nnh 2(n)(g) h2(n)与R8(n)的卷积y22nny 22(n )%=========内容3：谐振器分析======================== un=ones(1,256); %产生信号u(n) n=0:255;xsin=sin(0.014*n)+sin(0.4*n); %产生正弦信号A=[1,-1.8237,0.9801];B=[1/100.49,0,-1/100.49]; %系统差分方程系数向量B 和Ay31n=filter(B,A,un); %谐振器对u(n)的响应y31(n) y32n=filter(B,A,xsin); %谐振器对u(n)的响应y31(n) figure(3)n=0:length(y31n)-1;subplot(2,1,1);stem(n,y31n,'.'); title('(h) 谐振器对u(n)的响应y31n');xlabel('n');ylabel('y31(n)'); n=0:length(y32n)-1;subplot(2,1,2);stem(n,y32n,'.'); title('(i) 谐振器对正弦信号的响应y32n');xlabel('n');ylabel('y32(n)');050100150200250300(h) 谐振器对u(n)的响应y31nny 31(n )(i) 谐振器对正弦信号的响应y32nny 32(n )程序清单及波形显示：% DTMF 双频拨号信号产生6位电话号码 %clear all;clc;tm=[1,2,3,65;4,5,6,66;7,8,9,67;42,0,35,68]; % DTMF 信号代表的16个数N=205;K=[18,20,22,24,31,34,38,42];f1=[697,770,852,941]; % 行频率向量 f2=[1209,1336,1477,1633]; % 列频率向量 TN=input('键入6位电话号码= '); % 输入6位数字TNr=0; %接收端电话号码初值为零for l=1:6;d=fix(TN/10^(6-l)) TN=TN-d*10^(6-l); for p=1:4; for q=1:4;if tm(p,q)==abs(d); break,end % 检测码相符的列号q endif tm(p,q)==abs(d); break,end % 检测码相符的行号p endn=0:1023; % 为了发声，加长序列x = sin(2*pi*n*f1(p)/8000) + sin(2*pi*n*f2(q)/8000);% 构成双频信号sound(x,8000); % 发出声音 pause(0.1)% 接收检测端的程序X=goertzel(x(1:205),K+1); % 用Goertzel 算法计算八点DFT样本val = abs(X); % 列出八点DFT向量subplot(3,2,l);stem(K,val,'.');grid;xlabel('k');ylabel('|X(k)|') % 画出DFT(k)幅度axis([10 50 0 120])limit = 80; %for s=5:8;if val(s) > limit, break, end % 查找列号endfor r=1:4;if val(r) > limit, break, end % 查找行号endTNr=TNr+tm(r,s-4)*10^(6-l);enddisp('接收端检测到的号码为：') % 显示接收到的字符disp(TNr)显示结果：键入6位电话号码= 123456d = 1d = 2d = 3d = 4d = 5d = 6接收端检测到的号码为：123456050100|X (k )||X (k )|050100k|X (k )|050100k|X (k )|050100k|X (k )|050100k|X (k )|% DTMF 双频拨号信号产生8位电话号码 %clear all;clc;tm=[1,2,3,65;4,5,6,66;7,8,9,67;42,0,35,68]; % DTMF 信号代表的16个数 N=205;K=[18,20,22,24,31,34,38,42];f1=[697,770,852,941]; % 行频率向量 f2=[1209,1336,1477,1633]; % 列频率向量 TN=input('键入8位电话号码= '); % 输入8位数字TNr=0; %接收端电话号码初值为零for l=1:8;d=fix(TN/10^(8-l)) TN=TN-d*10^(8-l); for p=1:4; for q=1:4;if tm(p,q)==abs(d); break,end % 检测码相符的列号q endif tm(p,q)==abs(d); break,end % 检测码相符的行号p endn=0:1023; % 为了发声，加长序列 x = sin(2*pi*n*f1(p)/8000) + sin(2*pi*n*f2(q)/8000);% 构成双频信号 sound(x,8000); % 发出声音 pause(0.1)% 接收检测端的程序X=goertzel(x(1:205),K+1); % 用Goertzel 算法计算八点DFT 样本val = abs(X); % 列出八点DFT 向量 subplot(4,2,l);stem(K,val,'.');grid;xlabel('k');ylabel('|X(k)|') % 画出DFT(k)幅度axis([10 50 0 120])limit = 80; % for s=5:8;if val(s) > limit, break, end % 查找列号 endfor r=1:4;if val(r) > limit, break, end % 查找行号 endTNr=TNr+tm(r,s-4)*10^(8-l); enddisp('接收端检测到的号码为：') % 显示接收到的字符 disp(TNr) 显示结果：键入8位电话号码= 12345678 d = 1 d = 2 d = 3 d = 4 d = 5 d = 6 d = 7 d = 8接收端检测到的号码为：12345678|X (k )||X (k )||X (k )||X (k )||X (k )||X (k )|k|X (k )|k|X (k )|程序清单及波形显示： % 时域采样理论验证程序Tp=64/1000; %观察时间Tp=64微秒Fs=1000;T=1/Fs; M=Tp*Fs;n=0:M-1;t=n*T;A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;xat=A*exp(-alph*t).*sin(omega*t);Xk=T*fft(xat,M); %M点FFT[xat)]subplot(3,2,1); stem(n,xat,'.'); xlabel('n');ylabel('x1(n)'); title('(a) Fs=1000Hz');k=0:M-1;fk=k/Tp;subplot(3,2,2);plot(fk,abs(Xk));title('(a) T*FT[xa(nT)],Fs=1000Hz');xlabel('\omega/hz');ylabel('(H1(ejw))');axis([0,Fs,0,1.2*max(abs(Xk))]);Fs=300;T=1/Fs; M=Tp*Fs;n=0:M-1;t=n*T;A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;xat=A*exp(-alph*t).*sin(omega*t);Xk=T*fft(xat,M); %M点FFT[xat)]subplot(3,2,3); stem(n,xat,'.'); xlabel('n');ylabel('x2(n)'); title('(b)Fs=300Hz');k=0:M-1;fk=k/Tp;subplot(3,2,4);plot(fk,abs(Xk));title('(a) T*FT[xa(nT)],Fs=300Hz');xlabel('\omega/hz');ylabel('(H2(ejw))');axis([0,Fs,0,1.2*max(abs(Xk))]);Fs=200;T=1/Fs; M=Tp*Fs;n=0:M-1;t=n*T;A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;xat=A*exp(-alph*t).*sin(omega*t);Xk=T*fft(xat,M); %M点FFT[xat)]subplot(3,2,5); stem(n,xat,'.'); xlabel('n');ylabel('x3(n)'); title('(c) Fs=200Hz');k=0:M-1;fk=k/Tp;subplot(3,2,6);plot(fk,abs(Xk));title('(a) T*FT[xa(nT)],Fs=200Hz');xlabel('\omega/hz');ylabel('(H3(ejw))');axis([0,Fs,0,1.2*max(abs(Xk))])nx 1(n )(a) Fs=1000Hz(a) T*FT[xa(nT)],Fs=1000Hzω/hz(H 1(e j w ))nx 2(n )(b) Fs=300Hz(a) T*FT[xa(nT)],Fs=300Hzω/hz(H 2(e j w ))nx 3(n )(c) Fs=200Hz(a) T*FT[xa(nT)],Fs=200Hzω/hz(H 3(e j w ))%频域采样理论验证程序 clc;clear;close all; M=27;N=32;n=0:M; xa=0:(M/2); xb= ceil(M/2)-1:-1:0; xn=[xa,xb]; %产生M 长三角波序列x(n) Xk=fft(xn,1024); %1024点FFT[x(n)], 用于近似序列x(n)的TF X32k=fft(xn,32) ;%32点FFT[x(n)] x32n=ifft(X32k); %32点IFFT[X32(k)]得到x32(n) X16k=X32k(1:2:N); %隔点抽取X32k 得到X16(K) x16n=ifft(X16k,N/2); %16点IFFT[X16(k)]得到x16(n) subplot(3,2,2);stem(n,xn,'.'); title('(b) 三角波序列x(n)');xlabel('n');ylabel('x(n)');axis([0,32,0,20]) k=0:1023;wk=2*k/1024; subplot(3,2,1);plot(wk,abs(Xk));title('(a)FT[x(n)]'); xlabel('\omega/\pi');ylabel('|X(e^j^\omega)|');axis([0,1,0,200]) k=0:N/2-1; subplot(3,2,3);stem(k,abs(X16k),'.'); title('(c) 16点频域采样');xlabel('k');ylabel('|X_1_6(k)|');axis([0,8,0,200]) n1=0:N/2-1; subplot(3,2,4);stem(n1,x16n,'.'); title('(d) 16点IDFT[X_1_6(k)]');xlabel('n');ylabel('x_1_6(n)');axis([0,32,0,20]); k=0:N-1; subplot(3,2,5);stem(k,abs(X32k),'.'); title('(e) 32点频域采样');xlabel('k');ylabel('|X_3_2(k)|');axis([0,16,0,200]);n1=0:N-1; subplot(3,2,6);stem(n1,x32n,'.');box on title('(f) 32点IDFT[X_3_2(k)]');xlabel('n');ylabel('x_3_2(n)');axis([0,32,0,20])(b) 三角波序列x(n)nx (n )0100200(a)FT[x(n)]ω/π|X (e j ω)|(c) 16点频域采样k|X 16(k)|102030(d) 16点IDFT[X 16(k)]nx 16(n)(e) 32点频域采样k|X 32(k )|(f) 32点IDFT[X 32(k)]nx 32(n )程序清单及波形显示：% 用FFT 对信号作频谱分析 clear all;close all %实验内容(1)=================================================== x1n=[ones(1,4)]; %产生序列向量x1(n)=R4(n) M=8;xa=1:(M/2); xb=(M/2):-1:1; x2n=[xa,xb]; %产生长度为8的三角波序列x2(n) x3n=[xb,xa]; X1k8=fft(x1n,8); %计算x1n 的8点DFT X1k16=fft(x1n,16); %计算x1n 的16点DFT X2k8=fft(x2n,8); %计算x1n 的8点DFT X2k16=fft(x2n,16); %计算x1n 的16点DFT X3k8=fft(x3n,8); %计算x1n 的8点DFT X3k16=fft(x3n,16); %计算x1n 的16点DFT %以下绘制幅频特性曲线 subplot(1,2,1);stem(X1k8,'.'); %绘制8点DFT 的幅频特性图 title('(1a) 8点DFT[x_1(n)]');xlabel('ω/π');ylabel('幅度'); subplot(1,2,2);stem(X1k16,'.'); %绘制16点DFT 的幅频特性图 title('(1b)16点DFT[x_1(n)]');xlabel('ω/π');ylabel('幅度');figure(2)subplot(2,2,1);stem(X2k8,'.'); %绘制8点DFT 的幅频特性图 title('(2a) 8点DFT[x_2(n)]');xlabel('ω/π');ylabel('幅度');subplot(2,2,2);stem(X2k16,'.'); %绘制16点DFT 的幅频特性图 title('(2b)16点DFT[x_2(n)]');xlabel('ω/π');ylabel('幅度'); subplot(2,2,3);stem(X3k8,'.'); %绘制8点DFT 的幅频特性图 title('(3a) 8点DFT[x_3(n)]');xlabel('ω/π');ylabel('幅度');subplot(2,2,4);stem(X3k16,'.'); %绘制16点DFT 的幅频特性图 title('(3b)16点DFT[x_3(n)]');xlabel('ω/π');ylabel('幅度');-2.5-2-1.5-1-0.500.511.522.5(1a) 8点DFT[x 1(n)]ω/π幅度-2.5-2-1.5-1-0.500.511.522.5(1b)16点DFT[x 1(n)]ω/π幅度2468-4-2024(2a) 8点DFT[x 2(n)]ω/π幅度5101520-20-1001020(2b)16点DFT[x 2(n)]ω/π幅度2468-4-2024(3a) 8点DFT[x 3(n)]ω/π幅度-10-50510(3b)16点DFT[x 3(n)]ω/π幅度%实验内容(2) 周期序列谱分析==================================N=8;n=0:N-1; %FFT 的变换区间N=8 x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k8=fft(x4n); %计算x4n 的8点DFT X5k8=fft(x5n); %计算x5n 的8点DFT N=16;n=0:N-1; %FFT 的变换区间N=16 x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k16=fft(x4n); %计算x4n 的16点DFT X5k16=fft(x5n); %计算x5n 的16点DFT figure(3)subplot(2,2,1);stem(X4k8,'.'); %绘制8点DFT 的幅频特性图 title('(4a) 8点DFT[x_4(n)]');xlabel('ω/π');ylabel('幅度');subplot(2,2,3);stem(X4k16,'.'); %绘制16点DFT 的幅频特性图 title('(4b)16点DFT[x_4(n)]');xlabel('ω/π');ylabel('幅度');subplot(2,2,2);stem(X5k8,'.'); %绘制8点DFT 的幅频特性图 title('(5a) 8点DFT[x_5(n)]');xlabel('ω/π');ylabel('幅度');subplot(2,2,4);stem(X5k16,'.'); %绘制16点DFT 的幅频特性图title('(5b)16点DFT[x_5(n)]');xlabel('ω/π');ylabel('幅度');2468-1012x 10-15(4a) 8点DFT[x 4(n)]ω/π幅度5101520-4-2024x 10-15(4b)16点DFT[x 4(n)]ω/π幅度2468-4-2024(5a) 8点DFT[x 5(n)]ω/π幅度5101520-4-2024x 10-15(5b)16点DFT[x 5(n)]ω/π幅度%实验内容(3) 模拟周期信号谱分析=============================== figure(4) Fs=64;T=1/Fs; N=16;n=0:N-1; %FFT 的变换区间N=16 x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %对x6(t)16点采样 X6k16=fft(x6nT); %计算x6nT 的16点DFT X6k16=fftshift(X6k16); %将零频率移到频谱中心 Tp=N*T;F=1/Tp; %频率分辨率F k=-N/2:N/2-1;fk=k*F; %产生16点DFT 对应的采样点频率（以零频率为中心） subplot(3,1,1);stem(fk,abs(X6k16),'.');box on %绘制8点DFT 的幅频特性图 title('(6a) 16点|DFT[x_6(nT)]|');xlabel('f(Hz)');ylabel('幅度'); axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k16))])N=32;n=0:N-1; %FFT的变换区间N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %对x6(t)32点采样X6k32=fft(x6nT); %计算x6nT的32点DFTX6k32=fftshift(X6k32); %将零频率移到频谱中心Tp=N*T;F=1/Tp; %频率分辨率Fk=-N/2:N/2-1;fk=k*F; %产生16点DFT对应的采样点频率（以零频率为中心）subplot(3,1,2);stem(fk,abs(X6k32),'.');box on %绘制8点DFT的幅频特性图title('(6b) 32点|DFT[x_6(nT)]|');xlabel('f(Hz)');ylabel('幅度');axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k32))])N=64;n=0:N-1; %FFT的变换区间N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %对x6(t)64点采样X6k64=fft(x6nT); %计算x6nT的64点DFTX6k64=fftshift(X6k64); %将零频率移到频谱中心Tp=N*T;F=1/Tp; %频率分辨率Fk=-N/2:N/2-1;fk=k*F; %产生16点DFT对应的采样点频率（以零频率为中心）subplot(3,1,3);stem(fk,abs(X6k64),'.'); box on%绘制8点DFT的幅频特性图title('(6a) 64点|DFT[x_6(nT)]|');xlabel('f(Hz)');ylabel('幅度');axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k64))])-30-20-1001020300510(6a) 16点|DFT[x 6(nT)]|f(Hz)幅度010(6b) 32点|DFT[x 6(nT)]|f(Hz)幅度020(6a) 64点|DFT[x 6(nT)]|f(Hz)幅度程序清单及波形显示： clc;clear all;close all;fc1=250; fm1=15; fc2=500; fm2=50; fc3=1000; fm3=100;N=800; Fs=10000;Ts=1/Fs; n=[0:N-1];t=n*Ts;x11=cos(2*pi*fc1*t); x12=cos(2*pi*fm1*t); x1=x11.*x12; subplot(3,1,1);plot(t,x11,'g');plot(t,x12,'r');plot(t,x1,'b'); legend('载波','调制波 ','已调 ');xlabel('t/s');ylabel('波形')x=cos(2*pi*fc1*t).*cos(2*pi*fm1*t)+cos(2*pi*fc2*t).*cos(2*pi*fm2*t)+cos(2*pi*fc3*t).*cos(2*pi*fm3*t); subplot(3,1,2);plot(t,x);X=fft(x)subplot(3,1,3)k=[0:(N-1)/2]stem(k*2/N,abs(X(k+1))/max(abs(X(k+1))),'.');axis([0,0.3,0,1]);xlabel('\omeg a/\pi');ylabel('幅度');wp=[0.04,0.06];ws=[0.03,0.07];rp=0.1;rs=60;[N1,wp]=ellipord(wp,ws,rp,rs);[B,A]=ellip(N1,rp,rs,wp);y1=filter(B,A,x);figuresubplot(3,1,1);plot(t,x);xlabel('t/s');title('3路混合信号波形')subplot(3,1,2);[H1,w1]=freqz(B,A,N);plot(w1/pi,20*log10(abs(H1)));axis([0,0.5,-80,1]);xlabel('\omega/\pi');ylabel('|H(e^j\omega)|');title('中心频率为250Hz的频率响应');subplot(3,1,3);plot(t,y1);xlabel('t/s');ylabel('y1(t)');title('中心频率为250H的滤波信号')wp=[0.08,0.12];ws=[0.07,0.13];rp=0.1;rs=60;[N1,wp]=ellipord(wp,ws,rp,rs);[B,A]=ellip(N1,rp,rs,wp);y1=filter(B,A,x);figuresubplot(3,1,1);plot(t,x);xlabel('t/s');title('3路混合信号波形');subplot(3,1,2)[H1,w1]=freqz(B,A,N);plot(w1/pi,20*log10(abs(H1)));axis([0,0.5,-90,2]);xlabel('\omega/\pi');ylabel('|H(e^j\omega)|');title('中心频率为500Hz的频率响应')subplot(3,1,3);plot(t,y1);xlabel('t/s');ylabel('y1(t)');title('中心频率为500H的滤波信号')wp=[0.17,0.23];ws=[0.16,0.24];rp=0.1;rs=60;[N1,wp]=ellipord(wp,ws,rp,rs);[B,A]=ellip(N1,rp,rs,wp);y1=filter(B,A,x);figuresubplot(3,1,1);plot(t,x);xlabel('t/s');title('3路混合信号波形');subplot(3,1,2);[H1,w1]=freqz(B,A,N);plot(w1/pi,20*log10(abs(H1)));axis([0,0.5,-100,10]);xlabel('\omega/\pi');ylabel('|H(e^j\omega)|');title('中心频率为1000Hz的频率响应')subplot(3,1,3);plot(t,y1);xlabel('t/s');ylabel('y1(t)');title('中心频率为1000H的滤波信号')00.010.020.030.040.050.060.070.08-101t/s波形00.010.020.030.040.050.060.070.08-505ω/π幅度-505t/s3路混合信号波形-80-60-40-200ω/π|H (e j ω)|中心频率为250Hz 的频率响应-202t/sy 1(t )中心频率为250H 的滤波信号0.010.020.030.040.050.060.070.08-505t/s3路混合信号波形0.050.10.150.20.250.30.350.40.450.5-80-60-40-200ω/π|H (e j ω)|中心频率为500Hz 的频率响应0.010.020.030.040.050.060.070.08-202t/sy 1(t )中心频率为500H 的滤波信号0.010.020.030.040.050.060.070.08-505t/s3路混合信号波形0.050.10.150.20.250.30.350.40.450.5-100-500ω/π|H (e j ω)|中心频率为1000Hz 的频率响应0.010.020.030.040.050.060.070.08-202t/sy 1(t )中心频率为1000H 的滤波信号程序清单及波形显示： clc;clear;clear allN=1000; Fs=1000;T=1/Fs;Tp=N*T; t=0:T:(N-1)*T;fc=Fs/10;f0=fc/10; %载波频率fc=Fs/10，单频调制信号频率为f0=Fc/10; mt=cos(2*pi*f0*t); %产生单频正弦波调制信号mt ，频率为f0 ct=cos(2*pi*fc*t); %产生载波正弦波信号ct ，频率为fc xt=mt.*ct; %相乘产生单频调制信号xt nt=2*rand(1,N)-1; %产生随机噪声ntfp=150; fs=200;Rp=0.1;As=70; % 滤波器指标fb=[fp,fs];m=[0,1]; % 计算remezord 函数所需参数f,m,devdev=[10^(-As/20),(10^(Rp/20)-1)/(10^(Rp/20)+1)];[n,fo,mo,W]=remezord(fb,m,dev,Fs); % 确定remez 函数所需参数hn=remez(n,fo,mo,W); % 调用remez 函数进行设计,用于滤除噪声nt 中的低频成分yt=filter(hn,1,10*nt); %滤除随机噪声中低频成分，生成高通噪声yt xt=xt+yt; %噪声加信号 fst=fft(xt,N);k=0:N-1;f=k/Tp;subplot(2,1,1);plot(t,xt);grid;xlabel('t/s');ylabel('x(t)');axis([0,Tp/5,min(xt),max(xt)]);title('(a) 信号加噪声波形')subplot(2,1,2);plot(f,abs(fst)/max(abs(fst)));grid;title('(b) 信号加噪声的频谱')axis([0,Fs/2,0,1.2]);xlabel('f/Hz');ylabel('幅度')-10-50510t/sx (t )(a) 信号加噪声波形0501001502002503003504004505000.51(b) 信号加噪声的频谱f/Hz幅度%==调用xtg 产生信号xt, xt 长度N=1000,并显示xt 及其频谱,========= N=1000;xt=xtg(N);fp=120; fs=150;Rp=0.2;As=60;Fs=1000; % 输入给定指标% (1) 用窗函数法设计滤波器wc=(fp+fs)/Fs; %理想低通滤波器截止频率(关于pi归一化）B=2*pi*(fs-fp)/Fs; %过渡带宽度指标Nb=ceil(11*pi/B); %blackman窗的长度Nhn=fir1(Nb-1,wc,blackman(Nb));Hw=abs(fft(hn,1024)); % 求设计的滤波器频率特性ywt=fftfilt(hn,xt,N); %调用函数fftfilt对xt滤波subplot(4,1,1);plot(f,Hw);xlabel('f/Hz');ylabel('幅度'); title('(a)低通滤波器幅频特性');subplot(4,1,2);plot(t,ywt); title('(b)滤除噪声后的信号波形');xlabel('t/s');ylabel('ywt');% (2) 用等波纹最佳逼近法设计滤波器fb=[fp,fs];m=[1,0]; % 确定remezord函数所需参数f,m,dev dev=[(10^(Rp/20)-1)/(10^(Rp/20)+1),10^(-As/20)];[Ne,fo,mo,W]=remezord(fb,m,dev,Fs); % 确定remez函数所需参数hn=remez(Ne,fo,mo,W); % 调用remez函数进行设计Hw=abs(fft(hn,1024)); % 求设计的滤波器频率特性yet=fftfilt(hn,xt,N); % 调用函数fftfilt对xt滤波subplot(3,1,1);plot(f,Hw);xlabel('f/Hz');ylabel('幅度'); title('(c)低通滤波器幅频特性');subplot(3,1,2);plot(t,yet);title('(d)滤除噪声后的信号波形');xlabel('t/s');ylabel('yet');。