滤波器，窗函数以及切比雪夫逼近法

     
 滤波器，是指对输入信号起到滤波作用的系统，可分为模拟滤波器和数字滤波器两大类。其中，数字滤波器是数字信号处理的重要环节，具有可靠性好、精度高、灵活性大等优点，广泛地应用在语音和图像处理、模式识别、频谱分析等方面。数字滤波器按其单位样值响应的性质可分为有限冲激响应滤波器（FIR）和无限冲激响应滤波器（IIR）。

     
 本课题主要讨论了数字有限冲激响应滤波器（FIR）的实现，主要包括窗函数法和基于切比雪夫逼近法两种方法，分别对该两种方法进行了理论介绍，并进行了MATLAB仿真。并重点介绍了基于切比雪夫逼近法的滤波器的设计与实现，并使用了音乐WAV文件进行测试，通过切比雪夫逼近法滤波器对音频文件进行滤波，从时域和频域两方面分析其滤波效果。最后利用MATLAB具有强大的科学计算和图形显示这一优点,设计了一个简单的GUI界面，通过GUI界面，可以简单的进行系统操作。

滤波器，是指对输入信号起到滤波作用的系统。根据处理的信号不同，可分为模拟滤波器和数字滤波器两大类。若滤波器的输入、输出都是离散时间信号，则该滤波器的冲激响应也必然是离散的，我们称这样的滤波器为数字滤波器（Digital
Filter）。

数字滤波器是数字信号处理的重要环节，它实质是用有限精度算法实现的离散时间线性时不变系统，从而完成对信号进行滤波处理的功能。具有可靠性好、精度高、灵活性大等优点，广泛地应用在语音和图像处理、模式识别、频谱分析等方面。

数字滤波器是具有一定传输选择特性的数字信号处理装置，其输入输出均为数字信号，实质上是一个由有限精度算法实现的线性时不变（LTI）离散系统。它的基本工作原理是利用离散系统特性对系统输入信号进行加工、处理和变换，改变输入序列的频谱或信号波形，让有用频率的信号分量通过，抑制无用的信号分量输出。根据频率响应特性，数字滤波器可分为低通、高通、带通、带阻等类型，与模拟滤波器相比，除了具有数字信号处理的固有优点外，还有滤波精度高、稳定性好、灵活性强等优点。数字滤波器按其单位样值响应的性质可分为有限冲激响应滤波器（FIR）和无限冲激响应滤波器（IIR）。

因此，通过一种合适的方式来设计滤波器具有重要的意义。由于FIR滤波器具有跟多的优势，本课题将重点研究FIR滤波器的实现。FIR滤波器的实现主要有窗函数法和基于切比雪夫逼近法等多种实现方式。

窗函数法的基本思想是把频率响应通过IDTFT求得脉冲响应，然后利用加窗函数对它进行截断和平滑，以实现一个物理可实现且具有线性相位的 FIR
数字滤波器的设计目的。

切比雪夫逼近法，其目的是在所需要的区间内，使误差函数E(x)较均匀一致，并且通过合理地选择多项式，使E(x)的最大值达到最小。即，使最大误差最小化。

基于切比雪夫逼近法的滤波器的设计与仿真

数字滤波器的设计是现代数字信号处理的重要内容。常用的数字滤波器有FIR和IIR两种类型，两者比较而言。主要有如下几点区别：

从系统的幅频特性来看，IIR滤波器由于综合利用了系统的零极点，容易达到比较理想的设计效果；
而FIR滤波器由于只有零点，效果较IIR滤波器差。要达到与IIR滤波器相似的效果，往往要提高系统的阶数。

其次，从相位特性来看，用FIR滤波器可以得到线性相位数字滤波器，满足信号不失真传输的要求。 对于IIR滤波器而言，往往幅频特性越好，相位非线性就越严重。

再次，从系统稳定性来看，FIR滤波器由于没有极点，所以一定是稳定的；而IIR滤波器的稳定与否取决于其极点的位置。

最后，从设计方法来看，IIR滤波器的设计参照连续时间系统的传输函数进行，可以充分利用模拟滤波器的设计结果，但是要求设计者有一定的模拟滤波器的设计知识，而且必须保证在模拟滤波器中能够找到合适的滤波器原型作为设计基础；而FIR滤波器设计结果完全是根据系统频率进行，不需要设计者有其他滤波器的知识，设计方法比较简单。

综上所述，本课题决定使用FIR数字滤波器。

FIR滤波器设计的方法主要有窗函数法和切比雪夫逼近法等。本章将简单的介绍一下基于窗函数的FIR滤波器的设计与实现，然后就窗函数的缺陷，进一步介绍切比雪夫逼近法的滤波器的设计。

为了选取效果最佳的窗函数滤波器，我们这里统一设定a=0.5,Ts=1/1000,通过改变N理对比滤波器的性能，我们分别通过比较N=16,32,128来比较。

矩形窗

