00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #include "config.h"
00022
00023
00024 #include <stdio.h>
00025 #include <math.h>
00026 #include <string.h>
00027 #include <stdlib.h>
00028
00029 #ifdef HAVE_COMPLEX_H
00030 #include <complex.h>
00031 #endif
00032
00033 #include "nfft3.h"
00034 #include "nfft3util.h"
00035 #include "infft.h"
00036
00037 static void simple_test_nfsft(void)
00038 {
00039 const int N = 4;
00040 const int M = 8;
00041 nfsft_plan plan;
00042 int j, k, n;
00043
00044
00045 nfsft_precompute(N,1000.0,0U,0U);
00046
00047
00048
00049
00050
00051
00052
00053 nfsft_init_guru(&plan, N, M, NFSFT_MALLOC_X | NFSFT_MALLOC_F |
00054 NFSFT_MALLOC_F_HAT | NFSFT_NORMALIZED | NFSFT_PRESERVE_F_HAT,
00055 PRE_PHI_HUT | PRE_PSI | FFTW_INIT | FFT_OUT_OF_PLACE, 6);
00056
00057
00058 for (j = 0; j < plan.M_total; j++)
00059 {
00060 plan.x[2*j]= nfft_drand48() - K(0.5);
00061 plan.x[2*j+1]= K(0.5) * nfft_drand48();
00062 }
00063
00064
00065 nfsft_precompute_x(&plan);
00066
00067
00068 for (k = 0; k <= plan.N; k++)
00069 for (n = -k; n <= k; n++)
00070 plan.f_hat[NFSFT_INDEX(k,n,&plan)] =
00071 nfft_drand48() - K(0.5) + _Complex_I*(nfft_drand48() - K(0.5));
00072
00073
00074 nfsft_trafo_direct(&plan);
00075 printf("Vector f (NDSFT):\n");
00076 for (j = 0; j < plan.M_total; j++)
00077 printf("f[%+2d] = %+5.3" FE " %+5.3" FE "*I\n",j,
00078 creal(plan.f[j]), cimag(plan.f[j]));
00079
00080 printf("\n");
00081
00082
00083 printf("Vector f (NDSFT):\n");
00084 for (j = 0; j < plan.M_total; j++)
00085 printf("f[%+2d] = %+5.3" FE " %+5.3" FE "*I\n",j,
00086 creal(plan.f[j]), cimag(plan.f[j]));
00087
00088 printf("\n");
00089
00090
00091 nfsft_adjoint_direct(&plan);
00092 printf("Vector f_hat (NDSFT):\n");
00093 for (k = 0; k <= plan.N; k++)
00094 for (n = -k; n <= k; n++)
00095 fprintf(stdout,"f_hat[%+2d,%+2d] = %+5.3" FE " %+5.3" FE "*I\n",k,n,
00096 creal(plan.f_hat[NFSFT_INDEX(k,n,&plan)]),
00097 cimag(plan.f_hat[NFSFT_INDEX(k,n,&plan)]));
00098
00099 printf("\n");
00100
00101
00102 nfsft_adjoint(&plan);
00103 printf("Vector f_hat (NFSFT):\n");
00104 for (k = 0; k <= plan.N; k++)
00105 {
00106 for (n = -k; n <= k; n++)
00107 {
00108 fprintf(stdout,"f_hat[%+2d,%+2d] = %+5.3" FE " %+5.3" FE "*I\n",k,n,
00109 creal(plan.f_hat[NFSFT_INDEX(k,n,&plan)]),
00110 cimag(plan.f_hat[NFSFT_INDEX(k,n,&plan)]));
00111 }
00112 }
00113
00114
00115 nfsft_finalize(&plan);
00116
00117
00118 nfsft_forget();
00119 }
00120
00121 int main(void)
00122 {
00123 printf("Computing an NDSFT, an NFSFT, an adjoint NDSFT, and an adjoint NFSFT"
00124 "...\n\n");
00125 simple_test_nfsft();
00126 return EXIT_SUCCESS;
00127 }