LSSTApplications  10.0-2-g4f67435,11.0.rc2+1,11.0.rc2+12,11.0.rc2+3,11.0.rc2+4,11.0.rc2+5,11.0.rc2+6,11.0.rc2+7,11.0.rc2+8
LSSTDataManagementBasePackage
FourierTransform.cc
Go to the documentation of this file.
1 // -*- c++ -*-
2 /*
3  * LSST Data Management System
4  * Copyright 2008, 2009, 2010, 2011 LSST Corporation.
5  *
6  * This product includes software developed by the
7  * LSST Project (http://www.lsst.org/).
8  *
9  * This program is free software: you can redistribute it and/or modify
10  * it under the terms of the GNU General Public License as published by
11  * the Free Software Foundation, either version 3 of the License, or
12  * (at your option) any later version.
13  *
14  * This program is distributed in the hope that it will be useful,
15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17  * GNU General Public License for more details.
18  *
19  * You should have received a copy of the LSST License Statement and
20  * the GNU General Public License along with this program. If not,
21  * see <http://www.lsstcorp.org/LegalNotices/>.
22  */
23 #include "ndarray/fft/FFTWTraits.h"
25 
26 namespace ndarray {
27 
28 template <typename T, int N>
29 template <int M>
30 Array<typename FourierTransform<T,N>::ElementX,M,M>
33  return Array<ElementX,M,M>(external(xOwner.get(), shape, ROW_MAJOR, xOwner));
34 }
35 
36 template <typename T, int N>
37 template <int M>
40  Vector<int,M> kShape(shape);
43  return Array<ElementK,M,M>(external(kOwner.get(), kShape, ROW_MAJOR, kOwner));
44 }
45 
46 template <typename T, int N>
47 template <int M>
48 void
50  Vector<int,M> const & shape,
53 ) {
54  if (x.empty()) x = initializeX(shape);
55  if (k.empty()) k = initializeK(shape);
56  NDARRAY_ASSERT(x.getShape() == shape);
57  NDARRAY_ASSERT(std::equal(shape.begin(), shape.end()-1, k.getShape().begin()));
58 }
59 
60 template <typename T, int N>
63  Index const & shape,
66 ) {
67  initialize(shape,x,k);
68  return Ptr(
69  new FourierTransform(
71  N, shape.begin(), 1,
72  x.getData(), NULL, 1, 0,
73  k.getData(), NULL, 1, 0,
74  FFTW_MEASURE | FFTW_DESTROY_INPUT
75  ),
76  x.getManager(),
77  k.getManager()
78  )
79  );
80 }
81 
82 template <typename T, int N>
85  Index const & shape,
86  typename FourierTransform<T,N>::ArrayK & k,
88 ) {
89  initialize(shape,x,k);
90  return Ptr(
91  new FourierTransform(
93  N, shape.begin(), 1,
94  k.getData(), NULL, 1, 0,
95  x.getData(), NULL, 1, 0,
96  FFTW_MEASURE | FFTW_DESTROY_INPUT
97  ),
98  x.getManager(),
99  k.getManager()
100  )
101  );
102 }
103 
104 template <typename T, int N>
107  MultiplexIndex const & shape,
110 ) {
111  initialize(shape,x,k);
112  return Ptr(
113  new FourierTransform(
115  N, shape.begin()+1, shape[0],
116  x.getData(), NULL, 1, x.template getStride<0>(),
117  k.getData(), NULL, 1, k.template getStride<0>(),
118  FFTW_MEASURE | FFTW_DESTROY_INPUT
119  ),
120  x.getManager(),
121  k.getManager()
122  )
123  );
124 }
125 
126 template <typename T, int N>
129  MultiplexIndex const & shape,
132 ) {
133  initialize(shape,x,k);
134  return Ptr(
135  new FourierTransform(
137  N, shape.begin()+1, shape[0],
138  k.getData(), NULL, 1, k.template getStride<0>(),
139  x.getData(), NULL, 1, x.template getStride<0>(),
140  FFTW_MEASURE | FFTW_DESTROY_INPUT
141  ),
142  x.getManager(),
143  k.getManager()
144  )
145  );
146 }
147 
148 template <typename T, int N>
151  reinterpret_cast<typename detail::FFTWTraits<T>::Plan>(_plan)
152  );
153 }
154 
155 template <typename T, int N>
158  reinterpret_cast<typename detail::FFTWTraits<T>::Plan>(_plan)
159  );
160 }
161 
162 } // namespace ndarray
static Array< ElementX, M, M > initializeX(Vector< int, M > const &shape)
Create a new real-space array with the given real-space shape.
static void initialize(Vector< int, M > const &shape, Array< ElementX, M, M > &x, Array< ElementK, M, M > &k)
Initialize, as necessary, a pair of arrays with the given real-space shape.
static Ptr planMultiplexForward(MultiplexIndex const &shape, MultiplexArrayX &x, MultiplexArrayK &k)
Create a plan for forward-transforming a sequence of nested N-dimensional arrays. ...
static Ptr planMultiplexInverse(MultiplexIndex const &shape, MultiplexArrayK &k, MultiplexArrayX &x)
Create a plan for inverse-transforming a sequence of nested N-dimensional arrays. ...
Definitions for FourierTransform.
boost::shared_ptr< ElementX > OwnerX
detail::ExternalInitializer< T, N, Owner > external(T *data, Vector< int, N > const &shape, Vector< int, N > const &strides, Owner const &owner)
Create an expression that initializes an Array with externally allocated memory.
T product() const
Return the product of all elements.
Definition: Vector.h:216
#define NDARRAY_ASSERT(ARG)
Definition: ndarray_fwd.h:51
A fixed-size 1D array class.
Definition: Vector.h:93
Index getShape() const
Return a Vector of the sizes of all dimensions.
Definition: ArrayBase.h:136
iterator begin()
Return an iterator to the beginning of the Vector.
Definition: Vector.h:120
static Array< ElementK, M, M > initializeK(Vector< int, M > const &shape)
Create a new Fourier-space array with the given real-space shape.
void execute()
Execute the FFTW plan.
iterator end()
Return an iterator to the end of the Vector.
Definition: Vector.h:124
Element * getData() const
Return a raw pointer to the first element of the array.
Definition: ArrayBase.h:117
int x
boost::shared_ptr< FourierTransform > Ptr
bool empty() const
Return true if the first dimension has no elements.
A multidimensional strided array.
Definition: Array.h:47
static Ptr planForward(Index const &shape, ArrayX &x, ArrayK &k)
Create a plan for forward-transforming a single N-dimensional array.
< unspecified-expression-type > equal(ExpressionBase< Operand > const &operand, Scalar const &scalar)
Definition: operators.h:796
Traits classes that wrap FFTW in a template-friendly interface.
static Ptr planInverse(Index const &shape, ArrayK &k, ArrayX &x)
Create a plan for inverse-transforming a single N-dimensional array.
boost::shared_ptr< ElementK > OwnerK
Manager::Ptr getManager() const
Return the opaque object responsible for memory management.
Definition: ArrayBase.h:123
A wrapper for FFTW plans for fast Fourier transforms.