00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00024 #ifndef __NFFT3_H__
00025 #define __NFFT3_H__
00026
00028 #include <fftw3.h>
00029
00030
00031 #include "nfft3conf.h"
00032
00033
00034 extern void *nfft_malloc(size_t n);
00035 extern void nfft_free(void *p);
00036 extern void nfft_die(const char *s);
00037
00038
00039 typedef void *(*nfft_malloc_type_function) (size_t n);
00040 typedef void (*nfft_free_type_function) (void *p);
00041 typedef void (*nfft_die_type_function) (const char *errString);
00042 extern nfft_malloc_type_function nfft_malloc_hook;
00043 extern nfft_free_type_function nfft_free_hook;
00044 extern nfft_die_type_function nfft_die_hook;
00045
00047 #define MACRO_MV_PLAN(float_type) \
00048 int N_total; \
00050 int M_total; \
00051 \
00052 float_type *f_hat; \
00054 float_type *f; \
00056 void (*mv_trafo)(void*); \
00057 void (*mv_adjoint)(void*); \
00058
00059 typedef struct
00060 {
00061 MACRO_MV_PLAN(fftw_complex)
00062 } mv_plan_complex;
00063
00064 typedef struct
00065 {
00066 MACRO_MV_PLAN(double)
00067 } mv_plan_double;
00068
00070 #if defined(DIRAC_DELTA)
00071 #define PHI_HUT(k,d) 1.0
00072 #define PHI(x,d) (fabs((x))<10e-8)? 1.0 : 0.0
00073 #define WINDOW_HELP_INIT(d)
00074 #define WINDOW_HELP_FINALIZE
00075 #define WINDOW_HELP_ESTIMATE_m {ths->m = 0;}
00076 #elif defined(GAUSSIAN)
00077 #define PHI_HUT(k,d) ((double)exp(-(pow(PI*(k)/ths->n[d],2.0)*ths->b[d])))
00078 #define PHI(x,d) ((double)exp(-pow((x)*ths->n[d],2.0)/ ths->b[d])/sqrt(PI*ths->b[d]))
00079 #define WINDOW_HELP_INIT \
00080 { \
00081 int WINDOW_idx; \
00082 ths->b = (double*) nfft_malloc(ths->d*sizeof(double)); \
00083 for(WINDOW_idx=0; WINDOW_idx<ths->d; WINDOW_idx++) \
00084 ths->b[WINDOW_idx]=((double)2*ths->sigma[WINDOW_idx])/ \
00085 (2*ths->sigma[WINDOW_idx]-1)*(((double)ths->m) / PI); \
00086 }
00087 #define WINDOW_HELP_FINALIZE {nfft_free(ths->b);}
00088 #define WINDOW_HELP_ESTIMATE_m {ths->m =12;}
00089 #elif defined(B_SPLINE)
00090 #define PHI_HUT(k,d) ((double)(((k)==0)? 1.0/ths->n[(d)] : \
00091 pow(sin((k)*PI/ths->n[(d)])/((k)*PI/ths->n[(d)]),2*ths->m)/ths->n[(d)]))
00092 #define PHI(x,d) (nfft_bspline(2*ths->m,((x)*ths->n[(d)])+ \
00093 (double)ths->m,ths->spline_coeffs)/ths->n[(d)])
00094 #define WINDOW_HELP_INIT \
00095 { \
00096 ths->spline_coeffs= (double*)nfft_malloc(2*ths->m*sizeof(double)); \
00097 }
00098 #define WINDOW_HELP_FINALIZE {nfft_free(ths->spline_coeffs);}
00099 #define WINDOW_HELP_ESTIMATE_m {ths->m =11;}
00100 #elif defined(SINC_POWER)
00101 #define PHI_HUT(k,d) (nfft_bspline(2*ths->m,((double)2*ths->m*(k))/ \
00102 ((2*ths->sigma[(d)]-1)*ths->n[(d)]/ths->sigma[(d)])+ (double)ths->m, \
00103 ths->spline_coeffs))
00104 #define PHI(x,d) ((double)(ths->n[(d)]/ths->sigma[(d)]*(2*ths->sigma[(d)]-1)/\
00105 (2*ths->m)*pow(nfft_sinc(PI*ths->n[(d)]/ths->sigma[(d)]*(x)* \
00106 (2*ths->sigma[(d)]-1)/(2*ths->m)),2*ths->m)/ths->n[(d)]))
00107 #define WINDOW_HELP_INIT \
00108 { \
00109 ths->spline_coeffs= (double*)nfft_malloc(2*ths->m*sizeof(double)); \
00110 }
00111 #define WINDOW_HELP_FINALIZE {nfft_free(ths->spline_coeffs);}
00112 #define WINDOW_HELP_ESTIMATE_m {ths->m = 9;}
00113 #else
00114 #define PHI_HUT(k,d) ((double)nfft_i0( ths->m*sqrt(\
00115 pow((double)(ths->b[d]),2.0) - pow(2.0*PI*(k)/ths->n[d],2.0))))
00116 #define PHI(x,d) ((double)((pow((double)(ths->m),2.0)\
00117 -pow((x)*ths->n[d],2.0))>0)? \
00118 sinh(ths->b[d]*sqrt(pow((double)(ths->m),2.0)- \
00119 pow((x)*ths->n[d],2.0)))/(PI*sqrt(pow((double)(ths->m),2.0)- \
00120 pow((x)*ths->n[d],2.0))): (((pow((double)(ths->m),2.0)- \
00121 pow((x)*ths->n[d],2.0))<0)? sin(ths->b[d]* \
00122 sqrt(pow(ths->n[d]*(x),2.0)-pow((double)(ths->m),2.0)))/ \
00123 (PI*sqrt(pow(ths->n[d]*(x),2.0)-pow((double)(ths->m),2.0))):1.0))
00124 #define WINDOW_HELP_INIT \
00125 { \
00126 int WINDOW_idx; \
00127 ths->b = (double*) nfft_malloc(ths->d*sizeof(double)); \
00128 for(WINDOW_idx=0; WINDOW_idx<ths->d; WINDOW_idx++) \
00129 ths->b[WINDOW_idx] = ((double)PI*(2.0-1.0/ths->sigma[WINDOW_idx])); \
00130 }
00131 #define WINDOW_HELP_FINALIZE {nfft_free(ths->b);}
00132 #define WINDOW_HELP_ESTIMATE_m {ths->m = 6;}
00133 #endif
00134
00135
00136
00137
00138
00186
00187
00198 #define PRE_PHI_HUT (1U<< 0)
00199
00211 #define FG_PSI (1U<< 1)
00212
00228 #define PRE_LIN_PSI (1U<< 2)
00229
00241 #define PRE_FG_PSI (1U<< 3)
00242
00253 #define PRE_PSI (1U<< 4)
00254
00266 #define PRE_FULL_PSI (1U<< 5)
00267
00277 #define MALLOC_X (1U<< 6)
00278
00289 #define MALLOC_F_HAT (1U<< 7)
00290
00300 #define MALLOC_F (1U<< 8)
00301
00311 #define FFT_OUT_OF_PLACE (1U<< 9)
00312
00322 #define FFTW_INIT (1U<< 10)
00323
00337 #define PRE_ONE_PSI (PRE_LIN_PSI| PRE_FG_PSI| PRE_PSI| PRE_FULL_PSI)
00338
00340 typedef struct
00341 {
00342
00343 MACRO_MV_PLAN(fftw_complex)
00344
00345 int d;
00346 int *N;
00347 double *sigma;
00348 int *n;
00351 int n_total;
00352 int m;
00358 double *b;
00360 int K;
00364 unsigned nfft_flags;
00371 unsigned fftw_flags;
00375 double *x;
00378 double MEASURE_TIME_t[3];
00381
00382 fftw_plan my_fftw_plan1;
00383 fftw_plan my_fftw_plan2;
00385 double **c_phi_inv;
00388 double *psi;
00391 int *psi_index_g;
00393 int *psi_index_f;
00396 fftw_complex *g;
00399 fftw_complex *g_hat;
00402 fftw_complex *g1;
00403 fftw_complex *g2;
00405 double *spline_coeffs;
00408 } nfft_plan;
00409
00410
00418 void ndft_trafo(nfft_plan *ths);
00419
00427 void ndft_adjoint(nfft_plan *ths);
00428
00436 void nfft_trafo(nfft_plan *ths);
00437 void nfft_trafo_1d(nfft_plan *ths);
00438 void nfft_trafo_2d(nfft_plan *ths);
00439 void nfft_trafo_3d(nfft_plan *ths);
00440
00448 void nfft_adjoint(nfft_plan *ths);
00449 void nfft_adjoint_1d(nfft_plan *ths);
00450 void nfft_adjoint_2d(nfft_plan *ths);
00451 void nfft_adjoint_3d(nfft_plan *ths);
00452
00462 void nfft_init_1d(nfft_plan *ths, int N1, int M);
00463
00474 void nfft_init_2d(nfft_plan *ths, int N1, int N2, int M);
00475
00487 void nfft_init_3d(nfft_plan *ths, int N1, int N2, int N3, int M);
00488
00499 void nfft_init(nfft_plan *ths, int d, int *N, int M);
00500
00514 void nfft_init_advanced(nfft_plan *ths, int d, int *N, int M,
00515 unsigned nfft_flags_on, unsigned nfft_flags_off);
00516
00531 void nfft_init_guru(nfft_plan *ths, int d, int *N, int M, int *n,
00532 int m, unsigned nfft_flags, unsigned fftw_flags);
00533
00534
00547 void nfft_precompute_one_psi(nfft_plan *ths);
00548
00553 void nfft_precompute_full_psi(nfft_plan *ths);
00554
00559 void nfft_precompute_psi(nfft_plan *ths);
00560
00565 void nfft_precompute_lin_psi(nfft_plan *ths);
00566
00574 void nfft_check(nfft_plan *ths);
00575
00583 void nfft_finalize(nfft_plan *ths);
00584
00588
00589
00590
00591
00598 #ifdef HAVE_NFCT
00599
00601 typedef struct
00602 {
00603
00604 MACRO_MV_PLAN(double)
00605
00606 int d;
00607 int *N;
00608 int *n;
00609 double *sigma;
00610 int m;
00612 double nfct_full_psi_eps;
00613 double *b;
00615 unsigned nfct_flags;
00616 unsigned fftw_flags;
00618 double *x;
00620 double MEASURE_TIME_t[3];
00623 fftw_plan my_fftw_r2r_plan;
00624 fftw_r2r_kind *r2r_kind;
00626 double **c_phi_inv;
00627 double *psi;
00628 int size_psi;
00629 int *psi_index_g;
00630 int *psi_index_f;
00632 double *g;
00633 double *g_hat;
00634 double *g1;
00635 double *g2;
00637 double *spline_coeffs;
00639 } nfct_plan;
00640
00641
00651 void nfct_init_1d( nfct_plan *ths_plan, int N0, int M_total);
00652
00663 void nfct_init_2d( nfct_plan *ths_plan, int N0, int N1, int M_total);
00664
00676 void nfct_init_3d( nfct_plan *ths_plan, int N0, int N1, int N2, int M_total);
00677
00688 void nfct_init( nfct_plan *ths_plan, int d, int *N, int M_total);
00689
00704 void nfct_init_guru( nfct_plan *ths_plan, int d, int *N, int M_total, int *n,
00705 int m, unsigned nfct_flags, unsigned fftw_flags);
00706
00716 void nfct_precompute_psi( nfct_plan *ths_plan);
00717
00718
00727 void nfct_trafo( nfct_plan *ths_plan);
00728
00737 void ndct_trafo( nfct_plan *ths_plan);
00738
00747 void nfct_adjoint( nfct_plan *ths_plan);
00748
00757 void ndct_adjoint( nfct_plan *ths_plan);
00758
00766 void nfct_finalize( nfct_plan *ths_plan);
00767
00777 double nfct_phi_hut( nfct_plan *ths_plan, int k, int d);
00778
00788 double nfct_phi ( nfct_plan *ths_plan, double x, int d);
00789
00797 int nfct_fftw_2N( int n);
00798
00806 int nfct_fftw_2N_rev(int n);
00807
00808 #endif
00809
00810
00811
00812 #ifdef HAVE_NFST
00813
00815 typedef struct
00816 {
00817
00818 MACRO_MV_PLAN(double)
00819
00820 int d;
00821 int *N;
00822 int *n;
00823 double *sigma;
00824 int m;
00826 double nfst_full_psi_eps;
00827 double *b;
00829 unsigned nfst_flags;
00830 unsigned fftw_flags;
00832 double *x;
00834 double MEASURE_TIME_t[3];
00837 fftw_plan my_fftw_r2r_plan;
00838 fftw_r2r_kind *r2r_kind;
00840 double **c_phi_inv;
00841 double *psi;
00842 int size_psi;
00843 int *psi_index_g;
00844 int *psi_index_f;
00847 double *g;
00848 double *g_hat;
00849 double *g1;
00850 double *g2;
00852 double *spline_coeffs;
00854 } nfst_plan;
00855
00856
00866 void nfst_init_1d( nfst_plan *ths_plan, int N0, int M_total);
00867
00878 void nfst_init_2d( nfst_plan *ths_plan, int N0, int N1, int M_total);
00879
00891 void nfst_init_3d( nfst_plan *ths_plan, int N0, int N1, int N2, int M_total);
00892
00903 void nfst_init( nfst_plan *ths_plan, int d, int *N, int M_total);
00904
00917 void nfst_init_m( nfst_plan *ths_plan, int d, int *N, int M_total, int m);
00918
00933 void nfst_init_guru( nfst_plan *ths_plan, int d, int *N, int M_total, int *n,
00934 int m, unsigned nfst_flags, unsigned fftw_flags);
00935
00945 void nfst_precompute_psi( nfst_plan *ths_plan);
00946
00955 void nfst_trafo( nfst_plan *ths_plan);
00956
00965 void ndst_trafo( nfst_plan *ths_plan);
00966
00967
00968
00977 void nfst_adjoint( nfst_plan *ths_plan);
00978
00987 void ndst_adjoint( nfst_plan *ths_plan);
00988
00996 void nfst_finalize( nfst_plan *ths_plan);
00997
01004 void nfst_full_psi( nfst_plan *ths_plan, double eps);
01005
01015 double nfst_phi_hut( nfst_plan *ths_plan, int k, int d);
01016
01026 double nfst_phi ( nfst_plan *ths_plan, double x, int d);
01027
01035 int nfst_fftw_2N( int n);
01036
01044 int nfst_fftw_2N_rev( int n);
01045
01046 #endif
01047
01051
01052
01053
01054
01061 #ifdef HAVE_NNFFT
01062
01071 #define MALLOC_V (1U<< 11)
01072
01074 typedef struct
01075 {
01076
01077 MACRO_MV_PLAN(fftw_complex)
01078
01079 int d;
01080 double *sigma;
01081 double *a;
01082 int *N;
01083 int *N1;
01084 int *aN1;
01085 int m;
01086 double *b;
01087 int K;
01089 int aN1_total;
01091 nfft_plan *direct_plan;
01092 unsigned nnfft_flags;
01093 int *n;
01095 double *x;
01096 double *v;
01098 double *c_phi_inv;
01099 double *psi;
01100 int size_psi;
01101 int *psi_index_g;
01102 int *psi_index_f;
01103 fftw_complex *F;
01104
01105 double *spline_coeffs;
01107 } nnfft_plan;
01108
01109
01121 void nnfft_init(nnfft_plan *ths_plan, int d, int N_total, int M_total, int *N);
01122
01137 void nnfft_init_guru(nnfft_plan *ths_plan, int d, int N_total, int M_total,
01138 int *N, int *N1, int m, unsigned nnfft_flags);
01139
01151 void nndft_trafo(nnfft_plan *ths_plan);
01152
01164 void nndft_adjoint(nnfft_plan *ths_plan);
01165
01177 void nnfft_trafo(nnfft_plan *ths_plan);
01178
01190 void nnfft_adjoint(nnfft_plan *ths_plan);
01191
01203 void nnfft_precompute_lin_psi(nnfft_plan *ths_plan);
01204
01217 void nnfft_precompute_psi(nnfft_plan *ths_plan);
01218
01232 void nnfft_precompute_full_psi(nnfft_plan *ths_plan);
01233
01246 void nnfft_precompute_phi_hut(nnfft_plan *ths_plan);
01247
01255 void nnfft_finalize(nnfft_plan *ths_plan);
01256
01257 #endif
01258
01262
01263
01264
01265
01272 #ifdef HAVE_NSFFT
01273
01282 #define NSDFT (1U<< 12)
01283
01285 typedef struct
01286 {
01287 MACRO_MV_PLAN(fftw_complex)
01288
01289 int d;
01290 int J;
01294 int sigma;
01296 unsigned flags;
01298 int *index_sparse_to_full;
01301 int r_act_nfft_plan;
01302 nfft_plan *act_nfft_plan;
01303 nfft_plan *center_nfft_plan;
01305 fftw_plan* set_fftw_plan1;
01306 fftw_plan* set_fftw_plan2;
01308 nfft_plan *set_nfft_plan_1d;
01309 nfft_plan *set_nfft_plan_2d;
01311 double *x_transposed;
01312 double *x_102,*x_201,*x_120,*x_021;
01314 } nsfft_plan;
01315
01326 void nsdft_trafo(nsfft_plan *ths);
01327
01338 void nsdft_adjoint(nsfft_plan *ths);
01339
01351 void nsfft_trafo(nsfft_plan *ths);
01352
01364 void nsfft_adjoint(nsfft_plan *ths);
01365
01373 void nsfft_cp(nsfft_plan *ths, nfft_plan *ths_nfft);
01374
01382 void nsfft_init_random_nodes_coeffs(nsfft_plan *ths);
01383
01396 void nsfft_init(nsfft_plan *ths, int d, int J, int M, int m, unsigned flags);
01397
01405 void nsfft_finalize(nsfft_plan *ths);
01406
01407 #endif
01408
01412
01413
01414
01415
01420 #ifdef HAVE_MRI
01421
01425 typedef struct
01426 {
01427
01428 MACRO_MV_PLAN(fftw_complex)
01429
01430 nfft_plan plan;
01431
01432 int N3;
01433 double sigma3;
01434 double *t;
01435 double *w;
01436 } mri_inh_2d1d_plan;
01437
01441 typedef struct
01442 {
01443
01444 MACRO_MV_PLAN(fftw_complex)
01445
01446 nfft_plan plan;
01447
01448 int N3;
01449 double sigma3;
01450 double *t;
01451 double *w;
01452 } mri_inh_3d_plan;
01453
01454
01467 void mri_inh_2d1d_trafo(mri_inh_2d1d_plan *ths);
01468
01481 void mri_inh_2d1d_adjoint(mri_inh_2d1d_plan *ths);
01482
01496 void mri_inh_2d1d_init_guru(mri_inh_2d1d_plan *ths, int *N, int M, int *n,
01497 int m, double sigma, unsigned nfft_flags, unsigned fftw_flags);
01498
01506 void mri_inh_2d1d_finalize(mri_inh_2d1d_plan *ths);
01507
01520 void mri_inh_3d_trafo(mri_inh_3d_plan *ths);
01521
01534 void mri_inh_3d_adjoint(mri_inh_3d_plan *ths);
01535
01536 void mri_inh_3d_init_guru(mri_inh_3d_plan *ths, int *N, int M, int *n,
01537 int m, double sigma, unsigned nfft_flags, unsigned fftw_flags);
01538
01546 void mri_inh_3d_finalize(mri_inh_3d_plan *ths);
01547
01548 #endif
01549
01553
01554
01555
01556
01808 #ifdef HAVE_NFSFT
01809
01810
01811
01831 #define NFSFT_NORMALIZED (1U << 0)
01832
01843 #define NFSFT_USE_NDFT (1U << 1)
01844
01856 #define NFSFT_USE_DPT (1U << 2)
01857
01871 #define NFSFT_MALLOC_X (1U << 3)
01872
01886 #define NFSFT_MALLOC_F_HAT (1U << 5)
01887
01901 #define NFSFT_MALLOC_F (1U << 6)
01902
01913 #define NFSFT_PRESERVE_F_HAT (1U << 7)
01914
01926 #define NFSFT_PRESERVE_X (1U << 8)
01927
01938 #define NFSFT_PRESERVE_F (1U << 9)
01939
01949 #define NFSFT_DESTROY_F_HAT (1U << 10)
01950
01961 #define NFSFT_DESTROY_X (1U << 11)
01962
01972 #define NFSFT_DESTROY_F (1U << 12)
01973
01974
01975
01985 #define NFSFT_NO_DIRECT_ALGORITHM (1U << 13)
01986
01996 #define NFSFT_NO_FAST_ALGORITHM (1U << 14)
01997
02005 #define NFSFT_ZERO_F_HAT (1U << 16)
02006
02007
02008
02017 #define NFSFT_INDEX(k,n,plan) ((2*(plan)->N+2)*((plan)->N-n+1)+(plan)->N+k+1)
02018
02023 #define NFSFT_F_HAT_SIZE(N) ((2*N+2)*(2*N+2))
02024
02026 typedef struct
02027 {
02029 MACRO_MV_PLAN(fftw_complex)
02030
02031
02032 int N;
02033 double *x;
02040
02041
02043 int t;
02045 unsigned int flags;
02046 nfft_plan plan_nfft;
02047 fftw_complex *f_hat_intern;
02049 } nfsft_plan;
02050
02060 void nfsft_init(nfsft_plan *plan, int N, int M);
02061
02072 void nfsft_init_advanced(nfsft_plan* plan, int N, int M, unsigned int
02073 nfsft_flags);
02074
02086 void nfsft_init_guru(nfsft_plan *plan, int N, int M, unsigned int nfsft_flags,
02087 unsigned int nfft_flags, int nfft_cutoff);
02088
02102 void nfsft_precompute(int N, double kappa, unsigned int nfsft_flags,
02103 unsigned int fpt_flags);
02104
02110 void nfsft_forget(void);
02111
02123 void ndsft_trafo(nfsft_plan* plan);
02124
02137 void ndsft_adjoint(nfsft_plan* plan);
02138
02150 void nfsft_trafo(nfsft_plan* plan);
02151
02164 void nfsft_adjoint(nfsft_plan* plan);
02165
02173 void nfsft_finalize(nfsft_plan* plan);
02174
02175 void nfsft_precompute_x(nfsft_plan *plan);
02176
02177 #endif
02178
02182
02183
02184
02185
02194 #ifdef HAVE_FPT
02195
02196
02197 #define FPT_NO_FAST_ALGORITHM (1U << 2)
02198 #define FPT_NO_DIRECT_ALGORITHM (1U << 3)
02199 #define FPT_NO_STABILIZATION (1U << 0)
02202 #define FPT_PERSISTENT_DATA (1U << 4)
02204
02205 #define FPT_FUNCTION_VALUES (1U << 5)
02208 #define FPT_AL_SYMMETRY (1U << 6)
02210
02211 typedef struct fpt_set_s_ *fpt_set;
02228 fpt_set fpt_init(const int M, const int t, const unsigned int flags);
02229
02250 void fpt_precompute(fpt_set set, const int m, double *alpha, double *beta,
02251 double *gam, int k_start, const double threshold);
02252
02263 void dpt_trafo(fpt_set set, const int m, const fftw_complex *x, fftw_complex *y,
02264 const int k_end, const unsigned int flags);
02265
02276 void fpt_trafo(fpt_set set, const int m, const fftw_complex *x, fftw_complex *y,
02277 const int k_end, const unsigned int flags);
02278
02289 void dpt_transposed(fpt_set set, const int m, fftw_complex *x,
02290 fftw_complex *y, const int k_end, const unsigned int flags);
02291
02302 void fpt_transposed(fpt_set set, const int m, fftw_complex *x,
02303 fftw_complex *y, const int k_end, const unsigned int flags);
02304
02305 void fpt_finalize(fpt_set set);
02306
02307 #endif
02308
02312
02313
02314
02315
02326 #ifdef HAVE_NFSOFT
02327
02328
02348 #define NFSOFT_NORMALIZED (1U << 0)
02349
02360 #define NFSOFT_USE_NDFT (1U << 1)
02361
02372 #define NFSOFT_USE_DPT (1U << 2)
02373
02387 #define NFSOFT_MALLOC_X (1U << 3)
02388
02396 #define NFSOFT_REPRESENT (1U << 4)
02397
02398
02412 #define NFSOFT_MALLOC_F_HAT (1U << 5)
02413
02427 #define NFSOFT_MALLOC_F (1U << 6)
02428
02439 #define NFSOFT_PRESERVE_F_HAT (1U << 7)
02440
02452 #define NFSOFT_PRESERVE_X (1U << 8)
02453
02464 #define NFSOFT_PRESERVE_F (1U << 9)
02465
02475 #define NFSOFT_DESTROY_F_HAT (1U << 10)
02476
02487 #define NFSOFT_DESTROY_X (1U << 11)
02488
02498 #define NFSOFT_DESTROY_F (1U << 12)
02499
02508 #define NFSOFT_NO_STABILIZATION (1U << 13)
02509
02519 #define NFSOFT_CHOOSE_DPT (1U << 14)
02520
02531 #define NFSOFT_SOFT (1U << 15)
02532
02533
02541 #define NFSOFT_ZERO_F_HAT (1U << 16)
02542
02543
02544
02550 #define NFSOFT_INDEX(m,n,l,B) (((l)+((B)+1))+(2*(B)+2)*(((n)+((B)+1))+(2*(B)+2)*((m)+((B)+1))))
02551 #define NFSOFT_INDEX_TWO(m,n,l,B) ((B+1)*(B+1)+(B+1)*(B+1)*(m+B)-((m-1)*m*(2*m-1)+(B+1)*(B+2)*(2*B+3))/6)+(posN(n,m,B))+(l-MAX(ABS(m),ABS(n)))
02552 int posN(int n,int m, int B);
02553
02558 #define NFSOFT_F_HAT_SIZE(B) (((B)+1)*(4*((B)+1)*((B)+1)-1)/3)
02559
02561 typedef struct nfsoft_plan_
02562 {
02564 MACRO_MV_PLAN(fftw_complex)
02565
02566 double *x;
02568 fftw_complex *wig_coeffs;
02570 fftw_complex *cheby;
02572 fftw_complex *aux;
02576 int t;
02578 unsigned int flags;
02579 nfft_plan nfft_plan;
02580 fftw_plan fftw_plan;
02581 fpt_set fpt_set;
02583 int fpt_kappa;
02585 } nfsoft_plan;
02586
02587
02596 void nfsoft_precompute(nfsoft_plan *plan);
02597
02608 fpt_set SO3_single_fpt_init(int l, int k, int m, unsigned int flags, int kappa);
02609 void SO3_fpt(fftw_complex *coeffs, fpt_set set, int l, int k, int m, unsigned int nfsoft_flags);
02610 void SO3_fpt_transposed(fftw_complex *coeffs,fpt_set set,int l, int k, int m,unsigned int nfsoft_flags);
02611
02612
02622 void nfsoft_init(nfsoft_plan *plan, int N, int M);
02633 void nfsoft_init_advanced(nfsoft_plan *plan, int N, int M,unsigned int nfsoft_flags);
02647 void nfsoft_init_guru(nfsoft_plan *plan, int N, int M,unsigned int nfsoft_flags,unsigned int nfft_flags,int nfft_cutoff,int fpt_kappa);
02648
02660 void nfsoft_trafo(nfsoft_plan *plan_nfsoft);
02673 void nfsoft_adjoint(nfsoft_plan *plan_nfsoft);
02681 void nfsoft_finalize(nfsoft_plan *plan);
02682
02683
02684 #endif
02685
02690
02691
02692
02693
02698
02699
02706 #define LANDWEBER (1U<< 0)
02707
02714 #define STEEPEST_DESCENT (1U<< 1)
02715
02723 #define CGNR (1U<< 2)
02724
02732 #define CGNE (1U<< 3)
02733
02740 #define NORMS_FOR_LANDWEBER (1U<< 4)
02741
02748 #define PRECOMPUTE_WEIGHT (1U<< 5)
02749
02756 #define PRECOMPUTE_DAMP (1U<< 6)
02757
02758
02759 typedef struct
02760 {
02761 mv_plan_complex *mv;
02762 unsigned flags;
02764 double *w;
02765 double *w_hat;
02767 fftw_complex *y;
02769 fftw_complex *f_hat_iter;
02771 fftw_complex *r_iter;
02772 fftw_complex *z_hat_iter;
02774 fftw_complex *p_hat_iter;
02775 fftw_complex *v_iter;
02777 double alpha_iter;
02778 double beta_iter;
02780 double dot_r_iter;
02781 double dot_r_iter_old;
02782 double dot_z_hat_iter;
02784 double dot_z_hat_iter_old;
02785 double dot_p_hat_iter;
02787 double dot_v_iter;
02788 } solver_plan_complex;
02789
02790 void solver_init_advanced_complex(solver_plan_complex* ths, mv_plan_complex *mv, unsigned flags);
02791 void solver_init_complex(solver_plan_complex* ths, mv_plan_complex *mv);
02792 void solver_before_loop_complex(solver_plan_complex* ths);
02793 void solver_loop_one_step_complex(solver_plan_complex *ths);
02794 void solver_finalize_complex(solver_plan_complex *ths);
02795
02796 typedef struct
02797 {
02798 mv_plan_double *mv;
02799 unsigned flags;
02801 double *w;
02802 double *w_hat;
02804 double *y;
02806 double *f_hat_iter;
02808 double *r_iter;
02809 double *z_hat_iter;
02811 double *p_hat_iter;
02812 double *v_iter;
02814 double alpha_iter;
02815 double beta_iter;
02817 double dot_r_iter;
02818 double dot_r_iter_old;
02819 double dot_z_hat_iter;
02821 double dot_z_hat_iter_old;
02822 double dot_p_hat_iter;
02824 double dot_v_iter;
02825 } solver_plan_double;
02826
02827 void solver_init_advanced_double(solver_plan_double* ths, mv_plan_double *mv, unsigned flags);
02828 void solver_init_double(solver_plan_double* ths, mv_plan_double *mv);
02829 void solver_before_loop_double(solver_plan_double* ths);
02830 void solver_loop_one_step_double(solver_plan_double *ths);
02831 void solver_finalize_double(solver_plan_double *ths);
02832
02836 #endif
02837