FftTest.cpp

00001 /**
00002  *  @file       FftTest.cpp
00003  *  @brief      Calculate 256 point FFT with Hann window
00004  *              Cooley-Tukey algorithm Radix-2 Decimation-In-Time
00005  *  @date       2015.02.24
00006  *  @version    1.0.2
00007  */
00008 #include "FftTest.h"
00009 
00010 /** Constructor
00011 */
00012 FftTest::FftTest()
00013 {
00014     init();
00015 }
00016 
00017 /** Destructor
00018 */
00019 FftTest::~FftTest()
00020 {
00021 }
00022 
00023 /** Apply window fucntion
00024  *  @param      in_data[]   Input data
00025  *  @param      out_data[]  Output data
00026  */
00027 void FftTest::apply_window(float in_data[], float out_data[])
00028 {
00029     int i;
00030     for(i = 0; i < _dat_len; i++) {
00031         out_data[i] = in_data[i] * _window_array[i];
00032     }
00033 }
00034 
00035 /** Calculate FFT
00036  *  @param      in_data     Input data
00037  *  @param      out_re_data Output real data
00038  *  @param      out_im_data Output imaginary data
00039  */
00040 void FftTest::calc_fft(float in_data[], float out_re_data[], float out_im_data[])
00041 {
00042     int i, stage, phase, block, upper_idx, lower_idx;
00043     int32_t block_dist, butterfly_num;
00044     float w_re, w_im, u_re, u_im, temp_re, temp_im, temp_w_re, temp_w_im;
00045 
00046     /* Bit reversal */
00047     for (i = 0; i < _dat_len; i++) {
00048         out_re_data[i] = in_data[_br_idx_array[i]];
00049         out_im_data[i] = 0.0f;
00050     }
00051 
00052     /* Butterfly */
00053     for (stage = 1; stage <= _digit_len; stage++) {
00054         block_dist = 1 << stage;
00055         butterfly_num = block_dist >> 1;
00056 
00057         w_re = 1.0f;
00058         w_im = 0.0f;
00059         u_re = cosf(PI / butterfly_num);
00060         u_im = -sinf(PI / butterfly_num);
00061 
00062         for (phase = 0; phase < butterfly_num; phase++) {
00063             for (block = phase; block < _dat_len; block += block_dist) {
00064                 upper_idx = block;
00065                 lower_idx = upper_idx + butterfly_num;
00066                 temp_re = out_re_data[lower_idx] * w_re - out_im_data[lower_idx] * w_im;
00067                 temp_im = out_re_data[lower_idx] * w_im + out_im_data[lower_idx] * w_re;
00068                 out_re_data[lower_idx] = out_re_data[upper_idx] - temp_re;
00069                 out_im_data[lower_idx] = out_im_data[upper_idx] - temp_im;
00070                 out_re_data[upper_idx] += temp_re;
00071                 out_im_data[upper_idx] += temp_im;
00072             }
00073             temp_w_re = w_re * u_re - w_im * u_im;
00074             temp_w_im = w_re * u_im + w_im * u_re;
00075             w_re = temp_w_re;
00076             w_im = temp_w_im;
00077         }
00078     }
00079 }
00080 
00081 /** Calculate Power spectrum
00082  *  @param      re_data[]   Input real data
00083  *  @param      im_data[]   Input imaginary data
00084  *  @param      pow_data[]  Power data
00085  *  @param      len         Data length
00086  */
00087 void FftTest::calc_power(float re_data[], float im_data[], float pow_data[], int32_t len)
00088 {
00089     int i;
00090     int32_t d_len = (len > _dat_len)? _dat_len : len;
00091     for(i = 0; i < d_len; i++) {
00092         pow_data[i] = re_data[i] * re_data[i] + im_data[i] * im_data[i];
00093     }
00094 }
00095 
00096 /** Calculate Amplitude spectrum
00097  *  @param      pow_data[]  Input power data
00098  *  @param      amp_data[]  Output amplitude data
00099  *  @param      len         Data length
00100  */
00101 void FftTest::calc_amplitude(float pow_data[], float amp_data[], int32_t len)
00102 {
00103     int i;
00104     int32_t d_len = (len > _dat_len)? _dat_len : len;
00105     for(i = 0; i < d_len; i++) {
00106         amp_data[i] = sqrtf(pow_data[i]);
00107     }
00108 }
00109 
00110 /** Normalize Amplitude spectrum
00111  *  @param      pow_data[]  Input power data
00112  *  @param      amp_data[]  Output amplitude data
00113  *  @param      len         Data length
00114  */
00115 void FftTest::norm_amplitude(float amp_data[], int32_t len)
00116 {
00117     int i;
00118     int32_t d_len = (len > _dat_len)? _dat_len : len;
00119     amp_data[0] = amp_data[0] / (float)_dat_len;
00120     for(i = 1; i < d_len; i++) {
00121         amp_data[i] = amp_data[i] * SQRT_2F / (float)_dat_len;
00122     }
00123 }
00124 
00125 /** Initialization
00126  */
00127 void FftTest::init()
00128 {
00129     _digit_len = get_digit_len(DATA_LENGTH);
00130     _dat_len = (1 << _digit_len);
00131 
00132 #ifdef  _LARGE_RAM
00133     _window_array = new float[_dat_len];
00134     _br_idx_array = new int32_t[_dat_len];
00135     set_hann_window();
00136     set_bit_reversal();
00137 #endif  /* LARGE_RAM */
00138 }
00139 
00140 /** Get digit length
00141  *  @param      val         Source value
00142  *  @return                 Digit length
00143  */
00144 int32_t FftTest::get_digit_len(int32_t val)
00145 {
00146     int32_t ret_val = 0;
00147     while( val > 1 ) {
00148         ret_val++;
00149         val >>= 1;
00150     }
00151     return ret_val;
00152 }
00153 
00154 /** Calculate Hann window function
00155  */
00156 void FftTest::set_hann_window()
00157 {
00158 #ifdef  _LARGE_RAM
00159     int i;
00160     for(i = 0; i < _dat_len; i++) {
00161         _window_array[i] = (1.0f - cosf(2.0f * PI * i / (_dat_len - 1.0f))) / 2.0f;
00162     }
00163 #endif  /* LARGE_RAM */
00164 }
00165 
00166 /** Calculate bit reversal index
00167  */
00168 void FftTest::set_bit_reversal()
00169 {
00170 #ifdef  _LARGE_RAM
00171     int i, j;
00172     int32_t reversed_max_bit = (_dat_len >> 1);
00173 
00174     _br_idx_array[0] = 0;
00175     for (i = 1; i < _dat_len; i <<= 1) {
00176         for (j = 0; j < i; j++) {
00177             _br_idx_array[j + i] = _br_idx_array[j] + reversed_max_bit;
00178         }
00179         reversed_max_bit >>= 1;
00180     }
00181 #endif  /* LARGE_RAM */
00182 }
00183 
00184 #ifndef _LARGE_RAM
00185 /* Const */
00186 const float FftTest::_window_array[] = {
00187     0.0f,           0.000151774f,   0.0006070039f,  0.001365413f,   0.002426542f,
00188     0.003789745f,   0.005454196f,   0.007418883f,   0.009682614f,   0.01224402f,
00189     0.01510153f,    0.01825343f,    0.02169779f,    0.02543253f,    0.02945537f,
00190     0.03376389f,    0.03835545f,    0.04322727f,    0.0483764f,     0.05379971f,
00191     0.0594939f,     0.06545553f,    0.07168096f,    0.07816643f,    0.08490799f,
00192     0.09190154f,    0.09914286f,    0.1066275f,     0.114351f,  0.1223086f,
00193     0.1304955f, 0.1389068f, 0.1475372f, 0.1563817f, 0.1654347f, 0.1746908f, 0.1841445f, 0.1937899f, 0.2036212f, 0.2136324f,
00194     0.2238175f, 0.2341703f, 0.2446844f, 0.2553535f, 0.2661712f, 0.2771308f, 0.2882257f, 0.2994492f, 0.3107945f, 0.3222546f,
00195     0.3338226f, 0.3454915f, 0.3572542f, 0.3691036f, 0.3810324f, 0.3930334f, 0.4050995f, 0.4172231f, 0.4293969f, 0.4416136f,
00196     0.4538658f, 0.466146f,  0.4784467f, 0.4907605f, 0.50308f,   0.5153975f, 0.5277057f, 0.5399971f, 0.5522642f, 0.5644996f,
00197     0.5766958f, 0.5888455f, 0.6009412f, 0.6129757f, 0.6249415f, 0.6368315f, 0.6486384f, 0.6603551f, 0.6719745f, 0.6834894f,
00198     0.6948929f, 0.7061782f, 0.7173382f, 0.7283663f, 0.7392558f, 0.75f,      0.7605925f, 0.7710267f, 0.7812964f, 0.7913953f,
00199     0.8013173f, 0.8110564f, 0.8206066f, 0.8299623f, 0.8391176f, 0.848067f,  0.8568051f, 0.8653266f, 0.8736263f, 0.8816991f,
00200     0.8895403f, 0.897145f,  0.9045085f, 0.9116265f, 0.9184946f, 0.9251086f, 0.9314645f, 0.9375585f, 0.9433869f, 0.948946f,
00201     0.9542326f, 0.9592435f, 0.9639755f, 0.9684259f, 0.9725919f, 0.976471f,  0.9800608f, 0.9833592f, 0.9863641f, 0.9890738f,
00202     0.9914865f, 0.9936009f, 0.9954156f, 0.9969296f, 0.9981418f, 0.9990517f, 0.9996585f, 0.999962f,  0.999962f,  0.9996585f,
00203     0.9990517f, 0.9981418f, 0.9969296f, 0.9954156f, 0.9936009f, 0.9914865f, 0.9890738f, 0.9863641f, 0.9833592f, 0.9800608f,
00204     0.976471f,  0.9725919f, 0.9684259f, 0.9639755f, 0.9592435f, 0.9542326f, 0.948946f,  0.9433869f, 0.9375585f, 0.9314645f,
00205     0.9251086f, 0.9184946f, 0.9116265f, 0.9045085f, 0.897145f,  0.8895403f, 0.8816991f, 0.8736263f, 0.8653266f, 0.8568051f,
00206     0.848067f,  0.8391176f, 0.8299623f, 0.8206066f, 0.8110564f, 0.8013173f, 0.7913953f, 0.7812964f, 0.7710267f, 0.7605925f,
00207     0.75f,      0.7392558f, 0.7283663f, 0.7173382f, 0.7061782f, 0.6948929f, 0.6834894f, 0.6719745f, 0.6603551f, 0.6486384f,
00208     0.6368315f, 0.6249415f, 0.6129757f, 0.6009412f, 0.5888455f, 0.5766958f, 0.5644996f, 0.5522642f, 0.5399971f, 0.5277057f,
00209     0.5153975f, 0.50308f,   0.4907605f, 0.4784467f, 0.466146f,  0.4538658f, 0.4416136f, 0.4293969f, 0.4172231f, 0.4050995f,
00210     0.3930334f, 0.3810324f, 0.3691036f, 0.3572542f, 0.3454915f, 0.3338226f, 0.3222546f, 0.3107945f, 0.2994492f, 0.2882257f,
00211     0.2771308f, 0.2661712f, 0.2553535f, 0.2446844f, 0.2341703f, 0.2238175f, 0.2136324f, 0.2036212f, 0.1937899f, 0.1841445f,
00212     0.1746908f, 0.1654347f, 0.1563817f, 0.1475372f, 0.1389068f, 0.1304955f, 0.1223086f, 0.114351f,  0.1066275f, 0.09914286f,
00213     0.09190154f,    0.08490799f,    0.07816643f,    0.07168096f,    0.06545553f,
00214     0.0594939f,     0.05379971f,    0.0483764f,     0.04322727f,    0.03835545f,
00215     0.03376389f,    0.02945537f,    0.02543253f,    0.02169779f,    0.01825343f,
00216     0.01510153f,    0.01224402f,    0.009682614f,   0.007418883f,   0.005454196f,
00217     0.003789745f,   0.002426542f,   0.001365413f,   0.0006070039f,  0.000151774f,
00218     0.0f
00219 };
00220 
00221 const int FftTest::_br_idx_array[] = {
00222     0,      128,    64,     192,    32,     160,    96,     224,    16,     144,
00223     80,     208,    48,     176,    112,    240,    8,      136,    72,     200,
00224     40,     168,    104,    232,    24,     152,    88,     216,    56,     184,
00225     120,    248,    4,      132,    68,     196,    36,     164,    100,    228,
00226     20,     148,    84,     212,    52,     180,    116,    244,    12,     140,
00227     76,     204,    44,     172,    108,    236,    28,     156,    92,     220,
00228     60,     188,    124,    252,    2,      130,    66,     194,    34,     162,
00229     98,     226,    18,     146,    82,     210,    50,     178,    114,    242,
00230     10,     138,    74,     202,    42,     170,    106,    234,    26,     154,
00231     90,     218,    58,     186,    122,    250,    6,      134,    70,     198,
00232     38,     166,    102,    230,    22,     150,    86,     214,    54,     182,
00233     118,    246,    14,     142,    78,     206,    46,     174,    110,    238,
00234     30,     158,    94,     222,    62,     190,    126,    254,    1,      129,
00235     65,     193,    33,     161,    97,     225,    17,     145,    81,     209,
00236     49,     177,    113,    241,    9,      137,    73,     201,    41,     169,
00237     105,    233,    25,     153,    89,     217,    57,     185,    121,    249,
00238     5,      133,    69,     197,    37,     165,    101,    229,    21,     149,
00239     85,     213,    53,     181,    117,    245,    13,     141,    77,     205,
00240     45,     173,    109,    237,    29,     157,    93,     221,    61,     189,
00241     125,    253,    3,      131,    67,     195,    35,     163,    99,     227,
00242     19,     147,    83,     211,    51,     179,    115,    243,    11,     139,
00243     75,     203,    43,     171,    107,    235,    27,     155,    91,     219,
00244     59,     187,    123,    251,    7,      135,    71,     199,    39,     167,
00245     103,    231,    23,     151,    87,     215,    55,     183,    119,    247,
00246     15,     143,    79,     207,    47,     175,    111,    239,    31,     159,
00247     95,     223,    63,     191,    127,    255
00248 };
00249 #endif  /* LARGE_RAM */