1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2014 Jianwei Cui <thucjw@gmail.com>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #include "main.h"
11 #include <complex>
12 #include <cmath>
13 #include <Eigen/CXX11/Tensor>
14 
15 using Eigen::Tensor;
16 
17 template <int DataLayout>
test_1D_fft_ifft_invariant(int sequence_length)18 static void test_1D_fft_ifft_invariant(int sequence_length) {
19   Tensor<double, 1, DataLayout> tensor(sequence_length);
20   tensor.setRandom();
21 
22   array<int, 1> fft;
23   fft[0] = 0;
24 
25   Tensor<std::complex<double>, 1, DataLayout> tensor_after_fft;
26   Tensor<std::complex<double>, 1, DataLayout> tensor_after_fft_ifft;
27 
28   tensor_after_fft = tensor.template fft<Eigen::BothParts, Eigen::FFT_FORWARD>(fft);
29   tensor_after_fft_ifft = tensor_after_fft.template fft<Eigen::BothParts, Eigen::FFT_REVERSE>(fft);
30 
31   VERIFY_IS_EQUAL(tensor_after_fft.dimension(0), sequence_length);
32   VERIFY_IS_EQUAL(tensor_after_fft_ifft.dimension(0), sequence_length);
33 
34   for (int i = 0; i < sequence_length; ++i) {
35     VERIFY_IS_APPROX(static_cast<float>(tensor(i)), static_cast<float>(std::real(tensor_after_fft_ifft(i))));
36   }
37 }
38 
39 template <int DataLayout>
test_2D_fft_ifft_invariant(int dim0,int dim1)40 static void test_2D_fft_ifft_invariant(int dim0, int dim1) {
41   Tensor<double, 2, DataLayout> tensor(dim0, dim1);
42   tensor.setRandom();
43 
44   array<int, 2> fft;
45   fft[0] = 0;
46   fft[1] = 1;
47 
48   Tensor<std::complex<double>, 2, DataLayout> tensor_after_fft;
49   Tensor<std::complex<double>, 2, DataLayout> tensor_after_fft_ifft;
50 
51   tensor_after_fft = tensor.template fft<Eigen::BothParts, Eigen::FFT_FORWARD>(fft);
52   tensor_after_fft_ifft = tensor_after_fft.template fft<Eigen::BothParts, Eigen::FFT_REVERSE>(fft);
53 
54   VERIFY_IS_EQUAL(tensor_after_fft.dimension(0), dim0);
55   VERIFY_IS_EQUAL(tensor_after_fft.dimension(1), dim1);
56   VERIFY_IS_EQUAL(tensor_after_fft_ifft.dimension(0), dim0);
57   VERIFY_IS_EQUAL(tensor_after_fft_ifft.dimension(1), dim1);
58 
59   for (int i = 0; i < dim0; ++i) {
60     for (int j = 0; j < dim1; ++j) {
61       //std::cout << "[" << i << "][" << j << "]" <<  "  Original data: " << tensor(i,j) << " Transformed data:" << tensor_after_fft_ifft(i,j) << std::endl;
62       VERIFY_IS_APPROX(static_cast<float>(tensor(i,j)), static_cast<float>(std::real(tensor_after_fft_ifft(i,j))));
63     }
64   }
65 }
66 
67 template <int DataLayout>
test_3D_fft_ifft_invariant(int dim0,int dim1,int dim2)68 static void test_3D_fft_ifft_invariant(int dim0, int dim1, int dim2) {
69   Tensor<double, 3, DataLayout> tensor(dim0, dim1, dim2);
70   tensor.setRandom();
71 
72   array<int, 3> fft;
73   fft[0] = 0;
74   fft[1] = 1;
75   fft[2] = 2;
76 
77   Tensor<std::complex<double>, 3, DataLayout> tensor_after_fft;
78   Tensor<std::complex<double>, 3, DataLayout> tensor_after_fft_ifft;
79 
80   tensor_after_fft = tensor.template fft<Eigen::BothParts, Eigen::FFT_FORWARD>(fft);
81   tensor_after_fft_ifft = tensor_after_fft.template fft<Eigen::BothParts, Eigen::FFT_REVERSE>(fft);
82 
83   VERIFY_IS_EQUAL(tensor_after_fft.dimension(0), dim0);
84   VERIFY_IS_EQUAL(tensor_after_fft.dimension(1), dim1);
85   VERIFY_IS_EQUAL(tensor_after_fft.dimension(2), dim2);
86   VERIFY_IS_EQUAL(tensor_after_fft_ifft.dimension(0), dim0);
87   VERIFY_IS_EQUAL(tensor_after_fft_ifft.dimension(1), dim1);
88   VERIFY_IS_EQUAL(tensor_after_fft_ifft.dimension(2), dim2);
89 
90   for (int i = 0; i < dim0; ++i) {
91     for (int j = 0; j < dim1; ++j) {
92       for (int k = 0; k < dim2; ++k) {
93         VERIFY_IS_APPROX(static_cast<float>(tensor(i,j,k)), static_cast<float>(std::real(tensor_after_fft_ifft(i,j,k))));
94       }
95     }
96   }
97 }
98 
99 template <int DataLayout>
test_sub_fft_ifft_invariant(int dim0,int dim1,int dim2,int dim3)100 static void test_sub_fft_ifft_invariant(int dim0, int dim1, int dim2, int dim3) {
101   Tensor<double, 4, DataLayout> tensor(dim0, dim1, dim2, dim3);
102   tensor.setRandom();
103 
104   array<int, 2> fft;
105   fft[0] = 2;
106   fft[1] = 0;
107 
108   Tensor<std::complex<double>, 4, DataLayout> tensor_after_fft;
109   Tensor<double, 4, DataLayout> tensor_after_fft_ifft;
110 
111   tensor_after_fft = tensor.template fft<Eigen::BothParts, Eigen::FFT_FORWARD>(fft);
112   tensor_after_fft_ifft = tensor_after_fft.template fft<Eigen::RealPart, Eigen::FFT_REVERSE>(fft);
113 
114   VERIFY_IS_EQUAL(tensor_after_fft.dimension(0), dim0);
115   VERIFY_IS_EQUAL(tensor_after_fft.dimension(1), dim1);
116   VERIFY_IS_EQUAL(tensor_after_fft.dimension(2), dim2);
117   VERIFY_IS_EQUAL(tensor_after_fft.dimension(3), dim3);
118   VERIFY_IS_EQUAL(tensor_after_fft_ifft.dimension(0), dim0);
119   VERIFY_IS_EQUAL(tensor_after_fft_ifft.dimension(1), dim1);
120   VERIFY_IS_EQUAL(tensor_after_fft_ifft.dimension(2), dim2);
121   VERIFY_IS_EQUAL(tensor_after_fft_ifft.dimension(3), dim3);
122 
123   for (int i = 0; i < dim0; ++i) {
124     for (int j = 0; j < dim1; ++j) {
125       for (int k = 0; k < dim2; ++k) {
126         for (int l = 0; l < dim3; ++l) {
127           VERIFY_IS_APPROX(static_cast<float>(tensor(i,j,k,l)), static_cast<float>(tensor_after_fft_ifft(i,j,k,l)));
128         }
129       }
130     }
131   }
132 }
133 
test_cxx11_tensor_ifft()134 void test_cxx11_tensor_ifft() {
135   CALL_SUBTEST(test_1D_fft_ifft_invariant<ColMajor>(4));
136   CALL_SUBTEST(test_1D_fft_ifft_invariant<ColMajor>(16));
137   CALL_SUBTEST(test_1D_fft_ifft_invariant<ColMajor>(32));
138   CALL_SUBTEST(test_1D_fft_ifft_invariant<ColMajor>(1024*1024));
139 
140   CALL_SUBTEST(test_2D_fft_ifft_invariant<ColMajor>(4,4));
141   CALL_SUBTEST(test_2D_fft_ifft_invariant<ColMajor>(8,16));
142   CALL_SUBTEST(test_2D_fft_ifft_invariant<ColMajor>(16,32));
143   CALL_SUBTEST(test_2D_fft_ifft_invariant<ColMajor>(1024,1024));
144 
145   CALL_SUBTEST(test_3D_fft_ifft_invariant<ColMajor>(4,4,4));
146   CALL_SUBTEST(test_3D_fft_ifft_invariant<ColMajor>(8,16,32));
147   CALL_SUBTEST(test_3D_fft_ifft_invariant<ColMajor>(16,4,8));
148   CALL_SUBTEST(test_3D_fft_ifft_invariant<ColMajor>(256,256,256));
149 
150   CALL_SUBTEST(test_sub_fft_ifft_invariant<ColMajor>(4,4,4,4));
151   CALL_SUBTEST(test_sub_fft_ifft_invariant<ColMajor>(8,16,32,64));
152   CALL_SUBTEST(test_sub_fft_ifft_invariant<ColMajor>(16,4,8,12));
153   CALL_SUBTEST(test_sub_fft_ifft_invariant<ColMajor>(64,64,64,64));
154 }
155