Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Simple Discrete Fourier Transform
// Based on: http://nayuki.eigenstate.org/page/how-to-implement-the-discrete-fourier-transform
#include <iostream>
#include <vector>
#include <complex>
#include <cmath>
#include <cstdlib>
using namespace std;
vector<complex<double> > compute_dft(vector<complex<double> > in) {
int n = in.size();
vector<complex<double> > out(n);
// k-th harmonic
for (int k = 0; k < n; k++) {
double sumreal = 0;
double sumimag = 0;
// t-th sample
for (int t = 0; t < n; t++) {
// From Eulera's formula e^ix = cos x + i sin x
double arg = 2 * M_PI * t * k / n;
double cosinus = cos(arg);
double sinus = sin(arg);
sumreal += in[t].real() * cosinus + in[t].imag() * sinus;
sumimag += -in[t].real() * sinus + in[t].imag() * cosinus;
}
out[k] = complex<double>(sumreal, sumimag);
}
return out;
}
template<class T>
void getInput(T &input) {
bool inputOK;
cin.exceptions(ios_base::failbit); // Throw exceptions on fail
do {
try {
inputOK = true;
cin >> input;
} catch (const ios_base::failure& e) {
cout << "Bad input. Please retry." << endl;
inputOK = false;
cin.clear(); // Clear failbit
cin.ignore(); // Clear input
}
} while (!inputOK);
}
// You can check correctness in octave:
// $ octave
// > x = [ 1 + 3i, 2 + 2i, 3 + 1i ]
// > fft(x)
int main() {
int n; // Input samples count
vector<complex<double> > in, out;
cout << "How many samples?" << endl;
getInput(n);
in.resize(n);
for (int i = 0; i < n; i++) {
cout << "Reading sample " << i + 1 << "." << endl;
// Input in format: 1, (1,0), (1,3), etc.
getInput(in[i]);
}
out = compute_dft(in);
cout << scientific;
cout << "n = " << n << endl;
for (int i = 0; i < n; i++) {
cout << i + 1 << ". z = " << out[i] << ", |z| = " << abs(out[i])
<< endl; // Prints harmonic and module
}
return EXIT_SUCCESS;
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment