Converting tremor data in to frequency and filter data

11 views (last 30 days)
Can anyone help in apply fft and filtering to my IMU sensor data for detecting tremor
  4 Comments
Star Strider
Star Strider on 14 Sep 2022
I took each sample per 1sec..
Then the sampling frequency is 1 sample/sec or 1 Hz.
I have changed my code in my Answer and subsequent Comment to reflect that.

Sign in to comment.

Answers (1)

Star Strider
Star Strider on 13 Sep 2022
Edited: Star Strider on 14 Sep 2022
The fft plots are straightforward, however only the acceleration signals make sense. The ‘Gyr’ signals exhibit broadband noise.
How do you want to process them?
One approach —
T1 = readtable('https://www.mathworks.com/matlabcentral/answers/uploaded_files/1124040/labeled%20final%20data.csv')
T1 = 4502×7 table
xAcc yAcc zAcc xGyro yGyro zGyro label _____ _____ ____ _____ _____ _____ _____ 0.02 -0.06 0.99 0.12 -0.67 -0.24 1 0.02 -0.08 0.98 2.93 -2.56 -7.14 1 0.02 -0.09 0.98 7.32 -0.79 1.59 1 -0.01 -0.09 0.99 -3.97 -0.31 0.49 1 -0.01 -0.08 0.99 -3.54 -2.08 -1.71 1 -0.02 -0.09 0.99 -1.22 0.55 -1.95 1 -0.04 -0.08 0.99 -0.18 -0.24 0.12 1 -0.03 -0.09 0.99 0.49 -0.06 0.24 1 -0.04 -0.09 0.99 0.12 -0.43 0.12 1 -0.04 -0.09 0.99 0.24 -0.12 -0.12 1 -0.04 -0.09 0.99 0.24 0.37 0.12 1 -0.04 -0.09 0.99 0.12 0.18 -0.37 1 -0.04 -0.09 0.99 0.06 0.31 0.24 1 -0.04 -0.09 0.99 -0.31 0.06 0.18 1 -0.04 -0.08 0.99 -0.79 -0.24 0.06 1 -0.05 -0.09 0.99 0.12 -0.06 -0.18 1
Acc = T1{:,1:3}
Acc = 4502×3
0.0200 -0.0600 0.9900 0.0200 -0.0800 0.9800 0.0200 -0.0900 0.9800 -0.0100 -0.0900 0.9900 -0.0100 -0.0800 0.9900 -0.0200 -0.0900 0.9900 -0.0400 -0.0800 0.9900 -0.0300 -0.0900 0.9900 -0.0400 -0.0900 0.9900 -0.0400 -0.0900 0.9900
Gyr = T1{:,4:6};
Data = [Acc Gyr];
Fs = 1; % Sampling Frequency
Fn = Fs/2; % Nyquist Frequency
NrSp = size(Data,2); % Number Of Subplots
L = size(Data,1); % Length Of Data Vectors
t = linspace(0, L-1, L)/Fs; % Time Vector
NFFT = 2^nextpow2(L); % For Efficiency
FTData = fft(Data - mean(Data),NFFT)/L; % Fopurier Transform
Fv = linspace(0, 1, NFFT/2-1)*Fn; % Frequency Vector (For Plots)
Iv = 1:numel(Fv); % Index Vector (For Plots)
figure
sp = [1:2:NrSp 2:2:NrSp]; % Order 'subplot' Plots
for k = 1:NrSp
subplot(3,2,sp(k))
plot(Fv, abs(FTData(Iv,k))*2) % Plot Fourier Transforms
grid
xlabel('Frequency')
ylabel('Magnitude')
xlim([0 Fn/20])
% xlim([5 10])
title(sprintf('Column %d',k))
end
sgtitle('Fourier Transform of Data')
figure
sp = [1:2:NrSp 2:2:NrSp];
for k = 1:NrSp
subplot(3,2,sp(k))
plot(t, Data(:,k))
grid
ylim([-3 3])
xlabel('Time')
ylabel('Amplitude')
title(sprintf('Column %d',k))
end
sgtitle('Original Time-Domain Data')
Wp = [0.0001 0.005]/Fn; % Define Passband In Hz, Normalise To (0,pi)
Ws = Wp.*[0.95 1.05]; % Stopband (Normalised)
Rs = 50; % Stopband Ripple (Attenuation) dB
Rp = 1; % Passband Ripple dB
[n,Wn] = ellipord(Wp,Ws,Rp,Rs); % Calculate Filter Order
[z,p,k] = ellip(n,Rp,Rs,Wp); % Design Filter
[sos,g] = zp2sos(z,p,k); % Implement Filter As Second-Order-Section Representation
figure
freqz(sos, 2^16, Fs) % Filter Bode Plot
set(subplot(2,1,1), 'XLim',[0 0.01])
set(subplot(2,1,2), 'XLim',[0 0.01])
Data_Filt = filtfilt(sos,g,Data);
figure
for k = 1:NrSp
subplot(3,2,sp(k))
plot(t, Data_Filt(:,k))
grid
ylim([-1 1]*0.75)
xlabel('Time')
ylabel('Amplitude')
title(sprintf('Column %d',k))
end
sgtitle('Time-Domain Plots of Bandpass (0.0001 - 0.005 Hz) Filtered Data')
EDIT — (14 Sep 2022 at 11:40)
Changed ‘Fs’ and filter passbands to conform to the 1 Hz sampling frequency.
.
  12 Comments

Sign in to comment.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!