FIR低通滤波器的设计

1 技术要求

用窗函数法设计线性相位FIR低通滤波器。要求通带截止频率ωc=π/4，单位脉冲响应h(n)的长度N=23。绘出h(n)及其幅频响应特性曲线。

2 基本原理

2.1 FIR低通滤波器 

数字滤波器是数字信号处理学科的重要组成部分，应用非常广泛。数字滤波器，通常是一种算法，或是一种数字处理设备。它的功能是将一组输入的数字序列经过一定的运算后变成为另一组输出的数字序列。它的主要功能是对数字信号进行处理，保留数字信号中的有用成分，去除信号中的无用成分。

数字滤波器是在模拟滤波器的基础上发展起来的，但它们之间存在着一些重要差别。与模拟滤波器相比，数字滤波器具有精度高、稳定性好、设计灵活等优点。一般情况下数字滤波器是一个线性非移变系统。与模拟滤波器相同的是，数字滤波器也有低通、高通、带通、带阻之分。本次实验中，就是要设计数字低通滤波器。

低通滤波器是容许低于截止频率的信号通过， 但高于截止频率的信号不能通过的滤波器。

从结构上，数字滤波器可以分为递归型（IIR）数字滤波器和非递归型（FIR）数字滤波器。本次实验要求利用FIR设计线性相位的低通滤波器。FIR最大的特点之一就是能够做成严格的线性相位关系。所为线性相位，就是指滤波器对不同频率的正弦波产生的相位延迟与正弦波的频率呈线性关系。因而，在通过该滤波器后在滤波器通带内的所有信号频率成分，除了由相频特性决定的延迟外，可以全部保留。

窗函数法设计低通滤波器，窗函数法也称为傅里叶级数法。从单位取样响应的观点来看，就是使设计的滤波器的h(n)逼近理想的滤波器单位取样相应hd(n)。

设计思想:

先给定理想filter的频响Hd()，设计一个FIR的filter的频响为H()，使H()逼近Hd()。

设计过程：

先用傅氏反变换求出理想filter的单位抽样响应hd(n)，然后加时间窗w(n)对hd(n)截断，以求得FIR filter的单位抽样响应h(n)。

由此来看，窗函数的形状及长度的选择就尤为关键。由于长度已给定，利用矩形窗、三角窗、汉明窗、汉宁窗、布莱克曼窗和凯泽窗几种窗函数对窗函数形状不同进行对比。

2.2 图形用户界面GUI

图形用户界面(GUI)是指采用图形方式显示的计算机操作用户界面。与早期计算机使用的命令界面相比，图形界面对于用户来说在视觉上更易于接受。本次课程设计就是应用matlab的GUI制作简单的用户界面，并且基于该界面仿真常见的几种信号的模型。GUI界面的引入，使用户的操作更加便捷，也使仿真的结果呈现起来较为简便。

本次实验中，由于对比窗函数形状对低通滤波器幅频响应的影响，总共会有8个图形的显示。运行后，由于图形窗口较多，观察比较不方便，所以决定用图形用户界面GUI来实现图形的显示。

3 建立模型描述

3.1单位冲激响应及其幅频响应模块

根据实验指导书要求，不仅要求画出最终低通滤波器的幅频响应，还要画出单位冲激响应h(n)的图形，还有单位冲激响应在频域的图形。对单位冲激序列做离散傅里叶变换，得到其幅频响应。其中，离散傅里叶变换要编写相关函数m文件。理想单位冲激响应的幅频响应的曲线形状为抽样信号的形状。

3.2 不同窗函数低通幅频响应对比模块

在加窗处理的环节，窗函数长度已经确定，选择何种类型的窗函数就是影响最终低通滤波器幅频响应的性能的关键。本实验中，未规定窗函数的类型，所以在程序中，我决定将不同用窗函数截断后的低通滤波器的形状比较作为一个拓展功能。

Matlab中提供的窗函数有：

(1)矩形窗

(2)三角窗

