test_fft_operator.py 4.66 KB
Newer Older
Martin Reinecke's avatar
Martin Reinecke committed
1
2
3
4
5
6
7
8
9
10
11
12
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.
Theo Steininger's avatar
Theo Steininger committed
13
#
Martin Reinecke's avatar
Martin Reinecke committed
14
# Copyright(C) 2013-2018 Max-Planck-Society
Theo Steininger's avatar
Theo Steininger committed
15
16
17
#
# NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik
# and financially supported by the Studienstiftung des deutschen Volkes.
Martin Reinecke's avatar
Martin Reinecke committed
18

Martin Reinecke's avatar
Martin Reinecke committed
19
from __future__ import division
20

Martin Reinecke's avatar
Martin Reinecke committed
21
import unittest
Martin Reinecke's avatar
Martin Reinecke committed
22
23
from itertools import product
from test.common import expand
24

25
26
27
28
import nifty5 as ift
import numpy as np
from numpy.testing import assert_allclose

Theo Steininger's avatar
Theo Steininger committed
29

Martin Reinecke's avatar
Martin Reinecke committed
30
def _get_rtol(tp):
Theo Steininger's avatar
Theo Steininger committed
31
    if (tp == np.float64) or (tp == np.complex128):
Martin Reinecke's avatar
Martin Reinecke committed
32
33
34
        return 1e-10
    else:
        return 1e-5
Martin Reinecke's avatar
Martin Reinecke committed
35

Theo Steininger's avatar
Theo Steininger committed
36

37
class FFTOperatorTests(unittest.TestCase):
Martin Reinecke's avatar
Martin Reinecke committed
38
    @expand(product([16, ], [0.1, 1, 3.7],
Martin Reinecke's avatar
Martin Reinecke committed
39
40
                    [np.float64, np.float32, np.complex64, np.complex128]))
    def test_fft1D(self, dim1, d, itp):
Theo Steininger's avatar
Theo Steininger committed
41
        tol = _get_rtol(itp)
Martin Reinecke's avatar
Martin Reinecke committed
42
43
        a = ift.RGSpace(dim1, distances=d)
        b = ift.RGSpace(dim1, distances=1./(dim1*d), harmonic=True)
Martin Reinecke's avatar
Martin Reinecke committed
44
45
        np.random.seed(16)

Martin Reinecke's avatar
Martin Reinecke committed
46
        fft = ift.FFTOperator(domain=a, target=b)
Martin Reinecke's avatar
Martin Reinecke committed
47
48
49
        inp = ift.Field.from_random(domain=a, random_type='normal',
                                    std=7, mean=3, dtype=itp)
        out = fft.inverse_times(fft.times(inp))
Martin Reinecke's avatar
Martin Reinecke committed
50
        assert_allclose(inp.local_data, out.local_data, rtol=tol, atol=tol)
51

Martin Reinecke's avatar
Martin Reinecke committed
52
53
54
        a, b = b, a

        fft = ift.FFTOperator(domain=a, target=b)
Martin Reinecke's avatar
Martin Reinecke committed
55
56
        inp = ift.Field.from_random(domain=a, random_type='normal',
                                    std=7, mean=3, dtype=itp)
Martin Reinecke's avatar
Martin Reinecke committed
57
        out = fft.inverse_times(fft.times(inp))
Martin Reinecke's avatar
Martin Reinecke committed
58
        assert_allclose(inp.local_data, out.local_data, rtol=tol, atol=tol)
Martin Reinecke's avatar
Martin Reinecke committed
59

Martin Reinecke's avatar
Martin Reinecke committed
60
    @expand(product([12, 15], [9, 12], [0.1, 1, 3.7],
Martin Reinecke's avatar
stage1    
Martin Reinecke committed
61
                    [0.4, 1, 2.7],
Martin Reinecke's avatar
Martin Reinecke committed
62
63
                    [np.float64, np.float32, np.complex64, np.complex128]))
    def test_fft2D(self, dim1, dim2, d1, d2, itp):
Theo Steininger's avatar
Theo Steininger committed
64
        tol = _get_rtol(itp)
