zoukankan      html  css  js  c++  java
  • FastFourierTransform (FFT)

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

    FastFourierTransform.h

     1 #pragma once
     2 #include <stdio.h>
     3 #include <math.h>
     4 
     5 #ifndef INCLUDE_FASTTOURIERTRANSFORM
     6 #define INCLUDE_FASTTOURIERTRANSFORM
     7 
     8 /************************************************************************/
     9 /* CFastFourierTransform                                                */
    10 /************************************************************************/
    11 #define PI_2 6.283185F
    12 #define PI   3.1415925F
    13 class CFastFourierTransform
    14 {
    15 private:
    16     float* xre;
    17     float* xim;
    18     float* mag;
    19     float* fftSin;
    20     float* fftCos;
    21     int* fftBr;
    22     int ss;
    23     int ss2;
    24     int nu;
    25     int nu1;
    26 
    27     int BitRev(int j, int nu);
    28     void PrepareFFTTables();
    29 public:
    30     CFastFourierTransform(int pSampleSize);
    31     ~CFastFourierTransform(void);
    32 
    33     float* Calculate(float* pSample, size_t pSampleSize);
    34 };
    35 
    36 #endif

    FastFourierTransform.cpp

      1 #include "FastFourierTransform.h"
      2 
      3 /************************************************************************/
      4 /* CFastFourierTransform                                                */
      5 /************************************************************************/
      6 CFastFourierTransform::CFastFourierTransform(int pSampleSize)
      7 {
      8     xre = NULL;
      9     xim = NULL;
     10     mag = NULL;
     11     fftSin = NULL;
     12     fftCos = NULL;
     13     fftBr = NULL;
     14 
     15     ss = pSampleSize;
     16     ss2 = ss >> 1;
     17     nu = (int) (log((float)ss) / log((float)2));
     18     nu1 = nu - 1;
     19 
     20     xre = new float[ss]; // real part
     21     xim = new float[ss]; // image part
     22     mag = new float[ss2];
     23 
     24     PrepareFFTTables();
     25 }
     26 
     27 CFastFourierTransform::~CFastFourierTransform(void)
     28 {
     29     if(xre != NULL)
     30         delete [] xre;
     31 
     32     if(xim != NULL)
     33         delete [] xim;
     34 
     35     if(mag != NULL)
     36         delete [] mag;
     37 
     38     if(fftSin != NULL)
     39         delete [] fftSin;
     40 
     41     if(fftCos != NULL)
     42         delete [] fftCos;
     43 
     44     if(fftBr != NULL)
     45         delete [] fftBr;
     46 
     47     xre = NULL;
     48     xim = NULL;
     49     mag = NULL;
     50     fftSin = NULL;
     51     fftCos = NULL;
     52     fftBr = NULL;
     53 }
     54 
     55 void CFastFourierTransform::PrepareFFTTables()
     56 {
     57     int n2 = ss2;
     58     int nu1 = nu - 1;
     59 
     60     fftSin = new float[nu * n2];
     61     fftCos = new float[nu * n2];
     62 
     63     int k = 0;
     64     int x = 0;
     65     for (int l = 1; l <= nu; l++) {
     66         while (k < ss) {
     67             for (int i = 1; i <= n2; i++) {
     68                 float p = (float)BitRev(k >> nu1, nu);
     69                 float arg = (PI_2 * p) / (float) ss;
     70                 fftSin[x] = (float) sin(arg);
     71                 fftCos[x] = (float) cos(arg);
     72                 k++;
     73                 x++;
     74             }
     75 
     76             k += n2;
     77         }
     78 
     79         k = 0;
     80         nu1--;
     81         n2 >>= 1;
     82     }
     83 
     84     fftBr = new int[ss];
     85     for (k = 0; k < ss; k++)
     86         fftBr[k] = BitRev(k, nu);
     87 }
     88 
     89 int CFastFourierTransform::BitRev(int j, int nu) {
     90     int j1 = j;
     91     int k = 0;
     92     for (int i = 1; i <= nu; i++) {
     93         int j2 = j1 >> 1;
     94         k = ((k << 1) + j1) - (j2 << 1);
     95         j1 = j2;
     96     }
     97 
     98     return k;
     99 }
    100 
    101 float* CFastFourierTransform::Calculate(float* pSample, size_t pSampleSize) {
    102     int n2 = ss2;
    103     int nu1 = nu - 1;
    104     int wAps = pSampleSize / ss;
    105     size_t a = 0;
    106 
    107     for (size_t b = 0; a < pSampleSize; b++) {
    108         xre[b] = pSample[a];
    109         xim[b] = 0.0F;
    110         a += wAps;
    111     }
    112 
    113     int x = 0;
    114     for (int l = 1; l <= nu; l++) {
    115         for (int k = 0; k < ss; k += n2) {
    116             for (int i = 1; i <= n2; i++) {
    117                 float c = fftCos[x];
    118                 float s = fftSin[x];
    119                 int kn2 = k + n2;
    120                 float tr = xre[kn2] * c + xim[kn2] * s;
    121                 float ti = xim[kn2] * c - xre[kn2] * s;
    122                 xre[kn2] = xre[k] - tr;
    123                 xim[kn2] = xim[k] - ti;
    124                 xre[k] += tr;
    125                 xim[k] += ti;
    126                 k++;
    127                 x++;
    128             }
    129         }
    130 
    131         nu1--;
    132         n2 >>= 1;
    133     }
    134 
    135     for (int k = 0; k < ss; k++) {
    136         int r = fftBr[k];
    137         if (r > k) {
    138             float tr = xre[k];
    139             float ti = xim[k];
    140             xre[k] = xre[r];
    141             xim[k] = xim[r];
    142             xre[r] = tr;
    143             xim[r] = ti;
    144         }
    145     }
    146 
    147     mag[0] = (float) sqrt(xre[0] * xre[0] + xim[0] * xim[0]) / (float) ss;
    148     for (int i = 1; i < ss2; i++)
    149         mag[i] = (2.0F * (float) sqrt(xre[i] * xre[i] + xim[i] * xim[i])) / (float) ss;
    150 
    151     return mag;
    152 }
  • 相关阅读:
    spring框架
    自己来到这个世界的天数
    迭代器
    String
    mybatis-plus条件参数
    Linux常用命令
    web.xml
    log4j.properties
    springmvc.xml
    applicationContext.xml
  • 原文地址:https://www.cnblogs.com/crsky/p/5361941.html
Copyright © 2011-2022 走看看