(3)汉明窗

(4)汉宁窗

(5)布莱克曼窗

(6)凯泽窗

各窗函数参数，旁瓣峰值幅度，过渡带宽，阻带最小衰减如图一。

图1窗函数参数

利用窗函数设计低通滤波器，只需利用函数fir1(N,wn,窗函数(N+1))，就可求得滤波器的幅频响应和相频响应。

汉明窗的应用方式为用fir1(N,wn)函数得到低通滤波器，用[w,t]=freqz(b,1,512)，abs(w)求得幅频响应。

矩形窗的应用方式为用fir1(N,wn，boxcar(n))函数得到低通滤波器，[w,t]=freqz(b,1,512)，abs(w)求得幅频响应。

三角窗的应用方式为用fir1(N,wn，triang(n) )函数得到低通滤波器，[w,t]=freqz(b,1,512)，abs(w)求得幅频响应。

汉宁窗的应用为用fir1(N,wn，hanning(n))函数得到低通滤波器，[w,t]= freqz (b,1,512)，abs(w)求得幅频响应。

布莱克曼窗的方式为fir1(N,wn，blackman(n))函数得到低通滤波器[w,t]=freqz(b,1,512)，abs(w)求得幅频响应。

凯泽窗应用方式为用fir1(N,wnkaiser(n))函数得到低通滤波器，用[w,t]=freqz(b,1,512)，abs(w)求得幅频响应。  

3.3 截止频率调整模块

本次试验中，规定截止频率为0.25π，为了能直观的显示不同的截止频率，低通滤波器的形状会有什么不同，我在GUI上添加了一个可编辑文本框，把默认数值设定为0.25，同时将程序中对应的截止频率从一个常数变为从文本框中取得的数。

3.4界面的美化

为了制作出友好美观的界面，必须对波形的显示区域、幅度、频率等进行调整，改变曲线的粗细程度，对各个控件及整体界面进行颜色的设置，并利用其工具栏将各个控件对齐。

3.5 功能模块图

功能模块图如图2。

图2功能模块图

4 模块功能分析或源程序代码

4.1单位冲激响应及其幅频响应及其代码

单位冲激响应的绘制思路比较简单，就是将一条直线用stem函数绘制出来，其图形必须为23个点。其代码如下：

n=0:22;

x=n./n;

stem(n,x); title('h(n)');

axis([0,25,0,1.3]);

单位冲激响应的幅频响应要用到专门的函数m文件。该函数文件可以在主程序中调用多次，节省篇幅。函数m文件代码如下：

function xk=dft(xn,N)

n=[0:1:N-1]; 

k=n;

WN=exp(-j*2*pi/N);

nk=n'*k;

WNnk=WN.^nk;

xk=xn* WNnk;

4.2 不同窗函数低通幅频响应对比模块

窗函数不同，会对设计出的低通滤波器的性能不同。应尽量选取旁瓣小主瓣窄的窗函数。为了细致观察窗函数不同对设计的低通滤波器的不同影响，本次试验中设计了不同窗函数的低通滤波器的比较环节。各个窗函数设计低通滤波器的思路是相似的，只是其中的窗函数是不同的。

用矩形窗设计低通滤波器的代码：

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn,boxcar(N+1));

[w,t]=freqz(b,1,512);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0);

title('矩形窗');

用汉明窗设计低通滤波器的代码：

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn);

[w,t]=freqz(b,1,512);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0)

title('汉明窗');

用三角窗设计低通滤波器的代码：

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn,triang(N+1));

[w,t]=freqz(b,1,512);

axis([0,3.5,0,1.5]);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0)

axis([0,3.5,0,1.4]);

title('三角窗');

用汉宁窗设计低通滤波器的代码：

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn,hanning(N+1));

[w,t]=freqz(b,1,512);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0)

title('汉宁窗 ');

用布莱克曼窗设计低通滤波器的代码：

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn,blackman(N+1));

[w,t]=freqz(b,1,512);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0)