Martin Reinecke's avatar
Martin Reinecke committed
65
66
67
        a = ift.RGSpace([dim1, dim2], distances=[d1, d2])
        b = ift.RGSpace([dim1, dim2],
                        distances=[1./(dim1*d1), 1./(dim2*d2)], harmonic=True)
Martin Reinecke's avatar
Martin Reinecke committed
68

Martin Reinecke's avatar
Martin Reinecke committed
69
        fft = ift.FFTOperator(domain=a, target=b)
Martin Reinecke's avatar
Martin Reinecke committed
70
71
72
        inp = ift.Field.from_random(domain=a, random_type='normal',
                                    std=7, mean=3, dtype=itp)
        out = fft.inverse_times(fft.times(inp))
Martin Reinecke's avatar
Martin Reinecke committed
73
        assert_allclose(inp.local_data, out.local_data, rtol=tol, atol=tol)
74

Martin Reinecke's avatar
Martin Reinecke committed
75
        a, b = b, a
76

Martin Reinecke's avatar
Martin Reinecke committed
77
        fft = ift.FFTOperator(domain=a, target=b)
Martin Reinecke's avatar
Martin Reinecke committed
78
79
        inp = ift.Field.from_random(domain=a, random_type='normal',
                                    std=7, mean=3, dtype=itp)
Martin Reinecke's avatar
Martin Reinecke committed
80
        out = fft.inverse_times(fft.times(inp))
Martin Reinecke's avatar
Martin Reinecke committed
81
        assert_allclose(inp.local_data, out.local_data, rtol=tol, atol=tol)
82

Martin Reinecke's avatar
Martin Reinecke committed
83
    @expand(product([0, 1, 2],
Martin Reinecke's avatar
Martin Reinecke committed
84
85
                    [np.float64, np.float32, np.complex64, np.complex128]))
    def test_composed_fft(self, index, dtype):
86
        tol = _get_rtol(dtype)
Martin Reinecke's avatar
Martin Reinecke committed
87
88
89
        a = [a1, a2, a3] = [ift.RGSpace((32,)), ift.RGSpace((4, 4)),
                            ift.RGSpace((5, 6))]
        fft = ift.FFTOperator(domain=a, space=index)
90

Martin Reinecke's avatar
Martin Reinecke committed
91
92
        inp = ift.Field.from_random(domain=(a1, a2, a3), random_type='normal',
                                    std=7, mean=3, dtype=dtype)
Martin Reinecke's avatar
Martin Reinecke committed
93
        out = fft.inverse_times(fft.times(inp))
Martin Reinecke's avatar
Martin Reinecke committed
94
        assert_allclose(inp.local_data, out.local_data, rtol=tol, atol=tol)
Martin Reinecke's avatar
Martin Reinecke committed
95

Martin Reinecke's avatar
Martin Reinecke committed
96
    @expand(product([ift.RGSpace(128, distances=3.76, harmonic=True),
97
                     ift.RGSpace((15, 27), distances=(.7, .33), harmonic=True),
98
                     ift.RGSpace(73, distances=0.5643)],
Martin Reinecke's avatar
Martin Reinecke committed
99
100
101
                    [np.float64, np.float32, np.complex64, np.complex128]))
    def test_normalisation(self, space, tp):
        tol = 10 * _get_rtol(tp)
Martin Reinecke's avatar
Martin Reinecke committed
102
103
        cospace = space.get_default_codomain()
        fft = ift.FFTOperator(space, cospace)
Martin Reinecke's avatar
Martin Reinecke committed
104
105
106
        inp = ift.Field.from_random(domain=space, random_type='normal',
                                    std=1, mean=2, dtype=tp)
        out = fft.times(inp)
Martin Reinecke's avatar
Martin Reinecke committed
107
108
109
        fft2 = ift.FFTOperator(cospace, space)
        out2 = fft2.inverse_times(inp)
        zero_idx = tuple([0]*len(space.shape))
Martin Reinecke's avatar
Martin Reinecke committed
110
        assert_allclose(inp.to_global_data()[zero_idx], out.integrate(),
Martin Reinecke's avatar
Martin Reinecke committed
111
                        rtol=tol, atol=tol)
Martin Reinecke's avatar
Martin Reinecke committed
112
        assert_allclose(out.local_data, out2.local_data, rtol=tol, atol=tol)