ad_lp.m

最小均方自适应滤波方法举例

%用MATLAB设计自适应低通滤波器
%对噪声滤波
t=0:0.1:399.9
f=50
xs=10*sin(0.5*t)%XS是周期性信号sin(0.5*t)
figure
subplot(2,1,1)
plot(t,xs);grid;
ylabel('幅值');title('要提取的信号xs')
xn=randn(size(t)) %+cos(0.5*t)% XN是干扰信号
subplot(2,1,2)
plot(t,xn);grid;
ylabel('幅值');xlabel('时间');title('加入的噪声信号xn')
xt=xs+xn;%d是期望信号，长度与x相同
dt=xn%randn(size(t));dt是输入参考数组
u=0.001;% u是收敛因子
%,收敛速度和失调量是一对矛盾,要想得
%到较快的收敛速度,可选用大的μ值,这将导致较大的失调量
order=10    % N是FIR滤波器的长度（阶次）
w=[0 0 0 0 0 0 0 0 0 0]% w是估计的FIR滤波器输出
%M=length(xt) %M为权系数长度=2000 ???????????
M=length(dt)
y=zeros(1,M)%y是输出数组
w=zeros(1,order)%滤波器系数
e=zeros(1,M);
r=0.1  %变步长时的约束因子
for n=order:M     %lms算法10~2000
    x1=dt(n:-1:n-order+1);%抽样值
    y(n)=w*x1';
    e(n)=xt(n)-y(n);
    w=w+u*e(n)*x1;% w(n+1)=w(n)+u*e(n)*x(n) %固定步长
    %w=w+u/(r-x1*x1')*e(n)*x1;  %变步长
end
figure
subplot(2,1,1);
plot(t,xt);grid;title('含噪声的信号xt');ylabel('幅值');
subplot(2,1,2);
plot(t,y);grid;title('滤波器输出信号y');ylabel('幅值');
%subplot(3,1,3);
figure
subplot(2,1,1);
plot(t,e);grid;title('最终消噪后的信号e');ylabel('幅值');
subplot(2,1,2);
DD=var(e,1)     %求方差
plot(t,DD);grid;title('误差输出信号var(e)');xlabel('时间');ylabel('幅值');