图2-1 矩形窗频谱特性

   
图2-1，可以得知，虽然矩形窗函数FIR滤波器实现起来比较简单，但是随着N阶数的不断增加，其窗函数的旁瓣并没有明显的变化，能量也没有完全集中在主瓣上，应此基于矩形窗函数的FIR滤波器只能应用与简单的场合，而且随着阶数的增加，其系统新能并不能得到很的改善，应此，在实际信号处理中，不使用矩形窗函数的FIR滤波器。

汉宁窗

图2-2汉宁窗频谱特性

   
图2-2，可以得知，汉宁窗函数的性能相对于矩形窗函数有了明显的提升，但是当阶数较低的时候，其能力也并不完全集中在主瓣，当阶数较高的情况下，才能明显感觉到。但是当阶数较高的情况下，对设备成本有较大的影响。

汉明窗

图2-3汉明窗频谱特性

    汉明窗是在汉宁窗的基础改进的窗函数，有图可以看见，其能量基本集中在主瓣中，而且当阶数较低时，能量也基本集中在主瓣上，性能比汉宁窗和矩形窗更好。

布拉克曼窗

图2-4布拉克曼窗

布拉克曼窗函数其旁瓣衰减也是随着阶数的增加而增加，但是布拉克曼窗函数在硬件上实现较前三种窗函数而言，更加复杂。

以上就是窗函数的仿真。

   
从上面的仿真可以看到，用窗函数设计FIR滤波器的设计是用有限序列去代替无限长序列，这会引起误差，表现在频域就是通常所说的吉布斯(Gibbs)效应。下面，我们将重点介绍基于切比雪夫逼近法的FIR滤波器的设计与实现。

下面首先通过一个简单的测试进行测试，输入MATLAB代码如下所示：

f     =[0 0.25 0.5 1]; %给定频率轴分点;

a     =[1 1 0 0]; %给定在这些频率分点上理想的幅频响应;

weigh =[1 10]; %给定在这些频率分点上的加权;       

b=remez(10,f,a,weigh); %设计出切比雪夫最佳逼近滤波器;

测试信号为：

图2.5 测试信号及其频谱图

 

针对真个信号的特点，设计如下的切比雪夫逼近滤波器。

图2.6 基于切比雪夫逼近算法的滤波器频谱特性图

    滤波之后，可以看到如下的效果：

图2.7 滤波之后信号及其频谱

    通过上述仿真可以看到，设计的基于切比雪夫逼近法的滤波器具有良好的滤波效果。下面通过适当修改滤波器的参数，对一组音乐信号进行滤波。

这个就是remez算法的基本过程。

在MATLAB中通过如下代码实现：

%给定频率轴分点;

f     =[0 0.4  0.6 1];

%给定在这些频率分点上理想的幅频响应;

a     =[1 1 0 0];   

%给定在这些频率分点上的加权;

weigh =[1 10];  

%设计出切比雪夫最佳逼近滤波器;

b=remez(32,f,a,weigh);

运行MATLAB,可以得到通过切比雪夫逼近法得到的滤波器的幅频特性：

图2.8 切比雪夫逼近滤波器的抽样值和幅频特性

   
从上面的仿真可以看到，基于切比雪夫逼近法实现得到的滤波器具有跟好的幅频特性。下面我们将给予切比雪夫逼近法滤波器进行测试，对实际中带噪声的信息进行滤波处理。

   下面我们加入一个测试语音信号进行测试滤波器的性能。

假设输入的信号为：

图2.9 原始音乐信号的时域波形图

通过傅里叶变换得到原始音乐信号的傅里叶变换，在MATLAB中输入如下的代码：

y=fftshift(abs(fft(source)));

    运行MATLAB,可以看到信号的频谱为：

图2.10 原始音乐信号的频谱图

然后认为的加入噪声，加入噪声得到如下的仿真图：

图2.11 加入噪声后的信号时域图

图2.12 加入噪声的信号频谱图

从上面的频谱图可以看到，滤波之前的音乐信号，其频谱比较杂，在高频处存在较多的噪声，因此从听觉角度看，听起来会出现小的噪声干扰。

将此音乐信号通过滤波器，通过观察其频域信号就能看到滤波效果.

图2.13 滤波后的音乐信号

2.4 基于GUI的用户交互界面的设计

对于一个完整的系统，都具备一个GUI界面，MATLAB提供了功能强大的GUI界面设计模块，本节，我们将介绍使用MATLAB进行GUI界面的设计与实现。本系统的界面如下所示：

图2.14 GUI界面效果图

这里主要涉及到测试信号的滤波，还要打开外部WAV音频，然后对音频加入噪声，最后是滤波，这几个部分系统的主要功能模块。在GUI设计界面中首先搭好一个系统界面框架：

图2.15 GUI界面效果图

其中打开WAV文件对应的MATLAB代码如下所示：

function pushbutton3_Callback(hObject, eventdata, handles)

……….

 [filename,pathname]=uigetfile({ ...

    '*.*','All Files(*.*)';},...

    '选择文件');

