How to filter noise from the following time series data?

I am currently conducting vibrations tests to study the response of a cantilever beam with tip mass excited by harmonic base motion. The displacement response of the shaker armature and beam tip is measured using laser displacement sensors. Below is the following code I have created to post process the data and perform FFT on the time response of the beam and shaker.

close all

clear all

clc

freqs = 10:0.5:30;

orig_data = struct('freq',0,'time',0,'beam',0,'shaker',0);

ss_data = struct('freq',0,'time',0,'beam',0,'shaker',0);

fft_data = struct('freq',0,'f',0,'beam',0,'shaker',0);

rel_data = struct('freq',0,'time',0,'beam',0,'shaker',0);

for i = 1

% import displacement data

file_str = sprintf('f%d_01g_sine_dwell.csv',i); % file name

all_data = readmatrix(file_str,'NumHeaderLines',5,'ExpectedNumVariables',3); % import data

orig_data(i).time = all_data(:,1); % save time data

orig_data(i).beam = all_data(:,2); % save displacement data of beam measurement point (mm)

orig_data(i).shaker = all_data(:,3); % save displacement data of shaker armature measurement point (mm)

orig_data(i).freq = freqs(i); % save excitation frequency (Hz)

% plot original displacement data for shaker armature

figure

plot(orig_data(i).time,orig_data(i).shaker)

xlabel('Time (s)','fontweight','bold')

ylabel('Absolute displacement (mm)','fontweight','bold')

title('Original data for displacement of shaker armature')

set(gca,'FontSize',16)

% plot original displacement data for beam

figure

plot(orig_data(i).time,orig_data(i).beam)

xlabel('Time (s)','fontweight','bold')

ylabel('Absolute displacement (mm)','fontweight','bold')

title('Original data for the displacement of the beam')

set(gca,'FontSize',16)

% extract steady state data between 70 and 90 seconds

sst1 = 70;

sst2 = 90;

ssp1 = find(orig_data(i).time == sst1);

ssp2 = find(orig_data(i).time == sst2);

ss_data(i).time = orig_data(i).time(ssp1:ssp2);

ss_data(i).beam = orig_data(i).beam(ssp1:ssp2);

ss_data(i).shaker = orig_data(i).shaker(ssp1:ssp2);

ss_data(i).freq = freqs(i);

% detrend steady state displacement data to remove DC offset

ss_data(i).beam = detrend(ss_data(i).beam,1);

ss_data(i).shaker = detrend(ss_data(i).shaker,1);

% plot steady state displacement data for shaker armature

figure

plot(ss_data(i).time,ss_data(i).shaker)

xlabel('Time (s)','fontweight','bold')

ylabel('Absolute displacement (mm)','fontweight','bold')

title('Original data for displacement of shaker armature')

set(gca,'FontSize',16)

% plot steady state displacement data for beam

figure

plot(ss_data(i).time,ss_data(i).beam)

xlabel('Time (s)','fontweight','bold')

ylabel('Absolute displacement (mm)','fontweight','bold')

title('Original data for the displacement of the beam')

set(gca,'FontSize',16)

% perform FFT on time series data to determine the frequency spectra of the

% shaker displacement response

Fs = 1000;

T = 1/Fs; % sampling period

L = length(ss_data(i).time(1:end-1));

df = Fs/L;

w = hann(L); % create hanning window

Y = fft(w.*ss_data(i).shaker(1:end-1)); % perform fft on windowed signal

P2 = abs(Y/L); % normalise double sided spectrum

P1 = P2(1:L/2+1); % make single sided spectrum

P1(2:end-1) = 2*P1(2:end-1); % correct amplitude for single sided spectrum

Aw = length(w)/sum(w); % calculate window amplitude correction factor

P1(2:end-1) = Aw*P1(2:end-1); % correction factor for hanning window

f = (Fs*(0:(L/2))/L)';

fft_data(i).f = f;

fft_data(i).shaker = P1;

fft_data(i).freq = freqs(i);

figure

plot(fft_data(i).f,fft_data(i).shaker)

title('FFT of the shaker displacement')

xlabel('f (Hz)','fontweight','bold')

ylabel('Amplitude (mm)','fontweight','bold')

xlim([0 100])

set(gca,'FontSize',16)

set(gca, 'YScale', 'log')

% perform FFT on time series data to determine the frequency spectra of the

% beam displacement response

w = hann(L); % create hanning window

Y = fft(w.*ss_data(i).beam(1:end-1)); % perform fft on windowed signal

P2 = abs(Y/L); % normalise double sided spectrum

P1 = P2(1:L/2+1); % make single sided spectrum

P1(2:end-1) = 2*P1(2:end-1); % correct amplitude for single sided spectrum

Aw = length(w)/sum(w); % calculate window amplitude correction factor

P1(2:end-1) = Aw*P1(2:end-1); % correction factor for hanning window

f = (Fs*(0:(L/2))/L)';

fft_data(i).f = f;

fft_data(i).beam = P1;

fft_data(i).freq = freqs(i);

figure

plot(fft_data(i).f,fft_data(i).beam)

title('FFT of the beam displacement')

xlabel('f (Hz)','fontweight','bold')

ylabel('Amplitude (mm)|','fontweight','bold')

xlim([0 100])

set(gca,'FontSize',16)

set(gca, 'YScale', 'log')

end

Here are the displacement time series for both the shaker and beam.

Here are the frequency spectrums for both the shaker and beam.

The beam was excited by a harmonic excitation of 10Hz and 0.1g zero-to-peak amplitude.The fundamental frequency of 10Hz and the higher harmonics are evident, but the data appears to be corrupted by noise.

My question is about how I can filter the time series signals to remove this noise without losing the fundamental frequency and harmonics of the original signal. I will be very grateful for any guidance. Thank you.

Answer 1

[Star Strider]

The noise is broadband, so any frequency-selective filter will not work to eliminate the noise.

The best options are either wavelet denoising or the Signal Processing Toolbox sgolayfilt function, since either will work for this. An alternative is the smoothdata function with the 'sgolay' option. This may require the Signal Processing Toolbox, however the dependencies are not listed in the documentation, so I cannot determine that.

[Christopher Willett]

@Star Strider Thank you for your suggestions. I have the Signal Processing Toolbox so I will try both methods and hopefully I will have some success. Your help is much appreciated!

[Star Strider]

As always, my pleasure!