commit 7729d999239766f058c730418f4517db724a0cef
parent 5476e4e74b14a94e609a3af35ec7f447aaeb80a6
Author: fundamental <mark.d.mccurry@gmail.com>
Date: Sun, 14 Aug 2011 22:55:31 -0400
FFTW: R2C/C2R used over R2R interface
- Due to interleaving, this is marginally faster than old approach
- More importantly, it is easier to read
Diffstat:
5 files changed, 44 insertions(+), 112 deletions(-)
diff --git a/src/DSP/FFTwrapper.cpp b/src/DSP/FFTwrapper.cpp
@@ -22,23 +22,22 @@
#include <cmath>
#include <cassert>
+#include <cstring>
#include "FFTwrapper.h"
FFTwrapper::FFTwrapper(int fftsize_)
{
- fftsize = fftsize_;
- tmpfftdata1 = new fftw_real[fftsize];
- tmpfftdata2 = new fftw_real[fftsize];
- planfftw = fftw_plan_r2r_1d(fftsize,
- tmpfftdata1,
- tmpfftdata1,
- FFTW_R2HC,
- FFTW_ESTIMATE);
- planfftw_inv = fftw_plan_r2r_1d(fftsize,
- tmpfftdata2,
- tmpfftdata2,
- FFTW_HC2R,
+ fftsize = fftsize_;
+ time = new fftw_real[fftsize];
+ fft = new fftw_complex[fftsize];
+ planfftw = fftw_plan_dft_r2c_1d(fftsize,
+ time,
+ fft,
FFTW_ESTIMATE);
+ planfftw_inv = fftw_plan_dft_c2r_1d(fftsize,
+ fft,
+ time,
+ FFTW_ESTIMATE);
}
FFTwrapper::~FFTwrapper()
@@ -46,82 +45,34 @@ FFTwrapper::~FFTwrapper()
fftw_destroy_plan(planfftw);
fftw_destroy_plan(planfftw_inv);
- delete [] tmpfftdata1;
- delete [] tmpfftdata2;
+ delete [] time;
+ delete [] fft;
}
-/*
- * do the Fast Fourier Transform
- */
-void FFTwrapper::smps2freqs(const float *smps, FFTFREQS freqs)
+void FFTwrapper::smps2freqs(const float *smps, fft_t *freqs)
{
- for(int i = 0; i < fftsize; i++)
- tmpfftdata1[i] = smps[i];
- fftw_execute(planfftw);
- for(int i = 0; i < fftsize / 2; i++) {
- freqs.c[i] = tmpfftdata1[i];
- if(i != 0)
- freqs.s[i] = tmpfftdata1[fftsize - i];
- }
- tmpfftdata2[fftsize / 2] = 0.0;
-}
+ //Load data
+ for(int i = 0; i < fftsize; ++i)
+ time[i] = static_cast<double>(smps[i]);
-/*
- * do the Inverse Fast Fourier Transform
- */
-void FFTwrapper::freqs2smps(const FFTFREQS freqs, float *smps)
-{
- tmpfftdata2[fftsize / 2] = 0.0;
- for(int i = 0; i < fftsize / 2; i++) {
- tmpfftdata2[i] = freqs.c[i];
- if(i != 0)
- tmpfftdata2[fftsize - i] = freqs.s[i];
- }
- fftw_execute(planfftw_inv);
- for(int i = 0; i < fftsize; i++)
- smps[i] = tmpfftdata2[i];
-}
+ //DFT
+ fftw_execute(planfftw);
-//only OSCILLGEN SHOULD CALL THIS FOR NOW
-void FFTwrapper::smps2freqs(const float *smps, fft_t *freqs)
-{
- assert(fftsize==OSCIL_SIZE);
- FFTFREQS tmp;
- newFFTFREQS(&tmp, fftsize);
- smps2freqs(smps, tmp);
- for(int i = 0; i < fftsize / 2; ++i)
- freqs[i] = fft_t(tmp.c[i], tmp.s[i]);
- deleteFFTFREQS(&tmp);
+ //Grab data
+ memcpy((void*)freqs, (const void*)fft, fftsize*sizeof(double));
}
void FFTwrapper::freqs2smps(const fft_t *freqs, float *smps)
{
- assert(fftsize==OSCIL_SIZE);
- FFTFREQS tmp;
- newFFTFREQS(&tmp, fftsize);
- for(int i = 0; i < fftsize / 2; ++i) {
- tmp.c[i] = freqs[i].real();
- tmp.s[i] = freqs[i].imag();
- }
- freqs2smps(tmp, smps);
- deleteFFTFREQS(&tmp);
-}
+ //Load data
+ memcpy( (void*)fft, (const void*)freqs, fftsize*sizeof(double));
-void newFFTFREQS(FFTFREQS *f, int size)
-{
- f->c = new float[size];
- f->s = new float[size];
- for(int i = 0; i < size; i++) {
- f->c[i] = 0.0;
- f->s[i] = 0.0;
- }
-}
+ //IDFT
+ fftw_execute(planfftw_inv);
-void deleteFFTFREQS(FFTFREQS *f)
-{
- delete[] f->c;
- delete[] f->s;
- f->c = f->s = NULL;
+ //Grab data
+ for(int i = 0; i < fftsize; ++i)
+ smps[i] = static_cast<float>(time[i]);
}
void FFT_cleanup()
diff --git a/src/DSP/FFTwrapper.h b/src/DSP/FFTwrapper.h
@@ -22,15 +22,10 @@
#ifndef FFT_WRAPPER_H
#define FFT_WRAPPER_H
-
-#include "../globals.h"
-
#include <fftw3.h>
-typedef double fftw_real;
-typedef fftw_plan rfftw_plan;
-
#include <complex>
-typedef std::complex<float> fft_t;
+typedef double fftw_real;
+typedef std::complex<fftw_real> fft_t;
/**A wrapper for the FFTW library (Fast Fourier Transforms)*/
class FFTwrapper
@@ -44,18 +39,15 @@ class FFTwrapper
/**Convert Samples to Frequencies using Fourier Transform
* @param smps Pointer to Samples to be converted; has length fftsize_
* @param freqs Structure FFTFREQS which stores the frequencies*/
- void smps2freqs(const float *smps, FFTFREQS freqs);
- void freqs2smps(const FFTFREQS freqs, float *smps);
void smps2freqs(const float *smps, fft_t *freqs);
void freqs2smps(const fft_t *freqs, float *smps);
private:
- int fftsize;
- fftw_real *tmpfftdata1, *tmpfftdata2;
- rfftw_plan planfftw, planfftw_inv;
+ int fftsize;
+ fftw_real *time;
+ fftw_complex *fft;
+ fftw_plan planfftw, planfftw_inv;
};
-void newFFTFREQS(FFTFREQS *f, int size);
-void deleteFFTFREQS(FFTFREQS *f);
void FFT_cleanup();
#endif
diff --git a/src/Params/PADnoteParameters.cpp b/src/Params/PADnoteParameters.cpp
@@ -469,7 +469,7 @@ void PADnoteParameters::generatespectrum_bandwidthMode(float *spectrum,
float idfreq = i / (float)profilesize - 0.5;
idfreq *= ibw;
int spfreq = (int) (idfreq + ibasefreq);
- float fspfreq = fmod(idfreq + ibasefreq, 1.0);
+ float fspfreq = fmod((double)idfreq + ibasefreq, 1.0);
if(spfreq <= 0)
continue;
if(spfreq >= size - 1)
@@ -577,9 +577,8 @@ void PADnoteParameters::applyparameters(bool lockmutex)
samplemax = 1;
//prepare a BIG FFT stuff
- FFTwrapper *fft = new FFTwrapper(samplesize);
- FFTFREQS fftfreqs;
- newFFTFREQS(&fftfreqs, samplesize / 2);
+ FFTwrapper *fft = new FFTwrapper(samplesize);
+ fft_t *fftfreqs = new fft_t[samplesize / 2];
float adj[samplemax]; //this is used to compute frequency relation to the base frequency
for(int nsample = 0; nsample < samplemax; nsample++)
@@ -603,11 +602,8 @@ void PADnoteParameters::applyparameters(bool lockmutex)
newsample.smp = new float[samplesize + extra_samples];
newsample.smp[0] = 0.0;
- for(int i = 1; i < spectrumsize; i++) { //randomize the phases
- float phase = RND * 6.29;
- fftfreqs.c[i] = spectrum[i] * cos(phase);
- fftfreqs.s[i] = spectrum[i] * sin(phase);
- }
+ for(int i = 1; i < spectrumsize; i++) //randomize the phases
+ fftfreqs[i] = std::polar(spectrum[i], (float)RND * 6.29f);
fft->freqs2smps(fftfreqs, newsample.smp); //that's all; here is the only ifft for the whole sample; no windows are used ;-)
@@ -644,7 +640,7 @@ void PADnoteParameters::applyparameters(bool lockmutex)
newsample.smp = NULL;
}
delete (fft);
- deleteFFTFREQS(&fftfreqs);
+ delete[] fftfreqs;
//delete the additional samples that might exists and are not useful
if(lockmutex) {
diff --git a/src/Synth/OscilGen.cpp b/src/Synth/OscilGen.cpp
@@ -525,7 +525,7 @@ void OscilGen::spectrumadjust()
mag = 1.0;
break;
}
- oscilFFTfreqs[i] = std::polar(mag, phase);
+ oscilFFTfreqs[i] = std::polar<fftw_real>(mag, phase);
}
}
@@ -623,7 +623,7 @@ void OscilGen::prepare()
int k = i * (j + 1);
if(k >= OSCIL_SIZE / 2)
break;
- oscilFFTfreqs[k] += basefuncFFTfreqs[i] * std::polar(hmag[j], hphase[j] * k);
+ oscilFFTfreqs[k] += basefuncFFTfreqs[i] * std::polar<fftw_real>(hmag[j], hphase[j] * k);
}
}
}
@@ -727,7 +727,7 @@ void OscilGen::adaptiveharmonicpostprocess(fft_t *f, int size)
par = 1.0 - pow((1.0 - par), 1.5);
for(int i = 0; i < size; i++) {
- inf[i] = f[i] * par;
+ inf[i] = f[i] * double(par);
f[i] *= (1.0 - par);
}
@@ -846,7 +846,7 @@ short int OscilGen::get(float *smps, float freqHz, int resonance)
if((Prand > 64) && (freqHz >= 0.0) && (!ADvsPAD)) {
const float rnd = PI * pow((Prand - 64.0) / 64.0, 2.0);
for(int i = 1; i < nyquist - 1; i++) //to Nyquist only for AntiAliasing
- outoscilFFTfreqs[i] *= std::polar(1.0f, (float)(rnd * i * RND));
+ outoscilFFTfreqs[i] *= std::polar<fftw_real>(1.0f, (float)(rnd * i * RND));
}
//Harmonic Amplitude Randomness
diff --git a/src/globals.h b/src/globals.h
@@ -25,13 +25,6 @@
#ifndef GLOBALS_H
#define GLOBALS_H
-struct FFTFREQS {
- float *s, *c; //sine and cosine components
-};
-
-extern void newFFTFREQS(FFTFREQS *f, int size);
-extern void deleteFFTFREQS(FFTFREQS *f);
-
/**Sampling rate*/
extern int SAMPLE_RATE;
