数字信号处理实验报告二与三

------------------------------------------------------------------------------------------------

—————————————————————————————————————— 数字信号处理实验报告二与三 实验二 用FFT进行谱分析 一.实验目的： 1 进一步加深DFT算法原理和基本性质的理解（因为FFT只是DFT的一种快速算法，所以FFT的运算结果必然满足DFT的基本性质）。熟悉FFT程序结构及编程方法。 2 熟悉应用FFT对确定信号进行谱分析方法，熟悉FFT算法原理和FFT子程序的应用。 3 学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析误差及其原因，以便在实际中正确应有FFT。 二.实验内容： （1）用matlab编程产生并画出信号x1(n)、x2(n)、x3(n)、x4(n)、x5(n)。 x(n)?R(n)14 ?n?1,0?n?3 ? x2(n)??8?n,4?n?7 ?0,其它n? ?4 ?n,0?n?3 ? x(n)?n?3,4?n?7?3 ? ,其它n?0（ 2）用matlab编制FFT函数对上述信号进行频谱分析，并画出上------------------------------------------------------------------------------------------------

—————————————————————————————————————— 述信号谱图。 三.实验结果（1） 1.%This programm is to generate signal x1(n)=R4(n). k=-6:6; x=[zeros(1,6),ones(1,4),zeros(1,3)]; stem(k,x); （信号图如图1） title('图1'); 2.n=-5:1:10; x=(n+1).*(n>=0 & n<=3)+(8-n).*(n>=4 & n<=7)+0; stem(n,x)； title('图2'); 3.n=-5:10; x=(4-n).*(n>=0 & n<=3)+(n-3).*(n>=4 & n<=7); stem(n,x); title('图3'); 4.n=-10:10; x=cos(pi/4*n); stem(n,x); title('图4'); 5.n=-10:10; x=sin(pi/8*n); stem(n,x); title('图5'); 实验结果（2）: FFT算法 function y=myditfft(x) % y=myditfft(x) % 本程序对输入序列 x 实现DIT-FFT基2算法，点数取大于等于x长度的2的幂次 % x为给定时间序列 % y为x的离散傅立叶变换 m=nextpow2(x);N=2^m; % 求x的长度对应的2的最低幂次m if length(x)<N; ------------------------------------------------------------------------------------------------

—————————————————————————————————————— % 若x的长度不是2的幂，补零到2的整数幂 x=[x,zeros(1,N-length(x))]; end nxd=bin2dec(fliplr(dec2bin([1:N]-1,m)))+1; % 求1:2^m数列的倒序 y=x(nxd); % 将x倒序排列作为y的初始值 for mm=1:m; % 将DFT作m次基2分解,从左到右，对每次分解作DFT运算 Nmr=2^mm;u=1; % 旋转因子u初始化为WN^0=1 WN=exp(-i*2*pi/Nmr); % 本次分解的基本DFT因子WN=exp(-i*2*pi/Nmr) for j=1:Nmr/2; % 本次跨越间隔内的各次蝶形运算 for k=j:Nmr:N; % 本次蝶形运算的跨越间隔为Nmr=2^mm kp=k+Nmr/2; % 确定蝶形运算的对应单元下标 t=y(kp)*u; % 蝶形运算的乘积项 y(kp)=y(k)-t; % 蝶形运算 y(k)=y(k)+t; % 蝶形运算 end u=u*WN; % 修改旋转因子,多乘一个基本DFT因子WN end end 1. k=-6:6; x=[zeros(1,6),ones(1,4),zeros(1,3)]; y=myditfft(x); k=-6:9; stem(k,y); xlabel('m'); ylabel('X[M]'); X[M ------------------------------------------------------------------------------------------------

—————————————————————————————————————— ] title('FFT图'); 2. n=-5:1:10; x=(n+1).*(n>=0 & n<=3)+(8-n).*(n>=4 & n<=7)+0; y=myditfft(x); stem(n,y); xlabel('n'); ylabel('X[M]'); title('FFT图'); 3. n=-5:10; x=(4-n).*(n>=0 & n<=3)+(n-3).*(n>=4 & n<=7); y=myditfft(x); stem(n,y); xlabel('n'); ylabel('X[M]'); title('FFT3'); 4. n=-10:10; x=cos(pi/4*n); y=myditfft(x); n=-10:21; stem(n,y); xlabel('n'); ylabel('X[M]'); title('FFT4'); 5. n=-10:10; x=sin(pi/8*n); y=myditfft(x); n=-10:21; stem(n,y); xlabel('n'); ylabel('X[M]'); title('FFT5'); X[M] m n ------------------------------------------------------------------------------------------------

