test_smoothing_operator.py 2.7 KB
Newer Older
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
#
14
# Copyright(C) 2013-2019 Max-Planck-Society
Theo Steininger's avatar
Theo Steininger committed
15
#
16
# NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik.
17

18
import numpy as np
Philipp Arras's avatar
Philipp Arras committed
19
import pytest
20
21
from numpy.testing import assert_allclose

Martin Reinecke's avatar
5->6    
Martin Reinecke committed
22
import nifty6 as ift
Philipp Arras's avatar
Philipp Arras committed
23
24
25

from ..common import list2fixture

26
27
28
29
30
31
32

def _get_rtol(tp):
    if (tp == np.float64) or (tp == np.complex128):
        return 1e-10
    else:
        return 1e-5

Martin Reinecke's avatar
Martin Reinecke committed
33

Philipp Arras's avatar
Philipp Arras committed
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
space = list2fixture([ift.RGSpace(128)])
sigma = list2fixture([0., .5, 5.])
tp = list2fixture([np.float64, np.complex128])
pmp = pytest.mark.parametrize


def test_property(space, sigma):
    op = ift.HarmonicSmoothingOperator(space, sigma=sigma)
    if op.domain[0] != space:
        raise TypeError


def test_adjoint_times(space, sigma):
    op = ift.HarmonicSmoothingOperator(space, sigma=sigma)
    rand1 = ift.Field.from_random('normal', domain=space)
    rand2 = ift.Field.from_random('normal', domain=space)
Martin Reinecke's avatar
Martin Reinecke committed
50
51
    tt1 = rand1.s_vdot(op.times(rand2))
    tt2 = rand2.s_vdot(op.adjoint_times(rand1))
Philipp Arras's avatar
Philipp Arras committed
52
53
54
55
56
57
58
    assert_allclose(tt1, tt2)


def test_times(space, sigma):
    op = ift.HarmonicSmoothingOperator(space, sigma=sigma)
    fld = np.zeros(space.shape, dtype=np.float64)
    fld[0] = 1.
Martin Reinecke's avatar
Martin Reinecke committed
59
    rand1 = ift.Field.from_raw(space, fld)
Philipp Arras's avatar
Philipp Arras committed
60
    tt1 = op.times(rand1)
Martin Reinecke's avatar
Martin Reinecke committed
61
    assert_allclose(1, tt1.s_sum())
Philipp Arras's avatar
Philipp Arras committed
62
63
64
65
66
67
68
69
70
71
72


@pmp('sz', [128, 256])
@pmp('d', [1, 0.4])
def test_smooth_regular1(sz, d, sigma, tp):
    tol = _get_rtol(tp)
    sp = ift.RGSpace(sz, distances=d)
    smo = ift.HarmonicSmoothingOperator(sp, sigma=sigma)
    inp = ift.Field.from_random(
        domain=sp, random_type='normal', std=1, mean=4, dtype=tp)
    out = smo(inp)
Martin Reinecke's avatar
Martin Reinecke committed
73
    assert_allclose(inp.s_sum(), out.s_sum(), rtol=tol, atol=tol)
Philipp Arras's avatar
Philipp Arras committed
74
75
76
77
78
79
80
81
82
83
84
85
86


@pmp('sz1', [10, 15])
@pmp('sz2', [7, 10])
@pmp('d1', [1, 0.4])
@pmp('d2', [2, 0.3])
def test_smooth_regular2(sz1, sz2, d1, d2, sigma, tp):
    tol = _get_rtol(tp)
    sp = ift.RGSpace([sz1, sz2], distances=[d1, d2])
    smo = ift.HarmonicSmoothingOperator(sp, sigma=sigma)
    inp = ift.Field.from_random(
        domain=sp, random_type='normal', std=1, mean=4, dtype=tp)
    out = smo(inp)
Martin Reinecke's avatar
Martin Reinecke committed
87
    assert_allclose(inp.s_sum(), out.s_sum(), rtol=tol, atol=tol)