NFFT  3.3.0
nfsft_benchomp_createdataset.c
1 /*
2  * Copyright (c) 2002, 2015 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 /* $Id: simple_test.c 3372 2009-10-21 06:04:05Z skunis $ */
20 
21 #include <stdio.h>
22 #include <math.h>
23 #include <string.h>
24 #include <stdlib.h>
25 #include <complex.h>
26 
27 #include "config.h"
28 
29 #include "nfft3.h"
30 #include "infft.h"
31 
32 void nfsft_benchomp_createdataset(unsigned int trafo_adjoint, int N, int M)
33 {
34  int t, j, k, n;
35  R *x;
36  C *f, *f_hat;
37  int N_total = (2*N+2) * (2*N+2);
38  nfsft_plan ptemp;
39 
40  nfsft_init_guru(&ptemp, N, M, NFSFT_MALLOC_X | NFSFT_MALLOC_F |
41  NFSFT_MALLOC_F_HAT | NFSFT_NORMALIZED | NFSFT_PRESERVE_F_HAT,
42  PRE_PHI_HUT | PRE_PSI | FFTW_INIT | FFT_OUT_OF_PLACE, 6);
43 
44  x = (R*) nfft_malloc(2*M*sizeof(R));
45  f = (C*) nfft_malloc(M*sizeof(C));
46  f_hat = (C*) nfft_malloc(N_total*sizeof(C));
47 
48  /* init pseudo-random nodes */
49  for (j = 0; j < M; j++)
50  {
51  x[2*j]= X(drand48)() - K(0.5);
52  x[2*j+1]= K(0.5) * X(drand48)();
53  }
54 
55  if (trafo_adjoint==0)
56  {
57  for (k = 0; k <= N; k++)
58  for (n = -k; n <= k; n++)
59  nfft_vrand_unit_complex(f_hat+NFSFT_INDEX(k,n,&ptemp),1);
60  }
61  else
62  {
64  }
65 
66  printf("%d %d %d\n", trafo_adjoint, N, M);
67 
68  for (j=0; j < M; j++)
69  {
70  for (t=0; t < 2; t++)
71  printf("%.16e ", x[2*j+t]);
72  printf("\n");
73  }
74 
75  if (trafo_adjoint==0)
76  {
77  for (k = 0; k <= N; k++)
78  for (n = -k; n <= k; n++)
79  printf("%.16e %.16e\n", creal(f_hat[NFSFT_INDEX(k,n,&ptemp)]), cimag(f_hat[NFSFT_INDEX(k,n,&ptemp)]));
80  }
81  else
82  {
83  for (j=0; j < M; j++)
84  printf("%.16e %.16e\n", creal(f[j]), cimag(f[j]));
85  }
86 
87  nfft_free(x);
88  nfft_free(f);
89  nfft_free(f_hat);
90 }
91 
92 int main(int argc, char **argv)
93 {
94  int trafo_adjoint;
95  int N;
96  int M;
97 
98  if (argc < 4) {
99  fprintf(stderr, "usage: tr_adj N M\n");
100  return -1;
101  }
102 
103  trafo_adjoint = atoi(argv[1]);
104  if (trafo_adjoint < 0 && trafo_adjoint > 1)
105  trafo_adjoint = 1;
106 
107  N = atoi(argv[2]);
108  M = atoi(argv[3]);
109  fprintf(stderr, "tr_adj=%d, N=%d, M=%d\n", trafo_adjoint, N, M);
110 
111  nfsft_benchomp_createdataset(trafo_adjoint, N, M);
112 
113  return 0;
114 }
data structure for an NFSFT (nonequispaced fast spherical Fourier transform) plan with double precisi...
Definition: nfft3.h:576
void nfft_free(void *p)
void nfft_vrand_unit_complex(fftw_complex *x, const NFFT_INT n)
Inits a vector of random complex numbers in .
#define X(name)
Include header for C99 complex datatype.
Definition: fastsum.h:53
void * nfft_malloc(size_t n)