Discrete Fourier Transform
In the field of Mathematics, the Discrete Fourier Transform (DFT) is a way to convert a finite list of samples of a function (Complex Numbers) into a list of coefficients of a finite combination of complex sinusoids, ordered by their frequencies.
In Digital Signal Processing, the samples can be signals varying over time. For example, pressure of a sound wave, a radio signal etc.
The following is the basic fundamental equation behind the Discrete Fourier Transform:
In Digital Signal Processing, the samples can be signals varying over time. For example, pressure of a sound wave, a radio signal etc.
The following is the basic fundamental equation behind the Discrete Fourier Transform:
\[\color{silver} { \huge X(k) = \sum_{t=0}^{N-1}x(t)e^{-{i2\pi tk}/N}}\]
In this equation,
x(t) is the input sample/complex number
i is an imaginary number
N is the total number of input samples
t and k are both indices starting from 0 to N
e^(i * angle) is essentially an exponential that can be expressed as a sum of sinusoids according to Euler's Formula:
x(t) is the input sample/complex number
i is an imaginary number
N is the total number of input samples
t and k are both indices starting from 0 to N
e^(i * angle) is essentially an exponential that can be expressed as a sum of sinusoids according to Euler's Formula:
\[\color{silver} { \huge e^{i\theta} = cos(\theta) + isin(\theta)}\]
X(k) is the output complex number that encodes both the amplitude and phase of a sinusoidal component of function x(t) for each k cycle.
My program takes in two command line arguments, the number of samples and the name of the text file containing the complex number samples.
My algorithm is pretty straightforward following the equation directly. One loop to iterate through the k cycles and the second loop to add up all the complex sums.
My algorithm is pretty straightforward following the equation directly. One loop to iterate through the k cycles and the second loop to add up all the complex sums.
Sample Input (N = 8)1+1i -1+1i 1-1i -1-1i -3-1i 4-2i 0-3i 2-6i |
DFT Output3-12i 10.24+5.243i 3+2i 9.071+10.07i -5+4i 1.757-3.243i -9+6i -5.071-4.071i |
The program below takes about 235.013 millseconds to calculate the DFT for 1024 samples. This is not as efficient as we would want it to be. This is when the FFT algorithms solves the efficiency problem.
C++ Source Code
/*****************************************************************************/ /*! \Name Discrete Fourier Transform \Author Deepak Chennakkadan \Input <N> <FILENAME> <N>: Number of input Complex Numbers <FILENAME>: Text file containing the Complex Numbers \Output Prints the DFT Output */ /*****************************************************************************/ #include <fstream> #include <complex> #include <cstdlib> #include <cmath> #include <iostream> #include <iomanip> #include <string> #include <complex> #include <sys/time.h> using namespace std; /***************************************************************************/ /*! \Brief Overloaded operator to format Complex Number output */ /***************************************************************************/ ostream &operator<<(ostream &out, complex<double> c) { if (c.imag() < 0) { out << c.real() << c.imag() << "i"; } else { out << c.real() << "+" << c.imag() << "i"; } return out; } /***************************************************************************/ /*! \Brief Function to calculate elapsed program time */ /***************************************************************************/ float timedifference_msec(struct timeval t0, struct timeval t1) { return (t1.tv_sec - t0.tv_sec) * 1000.0f + (t1.tv_usec - t0.tv_usec) / 1000.0f; } /***************************************************************************/ /*! \Brief Main Function */ /***************************************************************************/ int main(int argc, char *argv[]) { struct timeval t0; struct timeval t1; float elapsed; gettimeofday(&t0, 0); int N; sscanf(argv[1], "%d", &N); const double PI = 4.0f*atan(1.0f); const double TWOPI = 2.0f*PI; complex<double>* inputs; inputs = new complex<double>[N]; ifstream f1(argv[2]); // Loop to read in the Complex Numbers from file for (int i = 0; i < N; i++) { if (f1) { char buf[50]; f1.getline(buf, 50); double d1 = 0, d2 = 0; sscanf(buf, "%lf %lf", &d1, &d2); complex<double> c1(d1, d2); inputs[i] = c1; } } // Loop for each k cycle for (int k = 0; k < N; k++) { double sumreal = 0; double sumimag = 0; int t; // Loop for calculating the sum for (t = 0; t < N; t++) { // e^((2*PI*k*t)/N) double angle = - TWOPI * t * k / N; // Real Sum after complex dot product sumreal += inputs[t].real() * cos(angle) - inputs[t].imag() * sin(angle); // Imaginary Sum after complex dot product sumimag += inputs[t].real() * sin(angle) + inputs[t].imag() * cos(angle); } complex<double> out(sumreal, sumimag); cout << out << endl; } delete[] inputs; // Time Stuff gettimeofday(&t1, 0); elapsed = timedifference_msec(t0, t1); printf("\nCode executed in %f milliseconds.\n", elapsed); return 0; }