my_a_nf.m

直接型自适应滤波器的源码

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 自适应算法实现 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


f1=90;  %信号频率 
N1=256;  % 每周期的采样点数
%fid = fopen('data.txt','wt');   
[b,a1]=ellip(4,0.5,60,16/2048);   %构造椭圆滤波器ellip(阶数，带内震荡幅度，带外衰减，截止频率)



phase_diff=0.001; %相位差
aaa= 0; %固定偏差
mycount=101;%计算100次求平均
mysum=0;
countN=100;
dertasum=0;
mycount=mycount-1;
N=7000;  %采样点数 N1*8个点 大约8个周期
% 构造三个不同频率的正弦信号的叠加作为试验信号 基波频率和各次协波频率定义
f2=3*f1;f3=5*f1;fs=N1*(f1+0.01);% 信号频率 采样频率 N1×（f1＋0.01）准同步采样
pn=0.000001;
%构造左路信号
u=randn(1,N);u=u*sqrt(pn);% 产生 零均值 幅值为0.001的白噪声
u1=u;      %加入直流噪声
%一弧度=pi/180=0.0174532925
faileft=2.00*pi/180.0; % 左路信号的基频对应的相角 
fai2=2*pi*rand(1) ;% 左路信号的二次协波对应的相角
fai3=2*pi*rand(1) ;% 左路信号的三次协波对应的相角
n=0:N-1;
t=2*pi*n/fs;
%左路信号
xl=sin(f1*t+faileft*ones(1,N))+0.1*sin(f2*t+fai2*ones(1,N))+0.01*sin(f3*t+fai3*ones(1,N))+u1;


u=randn(1,N);u=u*sqrt(pn);% 产生零均值 幅值为0.001的白噪声
u2=u;      %加入直流噪声
fairight=(2.0+phase_diff)*pi/180.0;% 右路信号的基频对应的相角
fai2=2*pi*rand(1) ;% 右路信号的二次协波对应的相角
fai3=2*pi*rand(1) ;% 右路信号的三次协波对应的相角
n=0:N-1;
t=2*pi*n/fs;
%右路信号

xl=filter(b,a1,xl);   %对信号进行滤波


for i=1:6000
    xl0(i)=xl(i+999);
end

r0=0.99;
rmax=0.9;
for i=1:6000
r(i)=0;
    a(i)=0;
    k(i)=0;
    q1(i)=0;
    q2(i)=0;
    n(i)=0;
    y(i)=0;
    f(i)=0;
end
p=[1000,0;0,1000];
w=[0,0];
yy=[0,0];
fi=[0,0];
r(1)=0.5;
for i=2:6000
    r(i)=r0*r(i-1)+(1-r0)*rmax;
end

l0=0.99;
lmax=1;
l(1)=0.95;
for i=2:6000
    l(i)=l0*l(i-1)+(1-l0)*lmax;
end
for i=3:5999
    aa=[a(i),0;0,a(i)*a(i)];
    yy=[y(i-1),y(i-2)];
    w1=w(1);
    w2=w(2);
    y(i)=xl0(i)+w1*a(i)*y(i-1)+w2*a(i)*a(i)*y(i-2);
    n(i)=y(i)-w1*y(i-1)-w2*y(i-2);
    fi=[xl0(i-1)-a(i)*n(i-1),xl0(i-2)-a(i)*a(i)*n(i-2)];
    q=p*fi';
    k=q'/(l(i)+fi*q);
    p=(p-q*k)/l(i);
    w=w+n(i)*k;
    

end
for i=1:2048
    sfxl0(i)=xl0(i+2048)-n(i+2048);

end

fl=fft(sfxl0);
omega=acos(w1/(2*sqrt(-w2)));
fre=omega*1024/pi