三，滤波器设计：

1、相关原理：

设计数字滤波器的任务就是寻求一个因果稳定的线性时不变系统，并使系统函数H(z)具有指定的频率特性。

数字滤波器从实现的网络结构或者从单位冲激响应分类，可以分成无限长单位冲激响应（IIR）数字滤波器和有限长单位冲激响应（FIR）数字滤波器。

数字滤波器频率响应的三个参数：

（1）    幅度平方响应：

（2）    相位响应          

（3）    群时延响应

IIR数字滤波器：

IIR数字滤波器的系统函数为 的有理分数，即

IIR数字滤波器的逼近问题就是求解滤波器的系数 和，使得在规定的物理意义上逼近所要求的特性的问题。如果是在s平面上逼近，就得到模拟滤波器，如果是在z平面上逼近，则得到数字滤波器。

FIR数字滤波器：

设FIR的单位脉冲响应h(n)为实数，长度为N，则其z变换和频率响应分别为

按频域采样定理FIR数字滤波器的传输函数H(z)和单位脉冲响应h(n)可由它的N个频域采样值H(k)唯一确定。

MATLAB中提供了几个函数，分别用于实现IIR滤波器和FIR滤波器。

（1）卷积函数conv

卷积函数conv的调用格式为           c=conv(a,b)

该格式可以计算两向量a和b的卷积，可以直接用于对有限长信号采用FIR滤波器的滤波。

（2）函数filter

函数filter的调用格式为          y=filter(b,a,x)

该格式采用数字滤波器对数据进行滤波，既可以用于IIR滤波器，也可以用于FIR滤波器。其中向量b和a分别表示系统函数的分子、分母多项式的系数，若a＝1，此时表示FIR滤波器，否则就是IIR滤波器。该函数是利用给出的向量b和a，对x中的数据进行滤波，结果放入向量y。

（3）函数fftfilt

函数fftfilt的调用格式为         y=fftfilt(b,x)

该格式是利用基于FFT的重叠相加法对数据进行滤波，这种频域滤波技术只对FIR滤波器有效。该函数是通过向量b描述的滤波器对x数据进行滤波。

关于用butter函数求系统函数分子与分母系数的几种形式。

[b,a]=butter(N,wc,'high'):设计N阶高通滤波器，wc为它的3dB边缘频率，以 为单位，故 。

[b,a]=butter(N,wc):当wc为具有两个元素的矢量wc=[w1,w2]时，它设计2N阶带通滤波器，3dB通带为，w的单位为 。

[b,a]=butter(N,wc,'stop'):若wc=[w1,w2]，则它设计2N阶带阻滤波器，3dB通带为 ，w的单位为。

如果在这个函数输入变元的最后，加一个变元“s”，表示设计的是模拟滤波器。这里不作讨论。

为了设计任意的选项巴特沃斯滤波器，必须知道阶数N和3dB边缘频率矢量wc。这可以直接利用信号处理工具箱中的buttord函数来计算。如果已知滤波器指标， ， 和 ，则调用格式为

[N,wc]=buttord(wp,ws,Rp,As)

对于不同类型的滤波器，参数wp和ws有一些限制：对于低通滤波器，wpws；对于带通滤波器，wp和ws分别为具有两个元素的矢量，wp=[wp1,wp2]和ws=[ws1,ws2]，并且 ws1

2、设计内容：

（1）滤波器示例：

在这里为了说明如何用MATLAB来实现滤波，特举出一个简单的函数信号滤波实例（对信号x(n)=sin( n/4)+5cos( n/2)进行滤波，信号长度为500点），从中了解滤波的实现过程。程序如下：

Wn=0.2*pi;

N=5;

[b,a]=butter(N,Wn/pi);

n=0:499;

x=sin(pi*n/4)+5*cos(pi*n/2);

X=fft(x,4096);

subplot(221);plot(x);title('滤波前信号的波形');

subplot(222);plot(X);title('滤波前信号的频谱');

y=filter(b,a,x);

Y=fft(y,4096);

subplot(223);plot(y);title('滤波后信号的波形');

subplot(224);plot(Y);title('滤波后信号的频谱');

在这里，是采用了butter命令，设计出一个巴特沃斯低通滤波器，从频谱图中可以很明显的看出来。下面，也就是本课题的主要内容，也都是运用到了butter函数，以便容易的得到系统函数的分子与分母系数，最终以此来实现信号的滤波。