title('布莱克曼窗');

用凯泽窗设计低通滤波器的代码：

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn,kaiser(N+1));

[w,t]=freqz(b,1,512);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0)

title('凯泽窗');

4.3 截止频率调整模块

为了研究低通滤波器的截止频率对滤波器的影响，本次实验设计了调整截止频率的模块，用以直观的显示不同的截止频率对低通滤波器的影响。为此，在控件中添加了一个可编辑文本框，用以输入要调整的截止频率。并且将从可编辑文本框获取的数值赋予D。但是要注意的是截止频率是个有实际意义的常数，不能够随意指定其数值。所以，在文本框内获取的是π的倍数。wc=0.25*π也变为wc=D*π。

另外，从文本框中获取的的值默认为字符型，但计算时用到的应为数值型，所以要进行类型的转换。其代码如下：

D=str2num(get(handles.edit1,'string'));

5 调试过程及结论

5.1 调试过程

在绘制h(n)时，本来打算绘制一条普通的直线，并用stem函数画出，但是直线中不含自变量，因此将其从y=1改为y=n./n。

在绘制h(n)的幅频响应时，我误认为dpt为matlab自带的函数，所以运行时报错，dpt未定义。因此查找资料，发现该函数为自己设定的函数文件，并且建立了该文件。

在绘制低通滤波器时，应用了freqz(b,1,512)函数，在m文件中运行正确。但在GUI中出现问题，后来发现，GUI中出现subplot函数后就会不能正常绘图。而freqz(b,1,512)函数，本身自带了同时绘制相频和幅频的命令，因此会出现错误。因此，我利用该函数求出幅频特性，但是利用plot绘制幅频特性。

在调整截止频率的模块，运行时报错，原因是忘记转换从文本框获取的数据类型。从文本框获取了字符型的数据，应当转换成字符型才能参与运算。

5.2 设计结果

下面将本次设计的结果按功能展示：

图1总体设计

图2单位冲激响应

图3单位冲激幅频响应

图4用矩形窗设计低通滤波器

图5用汉明窗设计低通滤波器

图6用三角窗设计低通滤波器

图7用汉宁窗设计低通滤波器

图8用布莱克曼窗设计低通滤波器

图9用凯泽窗设计低通滤波器

图10改变截止频率

6 思考题

线性相位满足的条件？线性相位FIR滤波器的特性及其应用领域？

答：(1)线性相位要求的条件是h(n)必须为实序列，且满足h(n)=h(N-n-1)。例如本设计中使用的单位冲激响应就符合实序列，且N=23，为奇数，h(n)=h(23-n-1)。

在此条件已满足的前提下共分四种情况。第一种是h(n)=h(N-n-1)，N为奇数；第二种是h(n)=h(N-n-1)，N为偶数；第三种是h(n)=-h(N-n-1)，N为奇数；第四种h(n)=-h(N-n-1)，N为偶数；由此易知本次设计所要求的为第一种情况。

(2)线性相位FIR滤波器的特性为：通过滤波器的信号中不同频率的正弦波产生的相位延迟与正弦波的频率呈线性关系。因而，在通过该滤波器后在滤波器通带内的所有信号频率成分，除了由相频特性决定的延迟外，可以全部保留。

正因为线性相位能够保留所有频率成分，且产生的相位延迟与频率呈线性关系，所以相位失真的产生是可以控制的。非线性相位滤波器由于产生的相位失真与频率不呈线性关系所以容易造成信号的色散。其应用领域为：语音信号、图像信号、视频信号要求线性相位，要求信道具有线性相位特性，因此该类信号在滤波、除噪时要用到FIR线性相位滤波器。

7 心得体会

本次课程设计要求在matlab软件中用窗函数法设计FIR线性相位低通滤波器。线性相位在实际之中应用的十分广泛。例如图像、语音信号要求信道具有线性相位的特性。

