gtext('d=0.01'), gtext('d=0.05'), gtext('d=0.1') 得到的如下结果图:
图2. 平均方差误差
图3.滤波器系数曲线
系数以时间常数的指数曲线收敛,δ越大,时间常数越小
3. RLS自适应滤波器的应用仿真
从噪声中提取信号
输入信号为:x(k)?A0cos(w0t??0)?b(k) 其中b(k)是附加的白噪声。
xr(k)?A0cos(w0t??0)
应用于RLS自适应滤波器的算法可描述如下:
?1cxx(k?1)x(k) 自适应增益行向量,大小(1,n); g(k)?T?11?x(k)cxx(k?1)x(k)e(k)?y(k)?hT(k?1)x(k) 先验误差
h(k)?h(k?1)?g(k)e(k) 自适应滤波器系数行向量,大小(1,n)
?1?1?1大小(1,n) cxx(k)?cxx(k?1)?g(k)xT(k)cxx(k?1) 输入信号x(k)的自相关转制矩阵,
y(k)?hT(k)x(k) 自适应滤波输出 所研究的滤波器阶数为200,采样周期等于1ms.
程序清单如下:
N=1000; n=200; k=12; Ts=1e-1
b=0.8*randn(1,N); for i=1:N
xr(1,i)=sin(k*2*pi*i/N); x(1,i)=xr(1,i)+b(i); end
Cxx=10000*eye(n); g=zeros(N,n); h=zeros(N,n); e=zeros(1,N); y=zeros(1,N); tr=zeros(1,N);
for i=n+1:N
g(i,:)=(Cxx*x(i-n+1:i)'./(1+x(i-n+1:i)*Cxx*x(i-n+1:i)'))';
e(1,i)=xr(i)-h(i-1,:)*x(i-n+1:i)'; h(i,:)=h(i-1,:)+e(1,i)*g(i,:); Cxx=Cxx-g(i,:)'*x(i-n+1:i)*Cxx; y(1,i)=h(i,:)*x(i-n+1:i)'; tr(1,i)=trace(Cxx); end
figure(1)
plot(0:N-n,x(1,n:N)),grid
title('x(k) input singnal in V') xlabel('Samples')
figure(2)
plot(0:N-n,xr(1,n:N),'r'),grid axis([0 800 -1.2 1.2])
title('xr(k) reference singnal in V') xlabel('Samples')
figure(3)
plot(0:N-n,e(1,n:N)),hold on
plot(0:N-n,y(1,n:N),'r'),hold on grid
title('e(k) error and y(k) output in V') xlabel('Samples')
gtext('e(k)'),gtext('y(k)')
figure(4)
plot(0:N-n,h(n:N,1)),hold on
plot(0:N-n,h(n:N,2),'r'),hold off grid
title('a(n-1) and a(n-2) coeffcients evolution') xlabel('Samples')
figure(5)
num1=fliplr(h(N,:)); sys1=tf(num1,1,Ts); bode(sys1),hold off
title('Synthesized filter') xlabel('Frequency in rad/s')
ylabel('Phase in degree;Module in dB')
figure(6)
semilogy(0:N-n,tr(n:N)),grid title('Cxx matrix trace') xlabel('Samples') 实验结果图如下:
图4,输入信号x(k)
图5 参考信号xr(k)
图6 误差e(k)和输出信号y(k)
图7.滤波器系数a(n-1)和a(n-1)变化曲线
系数的变化曲线在200步时有一个超调,这是由于h(k)向量为零,所以200步以后仅代表x值。获得的滤波器的传递函数也类似于LMS滤波器的传递函数,相应的预测也类似。它的中心频率调整为正弦信号频率,即75rad/s,如下图所示