（2）N阶高通滤波器的设计（在这里，以5阶为例，其中wc为其3dB边缘频率，以 为单位），程序设计如下：

x=wavread('ding.wav');

sound(x);

N=5;wc=0.3;

[b,a]=butter(N,wc,'high');

X=fft(x);

subplot(321);plot(x);title('滤波前信号的波形');

subplot(322);plot(X);title('滤波前信号的频谱');

y=filter(b,a,x);

Y=fft(y);

subplot(323);plot(y);title('IIR滤波后信号的波形');

subplot(324);plot(Y);title('IIR滤波后信号的频谱');

z=fftfilt(b,x);

Z=fft(z);

subplot(325);plot(z);title('FIR滤波后信号的波形');

subplot(326);plot(Z);title('FIR滤波后信号的频谱');

（3）N阶低通滤波器的设计（在这里，同样以5阶为例，其中wc为其3dB边缘频率，以 为单位），程序设计如下：

x=wavread('ding.wav');

sound(x);

N=5;wc=0.3;

[b,a]=butter(N,wc);

X=fft(x);

subplot(321);plot(x);title('滤波前信号的波形');

subplot(322);plot(X);title('滤波前信号的频谱');

y=filter(b,a,x);

Y=fft(y);

subplot(323);plot(y);title('IIR滤波后信号的波形');

subplot(324);plot(Y);title('IIR滤波后信号的频谱');

z=fftfilt(b,x);

Z=fft(z);

subplot(325);plot(z);title('FIR滤波后信号的波形');

subplot(326);plot(Z);title('FIR滤波后信号的频谱');

（4）2N阶带通滤波器的设计（在这里，以10阶为例，其中wc为其3dB边缘频率，以 为单位，wc=[w1,w2],w1 wc w2），程序设计如下：

x=wavread('ding.wav');

   sound(x);

N=5;wc=[0.3,0.6];

[b,a]=butter(N,wc);

X=fft(x);

subplot(321);plot(x);title('滤波前信号的波形');

subplot(322);plot(X);title('滤波前信号的频谱');

y=filter(b,a,x);

Y=fft(y);

subplot(323);plot(y);title('IIR滤波后信号的波形');

subplot(324);plot(Y);title('IIR滤波后信号的频谱');

z=fftfilt(b,x);

Z=fft(z);

subplot(325);plot(z);title('FIR滤波后信号的波形');

subplot(326);plot(Z);title('FIR滤波后信号的频谱');

 

（5）2N阶带阻滤波器的设计（在这里，以10阶为例，其中wc为其3dB边缘频率，以 为单位，wc=[w1,w2],w1 wc w2），程序设计如下：

x=wavread('ding.wav');

   sound(x);

N=5;wc=[0.2,0.7];

[b,a]=butter(N,wc,'stop');

X=fft(x);

subplot(321);plot(x);title('滤波前信号的波形');

subplot(322);plot(X);title('滤波前信号的频谱');

y=filter(b,a,x);

Y=fft(y);

subplot(323);plot(y);title('IIR滤波后信号的波形');

subplot(324);plot(Y);title('IIR滤波后信号的频谱');

z=fftfilt(b,x);

Z=fft(z);

subplot(325);plot(z);title('FIR滤波后信号的波形');

subplot(326);plot(Z);title('FIR滤波后信号的频谱');

（6）小结：以上几种滤波，我们都可以从信号滤波前后的波形图以及频谱图上看出变化。当然，也可以用sound()函数来播放滤波后的语音，从听觉上直接感受语音信号的变化，但由于人耳听力的限制，有些情况下我们是很难听出异同的。

同样，通过函数的调用，也可以将信号的频谱进行“分离观察”，如显出信号的幅值或相位。下面，通过改变系统函数的分子与分母系数比，来观察信号滤波前后的幅值与相位。并且使结果更加明显，使人耳得以很容易的辨听。

x=wavread('ding.wav');

   sound(x);

    b=100;a=5;

y=filter(b,a,x);

X=fft(x,4096);

subplot(221);plot(x);title('滤波前信号的波形');

subplot(222);plot(abs(X));title('滤波前信号的幅值');

Y=fft(y,4096);

