While the DFT transform above can be applied to any complex valued series, in practice for large series it can take considerable time to compute, the time taken being proportional to the square of the number on points in the series. A much faster algorithm has been developed by Cooley and Tukey around 1965 called the FFT (Fast Fourier Transform). The only requirement of the the most popular implementation of this algorithm (Radix-2 Cooley-Tukey) is that the number of points in the series be a power of 2. The computing time for the radix-2 FFT is proportional to

So for example a transform on 1024 points using the DFT takes 10 times longer than using the FFT, a significant speed increase. Note that in reality comparing speeds of various FFT routines is problematic, many of the reported timings have more to do with specific coding methods and their relationship to the hardware and operating system.

• The Fourier transform is linear, that is
  a f(t) + b g(t) ---> a F(f) + b G(f)

  a x_k + b y_k ---> a X_k + b Y_k

• Scaling relationship
  f(t / a) ---> a F(a f)
  f(a t) ---> F(f / a) / a

• Shifting
  f(t + a) ---> F(f) e^{-j 2 pi a f}

• Modulation
  f(t) e^{j 2 pi a t} ---> F(t - a)

• Duality
  X_k ---> (1/N) x_N-k

  Applying the DFT twice results in a scaled, time reversed version of the original series.

• The transform of a constant function is a DC value only.
  

• The transform of a delta function is a constant
  

• The transform of an infinite train of delta functions spaced by T is an infinite train of delta functions spaced by 1/T.

  

• The transform of a cos function is a positive delta at the appropriate positive and negative frequency.

  

• The transform of a sin function is a negative complex delta function at the appropriate positive frequency and a negative complex delta at the appropriate negative frequency.

  

• The transform of a square pulse is a sinc function
  

  More precisely, if f(t) = 1 for |t| < 0.5, and f(t) = 0 otherwise then F(f) = sin(pi f) / (pi f)

• Convolution
  f(t) x g(t) ---> F(f) G(f)

  F(f) x G(f) ---> f(t) g(t)

  x_k x y_k ---> N X_k Y_k

  x_k y_k ---> (1/N) X_k x Y_k

  Multiplication in one domain is equivalent to convolution in the other domain and visa versa. For example the transform of a truncated sin function are two delta functions convolved with a sinc function, a truncated sin function is a sin function multiplied by a square pulse.

• The transform of a triangular pulse is a sinc² function. This can be derived from first principles but is more easily derived by describing the triangular pulse as the convolution of two square pulses and using the convolution-multiplication relationship of the Fourier Transform.

The sampling theorem (often called "Shannons Sampling Theorem") states that a continuous signal must be discretely sampled at least twice the frequency of the highest frequency in the signal.

More precisely, a continuous function f(t) is completely defined by samples every 1/fs (fs is the sample frequency) if the frequency spectrum F(f) is zero for f > fs/2. fs/2 is called the Nyquist frequency and places the limit on the minimum sampling frequency when digitising a continuous sugnal.

If x(k) are the samples of f(t) every 1/fs then f(t) can be exactly reconstructed from these samples, if the sampling theorem has been satisfied, by

where

Normally the signal to be digitised would be appropriately filtered before sampling to remove higher frequency components. If the sampling frequency is not high enough the high frequency components will wrap around and appear in other locations in the discrete spectrum, thus corrupting it.

The key features and consequences of sampling a continuous signal can be shown graphically as follows.

Sample this signal with a sampling frequency f_s, time between samples is 1/f_s. This is equivalent to convolving in the frequency domain by delta function train with a spacing of f_s.

If the sampling frequency is too low the frequency spectrum overlaps, and become corrupted.

Another way to look at this is to consider a sine function sampled twice per period (Nyquist rate). There are other sinusoid functions of higher frequencies that would give exactly the same samples and thus can't be distinguished from the frequency of the original sinusoid.

/*
   Direct fourier transform
*/
int DFT(int dir,int m,double *x1,double *y1)
{
   long i,k;
   double arg;
   double cosarg,sinarg;
   double *x2=NULL,*y2=NULL;

   x2 = malloc(m*sizeof(double));
   y2 = malloc(m*sizeof(double));
   if (x2 == NULL || y2 == NULL)
      return(FALSE);

   for (i=0;i<m;i++) {
      x2[i] = 0;
      y2[i] = 0;
      arg = - dir * 2.0 * 3.141592654 * (double)i / (double)m;
      for (k=0;k<m;k++) {
         cosarg = cos(k * arg);
         sinarg = sin(k * arg);
         x2[i] += (x1[k] * cosarg - y1[k] * sinarg);
         y2[i] += (x1[k] * sinarg + y1[k] * cosarg);
      }
   }

   /* Copy the data back */
   if (dir == 1) {
      for (i=0;i<m;i++) {
         x1[i] = x2[i] / (double)m;
         y1[i] = y2[i] / (double)m;
      }
   } else {
      for (i=0;i<m;i++) {
         x1[i] = x2[i];
         y1[i] = y2[i];
      }
   }

   free(x2);
   free(y2);
   return(TRUE);
}

/*
   This computes an in-place complex-to-complex FFT 
   x and y are the real and imaginary arrays of 2^m points.
   dir =  1 gives forward transform
   dir = -1 gives reverse transform 
*/
short FFT(short int dir,long m,double *x,double *y)
{
   long n,i,i1,j,k,i2,l,l1,l2;
   double c1,c2,tx,ty,t1,t2,u1,u2,z;

   /* Calculate the number of points */
   n = 1;
   for (i=0;i<m;i++) 
      n *= 2;

   /* Do the bit reversal */
   i2 = n >> 1;
   j = 0;
   for (i=0;i<n-1;i++) {
      if (i < j) {
         tx = x[i];
         ty = y[i];
         x[i] = x[j];
         y[i] = y[j];
         x[j] = tx;
         y[j] = ty;
      }
      k = i2;
      while (k <= j) {
         j -= k;
         k >>= 1;
      }
      j += k;
   }

   /* Compute the FFT */
   c1 = -1.0; 
   c2 = 0.0;
   l2 = 1;
   for (l=0;l<m;l++) {
      l1 = l2;
      l2 <<= 1;
      u1 = 1.0; 
      u2 = 0.0;
      for (j=0;j<l1;j++) {
         for (i=j;i<n;i+=l2) {
            i1 = i + l1;
            t1 = u1 * x[i1] - u2 * y[i1];
            t2 = u1 * y[i1] + u2 * x[i1];
            x[i1] = x[i] - t1; 
            y[i1] = y[i] - t2;
            x[i] += t1;
            y[i] += t2;
         }
         z =  u1 * c1 - u2 * c2;
         u2 = u1 * c2 + u2 * c1;
         u1 = z;
      }
      c2 = sqrt((1.0 - c1) / 2.0);
      if (dir == 1) 
         c2 = -c2;
      c1 = sqrt((1.0 + c1) / 2.0);
   }

   /* Scaling for forward transform */
   if (dir == 1) {
      for (i=0;i<n;i++) {
         x[i] /= n;
         y[i] /= n;
      }
   }
   
   return(TRUE);
}