—————————————————————————————————————— XMn FFT4 X[M] nX[M] n 四.简要回答以下问题： ①在N=8时，x2(n)和x3(n)的幅频特性会相同吗？为什么？N=16呢？ 答：不相同。当N=8时，序列x1(n)和x2(n)中相同的元素值对应的 n 值是不同的，所乘的旋转因子的值也不同，因而得到的最终结果也是不同的。 同理，N=16时，所得的幅频特性也是不同的。 ②FFT在什么条件下也可以用来分析周期信号序列的频谱？如果正弦信号系统sin(2πf0k),f0=0.1Hz，用16点FFT来做DFT运算，得到的频谱是信号本身的真实谱吗？为什么？ 答：由于FFT算法对序列长度的要求是 N=2^M，M为正整数。所以，当周期信号序列一个周期的长度满足 N=2^M（M为正整数）的条件时，FFT可以用来分析周期信号的频谱。不是真实的频谱。因为序列的周期N=10不是 2 的整数次幂，所以不是真实的。 实验三：用双线性变换法设计IIR数字滤波器 一.实验目的： 1. 掌握用双线性变换法设计IIR DF的原理及具体设计方法，熟悉------------------------------------------------------------------------------------------------

—————————————————————————————————————— 用双线性变换法设计IIR DF的计算机编程。 2. 观察用双线性变换法设计的DF的频响特性，了解双线性变换法的特点。 3. 熟悉用双线性变换法设计BW和CB型DF的全过程。 二，实验原理： 为了克服冲激响应不变法产生的频率混叠现象，这是从S平面到Z平面的标准变换z＝esT的多值对应关系导致的，为了克服这一缺点，产生了双线性变换法。 12，阻带截频fS=3400Hz,衰耗ASmin=31dB,取样频率fs=8000Hz。画出其模拟滤波器及数字滤波器频响特性曲线。 程序：FS=8000; Fl=3000;Fh=3400; %通带、阻带截止频率 Rp=3;Rs=31; wp=Fl*2*pi; %临界频率采用角频率表示 ws=Fh*2*pi; %临界频率采用角频率表示 wp1=wp/FS; %求数字频率 ws1=ws/FS; %求数字频率 OmegaP=2*FS*tan(wp1/2);%频率预畸 OmegaS=2*FS*tan(ws1/2);%频率预畸 %选择滤波器的最小阶数 [n,Wn]=cheb1ord(OmegaP,OmegaS,Rp,Rs,'s'); %此处是代入经预畸变后获得的归一化模拟频率参数 ------------------------------------------------------------------------------------------------

—————————————————————————————————————— [bt,at]=cheby1(n,Rp,Wn,'s'); % 设计一个n阶的契比雪夫I型模拟滤波器 subplot(2,1,1); [h,w]=freqz(bt,at); plot(w*FS/(2*pi),abs(h)); grid; xlabel('频率/Hz');ylabel('幅值'); title('模拟滤波器频响特性'); [bz,az]=bilinear(bt,at,FS); %双线性变换为数字滤波器 subplot(2,1,2); [H,W] = freqz(bz,az); %求解数字滤波器的频率响应 plot(W*FS/(2*pi),abs(H));grid; xlabel('频率/Hz');ylabel('幅值'); title('数字滤波器频响特性'); 模拟滤波器频响特性 0.708 0.7079 0.707905001000150020002500 频率/Hz 数字滤波器频响特性 1 0.5 幅值幅值300035004000 3 ，阻带截频0050010001500fs=2000Hz2000。画出其模拟滤波器及数字滤波器频响特2500300035004000fS=400Hz,衰耗ASmin=15dB,取