subplot(223);plot(y);title('滤波后信号的波形');

subplot(224);plot(abs(Y));title('滤波后信号的幅值');

>> sound(y);

可以听到声音明显变得高亢了。从上面的波形与幅值（即幅频）图，也可看出，滤波后的幅值变成了滤波前的20倍。

>> figure,

subplot(211);plot(angle(X));title('滤波前信号相位');

subplot(212);plot(angle(Y));title('滤波后信号相位');

得图：

可以看到相位谱没什么变化。

（四）、界面设计：

直接用M文件编写GUI程序很繁琐，而使用GUIDE设计工具可以大大提高工作效率。GUIDE相当于一个控制面板，从中可以调用各种设计工具以辅助完成界面设计任务，例如控件的创建和布局、控件属性的编辑和菜单设计等。

使用GUIDE设计GUI程序的一般步骤如下：

1. 将所需控件从控件面板拖拽到GUIDE的设计区域；

2. 利用工具条中的工具（或相应的菜单和现场菜单），快速完成界面布局；

3. 设置控件的属性。尤其是tag属性，它是控件在程序内部的唯一标识；

4. 如果需要，打开菜单编辑器为界面添加菜单或现场菜单；

5. 保存设计。GUIDE默认把GUI程序保存为两个同名文件：一个是.fig文件，用来保存窗体布局和所有控件的界面信息；一个是.m文件，该文件的初始内容是GUIDE自动产生的程序框架，其中包括了各个控件回调函数的定义。该M文件与一般的M文件没有本质区别，但是鉴于它的特殊性，MATALAB把这类文件统称为GUI-M文件。保存完后GUI-M文件自动在编辑调试器中打开以供编辑。

6.      为每个回调函数添加代码以实现GUI程序的具体功能。这一步与一般函数文件的编辑调试过程相同。

设计过程及内容：

在MATLAB版面上，通过键入GUIDE弹出一个菜单栏进入gui制作界面（或者在File到new来进入gui）,从而开始应用界面的制作。

该界面主要实现了以下几个功能：

①打开wav格式的音频文件，并将该音频信号的值读取并赋予某一向量；

②播放音频文件，可以选择性的显示该音频信号的波形、频谱、幅值以及相位；

③对音频信号进行IIR与FIR的5阶固定滤波处理，可以选择性的显示滤波前后信号的波形、频谱、幅值以及相位，以及播放滤波后的声音。

界面如图所示：

通过该界面，可以方便用户进行语音信号的处理。

界面主程序见附件。

（五）、校验：

1、本设计圆满的完成了对语音信号的读取与打开，与课题的要求十分相符；

2、本设计也较好的完成了对语音信号的频谱分析，通过fft变换，得出了语音信号的频谱图；

3、在滤波这一块，课题主要是从巴特沃斯滤波器入手来设计滤波器，也从一方面基本实现了滤波；

4、初略的完成了界面的设计，但也存在相当的不足，只是很勉强的达到了打开语音文件、显示已定滤波前后的波形等图。

 

四、 结论：

语音信号处理是语音学与数字信号处理技术相结合的交叉学科，课题在这里不讨论语音学，而是将语音当做一种特殊的信号，即一种“复杂向量”来看待。也就是说，课题更多的还是体现了数字信号处理技术。

从课题的中心来看，课题是希望将数字信号处理技术应用于某一实际领域，这里就是指对语音的处理。作为存储于计算机中的语音信号，其本身就是离散化了的向量，我们只需将这些离散的量提取出来，就可以对其进行处理了。

下面介绍七种滤波方法的滤波结果

创建两个混合信号，便于更好测试滤波器效果。同时用七中滤波方法测试。
混合信号Mix_Signal_1 = 信号Signal_Original_1+白噪声。
http://bbs.armfly.com/attachment/Fid_2/2_58_fa6d1c834991ec4.png?6滤波器 及滤波效果对比(转载）" />
http://bbs.armfly.com/attachment/Fid_2/2_58_bbe8d85b202d1d7.png?8滤波器 及滤波效果对比(转载）" />



混合信号Mix_Signal_2 = 信号Signal_Original_2+白噪声。
http://bbs.armfly.com/attachment/Fid_2/2_58_bcc671609e379d9.png?5滤波器 及滤波效果对比(转载）" />
http://bbs.armfly.com/attachment/Fid_2/2_58_aad4cb3351b8ce2.png?7滤波器 及滤波效果对比(转载）" />