滤波器分为低通型、高通型、带通型和带阻型。而低通滤波器是最基本的滤波器。在本次设计中，基本的功能是比较容易实现的，因此增添了一些附加的功能。

指导书中规定了h(n)的长度及值，正好符合线性相位的第一类的第一种情况，所以加窗后的低通滤波器理论上具有线性相位。当然加窗后处理的低通滤波器由于阶段后出现的误差，其相位不一定为十分严格的直线，但在相当宽的一段频率中是符合严格的线性相位的。

在规定了h(n)及其长度后，窗函数的类型就成为影响低通滤波器的性能的主要因素了。因此本次设计最重要的一个附加功能就是加不同的窗函数后，所对应的低通滤波器的形状的对比。

通过对比发现，加不同的窗函数后，低通滤波器的衰减到0.9和衰减到0.1对应的频率的不同，和尾巴衰减的速率的不同。Matlab中自带有6种窗函数，因此制作了六种窗函数的对比。由于用脚本m文件画图，运行后生成六个图形的文件，不方便比较，于是制作了GUI界面用于显示绘制出的图形。

GUI界面的制作过程比较繁琐，并且subplot函数在坐标轴中绘图会出错所以，应用freqz函数绘图，不能够同时绘制相位图。由于指导书要求的只有幅频特性，因此将freqz函数与abs函数结合，绘制出低通滤波器的幅频响应。

通过本次设计，我了解了学习matlab这个强大的软件的方法，就是多多实际，注重自学。Matlab中的函数十分便捷，只需要用fir1函数就可以得到相应的低通滤波器。同时，只需要freqz函数就可以求出低通滤波器的幅频特性和相频特性。

8 参考文献

[1]刘泉.数字信号处理原理与实现.北京：电子工业出版社，2009

[2] 张志涌.精通matlab .北京：北京航空航天大学出版社，2011

附录

源程序

function varargout = chuanghanshu(varargin)

% CHUANGHANSHU M-file for chuanghanshu.fig

%      CHUANGHANSHU, by itself, creates a new CHUANGHANSHU or raises the existing

%      singleton*.

%

%      H = CHUANGHANSHU returns the handle to a new CHUANGHANSHU or the handle to

%      the existing singleton*.

%

%      CHUANGHANSHU('CALLBACK',hObject,eventData,handles,...) calls the local

%      function named CALLBACK in CHUANGHANSHU.M with the given input arguments.

%

%      CHUANGHANSHU('Property','Value',...) creates a new CHUANGHANSHU or raises the

%      existing singleton*.  Starting from the left, property value pairs are

%      applied to the GUI before chuanghanshu_OpeningFcn gets called.  An

%      unrecognized property name or invalid value makes property application

%      stop.  All inputs are passed to chuanghanshu_OpeningFcn via varargin.

%

%      *See GUI Options on GUIDE's Tools menu.  Choose "GUI allows only one

%      instance to run (singleton)".

%

% See also: GUIDE, GUIDATA, GUIHANDLES

% Edit the above text to modify the response to help chuanghanshu

% Last Modified by GUIDE v2.5 27-Jun-2013 10:14:02

% Begin initialization code - DO NOT EDIT

gui_Singleton = 1;

gui_State = struct('gui_Name',       mfilename, ...

                   'gui_Singleton',  gui_Singleton, ...

                   'gui_OpeningFcn', @chuanghanshu_OpeningFcn, ...

                   'gui_OutputFcn',  @chuanghanshu_OutputFcn, ...

                   'gui_LayoutFcn',  [] , ...

                   'gui_Callback',   []);

if nargin && ischar(varargin{1})

    gui_State.gui_Callback = str2func(varargin{1});

end

if nargout

    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});

else

    gui_mainfcn(gui_State, varargin{:});

end

% End initialization code - DO NOT EDIT

% --- Executes just before chuanghanshu is made visible.

function chuanghanshu_OpeningFcn(hObject, eventdata, handles, varargin)

% This function has no output args, see OutputFcn.

% hObject    handle to figure

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

