西安交通大学数字信号处理实验报告_频率采样型滤波器Word文档下载推荐.doc
《西安交通大学数字信号处理实验报告_频率采样型滤波器Word文档下载推荐.doc》由会员分享,可在线阅读,更多相关《西安交通大学数字信号处理实验报告_频率采样型滤波器Word文档下载推荐.doc(16页珍藏版)》请在冰点文库上搜索。
1.构造滤波器输入信号st=k=03sk(t),其中skt=Akcos(2πkf0t+ϕk),基波频率f0=50Hz,A0=0.5,A1=1,A2=0.5,A3=2,ϕ0=0,ϕ1=π2,ϕ0=π,ϕ3=-π2。
设时域信号s(t)的采样频率fs=Nf0,绘制出采样时刻从0到L-1的采样信号波形,其中采样点数为L=2N,确认时域信号采样正确。
2.对采样信号的第二个周期(n=N,N+1,⋯,L-1)进行离散傅里叶变换,画出幅频特性和相频特性图,观察并分析其特点。
3.设H(0)=1,H
(1)=e-jπ(N-1)N,H
(2)=e-j2π(N-1)N,H3=H4=⋯=H(13)=0,H14=-e-j14π(N-1)N,H(15)=-e-j15π(N-1)N,计算滤波器抽头系数h(n),n=0,1,…,N-1,画出该滤波器的频谱图,观察并分析其幅频特性和相频特性。
4.编程实现图1所示的频率采样型滤波器结构,其中r=0.999,H(k)取第3步中的值。
为了简化编程,梳妆滤波器可以调用CombFilter.m,谐振器可以调用Resonator2.m,使用helpCombFilter和helpResonator2查看如何配置参数。
将第1步生成的采样信号通过该滤波器,画出输出信号第二个周期(n=N,N+1,…,L-1)的时域波形和频谱,并与第2步的频谱进行对比,观察并分析二者的区别。
5.分别画出图1中前4路谐振器的输出信号第二个周期(n=N,N+1,…,L-1)的时域波形,观察并分析输出信号的特点。
6.将输入信号换成周期为N的冲激串,画出输出信号第二个周期(n=N,N+1,…,L-1)的幅频特性,并与第3步的滤波器幅频特性进行对比,观察并分析二者的关系。
7.思考题
(1)在第2步的幅频特性中,各次谐波的幅度与相应的时域信号幅度有什么关系?
(2)实验中为什么要观察第二个周期,如果直接观察第一个周期会怎么样?
(3)如果取r=0.95,观察会出现什么情况。
三、实验要求
1.按照实验内容编写Matlab程序,给出运行结果(图),并逐项进行分析讨论。
2.提交完整的Matlab源程序。
3.回答思考题,撰写实验报告。
四、实验结果与分析
输入信号是直流、一次波、二次波、三次波的叠加,采样频率为16f0,高于输入信号的2倍。
输入信号的时域图像如上。
由幅频特性图可知,在第0、1、2、3、13、14、15个点上有幅值,数值对称、相位反对称,其余点的幅度为0。
第0点表示基波分量;
第1点和第15点为直流分量,其幅值和是基波分量的2倍;
第2点和第14点为二次分量,其幅值和与基波分量相等;
第3点和第13点为三次分量,幅值和是基波分量的4倍。
同时,从FFT结果可发现,幅值为0的分量相位并不相同,这是由计算机计算精度造成的。
滤波器抽头系数h(n)可由H(k)作IFFT运算得出,第一行两幅图表示出H(k)的幅频特性和相频特性。
将H(k)内插后观察,发现H2为低通滤波器。
第四题输出图像:
信号通过滤波器后第二个周期的时域图像、FFT后的频谱
信号直接作FFT后的频谱
幅值:
与第2题对比,第0、1、2、14、15点不变,第3点和第13点的数值为0,这是因为三次分量在通过滤波器时被滤掉了,而滤波器在第0、1、2、14、15点的幅值为1,所以通过滤波器后的幅值没有发生变化。
相位:
直流分量的初始相位与滤波器H(0)相乘,得到新的相位;
一次分量的相位由基波初始相位0与H
(1)相位相乘而得;
其他分量的相位无法确定。
H(0)输出波形为直流分量,所以波形的幅度不变;
因为H(3)=0,所以H(3)输出波形也为0,幅度不变;
其他两路分别为基波和二次分量的波形。
输入信号是周期为N的冲激串时,输出信号的时域图像和频谱
滤波器幅频特性
冲激函数在频域上为1,通过滤波器后被滤掉k=3~13区域,所以其幅值和滤波器本身的幅值相等,保留了直流、基波、二次分量。
五、思考题
1.在第2步的幅频特性中,各次谐波的幅度与相应的时域信号幅度有什么关系?
频率为0时的幅度等于时域波形中直流的幅度,直流、二次、三次谐波的谱线的幅度为时域波形相应谐波对应幅度的一半。
因为正弦信号的频谱中有和两条谱线,幅度都为信号的一半,cos(t)=ejw+e-jw2。
2.实验中为什么要观察第二个周期,如果直接观察第一个周期会怎么样?
第二个周期信号与系统的单位脉冲响应卷积后,输出信号开始变得周期,不会发生失真;
直接观察第一个周期,信号没有完全进入系统,卷积结果不完整,输出会发生失真。
3.如果取r=0.95,观察会出现什么情况。
所有谐振器的极点从单位圆向内收缩一点,同时梳状滤波器的零点也移到r圆上,滤波器的幅频特性的幅度减小,通带内的仍然保持线性相位特性。
4.如何理解第3步与第6步在工程使用中的区别?
两者均为低通滤波器,但是第三步是进行频率采样后由内插函数得到传递函数H(z),从而得到系统频响H();
而题6中使用了冲激串,这是因为在理想条件下,冲激串的频响就是系统的频响,但是在实际生活中很难实现。
六、MATLAB程序代码
1.以下是MATLAB程序代码:
clc;
clear;
clf;
f0=50;
N=16;
fs=N*f0;
%时域信号s(t)的采样频率
T=1/fs;
%采样间隔
L=N*2;
%采样点数
t=0:
T:
(L-1)*T;
%构造输入信号
s0=0.5*cos(0);
s1=1*cos(2*pi*f0*t+pi/2);
s2=0.5*cos(2*pi*2*f0*t+pi);
s3=2*cos(2*pi*3*f0*t-pi/2);
s=s0+s1+s2+s3;
i=0:
1:
L-1;
stem(i,s(i+1));
gridon;
2.以下是MATLAB程序代码:
t=N*T:
s0=0.5*cos(2*pi*0*f0*t+0);
a=fft(s);
%FFT运算
a %输出a的各项值
y=abs(a);
z=angle(a);
N-1;
subplot(2,1,1);
stem(i,y(i+1));
subplot(2,1,2);
stem(i,z(i+1));
3.以下是MATLAB程序代码;
H=[1,exp(-1i*pi*(N-1)/N),exp(-j*2*pi*(N-1)/N),0,0,0,0,0,0,0,0,0,0,0,-exp(-j*14*pi*(N-1)/N),-exp(-j*15*pi*(N-1)/N)];
h=ifft(H,N);
h
H_abs=abs(H);
H_angle=angle(H);
subplot(2,2,1);
stem(i,H_abs(i+1));
title('
H幅频特性'
);
subplot(2,2,2);
stem(i,H_angle(i+1));
H相频特性'
w=0:
pi/200:
2*pi;
H2=zeros(1,length(w));
fori=1:
length(w)
forn=1:
N
e=exp(-j*w(i)*(n-1));
H2(i)=H2(i)+h(n)*e;
end
end
subplot(2,2,3);
plot(abs(H2));
H2幅频特性'
subplot(2,2,4);
plot(angle(H2));
H2相频特性'
4.以下是MATLAB程序代码;
%初始化数据
L=2*N;
s1=1*cos(2*pi*1*f0*t+pi/2);
x=s;
r=0.999;
y=CombFilter(x,N,r);
x=y;
z=zeros(1,48);
fori=0:
(N/2)
Order=i;
H=[1exp(-j*pi*(N-1)/N)exp((-j*2*pi*(N-1)/N))00000000000-exp((-j*14*pi*(N-1)/N))-exp((-j*15*pi*(N-1)/N))];
y=Resonator2(x,N,r,Order,H(i+1));
z=z+y;
y=z/N;
fori=N:
y1(i-N+1)=y(i+1);
%取出N到L-1的元素
%作时域图
subplot(2,2,1),stem(i,y1(i+1));
xlabel('
n'
ylabel('
振幅'
时域波形'
%作频谱
y=fft(y1);
y
y1=abs(y);
ang=angle(y);
subplot(2,2,3),stem(i,y1(i+1));
幅值'
幅频特性'
subplot(2,2,4),stem(i,ang(i+1));
相频特性'
5.以下是MATLAB程序代码;
(L-1)*T
%生成信号
%前4路谐振器的输出
y1=CombFilter(s,N,r);
H=[1,exp(-1i*pi*(N-1)/N),exp(-1i*2*pi*(N-1)/N),0,0,0,0,0,0,0,0,0,0,0,exp(-1i*14*pi*(N-1)/N),exp(-1i*15*pi*(N-1)/N)];
y_h0=Resonator2(y1,N,r,0,H
(1));
y_h0=y_h0/N;
y_h0=y_h0(N:
2*N-1);
y_h1=Resonator2(y1,N,r,1,H
(2));
y_h1=y_h1/N;
y_h1=y_h1(N:
y_h2=Resonator2(y1,N,r,2,H(3));
y_h2=y_h2/N;
y_h2=y_h2(N:
y_h3=Resonator2(y1,N,r,3,H(4));
y_h3=y_h3/N;
y_h3=y_h3(N:
%输出图像
k=0:
stem(k,y_h0);
H(0)路输出信号第二个周期'
stem(k,y_h1);
H
(1)路输出信号第二个周期'
stem(k,y_h2);
H
(2)路输出信号第二个周期'
stem(k,y_h3);
H(3)路输出信号第二个周期'
6.以下是MATLAB程序代码;
x=zeros(1,2*N);
x
(1)=1;
x(N+1)=1;
%冲激串
H=zeros(1,N);
h=zeros(1,N);
H=[1,exp(-j*pi*(N-1)/N),exp((-j*2*pi*(N-1)/N)),0,0,0,0,0,0,0,0,0,0,0,-exp((-j*14*pi*(N-1)/N)),-exp((-j*15*pi*(N-1)/N))];
%冲激串滤波后
y1=zeros;
fork=0:
N/2;
y1=y1+Resonator2(y,N,r,k,H(k+1));
end
y0=y1/N;
%经谐振器后
k=N:
2*N-1;
holdon;
stem(k,y0(k+1));
输出信号时域'
Y=zeros(1,N);
fork=N:
Y(k-N+1)=y0(k+1);
%输出信号第二个周期
F=fft(Y);
F
holdon
subplot(2,2,2)
stem(k,abs(F(k+1)))
FFT后幅频特性'
subplot(2,2,3)
stem(k,angle(F(k+1)))
FFT后相频特性'
15