The experiment using this program is introduced on "Interface" No.4, CQ publishing Co.,Ltd, 2015. 本プログラムを使った実験は,CQ出版社のインターフェース 2015年4月号で紹介しています.

Dependencies:   DSProcessingIO mbed

Committer:
CQpub0Mikami
Date:
Tue Jul 29 06:52:38 2014 +0000
Revision:
1:655f6a72e47e
Parent:
0:c740f515e0b7
ver.2

Who changed what in which revision?

UserRevisionLine numberNew contents of line
CQpub0Mikami 0:c740f515e0b7 1 //------------------------------------------------------------------------------
CQpub0Mikami 0:c740f515e0b7 2 // FFT class for real data
CQpub0Mikami 0:c740f515e0b7 3 // Copyright (c) 2014 MIKAMI, Naoki, 2014/06/29
CQpub0Mikami 0:c740f515e0b7 4 //------------------------------------------------------------------------------
CQpub0Mikami 0:c740f515e0b7 5
CQpub0Mikami 0:c740f515e0b7 6 #include "fftReal.hpp"
CQpub0Mikami 0:c740f515e0b7 7
CQpub0Mikami 0:c740f515e0b7 8 namespace Mikami
CQpub0Mikami 0:c740f515e0b7 9 {
CQpub0Mikami 0:c740f515e0b7 10 // Constructor
CQpub0Mikami 0:c740f515e0b7 11 FftReal::FftReal(uint16_t n) : N_FFT_(n)
CQpub0Mikami 0:c740f515e0b7 12 {
CQpub0Mikami 0:c740f515e0b7 13 wTable_ = new Complex[n/2];
CQpub0Mikami 0:c740f515e0b7 14 bTable_ = new uint16_t[n/2];
CQpub0Mikami 0:c740f515e0b7 15 u_ = new Complex[n];
CQpub0Mikami 0:c740f515e0b7 16
CQpub0Mikami 0:c740f515e0b7 17 // calculation of twiddle factor
CQpub0Mikami 0:c740f515e0b7 18 Complex arg = Complex(0, -6.283185f/N_FFT_);
CQpub0Mikami 0:c740f515e0b7 19 for (int k=0; k<N_FFT_/2; k++)
CQpub0Mikami 0:c740f515e0b7 20 wTable_[k] = exp(arg*(float)k);
CQpub0Mikami 0:c740f515e0b7 21
CQpub0Mikami 0:c740f515e0b7 22 // for bit reversal
CQpub0Mikami 0:c740f515e0b7 23 int nHalf = N_FFT_/2;
CQpub0Mikami 0:c740f515e0b7 24 bTable_[0] = 0;
CQpub0Mikami 0:c740f515e0b7 25 for (int m=1; m<N_FFT_/2; m*=2)
CQpub0Mikami 0:c740f515e0b7 26 {
CQpub0Mikami 0:c740f515e0b7 27 for (int k=0; k<m; k++)
CQpub0Mikami 0:c740f515e0b7 28 bTable_[k+m] = bTable_[k] + nHalf;
CQpub0Mikami 0:c740f515e0b7 29 nHalf = nHalf/2;
CQpub0Mikami 0:c740f515e0b7 30 }
CQpub0Mikami 0:c740f515e0b7 31 }
CQpub0Mikami 0:c740f515e0b7 32
CQpub0Mikami 0:c740f515e0b7 33 // Destructor
CQpub0Mikami 0:c740f515e0b7 34 FftReal::~FftReal()
CQpub0Mikami 0:c740f515e0b7 35 {
CQpub0Mikami 0:c740f515e0b7 36 delete[] wTable_;
CQpub0Mikami 0:c740f515e0b7 37 delete[] bTable_;
CQpub0Mikami 0:c740f515e0b7 38 delete[] u_;
CQpub0Mikami 0:c740f515e0b7 39 }
CQpub0Mikami 0:c740f515e0b7 40
CQpub0Mikami 0:c740f515e0b7 41 // Execute FFT
CQpub0Mikami 0:c740f515e0b7 42 void FftReal::Execute(const float x[], Complex y[])
CQpub0Mikami 0:c740f515e0b7 43 {
CQpub0Mikami 0:c740f515e0b7 44 for (int n=0; n<N_FFT_; n++) u_[n] = x[n];
CQpub0Mikami 0:c740f515e0b7 45
CQpub0Mikami 0:c740f515e0b7 46 // except for last stage
CQpub0Mikami 0:c740f515e0b7 47 uint16_t nHalf = N_FFT_/2;
CQpub0Mikami 0:c740f515e0b7 48 for (int stg=1; stg<N_FFT_/2; stg*=2)
CQpub0Mikami 0:c740f515e0b7 49 {
CQpub0Mikami 0:c740f515e0b7 50 uint16_t nHalf2 = nHalf*2;
CQpub0Mikami 0:c740f515e0b7 51 for (int kp=0; kp<N_FFT_; kp+=nHalf2)
CQpub0Mikami 0:c740f515e0b7 52 {
CQpub0Mikami 0:c740f515e0b7 53 uint16_t kx = 0;
CQpub0Mikami 0:c740f515e0b7 54 for (int k=kp; k<kp+nHalf; k++)
CQpub0Mikami 0:c740f515e0b7 55 {
CQpub0Mikami 0:c740f515e0b7 56 // Butterfly operation
CQpub0Mikami 0:c740f515e0b7 57 Complex uTmp = u_[k+nHalf];
CQpub0Mikami 0:c740f515e0b7 58 u_[k+nHalf] = (u_[k] - uTmp)*wTable_[kx];
CQpub0Mikami 0:c740f515e0b7 59 u_[k] = u_[k] + uTmp;
CQpub0Mikami 0:c740f515e0b7 60 kx = kx + stg;
CQpub0Mikami 0:c740f515e0b7 61 }
CQpub0Mikami 0:c740f515e0b7 62 }
CQpub0Mikami 0:c740f515e0b7 63 nHalf = nHalf/2;
CQpub0Mikami 0:c740f515e0b7 64 }
CQpub0Mikami 0:c740f515e0b7 65
CQpub0Mikami 0:c740f515e0b7 66 // Last stage
CQpub0Mikami 0:c740f515e0b7 67 y[0] = u_[0] + u_[1];
CQpub0Mikami 0:c740f515e0b7 68 y[N_FFT_/2] = u_[0] - u_[1];
CQpub0Mikami 0:c740f515e0b7 69 for (int k=2; k<N_FFT_; k+=2)
CQpub0Mikami 0:c740f515e0b7 70 u_[k] = u_[k] + u_[k+1];
CQpub0Mikami 0:c740f515e0b7 71
CQpub0Mikami 0:c740f515e0b7 72 // Reorder according to bit reversal
CQpub0Mikami 0:c740f515e0b7 73 for (int k=1; k<N_FFT_/2; k++)
CQpub0Mikami 0:c740f515e0b7 74 y[k] = u_[bTable_[k]];
CQpub0Mikami 0:c740f515e0b7 75 }
CQpub0Mikami 0:c740f515e0b7 76 }