一维数据序列滤波的matlab代码，
其实和之前做的图像滤波大同小异，
只是图像的噪声情况复杂得多，
而且是二维的。
做这个主要是手上有个心电的的mens传感器，
蓝牙把数据传过来做一个数据的100Hz左右的带通滤波，
用butterworht做个带通滤波，
再用c语言重构一下。
1 %**************************************************************************************** 2 % 3 % 创建两个信号Mix_Signal_1 和信号 Mix_Signal_2 4 % 5 %*************************************************************************************** 6 7 Fs = 1000; %采样率 8 N = 1000; %采样点数 9 n = 0:N-1; 10 t = 0:1/Fs:1-1/Fs; %时间序列 11 Signal_Original_1 =sin(2*pi*10*t)+sin(2*pi*20*t)+sin(2*pi*30*t); 12 Noise_White_1 = [0.3*randn(1,500), rand(1,500)]; %前500点高斯分部白噪声，后500点均匀分布白噪声 13 Mix_Signal_1 = Signal_Original_1 + Noise_White_1; %构造的混合信号 14 15 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)]; 16 Noise_White_2 = 0.5*randn(1,1000); %高斯白噪声 17 Mix_Signal_2 = Signal_Original_2 + Noise_White_2; %构造的混合信号
一、butterworth filter
1 %**************************************************************************************** 2 % 3 % 信号Mix_Signal_1 和 Mix_Signal_2 分别作巴特沃斯低通滤波。 4 % 5 %*************************************************************************************** 6 7 %混合信号 Mix_Signal_1 巴特沃斯低通滤波 8 figure(1); 9 Wc=2*50/Fs; %截止频率 50Hz 10 [b,a]=butter(4,Wc); 11 Signal_Filter=filter(b,a,Mix_Signal_1); 12 13 subplot(4,1,1); %Mix_Signal_1 原始信号 14 plot(Mix_Signal_1); 15 axis([0,1000,-4,4]); 16 title('原始信号 '); 17 18 subplot(4,1,2); %Mix_Signal_1 低通滤波滤波后信号 19 plot(Signal_Filter); 20 axis([0,1000,-4,4]); 21 title('巴特沃斯低通滤波后信号'); 22 23 %混合信号 Mix_Signal_2 巴特沃斯低通滤波 24 Wc=2*100/Fs; %截止频率 100Hz 25 [b,a]=butter(4,Wc); 26 Signal_Filter=filter(b,a,Mix_Signal_2); 27 28 subplot(4,1,3); %Mix_Signal_2 原始信号 29 plot(Mix_Signal_2); 30 axis([0,1000,-10,30]); 31 title('原始信号 '); 32 33 subplot(4,1,4); %Mix_Signal_2 低通滤波滤波后信号 34 plot(Signal_Filter); 35 axis([0,1000,-10,30]); 36 title('巴特沃斯低通滤波后信号');
二、FIR
1 %**************************************************************************************** 2 % 3 % 信号Mix_Signal_1 和 Mix_Signal_2 分别作FIR低通滤波。 4 % 5 %*************************************************************************************** 6 7 %混合信号 Mix_Signal_1 FIR低通滤波 8 figure(2); 9 F = [0:0.05:0.95]; 10 A = [1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] ; 11 b = firls(20,F,A); 12 Signal_Filter = filter(b,1,Mix_Signal_1); 13 14 subplot(4,1,1); %Mix_Signal_1 原始信号 15 plot(Mix_Signal_1); 16 axis([0,1000,-4,4]); 17 title('原始信号 '); 18 19 subplot(4,1,2); %Mix_Signal_1 FIR低通滤波滤波后信号 20 plot(Signal_Filter); 21 axis([0,1000,-5,5]); 22 title('FIR低通滤波后的信号'); 23 24 %混合信号 Mix_Signal_2 FIR低通滤波 25 F = [0:0.05:0.95]; 26 A = [1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0] ; 27 b = firls(20,F,A); 28 Signal_Filter = filter(b,1,Mix_Signal_2); 29 subplot(4,1,3); %Mix_Signal_2 原始信号 30 plot(Mix_Signal_2); 31 axis([0,1000,-10,30]); 32 title('原始信号 '); 33 34 subplot(4,1,4); %Mix_Signal_2 FIR低通滤波滤波后信号 35 plot(Signal_Filter); 36 axis([0,1000,-10,30]); 37 title('FIR低通滤波后的信号');
三、移动平均滤波
1 %**************************************************************************************** 2 % 3 % 信号Mix_Signal_1 和 Mix_Signal_2 分别作移动平均滤波 4 % 5 %*************************************************************************************** 6 7 %混合信号 Mix_Signal_1 移动平均滤波 8 figure(3); 9 b = [1 1 1 1 1 1]/6; 10 Signal_Filter = filter(b,1,Mix_Signal_1); 11 12 subplot(4,1,1); %Mix_Signal_1 原始信号 13 plot(Mix_Signal_1); 14 axis([0,1000,-4,4]); 15 title('原始信号 '); 16 17 subplot(4,1,2); %Mix_Signal_1 移动平均滤波后信号 18 plot(Signal_Filter); 19 axis([0,1000,-4,4]); 20 title('移动平均滤波后的信号'); 21 22 %混合信号 Mix_Signal_2 移动平均滤波 23 b = [1 1 1 1 1 1]/6; 24 Signal_Filter = filter(b,1,Mix_Signal_2); 25 subplot(4,1,3); %Mix_Signal_2 原始信号 26 plot(Mix_Signal_2); 27 axis([0,1000,-10,30]); 28 title('原始信号 '); 29 30 subplot(4,1,4); %Mix_Signal_2 移动平均滤波后信号 31 plot(Signal_Filter); 32 axis([0,1000,-10,30]); 33 title('移动平均滤波后的信号');
四、中值滤波
1 %**************************************************************************************** 2 % 3 % 信号Mix_Signal_1 和 Mix_Signal_2 分别作中值滤波 4 % 5 %*************************************************************************************** 6 7 %混合信号 Mix_Signal_1 中值滤波 8 figure(4); 9 Signal_Filter=medfilt1(Mix_Signal_1,10); 10 11 subplot(4,1,1); %Mix_Signal_1 原始信号 12 plot(Mix_Signal_1); 13 axis([0,1000,-5,5]); 14 title('原始信号 '); 15 16 subplot(4,1,2); %Mix_Signal_1 中值滤波后信号 17 plot(Signal_Filter); 18 axis([0,1000,-5,5]); 19 title('中值滤波后的信号'); 20 21 %混合信号 Mix_Signal_2 中值滤波 22 Signal_Filter=medfilt1(Mix_Signal_2,10); 23 subplot(4,1,3); %Mix_Signal_2 原始信号 24 plot(Mix_Signal_2); 25 axis([0,1000,-10,30]); 26 title('原始信号 '); 27 28 subplot(4,1,4); %Mix_Signal_2 中值滤波后信号 29 plot(Signal_Filter); 30 axis([0,1000,-10,30]); 31 title('中值滤波后的信号');
五、维纳滤波
1 %**************************************************************************************** 2 % 3 % 信号Mix_Signal_1 和 Mix_Signal_2 分别作维纳滤波 4 % 5 %*************************************************************************************** 6 7 %混合信号 Mix_Signal_1 维纳滤波 8 figure(5); 9 Rxx=xcorr(Mix_Signal_1,Mix_Signal_1); %得到混合信号的自相关函数 10 M=100; %维纳滤波器阶数 11 for i=1:M %得到混合信号的自相关矩阵 12 for j=1:M 13 rxx(i,j)=Rxx(abs(j-i)+N); 14 end 15 end 16 Rxy=xcorr(Mix_Signal_1,Signal_Original_1); %得到混合信号和原信号的互相关函数 17 for i=1:M 18 rxy(i)=Rxy(i+N-1); 19 end %得到混合信号和原信号的互相关向量 20 h = inv(rxx)*rxy'; %得到所要涉及的wiener滤波器系数 21 Signal_Filter=filter(h,1, Mix_Signal_1); %将输入信号通过维纳滤波器 22 23 subplot(4,1,1); %Mix_Signal_1 原始信号 24 plot(Mix_Signal_1); 25 axis([0,1000,-5,5]); 26 title('原始信号 '); 27 28 subplot(4,1,2); %Mix_Signal_1 维纳滤波后信号 29 plot(Signal_Filter); 30 axis([0,1000,-5,5]); 31 title('维纳滤波后的信号'); 32 33 %混合信号 Mix_Signal_2 维纳滤波 34 Rxx=xcorr(Mix_Signal_2,Mix_Signal_2); %得到混合信号的自相关函数 35 M=500; %维纳滤波器阶数 36 for i=1:M %得到混合信号的自相关矩阵 37 for j=1:M 38 rxx(i,j)=Rxx(abs(j-i)+N); 39 end 40 end 41 Rxy=xcorr(Mix_Signal_2,Signal_Original_2); %得到混合信号和原信号的互相关函数 42 for i=1:M 43 rxy(i)=Rxy(i+N-1); 44 end %得到混合信号和原信号的互相关向量 45 h=inv(rxx)*rxy'; %得到所要涉及的wiener滤波器系数 46 Signal_Filter=filter(h,1, Mix_Signal_2); %将输入信号通过维纳滤波器 47 48 subplot(4,1,3); %Mix_Signal_2 原始信号 49 plot(Mix_Signal_2); 50 axis([0,1000,-10,30]); 51 title('原始信号 '); 52 53 subplot(4,1,4); %Mix_Signal_2 维纳滤波后信号 54 plot(Signal_Filter); 55 axis([0,1000,-10,30]); 56 title('维纳滤波后的信号');
六、自适应滤波
1 %**************************************************************************************** 2 % 3 % 信号Mix_Signal_1 和 Mix_Signal_2 分别作自适应滤波 4 % 5 %*************************************************************************************** 6 7 %混合信号 Mix_Signal_1 自适应滤波 8 figure(6); 9 N=1000; %输入信号抽样点数N 10 k=100; %时域抽头LMS算法滤波器阶数 11 u=0.001; %步长因子 12 13 %设置初值 14 yn_1=zeros(1,N); %output signal 15 yn_1(1:k)=Mix_Signal_1(1:k); %将输入信号SignalAddNoise的前k个值作为输出yn_1的前k个值 16 w=zeros(1,k); %设置抽头加权初值 17 e=zeros(1,N); %误差信号 18 19 %用LMS算法迭代滤波 20 for i=(k+1):N 21 XN=Mix_Signal_1((i-k+1):(i)); 22 yn_1(i)=w*XN'; 23 e(i)=Signal_Original_1(i)-yn_1(i); 24 w=w+2*u*e(i)*XN; 25 end 26 27 subplot(4,1,1); 28 plot(Mix_Signal_1); %Mix_Signal_1 原始信号 29 axis([k+1,1000,-4,4]); 30 title('原始信号'); 31 32 subplot(4,1,2); 33 plot(yn_1); %Mix_Signal_1 自适应滤波后信号 34 axis([k+1,1000,-4,4]); 35 title('自适应滤波后信号'); 36 37 %混合信号 Mix_Signal_2 自适应滤波 38 N=1000; %输入信号抽样点数N 39 k=500; %时域抽头LMS算法滤波器阶数 40 u=0.000011; %步长因子 41 42 %设置初值 43 yn_1=zeros(1,N); %output signal 44 yn_1(1:k)=Mix_Signal_2(1:k); %将输入信号SignalAddNoise的前k个值作为输出yn_1的前k个值 45 w=zeros(1,k); %设置抽头加权初值 46 e=zeros(1,N); %误差信号 47 48 %用LMS算法迭代滤波 49 for i=(k+1):N 50 XN=Mix_Signal_2((i-k+1):(i)); 51 yn_1(i)=w*XN'; 52 e(i)=Signal_Original_2(i)-yn_1(i); 53 w=w+2*u*e(i)*XN; 54 end 55 56 subplot(4,1,3); 57 plot(Mix_Signal_2); %Mix_Signal_1 原始信号 58 axis([k+1,1000,-10,30]); 59 title('原始信号'); 60 61 subplot(4,1,4); 62 plot(yn_1); %Mix_Signal_1 自适应滤波后信号 63 axis([k+1,1000,-10,30]); 64 title('自适应滤波后信号');
七、小波滤波
1 %**************************************************************************************** 2 % 3 % 信号Mix_Signal_1 和 Mix_Signal_2 分别作小波滤波 4 % 5 %*************************************************************************************** 6 7 %混合信号 Mix_Signal_1 小波滤波 8 figure(7); 9 subplot(4,1,1); 10 plot(Mix_Signal_1); %Mix_Signal_1 原始信号 11 axis([0,1000,-5,5]); 12 title('原始信号 '); 13 14 subplot(4,1,2); 15 [xd,cxd,lxd] = wden(Mix_Signal_1,'sqtwolog','s','one',2,'db3'); 16 plot(xd); %Mix_Signal_1 小波滤波后信号 17 axis([0,1000,-5,5]); 18 title('小波滤波后信号 '); 19 20 %混合信号 Mix_Signal_2 小波滤波 21 subplot(4,1,3); 22 plot(Mix_Signal_2); %Mix_Signal_2 原始信号 23 axis([0,1000,-10,30]); 24 title('原始信号 '); 25 26 subplot(4,1,4); 27 [xd,cxd,lxd] = wden(Mix_Signal_2,'sqtwolog','h','sln',3,'db3'); 28 plot(xd); %Mix_Signal_2 小波滤波后信号 29 axis([0,1000,-10,30]); 30 title('小波滤波后信号 ');
