Search code examples
matlabmachine-learningfilteringsignal-processingbutterworth

How to align original and filtered signal for analysis?


I used butterworth bandpass filter for my acceleration data. Now i want to align/combine both original and smoothed signal for obervation/analysis. The issue is i am trying to align but it is not working. Probably i am doing wrong. I need guidance, how can align them. Below i am posting my output image and code i tried to use to align both signals. The output image i put here is consist of original signal and filtered signal done using butterworth bandpass filter. Thanks

delay = mean(grpdelay(B,A))
tt = Time(1:end-delay);
sn = norm(1:end-delay);

sf = Data_filtered;
sf(1:delay) = [];
plot(tt,sn)
hold on, plot(tt,sf,'-r','linewidth',1.5), hold off
title 'Gyro Filtration'
xlabel('Time (s)'), legend('Original Signal','Filtered Shifted Signal')

output original and filtered signals


Solution

  • I recommend that you filter with filtfilt() because filfilt() has exactly zero delay at all frequencies. It achieves this because it filters the signal forward and then backward. The time delays on the forward pass are cancelled out by the equal but opposite delays when it goes through backward. Note that because you are filtering twice, the final attenuation at each frequency is twice as much (if measured in dB) as what you get from a single pass. Therefore you get 6 dB attenuation at the original cutoff frequency, instead of the usual 3 dB. For most applications, including yours I believe, this is not a problem. If it is a problem, then you can adjust the cutoff of the initial filter definition, so that you get 3 dB attenuation at the desired cutoff after both passes.

    See attached code and figure for a comparison of delay with filt() and filtfilt(). The figure shows that the raw and filtfilt() signals both transition through y=0 at the exact same time, which reflects the absence of any delay with filtfilt(). The figure also shows that the high frequency noise in the raw signal has been removed by the filter.

    Comparison of raw signal and output from filter() and filtfilt(), for a 4th order bandpass Butterworth with Fhpc=0.1 Hz, Flpc=5 Hz.

    If for some reason you do not want to use filtfilt(), then a very good approximation to the delay of a lowpass Butterworth filter of order n, in the passband, is

    Delay = Kn/Fc (Delay in s, Fc=cutoff freq in Hz)

    where Kn is a constant that depends on the filter order, n:

    n=2,3,4,6,8

    Kn=0.225,0.318,0.416,0.615,0.816

    For more orders and the derivation, see my paper "A general solution for the time delay...", J Biomechanics (2007), https://pubmed.ncbi.nlm.nih.gov/16545388/.

    %filterCompare.m  WCR 20210118
    % Question on stack exchange
    clear;
    
    %Make an illustrative input signal
    Fs=200;                     %sampling frequency (Hz)
    N=800;                      %number of points
    Tdur=N/Fs;                  %duration (s)
    Time=(0:(N-1))/Fs;          %time vector
    x=square(Time*2*pi/Tdur)...     %square wave plus high frequency ripple
        +0.2*sin(2*pi*.25*Fs*Time);
    
    %Compute the filter coefficients
    order=4;
    Fhpc=0.1;                   %high pass cutoff frequency (Hz)
    Flpc=5;                     %low pass cutoff frequency (Hz)
    Wn1=[Fhpc,Flpc]*2/Fs;       %normalized cutoff frequencies
    [b,a]=butter(order,Wn1,'bandpass');  %compute filter coefficients
    
    %Filter the data
    y1=filter(b,a,x);         %standard (forward-only) filter
    y2=filtfilt(b,a,x);       %forward and backward filter
    
    %Plot the data
    subplot(2,1,1);
    plot(Time,x,'r.-',Time,y1,'g.-',Time,y2,'b.-');
    ylabel('Amplitude');
    legend('Raw Data','filter()','filtfilt()');
    grid on;
    subplot(2,1,2);
    plot(Time,x,'r.-',Time,y1,'g.-',Time,y2,'b.-');
    xlabel('Time (s)'); ylabel('Amplitude');
    legend('Raw Data','filter()','filtfilt()');
    axis([.45*Tdur .55*Tdur -inf inf]);
    grid on;