GNU Octave  3.8.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
oct-fftw.h
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 2001-2013 John W. Eaton
4 
5 This file is part of Octave.
6 
7 Octave is free software; you can redistribute it and/or modify it
8 under the terms of the GNU General Public License as published by the
9 Free Software Foundation; either version 3 of the License, or (at your
10 option) any later version.
11 
12 Octave is distributed in the hope that it will be useful, but WITHOUT
13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
15 for more details.
16 
17 You should have received a copy of the GNU General Public License
18 along with Octave; see the file COPYING. If not, see
19 <http://www.gnu.org/licenses/>.
20 
21 */
22 
23 #if !defined (octave_oct_fftw_h)
24 #define octave_oct_fftw_h 1
25 
26 #include <cstddef>
27 
28 #if defined (HAVE_FFTW3_H)
29 #include <fftw3.h>
30 #endif
31 
32 #include "oct-cmplx.h"
33 #include "dim-vector.h"
34 
35 #if defined (HAVE_FFTW)
36 
37 class
38 OCTAVE_API
40 {
41 protected:
42 
43  octave_fftw_planner (void);
44 
45 public:
46 
47  ~octave_fftw_planner (void);
48 
50  {
51  UNKNOWN = -1,
56  HYBRID
57  };
58 
59  static bool instance_ok (void);
60 
61  static fftw_plan
62  create_plan (int dir, const int rank, const dim_vector dims,
63  octave_idx_type howmany, octave_idx_type stride,
64  octave_idx_type dist, const Complex *in,
65  Complex *out)
66  {
67  static fftw_plan dummy;
68 
69  return instance_ok ()
70  ? instance->do_create_plan (dir, rank, dims, howmany, stride,
71  dist, in, out)
72  : dummy;
73  }
74 
75  static fftw_plan
76  create_plan (const int rank, const dim_vector dims,
77  octave_idx_type howmany, octave_idx_type stride,
78  octave_idx_type dist, const double *in, Complex *out)
79  {
80  static fftw_plan dummy;
81 
82  return instance_ok ()
83  ? instance->do_create_plan (rank, dims, howmany, stride, dist,
84  in, out)
85  : dummy;
86  }
87 
88  static FftwMethod method (void)
89  {
90  static FftwMethod dummy;
91 
92  return instance_ok () ? instance->do_method () : dummy;
93  }
94 
95  static FftwMethod method (FftwMethod _meth)
96  {
97  static FftwMethod dummy;
98 
99  return instance_ok () ? instance->do_method (_meth) : dummy;
100  }
101 
102 #if defined (HAVE_FFTW3F_THREADS)
103  static void threads (int _nthreads)
104  {
105  if (instance_ok () && _nthreads != threads ())
106  {
107  instance->nthreads = _nthreads;
108  fftw_plan_with_nthreads (_nthreads);
109  //Clear the current plans
110  instance->rplan = instance->plan[0] = instance->plan[1] = 0;
111  }
112  }
113 
114  static int threads ()
115  {
116  return instance_ok () ? instance->nthreads : 0;
117  }
118 #endif
119 
120 private:
121 
122  // No copying!
123 
125 
126  octave_fftw_planner& operator = (const octave_fftw_planner&);
127 
129 
130  static void cleanup_instance (void) { delete instance; instance = 0; }
131 
132  fftw_plan
133  do_create_plan (int dir, const int rank, const dim_vector dims,
134  octave_idx_type howmany, octave_idx_type stride,
135  octave_idx_type dist, const Complex *in,
136  Complex *out);
137 
138  fftw_plan
139  do_create_plan (const int rank, const dim_vector dims,
140  octave_idx_type howmany, octave_idx_type stride,
141  octave_idx_type dist, const double *in, Complex *out);
142 
143  FftwMethod do_method (void);
144 
145  FftwMethod do_method (FftwMethod _meth);
146 
148 
149  // FIXME: perhaps this should be split into two classes?
150 
151  // Plan for fft and ifft of complex values
152  fftw_plan plan[2];
153 
154  // dist
156 
157  // stride
159 
160  // rank
161  int r[2];
162 
163  // howmany
165 
166  // dims
168 
169  bool simd_align[2];
170  bool inplace[2];
171 
172  // Plan for fft of real values
173  fftw_plan rplan;
174 
175  // dist
177 
178  // stride
180 
181  // rank
182  int rr;
183 
184  // howmany
186 
187  // dims
189 
191 
192 #if defined (HAVE_FFTW3_THREADS)
193  //number of threads when compiled with Multi-threading support
194  int nthreads;
195 #endif
196 };
197 
198 class
199 OCTAVE_API
201 {
202 protected:
203 
205 
206 public:
207 
209 
211  {
212  UNKNOWN = -1,
217  HYBRID
218  };
219 
220  static bool instance_ok (void);
221 
222  static fftwf_plan
223  create_plan (int dir, const int rank, const dim_vector dims,
224  octave_idx_type howmany, octave_idx_type stride,
225  octave_idx_type dist, const FloatComplex *in,
226  FloatComplex *out)
227  {
228  static fftwf_plan dummy;
229 
230  return instance_ok ()
231  ? instance->do_create_plan (dir, rank, dims, howmany, stride,
232  dist, in, out)
233  : dummy;
234  }
235 
236  static fftwf_plan
237  create_plan (const int rank, const dim_vector dims,
238  octave_idx_type howmany, octave_idx_type stride,
239  octave_idx_type dist, const float *in, FloatComplex *out)
240  {
241  static fftwf_plan dummy;
242 
243  return instance_ok ()
244  ? instance->do_create_plan (rank, dims, howmany, stride, dist,
245  in, out)
246  : dummy;
247  }
248 
249  static FftwMethod method (void)
250  {
251  static FftwMethod dummy;
252 
253  return instance_ok () ? instance->do_method () : dummy;
254  }
255 
256  static FftwMethod method (FftwMethod _meth)
257  {
258  static FftwMethod dummy;
259 
260  return instance_ok () ? instance->do_method (_meth) : dummy;
261  }
262 
263 #if defined (HAVE_FFTW3F_THREADS)
264  static void threads (int _nthreads)
265  {
266  if (instance_ok () && _nthreads != threads ())
267  {
268  instance->nthreads = _nthreads;
269  fftwf_plan_with_nthreads (_nthreads);
270  //Clear the current plans
271  instance->rplan = instance->plan[0] = instance->plan[1] = 0;
272  }
273  }
274 
275  static int threads ()
276  {
277  return instance_ok () ? instance->nthreads : 0;
278  }
279 #endif
280 
281 private:
282 
283  // No copying!
284 
286 
288 
290 
291  static void cleanup_instance (void) { delete instance; instance = 0; }
292 
293  fftwf_plan
294  do_create_plan (int dir, const int rank, const dim_vector dims,
295  octave_idx_type howmany, octave_idx_type stride,
296  octave_idx_type dist, const FloatComplex *in,
297  FloatComplex *out);
298 
299  fftwf_plan
300  do_create_plan (const int rank, const dim_vector dims,
301  octave_idx_type howmany, octave_idx_type stride,
302  octave_idx_type dist, const float *in, FloatComplex *out);
303 
304  FftwMethod do_method (void);
305 
306  FftwMethod do_method (FftwMethod _meth);
307 
309 
310  // FIXME: perhaps this should be split into two classes?
311 
312  // Plan for fft and ifft of complex values
313  fftwf_plan plan[2];
314 
315  // dist
317 
318  // stride
320 
321  // rank
322  int r[2];
323 
324  // howmany
326 
327  // dims
329 
330  bool simd_align[2];
331  bool inplace[2];
332 
333  // Plan for fft of real values
334  fftwf_plan rplan;
335 
336  // dist
338 
339  // stride
341 
342  // rank
343  int rr;
344 
345  // howmany
347 
348  // dims
350 
352 
353 #if defined (HAVE_FFTW3F_THREADS)
354  //number of threads when compiled with Multi-threading support
355  int nthreads;
356 #endif
357 };
358 
359 class
360 OCTAVE_API
362 {
363 public:
364 
365  static int fft (const double *in, Complex *out, size_t npts,
366  size_t nsamples = 1, octave_idx_type stride = 1,
367  octave_idx_type dist = -1);
368  static int fft (const Complex *in, Complex *out, size_t npts,
369  size_t nsamples = 1, octave_idx_type stride = 1,
370  octave_idx_type dist = -1);
371  static int ifft (const Complex *in, Complex *out, size_t npts,
372  size_t nsamples = 1, octave_idx_type stride = 1,
373  octave_idx_type dist = -1);
374 
375  static int fftNd (const double*, Complex*, const int, const dim_vector &);
376  static int fftNd (const Complex*, Complex*, const int,
377  const dim_vector &);
378  static int ifftNd (const Complex*, Complex*, const int,
379  const dim_vector &);
380 
381  static int fft (const float *in, FloatComplex *out, size_t npts,
382  size_t nsamples = 1, octave_idx_type stride = 1,
383  octave_idx_type dist = -1);
384  static int fft (const FloatComplex *in, FloatComplex *out, size_t npts,
385  size_t nsamples = 1, octave_idx_type stride = 1,
386  octave_idx_type dist = -1);
387  static int ifft (const FloatComplex *in, FloatComplex *out, size_t npts,
388  size_t nsamples = 1, octave_idx_type stride = 1,
389  octave_idx_type dist = -1);
390 
391  static int fftNd (const float*, FloatComplex*, const int, const dim_vector &);
392  static int fftNd (const FloatComplex*, FloatComplex*, const int,
393  const dim_vector &);
394  static int ifftNd (const FloatComplex*, FloatComplex*, const int,
395  const dim_vector &);
396 
397 private:
398  octave_fftw (void);
399  octave_fftw (const octave_fftw&);
400  octave_fftw& operator = (const octave_fftw&);
401 };
402 
403 #endif
404 
405 #endif