From 841b7ecbd6c8406fc95b5c3e41ca9b7bd436788a Mon Sep 17 00:00:00 2001
From: Theo Steininger <theo.steininger@ultimanet.de>
Date: Sun, 9 Jul 2017 01:18:22 +0200
Subject: [PATCH] Fixed bug in _spec_to_rescaler and added stochastic tests for
power_analyze/power_synthesize
---
nifty/field.py | 7 +++++--
test/test_field.py | 43 +++++++++++++++++++++++++++++++++++++++++--
2 files changed, 46 insertions(+), 4 deletions(-)
diff --git a/nifty/field.py b/nifty/field.py
index 5b9355f43..961ec7fa2 100644
--- a/nifty/field.py
+++ b/nifty/field.py
@@ -672,8 +672,11 @@ class Field(Loggable, Versionable, object):
# Do this for every 'pindex-slice' in parallel using the 'slice(None)'s
local_pindex = pindex.get_local_data(copy=False)
- local_blow_up = [slice(None)]*len(self.shape)
- local_blow_up[self.domain_axes[power_space_index][0]] = local_pindex
+ local_blow_up = [slice(None)]*len(spec.shape)
+ # it is important to count from behind, since spec potentially grows
+ # with every iteration
+ index = self.domain_axes[power_space_index][0]-len(self.shape)
+ local_blow_up[index] = local_pindex
# here, the power_spectrum is distributed into the new shape
local_rescaler = spec[local_blow_up]
return local_rescaler
diff --git a/test/test_field.py b/test/test_field.py
index 6daaf99df..98249ec2e 100644
--- a/test/test_field.py
+++ b/test/test_field.py
@@ -20,12 +20,15 @@ import unittest
import numpy as np
from numpy.testing import assert_,\
- assert_almost_equal
+ assert_almost_equal,\
+ assert_allclose
from itertools import product
from nifty import Field,\
- RGSpace
+ RGSpace,\
+ LMSpace,\
+ PowerSpace
from d2o import distributed_data_object
@@ -80,3 +83,39 @@ class Test_Functionality(unittest.TestCase):
assert_almost_equal(a1.get_full_data(), a2.get_full_data())
assert_almost_equal(h1.get_full_data(), h3.get_full_data())
assert_almost_equal(a1.get_full_data(), a3.get_full_data())
+
+ @expand(product([RGSpace((8,), harmonic=True,
+ zerocenter=False),
+ RGSpace((8, 8), harmonic=True, distances=0.123,
+ zerocenter=True)],
+ [RGSpace((8,), harmonic=True,
+ zerocenter=False),
+ LMSpace(12)]))
+ def test_power_synthesize_analyze(self, space1, space2):
+ p1 = PowerSpace(space1)
+ spec1 = lambda k: 42/(1+k)**2
+ fp1 = Field(p1, val=spec1)
+
+ p2 = PowerSpace(space2)
+ spec2 = lambda k: 42/(1+k)**3
+ fp2 = Field(p2, val=spec2)
+
+ outer = np.outer(fp1.val.get_full_data(), fp2.val.get_full_data())
+ fp = Field((p1, p2), val=outer)
+
+ samples = 1000
+ ps1 = 0.
+ ps2 = 0.
+ for ii in xrange(samples):
+ sk = fp.power_synthesize(spaces=(0, 1), real_signal=True)
+
+ sp = sk.power_analyze(spaces=(0, 1), keep_phase_information=False)
+ ps1 += sp.sum(spaces=1)/fp2.sum()
+ ps2 += sp.sum(spaces=0)/fp1.sum()
+
+ assert_allclose(ps1.val.get_full_data()/samples,
+ fp1.val.get_full_data(),
+ rtol=0.1)
+ assert_allclose(ps2.val.get_full_data()/samples,
+ fp2.val.get_full_data(),
+ rtol=0.1)
--
GitLab