LAL  7.5.0.1-8083555
fftwf_wisdom.c
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2011 Joshua L. Willis
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 this program; if not, write to the Free Software
16  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
17  *
18  */
19 
20 /**
21  * \file
22  * \ingroup lal_fft
23  * \author Josh Willis
24  *
25  * \brief Utility to create single-precision FFTW wisdom files
26  *
27  * This utility is designed to serve as a replacement for the utility fftwf-wisdom that
28  * comes with FFTW. It will generate only the kind of plans that LAL supports, as well
29  * as a detailed, human readable description of those plans that is sent to stdout. This
30  * description can then be saved so that the user knows what transforms are included in
31  * the wisdom file (which is not itself human-readable).
32  *
33  * The transforms to be planned must be specified in an input file, with one transform
34  * per line. The format is:
35  *
36  * <type><direc><size>
37  *
38  * where:
39  *
40  * - <type> may be 'c' for a complex transform, or 'r' for a real transform
41  * - <direc> may be 'f' for a forward transform, or either 'r' or 'b' for a reverse or
42  * backward transform. The latter two are equivalent; the dual notation is
43  * designed to support both the 'backward' terminology of FFTW and the
44  * 'reverse' terminology of LAL.
45  * - <size> is the size of the desired transform.
46  *
47  * Both the <type> and <direc> specifiers are insensitive to case.
48  *
49  * The invocation must also specify an output file, where the wisdom itself will be
50  * written. You may redirect stdout to save the human-readable description. This should
51  * always be APPENDED to the human readable description of any existing wisdom files
52  * (system or otherwise) that were read in as part of the execution; there is no way
53  * for lal_fftwf_wisdom to print out the human readable description of any existing
54  * wisdom read in.
55  *
56  * So, as an example, suppose that there is an existing system-wide file in
57  * /etc/fftw/wisdomf, and that its human readable description is
58  * /etc/fftw/wisdomf_description. Then if there is presently no plan for a complex,
59  * reverse transform of size 1048576 (or 2^20), then put the line:
60  *
61  * cb1048576
62  *
63  * into the file input_plans (say) and run:
64  *
65  * cp /etc/fftw/wisdomf_description new_description
66  *
67  * lal_fftwf_wisdom -i input_plans -o new_wisdom -l 3 >> new_description
68  *
69  * When this has finished you can (as root) move new_wisdom to /etc/fftw/wisdomf and
70  * new_description to /etc/fftw/wisdomf_description. The same results could also be
71  * achieved by specifying any of cr1048576, CB1048576, etc, in input_plans.
72  *
73  * Aside from the options specifying input and output files (which are mandatory) there
74  * are three other optional arguments. They are:
75  *
76  * - -n or --no-system-wisdom This disables reading in /etc/fftw/wisdomf. Use this only
77  * if you do NOT want the wisdom there to be incorporated into the wisdom file you generate
78  * (for instance, if it is outdated because of hardware or FFTW library change). In such a
79  * case you should also not append to the existing human readable description, but replace it.
80  * - -w <FILE> or --wisdom=<FILE> Read in existing wisdom from <FILE>. Program exits if this
81  * option is given but the corresponding file is not readable. As with the system wisdom, you
82  * should append to a human readable description of the wisdom in this file; what is there
83  * already will not be descirbed in the output of lal_fftwf_wisdom.
84  * - -l <int> or --measurelvl=<int> The planning measure level used in creating the plans.
85  * Defaults to 3 (exhaustive) if not given. Levels 0 or 1 (estimate and measure, respectively)
86  * are possible, but not appropriate for a system-wide wisdom file, as plans of those levels
87  * can be efficiently generated on the fly in application code. Level 3 (exhaustive) is the
88  * highest level and would normally be expected to give the best performance, though for large
89  * transform sizes it can take hours or even days to create the plans. The human readable
90  * description of the generated wisdom will note the measure level at which it was created.
91  *
92  * \note All of the human readable descriptions will also indicate that the wisdom was created
93  * with FFTW_UNALIGNED. Users unfamiliar with this may safely ignore it; all LAL plans
94  * are presently created with this flag. The notation is there in case that should change
95  * at some point in the future, so that the record of the created wisdom preserves the
96  * distinction.
97  */
98 
99 #include <stdlib.h>
100 #include <stdio.h>
101 #include <fftw3.h>
102 #include <limits.h> /* For LINE_MAX */
103 
104 #include <lal/FileIO.h>
105 #include <lal/LALgetopt.h>
106 #include <lal/StringInput.h>
107 
108 /* prototypes */
109 void print_help(void);
110 int plan_problem(char type, char direc, UINT4 transform_size, int measurelvl);
111 
112 /** Print basic usage information about the program and exit */
113 void print_help(void)
114 {
115  fprintf(stderr,"\nlal_fftwf_wisdom -i INPUT_FILE -o OUTPUT_FILE [OPTIONS]:\n");
116  fprintf(stderr,"This program creates a wisdom file based on an input file specifying transform\n");
117  fprintf(stderr,"kinds and sizes. This is much like a simplified version of the fftwf-wisdom\n");
118  fprintf(stderr,"utility provided with FFTW. The transforms to be planned for must be specified\n");
119  fprintf(stderr,"through the input file, one transform per line, using the format:\n");
120  fprintf(stderr," <type><direction><size>\n");
121  fprintf(stderr,"where:\n");
122  fprintf(stderr," <type> is 'r' or 'c' for real or complex\n");
123  fprintf(stderr," <direction> is 'f' or either 'b' or 'r' for forward or backward/reverse\n");
124  fprintf(stderr," <size> is the size of the transform\n");
125  fprintf(stderr,"Backward and reverse are synonymous, and <type> and <direction> are not case-\n");
126  fprintf(stderr,"sensitive. For further details and warnings, see the doxygen documentation.\n");
127  fprintf(stderr,"\n");
128  fprintf(stderr,"Options and their behavior are:\n");
129  fprintf(stderr,"\n");
130  fprintf(stderr," --measurelvl=, -l <int> The measurelvl argument to plan creation calls. The\n");
131  fprintf(stderr," lowest level of zero corresponds to no planning,\n");
132  fprintf(stderr," whereas the highest level of 3 could take days, de-\n");
133  fprintf(stderr," pending on the transform size. Defaults to 3.\n");
134  fprintf(stderr," --input=, -i <file> Input file containing transforms to plan. Mandatory.\n");
135  fprintf(stderr," --output=, -o <file> File in which to save generated wisdom. Mandatory.\n");
136  fprintf(stderr," --wisdom=, -w <file> Existing wisdom file to import. Program exits if this\n");
137  fprintf(stderr," option is given but the file is unreadable. Optional.\n");
138  fprintf(stderr," --no-system-wisdom, -n By default, system wisdom file /etc/fftw/wisdomf is \n");
139  fprintf(stderr," imported. If this option is specified, that import is\n");
140  fprintf(stderr," disabled. If the system file is not present or not \n");
141  fprintf(stderr," readable, a warning is printed but execution continues.\n");
142  fprintf(stderr," --help, -h Print this help message and exit.\n");
143  fprintf(stderr,"In addition to writing the accumulated wisdom to the specified output file, the\n");
144  fprintf(stderr,"program also prints to stdout a human-readable description of the wisdom success-\n");
145  fprintf(stderr,"fully created during this invocation. That description should be appended to the\n");
146  fprintf(stderr,"corresponding description for any wisdom files (including system) read in.\n");
147  exit(EXIT_SUCCESS);
148 }
149 
150 
151 /**
152  * Function used only internally, to create an FFTW plan for a specified problem (thereby adding to wisdom)
153  */
154 int plan_problem(char type, /**< 'r' for real or 'c' for complex transform */
155  char direc, /**< 'f' for forward or 'b'/'r' for backward/reverse transform */
156  UINT4 transform_size, /**< Size of transform to plan */
157  int measurelvl) /**< Level of patience in planning (0 least, 3 most) */
158 {
159  fftwf_plan genericPlan;
160  void *indata, *outdata;
161  int fwdflag, planning_flags;
162 
163  fwdflag = ( (direc=='f') || (direc=='F') );
164 
165  /* We call FFTW routines directly, rather than through LAL, so that if XLAL planning routines
166  are changed to always read in wisdom, we can still toggle that behavior through the command line.
167  In case we ever allow for aligned memory, we allocate everything with fftwf_malloc().
168  */
169 
170  /* If we ever allow for aligned memory, this will have to toggle depending on input: */
171  planning_flags = FFTW_UNALIGNED;
172 
173  switch(measurelvl)
174  {
175  case 0:
176  planning_flags |= FFTW_ESTIMATE;
177  break;
178  default:
179  case 3:
180  planning_flags |= FFTW_EXHAUSTIVE;
181  /* Fall through: */
182 #if __GNUC__ >= 7 && !defined __INTEL_COMPILER
183  __attribute__ ((fallthrough));
184 #endif
185  case 2:
186  planning_flags |= FFTW_PATIENT;
187  /* Fall through */
188 #if __GNUC__ >= 7 && !defined __INTEL_COMPILER
189  __attribute__ ((fallthrough));
190 #endif
191  case 1:
192  planning_flags |= FFTW_MEASURE;
193  break;
194  }
195 
196  /* Ugly, but makes us 'locale' independent */
197 
198  if ( (type=='r') || (type=='R') )
199  {
200 
201  indata = (float *) fftwf_malloc(transform_size*sizeof(float));
202  outdata = (float *) fftwf_malloc(transform_size*sizeof(float));
203 
204  if ( (!indata) || (!outdata) )
205  {
206  if (indata) fftwf_free(indata);
207  if (outdata) fftwf_free(outdata);
208  return 1;
209  }
210 
211  genericPlan = fftwf_plan_r2r_1d(transform_size,indata,outdata,
212  (fwdflag ? FFTW_R2HC : FFTW_HC2R),
213  planning_flags);
214  if (!genericPlan)
215  {
216  fftwf_free(indata);
217  fftwf_free(outdata);
218  return 1;
219  }
220  else
221  {
222  fftwf_free(indata);
223  fftwf_free(outdata);
224  fftwf_destroy_plan(genericPlan);
225  return 0;
226  }
227  }
228  else
229  { /* type == 'c' */
230 
231  indata = (fftwf_complex *) fftwf_malloc(transform_size*sizeof(fftwf_complex));
232  outdata = (fftwf_complex *) fftwf_malloc(transform_size*sizeof(fftwf_complex));
233 
234  if ( (!indata) || (!outdata) )
235  {
236  if (indata) fftwf_free(indata);
237  if (outdata) fftwf_free(outdata);
238  return 1;
239  }
240 
241  genericPlan = fftwf_plan_dft_1d(transform_size,indata,outdata,
242  (fwdflag ? FFTW_FORWARD : FFTW_BACKWARD),
243  planning_flags);
244 
245  if (!genericPlan)
246  {
247  fftwf_free(indata);
248  fftwf_free(outdata);
249  return 1;
250  }
251  else
252  {
253  fftwf_free(indata);
254  fftwf_free(outdata);
255  fftwf_destroy_plan(genericPlan);
256  return 0;
257  }
258  }
259 }
260 
261 /**
262  * Main function
263  *
264  * Reads command line specifying input and output files, and optionally wisdom file,
265  * measure level, and flag preventing import of system wide wisdom, and then creates
266  * plans for each problem specified in the input file. Accumulated wisdom is written
267  * to the output file, and a human readable description of the plans successfully
268  * created is written to stdout. Any warnings or errors are written to stderr.
269  */
270 int main(int argc, char **argv)
271 {
272  static int measurelvl=3;
273  static int nosys=0;
274  UINT4 transform_size;
275  char input_line[LINE_MAX];
276  char type;
277  char direc;
278  FILE *infp=NULL, *outfp=NULL, *wisfp=NULL;
279  int optindex, optreturn, retval;
280 
281  static struct LALoption long_options[] =
282  {
283  /* Options setting flags */
284  {"no-system-wisdom",no_argument,&nosys,1},
285  /* Options specifying input/output */
286  {"input",required_argument,NULL,'i'},
287  {"output",required_argument,NULL,'o'},
288  {"wisdom",required_argument,NULL,'w'},
289  {"measurelvl",required_argument,NULL,'l'},
290  {"help",no_argument,NULL,'h'},
291  {0,0,0,0}
292  };
293 
294  while ( (optreturn = LALgetopt_long(argc,argv,"ni:o:w:l:h",long_options,&optindex)) != -1)
295  {
296  switch(optreturn)
297  {
298  case 0:
299  break; /* Everything done in setting flag */
300  case 'n':
301  nosys=1;
302  break;
303  case 'i':
304  infp = LALFopen(LALoptarg,"r+");
305  if (!infp)
306  {
307  fprintf(stderr,"Error: Could not open input file %s\n",LALoptarg);
308  if (outfp) LALFclose(outfp);
309  exit(EXIT_FAILURE);
310  }
311  break;
312  case 'o':
313  outfp = LALFopen(LALoptarg,"w+");
314  if (!outfp)
315  {
316  fprintf(stderr,"Error: Could not open output file %s\n",LALoptarg);
317  if (infp) LALFclose(infp);
318  exit(EXIT_FAILURE);
319  }
320  break;
321  case 'w':
322  wisfp = LALFopen(LALoptarg,"r+");
323  if (!wisfp)
324  {
325  fprintf(stderr,"Error: Could not open input wisdom file %s for reading\n",LALoptarg);
326  if (infp) LALFclose(infp);
327  if (outfp) LALFclose(outfp);
328  exit(EXIT_FAILURE);
329  }
330  else
331  {
332  retval = fftwf_import_wisdom_from_file(wisfp);
333  if (!retval)
334  {
335  /* Retval is zero if UNsuccessful */
336  fprintf(stderr,"Error: Could not read wisdom from input wisdom file %s\n",LALoptarg);
337  if (infp) LALFclose(infp);
338  if (outfp) LALFclose(outfp);
339  LALFclose(wisfp);
340  exit(EXIT_FAILURE);
341  }
342  LALFclose(wisfp);
343  }
344  fprintf(stderr,"Read in existing wisdom from file %s\n",LALoptarg);
345  break;
346  case 'l':
347  if ( sscanf(LALoptarg,"%d",&measurelvl) != 1)
348  {
349  fprintf(stderr,"Error: invalid measure level %s.\n",LALoptarg);
350  if (infp) LALFclose(infp);
351  if (outfp) LALFclose(outfp);
352  exit(EXIT_FAILURE);
353  }
354  if ( (measurelvl<0) || (measurelvl>3) )
355  {
356  fprintf(stderr,"Error: invalid measure level %d.\n",measurelvl);
357  if (infp) LALFclose(infp);
358  if (outfp) LALFclose(outfp);
359  }
360  break;
361  case 'h': /* Fall through */
362  case '?':
363  print_help();
364  break;
365  default:
366  exit(EXIT_FAILURE);
367  } /* switch(optreturn) */
368  } /* while(optreturn != -1) */
369 
370  /* Check to make sure mandatory options were given */
371 
372  if (!infp)
373  {
374  fprintf(stderr,"Error: You must specify an input file with -i <FILE> or --input=<FILE>\n");
375  if (outfp) LALFclose(outfp);
376  exit(EXIT_FAILURE);
377  }
378 
379  if (!outfp)
380  {
381  fprintf(stderr,"Error: You must specify an output file with -o <FILE> or --output=<FILE>\n");
382  if (infp) LALFclose(infp);
383  exit(EXIT_FAILURE);
384  }
385 
386  /* Only after processing all options do we know if we should read in system wisdom file */
387 
388  if (!nosys)
389  {
390  retval = fftwf_import_system_wisdom();
391  if (!retval)
392  {
393  /* Retval is zero if UNsuccessful */
394  fprintf(stderr,"Warning: Could not import system wisdom file /etc/fftw/wisdomf\n");
395  }
396  }
397  else
398  {
399  fprintf(stderr,"Skipped import of system wisdom file /etc/fftw/wisdomf\n");
400  }
401 
402 
403  /* Process the input file */
404 
405  while ( (fgets(input_line,LINE_MAX,infp) != NULL) )
406  {
407  if (sscanf(input_line,"%c%c%" LAL_UINT4_FORMAT, &type, &direc, &transform_size) == 3)
408  {
409  /* Yes, it's ugly, but we don't have to worry about locales: */
410  if ( !( (type=='r') || (type=='R') || (type=='c') || (type=='C') ) )
411  {
412  fprintf(stderr,"Error: Invalid type specifier %c; must be 'r' (real) or 'c' (complex). ",type);
413  fprintf(stderr,"Problem %c%c%" LAL_UINT4_FORMAT " will be skipped!\n", type, direc, transform_size);
414  }
415  else if ( !( (direc=='f') || (direc=='b') || (direc=='r') || (direc=='F') || (direc=='B') || (direc=='R') ) )
416  {
417  fprintf(stderr,"Error: Invalid direction specifier %c; must be 'f' (forward) or 'b'/'r' (backward/reverse). ",
418  direc);
419  fprintf(stderr,"Problem %c%c%" LAL_UINT4_FORMAT " will be skipped!\n",type,direc,transform_size);
420  }
421  else
422  {
423  retval = plan_problem(type,direc,transform_size,measurelvl);
424  if (retval)
425  {
426  fprintf(stderr,"Unable to create plan %c%c%" LAL_UINT4_FORMAT "; skipping!\n",
427  type,direc,transform_size);
428  }
429  else
430  {
431  fprintf(stdout,"Created single-precision %s %s plan, size %" LAL_UINT4_FORMAT
432  " with measure level %d and FFTW_UNALIGNED\n",
433  ( (type=='r') || (type=='R') ) ? "REAL4" : "COMPLEX8",
434  ( (direc=='f') || (direc=='F') ) ? "forward" : "reverse",
435  transform_size, measurelvl);
436  }
437  }
438  }
439  else
440  {
441  fprintf(stderr,"Error: Invalid problem specifier. Problem: %s will be skipped\n",input_line);
442  }
443  }
444 
445  fftwf_export_wisdom_to_file(outfp);
446  LALFclose(infp);
447  LALFclose(outfp);
448 
449  exit(EXIT_SUCCESS);
450 }
int LALgetopt_long(int argc, char *const *argv, const char *options, const struct LALoption *long_options, int *opt_index)
Definition: LALgetopt.c:178
char * LALoptarg
Definition: LALgetopt.c:64
#define no_argument
Definition: LALgetopt.h:100
#define required_argument
Definition: LALgetopt.h:101
#define fprintf
int main(int argc, char **argv)
Main function.
Definition: fftwf_wisdom.c:270
void print_help(void)
Print basic usage information about the program and exit.
Definition: fftwf_wisdom.c:113
int plan_problem(char type, char direc, UINT4 transform_size, int measurelvl)
Function used only internally, to create an FFTW plan for a specified problem (thereby adding to wisd...
Definition: fftwf_wisdom.c:154
#define __attribute__(x)
Definition: getdate.c:146
uint32_t UINT4
Four-byte unsigned integer.
#define LALFclose
Definition: LALStdio.h:51
#define LALFopen
Definition: LALStdio.h:50
#define LAL_UINT4_FORMAT
Definition: LALStdio.h:130