FastFourierTransform (FFT)

代码来自互联网开源，仅作为收集整理使用。

FastFourierTransform.h

 #pragma once

 #include <stdio.h>

 #include <math.h>

 #ifndef INCLUDE_FASTTOURIERTRANSFORM

 #define INCLUDE_FASTTOURIERTRANSFORM

 /************************************************************************/

 /* CFastFourierTransform                                                */

 /************************************************************************/

 #define PI_2 6.283185F

 #define PI   3.1415925F

 class CFastFourierTransform

 {

 private:

     float* xre;

     float* xim;

     float* mag;

     float* fftSin;

     float* fftCos;

     int* fftBr;

     int ss;

     int ss2;

     int nu;

     int nu1;

     int BitRev(int j, int nu);

     void PrepareFFTTables();

 public:

     CFastFourierTransform(int pSampleSize);

     ~CFastFourierTransform(void);

     float* Calculate(float* pSample, size_t pSampleSize);

 };

 #endif

FastFourierTransform.cpp

 #include "FastFourierTransform.h"

 /************************************************************************/

 /* CFastFourierTransform                                                */

 /************************************************************************/

 CFastFourierTransform::CFastFourierTransform(int pSampleSize)

 {

     xre = NULL;

     xim = NULL;

     mag = NULL;

     fftSin = NULL;

     fftCos = NULL;

     fftBr = NULL;

     ss = pSampleSize;

     ss2 = ss >> ;

     nu = (int) (log((float)ss) / log((float)));

     nu1 = nu - ;

     xre = new float[ss]; // real part

     xim = new float[ss]; // image part

     mag = new float[ss2];

     PrepareFFTTables();

 }

 CFastFourierTransform::~CFastFourierTransform(void)

 {

     if(xre != NULL)

         delete [] xre;

     if(xim != NULL)

         delete [] xim;

     if(mag != NULL)

         delete [] mag;

     if(fftSin != NULL)

         delete [] fftSin;

     if(fftCos != NULL)

         delete [] fftCos;

     if(fftBr != NULL)

         delete [] fftBr;

     xre = NULL;

     xim = NULL;

     mag = NULL;

     fftSin = NULL;

     fftCos = NULL;

     fftBr = NULL;

 }

 void CFastFourierTransform::PrepareFFTTables()

 {

     int n2 = ss2;

     int nu1 = nu - ;

     fftSin = new float[nu * n2];

     fftCos = new float[nu * n2];

     int k = ;

     int x = ;

     for (int l = ; l <= nu; l++) {

         while (k < ss) {

             for (int i = ; i <= n2; i++) {

                 float p = (float)BitRev(k >> nu1, nu);

                 float arg = (PI_2 * p) / (float) ss;

                 fftSin[x] = (float) sin(arg);

                 fftCos[x] = (float) cos(arg);

                 k++;

                 x++;

             }

             k += n2;

         }

         k = ;

         nu1--;

         n2 >>= ;

     }

     fftBr = new int[ss];

     for (k = ; k < ss; k++)

         fftBr[k] = BitRev(k, nu);

 }

 int CFastFourierTransform::BitRev(int j, int nu) {

     int j1 = j;

     int k = ;

     for (int i = ; i <= nu; i++) {

         int j2 = j1 >> ;

         k = ((k << ) + j1) - (j2 << );

         j1 = j2;

     }

     return k;

 }

 float* CFastFourierTransform::Calculate(float* pSample, size_t pSampleSize) {

     int n2 = ss2;

     int nu1 = nu - ;

     int wAps = pSampleSize / ss;

     size_t a = ;

     for (size_t b = ; a < pSampleSize; b++) {

         xre[b] = pSample[a];

         xim[b] = 0.0F;

         a += wAps;

     }

     int x = ;

     for (int l = ; l <= nu; l++) {

         for (int k = ; k < ss; k += n2) {

             for (int i = ; i <= n2; i++) {

                 float c = fftCos[x];

                 float s = fftSin[x];

                 int kn2 = k + n2;

                 float tr = xre[kn2] * c + xim[kn2] * s;

                 float ti = xim[kn2] * c - xre[kn2] * s;

                 xre[kn2] = xre[k] - tr;

                 xim[kn2] = xim[k] - ti;

                 xre[k] += tr;

                 xim[k] += ti;

                 k++;

                 x++;

             }

         }

         nu1--;

         n2 >>= ;

     }

     for (int k = ; k < ss; k++) {

         int r = fftBr[k];

         if (r > k) {

             float tr = xre[k];

             float ti = xim[k];

             xre[k] = xre[r];

             xim[k] = xim[r];

             xre[r] = tr;

             xim[r] = ti;

         }

     }

     mag[] = (float) sqrt(xre[] * xre[] + xim[] * xim[]) / (float) ss;

     for (int i = ; i < ss2; i++)

         mag[i] = (2.0F * (float) sqrt(xre[i] * xre[i] + xim[i] * xim[i])) / (float) ss;

     return mag;

 }