% varargin   command line arguments to chuanghanshu (see VARARGIN)

% Choose default command line output for chuanghanshu

handles.output = hObject;

% Update handles structure

guidata(hObject, handles);

% UIWAIT makes chuanghanshu wait for user response (see UIRESUME)

% uiwait(handles.figure1);

% --- Outputs from this function are returned to the command line.

function varargout = chuanghanshu_OutputFcn(hObject, eventdata, handles) 

% varargout  cell array for returning output args (see VARARGOUT);

% hObject    handle to figure

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

% Get default command line output from handles structure

varargout{1} = handles.output;

% --- Executes on button press in pushbutton1.

function pushbutton1_Callback(hObject, eventdata, handles)

% hObject    handle to pushbutton1 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn,boxcar(N+1));

[w,t]=freqz(b,1,512);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0);

title('矩形窗');

% --- Executes on button press in pushbutton2.

function pushbutton2_Callback(hObject, eventdata, handles)

% hObject    handle to pushbutton2 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn);

[w,t]=freqz(b,1,512);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0)

title('汉明窗');

% --- Executes on button press in pushbutton3.

function pushbutton3_Callback(hObject, eventdata, handles)

% hObject    handle to pushbutton3 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn,triang(N+1));

[w,t]=freqz(b,1,512);

axis([0,3.5,0,1.5]);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0)

axis([0,3.5,0,1.4]);

title('三角窗');

% --- Executes on button press in pushbutton4.

function pushbutton4_Callback(hObject, eventdata, handles)

% hObject    handle to pushbutton4 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn,hanning(N+1));

[w,t]=freqz(b,1,512);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0)

title('汉宁窗 ');

% --- Executes on button press in pushbutton5.

function pushbutton5_Callback(hObject, eventdata, handles)

% hObject    handle to pushbutton5 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn,blackman(N+1));

[w,t]=freqz(b,1,512);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0)

title('布莱克曼窗');

% --- Executes on button press in pushbutton6.

function pushbutton6_Callback(hObject, eventdata, handles)

% hObject    handle to pushbutton6 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

D=str2num(get(handles.edit1,'string'));

wn=D*pi;

N=23;

b=fir1(N,wn,kaiser(N+1));

[w,t]=freqz(b,1,512);

QX=plot(t,abs(w));

set(QX,'LineWidth',2.0)

title('凯泽窗');

% --- Executes during object creation, after setting all properties.

function axes1_CreateFcn(hObject, eventdata, handles)

% hObject    handle to axes1 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    empty - handles not created until after all CreateFcns called

% Hint: place code in OpeningFcn to populate axes1

function pushbutton13_Callback(hObject, eventdata, handles)

% hObject    handle to pushbutton13 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

n=0:22;

x=n./n;

stem(n,x); title('h(n)');

axis([0,25,0,1.3]);

% --- Executes on button press in pushbutton14.

function pushbutton14_Callback(hObject, eventdata, handles)

% hObject    handle to pushbutton14 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

D=str2num(get(handles.edit1,'string'));

N=23;

n=0:1:N-1;

Wc=D*pi;

hd=ideal_lp(Wc,N);

QX=plot(n,hd);

set(QX,'LineWidth',2.0)

title('h(n)的频率响应');

function edit1_Callback(hObject, eventdata, handles)

% hObject    handle to edit1 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    structure with handles and user data (see GUIDATA)

% Hints: get(hObject,'String') returns contents of edit1 as text

%        str2double(get(hObject,'String')) returns contents of edit1 as a double

% --- Executes during object creation, after setting all properties.

function edit1_CreateFcn(hObject, eventdata, handles)

% hObject    handle to edit1 (see GCBO)

% eventdata  reserved - to be defined in a future version of MATLAB

% handles    empty - handles not created until after all CreateFcns called

% Hint: edit controls usually have a white background on Windows.

%       See ISPC and COMPUTER.

if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))

    set(hObject,'BackgroundColor','white');

end