How can I keep the values passed through the fft the same?

So I have a list of 500000 data points taken from an accelerometer with a sampling frequency of 100000 Hz for 5 seconds and I wish to pass it through an FFT to see it's frequency response and later calculate it's PSD. I've done some reading around and managed to write my own code and this is what I came up with. This is intended to be a vibration analysis so that we can find the natural frequency of the metal pipe we are testing.

Fs = 100000;%sampling frequency

DATA = xlsread('FFTData 1.xlsx');%reading data in from excel file

t = DATA(:,1);%time values from DATA file in Second's

x = DATA(:,4);%Acceleration values read from DATA file in G's

L=length(x);%the number of values read from the data file (500000)

NFFT = 2^nextpow2(L);%better fft makes zero padding.

Y = fft((x - mean(x)),NFFT);%fourier transform of the signal subracting DC Bias Voltage

%Also need to specify NFFT, if we don't, fft() defaults to 512 points

f = Fs/2*linspace(0,1,NFFT/2+1);%span of frequency we want to run through, all the way up to NFFT %value

figure(1);

plot(t,x);%plotting time vs. acceleration

title('Drop Shock Test 1');

xlabel('Time (s)');

ylabel('Acceleration (G''s)');

figure(2);

plot(f,abs((Y(1:NFFT/2+1)))/L);%plotting frequency vs. acceleration. FFT is in Volts*Seconds.

title('Drop Shock Test 1 FFT');

xlabel('Frequency (Hz)');

ylabel('Acceleration (G''s)');

h = spectrum.welch;    % Create a Welch spectral estimator. 

Hpsd = psd(h,x,'Fs',Fs);             % Calculate the PSD 

figure(3);

plot(Hpsd)

My question is, is there a way to get the same y-axis values seen in figure (1), to appear on the y-axis in figure (2) in the frequency domain? I'm also new to asking questions on here (first timer) and I'm not very sure how to attach a file with data points or pictures.

Thanks, Brandon Deal

Answer 1

[Matt J]

You can use ylim() to control the y-axis limits. However, I don't see how it makes sense to have them the same. Since one is a time-domain plot and the other a frequency domain plot, the quantities plotted are completely different in nature. Why expect them to have the same range of values?

[Brandon Deal]

So when passing these units through an FFT (Acceleration in G's) the outcome on the y-axis should be G*s (G's times seconds) correct? I'm trying to analyze my output plots, but I am unsure if they are correctly displaying what it should. I'm going to attach links in a few moments to the graphs I have plotted.

[Brandon Deal]

Here is a link to the plots. I am mostly confused of the high peaks very near 0 Hz. I know at 0 Hz I had the largest peak (due to the DC bias voltage) and I accounted for that and removed it. However, why is my Y values so small in the frequency domain? Any reference, or guidance, you could give me would be awesome.

https://picasaweb.google.com/109805715000182735269/DropShockData?authuser=0&feat=directlink

[Iain]

Er no, the units of the y axis, I believe, are the same.

The range of values is likely to be massively different.

A short-lived high intensity pulse in the time domain gives a very broad, much lower level signature in the frequency domain. A long-lived low intensity pulse in the time domain gives a higher level, narrower signature in the frequency domain.

[Matt J]

So when passing these units through an FFT (Acceleration in G's) the outcome on the y-axis should be G*s (G's times seconds) correct?

That's true of continuous Fourier Transforms. The FFT, however, is not an integral and simply multiplies all the input data by dimensionless exponentials and sums. There is therefore no change in units, unless you would multiply the FFT output by the time sampling interval.

Dividing by L is not appropriate if you're trying to have the same scale as the continuous FT.

[Brandon Deal]

Hmmm interesting. I thought I read somewhere that when passing data through an fft according to theory if you input Volts, then the output should be Volts*seconds. Also, is there any error you could see in my coding? I realized a huge difference in these lines maybe someone could explain.

In this line if I divide by the Length (L) here I get some HUGE numbers on the y-axis fft spectrum.

Y = fft((x - mean(x)),NFFT);%fourier transform of the signal subracting DC Bias Voltage

But when I divided by L in this line of code I got the low values needed in my spectrum. Is it because in the earlier line of code I had complex values when could have messed with my data?

plot(f,abs((Y(1:NFFT/2+1)))/L);%plotting frequency vs. acceleration. FFT is in Volts*Seconds.

Also here is the link to the theory showing output of the fft is Volts*seconds. If you could explain, or reference, why this shouldn't be true, it would greatly be appreciated. Thanks!

http://physics.stackexchange.com/questions/15073/how-does-the-fourier-transform-invert-units

[Matt J]

Also here is the link to the theory showing output of the fft is Volts*seconds. If you could explain, or reference, why this shouldn't be true, it would greatly be appreciated. Thanks!

Your link is talking exclusively about continuous Fourier Transforms. I see only integrals there, not sums.

[Brandon Deal]

That's true of continuous Fourier Transforms. The FFT, however, is not an integral and simply multiplies all the input data by dimensionless exponentials and sums. There is therefore no change in units, unless you would multiply the FFT output by the time sampling interval... Dividing by L is not appropriate if you're trying to have the same scale as the continuous FT.

@Matt So you are saying don't divide anything by L? I should just have the abs(Y) plotted and that would be my complete fft?

And yes I see what you mean about integrals and sums, makes a lot more sense now!

[Matt J]

I'm saying I don't know what Y is trying to measure. If it is supposed to be approximating the continuous Fourier transform of your data, you would do this

 dt=t(2)-t(1); %time sampling interval

 Y = fft((x - mean(x)),NFFT) * dt;

However, only you know what you want. There are many applications of FFTs besides approximating continuous Fourier spectra.

[Brandon Deal]

@Matt Okay and finally my last question and I should be done. For my PSD graph I want to make sure I'm understanding my Y-axis units since I've seen different things. Mostly I've seen that when using the spectrum.welch and then passing through a PSD we get our units as db/Hz. I've also seen that when performing a PSD I can square my fft output, divide by it's length, and then pass through the PSD command. What's the difference between these two essentially?