if isequal([filename,pathname],[0,0])

    return

else

    pic = fullfile(pathname,filename);

    [source2,Fs,nbits]= wavread(pic);

end

……….

通过执行这个语句，在GUI界面上，就是自动弹出对话框，然后选择需要打开的WAV文件。

给加入的WAV文件加入噪声，其对应的MATLAB代码如下所示：

function pushbutton6_Callback(hObject, eventdata, handles)

……….

%添加噪声

sums = zeros(1,length(source2));

for i=1:60

n=1:((80000+1000*i)*pi-1)/length(source2):(80000+1000*i)*pi;

noise = 0.3*rand(1)*sin(n);

sums   = sums + noise(1:length(source2));

end

source_noise = source2 + sums';

source_noise = AWGN(source_noise,20);

……….

通过执行这个语句，在GUI界面上，就会自动对系统添加噪声。

最后一个主要功能模块是对带有噪声的信号进行滤波，其对应的MATLAB代码：

function pushbutton4_Callback(hObject, eventdata, handles)

……….

R1T=get(handles.edit1,'string')

if isempty(str2num(R1T))

        warndlg('Input is wrong,please input numeric type')

else

        R1T = str2num(R1T);

end

R2T=get(handles.edit2,'string')

if isempty(str2num(R2T))

        warndlg('Input is wrong,please input numeric type')

else

        R2T = str2num(R2T);

end

R3T=get(handles.edit3,'string')

if isempty(str2num(R3T))

        warndlg('Input is wrong,please input numeric type')

else

        R3T = str2num(R3T);

end

f     =[0 R1T R2T 1];    %给定频率轴分点;

a     =[1 1 0 0];      %给定在这些频率分点上理想的幅频响应;

weigh =[1 10];         %给定在这些频率分点上的加权;

b=remez(R3T,f,a,weigh); %设计出切比雪夫最佳逼近滤波器;

[h,w]=freqz(b,1,512,1);%数字滤波器的频率响应;

h=abs(h);%绝对值;

h=20*log10(h);

source2S=filter(b,1,source_noise);    %滤波正常实现

……….

通过执行这个语句，在GUI界面上，就会显示滤波之后的效果。

 

附录：

clc;

clear;

close all;

select = 2;%1：进行普通信号测试；2：音乐信号的测试

if select == 1

%普通信号测试

source = func_test_signal();

figure;

subplot(121);plot(source):title('原始音乐信号');

f     =[0 0.25 0.5 1]; %给定频率轴分点;

a     =[1 1 0 0];     %给定在这些频率分点上理想的幅频响应;

weigh =[1 10];         %给定在这些频率分点上的加权;

y=fftshift(abs(fft(source)));

subplot(122);plot(y);title('原始信号的频谱');

figure;

b=remez(10,f,a,weigh);%设计出切比雪夫最佳逼近滤波器;

[h,w]=freqz(b,1,512,1);%数字滤波器的频率响应;

h=abs(h);%绝对值;

h=20*log10(h);

subplot(121)%改置定位坐标系;

stem(b,'.');

grid;

title('切比雪夫逼近滤波器的抽样值');

subplot(122);

plot(w,h);%生成参数方程的图形;

grid;

title('切比雪夫逼近滤波器幅频特性(dB)');

source2=filter(b,1,source);    %滤波正常实现

figure;

subplot(121);plot(source2):title('通过切比雪夫逼近法滤波后的音乐信号');

y2=fftshift(abs(fft(source2)));

subplot(122);plot(y2);title('滤波后的信号频谱');

end

if select == 2

%音乐信号测试

[source,Fs,nbits] = wavread('source.wav');

figure;

plot(source):title('原始音乐信号');

f     =[0 0.4 0.6 1]; %给定频率轴分点;

a     =[1 1 0 0];     %给定在这些频率分点上理想的幅频响应;

weigh =[1 10];         %给定在这些频率分点上的加权;

y=fftshift(abs(fft(source)));

figure;

plot(y):title('原始信号的频谱');

sound(source,Fs);

figure;

b=remez(32,f,a,weigh);%设计出切比雪夫最佳逼近滤波器;

[h,w]=freqz(b,1,512,1);%数字滤波器的频率响应;

h=abs(h);%绝对值;

h=20*log10(h);

subplot(121)%改置定位坐标系;

stem(b,'.');

grid;

title('切比雪夫逼近滤波器的抽样值');

subplot(122);

plot(w,h);%生成参数方程的图形;

grid;

title('切比雪夫逼近滤波器幅频特性(dB)');

source2=filter(b,1,source);    %滤波正常实现

figure;

plot(source2);title('通过切比雪夫逼近法滤波后的音乐信号');

y2=fftshift(abs(fft(source2)));

figure;

plot(y2);title('滤波后的信号频谱');

sound(source2,Fs);

end

 

 

• « 上一篇：智慧城市发展情况及未来
• » 下一篇：前端三剑客（入门）