GNU Octave  4.0.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
fftw.cc
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 2006-2015 David Bateman
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 #ifdef HAVE_CONFIG_H
24 #include <config.h>
25 #endif
26 
27 #include <algorithm>
28 
29 #include "oct-fftw.h"
30 
31 #include "defun-dld.h"
32 #include "error.h"
33 #include "ov.h"
34 
35 #include "gripes.h"
36 
37 DEFUN_DLD (fftw, args, ,
38  "-*- texinfo -*-\n\
39 @deftypefn {Loadable Function} {@var{method} =} fftw (\"planner\")\n\
40 @deftypefnx {Loadable Function} {} fftw (\"planner\", @var{method})\n\
41 @deftypefnx {Loadable Function} {@var{wisdom} =} fftw (\"dwisdom\")\n\
42 @deftypefnx {Loadable Function} {} fftw (\"dwisdom\", @var{wisdom})\n\
43 @deftypefnx {Loadable Function} {} fftw (\"threads\", @var{nthreads})\n\
44 @deftypefnx {Loadable Function} {@var{nthreads} =} fftw (\"threads\")\n\
45 \n\
46 Manage @sc{fftw} wisdom data.\n\
47 \n\
48 Wisdom data can be used to significantly accelerate the calculation of the\n\
49 FFTs, but implies an initial cost in its calculation. When the @sc{fftw}\n\
50 libraries are initialized, they read a system wide wisdom file (typically in\n\
51 @file{/etc/fftw/wisdom}), allowing wisdom to be shared between applications\n\
52 other than Octave. Alternatively, the @code{fftw} function can be used to\n\
53 import wisdom. For example,\n\
54 \n\
55 @example\n\
56 @var{wisdom} = fftw (\"dwisdom\")\n\
57 @end example\n\
58 \n\
59 @noindent\n\
60 will save the existing wisdom used by Octave to the string @var{wisdom}.\n\
61 This string can then be saved to a file and restored using the @code{save}\n\
62 and @code{load} commands respectively. This existing wisdom can be\n\
63 re-imported as follows\n\
64 \n\
65 @example\n\
66 fftw (\"dwisdom\", @var{wisdom})\n\
67 @end example\n\
68 \n\
69 If @var{wisdom} is an empty string, then the wisdom used is cleared.\n\
70 \n\
71 During the calculation of Fourier transforms further wisdom is generated.\n\
72 The fashion in which this wisdom is generated is also controlled by\n\
73 the @code{fftw} function. There are five different manners in which the\n\
74 wisdom can be treated:\n\
75 \n\
76 @table @asis\n\
77 @item @qcode{\"estimate\"}\n\
78 Specifies that no run-time measurement of the optimal means of\n\
79 calculating a particular is performed, and a simple heuristic is used\n\
80 to pick a (probably sub-optimal) plan. The advantage of this method is\n\
81 that there is little or no overhead in the generation of the plan, which\n\
82 is appropriate for a Fourier transform that will be calculated once.\n\
83 \n\
84 @item @qcode{\"measure\"}\n\
85 In this case a range of algorithms to perform the transform is considered\n\
86 and the best is selected based on their execution time.\n\
87 \n\
88 @item @qcode{\"patient\"}\n\
89 Similar to @qcode{\"measure\"}, but a wider range of algorithms is\n\
90 considered.\n\
91 \n\
92 @item @qcode{\"exhaustive\"}\n\
93 Like @qcode{\"measure\"}, but all possible algorithms that may be used to\n\
94 treat the transform are considered.\n\
95 \n\
96 @item @qcode{\"hybrid\"}\n\
97 As run-time measurement of the algorithm can be expensive, this is a\n\
98 compromise where @qcode{\"measure\"} is used for transforms up to the size\n\
99 of 8192 and beyond that the @qcode{\"estimate\"} method is used.\n\
100 @end table\n\
101 \n\
102 The default method is @qcode{\"estimate\"}. The current method can\n\
103 be queried with\n\
104 \n\
105 @example\n\
106 @var{method} = fftw (\"planner\")\n\
107 @end example\n\
108 \n\
109 @noindent\n\
110 or set by using\n\
111 \n\
112 @example\n\
113 fftw (\"planner\", @var{method})\n\
114 @end example\n\
115 \n\
116 Note that calculated wisdom will be lost when restarting Octave. However,\n\
117 the wisdom data can be reloaded if it is saved to a file as described\n\
118 above. Saved wisdom files should not be used on different platforms since\n\
119 they will not be efficient and the point of calculating the wisdom is lost.\n\
120 \n\
121 The number of threads used for computing the plans and executing the\n\
122 transforms can be set with\n\
123 \n\
124 @example\n\
125 fftw (\"threads\", @var{NTHREADS})\n\
126 @end example\n\
127 \n\
128 Note that octave must be compiled with multi-threaded @sc{fftw} support for\n\
129 this feature. The number of processors available to the current process is\n\
130 used per default.\n\
131 \n\
132 @seealso{fft, ifft, fft2, ifft2, fftn, ifftn}\n\
133 @end deftypefn")
134 {
135  octave_value retval;
136 
137  int nargin = args.length ();
138 
139  if (nargin < 1 || nargin > 2)
140  {
141  print_usage ();
142  return retval;
143  }
144 
145 #if defined (HAVE_FFTW)
146  if (args(0).is_string ())
147  {
148  std::string arg0 = args(0).string_value ();
149 
150  if (arg0 == "planner")
151  {
152  if (nargin == 2) //planner setter
153  {
154  if (args(1).is_string ())
155  {
156  // Use STL function to convert to lower case
157  std::transform (arg0.begin (), arg0.end (), arg0.begin (),
158  tolower);
159  std::string arg1 = args(1).string_value ();
160  if (!error_state)
161  {
162  std::transform (arg1.begin (), arg1.end (),
163  arg1.begin (), tolower);
168 
169  if (arg1 == "estimate")
170  {
173  }
174  else if (arg1 == "measure")
175  {
178  }
179  else if (arg1 == "patient")
180  {
183  }
184  else if (arg1 == "exhaustive")
185  {
188  }
189  else if (arg1 == "hybrid")
190  {
193  }
194  else
195  error ("fftw: unrecognized planner METHOD");
196 
197  if (!error_state)
198  {
199  meth = octave_fftw_planner::method (meth);
201 
202  if (meth == octave_fftw_planner::MEASURE)
203  retval = octave_value ("measure");
204  else if (meth == octave_fftw_planner::PATIENT)
205  retval = octave_value ("patient");
206  else if (meth == octave_fftw_planner::EXHAUSTIVE)
207  retval = octave_value ("exhaustive");
208  else if (meth == octave_fftw_planner::HYBRID)
209  retval = octave_value ("hybrid");
210  else
211  retval = octave_value ("estimate");
212  }
213  }
214  }
215  else
216  error ("fftw: planner expects a string value as METHOD");
217  }
218  else //planner getter
219  {
222 
223  if (meth == octave_fftw_planner::MEASURE)
224  retval = octave_value ("measure");
225  else if (meth == octave_fftw_planner::PATIENT)
226  retval = octave_value ("patient");
227  else if (meth == octave_fftw_planner::EXHAUSTIVE)
228  retval = octave_value ("exhaustive");
229  else if (meth == octave_fftw_planner::HYBRID)
230  retval = octave_value ("hybrid");
231  else
232  retval = octave_value ("estimate");
233  }
234  }
235  else if (arg0 == "dwisdom")
236  {
237  if (nargin == 2) //dwisdom setter
238  {
239  if (args(1).is_string ())
240  {
241  // Use STL function to convert to lower case
242  std::transform (arg0.begin (), arg0.end (), arg0.begin (),
243  tolower);
244  std::string arg1 = args(1).string_value ();
245  if (!error_state)
246  {
247  char *str = fftw_export_wisdom_to_string ();
248 
249  if (arg1.length () < 1)
250  fftw_forget_wisdom ();
251  else if (! fftw_import_wisdom_from_string (arg1.c_str ()))
252  error ("fftw: could not import supplied WISDOM");
253 
254  if (!error_state)
255  retval = octave_value (std::string (str));
256 
257  free (str);
258  }
259  }
260  }
261  else //dwisdom getter
262  {
263  char *str = fftw_export_wisdom_to_string ();
264  retval = octave_value (std::string (str));
265  free (str);
266  }
267  }
268  else if (arg0 == "swisdom")
269  {
270  //swisdom uses fftwf_ functions (float), dwisdom fftw_ (real)
271  if (nargin == 2) //swisdom setter
272  {
273  if (args(1).is_string ())
274  {
275  // Use STL function to convert to lower case
276  std::transform (arg0.begin (), arg0.end (), arg0.begin (),
277  tolower);
278  std::string arg1 = args(1).string_value ();
279  if (!error_state)
280  {
281  char *str = fftwf_export_wisdom_to_string ();
282 
283  if (arg1.length () < 1)
284  fftwf_forget_wisdom ();
285  else if (! fftwf_import_wisdom_from_string (arg1.c_str ()))
286  error ("fftw: could not import supplied WISDOM");
287 
288  if (!error_state)
289  retval = octave_value (std::string (str));
290 
291  free (str);
292  }
293  }
294  }
295  else //swisdom getter
296  {
297  char *str = fftwf_export_wisdom_to_string ();
298  retval = octave_value (std::string (str));
299  free (str);
300  }
301  }
302  else if (arg0 == "threads")
303  {
304  if (nargin == 2) //threads setter
305  {
306  if (args(1).is_real_scalar ())
307  {
308  int nthreads = args(1).int_value();
309  if (nthreads >= 1)
310  {
311 #if defined (HAVE_FFTW3_THREADS)
312  octave_fftw_planner::threads (nthreads);
313 #else
314  gripe_disabled_feature ("fftw", "multithreaded FFTW");
315 #endif
316 #if defined (HAVE_FFTW3F_THREADS)
318 #else
319  gripe_disabled_feature ("fftw", "multithreaded FFTW");
320 #endif
321  }
322  else
323  error ("fftw: number of threads must be >=1");
324  }
325  else
326  error ("fftw: setting threads needs one integer argument");
327  }
328  else //threads getter
329 #if defined (HAVE_FFTW3_THREADS)
331 #else
332  retval = 1;
333 #endif
334  }
335  else
336  error ("fftw: unrecognized argument");
337  }
338  else
339  error ("fftw: unrecognized argument");
340 #else
341 
342  gripe_disabled_feature ("fftw", "the FFTW3 planner");
343 
344 #endif
345 
346  return retval;
347 }
348 
349 /*
350 %!testif HAVE_FFTW
351 %! def_method = fftw ("planner");
352 %! unwind_protect
353 %! method = "estimate";
354 %! fftw ("planner", method);
355 %! assert (fftw ("planner"), method);
356 %! method = "measure";
357 %! fftw ("planner", method);
358 %! assert (fftw ("planner"), method);
359 %! method = "patient";
360 %! fftw ("planner", method);
361 %! assert (fftw ("planner"), method);
362 %! method = "exhaustive";
363 %! fftw ("planner", method);
364 %! assert (fftw ("planner"), method);
365 %! method = "hybrid";
366 %! fftw ("planner", method);
367 %! assert (fftw ("planner"), method);
368 %! unwind_protect_cleanup
369 %! fftw ("planner", def_method);
370 %! end_unwind_protect
371 
372 %!testif HAVE_FFTW
373 %! def_dwisdom = fftw ("dwisdom");
374 %! def_swisdom = fftw ("swisdom");
375 %! unwind_protect
376 %! wisdom = fftw ("dwisdom");
377 %! assert (ischar (wisdom));
378 %! fftw ("dwisdom", wisdom);
379 %! assert (fftw ("dwisdom"), wisdom);
380 %! wisdom = fftw ("swisdom");
381 %! assert (ischar (wisdom));
382 %! fftw ("swisdom", wisdom);
383 %! assert (fftw ("swisdom"), wisdom);
384 %! unwind_protect_cleanup
385 %! fftw ("dwisdom", def_dwisdom);
386 %! fftw ("swisdom", def_swisdom);
387 %! end_unwind_protect
388 
389 %!testif HAVE_FFTW3_THREADS
390 %! n = fftw ("threads");
391 %! unwind_protect
392 %! fftw ("threads", 3);
393 %! assert (fftw ("threads"), 3);
394 %! unwind_protect_cleanup
395 %! fftw ("threads", n);
396 %! end_unwind_protect
397 
398 %!error <Invalid call to fftw> fftw ();
399 %!error <Invalid call to fftw> fftw ("planner", "estimate", "measure");
400 %!error fftw (3);
401 %!error fftw ("invalid");
402 %!error fftw ("planner", "invalid");
403 %!error fftw ("planner", 2);
404 %!error fftw ("dwisdom", "invalid");
405 %!error fftw ("swisdom", "invalid");
406 %!error fftw ("threads", "invalid");
407 %!error fftw ("threads", -3);
408  */
OCTINTERP_API void print_usage(void)
Definition: defun.cc:51
void error(const char *fmt,...)
Definition: error.cc:476
static FftwMethod method(void)
Definition: oct-fftw.h:88
static FftwMethod method(void)
Definition: oct-fftw.h:249
std::string string_value(bool force=false) const
Definition: ov.h:897
void gripe_disabled_feature(const std::string &func, const std::string &feature, const std::string &pkg)
Definition: gripes.cc:242
int error_state
Definition: error.cc:101
octave_idx_type length(void) const
Definition: ov.cc:1525
void free(void *)
static int threads()
Definition: oct-fftw.h:114
ColumnVector transform(const Matrix &m, double x, double y, double z)
Definition: graphics.cc:5259
#define DEFUN_DLD(name, args_name, nargout_name, doc)
Definition: defun-dld.h:59
static int threads()
Definition: oct-fftw.h:275
return octave_value(v1.char_array_value().concat(v2.char_array_value(), ra_idx),((a1.is_sq_string()||a2.is_sq_string())? '\'': '"'))