1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22 package org.hipparchus.transform;
23
24 import java.util.Random;
25
26 import org.hipparchus.analysis.UnivariateFunction;
27 import org.hipparchus.analysis.function.Sin;
28 import org.hipparchus.analysis.function.Sinc;
29 import org.hipparchus.complex.Complex;
30 import org.hipparchus.exception.MathIllegalArgumentException;
31 import org.hipparchus.util.FastMath;
32 import org.junit.Assert;
33 import org.junit.Test;
34
35
36
37
38
39
40
41
42 public final class FastFourierTransformerTest {
43
44 private final static long SEED = 20110111L;
45
46
47
48
49
50 @Test
51 public void testTransformComplexSizeNotAPowerOfTwo() {
52 final int n = 127;
53 final Complex[] x = createComplexData(n);
54 final DftNormalization[] norm;
55 norm = DftNormalization.values();
56 final TransformType[] type;
57 type = TransformType.values();
58 for (int i = 0; i < norm.length; i++) {
59 for (int j = 0; j < type.length; j++) {
60 final FastFourierTransformer fft;
61 fft = new FastFourierTransformer(norm[i]);
62 try {
63 fft.transform(x, type[j]);
64 Assert.fail(norm[i] + ", " + type[j] +
65 ": MathIllegalArgumentException was expected");
66 } catch (MathIllegalArgumentException e) {
67
68 }
69 }
70 }
71 }
72
73 @Test
74 public void testTransformRealSizeNotAPowerOfTwo() {
75 final int n = 127;
76 final double[] x = createRealData(n);
77 final DftNormalization[] norm;
78 norm = DftNormalization.values();
79 final TransformType[] type;
80 type = TransformType.values();
81 for (int i = 0; i < norm.length; i++) {
82 for (int j = 0; j < type.length; j++) {
83 final FastFourierTransformer fft;
84 fft = new FastFourierTransformer(norm[i]);
85 try {
86 fft.transform(x, type[j]);
87 Assert.fail(norm[i] + ", " + type[j] +
88 ": MathIllegalArgumentException was expected");
89 } catch (MathIllegalArgumentException e) {
90
91 }
92 }
93 }
94 }
95
96 @Test
97 public void testTransformFunctionSizeNotAPowerOfTwo() {
98 final int n = 127;
99 final UnivariateFunction f = new Sin();
100 final DftNormalization[] norm;
101 norm = DftNormalization.values();
102 final TransformType[] type;
103 type = TransformType.values();
104 for (int i = 0; i < norm.length; i++) {
105 for (int j = 0; j < type.length; j++) {
106 final FastFourierTransformer fft;
107 fft = new FastFourierTransformer(norm[i]);
108 try {
109 fft.transform(f, 0.0, Math.PI, n, type[j]);
110 Assert.fail(norm[i] + ", " + type[j] +
111 ": MathIllegalArgumentException was expected");
112 } catch (MathIllegalArgumentException e) {
113
114 }
115 }
116 }
117 }
118
119 @Test
120 public void testTransformFunctionNotStrictlyPositiveNumberOfSamples() {
121 final int n = -128;
122 final UnivariateFunction f = new Sin();
123 final DftNormalization[] norm;
124 norm = DftNormalization.values();
125 final TransformType[] type;
126 type = TransformType.values();
127 for (int i = 0; i < norm.length; i++) {
128 for (int j = 0; j < type.length; j++) {
129 final FastFourierTransformer fft;
130 fft = new FastFourierTransformer(norm[i]);
131 try {
132 fft.transform(f, 0.0, Math.PI, n, type[j]);
133 fft.transform(f, 0.0, Math.PI, n, type[j]);
134 Assert.fail(norm[i] + ", " + type[j] +
135 ": MathIllegalArgumentException was expected");
136 } catch (MathIllegalArgumentException e) {
137
138 }
139 }
140 }
141 }
142
143 @Test
144 public void testTransformFunctionInvalidBounds() {
145 final int n = 128;
146 final UnivariateFunction f = new Sin();
147 final DftNormalization[] norm;
148 norm = DftNormalization.values();
149 final TransformType[] type;
150 type = TransformType.values();
151 for (int i = 0; i < norm.length; i++) {
152 for (int j = 0; j < type.length; j++) {
153 final FastFourierTransformer fft;
154 fft = new FastFourierTransformer(norm[i]);
155 try {
156 fft.transform(f, Math.PI, 0.0, n, type[j]);
157 Assert.fail(norm[i] + ", " + type[j] +
158 ": MathIllegalArgumentException was expected");
159 } catch (MathIllegalArgumentException e) {
160
161 }
162 }
163 }
164 }
165
166 @Test
167 public void testTransformOnePoint() {
168 double[][] data = new double[][] { { 1.0 }, { 2.0} };
169 FastFourierTransformer.transformInPlace(data,
170 DftNormalization.STANDARD,
171 TransformType.FORWARD);
172 Assert.assertEquals(1.0, data[0][0], 1.0e-15);
173 Assert.assertEquals(2.0, data[1][0], 1.0e-15);
174 }
175
176
177
178
179
180 private static Complex[] createComplexData(final int n) {
181 final Random random = new Random(SEED);
182 final Complex[] data = new Complex[n];
183 for (int i = 0; i < n; i++) {
184 final double re = 2.0 * random.nextDouble() - 1.0;
185 final double im = 2.0 * random.nextDouble() - 1.0;
186 data[i] = new Complex(re, im);
187 }
188 return data;
189 }
190
191 private static double[] createRealData(final int n) {
192 final Random random = new Random(SEED);
193 final double[] data = new double[n];
194 for (int i = 0; i < n; i++) {
195 data[i] = 2.0 * random.nextDouble() - 1.0;
196 }
197 return data;
198 }
199
200
201 private static Complex[] dft(final Complex[] x, final int sgn) {
202 final int n = x.length;
203 final double[] cos = new double[n];
204 final double[] sin = new double[n];
205 final Complex[] y = new Complex[n];
206 for (int i = 0; i < n; i++) {
207 final double arg = 2.0 * FastMath.PI * i / n;
208 cos[i] = FastMath.cos(arg);
209 sin[i] = FastMath.sin(arg);
210 }
211 for (int i = 0; i < n; i++) {
212 double yr = 0.0;
213 double yi = 0.0;
214 for (int j = 0; j < n; j++) {
215 final int index = (i * j) % n;
216 final double c = cos[index];
217 final double s = sin[index];
218 final double xr = x[j].getReal();
219 final double xi = x[j].getImaginary();
220 yr += c * xr - sgn * s * xi;
221 yi += sgn * s * xr + c * xi;
222 }
223 y[i] = new Complex(yr, yi);
224 }
225 return y;
226 }
227
228 private static void doTestTransformComplex(final int n, final double tol,
229 final DftNormalization normalization,
230 final TransformType type) {
231 final FastFourierTransformer fft;
232 fft = new FastFourierTransformer(normalization);
233 final Complex[] x = createComplexData(n);
234 final Complex[] expected;
235 final double s;
236 if (type==TransformType.FORWARD) {
237 expected = dft(x, -1);
238 if (normalization == DftNormalization.STANDARD){
239 s = 1.0;
240 } else {
241 s = 1.0 / FastMath.sqrt(n);
242 }
243 } else {
244 expected = dft(x, 1);
245 if (normalization == DftNormalization.STANDARD) {
246 s = 1.0 / n;
247 } else {
248 s = 1.0 / FastMath.sqrt(n);
249 }
250 }
251 final Complex[] actual = fft.transform(x, type);
252 for (int i = 0; i < n; i++) {
253 final String msg;
254 msg = String.format("%s, %s, %d, %d", normalization, type, n, i);
255 final double re = s * expected[i].getReal();
256 Assert.assertEquals(msg, re, actual[i].getReal(),
257 tol * FastMath.abs(re));
258 final double im = s * expected[i].getImaginary();
259 Assert.assertEquals(msg, im, actual[i].getImaginary(), tol *
260 FastMath.abs(re));
261 }
262 }
263
264 private static void doTestTransformReal(final int n, final double tol,
265 final DftNormalization normalization,
266 final TransformType type) {
267 final FastFourierTransformer fft;
268 fft = new FastFourierTransformer(normalization);
269 final double[] x = createRealData(n);
270 final Complex[] xc = new Complex[n];
271 for (int i = 0; i < n; i++) {
272 xc[i] = new Complex(x[i], 0.0);
273 }
274 final Complex[] expected;
275 final double s;
276 if (type == TransformType.FORWARD) {
277 expected = dft(xc, -1);
278 if (normalization == DftNormalization.STANDARD) {
279 s = 1.0;
280 } else {
281 s = 1.0 / FastMath.sqrt(n);
282 }
283 } else {
284 expected = dft(xc, 1);
285 if (normalization == DftNormalization.STANDARD) {
286 s = 1.0 / n;
287 } else {
288 s = 1.0 / FastMath.sqrt(n);
289 }
290 }
291 final Complex[] actual = fft.transform(x, type);
292 for (int i = 0; i < n; i++) {
293 final String msg;
294 msg = String.format("%s, %s, %d, %d", normalization, type, n, i);
295 final double re = s * expected[i].getReal();
296 Assert.assertEquals(msg, re, actual[i].getReal(),
297 tol * FastMath.abs(re));
298 final double im = s * expected[i].getImaginary();
299 Assert.assertEquals(msg, im, actual[i].getImaginary(), tol *
300 FastMath.abs(re));
301 }
302 }
303
304 private static void doTestTransformFunction(final UnivariateFunction f,
305 final double min, final double max, int n, final double tol,
306 final DftNormalization normalization,
307 final TransformType type) {
308 final FastFourierTransformer fft;
309 fft = new FastFourierTransformer(normalization);
310 final Complex[] x = new Complex[n];
311 for (int i = 0; i < n; i++) {
312 final double t = min + i * (max - min) / n;
313 x[i] = new Complex(f.value(t));
314 }
315 final Complex[] expected;
316 final double s;
317 if (type == TransformType.FORWARD) {
318 expected = dft(x, -1);
319 if (normalization == DftNormalization.STANDARD) {
320 s = 1.0;
321 } else {
322 s = 1.0 / FastMath.sqrt(n);
323 }
324 } else {
325 expected = dft(x, 1);
326 if (normalization == DftNormalization.STANDARD) {
327 s = 1.0 / n;
328 } else {
329 s = 1.0 / FastMath.sqrt(n);
330 }
331 }
332 final Complex[] actual = fft.transform(f, min, max, n, type);
333 for (int i = 0; i < n; i++) {
334 final String msg = String.format("%d, %d", n, i);
335 final double re = s * expected[i].getReal();
336 Assert.assertEquals(msg, re, actual[i].getReal(),
337 tol * FastMath.abs(re));
338 final double im = s * expected[i].getImaginary();
339 Assert.assertEquals(msg, im, actual[i].getImaginary(), tol *
340 FastMath.abs(re));
341 }
342 }
343
344
345
346
347
348 @Test
349 public void testTransformComplex() {
350 final DftNormalization[] norm;
351 norm = DftNormalization.values();
352 final TransformType[] type;
353 type = TransformType.values();
354 for (int i = 0; i < norm.length; i++) {
355 for (int j = 0; j < type.length; j++) {
356 doTestTransformComplex(2, 1.0E-15, norm[i], type[j]);
357 doTestTransformComplex(4, 1.0E-14, norm[i], type[j]);
358 doTestTransformComplex(8, 1.0E-14, norm[i], type[j]);
359 doTestTransformComplex(16, 1.0E-13, norm[i], type[j]);
360 doTestTransformComplex(32, 1.0E-13, norm[i], type[j]);
361 doTestTransformComplex(64, 1.0E-12, norm[i], type[j]);
362 doTestTransformComplex(128, 1.0E-12, norm[i], type[j]);
363 }
364 }
365 }
366
367 @Test
368 public void testStandardTransformReal() {
369 final DftNormalization[] norm;
370 norm = DftNormalization.values();
371 final TransformType[] type;
372 type = TransformType.values();
373 for (int i = 0; i < norm.length; i++) {
374 for (int j = 0; j < type.length; j++) {
375 doTestTransformReal(2, 1.0E-15, norm[i], type[j]);
376 doTestTransformReal(4, 1.0E-14, norm[i], type[j]);
377 doTestTransformReal(8, 1.0E-14, norm[i], type[j]);
378 doTestTransformReal(16, 1.0E-13, norm[i], type[j]);
379 doTestTransformReal(32, 1.0E-13, norm[i], type[j]);
380 doTestTransformReal(64, 1.0E-13, norm[i], type[j]);
381 doTestTransformReal(128, 1.0E-11, norm[i], type[j]);
382 }
383 }
384 }
385
386 @Test
387 public void testStandardTransformFunction() {
388 final UnivariateFunction f = new Sinc();
389 final double min = -FastMath.PI;
390 final double max = FastMath.PI;
391 final DftNormalization[] norm;
392 norm = DftNormalization.values();
393 final TransformType[] type;
394 type = TransformType.values();
395 for (int i = 0; i < norm.length; i++) {
396 for (int j = 0; j < type.length; j++) {
397 doTestTransformFunction(f, min, max, 2, 1.0E-15, norm[i], type[j]);
398 doTestTransformFunction(f, min, max, 4, 1.0E-14, norm[i], type[j]);
399 doTestTransformFunction(f, min, max, 8, 1.0E-14, norm[i], type[j]);
400 doTestTransformFunction(f, min, max, 16, 1.0E-13, norm[i], type[j]);
401 doTestTransformFunction(f, min, max, 32, 1.0E-13, norm[i], type[j]);
402 doTestTransformFunction(f, min, max, 64, 1.0E-12, norm[i], type[j]);
403 doTestTransformFunction(f, min, max, 128, 1.0E-11, norm[i], type[j]);
404 }
405 }
406 }
407
408
409
410
411
412
413
414
415 @Test
416 public void testAdHocData() {
417 FastFourierTransformer transformer;
418 transformer = new FastFourierTransformer(DftNormalization.STANDARD);
419 Complex[] result; double tolerance = 1E-12;
420
421 double[] x = {1.3, 2.4, 1.7, 4.1, 2.9, 1.7, 5.1, 2.7};
422 Complex[] y = {
423 new Complex(21.9, 0.0),
424 new Complex(-2.09497474683058, 1.91507575950825),
425 new Complex(-2.6, 2.7),
426 new Complex(-1.10502525316942, -4.88492424049175),
427 new Complex(0.1, 0.0),
428 new Complex(-1.10502525316942, 4.88492424049175),
429 new Complex(-2.6, -2.7),
430 new Complex(-2.09497474683058, -1.91507575950825)};
431
432 result = transformer.transform(x, TransformType.FORWARD);
433 for (int i = 0; i < result.length; i++) {
434 Assert.assertEquals(y[i].getReal(), result[i].getReal(), tolerance);
435 Assert.assertEquals(y[i].getImaginary(), result[i].getImaginary(), tolerance);
436 }
437
438 result = transformer.transform(y, TransformType.INVERSE);
439 for (int i = 0; i < result.length; i++) {
440 Assert.assertEquals(x[i], result[i].getReal(), tolerance);
441 Assert.assertEquals(0.0, result[i].getImaginary(), tolerance);
442 }
443
444 double[] x2 = {10.4, 21.6, 40.8, 13.6, 23.2, 32.8, 13.6, 19.2};
445 TransformUtils.scaleArray(x2, 1.0 / FastMath.sqrt(x2.length));
446 Complex[] y2 = y;
447
448 transformer = new FastFourierTransformer(DftNormalization.UNITARY);
449 result = transformer.transform(y2, TransformType.FORWARD);
450 for (int i = 0; i < result.length; i++) {
451 Assert.assertEquals(x2[i], result[i].getReal(), tolerance);
452 Assert.assertEquals(0.0, result[i].getImaginary(), tolerance);
453 }
454
455 result = transformer.transform(x2, TransformType.INVERSE);
456 for (int i = 0; i < result.length; i++) {
457 Assert.assertEquals(y2[i].getReal(), result[i].getReal(), tolerance);
458 Assert.assertEquals(y2[i].getImaginary(), result[i].getImaginary(), tolerance);
459 }
460 }
461
462
463
464
465 @Test
466 public void testSinFunction() {
467 UnivariateFunction f = new Sin();
468 FastFourierTransformer transformer;
469 transformer = new FastFourierTransformer(DftNormalization.STANDARD);
470 Complex[] result; int N = 1 << 8;
471 double min, max, tolerance = 1E-12;
472
473 min = 0.0; max = 2.0 * FastMath.PI;
474 result = transformer.transform(f, min, max, N, TransformType.FORWARD);
475 Assert.assertEquals(0.0, result[1].getReal(), tolerance);
476 Assert.assertEquals(-(N >> 1), result[1].getImaginary(), tolerance);
477 Assert.assertEquals(0.0, result[N-1].getReal(), tolerance);
478 Assert.assertEquals(N >> 1, result[N-1].getImaginary(), tolerance);
479 for (int i = 0; i < N-1; i += (i == 0 ? 2 : 1)) {
480 Assert.assertEquals(0.0, result[i].getReal(), tolerance);
481 Assert.assertEquals(0.0, result[i].getImaginary(), tolerance);
482 }
483
484 min = -FastMath.PI; max = FastMath.PI;
485 result = transformer.transform(f, min, max, N, TransformType.INVERSE);
486 Assert.assertEquals(0.0, result[1].getReal(), tolerance);
487 Assert.assertEquals(-0.5, result[1].getImaginary(), tolerance);
488 Assert.assertEquals(0.0, result[N-1].getReal(), tolerance);
489 Assert.assertEquals(0.5, result[N-1].getImaginary(), tolerance);
490 for (int i = 0; i < N-1; i += (i == 0 ? 2 : 1)) {
491 Assert.assertEquals(0.0, result[i].getReal(), tolerance);
492 Assert.assertEquals(0.0, result[i].getImaginary(), tolerance);
493 }
494 }
495
496 }