Accelerometer data to obtain tremor frequency

Hi all,

I have tri-axial accelerometer data which has 6 channels. the first 3 are XYZ of right and next 3 are XYZ of left. I want to know the tremor frequency from this. this has been sample at Fs=333. what steps should i do to get my output in frequency domain. kindly help me out as im very new to MATLAB.

Thanks in advance!

Attached is the tremor file obtained

[Paras Shaikh]

If anyone of u know about this error? As i cant find my algorithm name in this model.h Continuosly finding this error

Answer 1

[Walter Roberson]

triax = readmatrix('https://www.mathworks.com/matlabcentral/answers/uploaded_files/587501/anj%20off%20rt.txt');
Fs = 333;

nsamp = size(triax,1);
time = (0:nsamp-1) / Fs;
first_half = 1:ceil(nsamp/2);
nhalf = length(first_half)

nhalf = 2509

freqs = linspace(0, Fs/2, nhalf);

XL = triax(:,1);
XL_fft = fft(XL - mean(XL));

subplot(2,1,1);
plot(time, XL); xlabel('s'); ylabel('mag')
subplot(2,1,2);

plot(freqs, abs(XL_fft(first_half))); xlabel('Hz'); ylabel('mag');

I do not promise that the frequencies are correct. Read the first example of the documentation for fft() to see how to create the frequencies.

[Hashvi Sai]

thanks a lot for responding and will definitely check

[Paras Shaikh]

Hey . I have used this fft code. Worked well . Can u plz mention about filtering of obtained frequency. Waiting for reply

[Walter Roberson]

@Paras Shaikh

Sorry, I do not have any experience with filtering the frequencies in this context.

@Star Strider has far more filtering experience than I do.

[Star Strider]

@Paras Shaikh — What do you want to filter? What version of MATLAB do you have? Do you have the Signal Processing Toolbox?

@Walter Roberson — Thank you!

[Paras Shaikh]

I have matlab 2016a I want to filter out my imu data(acc+gyro) But i dont know.. i have to filter my data before apply fft or after. Actually i have no idea about it.. its the part of my project . So how to do?

[Star Strider]

This uses the existing data to present and filter the highest peak (about 7.5 Hz) as well as filter out the D-C offset.

A prototype bandpass filter that should work in R2016a would be implemented as described here —

Data = readmatrix('https://www.mathworks.com/matlabcentral/answers/uploaded_files/587501/anj%20off%20rt.txt')

Data = 5018×6

43 47 45 34 42 40 42 46 46 34 41 41 45 47 47 34 39 42 43 45 47 35 40 42 43 49 47 34 38 41 43 47 47 35 38 42 46 47 47 36 38 43 44 47 48 36 39 41 45 46 48 37 39 42 42 47 47 37 39 40

Fs = 333;                                                   % 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 50])
%     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([20 60])
    xlabel('Time')
    ylabel('Amplitude')
    title(sprintf('Column %d',k))
end
sgtitle('Original Time-Domain Data')

Wp = [7.25 7.75]/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 15])
set(subplot(2,1,2), 'XLim',[0 15])

Data_Filt = filtfilt(sos,g,Data);

figure
for k = 1:NrSp
    subplot(3,2,sp(k))
    plot(t, Data_Filt(:,k))
    grid
    ylim([-10 10])
    xlabel('Time')
    ylabel('Amplitude')
    title(sprintf('Column %d',k))
end
sgtitle('Time-Domain Plots of Bandpass (7.25-7.75 Hz) Filtered Data')

Make appropriate changes to the code and filter characteristics to work with your data.

EDIT — Minor appearance tweaks.

.

[Walter Roberson]

@Paras Shaikh

What do you want the filtering to accomplish ?

Is your X Y Z acceleration data or is it position data? Your use of fft() on the row data implies that it is position data.

[Paras Shaikh]

Hey this is my data

[Paras Shaikh]

I have applied fft on accelerometer data (xyz)

[Star Strider]

@Paras Shaikh — Post this as a new question, along with what you want to do.

Post a link to it here.

[Paras Shaikh]

Star srider if u can do this coding to my data.. that will be great help.. ok im trying to post

[Paras Shaikh]

https://www.mathworks.com/matlabcentral/answers/1804590-converting-tremor-data-in-to-frequency-and-filter-data