1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
|
PROGRAM test
use, intrinsic :: iso_c_binding
IMPLICIT none
INCLUDE 'fftw3.f03'
integer, parameter :: npoints = 1024, ncomplex = npoints / 2 + 1
real(C_DOUBLE), parameter :: start = 0, end = 1
integer :: i
real(C_DOUBLE) :: t, &
delta = (end - start) / real(npoints - 1, C_DOUBLE)
type(C_PTR) :: p_real, p_complex, plan
real(C_DOUBLE), pointer :: arr_real(:)
complex(C_DOUBLE_COMPLEX), pointer :: arr_complex(:)
p_real = fftw_alloc_real(int(npoints, C_SIZE_T))
p_complex = fftw_alloc_complex(int(ncomplex, C_SIZE_T))
call c_f_pointer(p_real, arr_real, [npoints])
call c_f_pointer(p_complex, arr_complex, [ncomplex])
plan = fftw_plan_dft_r2c_1d(int(npoints, C_INT), arr_real, &
arr_complex, FFTW_ESTIMATE)
DO i = 1, npoints
t = (i - 1) * delta + start
arr_real(i) = f1(t)
! write(*,*) t, " ", arr_real(i)
END DO
call fftw_execute_dft_r2c(plan, arr_real, arr_complex)
DO i = 1, ncomplex
write(*,*) abs(arr_complex(i))
END DO
call fftw_free(p_real)
call fftw_free(p_complex)
CONTAINS
PURE real(C_DOUBLE) FUNCTION f1(t)
real(C_DOUBLE), intent(in) :: t
real(C_DOUBLE), parameter :: pi = acos(-1.0)
f1 = sin(2 * pi * t * 200) + 2 * sin(2 * pi * t * 400)
END FUNCTION f1
END PROGRAM test
|