1.巴特沃斯低通滤波器去噪
巴特沃斯滤波器适合用于信号和噪声没有重叠的情况下。下图是巴特沃斯对两个信号的滤波效果。
http://bbs.armfly.com/attachment/Fid_2/2_58_a5aba2264d4ef7b.png?11滤波器 及滤波效果对比(转载）" />

从图上可以看出巴特沃斯低通滤波器对信号一的滤波效果还是可以的，主要是因为有效的信号最高频率才30Hz,本程序将50Hz以上的信号全部滤除，通过的频率成分中仍然是有白噪声的。
对于信号二，滤波后的信号与没有加噪声的信号相比就有失真了，上升沿和下降沿的高频信号被滤除了。
2.FIR低通滤波器去噪
情况同巴特沃斯低通滤波器相似。滤波后的效果如下：
http://bbs.armfly.com/attachment/Fid_2/2_58_0b0c108b60e5bfb.png?10滤波器 及滤波效果对比(转载）" />

3. 移动平均滤波去噪
滤波效果如下：
http://bbs.armfly.com/attachment/Fid_2/2_58_4885c9cd30e49a4.png?10滤波器 及滤波效果对比(转载）" />

4. 中值滤波去噪
http://bbs.armfly.com/attachment/Fid_2/2_58_f9d7d4a94f79e18.png?10滤波器 及滤波效果对比(转载）" />

从上图可以看出，无论是对信号一还是对信号二，中值滤波的滤波效果都是很不错，特备是对于信号二，上升沿和下降失真比较的小。5. 维纳滤波去噪

维纳滤波器属于现代滤波器，传统的滤波器只能滤除信号和干扰频带没有重叠的情况，当信号和干扰频带有重叠的时候传统滤波器将无能为力，这时就需要用到现代滤波器，现代滤波器利用信号和干扰的统计特征（如自相关函数、功率谱等）导出一套最佳估值算法，然后用硬件或软件予以实现。
维纳滤波是以均方误差最小（LMS（Least MeanSquare）为准则的，它根据过去观测值和当前观测值来估计信号的当前值，因此它的解形式是系统的传递函数或单位脉冲响应。
均方误差为：
http://bbs.armfly.com/attachment/Fid_2/2_58_fb11735b6cbccca.png?2滤波器 及滤波效果对比(转载）" />

维纳－霍夫（Wiener－Hopf）方程最小均方误差下的解为：
http://bbs.armfly.com/attachment/Fid_2/2_58_102e47dcbb194cf.png?7滤波器 及滤波效果对比(转载）" />

由于理解不深，对于信号二，没有什么滤波效果
http://bbs.armfly.com/attachment/Fid_2/2_58_8c9515b19fa120a.png?11滤波器 及滤波效果对比(转载）" />

6. 自适应滤波去噪
维纳滤波器参数是固定的，适合于平稳随机信号。卡尔曼滤波器参数是时变的，适合于非平稳随机信号。然而，只有在信号和噪声的统计特性先验已知的情况下，这两种滤波技术才能获得最优滤波。在实际应用中，常常无法得到信号和噪声统计特性的先验知识。在这种情况下，自适应滤波技术能够获得极佳的滤波性能，因而具有很好的应用价值。
自适应滤波的滤波效果如下：
http://bbs.armfly.com/attachment/Fid_2/2_58_dd2d3c8273f8dc4.png?10滤波器 及滤波效果对比(转载）" />

本程序是基于LMS算法的自适应滤波，从上图可以看出，滤波效果也是很不错的，特别是对于信号二，上升沿有失真，下降沿保持还可以，最要的是得到的波形十分的平滑。由此可见自适应滤波极具使用价值。7. 小波去噪

首先看一下小波的去噪效果。
http://bbs.armfly.com/attachment/Fid_2/2_58_f3636cd8914fb12.png?10滤波器 及滤波效果对比(转载）" />

对于信号二，小波的去噪效果非常不错，虽然得到波形不是很平滑，但是上升沿和下降沿保持的非常高，基本可以看到棱角.
%****************************************************************************************
%  
%                      创建两个信号Mix_Signal_1 和信号 Mix_Signal_2 
%
%***************************************************************************************

