FFT on vibration analysis

Hi. My task is to analyze acceleration signal for driving comfort. I have to multiply this acceleration signal with appropriate weighting function using values from a chart. This chart has values of weighting function values for "different third octave frequency". So, for this I am converting acceleration signal(71376 values in one column vector)from time to frequency domain by FFT.my code:

% x is the my acceleration signal

 N=size(x); % length of the vector
 Fs = 1230.77; % sampling frequency in Hz

% carrying out FFT on acceleration data%

 xdft = fft(x-mean(x))/N;

%just retaining the positive frequencies%

 xdft = xdft(1:length(x)/2+1);

%generating suitable frequency scale%

 freq = linspace(0,Fs/2,length(x)/2+1);
 plot(freq,abs(xdft));
 xlabel('Hz'); ylabel('Magnitude');

Now this works for me. I get picks at non zero frequency value in my graph.But I do not know if it's giving me correct transformation form time domain to frequency domain. Can somebody guide me how to check it?

Second thing is I need "one third octave frequency". So how can i get it form the frequency scale which I am generating in "freq = linspace(0,Fs/2,length(x)/2+1);" ?

Next thing is when I make fft in "xdft = fft(x-mean(x))/N;" , subtracting Mean form each value means detrending or centering data? Thanks in advance.

Answer 1

[Walter Roberson]

For any constant C, x+C is the same taking x and adding to it the signal C with frequency 0 HZ. Thus, fft(x-mean(x)) affects only the 0 Hz entry of the fft, the first element in the matrix returned by fft. It is therefore, in your terms, centering and not de-trending.

Answer 2

[Brad Humphreys]

Dharmresh,

Go about it this way:

1) Use Pwelch to calulate the Power Spectral density of the acceleration. (Dentrend the time domain data as you have done above first). In the simplest of terms, this is going to give you a scaled fft. The units will be g^2/hz.

2) Next you want to use something called parseval's theroem to convert to unit g (actually grms). A quick look on the file exchange and I don't see an example of how to do this, but try wiki. Alls parseval's tells you to do is to break the PSD up into frequency blocks/bands, integrate the area underneath the PSD in each of those blocks, take the square root of the area, and finally multiply by the frequency width of the band. Thats how grms^2/Hz becomes g.

3) When you section the PSD up into bands, you use the octave 1/3 ocatve bands edges. I have not checked this particular function, but I see an entry on the file exchange at: http://www.mathworks.com/matlabcentral/fileexchange/19033-nth-octave-frequency-bands

4) Then you can apply your scaling factors.

Answer 3

[Stephen Hanly]

Dharmesh, here's a link to a blog I wrote that has a lot of vibration analysis examples and where all the MATLAB code can be downloaded that computes FFT, PSD (power spectral density), and/or spectrogram.

http://blog.mide.com/vibration-analysis-fft-psd-and-spectrogram

Let me know if you have any questions!

[Christian Keelor]

I have some data that I would like to read from an excel file. Time(s) and acceleration m/(s^2). Need to turn time(s)to Hz and m/(s^2) to g^2/Hz. I am not sure how to approach this. Could you help?

[Stephen Hanly]

Christian,

You'll need to convert the acceleration to gs (divide by the gravity constant = 9.81 m/s^2). Then you need to calculate a power spectral density (PSD). This is done first by computing the FFT, then multiplying each frequency component by its complex conjugate.

MATLAB has a nice example here: https://www.mathworks.com/help/signal/ug/power-spectral-density-estimates-using-fft.html

And the code that is included in my vibration analysis blog also shows how this is done on some real world examples.