| /***************************************************************************** |
| |
| FFTRealFixLen.hpp |
| Copyright (c) 2005 Laurent de Soras |
| |
| --- Legal stuff --- |
| |
| This library is free software; you can redistribute it and/or |
| modify it under the terms of the GNU Lesser General Public |
| License as published by the Free Software Foundation; either |
| version 2.1 of the License, or (at your option) any later version. |
| |
| This library is distributed in the hope that it will be useful, |
| but WITHOUT ANY WARRANTY; without even the implied warranty of |
| MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
| Lesser General Public License for more details. |
| |
| You should have received a copy of the GNU Lesser General Public |
| License along with this library; if not, write to the Free Software |
| Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA |
| |
| *Tab=3***********************************************************************/ |
| |
| |
| |
| #if defined (FFTRealFixLen_CURRENT_CODEHEADER) |
| #error Recursive inclusion of FFTRealFixLen code header. |
| #endif |
| #define FFTRealFixLen_CURRENT_CODEHEADER |
| |
| #if ! defined (FFTRealFixLen_CODEHEADER_INCLUDED) |
| #define FFTRealFixLen_CODEHEADER_INCLUDED |
| |
| |
| |
| /*\\\ INCLUDE FILES \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/ |
| |
| #include "def.h" |
| #include "FFTRealPassDirect.h" |
| #include "FFTRealPassInverse.h" |
| #include "FFTRealSelect.h" |
| |
| #include <cassert> |
| #include <cmath> |
| |
| namespace std { } |
| |
| |
| |
| /*\\\ PUBLIC \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/ |
| |
| |
| |
| template <int LL2> |
| FFTRealFixLen <LL2>::FFTRealFixLen () |
| : _buffer (FFT_LEN) |
| , _br_data (BR_ARR_SIZE) |
| , _trigo_data (TRIGO_TABLE_ARR_SIZE) |
| , _trigo_osc () |
| { |
| build_br_lut (); |
| build_trigo_lut (); |
| build_trigo_osc (); |
| } |
| |
| |
| |
| template <int LL2> |
| long FFTRealFixLen <LL2>::get_length () const |
| { |
| return (FFT_LEN); |
| } |
| |
| |
| |
| // General case |
| template <int LL2> |
| void FFTRealFixLen <LL2>::do_fft (DataType f [], const DataType x []) |
| { |
| assert (f != 0); |
| assert (x != 0); |
| assert (x != f); |
| assert (FFT_LEN_L2 >= 3); |
| |
| // Do the transform in several passes |
| const DataType * cos_ptr = &_trigo_data [0]; |
| const long * br_ptr = &_br_data [0]; |
| |
| FFTRealPassDirect <FFT_LEN_L2 - 1>::process ( |
| FFT_LEN, |
| f, |
| &_buffer [0], |
| x, |
| cos_ptr, |
| TRIGO_TABLE_ARR_SIZE, |
| br_ptr, |
| &_trigo_osc [0] |
| ); |
| } |
| |
| // 4-point FFT |
| template <> |
| void FFTRealFixLen <2>::do_fft (DataType f [], const DataType x []) |
| { |
| assert (f != 0); |
| assert (x != 0); |
| assert (x != f); |
| |
| f [1] = x [0] - x [2]; |
| f [3] = x [1] - x [3]; |
| |
| const DataType b_0 = x [0] + x [2]; |
| const DataType b_2 = x [1] + x [3]; |
| |
| f [0] = b_0 + b_2; |
| f [2] = b_0 - b_2; |
| } |
| |
| // 2-point FFT |
| template <> |
| void FFTRealFixLen <1>::do_fft (DataType f [], const DataType x []) |
| { |
| assert (f != 0); |
| assert (x != 0); |
| assert (x != f); |
| |
| f [0] = x [0] + x [1]; |
| f [1] = x [0] - x [1]; |
| } |
| |
| // 1-point FFT |
| template <> |
| void FFTRealFixLen <0>::do_fft (DataType f [], const DataType x []) |
| { |
| assert (f != 0); |
| assert (x != 0); |
| |
| f [0] = x [0]; |
| } |
| |
| |
| |
| // General case |
| template <int LL2> |
| void FFTRealFixLen <LL2>::do_ifft (const DataType f [], DataType x []) |
| { |
| assert (f != 0); |
| assert (x != 0); |
| assert (x != f); |
| assert (FFT_LEN_L2 >= 3); |
| |
| // Do the transform in several passes |
| DataType * s_ptr = |
| FFTRealSelect <FFT_LEN_L2 & 1>::sel_bin (&_buffer [0], x); |
| DataType * d_ptr = |
| FFTRealSelect <FFT_LEN_L2 & 1>::sel_bin (x, &_buffer [0]); |
| const DataType * cos_ptr = &_trigo_data [0]; |
| const long * br_ptr = &_br_data [0]; |
| |
| FFTRealPassInverse <FFT_LEN_L2 - 1>::process ( |
| FFT_LEN, |
| d_ptr, |
| s_ptr, |
| f, |
| cos_ptr, |
| TRIGO_TABLE_ARR_SIZE, |
| br_ptr, |
| &_trigo_osc [0] |
| ); |
| } |
| |
| // 4-point IFFT |
| template <> |
| void FFTRealFixLen <2>::do_ifft (const DataType f [], DataType x []) |
| { |
| assert (f != 0); |
| assert (x != 0); |
| assert (x != f); |
| |
| const DataType b_0 = f [0] + f [2]; |
| const DataType b_2 = f [0] - f [2]; |
| |
| x [0] = b_0 + f [1] * 2; |
| x [2] = b_0 - f [1] * 2; |
| x [1] = b_2 + f [3] * 2; |
| x [3] = b_2 - f [3] * 2; |
| } |
| |
| // 2-point IFFT |
| template <> |
| void FFTRealFixLen <1>::do_ifft (const DataType f [], DataType x []) |
| { |
| assert (f != 0); |
| assert (x != 0); |
| assert (x != f); |
| |
| x [0] = f [0] + f [1]; |
| x [1] = f [0] - f [1]; |
| } |
| |
| // 1-point IFFT |
| template <> |
| void FFTRealFixLen <0>::do_ifft (const DataType f [], DataType x []) |
| { |
| assert (f != 0); |
| assert (x != 0); |
| assert (x != f); |
| |
| x [0] = f [0]; |
| } |
| |
| |
| |
| |
| template <int LL2> |
| void FFTRealFixLen <LL2>::rescale (DataType x []) const |
| { |
| assert (x != 0); |
| |
| const DataType mul = DataType (1.0 / FFT_LEN); |
| |
| if (FFT_LEN < 4) |
| { |
| long i = FFT_LEN - 1; |
| do |
| { |
| x [i] *= mul; |
| --i; |
| } |
| while (i >= 0); |
| } |
| |
| else |
| { |
| assert ((FFT_LEN & 3) == 0); |
| |
| // Could be optimized with SIMD instruction sets (needs alignment check) |
| long i = FFT_LEN - 4; |
| do |
| { |
| x [i + 0] *= mul; |
| x [i + 1] *= mul; |
| x [i + 2] *= mul; |
| x [i + 3] *= mul; |
| i -= 4; |
| } |
| while (i >= 0); |
| } |
| } |
| |
| |
| |
| /*\\\ PROTECTED \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/ |
| |
| |
| |
| /*\\\ PRIVATE \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/ |
| |
| |
| |
| template <int LL2> |
| void FFTRealFixLen <LL2>::build_br_lut () |
| { |
| _br_data [0] = 0; |
| for (long cnt = 1; cnt < BR_ARR_SIZE; ++cnt) |
| { |
| long index = cnt << 2; |
| long br_index = 0; |
| |
| int bit_cnt = FFT_LEN_L2; |
| do |
| { |
| br_index <<= 1; |
| br_index += (index & 1); |
| index >>= 1; |
| |
| -- bit_cnt; |
| } |
| while (bit_cnt > 0); |
| |
| _br_data [cnt] = br_index; |
| } |
| } |
| |
| |
| |
| template <int LL2> |
| void FFTRealFixLen <LL2>::build_trigo_lut () |
| { |
| const double mul = (0.5 * PI) / TRIGO_TABLE_ARR_SIZE; |
| for (long i = 0; i < TRIGO_TABLE_ARR_SIZE; ++ i) |
| { |
| using namespace std; |
| |
| _trigo_data [i] = DataType (cos (i * mul)); |
| } |
| } |
| |
| |
| |
| template <int LL2> |
| void FFTRealFixLen <LL2>::build_trigo_osc () |
| { |
| for (int i = 0; i < NBR_TRIGO_OSC; ++i) |
| { |
| OscType & osc = _trigo_osc [i]; |
| |
| const long len = static_cast <long> (TRIGO_TABLE_ARR_SIZE) << (i + 1); |
| const double mul = (0.5 * PI) / len; |
| osc.set_step (mul); |
| } |
| } |
| |
| |
| |
| #endif // FFTRealFixLen_CODEHEADER_INCLUDED |
| |
| #undef FFTRealFixLen_CURRENT_CODEHEADER |
| |
| |
| |
| /*\\\ EOF \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/ |