00001 #include <stdio.h>
00002 #include <math.h>
00003 #include <string.h>
00004 #include <stdlib.h>
00005 #include "util.h"
00006 #include "nfft3.h"
00007
00015 void simple_test_infft_1d(int N, int M, int iter)
00016 {
00017 int k,l;
00018 nfft_plan p;
00019 infft_plan ip;
00022 nfft_init_1d(&p, N, M);
00023
00025 nfft_vrand_shifted_unit_double(p.x,p.M_total);
00026
00028 if(p.nfft_flags & PRE_ONE_PSI)
00029 nfft_precompute_one_psi(&p);
00030
00032 infft_init(&ip,&p);
00033
00035 nfft_vrand_unit_complex(ip.y,p.M_total);
00036 nfft_vpr_complex(ip.y,p.M_total,"Given data, vector y");
00037
00039 for(k=0;k<p.N_total;k++)
00040 ip.f_hat_iter[k]=0;
00041
00042 nfft_vpr_complex(ip.f_hat_iter,p.N_total,"Initial guess, vector f_hat_iter");
00043
00044 NFFT_SWAP_complex(ip.f_hat_iter,p.f_hat);
00045 nfft_trafo(&p);
00046 nfft_vpr_complex(p.f,p.M_total,"Data fit, vector f");
00047 NFFT_SWAP_complex(ip.f_hat_iter,p.f_hat);
00048
00049 infft_before_loop(&ip);
00050 printf("\n Residual r=%e\n",ip.dot_r_iter);
00051
00052 for(l=0;l<iter;l++)
00053 {
00054 printf("\n********** Iteration l=%d **********\n",l);
00055 infft_loop_one_step(&ip);
00056 nfft_vpr_complex(ip.f_hat_iter,p.N_total,
00057 "Approximate solution, vector f_hat_iter");
00058
00059 NFFT_SWAP_complex(ip.f_hat_iter,p.f_hat);
00060 nfft_trafo(&p);
00061 nfft_vpr_complex(p.f,p.M_total,"Data fit, vector f");
00062 NFFT_SWAP_complex(ip.f_hat_iter,p.f_hat);
00063
00064 printf("\n Residual r=%e\n",ip.dot_r_iter);
00065 }
00066
00067 infft_finalize(&ip);
00068 nfft_finalize(&p);
00069 }
00070
00072 int main()
00073 {
00074 printf("\n Computing a one dimensional inverse nfft\n");
00075 simple_test_infft_1d(8,4,5);
00076
00077 return 1;
00078 }
00079
