discrete low pass in an array

I have an array that is the power measurement. It is recorded each second, but sometimes the sample rate is higher. I want to apply a low-pass filter to that array, to eliminate noise. The signal is an already "discretized" array, and the sample time is a value per second, but it is not always constant. I want to be able to modify some parameters to check the denoise level applied. Is there some code, function or algorithm available to apply this low-pass filter?

P.S. I attached a file 'data.mat' with data.

Answer 1

[Star Strider]

You willl first need to use the resample function to create uniform sampling intervals for your data. Then there are any number of ways you can design a filter for it.

Your ‘test_data.mat’ file contains a (1000x1) double array that I assume is your signal (I have not examined it). If you want help with it, you need to provide a matching array of sampling times.

[Neko Benítez]

Here you have the "new" test_data, in which I attached the timestamp for each value. In this example, there are 122 cases in which the timestamp is greater than 1. Hope you can help me with this. Thank you again.

[Star Strider]

My pleaure.

Try this:

D = load('test_Data.mat');
t = D.test_data(:,1);
s = D.test_data(:,2);

Fs = 1;                                                             % Sampling Frequency
Ts = 1;                                                             % Sampling Interval
[sr,tr] = resample(s, t, Fs);                                       % Resample, Return Resampled Signal & New Time Vector
sre = sr(1:end-2);                                                  % Eliminate End Transient
tre = tr(1:end-2);                                                  % Eliminate End Transient

figure
plot(t, s)
hold on
plot(tre, sre, '--')
hold off
grid
legend('Original Signal', 'Resampled Signal')

L = numel(t);                                                       % Signal Length
Fn = Fs/2;                                                          % Nyquist Frequency
sm = sre - mean(sre);                                               % Subtract Mean
FTs = fft(sm)/L;                                                    % Scaled Fourier Transform
Fv = linspace(0, 1, fix(L/2)+1)*Fn;                                 % Frequency Vector
Iv = 1:numel(Fv);                                                   % Index Vector

figure
plot(Fv, abs(FTs(Iv))*2)
grid
title('Fourier Transform')

Wp = [0.05]/Fn;                                                     % Passband Frequency (Normalised)
Ws = [0.09]/Fn;                                                     % Stopband Frequency (Normalised)
Rp =  1;                                                            % Passband Ripple
Rs = 60;                                                            % Passband Ripple (Attenuation)
[n,Wp] = ellipord(Wp,Ws,Rp,Rs);                                     % Elliptic Order Calculation
[z,p,k] = ellip(n,Rp,Rs,Wp,'low');                                  % Elliptic Filter Design: Zero-Pole-Gain 
[sos,g] = zp2sos(z,p,k);                                            % Second-Order Section For Stability

figure
freqz(sos, 2^16, Fs)                                                % Filter Bode Plot

sre_filt = filtfilt(sos, g, sre);                                   % Filter Signal

figure
subplot(2,1,1)
plot(tre, sre)
grid
title('Resampled Signal')
subplot(2,1,2)
plot(tre, sre_filt, '-')
grid
title('Filtered Resampled Signal')

Experiment with the filter passband and stopband to get the result you want. Remember that for a lowpass filter, the stopband must be greater than the passband. Use the Fourier transform plot to help choose the correct frequencies for the filter. I will of course help as necessary.

[Neko Benítez]

Thank you very much for your code. I really appreciate it. I tried to understand what you are doing, and I have some doubts. But before I deepen with those doubts, I have problems with the function 'resample'. I guess I need a particular Toolbox or similar. Is that right?

Undefined function 'resample' for input arguments of type 'double'.

[Star Strider]

Yes. You need the Signal Processing Toolbox to use resample. It has advantages over interp1 in that it uses an anti-aliasing filter, so it is preferable for signal processing applications.

You can run my code without the resample function. To do so,

replace this code:

[sr,tr] = resample(s, t, Fs);                                       % Resample, Return Resampled Signal & New Time Vector
sre = sr(1:end-2);                                                  % Eliminate End Transient
tre = tr(1:end-2);                                                  % Eliminate End Transient

with this code:

tre = linspace(min(t), max(t), numel(t)*mean(diff(t)));             % New Time Vector
[Ut,ia,ic] = unique(t, 'stable');
sre = interp1(Ut, s(ia), tre);                                      % Resampled Signal

to get the same essential result.

That worked when I ran it.

[Neko Benítez]

Thank you very much again for the help. It was sorted out the problem with the 'resample' function. However, when I run the code, there is another function that also needs the Signal Processing Toolbox: 'ellipord'.

Do you suggest another option instead of 'ellipord'?

Thank you very much again!

[Star Strider]

My pleasure.

You need the Signal Processing Toolbox to do signal processing in MATLAB.

You can possibly design your own linear-phase FIR filter that would emulate the IIR elliptic filter I designed here to use with the filter function. There are several such available online. One that looks promising is: TFilter - Free online FIR filter design. Experiment with different passband and stopband combinations until you get the result you want. Note that FIR filters are generally ‘long’ so usually 40 to 60 elements to get a decent result.

[Neko Benítez]

I don't know how to use that link to design the FIR filter. Do I also need the Signal Processing Toolbox to implement the FIR filter? The link you sent me does not clarify how to do it.

[Star Strider]

There are several websites that can likely design the FIR filter for you (I linked to one of them), or you can design it yourself.

If you are going to do any significant amount of signal processing, I definitely recommend that you get the Signal Processing Toolbox.

If you design a linear-phase FIR filter, the MATLAB filter function is all you need.