1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
|
// This file is part of Eigen, a lightweight C++ template library
// for linear algebra.
//
// Copyright (C) 2009 Mark Borgerding mark a borgerding net
//
// Eigen 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 3 of the License, or (at your option) any later version.
//
// Alternatively, you can redistribute it and/or
// modify it under the terms of the GNU General Public License as
// published by the Free Software Foundation; either version 2 of
// the License, or (at your option) any later version.
//
// Eigen 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 or the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public
// License and a copy of the GNU General Public License along with
// Eigen. If not, see <http://www.gnu.org/licenses/>.
#ifndef EIGEN_FFT_H
#define EIGEN_FFT_H
// ei_kissfft_impl: small, free, reasonably efficient default, derived from kissfft
#include "src/FFT/ei_kissfft_impl.h"
#define DEFAULT_FFT_IMPL ei_kissfft_impl
// FFTW: faster, GPL -- incompatible with Eigen in LGPL form, bigger code size
#ifdef FFTW_ESTIMATE // definition of FFTW_ESTIMATE indicates the caller has included fftw3.h, we can use FFTW routines
#include "src/FFT/ei_fftw_impl.h"
#undef DEFAULT_FFT_IMPL
#define DEFAULT_FFT_IMPL ei_fftw_impl
#endif
// intel Math Kernel Library: fastest, commercial -- incompatible with Eigen in GPL form
#ifdef _MKL_DFTI_H_ // mkl_dfti.h has been included, we can use MKL FFT routines
// TODO
// #include "src/FFT/ei_imkl_impl.h"
// #define DEFAULT_FFT_IMPL ei_imkl_impl
#endif
namespace Eigen {
template <typename _Scalar,
typename _Traits=DEFAULT_FFT_IMPL<_Scalar>
>
class FFT
{
public:
typedef _Traits traits_type;
typedef typename traits_type::Scalar Scalar;
typedef typename traits_type::Complex Complex;
FFT(const traits_type & traits=traits_type() ) :m_traits(traits) { }
template <typename _Input>
void fwd( Complex * dst, const _Input * src, int nfft)
{
m_traits.fwd(dst,src,nfft);
}
template <typename _Input>
void fwd( std::vector<Complex> & dst, const std::vector<_Input> & src)
{
dst.resize( src.size() );
fwd( &dst[0],&src[0],src.size() );
}
template<typename InputDerived, typename ComplexDerived>
void fwd( MatrixBase<ComplexDerived> & dst, const MatrixBase<InputDerived> & src)
{
EIGEN_STATIC_ASSERT_VECTOR_ONLY(InputDerived)
EIGEN_STATIC_ASSERT_VECTOR_ONLY(ComplexDerived)
EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(ComplexDerived,InputDerived) // size at compile-time
EIGEN_STATIC_ASSERT((ei_is_same_type<typename ComplexDerived::Scalar, Complex>::ret),
YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY)
EIGEN_STATIC_ASSERT(int(InputDerived::Flags)&int(ComplexDerived::Flags)&DirectAccessBit,
THIS_METHOD_IS_ONLY_FOR_EXPRESSIONS_WITH_DIRECT_MEMORY_ACCESS_SUCH_AS_MAP_OR_PLAIN_MATRICES)
dst.derived().resize( src.size() );
fwd( &dst[0],&src[0],src.size() );
}
template <typename _Output>
void inv( _Output * dst, const Complex * src, int nfft)
{
m_traits.inv( dst,src,nfft );
}
template <typename _Output>
void inv( std::vector<_Output> & dst, const std::vector<Complex> & src)
{
dst.resize( src.size() );
inv( &dst[0],&src[0],src.size() );
}
template<typename OutputDerived, typename ComplexDerived>
void inv( MatrixBase<OutputDerived> & dst, const MatrixBase<ComplexDerived> & src)
{
EIGEN_STATIC_ASSERT_VECTOR_ONLY(OutputDerived)
EIGEN_STATIC_ASSERT_VECTOR_ONLY(ComplexDerived)
EIGEN_STATIC_ASSERT_SAME_VECTOR_SIZE(ComplexDerived,OutputDerived) // size at compile-time
EIGEN_STATIC_ASSERT((ei_is_same_type<typename ComplexDerived::Scalar, Complex>::ret),
YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY)
EIGEN_STATIC_ASSERT(int(OutputDerived::Flags)&int(ComplexDerived::Flags)&DirectAccessBit,
THIS_METHOD_IS_ONLY_FOR_EXPRESSIONS_WITH_DIRECT_MEMORY_ACCESS_SUCH_AS_MAP_OR_PLAIN_MATRICES)
dst.derived().resize( src.size() );
inv( &dst[0],&src[0],src.size() );
}
// TODO: multi-dimensional FFTs
// TODO: handle Eigen MatrixBase
// ---> i added fwd and inv specializations above + unit test, is this enough? (bjacob)
traits_type & traits() {return m_traits;}
private:
traits_type m_traits;
};
#undef DEFAULT_FFT_IMPL
}
#endif
|