LAL  7.5.0.1-ec27e42
ComplexFFT.h
Go to the documentation of this file.
1 /*
2 * Copyright (C) 2007 Duncan Brown, Jolien Creighton
3 *
4 * This program is free software; you can redistribute it and/or modify
5 * it under the terms of the GNU General Public License as published by
6 * the Free Software Foundation; either version 2 of the License, or
7 * (at your option) any later version.
8 *
9 * This program is distributed in the hope that it will be useful,
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 * GNU General Public License for more details.
13 *
14 * You should have received a copy of the GNU General Public License
15 * along with with program; see the file COPYING. If not, write to the
16 * Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
17 * MA 02110-1301 USA
18 */
19 
20 #ifndef _COMPLEXFFT_H
21 #define _COMPLEXFFT_H
22 
23 #include <lal/LALDatatypes.h>
24 
25 #if defined(__cplusplus)
26 extern "C" {
27 #elif 0
28 } /* so that editors will match preceding brace */
29 #endif
30 
31 /**
32  * \defgroup ComplexFFT_h Header ComplexFFT_h
33  * \ingroup lal_fft
34  *
35  * \brief Performs complex-to-complex FFTs.
36  *
37  * ### Synopsis ###
38  *
39  * \code
40  * #include <lal/ComplexFFT.h>
41  * \endcode
42  *
43  * Perform complex-to-complex fast Fourier transforms of vectors using the
44  * package FFTW \cite fj_1998 .
45  *
46  */
47 /** @{ */
48 
49 /** \name Error Codes */
50 /**@{*/
51 #define COMPLEXFFTH_ENULL 1 /**< Null pointer */
52 #define COMPLEXFFTH_ENNUL 2 /**< Non-null pointer */
53 #define COMPLEXFFTH_ESIZE 4 /**< Invalid input size */
54 #define COMPLEXFFTH_ESZMM 8 /**< Size mismatch */
55 #define COMPLEXFFTH_ESLEN 16 /**< Invalid/mismatched sequence lengths */
56 #define COMPLEXFFTH_ESAME 32 /**< Input/Output data vectors are the same */
57 #define COMPLEXFFTH_EALOC 64 /**< Memory allocation failed */
58 #define COMPLEXFFTH_EFFTW 128 /**< Error in FFTW */
59 #define COMPLEXFFTH_ESNGL 256 /**< FFTW library is not single-precision */
60 #define COMPLEXFFTH_EINTL 512 /**< Error in Intel FFT library */
61 #define COMPLEXFFTH_ESIGN 1024 /**< Unknown sign of transform in plan */
62 /**@}*/
63 
64 /** \cond DONT_DOXYGEN */
65 #define COMPLEXFFTH_MSGENULL "Null pointer"
66 #define COMPLEXFFTH_MSGENNUL "Non-null pointer"
67 #define COMPLEXFFTH_MSGESIZE "Invalid input size"
68 #define COMPLEXFFTH_MSGESZMM "Size mismatch"
69 #define COMPLEXFFTH_MSGESLEN "Invalid/mismatched sequence lengths"
70 #define COMPLEXFFTH_MSGESAME "Input/Output data vectors are the same"
71 #define COMPLEXFFTH_MSGEALOC "Memory allocation failed"
72 #define COMPLEXFFTH_MSGEFFTW "Error in FFTW"
73 #define COMPLEXFFTH_MSGESNGL "FFTW library is not single-precision"
74 #define COMPLEXFFTH_MSGEINTL "Error in Intel FFT library"
75 #define COMPLEXFFTH_MSGESIGN "Unknown sign of transform in plan"
76 /** \endcond */
77 
78 /** Plan to perform FFT of COMPLEX8 data */
79 typedef struct tagCOMPLEX8FFTPlan COMPLEX8FFTPlan;
80 /** Plan to perform FFT of COMPLEX16 data */
81 typedef struct tagCOMPLEX16FFTPlan COMPLEX16FFTPlan;
82 #define tagComplexFFTPlan tagCOMPLEX8FFTPlan
83 #define ComplexFFTPlan COMPLEX8FFTPlan
84 
85 #ifdef SWIG /* SWIG interface directives */
86 SWIGLAL(VIEWIN_ARRAYS(COMPLEX8Vector, output));
87 SWIGLAL(VIEWIN_ARRAYS(COMPLEX16Vector, output));
88 #endif /* SWIG */
89 
90 /*
91  *
92  * XLAL COMPLEX8 functions
93  *
94  */
95 
96 /**
97  * Returns a new COMPLEX8FFTPlan
98  * A COMPLEX8FFTPlan is required to perform a FFT that involves complex data.
99  * A different plan is required for each size of the complex data vectors
100  * and for each direction of transform (forward or reverse).
101  * A forward transform performs
102  * \f[Z[k] = \sum_{j=0}^{N-1} e^{-2\pi ijk/N}\,z[j]\f]
103  * where N, the size of the transform, is the length of the vectors z and Z.
104  * A reverse transform performs
105  * \f[w[j] = \sum_{k=0}^{N-1} e^{+2\pi ijk/N}\,W[k]\f]
106  * where N, the size of the transform, is the length of the vectors w and W.
107  *
108  * @note
109  * The reverse transform of the forward transform of some data is
110  * equal to N times the original data (we therefore call it a "reverse"
111  * transform rather than an "inverse" transform).
112  *
113  * @param[in] size The number of points in the complex data.
114  * @param[in] fwdflg Set non-zero for a forward FFT plan;
115  * otherwise create a reverse plan
116  * @param[in] measurelvl Measurement level for plan creation:
117  * - 0: no measurement, just estimate the plan;
118  * - 1: measure the best plan;
119  * - 2: perform a lengthy measurement of the best plan;
120  * - 3: perform an exhasutive measurement of the best plan.
121  * @return A pointer to an allocated \c COMPLEX8FFTPlan structure is returned
122  * upon successful completion. Otherwise, a \c NULL pointer is returned
123  * and \c xlalErrno is set to indicate the error.
124  * @par Errors:
125  * The \c XLALCreateCOMPLEX8Plan() function shall fail if:
126  * - [\c XLAL_EBADLEN] The size of the requested plan is 0.
127  * - [\c XLAL_ENOMEM] Insufficient storage space is available.
128  * - [\c XLAL_EFAILED] The call to the underlying FFTW routine failed.
129  * .
130  */
131 COMPLEX8FFTPlan * XLALCreateCOMPLEX8FFTPlan( UINT4 size, int fwdflg, int measurelvl );
132 
133 /**
134  * Returns a new COMPLEX8FFTPlan for a forward transform
135  *
136  * A COMPLEX8FFTPlan is required to perform a FFT that involves complex data.
137  * A different plan is required for each size of the complex data vectors.
138  * A forward transform performs
139  * \f[Z[k] = \sum_{j=0}^{N-1} e^{-2\pi ijk/N}\,z[j]\f]
140  * where N, the size of the transform, is the length of the vector z and Z.
141  *
142  * @param[in] size The number of points in the complex data.
143  * @param[in] measurelvl Measurement level for plan creation:
144  * - 0: no measurement, just estimate the plan;
145  * - 1: measure the best plan;
146  * - 2: perform a lengthy measurement of the best plan;
147  * - 3: perform an exhasutive measurement of the best plan.
148  * @return A pointer to an allocated \c COMPLEX8FFTPlan structure is returned
149  * upon successful completion. Otherwise, a \c NULL pointer is returned
150  * and \c xlalErrno is set to indicate the error.
151  * @par Errors:
152  * The \c XLALCreateForwardCOMPLEX8Plan() function shall fail if:
153  * - [\c XLAL_EBADLEN] The size of the requested plan is 0.
154  * - [\c XLAL_ENOMEM] Insufficient storage space is available.
155  * - [\c XLAL_EFAILED] The call to the underlying FFTW routine failed.
156  * .
157  */
158 COMPLEX8FFTPlan * XLALCreateForwardCOMPLEX8FFTPlan( UINT4 size, int measurelvl );
159 
160 /**
161  * Returns a new COMPLEX8FFTPlan for a reverse transform
162  *
163  * A COMPLEX8FFTPlan is required to perform a FFT that involves complex data.
164  * A reverse transform performs
165  * \f[w[j] = \sum_{k=0}^{N-1} e^{+2\pi ijk/N}\,W[k]\f]
166  * where N, the size of the transform, is the length of the vectors w and W.
167  *
168  * @note
169  * The reverse transform of the forward transform of some data is
170  * equal to N times the original data (we therefore call it a "reverse"
171  * transform rather than an "inverse" transform).
172  *
173  * @param[in] size The number of points in the complex data.
174  * @param[in] measurelvl Measurement level for plan creation:
175  * - 0: no measurement, just estimate the plan;
176  * - 1: measure the best plan;
177  * - 2: perform a lengthy measurement of the best plan;
178  * - 3: perform an exhasutive measurement of the best plan.
179  * @return A pointer to an allocated \c COMPLEX8FFTPlan structure is returned
180  * upon successful completion. Otherwise, a \c NULL pointer is returned
181  * and \c xlalErrno is set to indicate the error.
182  * @par Errors:
183  * The \c XLALCreateReverseCOMPLEX8Plan() function shall fail if:
184  * - [\c XLAL_EBADLEN] The size of the requested plan is 0.
185  * - [\c XLAL_ENOMEM] Insufficient storage space is available.
186  * - [\c XLAL_EFAILED] The call to the underlying FFTW routine failed.
187  * .
188  */
189 COMPLEX8FFTPlan * XLALCreateReverseCOMPLEX8FFTPlan( UINT4 size, int measurelvl );
190 
191 /**
192  * Destroys a COMPLEX8FFTPlan
193  * @param[in] plan A pointer to the COMPLEX8FFTPlan to be destroyed.
194  */
195 void XLALDestroyCOMPLEX8FFTPlan( COMPLEX8FFTPlan *plan );
196 
197 /**
198  * Perform a COMPLEX8Vector to COMPLEX8Vector FFT
199  *
200  * This routine computes
201  * \f[Z[k] = \sum_{j=0}^{N-1} e^{\mp2\pi ijk/N}\,z[j],\f]
202  * and where the minus sign is used if a forward plan is provided as the argument
203  * and the plus sign is used if a reverse plan is provided as the argument;
204  * here N is the length of the input and output vectors z and Z.
205  *
206  * @param[out] output The complex output data vector Z of length N
207  * @param[in] input The input complex data vector z of length N
208  * @param[in] plan The FFT plan to use for the transform
209  * @note
210  * The input and output vectors must be distinct.
211  * @return 0 upon successful completion or non-zero upon failure.
212  * @par Errors:
213  * The \c XLALCOMPLEX8VectorFFT() function shall fail if:
214  * - [\c XLAL_EFAULT] A \c NULL pointer is provided as one of the arguments.
215  * - [\c XLAL_EINVAL] A argument is invalid or the input and output data
216  * vectors are the same.
217  * - [\c XLAL_EBADLEN] The input vector, output vector, and plan size are
218  * incompatible.
219  * - [\c XLAL_ENOMEM] Insufficient storage space is available.
220  * .
221  */
222 int XLALCOMPLEX8VectorFFT( COMPLEX8Vector * _LAL_RESTRICT_ output, const COMPLEX8Vector * _LAL_RESTRICT_ input, const COMPLEX8FFTPlan *plan );
223 
224 /*
225  *
226  * XLAL COMPLEX16 functions
227  *
228  */
229 
230 /**
231  * Returns a new COMPLEX16FFTPlan
232  *
233  * A COMPLEX16FFTPlan is required to perform a FFT that involves complex data.
234  * A different plan is required for each size of the complex data vectors
235  * and for each direction of transform (forward or reverse).
236  * A forward transform performs
237  * \f[Z[k] = \sum_{j=0}^{N-1} e^{-2\pi ijk/N}\,z[j]\f]
238  * where N, the size of the transform, is the length of the vectors z and Z.
239  * A reverse transform performs
240  * \f[w[j] = \sum_{k=0}^{N-1} e^{+2\pi ijk/N}\,W[k]\f]
241  * where N, the size of the transform, is the length of the vectors w and W.
242  *
243  * @note
244  * The reverse transform of the forward transform of some data is
245  * equal to N times the original data (we therefore call it a "reverse"
246  * transform rather than an "inverse" transform).
247  *
248  * @param[in] size The number of points in the complex data.
249  * @param[in] fwdflg Set non-zero for a forward FFT plan;
250  * otherwise create a reverse plan
251  * @param[in] measurelvl Measurement level for plan creation:
252  * - 0: no measurement, just estimate the plan;
253  * - 1: measure the best plan;
254  * - 2: perform a lengthy measurement of the best plan;
255  * - 3: perform an exhasutive measurement of the best plan.
256  * @return A pointer to an allocated \c COMPLEX16FFTPlan structure is returned
257  * upon successful completion. Otherwise, a \c NULL pointer is returned
258  * and \c xlalErrno is set to indicate the error.
259  * @par Errors:
260  * The \c XLALCreateCOMPLEX16Plan() function shall fail if:
261  * - [\c XLAL_EBADLEN] The size of the requested plan is 0.
262  * - [\c XLAL_ENOMEM] Insufficient storage space is available.
263  * - [\c XLAL_EFAILED] The call to the underlying FFTW routine failed.
264  * .
265  */
266 COMPLEX16FFTPlan * XLALCreateCOMPLEX16FFTPlan( UINT4 size, int fwdflg, int measurelvl );
267 
268 /**
269  * Returns a new COMPLEX16FFTPlan for a forward transform
270  *
271  * A COMPLEX16FFTPlan is required to perform a FFT that involves complex data.
272  * A different plan is required for each size of the complex data vectors.
273  * A forward transform performs
274  * \f[Z[k] = \sum_{j=0}^{N-1} e^{-2\pi ijk/N}\,z[j]\f]
275  * where N, the size of the transform, is the length of the vector z and Z.
276  *
277  * @param[in] size The number of points in the complex data.
278  * @param[in] measurelvl Measurement level for plan creation:
279  * - 0: no measurement, just estimate the plan;
280  * - 1: measure the best plan;
281  * - 2: perform a lengthy measurement of the best plan;
282  * - 3: perform an exhasutive measurement of the best plan.
283  * @return A pointer to an allocated \c COMPLEX16FFTPlan structure is returned
284  * upon successful completion. Otherwise, a \c NULL pointer is returned
285  * and \c xlalErrno is set to indicate the error.
286  * @par Errors:
287  * The \c XLALCreateForwardCOMPLEX16Plan() function shall fail if:
288  * - [\c XLAL_EBADLEN] The size of the requested plan is 0.
289  * - [\c XLAL_ENOMEM] Insufficient storage space is available.
290  * - [\c XLAL_EFAILED] The call to the underlying FFTW routine failed.
291  * .
292  */
293 COMPLEX16FFTPlan * XLALCreateForwardCOMPLEX16FFTPlan( UINT4 size, int measurelvl );
294 
295 /**
296  * Returns a new COMPLEX16FFTPlan for a reverse transform
297  *
298  * A COMPLEX16FFTPlan is required to perform a FFT that involves complex data.
299  * A reverse transform performs
300  * \f[w[j] = \sum_{k=0}^{N-1} e^{+2\pi ijk/N}\,W[k]\f]
301  * where N, the size of the transform, is the length of the vectors w and W.
302  *
303  * @note
304  * The reverse transform of the forward transform of some data is
305  * equal to N times the original data (we therefore call it a "reverse"
306  * transform rather than an "inverse" transform).
307  *
308  * @param[in] size The number of points in the complex data.
309  * @param[in] measurelvl Measurement level for plan creation:
310  * - 0: no measurement, just estimate the plan;
311  * - 1: measure the best plan;
312  * - 2: perform a lengthy measurement of the best plan;
313  * - 3: perform an exhasutive measurement of the best plan.
314  * @return A pointer to an allocated \c COMPLEX16FFTPlan structure is returned
315  * upon successful completion. Otherwise, a \c NULL pointer is returned
316  * and \c xlalErrno is set to indicate the error.
317  * @par Errors:
318  * The \c XLALCreateReverseCOMPLEX16Plan() function shall fail if:
319  * - [\c XLAL_EBADLEN] The size of the requested plan is 0.
320  * - [\c XLAL_ENOMEM] Insufficient storage space is available.
321  * - [\c XLAL_EFAILED] The call to the underlying FFTW routine failed.
322  * .
323  */
324 COMPLEX16FFTPlan * XLALCreateReverseCOMPLEX16FFTPlan( UINT4 size, int measurelvl );
325 
326 /**
327  * Destroys a COMPLEX16FFTPlan
328  * @param[in] plan A pointer to the COMPLEX16FFTPlan to be destroyed.
329  */
330 void XLALDestroyCOMPLEX16FFTPlan( COMPLEX16FFTPlan *plan );
331 
332 /**
333  * Perform a COMPLEX16Vector to COMPLEX16Vector FFT
334  *
335  * This routine computes
336  * \f[Z[k] = \sum_{j=0}^{N-1} e^{\mp2\pi ijk/N}\,z[j],\f]
337  * and where the minus sign is used if a forward plan is provided as the argument
338  * and the plus sign is used if a reverse plan is provided as the argument;
339  * here N is the length of the input and output vectors z and Z.
340  *
341  * @param[out] output The complex output data vector Z of length N
342  * @param[in] input The input complex data vector z of length N
343  * @param[in] plan The FFT plan to use for the transform
344  * @note
345  * The input and output vectors must be distinct.
346  * @return 0 upon successful completion or non-zero upon failure.
347  * @par Errors:
348  * The \c XLALCOMPLEX16VectorFFT() function shall fail if:
349  * - [\c XLAL_EFAULT] A \c NULL pointer is provided as one of the arguments.
350  * - [\c XLAL_EINVAL] A argument is invalid or the input and output data
351  * vectors are the same.
352  * - [\c XLAL_EBADLEN] The input vector, output vector, and plan size are
353  * incompatible.
354  * - [\c XLAL_ENOMEM] Insufficient storage space is available.
355  * .
356  */
357 int XLALCOMPLEX16VectorFFT( COMPLEX16Vector * _LAL_RESTRICT_ output, const COMPLEX16Vector * _LAL_RESTRICT_ input, const COMPLEX16FFTPlan *plan );
358 
359 /** @} */
360 
361 #if 0
362 { /* so that editors will match succeeding brace */
363 #elif defined(__cplusplus)
364 }
365 #endif
366 
367 #endif /* _COMPLEXFFT_H */
#define _LAL_RESTRICT_
Definition: LALStddef.h:35
COMPLEX8FFTPlan * XLALCreateCOMPLEX8FFTPlan(UINT4 size, int fwdflg, int measurelvl)
Returns a new COMPLEX8FFTPlan A COMPLEX8FFTPlan is required to perform a FFT that involves complex da...
COMPLEX8FFTPlan * XLALCreateReverseCOMPLEX8FFTPlan(UINT4 size, int measurelvl)
Returns a new COMPLEX8FFTPlan for a reverse transform.
int XLALCOMPLEX8VectorFFT(COMPLEX8Vector *_LAL_RESTRICT_ output, const COMPLEX8Vector *_LAL_RESTRICT_ input, const COMPLEX8FFTPlan *plan)
Perform a COMPLEX8Vector to COMPLEX8Vector FFT.
void XLALDestroyCOMPLEX16FFTPlan(COMPLEX16FFTPlan *plan)
Destroys a COMPLEX16FFTPlan.
COMPLEX8FFTPlan * XLALCreateForwardCOMPLEX8FFTPlan(UINT4 size, int measurelvl)
Returns a new COMPLEX8FFTPlan for a forward transform.
COMPLEX16FFTPlan * XLALCreateForwardCOMPLEX16FFTPlan(UINT4 size, int measurelvl)
Returns a new COMPLEX16FFTPlan for a forward transform.
int XLALCOMPLEX16VectorFFT(COMPLEX16Vector *_LAL_RESTRICT_ output, const COMPLEX16Vector *_LAL_RESTRICT_ input, const COMPLEX16FFTPlan *plan)
Perform a COMPLEX16Vector to COMPLEX16Vector FFT.
void XLALDestroyCOMPLEX8FFTPlan(COMPLEX8FFTPlan *plan)
Destroys a COMPLEX8FFTPlan.
COMPLEX16FFTPlan * XLALCreateCOMPLEX16FFTPlan(UINT4 size, int fwdflg, int measurelvl)
Returns a new COMPLEX16FFTPlan.
COMPLEX16FFTPlan * XLALCreateReverseCOMPLEX16FFTPlan(UINT4 size, int measurelvl)
Returns a new COMPLEX16FFTPlan for a reverse transform.
uint32_t UINT4
Four-byte unsigned integer.
Vector of type COMPLEX16, see DATATYPE-Vector types for more details.
Definition: LALDatatypes.h:172
Vector of type COMPLEX8, see DATATYPE-Vector types for more details.
Definition: LALDatatypes.h:163
Plan to perform an FFT of COMPLEX16 data.
Definition: ComplexFFT.c:137
fftw_plan plan
the FFTW plan
Definition: ComplexFFT.c:140
UINT4 size
length of the complex data vector for this plan
Definition: ComplexFFT.c:139
Plan to perform an FFT of COMPLEX8 data.
Definition: ComplexFFT.c:126
void output(int gps_sec, int output_type)
Definition: tconvert.c:440