Fs = 1000;                                                                        %采样率
N  = 1000;                                                                        %采样点数
n  = 0:N-1;
t   = 0:1/Fs:1-1/Fs;                                                            %时间序列 
Signal_Original_1 =sin(2*pi*10*t)+sin(2*pi*20*t)+sin(2*pi*30*t); 
Noise_White_1    = [0.3*randn(1,500), rand(1,500)];           %前500点高斯分部白噪声，后500点均匀分布白噪声
Mix_Signal_1   = Signal_Original_1 + Noise_White_1;        %构造的混合信号

Signal_Original_2  =  [zeros(1,100), 20*ones(1,20), -2*ones(1,30), 5*ones(1,80), -5*ones(1,30), 9*ones(1,140), -4*ones(1,40), 3*ones(1,220), 12*ones(1,100), 5*ones(1,20), 25*ones(1,30), 7 *ones(1,190)]; 
Noise_White_2     =  0.5*randn(1,1000);                                 %高斯白噪声
Mix_Signal_2        =  Signal_Original_2 + Noise_White_2;      %构造的混合信号

%****************************************************************************************
%  
%                信号Mix_Signal_1 和 Mix_Signal_2  分别作巴特沃斯低通滤波。
%
%***************************************************************************************

%混合信号 Mix_Signal_1  巴特沃斯低通滤波
figure(1);
Wc=2*50/Fs;                                          %截止频率 50Hz
[b,a]=butter(4,Wc);
Signal_Filter=filter(b,a,Mix_Signal_1);

subplot(4,1,1);                                        %Mix_Signal_1 原始信号                 
plot(Mix_Signal_1);
axis([0,1000,-4,4]);
title('原始信号 ');

subplot(4,1,2);                                        %Mix_Signal_1 低通滤波滤波后信号  
plot(Signal_Filter);
axis([0,1000,-4,4]);
title('巴特沃斯低通滤波后信号');

%混合信号 Mix_Signal_2  巴特沃斯低通滤波
Wc=2*100/Fs;                                          %截止频率 100Hz
[b,a]=butter(4,Wc);
Signal_Filter=filter(b,a,Mix_Signal_2);

subplot(4,1,3);                                        %Mix_Signal_2 原始信号                 
plot(Mix_Signal_2);
axis([0,1000,-10,30]);
title('原始信号 ');

subplot(4,1,4);                                       %Mix_Signal_2 低通滤波滤波后信号  
plot(Signal_Filter);
axis([0,1000,-10,30]);
title('巴特沃斯低通滤波后信号');

%****************************************************************************************
%  
%                信号Mix_Signal_1 和 Mix_Signal_2  分别作FIR低通滤波。
%
%***************************************************************************************

%混合信号 Mix_Signal_1  FIR低通滤波
figure(2);
F   =  [0:0.05:0.95]; 
A  =  [1    1      0     0     0    0      0     0     0    0     0     0     0     0     0     0    0   0   0   0] ;
b  =  firls(20,F,A);
Signal_Filter = filter(b,1,Mix_Signal_1);

subplot(4,1,1);                                          %Mix_Signal_1 原始信号                 
plot(Mix_Signal_1);
axis([0,1000,-4,4]);
title('原始信号 ');

subplot(4,1,2);                                          %Mix_Signal_1 FIR低通滤波滤波后信号  
plot(Signal_Filter);
axis([0,1000,-5,5]);
title('FIR低通滤波后的信号');

%混合信号 Mix_Signal_2  FIR低通滤波
F   =  [0:0.05:0.95]; 
A  =  [1    1      1     1     1    0      0    0     0    0     0     0     0     0     0     0    0   0   0   0] ;
b  =  firls(20,F,A);
Signal_Filter = filter(b,1,Mix_Signal_2);
subplot(4,1,3);                                          %Mix_Signal_2 原始信号                 
plot(Mix_Signal_2);
axis([0,1000,-10,30]);
title('原始信号 ');

subplot(4,1,4);                                          %Mix_Signal_2 FIR低通滤波滤波后信号  
plot(Signal_Filter);
axis([0,1000,-10,30]);
title('FIR低通滤波后的信号');

%****************************************************************************************
%  
%                信号Mix_Signal_1 和 Mix_Signal_2  分别作移动平均滤波
%
%***************************************************************************************

