Mercurial > hg > constant-q-cpp

comparison test/TestFFT.cpp @ 131:6b13f9c694a8

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Start introducing the tests
author Chris Cannam <c.cannam@qmul.ac.uk>
date Mon, 19 May 2014 10:21:51 +0100
parents
children
comparison
equal deleted inserted replaced
130:1e33f719dde1 131:6b13f9c694a8
1 /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*- vi:set ts=8 sts=4 sw=4: */
2
3 #include "dsp/FFT.h"
4
5 #define BOOST_TEST_DYN_LINK
6 #define BOOST_TEST_MAIN
7
8 #include <boost/test/unit_test.hpp>
9
10 #include <stdexcept>
11
12 BOOST_AUTO_TEST_SUITE(TestFFT)
13
14 #define COMPARE_CONST(a, n) \
15 for (int cmp_i = 0; cmp_i < (int)(sizeof(a)/sizeof(a[0])); ++cmp_i) { \
16 BOOST_CHECK_SMALL(a[cmp_i] - n, 1e-14); \
17 }
18
19 #define COMPARE_ARRAY(a, b) \
20 for (int cmp_i = 0; cmp_i < (int)(sizeof(a)/sizeof(a[0])); ++cmp_i) { \
21 BOOST_CHECK_SMALL(a[cmp_i] - b[cmp_i], 1e-14); \
22 }
23
24 //!!! need at least one test with complex time-domain signal
25
26 BOOST_AUTO_TEST_CASE(forwardArrayBounds)
27 {
28 // initialise bins to something recognisable, so we can tell if
29 // they haven't been written; and allocate the inputs on the heap
30 // so that, if running under valgrind, we get warnings about
31 // overruns
32 double *in = new double[4];
33 in[0] = 1;
34 in[1] = 1;
35 in[2] = -1;
36 in[3] = -1;
37 double re[] = { 999, 999, 999, 999, 999, 999 };
38 double im[] = { 999, 999, 999, 999, 999, 999 };
39 FFT(4).process(false, in, 0, re+1, im+1);
40 // And check we haven't overrun the arrays
41 BOOST_CHECK_EQUAL(re[0], 999.0);
42 BOOST_CHECK_EQUAL(im[0], 999.0);
43 BOOST_CHECK_EQUAL(re[5], 999.0);
44 BOOST_CHECK_EQUAL(im[5], 999.0);
45 delete[] in;
46 }
47
48 BOOST_AUTO_TEST_CASE(r_forwardArrayBounds)
49 {
50 // initialise bins to something recognisable, so we can tell if
51 // they haven't been written; and allocate the inputs on the heap
52 // so that, if running under valgrind, we get warnings about
53 // overruns
54 double *in = new double[4];
55 in[0] = 1;
56 in[1] = 1;
57 in[2] = -1;
58 in[3] = -1;
59 double re[] = { 999, 999, 999, 999, 999, 999 };
60 double im[] = { 999, 999, 999, 999, 999, 999 };
61 FFTReal(4).forward(in, re+1, im+1);
62 // And check we haven't overrun the arrays
63 BOOST_CHECK_EQUAL(re[0], 999.0);
64 BOOST_CHECK_EQUAL(im[0], 999.0);
65 BOOST_CHECK_EQUAL(re[5], 999.0);
66 BOOST_CHECK_EQUAL(im[5], 999.0);
67 delete[] in;
68 }
69
70 BOOST_AUTO_TEST_CASE(inverseArrayBounds)
71 {
72 // initialise bins to something recognisable, so we can tell if
73 // they haven't been written; and allocate the inputs on the heap
74 // so that, if running under valgrind, we get warnings about
75 // overruns
76 double *re = new double[4];
77 double *im = new double[4];
78 re[0] = 0;
79 re[1] = 1;
80 re[2] = 0;
81 re[3] = 1;
82 im[0] = 0;
83 im[1] = -2;
84 im[2] = 0;
85 im[3] = 2;
86 double outre[] = { 999, 999, 999, 999, 999, 999 };
87 double outim[] = { 999, 999, 999, 999, 999, 999 };
88 FFT(4).process(true, re, im, outre+1, outim+1);
89 // And check we haven't overrun the arrays
90 BOOST_CHECK_EQUAL(outre[0], 999.0);
91 BOOST_CHECK_EQUAL(outim[0], 999.0);
92 BOOST_CHECK_EQUAL(outre[5], 999.0);
93 BOOST_CHECK_EQUAL(outim[5], 999.0);
94 delete[] re;
95 delete[] im;
96 }
97
98 BOOST_AUTO_TEST_CASE(r_inverseArrayBounds)
99 {
100 // initialise bins to something recognisable, so we can tell if
101 // they haven't been written; and allocate the inputs on the heap
102 // so that, if running under valgrind, we get warnings about
103 // overruns
104 double *re = new double[3];
105 double *im = new double[3];
106 re[0] = 0;
107 re[1] = 1;
108 re[2] = 0;
109 im[0] = 0;
110 im[1] = -2;
111 im[2] = 0;
112 double outre[] = { 999, 999, 999, 999, 999, 999 };
113 FFTReal(4).inverse(re, im, outre+1);
114 // And check we haven't overrun the arrays
115 BOOST_CHECK_EQUAL(outre[0], 999.0);
116 BOOST_CHECK_EQUAL(outre[5], 999.0);
117 delete[] re;
118 delete[] im;
119 }
120
121 BOOST_AUTO_TEST_CASE(dc)
122 {
123 // DC-only signal. The DC bin is purely real
124 double in[] = { 1, 1, 1, 1 };
125 double re[] = { 999, 999, 999, 999 };
126 double im[] = { 999, 999, 999, 999 };
127 FFT(4).process(false, in, 0, re, im);
128 BOOST_CHECK_EQUAL(re[0], 4.0);
129 BOOST_CHECK_EQUAL(re[1], 0.0);
130 BOOST_CHECK_EQUAL(re[2], 0.0);
131 BOOST_CHECK_EQUAL(re[3], 0.0);
132 COMPARE_CONST(im, 0.0);
133 double back[4];
134 double backim[4];
135 FFT(4).process(true, re, im, back, backim);
136 COMPARE_ARRAY(back, in);
137 COMPARE_CONST(backim, 0.0);
138 }
139
140 BOOST_AUTO_TEST_CASE(r_dc)
141 {
142 // DC-only signal. The DC bin is purely real
143 double in[] = { 1, 1, 1, 1 };
144 double re[] = { 999, 999, 999, 999 };
145 double im[] = { 999, 999, 999, 999 };
146 FFTReal(4).forward(in, re, im);
147 BOOST_CHECK_EQUAL(re[0], 4.0);
148 BOOST_CHECK_EQUAL(re[1], 0.0);
149 BOOST_CHECK_EQUAL(re[2], 0.0);
150 BOOST_CHECK_EQUAL(re[3], 0.0);
151 COMPARE_CONST(im, 0.0);
152 double back[4];
153 // check conjugates are reconstructed
154 re[3] = 999;
155 im[3] = 999;
156 FFTReal(4).inverse(re, im, back);
157 COMPARE_ARRAY(back, in);
158 }
159
160 BOOST_AUTO_TEST_CASE(c_dc)
161 {
162 // DC-only signal. The DC bin is purely real
163 double rin[] = { 1, 1, 1, 1 };
164 double iin[] = { 1, 1, 1, 1 };
165 double re[] = { 999, 999, 999, 999 };
166 double im[] = { 999, 999, 999, 999 };
167 FFT(4).process(false, rin, iin, re, im);
168 BOOST_CHECK_EQUAL(re[0], 4.0);
169 BOOST_CHECK_EQUAL(re[1], 0.0);
170 BOOST_CHECK_EQUAL(re[2], 0.0);
171 BOOST_CHECK_EQUAL(re[3], 0.0);
172 BOOST_CHECK_EQUAL(im[0], 4.0);
173 BOOST_CHECK_EQUAL(im[1], 0.0);
174 BOOST_CHECK_EQUAL(im[2], 0.0);
175 BOOST_CHECK_EQUAL(im[3], 0.0);
176 double back[4];
177 double backim[4];
178 FFT(4).process(true, re, im, back, backim);
179 COMPARE_ARRAY(back, rin);
180 COMPARE_ARRAY(backim, iin);
181 }
182
183 BOOST_AUTO_TEST_CASE(sine)
184 {
185 // Sine. Output is purely imaginary
186 double in[] = { 0, 1, 0, -1 };
187 double re[] = { 999, 999, 999, 999 };
188 double im[] = { 999, 999, 999, 999 };
189 FFT(4).process(false, in, 0, re, im);
190 COMPARE_CONST(re, 0.0);
191 BOOST_CHECK_EQUAL(im[0], 0.0);
192 BOOST_CHECK_EQUAL(im[1], -2.0);
193 BOOST_CHECK_EQUAL(im[2], 0.0);
194 BOOST_CHECK_EQUAL(im[3], 2.0);
195 double back[4];
196 double backim[4];
197 FFT(4).process(true, re, im, back, backim);
198 COMPARE_ARRAY(back, in);
199 COMPARE_CONST(backim, 0.0);
200 }
201
202 BOOST_AUTO_TEST_CASE(r_sine)
203 {
204 // Sine. Output is purely imaginary
205 double in[] = { 0, 1, 0, -1 };
206 double re[] = { 999, 999, 999, 999 };
207 double im[] = { 999, 999, 999, 999 };
208 FFTReal(4).forward(in, re, im);
209 COMPARE_CONST(re, 0.0);
210 BOOST_CHECK_EQUAL(im[0], 0.0);
211 BOOST_CHECK_EQUAL(im[1], -2.0);
212 BOOST_CHECK_EQUAL(im[2], 0.0);
213 BOOST_CHECK_EQUAL(im[3], 2.0);
214 double back[4];
215 // check conjugates are reconstructed
216 re[3] = 999;
217 im[3] = 999;
218 FFTReal(4).inverse(re, im, back);
219 COMPARE_ARRAY(back, in);
220 }
221
222 BOOST_AUTO_TEST_CASE(cosine)
223 {
224 // Cosine. Output is purely real
225 double in[] = { 1, 0, -1, 0 };
226 double re[] = { 999, 999, 999, 999 };
227 double im[] = { 999, 999, 999, 999 };
228 FFT(4).process(false, in, 0, re, im);
229 BOOST_CHECK_EQUAL(re[0], 0.0);
230 BOOST_CHECK_EQUAL(re[1], 2.0);
231 BOOST_CHECK_EQUAL(re[2], 0.0);
232 BOOST_CHECK_EQUAL(re[3], 2.0);
233 COMPARE_CONST(im, 0.0);
234 double back[4];
235 double backim[4];
236 FFT(4).process(true, re, im, back, backim);
237 COMPARE_ARRAY(back, in);
238 COMPARE_CONST(backim, 0.0);
239 }
240
241 BOOST_AUTO_TEST_CASE(r_cosine)
242 {
243 // Cosine. Output is purely real
244 double in[] = { 1, 0, -1, 0 };
245 double re[] = { 999, 999, 999, 999 };
246 double im[] = { 999, 999, 999, 999 };
247 FFTReal(4).forward(in, re, im);
248 BOOST_CHECK_EQUAL(re[0], 0.0);
249 BOOST_CHECK_EQUAL(re[1], 2.0);
250 BOOST_CHECK_EQUAL(re[2], 0.0);
251 BOOST_CHECK_EQUAL(re[3], 2.0);
252 COMPARE_CONST(im, 0.0);
253 double back[4];
254 // check conjugates are reconstructed
255 re[3] = 999;
256 im[3] = 999;
257 FFTReal(4).inverse(re, im, back);
258 COMPARE_ARRAY(back, in);
259 }
260
261 BOOST_AUTO_TEST_CASE(c_cosine)
262 {
263 // Cosine. Output is purely real
264 double rin[] = { 1, 0, -1, 0 };
265 double iin[] = { 1, 0, -1, 0 };
266 double re[] = { 999, 999, 999, 999 };
267 double im[] = { 999, 999, 999, 999 };
268 FFT(4).process(false, rin, iin, re, im);
269 BOOST_CHECK_EQUAL(re[0], 0.0);
270 BOOST_CHECK_EQUAL(re[1], 2.0);
271 BOOST_CHECK_EQUAL(re[2], 0.0);
272 BOOST_CHECK_EQUAL(re[3], 2.0);
273 BOOST_CHECK_EQUAL(im[0], 0.0);
274 BOOST_CHECK_EQUAL(im[1], 2.0);
275 BOOST_CHECK_EQUAL(im[2], 0.0);
276 BOOST_CHECK_EQUAL(im[3], 2.0);
277 double back[4];
278 double backim[4];
279 FFT(4).process(true, re, im, back, backim);
280 COMPARE_ARRAY(back, rin);
281 COMPARE_ARRAY(backim, iin);
282 }
283
284 BOOST_AUTO_TEST_CASE(sineCosine)
285 {
286 // Sine and cosine mixed
287 double in[] = { 0.5, 1, -0.5, -1 };
288 double re[] = { 999, 999, 999, 999 };
289 double im[] = { 999, 999, 999, 999 };
290 FFT(4).process(false, in, 0, re, im);
291 BOOST_CHECK_EQUAL(re[0], 0.0);
292 BOOST_CHECK_CLOSE(re[1], 1.0, 1e-12);
293 BOOST_CHECK_EQUAL(re[2], 0.0);
294 BOOST_CHECK_CLOSE(re[3], 1.0, 1e-12);
295 BOOST_CHECK_EQUAL(im[0], 0.0);
296 BOOST_CHECK_CLOSE(im[1], -2.0, 1e-12);
297 BOOST_CHECK_EQUAL(im[2], 0.0);
298 BOOST_CHECK_CLOSE(im[3], 2.0, 1e-12);
299 double back[4];
300 double backim[4];
301 FFT(4).process(true, re, im, back, backim);
302 COMPARE_ARRAY(back, in);
303 COMPARE_CONST(backim, 0.0);
304 }
305
306 BOOST_AUTO_TEST_CASE(r_sineCosine)
307 {
308 // Sine and cosine mixed
309 double in[] = { 0.5, 1, -0.5, -1 };
310 double re[] = { 999, 999, 999, 999 };
311 double im[] = { 999, 999, 999, 999 };
312 FFTReal(4).forward(in, re, im);
313 BOOST_CHECK_EQUAL(re[0], 0.0);
314 BOOST_CHECK_CLOSE(re[1], 1.0, 1e-12);
315 BOOST_CHECK_EQUAL(re[2], 0.0);
316 BOOST_CHECK_CLOSE(re[3], 1.0, 1e-12);
317 BOOST_CHECK_EQUAL(im[0], 0.0);
318 BOOST_CHECK_CLOSE(im[1], -2.0, 1e-12);
319 BOOST_CHECK_EQUAL(im[2], 0.0);
320 BOOST_CHECK_CLOSE(im[3], 2.0, 1e-12);
321 double back[4];
322 // check conjugates are reconstructed
323 re[3] = 999;
324 im[3] = 999;
325 FFTReal(4).inverse(re, im, back);
326 COMPARE_ARRAY(back, in);
327 }
328
329 BOOST_AUTO_TEST_CASE(c_sineCosine)
330 {
331 double rin[] = { 1, 0, -1, 0 };
332 double iin[] = { 0, 1, 0, -1 };
333 double re[] = { 999, 999, 999, 999 };
334 double im[] = { 999, 999, 999, 999 };
335 FFT(4).process(false, rin, iin, re, im);
336 BOOST_CHECK_EQUAL(re[0], 0.0);
337 BOOST_CHECK_EQUAL(re[1], 4.0);
338 BOOST_CHECK_EQUAL(re[2], 0.0);
339 BOOST_CHECK_EQUAL(re[3], 0.0);
340 COMPARE_CONST(im, 0.0);
341 double back[4];
342 double backim[4];
343 FFT(4).process(true, re, im, back, backim);
344 COMPARE_ARRAY(back, rin);
345 COMPARE_ARRAY(backim, iin);
346 }
347
348 BOOST_AUTO_TEST_CASE(nyquist)
349 {
350 double in[] = { 1, -1, 1, -1 };
351 double re[] = { 999, 999, 999, 999 };
352 double im[] = { 999, 999, 999, 999 };
353 FFT(4).process(false, in, 0, re, im);
354 BOOST_CHECK_EQUAL(re[0], 0.0);
355 BOOST_CHECK_EQUAL(re[1], 0.0);
356 BOOST_CHECK_EQUAL(re[2], 4.0);
357 BOOST_CHECK_EQUAL(re[3], 0.0);
358 COMPARE_CONST(im, 0.0);
359 double back[4];
360 double backim[4];
361 FFT(4).process(true, re, im, back, backim);
362 COMPARE_ARRAY(back, in);
363 COMPARE_CONST(backim, 0.0);
364 }
365
366 BOOST_AUTO_TEST_CASE(r_nyquist)
367 {
368 double in[] = { 1, -1, 1, -1 };
369 double re[] = { 999, 999, 999, 999 };
370 double im[] = { 999, 999, 999, 999 };
371 FFTReal(4).forward(in, re, im);
372 BOOST_CHECK_EQUAL(re[0], 0.0);
373 BOOST_CHECK_EQUAL(re[1], 0.0);
374 BOOST_CHECK_EQUAL(re[2], 4.0);
375 BOOST_CHECK_EQUAL(re[3], 0.0);
376 COMPARE_CONST(im, 0.0);
377 double back[4];
378 // check conjugates are reconstructed
379 re[3] = 999;
380 im[3] = 999;
381 FFTReal(4).inverse(re, im, back);
382 COMPARE_ARRAY(back, in);
383 }
384
385 BOOST_AUTO_TEST_CASE(dirac)
386 {
387 double in[] = { 1, 0, 0, 0 };
388 double re[] = { 999, 999, 999, 999 };
389 double im[] = { 999, 999, 999, 999 };
390 FFT(4).process(false, in, 0, re, im);
391 BOOST_CHECK_EQUAL(re[0], 1.0);
392 BOOST_CHECK_EQUAL(re[1], 1.0);
393 BOOST_CHECK_EQUAL(re[2], 1.0);
394 BOOST_CHECK_EQUAL(re[3], 1.0);
395 COMPARE_CONST(im, 0.0);
396 double back[4];
397 double backim[4];
398 FFT(4).process(true, re, im, back, backim);
399 COMPARE_ARRAY(back, in);
400 COMPARE_CONST(backim, 0.0);
401 }
402
403 BOOST_AUTO_TEST_CASE(r_dirac)
404 {
405 double in[] = { 1, 0, 0, 0 };
406 double re[] = { 999, 999, 999, 999 };
407 double im[] = { 999, 999, 999, 999 };
408 FFTReal(4).forward(in, re, im);
409 BOOST_CHECK_EQUAL(re[0], 1.0);
410 BOOST_CHECK_EQUAL(re[1], 1.0);
411 BOOST_CHECK_EQUAL(re[2], 1.0);
412 BOOST_CHECK_EQUAL(re[3], 1.0);
413 COMPARE_CONST(im, 0.0);
414 double back[4];
415 // check conjugates are reconstructed
416 re[3] = 999;
417 im[3] = 999;
418 FFTReal(4).inverse(re, im, back);
419 COMPARE_ARRAY(back, in);
420 }
421
422 BOOST_AUTO_TEST_CASE(sizes)
423 {
424 // Complex supports any size. A single test with an odd size
425 // will do here, without getting too much into our expectations
426 // about supported butterflies etc
427
428 double in[] = { 1, 1, 1 };
429 double re[] = { 999, 999, 999 };
430 double im[] = { 999, 999, 999 };
431 FFT(3).process(false, in, 0, re, im);
432 BOOST_CHECK_EQUAL(re[0], 3.0);
433 BOOST_CHECK_EQUAL(re[1], 0.0);
434 BOOST_CHECK_EQUAL(re[2], 0.0);
435 COMPARE_CONST(im, 0.0);
436 double back[3];
437 double backim[3];
438 FFT(3).process(true, re, im, back, backim);
439 COMPARE_ARRAY(back, in);
440 COMPARE_CONST(backim, 0.0);
441 }
442
443 BOOST_AUTO_TEST_CASE(r_sizes)
444 {
445 // Real supports any even size, but not odd ones
446
447 BOOST_CHECK_THROW(FFTReal r(3), std::invalid_argument);
448
449 double in[] = { 1, 1, 1, 1, 1, 1 };
450 double re[] = { 999, 999, 999, 999, 999, 999 };
451 double im[] = { 999, 999, 999, 999, 999, 999 };
452 FFTReal(6).forward(in, re, im);
453 BOOST_CHECK_EQUAL(re[0], 6.0);
454 BOOST_CHECK_EQUAL(re[1], 0.0);
455 BOOST_CHECK_EQUAL(re[2], 0.0);
456 BOOST_CHECK_EQUAL(re[3], 0.0);
457 BOOST_CHECK_EQUAL(re[4], 0.0);
458 BOOST_CHECK_EQUAL(re[5], 0.0);
459 COMPARE_CONST(im, 0.0);
460 double back[6];
461 FFTReal(6).inverse(re, im, back);
462 COMPARE_ARRAY(back, in);
463 }
464
465 BOOST_AUTO_TEST_SUITE_END()
466