NFFT 3.5.3alpha
nfsft/simple_test.c
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
19/* standard headers */
20#include <stdio.h>
21#include <math.h>
22#include <string.h>
23#include <stdlib.h>
24/* It is important to include complex.h before nfft3.h. */
25#include <complex.h>
26
27#include "nfft3.h" /* NFFT3 header */
28
29#define __FES__ "E"
30#define K(x) ((double) x)
31
32static void simple_test_nfsft(void)
33{
34 const int N = 4; /* bandwidth/maximum degree */
35 const int M = 8; /* number of nodes */
36 nfsft_plan plan; /* transform plan */
37 int j, k, n; /* loop variables */
38
39 /* precomputation (for fast polynomial transform) */
40 nfsft_precompute(N,1000.0,0U,0U);
41
42 /* Initialize transform plan using the guru interface. All input and output
43 * arrays are allocated by nfsft_init_guru(). Computations are performed with
44 * respect to L^2-normalized spherical harmonics Y_k^n. The array of spherical
45 * Fourier coefficients is preserved during transformations. The NFFT uses a
46 * cut-off parameter m = 6. See the NFFT 3 manual for details.
47 */
48 nfsft_init_guru(&plan, N, M, NFSFT_MALLOC_X | NFSFT_MALLOC_F |
51
52 /* pseudo-random nodes */
53 for (j = 0; j < plan.M_total; j++)
54 {
55 plan.x[2*j]= nfft_drand48() - K(0.5);
56 plan.x[2*j+1]= K(0.5) * nfft_drand48();
57 }
58
59 /* precomputation (for NFFT, node-dependent) */
60 nfsft_precompute_x(&plan);
61
62 /* pseudo-random Fourier coefficients */
63 for (k = 0; k <= plan.N; k++)
64 for (n = -k; n <= k; n++)
65 plan.f_hat[NFSFT_INDEX(k,n,&plan)] =
66 nfft_drand48() - K(0.5) + _Complex_I*(nfft_drand48() - K(0.5));
67
68 /* Direct transformation, display result. */
69 nfsft_trafo_direct(&plan);
70 printf("Vector f (NDSFT):\n");
71 for (j = 0; j < plan.M_total; j++)
72 printf("f[%+2d] = %+5.3" __FES__ " %+5.3" __FES__ "*I\n",j,
73 creal(plan.f[j]), cimag(plan.f[j]));
74
75 printf("\n");
76
77 /* Fast approximate transformation, display result. */
78 nfsft_trafo(&plan);
79 printf("Vector f (NFSFT):\n");
80 for (j = 0; j < plan.M_total; j++)
81 printf("f[%+2d] = %+5.3" __FES__ " %+5.3" __FES__ "*I\n",j,
82 creal(plan.f[j]), cimag(plan.f[j]));
83
84 printf("\n");
85
86 /* Direct adjoint transformation, display result. */
87 nfsft_adjoint_direct(&plan);
88 printf("Vector f_hat (NDSFT):\n");
89 for (k = 0; k <= plan.N; k++)
90 for (n = -k; n <= k; n++)
91 fprintf(stdout,"f_hat[%+2d,%+2d] = %+5.3" __FES__ " %+5.3" __FES__ "*I\n",k,n,
92 creal(plan.f_hat[NFSFT_INDEX(k,n,&plan)]),
93 cimag(plan.f_hat[NFSFT_INDEX(k,n,&plan)]));
94
95 printf("\n");
96
97 /* Fast approximate adjoint transformation, display result. */
98 nfsft_adjoint(&plan);
99 printf("Vector f_hat (NFSFT):\n");
100 for (k = 0; k <= plan.N; k++)
101 {
102 for (n = -k; n <= k; n++)
103 {
104 fprintf(stdout,"f_hat[%+2d,%+2d] = %+5.3" __FES__ " %+5.3" __FES__ "*I\n",k,n,
105 creal(plan.f_hat[NFSFT_INDEX(k,n,&plan)]),
106 cimag(plan.f_hat[NFSFT_INDEX(k,n,&plan)]));
107 }
108 }
109
110 /* Finalize the plan. */
111 nfsft_finalize(&plan);
112
113 /* Destroy data precomputed for fast polynomial transform. */
114 nfsft_forget();
115}
116
117int main(void)
118{
119 printf("Computing an NDSFT, an NFSFT, an adjoint NDSFT, and an adjoint NFSFT"
120 "...\n\n");
121 simple_test_nfsft();
122 return EXIT_SUCCESS;
123}
#define FFT_OUT_OF_PLACE
Definition nfft3.h:196
#define PRE_PSI
Definition nfft3.h:191
#define FFTW_INIT
Definition nfft3.h:197
#define PRE_PHI_HUT
Definition nfft3.h:187
#define NFSFT_MALLOC_X
Definition nfft3.h:588
void nfsft_forget(void)
Definition nfsft.c:579
void nfsft_trafo(nfsft_plan *plan)
Definition nfsft.c:1068
#define NFSFT_NORMALIZED
Definition nfft3.h:585
void nfsft_adjoint(nfsft_plan *plan)
Definition nfsft.c:1316
#define NFSFT_INDEX(k, n, plan)
Definition nfft3.h:605
void nfsft_finalize(nfsft_plan *plan)
Definition nfsft.c:625
#define NFSFT_MALLOC_F_HAT
Definition nfft3.h:589
void nfsft_precompute(int N, double kappa, unsigned int nfsft_flags, unsigned int fpt_flags)
Definition nfsft.c:376
#define NFSFT_PRESERVE_F_HAT
Definition nfft3.h:591
#define NFSFT_MALLOC_F
Definition nfft3.h:590
Header file for the nfft3 library.
data structure for an NFSFT (nonequispaced fast spherical Fourier transform) plan with double precisi...
Definition nfft3.h:581
NFFT_INT M_total
Total number of samples.
Definition nfft3.h:581
int N
the bandwidth
Definition nfft3.h:581
fftw_complex * f
Samples.
Definition nfft3.h:581
double * x
the nodes for ,
Definition nfft3.h:581
fftw_complex * f_hat
Fourier coefficients.
Definition nfft3.h:581