%混合信号 Mix_Signal_1  移动平均滤波
figure(3);
b  =  [1 1 1 1 1 1]/6;
Signal_Filter = filter(b,1,Mix_Signal_1);

subplot(4,1,1);                                          %Mix_Signal_1 原始信号                 
plot(Mix_Signal_1);
axis([0,1000,-4,4]);
title('原始信号 ');

subplot(4,1,2);                                          %Mix_Signal_1 移动平均滤波后信号  
plot(Signal_Filter);
axis([0,1000,-4,4]);
title('移动平均滤波后的信号');

%混合信号 Mix_Signal_2  移动平均滤波
b  =  [1 1 1 1 1 1]/6;
Signal_Filter = filter(b,1,Mix_Signal_2);
subplot(4,1,3);                                          %Mix_Signal_2 原始信号                 
plot(Mix_Signal_2);
axis([0,1000,-10,30]);
title('原始信号 ');

subplot(4,1,4);                                          %Mix_Signal_2 移动平均滤波后信号  
plot(Signal_Filter);
axis([0,1000,-10,30]);
title('移动平均滤波后的信号');

%****************************************************************************************
%  
%                信号Mix_Signal_1 和 Mix_Signal_2  分别作中值滤波
%
%***************************************************************************************

%混合信号 Mix_Signal_1  中值滤波
figure(4);
Signal_Filter=medfilt1(Mix_Signal_1,10);

subplot(4,1,1);                                          %Mix_Signal_1 原始信号                 
plot(Mix_Signal_1);
axis([0,1000,-5,5]);
title('原始信号 ');

subplot(4,1,2);                                          %Mix_Signal_1 中值滤波后信号  
plot(Signal_Filter);
axis([0,1000,-5,5]);
title('中值滤波后的信号');

%混合信号 Mix_Signal_2  中值滤波
Signal_Filter=medfilt1(Mix_Signal_2,10);
subplot(4,1,3);                                          %Mix_Signal_2 原始信号                 
plot(Mix_Signal_2);
axis([0,1000,-10,30]);
title('原始信号 ');

subplot(4,1,4);                                          %Mix_Signal_2 中值滤波后信号  
plot(Signal_Filter);
axis([0,1000,-10,30]);
title('中值滤波后的信号');

%****************************************************************************************
%  
%                信号Mix_Signal_1 和 Mix_Signal_2  分别作维纳滤波
%
%***************************************************************************************

%混合信号 Mix_Signal_1  维纳滤波
figure(5);
Rxx=xcorr(Mix_Signal_1,Mix_Signal_1);              %得到混合信号的自相关函数
M=100;                                                             %维纳滤波器阶数
for i=1:M                                                           %得到混合信号的自相关矩阵
    for j=1:M
        rxx(i,j)=Rxx(abs(j-i)+N);
    end
end
Rxy=xcorr(Mix_Signal_1,Signal_Original_1);       %得到混合信号和原信号的互相关函数
for i=1:M
    rxy(i)=Rxy(i+N-1);
end                                                                  %得到混合信号和原信号的互相关向量
h = inv(rxx)*rxy';                                               %得到所要涉及的wiener滤波器系数
Signal_Filter=filter(h,1, Mix_Signal_1);               %将输入信号通过维纳滤波器

subplot(4,1,1);                                                   %Mix_Signal_1 原始信号                 
plot(Mix_Signal_1);
axis([0,1000,-5,5]);
title('原始信号 ');

subplot(4,1,2);                                                   %Mix_Signal_1 维纳滤波后信号  
plot(Signal_Filter);
axis([0,1000,-5,5]);
title('维纳滤波后的信号');

%混合信号 Mix_Signal_2  维纳滤波
Rxx=xcorr(Mix_Signal_2,Mix_Signal_2);              %得到混合信号的自相关函数
M=500;                                                             %维纳滤波器阶数
for i=1:M                                                           %得到混合信号的自相关矩阵
    for j=1:M
        rxx(i,j)=Rxx(abs(j-i)+N);
    end
end
Rxy=xcorr(Mix_Signal_2,Signal_Original_2);       %得到混合信号和原信号的互相关函数
for i=1:M
    rxy(i)=Rxy(i+N-1);
