七种滤波方法的matlab实现和测试

创建两个混合信号,便于更好测试滤波器效果。同时用七中滤波方法测试。
混合信号Mix_Signal_1 = 信号Signal_Original_1+白噪声。
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸



混合信号Mix_Signal_2 = 信号Signal_Original_2+白噪声。
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸



1.巴特沃斯低通滤波器去噪
巴特沃斯滤波器适合用于信号和噪声没有重叠的情况下。下图是巴特沃斯对两个信号的滤波效果。
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸

从图上可以看出巴特沃斯低通滤波器对信号一的滤波效果还是可以的,主要是因为有效的信号最高频率才30Hz,本程序将50Hz以上的信号全部滤除,通过的频率成分中仍然是有白噪声的。

对于信号二,滤波后的信号与没有加噪声的信号相比就有失真了,上升沿和下降沿的高频信号被滤除了。


2.FIR低通滤波器去噪
情况同巴特沃斯低通滤波器相似。滤波后的效果如下:
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸

3. 移动平均滤波去噪
滤波效果如下:
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸

4. 中值滤波去噪
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸

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


5. 维纳滤波去噪

维纳滤波器属于现代滤波器,传统的滤波器只能滤除信号和干扰频带没有重叠的情况,当信号和干扰频带有重叠的时候传统滤波器将无能为力,这时就需要用到现代滤波器,现代滤波器利用信号和干扰的统计特征(如自相关函数、功率谱等)导出一套最佳估值算法,然后用硬件或软件予以实现。
维纳滤波是以均方误差最小( LMS(Least MeanSquare) 为准则的,它根据过去观测值和当前观测值来估计信号的当前值,因此它的解形式是系统的传递函数或单位脉冲响应。
均方误差为:
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸

维纳-霍夫(Wiener-Hopf)方程最小均方误差下的解为:
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸

由于理解不深,对于信号二,没有什么滤波效果
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸

6. 自适应滤波去噪
维纳滤波器参数是固定的,适合于平稳随机信号。卡尔曼滤波器参数是时变的,适合于非平稳随机信号。然而,只有在信号和噪声的统计特性先验已知的情况下,这两种滤波技术才能获得最优滤波。在实际应用中,常常无法得到信号和噪声统计特性的先验知识。在这种情况下,自适应滤波技术能够获得极佳的滤波性能,因而具有很好的应用价值。
自适应滤波的滤波效果如下:
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸

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


7. 小波去噪

首先看一下小波的去噪效果。
七中滤波方法测试matlab实现(转) - 夏衡 - 槿の伊甸

   
   
  1. %****************************************************************************************
  2. %
  3. % 创建两个信号Mix_Signal_1 和信号 Mix_Signal_2
  4. %
  5. %***************************************************************************************
  6. Fs = 1000; %采样率
  7. N = 1000; %采样点数
  8. n = 0:N -1;
  9. t = 0: 1/Fs: 1 -1/Fs; %时间序列
  10. Signal_Original_1 =sin( 2*pi* 10*t)+sin( 2*pi* 20*t)+sin( 2*pi* 30*t);
  11. Noise_White_1 = [ 0.3*randn( 1, 500), rand( 1, 500)]; %前 500点高斯分部白噪声,后 500点均匀分布白噪声
  12. Mix_Signal_1 = Signal_Original_1 + Noise_White_1; %构造的混合信号
  13. 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)];
  14. Noise_White_2 = 0.5*randn( 1, 1000); %高斯白噪声
  15. Mix_Signal_2 = Signal_Original_2 + Noise_White_2; %构造的混合信号
  16. %****************************************************************************************
  17. %
  18. % 信号Mix_Signal_1 和 Mix_Signal_2 分别作巴特沃斯低通滤波。
  19. %
  20. %***************************************************************************************
  21. %混合信号 Mix_Signal_1 巴特沃斯低通滤波
  22. figure( 1);
  23. Wc= 2* 50/Fs; %截止频率 50Hz
  24. [b,a]=butter( 4,Wc);
  25. Signal_Filter=filter(b,a,Mix_Signal_1);
  26. subplot( 4, 1, 1); %Mix_Signal_1 原始信号
  27. plot(Mix_Signal_1);
  28. axis([ 0, 1000, -4, 4]);
  29. title( '原始信号 ');
  30. subplot( 4, 1, 2); %Mix_Signal_1 低通滤波滤波后信号
  31. plot(Signal_Filter);
  32. axis([ 0, 1000, -4, 4]);
  33. title( '巴特沃斯低通滤波后信号');
  34. %混合信号 Mix_Signal_2 巴特沃斯低通滤波
  35. Wc= 2* 100/Fs; %截止频率 100Hz
  36. [b,a]=butter( 4,Wc);
  37. Signal_Filter=filter(b,a,Mix_Signal_2);
  38. subplot( 4, 1, 3); %Mix_Signal_2 原始信号
  39. plot(Mix_Signal_2);
  40. axis([ 0, 1000, -10, 30]);
  41. title( '原始信号 ');
  42. subplot( 4, 1, 4); %Mix_Signal_2 低通滤波滤波后信号
  43. plot(Signal_Filter);
  44. axis([ 0, 1000, -10, 30]);
  45. title( '巴特沃斯低通滤波后信号');
  46. %****************************************************************************************
  47. %
  48. % 信号Mix_Signal_1 和 Mix_Signal_2 分别作FIR低通滤波。
  49. %
  50. %***************************************************************************************
  51. %混合信号 Mix_Signal_1 FIR低通滤波
  52. figure( 2);
  53. F = [ 0: 0.05: 0.95];
  54. A = [ 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] ;
  55. b = firls( 20,F,A);
  56. Signal_Filter = filter(b, 1,Mix_Signal_1);
  57. subplot( 4, 1, 1); %Mix_Signal_1 原始信号
  58. plot(Mix_Signal_1);
  59. axis([ 0, 1000, -4, 4]);
  60. title( '原始信号 ');
  61. subplot( 4, 1, 2); %Mix_Signal_1 FIR低通滤波滤波后信号
  62. plot(Signal_Filter);
  63. axis([ 0, 1000, -5, 5]);
  64. title( 'FIR低通滤波后的信号');
  65. %混合信号 Mix_Signal_2 FIR低通滤波
  66. F = [ 0: 0.05: 0.95];
  67. A = [ 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] ;
  68. b = firls( 20,F,A);
  69. Signal_Filter = filter(b, 1,Mix_Signal_2);
  70. subplot( 4, 1, 3); %Mix_Signal_2 原始信号
  71. plot(Mix_Signal_2);
  72. axis([ 0, 1000, -10, 30]);
  73. title( '原始信号 ');
  74. subplot( 4, 1, 4); %Mix_Signal_2 FIR低通滤波滤波后信号
  75. plot(Signal_Filter);
  76. axis([ 0, 1000, -10, 30]);
  77. title( 'FIR低通滤波后的信号');
  78. %****************************************************************************************
  79. %
  80. % 信号Mix_Signal_1 和 Mix_Signal_2 分别作移动平均滤波
  81. %
  82. %***************************************************************************************
  83. %混合信号 Mix_Signal_1 移动平均滤波
  84. figure( 3);
  85. b = [ 1 1 1 1 1 1]/ 6;
  86. Signal_Filter = filter(b, 1,Mix_Signal_1);
  87. subplot( 4, 1, 1); %Mix_Signal_1 原始信号
  88. plot(Mix_Signal_1);
  89. axis([ 0, 1000, -4, 4]);
  90. title( '原始信号 ');
  91. subplot( 4, 1, 2); %Mix_Signal_1 移动平均滤波后信号
  92. plot(Signal_Filter);
  93. axis([ 0, 1000, -4, 4]);
  94. title( '移动平均滤波后的信号');
  95. %混合信号 Mix_Signal_2 移动平均滤波
  96. b = [ 1 1 1 1 1 1]/ 6;
  97. Signal_Filter = filter(b, 1,Mix_Signal_2);
  98. subplot( 4, 1, 3); %Mix_Signal_2 原始信号
  99. plot(Mix_Signal_2);
  100. axis([ 0, 1000, -10, 30]);
  101. title( '原始信号 ');
  102. subplot( 4, 1, 4); %Mix_Signal_2 移动平均滤波后信号
  103. plot(Signal_Filter);
  104. axis([ 0, 1000, -10, 30]);
  105. title( '移动平均滤波后的信号');
  106. %****************************************************************************************
  107. %
  108. % 信号Mix_Signal_1 和 Mix_Signal_2 分别作中值滤波
  109. %
  110. %***************************************************************************************
  111. %混合信号 Mix_Signal_1 中值滤波
  112. figure( 4);
  113. Signal_Filter=medfilt1(Mix_Signal_1, 10);
  114. subplot( 4, 1, 1); %Mix_Signal_1 原始信号
  115. plot(Mix_Signal_1);
  116. axis([ 0, 1000, -5, 5]);
  117. title( '原始信号 ');
  118. subplot( 4, 1, 2); %Mix_Signal_1 中值滤波后信号
  119. plot(Signal_Filter);
  120. axis([ 0, 1000, -5, 5]);
  121. title( '中值滤波后的信号');
  122. %混合信号 Mix_Signal_2 中值滤波
  123. Signal_Filter=medfilt1(Mix_Signal_2, 10);
  124. subplot( 4, 1, 3); %Mix_Signal_2 原始信号
  125. plot(Mix_Signal_2);
  126. axis([ 0, 1000, -10, 30]);
  127. title( '原始信号 ');
  128. subplot( 4, 1, 4); %Mix_Signal_2 中值滤波后信号
  129. plot(Signal_Filter);
  130. axis([ 0, 1000, -10, 30]);
  131. title( '中值滤波后的信号');
  132. %****************************************************************************************
  133. %
  134. % 信号Mix_Signal_1 和 Mix_Signal_2 分别作维纳滤波
  135. %
  136. %***************************************************************************************
  137. %混合信号 Mix_Signal_1 维纳滤波
  138. figure( 5);
  139. Rxx=xcorr(Mix_Signal_1,Mix_Signal_1); %得到混合信号的自相关函数
  140. M= 100; %维纳滤波器阶数
  141. for i= 1:M %得到混合信号的自相关矩阵
  142. for j= 1:M
  143. rxx(i,j)=Rxx(abs(j-i)+N);
  144. end
  145. end
  146. Rxy=xcorr(Mix_Signal_1,Signal_Original_1); %得到混合信号和原信号的互相关函数
  147. for i= 1:M
  148. rxy(i)=Rxy(i+N -1);
  149. end %得到混合信号和原信号的互相关向量
  150. h = inv(rxx)*rxy '; %得到所要涉及的wiener滤波器系数
  151. Signal_Filter=filter(h,1, Mix_Signal_1); %将输入信号通过维纳滤波器
  152. subplot(4,1,1); %Mix_Signal_1 原始信号
  153. plot(Mix_Signal_1);
  154. axis([0,1000,-5,5]);
  155. title('原始信号 ');
  156. subplot(4,1,2); %Mix_Signal_1 维纳滤波后信号
  157. plot(Signal_Filter);
  158. axis([0,1000,-5,5]);
  159. title('维纳滤波后的信号 ');
  160. %混合信号 Mix_Signal_2 维纳滤波
  161. Rxx=xcorr(Mix_Signal_2,Mix_Signal_2); %得到混合信号的自相关函数
  162. M=500; %维纳滤波器阶数
  163. for i=1:M %得到混合信号的自相关矩阵
  164. for j=1:M
  165. rxx(i,j)=Rxx(abs(j-i)+N);
  166. end
  167. end
  168. Rxy=xcorr(Mix_Signal_2,Signal_Original_2); %得到混合信号和原信号的互相关函数
  169. for i=1:M
  170. rxy(i)=Rxy(i+N-1);
  171. end %得到混合信号和原信号的互相关向量
  172. h=inv(rxx)*rxy'; %得到所要涉及的wiener滤波器系数
  173. Signal_Filter=filter(h, 1, Mix_Signal_2); %将输入信号通过维纳滤波器
  174. subplot( 4, 1, 3); %Mix_Signal_2 原始信号
  175. plot(Mix_Signal_2);
  176. axis([ 0, 1000, -10, 30]);
  177. title( '原始信号 ');
  178. subplot( 4, 1, 4); %Mix_Signal_2 维纳滤波后信号
  179. plot(Signal_Filter);
  180. axis([ 0, 1000, -10, 30]);
  181. title( '维纳滤波后的信号');
  182. %****************************************************************************************
  183. %
  184. % 信号Mix_Signal_1 和 Mix_Signal_2 分别作自适应滤波
  185. %
  186. %***************************************************************************************
  187. %混合信号 Mix_Signal_1 自适应滤波
  188. figure( 6);
  189. N= 1000; %输入信号抽样点数N
  190. k= 100; %时域抽头LMS算法滤波器阶数
  191. u= 0.001; %步长因子
  192. %设置初值
  193. yn_1=zeros( 1,N); %output signal
  194. yn_1( 1:k)=Mix_Signal_1( 1:k); %将输入信号SignalAddNoise的前k个值作为输出yn_1的前k个值
  195. w=zeros( 1,k); %设置抽头加权初值
  196. e=zeros( 1,N); %误差信号
  197. %用LMS算法迭代滤波
  198. for i=(k+ 1):N
  199. XN=Mix_Signal_1((i-k+ 1):(i));
  200. yn_1(i)=w*XN ';
  201. e(i)=Signal_Original_1(i)-yn_1(i);
  202. w=w+2*u*e(i)*XN;
  203. end
  204. subplot(4,1,1);
  205. plot(Mix_Signal_1); %Mix_Signal_1 原始信号
  206. axis([k+1,1000,-4,4]);
  207. title('原始信号 ');
  208. subplot(4,1,2);
  209. plot(yn_1); %Mix_Signal_1 自适应滤波后信号
  210. axis([k+1,1000,-4,4]);
  211. title('自适应滤波后信号 ');
  212. %混合信号 Mix_Signal_2 自适应滤波
  213. N=1000; %输入信号抽样点数N
  214. k=500; %时域抽头LMS算法滤波器阶数
  215. u=0.000011; %步长因子
  216. %设置初值
  217. yn_1=zeros(1,N); %output signal
  218. yn_1(1:k)=Mix_Signal_2(1:k); %将输入信号SignalAddNoise的前k个值作为输出yn_1的前k个值
  219. w=zeros(1,k); %设置抽头加权初值
  220. e=zeros(1,N); %误差信号
  221. %用LMS算法迭代滤波
  222. for i=(k+1):N
  223. XN=Mix_Signal_2((i-k+1):(i));
  224. yn_1(i)=w*XN';
  225. e(i)=Signal_Original_2(i)-yn_1(i);
  226. w=w+ 2*u*e(i)*XN;
  227. end
  228. subplot( 4, 1, 3);
  229. plot(Mix_Signal_2); %Mix_Signal_1 原始信号
  230. axis([k+ 1, 1000, -10, 30]);
  231. title( '原始信号');
  232. subplot( 4, 1, 4);
  233. plot(yn_1); %Mix_Signal_1 自适应滤波后信号
  234. axis([k+ 1, 1000, -10, 30]);
  235. title( '自适应滤波后信号');
  236. %****************************************************************************************
  237. %
  238. % 信号Mix_Signal_1 和 Mix_Signal_2 分别作小波滤波
  239. %
  240. %***************************************************************************************
  241. %混合信号 Mix_Signal_1 小波滤波
  242. figure( 7);
  243. subplot( 4, 1, 1);
  244. plot(Mix_Signal_1); %Mix_Signal_1 原始信号
  245. axis([ 0, 1000, -5, 5]);
  246. title( '原始信号 ');
  247. subplot( 4, 1, 2);
  248. [xd,cxd,lxd] = wden(Mix_Signal_1, 'sqtwolog', 's', 'one', 2, 'db3');
  249. plot(xd); %Mix_Signal_1 小波滤波后信号
  250. axis([ 0, 1000, -5, 5]);
  251. title( '小波滤波后信号 ');
  252. %混合信号 Mix_Signal_2 小波滤波
  253. subplot( 4, 1, 3);
  254. plot(Mix_Signal_2); %Mix_Signal_2 原始信号
  255. axis([ 0, 1000, -10, 30]);
  256. title( '原始信号 ');
  257. subplot( 4, 1, 4);
  258. [xd,cxd,lxd] = wden(Mix_Signal_2, 'sqtwolog', 'h', 'sln', 3, 'db3');
  259. plot(xd); %Mix_Signal_2 小波滤波后信号
  260. axis([ 0, 1000, -10, 30]);
  261. title( '小波滤波后信号 ');



  • 0
    点赞
  • 7
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值