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