end                                                                  %得到混合信号和原信号的互相关向量
h=inv(rxx)*rxy';                                               %得到所要涉及的wiener滤波器系数
Signal_Filter=filter(h,1, Mix_Signal_2);             %将输入信号通过维纳滤波器

subplot(4,1,3);                                                  %Mix_Signal_2 原始信号                 
plot(Mix_Signal_2);
axis([0,1000,-10,30]);
title('原始信号 ');

subplot(4,1,4);                                                  %Mix_Signal_2 维纳滤波后信号  
plot(Signal_Filter);
axis([0,1000,-10,30]);
title('维纳滤波后的信号');

%****************************************************************************************
%  
%                信号Mix_Signal_1 和 Mix_Signal_2  分别作自适应滤波
%
%***************************************************************************************

%混合信号 Mix_Signal_1 自适应滤波
figure(6);
N=1000;                                             %输入信号抽样点数N
k=100;                                                  %时域抽头LMS算法滤波器阶数
u=0.001;                                             %步长因子

%设置初值
yn_1=zeros(1,N);                                  %output signal
yn_1(1:k)=Mix_Signal_1(1:k);                 %将输入信号SignalAddNoise的前k个值作为输出yn_1的前k个值
w=zeros(1,k);                                        %设置抽头加权初值
e=zeros(1,N);                                        %误差信号

%用LMS算法迭代滤波
for i=(k+1):N
        XN=Mix_Signal_1((i-k+1):(i));
        yn_1(i)=w*XN';
        e(i)=Signal_Original_1(i)-yn_1(i);
        w=w+2*u*e(i)*XN;
end

subplot(4,1,1);
plot(Mix_Signal_1);                               %Mix_Signal_1 原始信号
axis([k+1,1000,-4,4]);
title('原始信号');

subplot(4,1,2);
plot(yn_1);                                            %Mix_Signal_1 自适应滤波后信号
axis([k+1,1000,-4,4]);
title('自适应滤波后信号');

%混合信号 Mix_Signal_2 自适应滤波
N=1000;                                             %输入信号抽样点数N
k=500;                                                %时域抽头LMS算法滤波器阶数
u=0.000011;                                        %步长因子

%设置初值
yn_1=zeros(1,N);                                   %output signal
yn_1(1:k)=Mix_Signal_2(1:k);                  %将输入信号SignalAddNoise的前k个值作为输出yn_1的前k个值
w=zeros(1,k);                                        %设置抽头加权初值
e=zeros(1,N);                                        %误差信号

%用LMS算法迭代滤波
for i=(k+1):N
        XN=Mix_Signal_2((i-k+1):(i));
        yn_1(i)=w*XN';
        e(i)=Signal_Original_2(i)-yn_1(i);
        w=w+2*u*e(i)*XN;
end

subplot(4,1,3);
plot(Mix_Signal_2);                               %Mix_Signal_1 原始信号
axis([k+1,1000,-10,30]);
title('原始信号');

subplot(4,1,4);
plot(yn_1);                                            %Mix_Signal_1 自适应滤波后信号
axis([k+1,1000,-10,30]);
title('自适应滤波后信号');

%****************************************************************************************
%  
%                信号Mix_Signal_1 和 Mix_Signal_2  分别作小波滤波
%
%***************************************************************************************

%混合信号 Mix_Signal_1  小波滤波
figure(7);
subplot(4,1,1);
plot(Mix_Signal_1);                                 %Mix_Signal_1 原始信号
axis([0,1000,-5,5]);
title('原始信号 ');

subplot(4,1,2);
[xd,cxd,lxd] = wden(Mix_Signal_1,'sqtwolog','s','one',2,'db3');
plot(xd);                                                 %Mix_Signal_1 小波滤波后信号
axis([0,1000,-5,5]);
title('小波滤波后信号 ');

%混合信号 Mix_Signal_2  小波滤波
subplot(4,1,3);
plot(Mix_Signal_2);                                 %Mix_Signal_2 原始信号
axis([0,1000,-10,30]);
title('原始信号 ');

subplot(4,1,4);
[xd,cxd,lxd] = wden(Mix_Signal_2,'sqtwolog','h','sln',3,'db3');
plot(xd);                                                %Mix_Signal_2 小波滤波后信号
axis([0,1000,-10,30]);
title('小波滤波后信号 ');