Digital filter with lowpass, then filtfilt: different output

Hello world!

Here I am again dealing with digital filters; I am using the latest MATLAB version (R2020a) with DSP system toolbox.

I am manipulating accelerometers data that need to be filtered by a lowpass filter: satisfactory results have been achieved with function lowpass. Since next acquistions have to pass through the very same filter, I extracted the filter built in lowpass writing

[dataflt, Hd] = lowpass(data, cut_off, fs, 'ImpulseResponse', 'fir');

Subsequent steps would be to apply Hd via filtfilt function to all the other acquisitions: as stated in past MathWorks questions, lowpass inherently uses filtfilt function. In order to double-check this statement, I filtered the data used in lowpass with filter Hd writing

dataflt = filtfilt(Hd,data);

I plotted the two outputs to compare them (in the time domain, figure below): it can be concluded that discrepancies arise at the beginning and at the end of the time history.

Why does such a difference exist? Am I neglecting anything fundamental? Is there any way to solve this problem?

I would prefer to use filtfilt (or similar) instead of lowpass because it seems to me the faster way to filter data, given that the filter is precalculated: is this statement true?

Thank you in advance for your wothwhile help.

Riccardo

Answer 1

[Star Strider]

Occasionally, filter transients can appear, such as those you are seeing. The way I deal with those (to my satisfaction at least) is to add a vector of ones multiplied by the initial value of the vector to the beginning (and similarly at the end if necessary) to the vector. I then remove these elements from the filtered signal. The length of the added elements can vary, I usually use 10 samples. Experiment to get the result you want.

[Riccardo Sorvillo]

First of all, thank you for the fast and clear answer.

What if instead I use lowpass function iteratively for each acquistion? I would avoid the removal of elements and further manipulations of my data. The code would be something like:

% N_acq = number of acquisitions
for ii = 1:N_acq
    datafltr = lowpass(data(:,ii), cut_off, fs, 'ImpulseResponse', 'fir');
end

Is it guaranteed that the same filter would be used at each iteration? And what about computational time?

Thank you for your help.

[Star Strider]

As always, my pleasure!

To use the procedure I discussed, it would not be necessary to remove any elements of the data. Add the ones vectors multiplied by the beginning (and end if necessary) values of your data, then remove those elements from the filtered data. The data itself remains intact. That is simply a way of dealing with the transients.

I believe lowpass (and related functions) use filtfilt, although they may also incorporate a way of setting the initial conditions of filtfilt to deal with the transients. (The leading and trailing vectors approach is how I deal with them.) The problem with running lowpass (or the others) in the loop as posted is that lowpass designs and implements the filter again each time it is called, and that takes time. I would design the filter once, and then use it subsequently for each signal.

The filtfilt function can also filter arrays, the constraint being that the signal lengths all have to be equal in order to create the arrays. That eliminates the need for the loop. (See the documentationn section on x for details.)

There are also ways to design filters using command-line functions. I like the filters that lowpass and the others produce, and I have learned from those designs. I do not always use them, preferring to design my own filters in some situations.

[Riccardo Sorvillo]

Thank you very much for your valuable help: your explaination was clear and solved my problem.

[Star Strider]

As always, my pleasure!

Thank you! I do my best!

Answer 2

[Paul]

lowpass calls filtfilt (or executes code equivalent to filtfilt) only if the underlying filter is IIR. If it's FIR, as in the question by @Riccardo Sorvillo, then lowpass "compensates for the delay introduced by the filter" by using the known delay of a FIR filter that depends on its order.

Example:

rng(100);
fs = 1e3;
t = 0:1/fs:1;

x = [1 2]*sin(2*pi*[50 250]'.*t) + randn(size(t))/10;

% lowpass() with IIR filter, same result as filtfilt()
[ylp,d] = lowpass(x,200,fs,'ImpulseResponse','iir','Steepness',0.5);
yff = filtfilt(d,x);

figure
plot(t,ylp-yff)

% lowpass() with FIR filter, different result than filtfilt()
[ylp,d] = lowpass(x,150,fs);
yff     = filtfilt(d,x);
figure
plot(t,ylp-yff),grid