NFFT  3.4.1
infft.h
1 /*
2  * Copyright (c) 2002, 2017 Jens Keiner, Stefan Kunis, Daniel Potts
3  *
4  * This program is free software; you can redistribute it and/or modify it under
5  * the terms of the GNU General Public License as published by the Free Software
6  * Foundation; either version 2 of the License, or (at your option) any later
7  * version.
8  *
9  * This program is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
11  * FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
12  * details.
13  *
14  * You should have received a copy of the GNU General Public License along with
15  * this program; if not, write to the Free Software Foundation, Inc., 51
16  * Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
17  */
18 
22 #ifndef __INFFT_H__
23 #define __INFFT_H__
24 
25 #include "config.h"
26 
27 #include <math.h>
28 #include <float.h>
29 #ifdef HAVE_COMPLEX_H
30 #include <complex.h>
31 #endif
32 #include <stdio.h>
33 #include <string.h>
34 
35 #include <stdlib.h> /* size_t */
36 #include <stdarg.h> /* va_list */
37 #include <stddef.h> /* ptrdiff_t */
38 
39 #if HAVE_SYS_TYPES_H
40 #include <sys/types.h>
41 #endif
42 
43 #if HAVE_STDINT_H
44 #include <stdint.h> /* uintptr_t, maybe */
45 #endif
46 
47 #if HAVE_INTTYPES_H
48 #include <inttypes.h> /* uintptr_t, maybe */
49 #endif
50 
51 #include <fftw3.h>
52 
53 #include "ticks.h"
54 
66 /* Determine precision and name-mangling scheme. */
67 #define CONCAT(prefix, name) prefix ## name
68 #if defined(NFFT_SINGLE)
69 typedef float R;
70 typedef float _Complex C;
71 #define Y(name) CONCAT(nfftf_,name)
72 #define FFTW(name) CONCAT(fftwf_,name)
73 #define NFFT(name) CONCAT(nfftf_,name)
74 #define NFCT(name) CONCAT(nfctf_,name)
75 #define NFST(name) CONCAT(nfstf_,name)
76 #define NFSFT(name) CONCAT(nfsftf_,name)
77 #define SOLVER(name) CONCAT(solverf_,name)
78 #elif defined(NFFT_LDOUBLE)
79 typedef long double R;
80 typedef long double _Complex C;
81 #define Y(name) CONCAT(nfftl_,name)
82 #define FFTW(name) CONCAT(fftwl_,name)
83 #define NFFT(name) CONCAT(nfftl_,name)
84 #define NFCT(name) CONCAT(nfctl_,name)
85 #define NFST(name) CONCAT(nfstl_,name)
86 #define NFSFT(name) CONCAT(nfsftl_,name)
87 #define SOLVER(name) CONCAT(solverl_,name)
88 #else
89 typedef double R;
90 typedef double _Complex C;
91 #define Y(name) CONCAT(nfft_,name)
92 #define FFTW(name) CONCAT(fftw_,name)
93 #define NFFT(name) CONCAT(nfft_,name)
94 #define NFCT(name) CONCAT(nfct_,name)
95 #define NFST(name) CONCAT(nfst_,name)
96 #define NFSFT(name) CONCAT(nfsft_,name)
97 #define SOLVER(name) CONCAT(solver_,name)
98 #endif
99 #define X(name) Y(name)
100 
101 #define STRINGIZEx(x) #x
102 #define STRINGIZE(x) STRINGIZEx(x)
103 
104 #ifdef NFFT_LDOUBLE
105 # define K(x) ((R) x##L)
106 #else
107 # define K(x) ((R) x)
108 #endif
109 #define DK(name, value) const R name = K(value)
110 
111 #if defined __CYGWIN32__ && !defined __CYGWIN__
112  /* For backwards compatibility with Cygwin b19 and
113  earlier, we define __CYGWIN__ here, so that
114  we can rely on checking just for that macro. */
115 # define __CYGWIN__ __CYGWIN32__
116 #endif
117 
118 /* Integral type large enough to contain a stride (what ``int'' should have been
119  * in the first place) */
120 typedef ptrdiff_t INT;
121 
122 #define KPI K(3.1415926535897932384626433832795028841971693993751)
123 #define K2PI K(6.2831853071795864769252867665590057683943387987502)
124 #define K4PI K(12.5663706143591729538505735331180115367886775975004)
125 #define KE K(2.7182818284590452353602874713526624977572470937000)
126 
127 #define IF(x,a,b) ((x)?(a):(b))
128 #define MIN(a,b) (((a)<(b))?(a):(b))
129 #define MAX(a,b) (((a)>(b))?(a):(b))
130 #define ABS(x) (((x)>K(0.0))?(x):(-(x)))
131 #define SIGN(a) (((a)>=0)?1:-1)
132 #define SIGN(a) (((a)>=0)?1:-1)
133 #define SIGNF(a) IF((a)<K(0.0),K(-1.0),K(1.0))
134 
135 /* Size of array. */
136 #define SIZE(x) sizeof(x)/sizeof(x[0])
137 
139 #define CSWAP(x,y) {C* NFFT_SWAP_temp__; \
140  NFFT_SWAP_temp__=(x); (x)=(y); (y)=NFFT_SWAP_temp__;}
141 
143 #define RSWAP(x,y) {R* NFFT_SWAP_temp__; NFFT_SWAP_temp__=(x); \
144  (x)=(y); (y)=NFFT_SWAP_temp__;}
145 
146 /* macros for window functions */
147 
148 #if defined(DIRAC_DELTA)
149  #define PHI_HUT(n,k,d) K(1.0)
150  #define PHI(n,x,d) IF(FABS((x)) < K(10E-8),K(1.0),K(0.0))
151  #define WINDOW_HELP_INIT(d)
152  #define WINDOW_HELP_FINALIZE
153  #define WINDOW_HELP_ESTIMATE_m 0
154 #elif defined(GAUSSIAN)
155  #define PHI_HUT(n,k,d) ((R)EXP(-(POW(KPI*(k)/n,K(2.0))*ths->b[d])))
156  #define PHI(n,x,d) ((R)EXP(-POW((x)*((R)n),K(2.0)) / \
157  ths->b[d])/SQRT(KPI*ths->b[d]))
158  #define WINDOW_HELP_INIT \
159  { \
160  int WINDOW_idx; \
161  ths->b = (R*) Y(malloc)(ths->d*sizeof(R)); \
162  for (WINDOW_idx = 0; WINDOW_idx < ths->d; WINDOW_idx++) \
163  ths->b[WINDOW_idx]=(K(2.0)*ths->sigma[WINDOW_idx]) / \
164  (K(2.0)*ths->sigma[WINDOW_idx] - K(1.0)) * (((R)ths->m) / KPI); \
165  }
166  #define WINDOW_HELP_FINALIZE {Y(free)(ths->b);}
167 #if defined(NFFT_LDOUBLE)
168  #define WINDOW_HELP_ESTIMATE_m 17
169 #elif defined(NFFT_SINGLE)
170  #define WINDOW_HELP_ESTIMATE_m 5
171 #else
172  #define WINDOW_HELP_ESTIMATE_m 13
173 #endif
174 #elif defined(B_SPLINE)
175  #define PHI_HUT(n,k,d) ((R)(((k) == 0) ? K(1.0) / n : \
176  POW(SIN((k) * KPI / n) / ((k) * KPI / n), \
177  K(2.0) * ths->m)/n))
178  #define PHI(n,x,d) (Y(bsplines)(2*ths->m,((x)*n) + \
179  (R)ths->m) / n)
180  #define WINDOW_HELP_INIT
181  #define WINDOW_HELP_FINALIZE
182 #if defined(NFFT_LDOUBLE)
183  #define WINDOW_HELP_ESTIMATE_m 11
184 #elif defined(NFFT_SINGLE)
185  #define WINDOW_HELP_ESTIMATE_m 11
186 #else
187  #define WINDOW_HELP_ESTIMATE_m 11
188 #endif
189 #elif defined(SINC_POWER)
190  #define PHI_HUT(n,k,d) (Y(bsplines)(2 * ths->m, (K(2.0) * ths->m*(k)) / \
191  ((K(2.0) * ths->sigma[(d)] - 1) * n / \
192  ths->sigma[(d)]) + (R)ths->m))
193  #define PHI(n,x,d) ((R)(n / ths->sigma[(d)] * \
194  (K(2.0) * ths->sigma[(d)] - K(1.0))/ (K(2.0)*ths->m) * \
195  POW(Y(sinc)(KPI * n / ths->sigma[(d)] * (x) * \
196  (K(2.0) * ths->sigma[(d)] - K(1.0)) / (K(2.0)*ths->m)) , 2*ths->m) / \
197  n))
198  #define WINDOW_HELP_INIT
199  #define WINDOW_HELP_FINALIZE
200 #if defined(NFFT_LDOUBLE)
201  #define WINDOW_HELP_ESTIMATE_m 13
202 #elif defined(NFFT_SINGLE)
203  #define WINDOW_HELP_ESTIMATE_m 11
204 #else
205  #define WINDOW_HELP_ESTIMATE_m 11
206 #endif
207 #else /* Kaiser-Bessel is the default. */
208  #define PHI_HUT(n,k,d) (Y(bessel_i0)((R)(ths->m) * SQRT(ths->b[d] * ths->b[d] - (K(2.0) * KPI * (R)(k) / (R)(n)) * (K(2.0) * KPI * (R)(k) / (R)(n)))))
209  #define PHI(n,x,d) ( (((R)(ths->m) * (R)(ths->m) - (x) * (R)(n) * (x) * (R)(n)) > K(0.0)) \
210  ? SINH(ths->b[d] * SQRT((R)(ths->m) * (R)(ths->m) - (x) * (R)(n) * (x) * (R)(n))) \
211  / (KPI * SQRT((R)(ths->m) * (R)(ths->m) - (x) * (R)(n) * (x) * (R)(n))) \
212  : ((((R)(ths->m) * (R)(ths->m) - (x) * (R)(n) * (x) * (R)(n)) < K(0.0)) \
213  ? SIN(ths->b[d] * SQRT((x) * (R)(n) * (x) * (R)(n) - (R)(ths->m) * (R)(ths->m))) \
214  / (KPI * SQRT((x) * (R)(n) * (x) * (R)(n) - (R)(ths->m) * (R)(ths->m))) \
215  : ths->b[d] / KPI))
216  #define WINDOW_HELP_INIT \
217  { \
218  int WINDOW_idx; \
219  ths->b = (R*) Y(malloc)((size_t)(ths->d) * sizeof(R)); \
220  for (WINDOW_idx = 0; WINDOW_idx < ths->d; WINDOW_idx++) \
221  ths->b[WINDOW_idx] = (KPI * (K(2.0) - K(1.0) / ths->sigma[WINDOW_idx])); \
222  }
223  #define WINDOW_HELP_FINALIZE {Y(free)(ths->b);}
224  #if defined(NFFT_LDOUBLE)
225  #define WINDOW_HELP_ESTIMATE_m 9
226  #elif defined(NFFT_SINGLE)
227  #define WINDOW_HELP_ESTIMATE_m 4
228  #else
229  #define WINDOW_HELP_ESTIMATE_m 8
230  #endif
231 #endif
232 
233 /* window.c */
234 INT Y(m2K)(const INT m);
235 
236 #if defined(NFFT_LDOUBLE)
237 #if HAVE_DECL_COPYSIGNL == 0
238 extern long double copysignl(long double, long double);
239 #endif
240 #if HAVE_DECL_NEXTAFTERL == 0
241 extern long double nextafterl(long double, long double);
242 #endif
243 #if HAVE_DECL_NANL == 0
244 extern long double nanl(const char *tag);
245 #endif
246 #if HAVE_DECL_CEILL == 0
247 extern long double ceill(long double);
248 #endif
249 #if HAVE_DECL_FLOORL == 0
250 extern long double floorl(long double);
251 #endif
252 #if HAVE_DECL_NEARBYINTL == 0
253 extern long double nearbyintl(long double);
254 #endif
255 #if HAVE_DECL_RINTL == 0
256 extern long double rintl(long double);
257 #endif
258 #if HAVE_DECL_ROUNDL == 0
259 extern long double roundl(long double);
260 #endif
261 #if HAVE_DECL_LRINTL == 0
262 extern long int lrintl(long double);
263 #endif
264 #if HAVE_DECL_LROUNDL == 0
265 extern long int lroundl(long double);
266 #endif
267 #if HAVE_DECL_LLRINTL == 0
268 extern long long int llrintl(long double);
269 #endif
270 #if HAVE_DECL_LLROUNDL == 0
271 extern long long int llroundl(long double);
272 #endif
273 #if HAVE_DECL_TRUNCL == 0
274 extern long double truncl(long double);
275 #endif
276 #if HAVE_DECL_FMODL == 0
277 extern long double fmodl(long double, long double);
278 #endif
279 #if HAVE_DECL_REMAINDERL == 0
280 extern long double remainderl(long double, long double);
281 #endif
282 #if HAVE_DECL_REMQUOL == 0
283 extern long double remquol(long double x, long double y, int *);
284 #endif
285 #if HAVE_DECL_FDIML == 0
286 extern long double fdiml(long double, long double);
287 #endif
288 #if HAVE_DECL_FMAXL == 0
289 extern long double fmaxl(long double, long double);
290 #endif
291 #if HAVE_DECL_FMINL == 0
292 extern long double fminl(long double, long double);
293 #endif
294 #if HAVE_DECL_FMAL == 0
295 extern long double fmal(long double x, long double y, long double z);
296 #endif
297 #if HAVE_DECL_FABSL == 0
298 extern long double fabsl(long double);
299 #endif
300 #if HAVE_DECL_SQRTL == 0
301 extern long double sqrtl(long double);
302 #endif
303 #if HAVE_DECL_CBRTL == 0
304 extern long double cbrtl(long double);
305 #endif
306 #if HAVE_DECL_HYPOTL == 0
307 extern long double hypotl(long double, long double);
308 #endif
309 #if HAVE_DECL_EXPL == 0
310 extern long double expl(long double);
311 #endif
312 #if HAVE_DECL_EXP2L == 0
313 extern long double exp2l(long double);
314 #endif
315 #if HAVE_DECL_EXPM1L == 0
316 extern long double expm1l(long double);
317 #endif
318 #if HAVE_DECL_LOGL == 0
319 extern long double logl(long double);
320 #endif
321 #if HAVE_DECL_LOG2L == 0
322 extern long double log2l(long double);
323 #endif
324 #if HAVE_DECL_LOG10L == 0
325 extern long double log10l(long double);
326 #endif
327 #if HAVE_DECL_LOG1PL == 0
328 extern long double log1pl(long double);
329 #endif
330 #if HAVE_DECL_LOGBL == 0
331 extern long double logbl(long double);
332 #endif
333 #if HAVE_DECL_ILOGBL == 0
334 extern int ilogbl(long double);
335 #endif
336 #if HAVE_DECL_MODFL == 0
337 extern long double modfl(long double, long double *);
338 #endif
339 #if HAVE_DECL_FREXPL == 0
340 extern long double frexpl(long double, int *);
341 #endif
342 #if HAVE_DECL_LDEXPL == 0
343 extern long double ldexpl(long double, int);
344 #endif
345 #if HAVE_DECL_SCALBNL == 0
346 extern long double scalbnl(long double, int);
347 #endif
348 #if HAVE_DECL_SCALBLNL == 0
349 extern long double scalblnl(long double, long int);
350 #endif
351 #if HAVE_DECL_POWL == 0
352 extern long double powl(long double, long double);
353 #endif
354 #if HAVE_DECL_COSL == 0
355 extern long double cosl(long double);
356 #endif
357 #if HAVE_DECL_SINL == 0
358 extern long double sinl(long double);
359 #endif
360 #if HAVE_DECL_TANL == 0
361 extern long double tanl(long double);
362 #endif
363 #if HAVE_DECL_COSHL == 0
364 extern long double coshl(long double);
365 #endif
366 #if HAVE_DECL_SINHL == 0
367 extern long double sinhl(long double);
368 #endif
369 #if HAVE_DECL_TANHL == 0
370 extern long double tanhl(long double);
371 #endif
372 #if HAVE_DECL_ACOSL == 0
373 extern long double acosl(long double);
374 #endif
375 #if HAVE_DECL_ASINL == 0
376 extern long double asinl(long double);
377 #endif
378 #if HAVE_DECL_ATANL == 0
379 extern long double atanl(long double);
380 #endif
381 #if HAVE_DECL_ATAN2L == 0
382 extern long double atan2l(long double, long double);
383 #endif
384 #if HAVE_DECL_ACOSHL == 0
385 extern long double acoshl(long double);
386 #endif
387 #if HAVE_DECL_ASINHL == 0
388 extern long double asinhl(long double);
389 #endif
390 #if HAVE_DECL_ATANHL == 0
391 extern long double atanhl(long double);
392 #endif
393 #if HAVE_DECL_TGAMMAL == 0
394 extern long double tgammal(long double);
395 #endif
396 #if HAVE_DECL_LGAMMAL == 0
397 extern long double lgammal(long double);
398 #endif
399 #if HAVE_DECL_J0L == 0
400 extern long double j0l(long double);
401 #endif
402 #if HAVE_DECL_J1L == 0
403 extern long double j1l(long double);
404 #endif
405 #if HAVE_DECL_JNL == 0
406 extern long double jnl(int, long double);
407 #endif
408 #if HAVE_DECL_Y0L == 0
409 extern long double y0l(long double);
410 #endif
411 #if HAVE_DECL_Y1L == 0
412 extern long double y1l(long double);
413 #endif
414 #if HAVE_DECL_YNL == 0
415 extern long double ynl(int, long double);
416 #endif
417 #if HAVE_DECL_ERFL == 0
418 extern long double erfl(long double);
419 #endif
420 #if HAVE_DECL_ERFCL == 0
421 extern long double erfcl(long double);
422 #endif
423 #if HAVE_DECL_CREALL == 0
424 extern long double creall(long double _Complex z);
425 #endif
426 #if HAVE_DECL_CIMAGL == 0
427 extern long double cimagl(long double _Complex z);
428 #endif
429 #if HAVE_DECL_CABSL == 0
430 extern long double cabsl(long double _Complex z);
431 #endif
432 #if HAVE_DECL_CARGL == 0
433 extern long double cargl(long double _Complex z);
434 #endif
435 #if HAVE_DECL_CONJL == 0
436 extern long double _Complex conjl(long double _Complex z);
437 #endif
438 #if HAVE_DECL_CPROJL == 0
439 extern long double _Complex cprojl(long double _Complex z);
440 #endif
441 #if HAVE_DECL_CSQRTL == 0
442 extern long double _Complex csqrtl(long double _Complex z);
443 #endif
444 #if HAVE_DECL_CEXPL == 0
445 extern long double _Complex cexpl(long double _Complex z);
446 #endif
447 #if HAVE_DECL_CLOGL == 0
448 extern long double _Complex clogl(long double _Complex z);
449 #endif
450 #if HAVE_DECL_CPOWL == 0
451 extern long double _Complex cpowl(long double _Complex z, long double _Complex w);
452 #endif
453 #if HAVE_DECL_CSINL == 0
454 extern long double _Complex csinl(long double _Complex z);
455 #endif
456 #if HAVE_DECL_CCOSL == 0
457 extern long double _Complex ccosl(long double _Complex z);
458 #endif
459 #if HAVE_DECL_CTANL == 0
460 extern long double _Complex ctanl(long double _Complex z);
461 #endif
462 #if HAVE_DECL_CASINL == 0
463 extern long double _Complex casinl(long double _Complex z);
464 #endif
465 #if HAVE_DECL_CACOSL == 0
466 extern long double _Complex cacosl(long double _Complex z);
467 #endif
468 #if HAVE_DECL_CATANL == 0
469 extern long double _Complex catanl(long double _Complex z);
470 #endif
471 #if HAVE_DECL_CSINHL == 0
472 extern long double _Complex csinhl(long double _Complex z);
473 #endif
474 #if HAVE_DECL_CCOSHL == 0
475 extern long double _Complex ccoshl(long double _Complex z);
476 #endif
477 #if HAVE_DECL_CTANHL == 0
478 extern long double _Complex ctanhl(long double _Complex z);
479 #endif
480 #if HAVE_DECL_CASINHL == 0
481 extern long double _Complex casinhl(long double _Complex z);
482 #endif
483 #if HAVE_DECL_CACOSHL == 0
484 extern long double _Complex cacoshl(long double _Complex z);
485 #endif
486 #if HAVE_DECL_CATANHL == 0
487 extern long double _Complex catanhl(long double _Complex z);
488 #endif
489 #define COPYSIGN copysignl
490 #define NEXTAFTER nextafterl
491 #define MKNAN nanl
492 #define CEIL ceill
493 #define FLOOR floorl
494 #define NEARBYINT nearbyintl
495 #define RINT rintl
496 #define ROUND roundl
497 #define LRINT lrintl
498 #define LROUND lroundl
499 #define LLRINT llrintl
500 #define LLROUND llroundl
501 #define TRUNC truncl
502 #define FMOD fmodl
503 #define REMAINDER remainderl
504 #define REMQUO remquol
505 #define FDIM fdiml
506 #define FMAX fmaxl
507 #define FMIN fminl
508 #define FFMA fmal
509 #define FABS fabsl
510 #define SQRT sqrtl
511 #define CBRT cbrtl
512 #define HYPOT hypotl
513 #define EXP expl
514 #define EXP2 exp2l
515 #define EXPM1 expm1l
516 #define LOG logl
517 #define LOG2 log2l
518 #define LOG10 log10l
519 #define LOG1P log1pl
520 #define LOGB logbl
521 #define ILOGB ilogbl
522 #define MODF modfl
523 #define FREXP frexpl
524 #define LDEXP ldexpl
525 #define SCALBN scalbnl
526 #define SCALBLN scalblnl
527 #define POW powl
528 #define COS cosl
529 #define SIN sinl
530 #define TAN tanl
531 #define COSH coshl
532 #define SINH sinhl
533 #define TANH tanhl
534 #define ACOS acosl
535 #define ASIN asinl
536 #define ATAN atanl
537 #define ATAN2 atan2l
538 #define ACOSH acoshl
539 #define ASINH asinhl
540 #define ATANH atanhl
541 #define TGAMMA tgammal
542 #define LGAMMA lgammal
543 #define J0 j0l
544 #define J1 j1l
545 #define JN jnl
546 #define Y0 y0l
547 #define Y1 y1l
548 #define YN ynl
549 #define ERF erfl
550 #define ERFC erfcl
551 #define CREAL creall
552 #define CIMAG cimagl
553 #define CABS cabsl
554 #define CARG cargl
555 #define CONJ conjl
556 #define CPROJ cprojl
557 #define CSQRT csqrtl
558 #define CEXP cexpl
559 #define CLOG clogl
560 #define CPOW cpowl
561 #define CSIN csinl
562 #define CCOS ccosl
563 #define CTAN ctanl
564 #define CASIN casinl
565 #define CACOS cacosl
566 #define CATAN catanl
567 #define CSINH csinhl
568 #define CCOSH ccoshl
569 #define CTANH ctanhl
570 #define CASINH casinhl
571 #define CACOSH cacoshl
572 #define CATANH catanhl
573 #elif defined(NFFT_SINGLE)
574 #if HAVE_DECL_COPYSIGNF == 0
575 extern float copysignf(float, float);
576 #endif
577 #if HAVE_DECL_NEXTAFTERF == 0
578 extern float nextafterf(float, float);
579 #endif
580 #if HAVE_DECL_NANF == 0
581 extern float nanf(const char *tag);
582 #endif
583 #if HAVE_DECL_CEILF == 0
584 extern float ceilf(float);
585 #endif
586 #if HAVE_DECL_FLOORF == 0
587 extern float floorf(float);
588 #endif
589 #if HAVE_DECL_NEARBYINTF == 0
590 extern float nearbyintf(float);
591 #endif
592 #if HAVE_DECL_RINTF == 0
593 extern float rintf(float);
594 #endif
595 #if HAVE_DECL_ROUNDF == 0
596 extern float roundf(float);
597 #endif
598 #if HAVE_DECL_LRINTF == 0
599 extern long int lrintf(float);
600 #endif
601 #if HAVE_DECL_LROUNDF == 0
602 extern long int lroundf(float);
603 #endif
604 #if HAVE_DECL_LLRINTF == 0
605 extern long long int llrintf(float);
606 #endif
607 #if HAVE_DECL_LLROUNDF == 0
608 extern long long int llroundf(float);
609 #endif
610 #if HAVE_DECL_TRUNCF == 0
611 extern float truncf(float);
612 #endif
613 #if HAVE_DECL_FMODF == 0
614 extern float fmodf(float, float);
615 #endif
616 #if HAVE_DECL_REMAINDERF == 0
617 extern float remainderf(float, float);
618 #endif
619 #if HAVE_DECL_REMQUOF == 0
620 extern float remquof(float x, float y, int *);
621 #endif
622 #if HAVE_DECL_FDIMF == 0
623 extern float fdimf(float, float);
624 #endif
625 #if HAVE_DECL_FMAXF == 0
626 extern float fmaxf(float, float);
627 #endif
628 #if HAVE_DECL_FMINF == 0
629 extern float fminf(float, float);
630 #endif
631 #if HAVE_DECL_FMAF == 0
632 extern float fmaf(float x, float y, float z);
633 #endif
634 #if HAVE_DECL_FABSF == 0
635 extern float fabsf(float);
636 #endif
637 #if HAVE_DECL_SQRTF == 0
638 extern float sqrtf(float);
639 #endif
640 #if HAVE_DECL_CBRTF == 0
641 extern float cbrtf(float);
642 #endif
643 #if HAVE_DECL_HYPOTF == 0
644 extern float hypotf(float, float);
645 #endif
646 #if HAVE_DECL_EXPF == 0
647 extern float expf(float);
648 #endif
649 #if HAVE_DECL_EXP2F == 0
650 extern float exp2f(float);
651 #endif
652 #if HAVE_DECL_EXPM1F == 0
653 extern float expm1f(float);
654 #endif
655 #if HAVE_DECL_LOGF == 0
656 extern float logf(float);
657 #endif
658 #if HAVE_DECL_LOG2F == 0
659 extern float log2f(float);
660 #endif
661 #if HAVE_DECL_LOG10F == 0
662 extern float log10f(float);
663 #endif
664 #if HAVE_DECL_LOG1PF == 0
665 extern float log1pf(float);
666 #endif
667 #if HAVE_DECL_LOGBF == 0
668 extern float logbf(float);
669 #endif
670 #if HAVE_DECL_ILOGBF == 0
671 extern int ilogbf(float);
672 #endif
673 #if HAVE_DECL_MODFF == 0
674 extern float modff(float, float *);
675 #endif
676 #if HAVE_DECL_FREXPF == 0
677 extern float frexpf(float, int *);
678 #endif
679 #if HAVE_DECL_LDEXPF == 0
680 extern float ldexpf(float, int);
681 #endif
682 #if HAVE_DECL_SCALBNF == 0
683 extern float scalbnf(float, int);
684 #endif
685 #if HAVE_DECL_SCALBLNF == 0
686 extern float scalblnf(float, long int);
687 #endif
688 #if HAVE_DECL_POWF == 0
689 extern float powf(float, float);
690 #endif
691 #if HAVE_DECL_COSF == 0
692 extern float cosf(float);
693 #endif
694 #if HAVE_DECL_SINF == 0
695 extern float sinf(float);
696 #endif
697 #if HAVE_DECL_TANF == 0
698 extern float tanf(float);
699 #endif
700 #if HAVE_DECL_COSHF == 0
701 extern float coshf(float);
702 #endif
703 #if HAVE_DECL_SINHF == 0
704 extern float sinhf(float);
705 #endif
706 #if HAVE_DECL_TANHF == 0
707 extern float tanhf(float);
708 #endif
709 #if HAVE_DECL_ACOSF == 0
710 extern float acosf(float);
711 #endif
712 #if HAVE_DECL_ASINF == 0
713 extern float asinf(float);
714 #endif
715 #if HAVE_DECL_ATANF == 0
716 extern float atanf(float);
717 #endif
718 #if HAVE_DECL_ATAN2F == 0
719 extern float atan2f(float, float);
720 #endif
721 #if HAVE_DECL_ACOSHF == 0
722 extern float acoshf(float);
723 #endif
724 #if HAVE_DECL_ASINHF == 0
725 extern float asinhf(float);
726 #endif
727 #if HAVE_DECL_ATANHF == 0
728 extern float atanhf(float);
729 #endif
730 #if HAVE_DECL_TGAMMAF == 0
731 extern float tgammaf(float);
732 #endif
733 #if HAVE_DECL_LGAMMAF == 0
734 extern float lgammaf(float);
735 #endif
736 #if HAVE_DECL_J0F == 0
737 extern float j0f(float);
738 #endif
739 #if HAVE_DECL_J1F == 0
740 extern float j1f(float);
741 #endif
742 #if HAVE_DECL_JNF == 0
743 extern float jnf(int, float);
744 #endif
745 #if HAVE_DECL_Y0F == 0
746 extern float y0f(float);
747 #endif
748 #if HAVE_DECL_Y1F == 0
749 extern float y1f(float);
750 #endif
751 #if HAVE_DECL_YNF == 0
752 extern float ynf(int, float);
753 #endif
754 #if HAVE_DECL_ERFF == 0
755 extern float erff(float);
756 #endif
757 #if HAVE_DECL_ERFCF == 0
758 extern float erfcf(float);
759 #endif
760 #if HAVE_DECL_CREALF == 0
761 extern float crealf(float _Complex z);
762 #endif
763 #if HAVE_DECL_CIMAGF == 0
764 extern float cimagf(float _Complex z);
765 #endif
766 #if HAVE_DECL_CABSF == 0
767 extern float cabsf(float _Complex z);
768 #endif
769 #if HAVE_DECL_CARGF == 0
770 extern float cargf(float _Complex z);
771 #endif
772 #if HAVE_DECL_CONJF == 0
773 extern float _Complex conjf(float _Complex z);
774 #endif
775 #if HAVE_DECL_CPROJF == 0
776 extern float _Complex cprojf(float _Complex z);
777 #endif
778 #if HAVE_DECL_CSQRTF == 0
779 extern float _Complex csqrtf(float _Complex z);
780 #endif
781 #if HAVE_DECL_CEXPF == 0
782 extern float _Complex cexpf(float _Complex z);
783 #endif
784 #if HAVE_DECL_CLOGF == 0
785 extern float _Complex clogf(float _Complex z);
786 #endif
787 #if HAVE_DECL_CPOWF == 0
788 extern float _Complex cpowf(float _Complex z, float _Complex w);
789 #endif
790 #if HAVE_DECL_CSINF == 0
791 extern float _Complex csinf(float _Complex z);
792 #endif
793 #if HAVE_DECL_CCOSF == 0
794 extern float _Complex ccosf(float _Complex z);
795 #endif
796 #if HAVE_DECL_CTANF == 0
797 extern float _Complex ctanf(float _Complex z);
798 #endif
799 #if HAVE_DECL_CASINF == 0
800 extern float _Complex casinf(float _Complex z);
801 #endif
802 #if HAVE_DECL_CACOSF == 0
803 extern float _Complex cacosf(float _Complex z);
804 #endif
805 #if HAVE_DECL_CATANF == 0
806 extern float _Complex catanf(float _Complex z);
807 #endif
808 #if HAVE_DECL_CSINHF == 0
809 extern float _Complex csinhf(float _Complex z);
810 #endif
811 #if HAVE_DECL_CCOSHF == 0
812 extern float _Complex ccoshf(float _Complex z);
813 #endif
814 #if HAVE_DECL_CTANHF == 0
815 extern float _Complex ctanhf(float _Complex z);
816 #endif
817 #if HAVE_DECL_CASINHF == 0
818 extern float _Complex casinhf(float _Complex z);
819 #endif
820 #if HAVE_DECL_CACOSHF == 0
821 extern float _Complex cacoshf(float _Complex z);
822 #endif
823 #if HAVE_DECL_CATANHF == 0
824 extern float _Complex catanhf(float _Complex z);
825 #endif
826 #define COPYSIGN copysignf
827 #define NEXTAFTER nextafterf
828 #define MKNAN nanf
829 #define CEIL ceilf
830 #define FLOOR floorf
831 #define NEARBYINT nearbyintf
832 #define RINT rintf
833 #define ROUND roundf
834 #define LRINT lrintf
835 #define LROUND lroundf
836 #define LLRINT llrintf
837 #define LLROUND llroundf
838 #define TRUNC truncf
839 #define FMOD fmodf
840 #define REMAINDER remainderf
841 #define REMQUO remquof
842 #define FDIM fdimf
843 #define FMAX fmaxf
844 #define FMIN fminf
845 #define FFMA fmaf
846 #define FABS fabsf
847 #define SQRT sqrtf
848 #define CBRT cbrtf
849 #define HYPOT hypotf
850 #define EXP expf
851 #define EXP2 exp2f
852 #define EXPM1 expm1f
853 #define LOG logf
854 #define LOG2 log2f
855 #define LOG10 log10f
856 #define LOG1P log1pf
857 #define LOGB logbf
858 #define ILOGB ilogbf
859 #define MODF modff
860 #define FREXP frexpf
861 #define LDEXP ldexpf
862 #define SCALBN scalbnf
863 #define SCALBLN scalblnf
864 #define POW powf
865 #define COS cosf
866 #define SIN sinf
867 #define TAN tanf
868 #define COSH coshf
869 #define SINH sinhf
870 #define TANH tanhf
871 #define ACOS acosf
872 #define ASIN asinf
873 #define ATAN atanf
874 #define ATAN2 atan2f
875 #define ACOSH acoshf
876 #define ASINH asinhf
877 #define ATANH atanhf
878 #define TGAMMA tgammaf
879 #define LGAMMA lgammaf
880 #define J0 j0f
881 #define J1 j1f
882 #define JN jnf
883 #define Y0 y0f
884 #define Y1 y1f
885 #define YN ynf
886 #define ERF erff
887 #define ERFC erfcf
888 #define CREAL crealf
889 #define CIMAG cimagf
890 #define CABS cabsf
891 #define CARG cargf
892 #define CONJ conjf
893 #define CPROJ cprojf
894 #define CSQRT csqrtf
895 #define CEXP cexpf
896 #define CLOG clogf
897 #define CPOW cpowf
898 #define CSIN csinf
899 #define CCOS ccosf
900 #define CTAN ctanf
901 #define CASIN casinf
902 #define CACOS cacosf
903 #define CATAN catanf
904 #define CSINH csinhf
905 #define CCOSH ccoshf
906 #define CTANH ctanhf
907 #define CASINH casinhf
908 #define CACOSH cacoshf
909 #define CATANH catanhf
910 #else
911 #if HAVE_DECL_COPYSIGN == 0
912 extern double copysign(double, double);
913 #endif
914 #if HAVE_DECL_NEXTAFTER == 0
915 extern double nextafter(double, double);
916 #endif
917 #if HAVE_DECL_NAN == 0
918 extern double nan(const char *tag);
919 #endif
920 #if HAVE_DECL_CEIL == 0
921 extern double ceil(double);
922 #endif
923 #if HAVE_DECL_FLOOR == 0
924 extern double floor(double);
925 #endif
926 #if HAVE_DECL_NEARBYINT == 0
927 extern double nearbyint(double);
928 #endif
929 #if HAVE_DECL_RINT == 0
930 extern double rint(double);
931 #endif
932 #if HAVE_DECL_ROUND == 0
933 extern double round(double);
934 #endif
935 #if HAVE_DECL_LRINT == 0
936 extern long int lrint(double);
937 #endif
938 #if HAVE_DECL_LROUND == 0
939 extern long int lround(double);
940 #endif
941 #if HAVE_DECL_LLRINT == 0
942 extern long long int llrint(double);
943 #endif
944 #if HAVE_DECL_LLROUND == 0
945 extern long long int llround(double);
946 #endif
947 #if HAVE_DECL_TRUNC == 0
948 extern double trunc(double);
949 #endif
950 #if HAVE_DECL_FMOD == 0
951 extern double fmod(double, double);
952 #endif
953 #if HAVE_DECL_REMAINDER == 0
954 extern double remainder(double, double);
955 #endif
956 #if HAVE_DECL_REMQUO == 0
957 extern double remquo(double x, double y, int *);
958 #endif
959 #if HAVE_DECL_FDIM == 0
960 extern double fdim(double, double);
961 #endif
962 #if HAVE_DECL_FMAX == 0
963 extern double fmax(double, double);
964 #endif
965 #if HAVE_DECL_FMIN == 0
966 extern double fmin(double, double);
967 #endif
968 #if HAVE_DECL_FMA == 0
969 extern double fma(double x, double y, double z);
970 #endif
971 #if HAVE_DECL_FABS == 0
972 extern double fabs(double);
973 #endif
974 #if HAVE_DECL_SQRT == 0
975 extern double sqrt(double);
976 #endif
977 #if HAVE_DECL_CBRT == 0
978 extern double cbrt(double);
979 #endif
980 #if HAVE_DECL_HYPOT == 0
981 extern double hypot(double, double);
982 #endif
983 #if HAVE_DECL_EXP == 0
984 extern double exp(double);
985 #endif
986 #if HAVE_DECL_EXP2 == 0
987 extern double exp2(double);
988 #endif
989 #if HAVE_DECL_EXPM1 == 0
990 extern double expm1(double);
991 #endif
992 #if HAVE_DECL_LOG == 0
993 extern double log(double);
994 #endif
995 #if HAVE_DECL_LOG2 == 0
996 extern double log2(double);
997 #endif
998 #if HAVE_DECL_LOG10 == 0
999 extern double log10(double);
1000 #endif
1001 #if HAVE_DECL_LOG1P == 0
1002 extern double log1p(double);
1003 #endif
1004 #if HAVE_DECL_LOGB == 0
1005 extern double logb(double);
1006 #endif
1007 #if HAVE_DECL_ILOGB == 0
1008 extern int ilogb(double);
1009 #endif
1010 #if HAVE_DECL_MODF == 0
1011 extern double modf(double, double *);
1012 #endif
1013 #if HAVE_DECL_FREXP == 0
1014 extern double frexp(double, int *);
1015 #endif
1016 #if HAVE_DECL_LDEXP == 0
1017 extern double ldexp(double, int);
1018 #endif
1019 #if HAVE_DECL_SCALBN == 0
1020 extern double scalbn(double, int);
1021 #endif
1022 #if HAVE_DECL_SCALBLN == 0
1023 extern double scalbln(double, long int);
1024 #endif
1025 #if HAVE_DECL_POW == 0
1026 extern double pow(double, double);
1027 #endif
1028 #if HAVE_DECL_COS == 0
1029 extern double cos(double);
1030 #endif
1031 #if HAVE_DECL_SIN == 0
1032 extern double sin(double);
1033 #endif
1034 #if HAVE_DECL_TAN == 0
1035 extern double tan(double);
1036 #endif
1037 #if HAVE_DECL_COSH == 0
1038 extern double cosh(double);
1039 #endif
1040 #if HAVE_DECL_SINH == 0
1041 extern double sinh(double);
1042 #endif
1043 #if HAVE_DECL_TANH == 0
1044 extern double tanh(double);
1045 #endif
1046 #if HAVE_DECL_ACOS == 0
1047 extern double acos(double);
1048 #endif
1049 #if HAVE_DECL_ASIN == 0
1050 extern double asin(double);
1051 #endif
1052 #if HAVE_DECL_ATAN == 0
1053 extern double atan(double);
1054 #endif
1055 #if HAVE_DECL_ATAN2 == 0
1056 extern double atan2(double, double);
1057 #endif
1058 #if HAVE_DECL_ACOSH == 0
1059 extern double acosh(double);
1060 #endif
1061 #if HAVE_DECL_ASINH == 0
1062 extern double asinh(double);
1063 #endif
1064 #if HAVE_DECL_ATANH == 0
1065 extern double atanh(double);
1066 #endif
1067 #if HAVE_DECL_TGAMMA == 0
1068 extern double tgamma(double);
1069 #endif
1070 #if HAVE_DECL_LGAMMA == 0
1071 extern double lgamma(double);
1072 #endif
1073 #if HAVE_DECL_J0 == 0
1074 extern double j0(double);
1075 #endif
1076 #if HAVE_DECL_J1 == 0
1077 extern double j1(double);
1078 #endif
1079 #if HAVE_DECL_JN == 0
1080 extern double jn(int, double);
1081 #endif
1082 #if HAVE_DECL_Y0 == 0
1083 extern double y0(double);
1084 #endif
1085 #if HAVE_DECL_Y1 == 0
1086 extern double y1(double);
1087 #endif
1088 #if HAVE_DECL_YN == 0
1089 extern double yn(int, double);
1090 #endif
1091 #if HAVE_DECL_ERF == 0
1092 extern double erf(double);
1093 #endif
1094 #if HAVE_DECL_ERFC == 0
1095 extern double erfc(double);
1096 #endif
1097 #if HAVE_DECL_CREAL == 0
1098 extern double creal(double _Complex z);
1099 #endif
1100 #if HAVE_DECL_CIMAG == 0
1101 extern double cimag(double _Complex z);
1102 #endif
1103 #if HAVE_DECL_CABS == 0
1104 extern double cabs(double _Complex z);
1105 #endif
1106 #if HAVE_DECL_CARG == 0
1107 extern double carg(double _Complex z);
1108 #endif
1109 #if HAVE_DECL_CONJ == 0
1110 extern double _Complex conj(double _Complex z);
1111 #endif
1112 #if HAVE_DECL_CPROJ == 0
1113 extern double _Complex cproj(double _Complex z);
1114 #endif
1115 #if HAVE_DECL_CSQRT == 0
1116 extern double _Complex csqrt(double _Complex z);
1117 #endif
1118 #if HAVE_DECL_CEXP == 0
1119 extern double _Complex cexp(double _Complex z);
1120 #endif
1121 #if HAVE_DECL_CLOG == 0
1122 extern double _Complex clog(double _Complex z);
1123 #endif
1124 #if HAVE_DECL_CPOW == 0
1125 extern double _Complex cpow(double _Complex z, double _Complex w);
1126 #endif
1127 #if HAVE_DECL_CSIN == 0
1128 extern double _Complex csin(double _Complex z);
1129 #endif
1130 #if HAVE_DECL_CCOS == 0
1131 extern double _Complex ccos(double _Complex z);
1132 #endif
1133 #if HAVE_DECL_CTAN == 0
1134 extern double _Complex ctan(double _Complex z);
1135 #endif
1136 #if HAVE_DECL_CASIN == 0
1137 extern double _Complex casin(double _Complex z);
1138 #endif
1139 #if HAVE_DECL_CACOS == 0
1140 extern double _Complex cacos(double _Complex z);
1141 #endif
1142 #if HAVE_DECL_CATAN == 0
1143 extern double _Complex catan(double _Complex z);
1144 #endif
1145 #if HAVE_DECL_CSINH == 0
1146 extern double _Complex csinh(double _Complex z);
1147 #endif
1148 #if HAVE_DECL_CCOSH == 0
1149 extern double _Complex ccosh(double _Complex z);
1150 #endif
1151 #if HAVE_DECL_CTANH == 0
1152 extern double _Complex ctanh(double _Complex z);
1153 #endif
1154 #if HAVE_DECL_CASINH == 0
1155 extern double _Complex casinh(double _Complex z);
1156 #endif
1157 #if HAVE_DECL_CACOSH == 0
1158 extern double _Complex cacosh(double _Complex z);
1159 #endif
1160 #if HAVE_DECL_CATANH == 0
1161 extern double _Complex catanh(double _Complex z);
1162 #endif
1163 #define COPYSIGN copysign
1164 #define NEXTAFTER nextafter
1165 #define MKNAN nan
1166 #define CEIL ceil
1167 #define FLOOR floor
1168 #define NEARBYINT nearbyint
1169 #define RINT rint
1170 #define ROUND round
1171 #define LRINT lrint
1172 #define LROUND lround
1173 #define LLRINT llrint
1174 #define LLROUND llround
1175 #define TRUNC trunc
1176 #define FMOD fmod
1177 #define REMAINDER remainder
1178 #define REMQUO remquo
1179 #define FDIM fdim
1180 #define FMAX fmax
1181 #define FMIN fmin
1182 #define FFMA fma
1183 #define FABS fabs
1184 #define SQRT sqrt
1185 #define CBRT cbrt
1186 #define HYPOT hypot
1187 #define EXP exp
1188 #define EXP2 exp2
1189 #define EXPM1 expm1
1190 #define LOG log
1191 #define LOG2 log2
1192 #define LOG10 log10
1193 #define LOG1P log1p
1194 #define LOGB logb
1195 #define ILOGB ilogb
1196 #define MODF modf
1197 #define FREXP frexp
1198 #define LDEXP ldexp
1199 #define SCALBN scalbn
1200 #define SCALBLN scalbln
1201 #define POW pow
1202 #define COS cos
1203 #define SIN sin
1204 #define TAN tan
1205 #define COSH cosh
1206 #define SINH sinh
1207 #define TANH tanh
1208 #define ACOS acos
1209 #define ASIN asin
1210 #define ATAN atan
1211 #define ATAN2 atan2
1212 #define ACOSH acosh
1213 #define ASINH asinh
1214 #define ATANH atanh
1215 #define TGAMMA tgamma
1216 #define LGAMMA lgamma
1217 #define J0 j0
1218 #define J1 j1
1219 #define JN jn
1220 #define Y0 y0
1221 #define Y1 y1
1222 #define YN yn
1223 #define ERF erf
1224 #define ERFC erfc
1225 #define CREAL creal
1226 #define CIMAG cimag
1227 #define CABS cabs
1228 #define CARG carg
1229 #define CONJ conj
1230 #define CPROJ cproj
1231 #define CSQRT csqrt
1232 #define CEXP cexp
1233 #define CLOG clog
1234 #define CPOW cpow
1235 #define CSIN csin
1236 #define CCOS ccos
1237 #define CTAN ctan
1238 #define CASIN casin
1239 #define CACOS cacos
1240 #define CATAN catan
1241 #define CSINH csinh
1242 #define CCOSH ccosh
1243 #define CTANH ctanh
1244 #define CASINH casinh
1245 #define CACOSH cacosh
1246 #define CATANH catanh
1247 #endif
1248 
1249 #if defined(NFFT_LDOUBLE)
1250  #define MANT_DIG LDBL_MANT_DIG
1251  #define MIN_EXP LDBL_MIN_EXP
1252  #define MAX_EXP LDBL_MAX_EXP
1253 #elif defined(NFFT_SINGLE)
1254  #define MANT_DIG FLT_MANT_DIG
1255  #define MIN_EXP FLT_MIN_EXP
1256  #define MAX_EXP FLT_MAX_EXP
1257 #else
1258  #define MANT_DIG DBL_MANT_DIG
1259  #define MIN_EXP DBL_MIN_EXP
1260  #define MAX_EXP DBL_MAX_EXP
1261 #endif
1262 
1263 #if defined(FLT_ROUND)
1264  #if FLT_ROUND != -1
1265  #define FLTROUND 1.0
1266  #else
1267  #define FLTROUND 0.0
1268  #endif
1269 #else
1270  #define FLTROUND 0.0
1271 #endif
1272 
1273 #if HAVE_DECL_DRAND48 == 0
1274  extern double drand48(void);
1275 #endif
1276 #if HAVE_DECL_SRAND48 == 0
1277  extern void srand48(long int);
1278 #endif
1279 #define R_RADIX FLT_RADIX
1280 #define II _Complex_I
1281 
1282 /* format strings */
1283 #if defined(NFFT_LDOUBLE)
1284 # define __FGS__ "Lg"
1285 # define __FES__ "LE"
1286 # define __FE__ "% 36.32LE"
1287 # define __FI__ "%Lf"
1288 # define __FIS__ "Lf"
1289 # define __FR__ "%Le"
1290 #elif defined(NFFT_SINGLE)
1291 # define __FGS__ "g"
1292 # define __FES__ "E"
1293 # define __FE__ "% 12.8E"
1294 # define __FI__ "%f"
1295 # define __FIS__ "f"
1296 # define __FR__ "%e"
1297 #else
1298 # define __FGS__ "lg"
1299 # define __FES__ "lE"
1300 # define __FE__ "% 20.16lE"
1301 # define __FI__ "%lf"
1302 # define __FIS__ "lf"
1303 # define __FR__ "%le"
1304 #endif
1305 
1306 #define TRUE 1
1307 #define FALSE 0
1308 
1309 #if defined(_WIN32) || defined(_WIN64)
1310 # define __D__ "%Id"
1311 #else
1312 # define __D__ "%td"
1313 #endif
1314 
1316 #define UNUSED(x) (void)x
1317 
1318 #ifdef HAVE_ALLOCA
1319  /* Use alloca if available. */
1320  #ifndef alloca
1321  #ifdef __GNUC__
1322  /* No alloca defined but can use GCC's builtin version. */
1323  #define alloca __builtin_alloca
1324  #else
1325  /* No alloca defined and not using GCC. */
1326  #ifdef _MSC_VER
1327  /* Using Microsoft's C compiler. Include header file and use _alloca
1328  * defined therein. */
1329  #include <malloc.h>
1330  #define alloca _alloca
1331  #else
1332  /* Also not using Microsoft's C compiler. */
1333  #if HAVE_ALLOCA_H
1334  /* Alloca header is available. */
1335  #include <alloca.h>
1336  #else
1337  /* No alloca header available. */
1338  #ifdef _AIX
1339  /* We're using the AIX C compiler. Use pragma. */
1340  #pragma alloca
1341  #else
1342  /* Not using AIX compiler. */
1343  #ifndef alloca /* HP's cc +Olibcalls predefines alloca. */
1344  void *alloca(size_t);
1345  #endif
1346  #endif
1347  #endif
1348  #endif
1349  #endif
1350  #endif
1351  /* So we have alloca. */
1352  #define STACK_MALLOC(T, p, x) p = (T)alloca(x)
1353  #define STACK_FREE(x) /* Nothing. Cleanup done automatically. */
1354 #else /* ! HAVE_ALLOCA */
1355  /* Use malloc instead of alloca. So we allocate memory on the heap instead of
1356  * on the stack which is slower. */
1357  #define STACK_MALLOC(T, p, x) p = (T)Y(malloc)(x)
1358  #define STACK_FREE(x) Y(free)(x)
1359 #endif /* ! HAVE_ALLOCA */
1360 
1362 R Y(elapsed_seconds)(ticks t1, ticks t0);
1363 
1365 #define UNUSED(x) (void)x
1366 
1373 #ifdef MEASURE_TIME
1374  int MEASURE_TIME_r;
1375  double MEASURE_TIME_tt;
1376  ticks MEASURE_TIME_t0, MEASURE_TIME_t1;
1377 
1378 #define TIC(a) \
1379  ths->MEASURE_TIME_t[(a)]=0; \
1380  MEASURE_TIME_r=0; \
1381  /* DISABLED LOOP due to code blocks causing segfault when repeatedly run */ \
1382  /*while(ths->MEASURE_TIME_t[(a)]<0.01)*/ \
1383  { \
1384  MEASURE_TIME_r++; \
1385  MEASURE_TIME_t0 = getticks(); \
1386 
1387 /* THE MEASURED FUNCTION IS CALLED REPEATEDLY */
1388 
1389 #define TOC(a) \
1390  MEASURE_TIME_t1 = getticks(); \
1391  MEASURE_TIME_tt = Y(elapsed_seconds)(MEASURE_TIME_t1,MEASURE_TIME_t0);\
1392  ths->MEASURE_TIME_t[(a)]+=MEASURE_TIME_tt; \
1393  } \
1394  ths->MEASURE_TIME_t[(a)]/=MEASURE_TIME_r; \
1395 
1396 #else
1397 #define TIC(a)
1398 #define TOC(a)
1399 #endif
1400 
1401 #ifdef MEASURE_TIME_FFTW
1402 #define TIC_FFTW(a) TIC(a)
1403 #define TOC_FFTW(a) TOC(a)
1404 #else
1405 #define TIC_FFTW(a)
1406 #define TOC_FFTW(a)
1407 #endif
1408 
1409 /* sinc.c: */
1410 
1411 /* Sinus cardinalis. */
1412 R Y(sinc)(R x);
1413 
1414 /* lambda.c: */
1415 
1416 /* lambda(z, eps) = gamma(z + eps) / gamma(z + 1) */
1417 R Y(lambda)(R z, R eps);
1418 
1419 /* lambda2(mu, nu) = sqrt(gamma(mu + nu + 1) / (gamma(mu + 1) * gamma(nu + 1))) */
1420 R Y(lambda2)(R mu, R nu);
1421 
1422 /* bessel_i0.c: */
1423 R Y(bessel_i0)(R x);
1424 
1425 /* bspline.c: */
1426 R Y(bsplines)(const INT, const R x);
1427 
1428 /* float.c: */
1429 typedef enum {NFFT_SAFE__MIN = 1, NFFT_BASE = 2,
1430  NFFT_PRECISION = 3, NFFT_MANT_DIG = 4, NFFT_FLTROUND = 5, NFFT_E_MIN = 6,
1431  NFFT_R_MIN = 7, NFFT_E_MAX = 8, NFFT_R_MAX = 9} float_property;
1432 
1433 R Y(float_property)(float_property);
1434 R Y(prod_real)(R *vec, INT d);
1435 
1436 /* int.c: */
1437 INT Y(log2i)(const INT m);
1438 void Y(next_power_of_2_exp)(const INT N, INT *N2, INT *t);
1439 void Y(next_power_of_2_exp_int)(const int N, int *N2, int *t);
1440 
1441 /* error.c: */
1442 /* not used */ R Y(error_l_infty_double)(const R *x, const R *y, const INT n);
1443 /* not used */ R Y(error_l_infty_1_double)(const R *x, const R *y, const INT n, const R *z,
1444  const INT m);
1445 R Y(error_l_2_complex)(const C *x, const C *y, const INT n);
1446 /* not used */ R Y(error_l_2_double)(const R *x, const R *y, const INT n);
1447 
1448 /* sort.c: */
1449 void Y(sort_node_indices_radix_msdf)(INT n, INT *keys0, INT *keys1, INT rhigh);
1450 void Y(sort_node_indices_radix_lsdf)(INT n, INT *keys0, INT *keys1, INT rhigh);
1451 
1452 /* assert.c */
1453 void Y(assertion_failed)(const char *s, int line, const char *file);
1454 
1455 /* vector1.c */
1457 R Y(dot_double)(R *x, INT n);
1459 R Y(dot_w_complex)(C *x, R *w, INT n);
1461 R Y(dot_w_double)(R *x, R *w, INT n);
1463 R Y(dot_w_w2_complex)(C *x, R *w, R *w2, INT n);
1465 R Y(dot_w2_complex)(C *x, R *w2, INT n);
1466 
1467 /* vector2.c */
1469 void Y(cp_complex)(C *x, C *y, INT n);
1471 void Y(cp_double)(R *x, R *y, INT n);
1473 void Y(cp_a_complex)(C *x, R a, C *y, INT n);
1475 void Y(cp_a_double)(R *x, R a, R *y, INT n);
1477 void Y(cp_w_complex)(C *x, R *w, C *y, INT n);
1479 void Y(cp_w_double)(R *x, R *w, R *y, INT n);
1480 
1481 /* vector3.c */
1483 void Y(upd_axpy_double)(R *x, R a, R *y, INT n);
1485 void Y(upd_xpay_complex)(C *x, R a, C *y, INT n);
1487 void Y(upd_xpay_double)(R *x, R a, R *y, INT n);
1489 void Y(upd_axpby_complex)(C *x, R a, C *y, R b, INT n);
1491 void Y(upd_axpby_double)(R *x, R a, R *y, R b, INT n);
1493 void Y(upd_xpawy_complex)(C *x, R a, R *w, C *y, INT n);
1495 void Y(upd_xpawy_double)(R *x, R a, R *w, R *y, INT n);
1497 void Y(upd_axpwy_complex)(C *x, R a, R *w, C *y, INT n);
1499 void Y(upd_axpwy_double)(R *x, R a, R *w, R *y, INT n);
1500 
1501 /* voronoi.c */
1502 void Y(voronoi_weights_1d)(R *w, R *x, const INT M);
1503 
1504 /* damp.c */
1509 R Y(modified_fejer)(const INT N, const INT kk);
1511 R Y(modified_jackson2)(const INT N, const INT kk);
1513 R Y(modified_jackson4)(const INT N, const INT kk);
1515 R Y(modified_sobolev)(const R mu, const INT kk);
1517 R Y(modified_multiquadric)(const R mu, const R c, const INT kk);
1518 
1519 /* always check */
1520 #define CK(ex) \
1521  (void)((ex) || (Y(assertion_failed)(#ex, __LINE__, __FILE__), 0))
1522 
1523 #ifdef NFFT_DEBUG
1524  /* check only if debug enabled */
1525  #define A(ex) \
1526  (void)((ex) || (Y(assertion_failed)(#ex, __LINE__, __FILE__), 0))
1527 #else
1528  #define A(ex) /* nothing */
1529 #endif
1530 
1534 #endif
Generated on Thu Apr 12 2018 21:00:49 for NFFT by   doxygen 1.8.6