FFT: scaling for correct amplitude

What is the general method for which I can find the correct amplitude of a signal's Fourier Transform using FFT()? Here's my code a rectangular pulse as an example...

fs = 4000; t = -25:1/fs:25;  % the sampling frequency and time vector 
g = rectpuls(t, 20);        % g the rectangular pulse
figure; plot(t, g, 'LineWidth', 2); axis([-25 25 -0.5 1.5]); hold on   
xlabel("Time (seconds)"); ylabel("Amplitude"); 
title("Rectangular Pulses, with {\tau} = 10, 20, 40");
N = 2^nextpow2(length(g)); % the # of samples
f = fs*(-N/2:N/2 - 1)/N; % the frequency vector
G = fftshift(fft(g, N));    % G is the result of utilizing FFT on g, then using FFT shift
figure; 
subplot(2,1,1); plot(f, abs(G)/fs, 'LineWidth', 2); 
axis([-0.2*pi 0.2*pi 0 tau(i)]); hold on % 1st subplot - the magnitude spectra
xlabel("Frequency (Hz)"); ylabel("Experimental Amplitude"); 
title(["Magnitude Spectra of Rectangular Pulse, {\tau} = ", num2str(tau(i))]); xticks([-0.6:1/40:0.6])
subplot(2,1,2); plot(f, angle(G)*180/pi, 'LineWidth', 2); 
axis([-0.2*pi 0.2*pi 0 200]); hold on   % 2nd subplot - the angle spectra
xlabel("Frequency (Hz)"); ylabel("Angle in Degrees"); 
title(["Experimental Phase Spectra of Rectangular Pulse, {\tau} = ", num2str(tau(i))]); yticks([180]);

The most important line of code that's bugging me most is this one...

subplot(2,1,1); plot(f, abs(G)/fs, 'LineWidth', 2); 

Now the general method is you divide the magnitude by the length of the original signal's vector like so "abs(G)/length(g)"; however, after applying FFT and FFTshift to a rectangular pulse, if I divide by the length of the original signal's vector, the magnitude is off by a huge factor. It's only correct if I divide by the sampling frequency fs.

So this gets me thinking... which signals (i.e.: rectangular pulse, triangular pulse, sine, sinc, exp()) do I divide the original signal's vector or by the sampling frequency?

Answer 1

[David Goodmanson]

Hi Eric,

It's not the specific details of the signal so much as it is the context. If you are transforming what is basically a continuous wave, say the sum of a few trig fuctions of various frequencies, and you want to know the amplitude of each component, then you divide by N, the number of points in the array as you said. The reasoning is, suppose the time array has index m and array spacing delt, and you have a signal A*exp(2*pi*i*f0*m*delt), i.e. frequency f0 and amplitude A. For the set of frequencies 'fn' in the frequency domain, each amplitude 'an' is

an = sum{m=0,N-1}  A*exp(2*pi*i*f0*m*delt).*exp(-2*pi*i*fn*m*delt)

When fn matches f0, there is a sum of N identical terms of size A which produces N*A. So you have to divide by N to get back to an = A.

[the foregoing ignored some inessential off-by-one indexing issues].

On the other hand, suppose the time domain function is not CW but dies out on each end, such as a gaussian. The die-off takes the cooperation of a lot of closely spaced frequency components, meaning that the fft being used to approximate a continuous integral:

Int g(t)*exp(-2*pi*i*f*t) dt = g(f)

The fft just does an array sum and doesn't know anything about delt, so you need to multiply the sum by dt ~~ delt to approximate the integal. And since delt = 1/fs, you divide by fs.

It's common to use fs as the starting point, with delt = 1/fs and (in the frequency domain) delf = fs/N being derived quantiies. That works, but in my opinion it can sometimes obscure the fundamentals of what is going on.