1/*****************************************************************************
2
3 FFTRealFixLen.hpp
4 Copyright (c) 2005 Laurent de Soras
5
6--- Legal stuff ---
7
8This library is free software; you can redistribute it and/or
9modify it under the terms of the GNU Lesser General Public
10License as published by the Free Software Foundation; either
11version 2.1 of the License, or (at your option) any later version.
12
13This library is distributed in the hope that it will be useful,
14but WITHOUT ANY WARRANTY; without even the implied warranty of
15MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16Lesser General Public License for more details.
17
18You should have received a copy of the GNU Lesser General Public
19License along with this library; if not, write to the Free Software
20Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
21
22*Tab=3***********************************************************************/
23
24
25
26#if defined (FFTRealFixLen_CURRENT_CODEHEADER)
27 #error Recursive inclusion of FFTRealFixLen code header.
28#endif
29#define FFTRealFixLen_CURRENT_CODEHEADER
30
31#if ! defined (FFTRealFixLen_CODEHEADER_INCLUDED)
32#define FFTRealFixLen_CODEHEADER_INCLUDED
33
34
35
36/*\\\ INCLUDE FILES \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/
37
38#include "def.h"
39#include "FFTRealPassDirect.h"
40#include "FFTRealPassInverse.h"
41#include "FFTRealSelect.h"
42
43#include <cassert>
44#include <cmath>
45
46namespace std { }
47
48
49
50/*\\\ PUBLIC \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/
51
52
53
54template <int LL2>
55FFTRealFixLen <LL2>::FFTRealFixLen ()
56: _buffer (FFT_LEN)
57, _br_data (BR_ARR_SIZE)
58, _trigo_data (TRIGO_TABLE_ARR_SIZE)
59, _trigo_osc ()
60{
61 build_br_lut ();
62 build_trigo_lut ();
63 build_trigo_osc ();
64}
65
66
67
68template <int LL2>
69long FFTRealFixLen <LL2>::get_length () const
70{
71 return (FFT_LEN);
72}
73
74
75
76// General case
77template <int LL2>
78void FFTRealFixLen <LL2>::do_fft (DataType f [], const DataType x [])
79{
80 assert (f != 0);
81 assert (x != 0);
82 assert (x != f);
83 assert (FFT_LEN_L2 >= 3);
84
85 // Do the transform in several passes
86 const DataType * cos_ptr = &_trigo_data [0];
87 const long * br_ptr = &_br_data [0];
88
89 FFTRealPassDirect <FFT_LEN_L2 - 1>::process (
90 FFT_LEN,
91 f,
92 &_buffer [0],
93 x,
94 cos_ptr,
95 TRIGO_TABLE_ARR_SIZE,
96 br_ptr,
97 &_trigo_osc [0]
98 );
99}
100
101// 4-point FFT
102template <>
103void FFTRealFixLen <2>::do_fft (DataType f [], const DataType x [])
104{
105 assert (f != 0);
106 assert (x != 0);
107 assert (x != f);
108
109 f [1] = x [0] - x [2];
110 f [3] = x [1] - x [3];
111
112 const DataType b_0 = x [0] + x [2];
113 const DataType b_2 = x [1] + x [3];
114
115 f [0] = b_0 + b_2;
116 f [2] = b_0 - b_2;
117}
118
119// 2-point FFT
120template <>
121void FFTRealFixLen <1>::do_fft (DataType f [], const DataType x [])
122{
123 assert (f != 0);
124 assert (x != 0);
125 assert (x != f);
126
127 f [0] = x [0] + x [1];
128 f [1] = x [0] - x [1];
129}
130
131// 1-point FFT
132template <>
133void FFTRealFixLen <0>::do_fft (DataType f [], const DataType x [])
134{
135 assert (f != 0);
136 assert (x != 0);
137
138 f [0] = x [0];
139}
140
141
142
143// General case
144template <int LL2>
145void FFTRealFixLen <LL2>::do_ifft (const DataType f [], DataType x [])
146{
147 assert (f != 0);
148 assert (x != 0);
149 assert (x != f);
150 assert (FFT_LEN_L2 >= 3);
151
152 // Do the transform in several passes
153 DataType * s_ptr =
154 FFTRealSelect <FFT_LEN_L2 & 1>::sel_bin (&_buffer [0], x);
155 DataType * d_ptr =
156 FFTRealSelect <FFT_LEN_L2 & 1>::sel_bin (x, &_buffer [0]);
157 const DataType * cos_ptr = &_trigo_data [0];
158 const long * br_ptr = &_br_data [0];
159
160 FFTRealPassInverse <FFT_LEN_L2 - 1>::process (
161 FFT_LEN,
162 d_ptr,
163 s_ptr,
164 f,
165 cos_ptr,
166 TRIGO_TABLE_ARR_SIZE,
167 br_ptr,
168 &_trigo_osc [0]
169 );
170}
171
172// 4-point IFFT
173template <>
174void FFTRealFixLen <2>::do_ifft (const DataType f [], DataType x [])
175{
176 assert (f != 0);
177 assert (x != 0);
178 assert (x != f);
179
180 const DataType b_0 = f [0] + f [2];
181 const DataType b_2 = f [0] - f [2];
182
183 x [0] = b_0 + f [1] * 2;
184 x [2] = b_0 - f [1] * 2;
185 x [1] = b_2 + f [3] * 2;
186 x [3] = b_2 - f [3] * 2;
187}
188
189// 2-point IFFT
190template <>
191void FFTRealFixLen <1>::do_ifft (const DataType f [], DataType x [])
192{
193 assert (f != 0);
194 assert (x != 0);
195 assert (x != f);
196
197 x [0] = f [0] + f [1];
198 x [1] = f [0] - f [1];
199}
200
201// 1-point IFFT
202template <>
203void FFTRealFixLen <0>::do_ifft (const DataType f [], DataType x [])
204{
205 assert (f != 0);
206 assert (x != 0);
207 assert (x != f);
208
209 x [0] = f [0];
210}
211
212
213
214
215template <int LL2>
216void FFTRealFixLen <LL2>::rescale (DataType x []) const
217{
218 assert (x != 0);
219
220 const DataType mul = DataType (1.0 / FFT_LEN);
221
222 if (FFT_LEN < 4)
223 {
224 long i = FFT_LEN - 1;
225 do
226 {
227 x [i] *= mul;
228 --i;
229 }
230 while (i >= 0);
231 }
232
233 else
234 {
235 assert ((FFT_LEN & 3) == 0);
236
237 // Could be optimized with SIMD instruction sets (needs alignment check)
238 long i = FFT_LEN - 4;
239 do
240 {
241 x [i + 0] *= mul;
242 x [i + 1] *= mul;
243 x [i + 2] *= mul;
244 x [i + 3] *= mul;
245 i -= 4;
246 }
247 while (i >= 0);
248 }
249}
250
251
252
253/*\\\ PROTECTED \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/
254
255
256
257/*\\\ PRIVATE \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/
258
259
260
261template <int LL2>
262void FFTRealFixLen <LL2>::build_br_lut ()
263{
264 _br_data [0] = 0;
265 for (long cnt = 1; cnt < BR_ARR_SIZE; ++cnt)
266 {
267 long index = cnt << 2;
268 long br_index = 0;
269
270 int bit_cnt = FFT_LEN_L2;
271 do
272 {
273 br_index <<= 1;
274 br_index += (index & 1);
275 index >>= 1;
276
277 -- bit_cnt;
278 }
279 while (bit_cnt > 0);
280
281 _br_data [cnt] = br_index;
282 }
283}
284
285
286
287template <int LL2>
288void FFTRealFixLen <LL2>::build_trigo_lut ()
289{
290 const double mul = (0.5 * PI) / TRIGO_TABLE_ARR_SIZE;
291 for (long i = 0; i < TRIGO_TABLE_ARR_SIZE; ++ i)
292 {
293 using namespace std;
294
295 _trigo_data [i] = DataType (cos (x: i * mul));
296 }
297}
298
299
300
301template <int LL2>
302void FFTRealFixLen <LL2>::build_trigo_osc ()
303{
304 for (int i = 0; i < NBR_TRIGO_OSC; ++i)
305 {
306 OscType & osc = _trigo_osc [i];
307
308 const long len = static_cast <long> (TRIGO_TABLE_ARR_SIZE) << (i + 1);
309 const double mul = (0.5 * PI) / len;
310 osc.set_step (mul);
311 }
312}
313
314
315
316#endif // FFTRealFixLen_CODEHEADER_INCLUDED
317
318#undef FFTRealFixLen_CURRENT_CODEHEADER
319
320
321
322/*\\\ EOF \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\*/
323

source code of qtmultimedia/examples/multimedia/spectrum/3rdparty/fftreal/FFTRealFixLen.hpp