NFFT Logo 3.1.1 API Reference

nfft3.h

Go to the documentation of this file.
00001 /*
00002  * Copyright (c) 2002, 2009 Jens Keiner, Stefan Kunis, Daniel Potts
00003  *
00004  * This program is free software; you can redistribute it and/or modify it under
00005  * the terms of the GNU General Public License as published by the Free Software
00006  * Foundation; either version 2 of the License, or (at your option) any later
00007  * version.
00008  *
00009  * This program is distributed in the hope that it will be useful, but WITHOUT
00010  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
00011  * FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
00012  * details.
00013  *
00014  * You should have received a copy of the GNU General Public License along with
00015  * this program; if not, write to the Free Software Foundation, Inc., 51
00016  * Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
00017  */
00018 
00019 /* $Id: nfft3.h 3198 2009-05-27 14:16:50Z keiner $ */
00020 
00024 #ifndef __NFFT3_H__
00025 #define __NFFT3_H__
00026 
00028 #include <fftw3.h>
00029 
00030 /* config header */
00031 #include "nfft3conf.h"
00032 
00033 /* Malloc and free functions */
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 /* Malloc and free hooks */
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 /* Kaiser-Bessel is the default. */
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 /* Planner flags, i.e. constant symbols for precomputation and memory usage */
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   /* api */
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   /* internal*/
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   /* api */
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   /* api */
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   /* api */
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   /* api */
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   /* api */
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 /* Planner flags */
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 /* Precomputation flags */
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   /* Public members */
02032   int N;                              
02033   double *x;                          
02040   /* Private members */
02041   /*int NPT;*/                        
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 /* Flags for fpt_init() */
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 /* Flags for fpt_trafo(), dpt_transposed(), fpt_trafo(), fpt_transposed() */
02205 #define FPT_FUNCTION_VALUES   (1U << 5) 
02208 #define FPT_AL_SYMMETRY       (1U << 6) 
02210 /* Data structures */
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 /* Planner flags */
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 /* Helper macros*/
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 /* Planner flags, i.e. constant symbols for methods */
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 /* nfft3.h */

Generated on 17 Aug 2009 by Doxygen 1.5.3