compute output, error and coefficients using recursive least squares (rls) algorithm -凯发k8网页登录

use a recursive least squares (rls) filter to identify an unknown system modeled with a lowpass fir filter. compare the frequency responses of the unknown and estimated systems.

filt = dsp.firfilter('numerator',fir1(10,.25));
x = randn(1000,1);
n = 0.01*randn(1000,1);
d = filt(x)   n;

rls = dsp.rlsfilter(11, 'forgettingfactor', 0.98);
[y,e] = rls(x,d);
w = rls.coefficients;

plot(1:1000, [d,y,e]);
title('system identification of an fir filter');
legend('desired', 'output', 'error');
xlabel('time index'); 
ylabel('signal value');

stem([filt.numerator; w].');
legend('actual','estimated');
xlabel('coefficient #'); 
ylabel('coefficient value');

this example demonstrates the rls adaptive algorithm using the inverse system identification model shown here.

cascading the adaptive filter with an unknown filter causes the adaptive filter to converge to a solution that is the inverse of the unknown system.

if the transfer function of the unknown system and the adaptive filter are h(z) and g(z), respectively, the error measured between the desired signal and the signal from the cascaded system reaches its minimum when g(z)×h(z) = 1. for this relation to be true, g(z) must equal 1/h(z), the inverse of the transfer function of the unknown system.

to demonstrate that this is true, create a signal s to input to the cascaded filter pair.

s = randn(3000,1);

in the cascaded filters case, the unknown filter results in a delay in the signal arriving at the summation point after both filters. to prevent the adaptive filter from trying to adapt to a signal it has not yet seen (equivalent to predicting the future), delay the desired signal by 12 samples, which is the order of the unknown system.

generally, you do not know the order of the system you are trying to identify. in that case, delay the desired signal by number of samples equal to half the order of the adaptive filter. delaying the input requires prepending 12 zero-value samples to the input s.

delay = zeros(12,1);
d = [delay; s(1:2988)]; % concatenate the delay and the signal.

you have to keep the desired signal vector d the same length as x, so adjust the signal element count to allow for the delay samples.

although not generally the case, for this example you know the order of the unknown filter, so add a delay equal to the order of the unknown filter.

for the unknown system, use a lowpass, 12th-order fir filter.

filt = dsp.firfilter;
filt.numerator = fir1(12,0.55,'low');

filtering s provides the input data signal for the adaptive algorithm function.

x = filt(s);

to use the rls algorithm, create a dsp.rlsfilter object and set its length, forgettingfactor, and initialinversecovariance properties.

for more information about the input conditions to prepare the rls algorithm object, refer to dsp.rlsfilter.

p0 = 2 * eye(13);
lambda = 0.99;
rls = dsp.rlsfilter(13,'forgettingfactor',lambda,...
   'initialinversecovariance',p0);

this example seeks to develop an inverse solution, you need to be careful about which signal carries the data and which is the desired signal.

earlier examples of adaptive filters use the filtered noise as the desired signal. in this case, the filtered noise (x) carries the unknown system's information. with gaussian distribution and variance of 1, the unfiltered noise d is the desired signal. the code to run this adaptive filter is:

[y,e] = rls(x,d);

where y returns the filtered output and e contains the error signal as the filter adapts to find the inverse of the unknown system.

obtain the estimated coefficients of the rls filter.

b = rls.coefficients;

view the frequency response of the adapted rls filter (inverse system, g(z)) using freqz. the inverse system looks like a highpass filter with linear phase.

freqz(b,1)

view the frequency response of the unknown system, h(z). the response is that of a lowpass filter with a cutoff frequency of 0.55.

freqz(filt.numerator,1)

the result of the cascade of the unknown system and the adapted filter is a compensated system with an extended cutoff frequency of 0.8.

overallcoeffs = conv(filt.numerator,b);
freqz(overallcoeffs,1)

cancel additive noise n added to an unknown system using an rls filter. the rls filter adapts its coefficients until its transfer function matches the transfer function of the unknown system as closely as possible. the difference between the output of the adaptive filter and the output of the unknown system is the error signal e, which represents the additive white noise. minimizing this error signal is the objective of the adaptive filter.

framesize = 100;
niter = 10;
rls = dsp.rlsfilter('length',11,...
    'method','householder rls');
filt = dsp.firfilter('numerator',...
    fir1(10,[.5,.75]));
sinewave = dsp.sinewave('frequency',0.01,...
    'samplerate',1,...
    'samplesperframe',framesize);
scope = timescope('layoutdimensions',[2 1],...
    'numinputports',2, ...
    'timeunits','seconds',...
    'ylimits',[-2.5 2.5], ...
    'bufferlength',2*framesize*niter,...
    'activedisplay',1,...
    'showlegend',true,...
    'channelnames',{'noisy signal'},...
    'activedisplay',2,...
    'showlegend',true,...
    'channelnames',{'error signal'});
for k = 1:niter
    x = randn(framesize,1);
    d = filt(x)   sinewave();
    [y,e] = rls(x,d);
    w = rls.coefficients;
    scope(d,e)
end
release(scope)