Commit 7453ea72 authored by Philipp Haim's avatar Philipp Haim

Again merge branch 'normalized_amplitudes_pp' into normalized_amplitudes_DomT

parents 94309404 5e871700
Pipeline #63674 failed with stages
in 5 minutes and 59 seconds
......@@ -28,16 +28,16 @@ def _default_pspace(dom):
if __name__ == '__main__':
np.random.seed(42)
fa = ift.CorrelatedFieldMaker()
fa = ift.CorrelatedFieldMaker.make(10, 0.1, '')
n_samps = 20
slope_means = [-2, -3]
fa.add_fluctuations(ift.RGSpace(128, 0.1), 10, 2, 1, 1e-6,
2, 1e-6, slope_means[0], 0.2, 'spatial')
# fa.add_fluctuations(_default_pspace(ift.RGSpace((128, 64))), 10, 2, 1,
# 1e-6, 2, 1e-6, slope_means[0], 0.2, 'spatial')
fa.add_fluctuations(ift.RGSpace(32), 10, 5, 1, 1e-6, 2,
fa.add_fluctuations(ift.RGSpace(32), 3, 5, 1, 1e-6, 2,
1e-6, slope_means[1], 1, 'freq')
correlated_field = fa.finalize(10, 0.1, '')
correlated_field = fa.finalize()
amplitudes = fa.amplitudes
plt.style.use('seaborn-notebook')
......
......@@ -57,10 +57,10 @@ if __name__ == '__main__':
position_space = ift.RGSpace([128, 128])
cfmaker = ift.CorrelatedFieldMaker()
cfmaker = ift.CorrelatedFieldMaker.make(1e-3, 1e-6, '')
cfmaker.add_fluctuations(position_space,
1, 1e-2, 1, .5, .1, .5, -3, 0.5, '')
correlated_field = cfmaker.finalize(1e-3, 1e-6, '')
1., 1e-2, 1, .5, .1, .5, -3, 0.5, '')
correlated_field = cfmaker.finalize()
A = cfmaker.amplitudes[0]
# Apply a nonlinearity
......
......@@ -56,7 +56,7 @@ def radial_los(n_los):
if __name__ == '__main__':
np.random.seed(45)
np.random.seed(43)
# Choose between random line-of-sight response (mode=0) and radial lines
# of sight (mode=1)
......@@ -71,12 +71,12 @@ if __name__ == '__main__':
sp1 = ift.RGSpace(npix1)
sp2 = ift.RGSpace(npix2)
cfmaker = ift.CorrelatedFieldMaker()
cfmaker = ift.CorrelatedFieldMaker.make(1e-2, 1e-6, '')
amp1 = 0.5
cfmaker.add_fluctuations(sp1, amp1, 1e-2, 1, .1, .01, .5, -2, 1., 'amp1')
cfmaker.add_fluctuations(sp2, np.sqrt(1. - amp1**2), 1e-2, 1, .1, .01, .5,
cfmaker.add_fluctuations(sp1, 0.1, 1e-2, 1, .1, .01, .5, -2, 1., 'amp1')
cfmaker.add_fluctuations(sp2, 0.1, 1e-2, 1, .1, .01, .5,
-1.5, .5, 'amp2')
correlated_field = cfmaker.finalize(1e-3, 1e-6, '')
correlated_field = cfmaker.finalize()
A1 = cfmaker.amplitudes[0]
A2 = cfmaker.amplitudes[1]
......
......@@ -11,7 +11,7 @@ def testAmplitudesConsistency(seed, sspace):
sc.add(op(s.extract(op.domain)))
return sc.mean.to_global_data(), sc.var.sqrt().to_global_data()
np.random.seed(seed)
offset_std = 30
offset_std = .1
intergated_fluct_std0 = .003
intergated_fluct_std1 = 0.1
......@@ -20,12 +20,12 @@ def testAmplitudesConsistency(seed, sspace):
fsspace = ift.RGSpace((12,), (0.4,))
fa = ift.CorrelatedFieldMaker()
fa = ift.CorrelatedFieldMaker.make(offset_std, 1E-8, '')
fa.add_fluctuations(sspace, intergated_fluct_std0, 1E-8, 1.1, 2., 2.1, .5,
-2, 1., 'spatial')
fa.add_fluctuations(fsspace, intergated_fluct_std1, 1E-8, 3.1, 1., .5, .1,
-4, 1., 'freq')
op = fa.finalize(offset_std, 1E-8, '')
op = fa.finalize()
samples = [ift.from_random('normal',op.domain) for _ in range(nsam)]
tot_flm, _ = stats(fa.total_fluctuation,samples)
......@@ -44,13 +44,6 @@ def testAmplitudesConsistency(seed, sspace):
sl_fluct_space = fa.slice_fluctuation_realized(sams,0)
sl_fluct_freq = fa.slice_fluctuation_realized(sams,1)
np.testing.assert_allclose(offset_std, zm_std_mean, rtol=0.5)
np.testing.assert_allclose(intergated_fluct_std0, fluct_space, rtol=0.5)
np.testing.assert_allclose(intergated_fluct_std1, fluct_freq, rtol=0.5)
np.testing.assert_allclose(tot_flm, fluct_total, rtol=0.5)
np.testing.assert_allclose(slice_fluct_std0, sl_fluct_space, rtol=0.5)
np.testing.assert_allclose(slice_fluct_std1, sl_fluct_freq, rtol=0.5)
print("Expected offset Std: " + str(offset_std))
print("Estimated offset Std: " + str(zm_std_mean))
......@@ -73,14 +66,23 @@ def testAmplitudesConsistency(seed, sspace):
print("Expected total fluct. Std: " + str(tot_flm))
print("Estimated total fluct. Std: " + str(fluct_total))
fa = ift.CorrelatedFieldMaker()
np.testing.assert_allclose(offset_std, zm_std_mean, rtol=0.5)
np.testing.assert_allclose(intergated_fluct_std0, fluct_space, rtol=0.5)
np.testing.assert_allclose(intergated_fluct_std1, fluct_freq, rtol=0.5)
np.testing.assert_allclose(tot_flm, fluct_total, rtol=0.5)
np.testing.assert_allclose(slice_fluct_std0, sl_fluct_space, rtol=0.5)
np.testing.assert_allclose(slice_fluct_std1, sl_fluct_freq, rtol=0.5)
fa = ift.CorrelatedFieldMaker.make(offset_std, .1, '')
fa.add_fluctuations(fsspace, intergated_fluct_std1, 1., 3.1, 1., .5, .1,
-4, 1., 'freq')
m = 3.
x = fa.moment_slice_to_average(m)
fa.add_fluctuations(sspace, x, 1.5, 1.1, 2., 2.1, .5,
-2, 1., 'spatial', 0)
op = fa.finalize(offset_std, .1, '')
op = fa.finalize()
em, estd = stats(fa.slice_fluctuation(0),samples)
print("Forced slice fluct. space Std: "+str(m))
print("Expected slice fluct. Std: " + str(em))
......
......@@ -4,9 +4,9 @@ np.random.seed(42)
sspace = ift.RGSpace((128,))
fa = ift.CorrelatedFieldMaker()
fa = ift.CorrelatedFieldMaker.make(10, 0.1, '')
fa.add_fluctuations(sspace, 10, 2, 1, 1e-6, 2, 1e-6, -2, 1e-6, 'spatial')
op = fa.finalize(10, 0.1, '')
op = fa.finalize()
A = fa.amplitudes[0]
cstpos = ift.from_random('normal', op.domain)
......
......@@ -317,11 +317,24 @@ class _Amplitude(Operator):
class CorrelatedFieldMaker:
def __init__(self):
def __init__(self, amplitude_offset, prefix):
self._a = []
self._spaces = []
self._azm = None
self._position_spaces = []
self._azm = amplitude_offset
self._prefix = prefix
@staticmethod
def make(offset_amplitude_mean, offset_amplitude_stddev, prefix):
offset_amplitude_stddev = float(offset_amplitude_stddev)
offset_amplitude_mean = float(offset_amplitude_mean)
assert offset_amplitude_stddev > 0
assert offset_amplitude_mean > 0
zm = _LognormalMomentMatching(offset_amplitude_mean,
offset_amplitude_stddev,
prefix + 'zeromode')
return CorrelatedFieldMaker(zm, prefix)
def add_fluctuations(self,
position_space,
......@@ -342,6 +355,8 @@ class CorrelatedFieldMaker:
power_space = makeDomain(power_space)
prefix = str(prefix)
#assert isinstance(position_space[space], (RGSpace, HPSpace, GLSpace)
#NOTE alternative to get auxilliary domain
#auxdom = ContractionOperator(position_space, space).domain
auxdom = makeDomain(tuple(dom for i, dom in enumerate(position_space)
if i != space))
......@@ -349,6 +364,8 @@ class CorrelatedFieldMaker:
fluctuations_stddev,
prefix + 'fluctuations',
auxdom)
#FIXME How should this work on domain tuples?
#fluct = fluct*self._azm.one_over()
flex = _LognormalMomentMatching(flexibility_mean, flexibility_stddev,
prefix + 'flexibility',
auxdom)
......@@ -381,7 +398,6 @@ class CorrelatedFieldMaker:
zeroind += (0,)*len(dd[self._spaces[i]].shape)
zeroind += (slice(None),)*(len(dd)-self._spaces[i]-1)
#tuple(zeroind = zeroind + (slice(None),)*space + len(dd.shape)*(0,)
foo = np.ones(hspace.shape)
foo[zeroind] = 0
......@@ -420,27 +436,16 @@ class CorrelatedFieldMaker:
return ht(azm*(pd @ a)*ducktape(hspace, None, prefix + 'xi'))
def finalize(self,
offset_amplitude_mean,
offset_amplitude_stddev,
prefix='',
offset=None,
prior_info=100):
"""
offset vs zeromode: volume factor
"""
prefix = str(prefix)
if offset is not None:
raise NotImplementedError
offset = float(offset)
auxdom = makeDomain((dom for i, domT in enumerate(self._position_spaces)
for j, dom in enumerate(domT)
if j != self._spaces[i]))
azm = _LognormalMomentMatching(offset_amplitude_mean,
offset_amplitude_stddev,
prefix + 'zeromode',
auxdom)
op = self.finalize_from_op(azm, prefix)
op = self.finalize_from_op(self._azm, self._prefix)
if prior_info > 0:
from ..sugar import from_random
samps = [
......@@ -469,12 +474,13 @@ class CorrelatedFieldMaker:
def moment_slice_to_average(self, fluctuations_slice_mean, nsamples=1000):
fluctuations_slice_mean = float(fluctuations_slice_mean)
assert fluctuations_slice_mean > 0
from ..sugar import from_random
scm = 1.
for a in self._a:
m, std = a.fluctuation_amplitude.mean, a.fluctuation_amplitude.std
mu, sig = _lognormal_moments(m, std)
flm = np.exp(mu + sig*np.random.normal(size=nsamples))
scm *= flm**2 + 1.
op = a.fluctuation_amplitude
res= np.array([op(from_random('normal',op.domain)).to_global_data()
for _ in range(nsamples)])
scm *= res**2 + 1.
return fluctuations_slice_mean/np.mean(np.sqrt(scm))
@property
......@@ -491,12 +497,12 @@ class CorrelatedFieldMaker:
if len(self._a) == 0:
raise NotImplementedError
if len(self._a) == 1:
return self._a[0].fluctuation_amplitude
return self.average_fluctuation(0)
q = 1.
for a in self._a:
fl = a.fluctuation_amplitude
q = q*(Adder(full(fl.target, 1.)) @ fl**2)
return (Adder(full(q.target, -1.)) @ q).sqrt()
return (Adder(full(q.target, -1.)) @ q).sqrt() * self._azm
def slice_fluctuation(self, space):
"""Returns operator which acts on prior or posterior samples"""
......@@ -504,7 +510,7 @@ class CorrelatedFieldMaker:
raise NotImplementedError
assert space < len(self._a)
if len(self._a) == 1:
return self._a[0].fluctuation_amplitude
return self.average_fluctuation(0)
q = 1.
for j in range(len(self._a)):
fl = self._a[j].fluctuation_amplitude
......@@ -512,7 +518,7 @@ class CorrelatedFieldMaker:
q = q*fl**2
else:
q = q*(Adder(full(fl.target, 1.)) @ fl**2)
return q.sqrt()
return q.sqrt() * self._azm
def average_fluctuation(self, space):
"""Returns operator which acts on prior or posterior samples"""
......@@ -520,8 +526,8 @@ class CorrelatedFieldMaker:
raise NotImplementedError
assert space < len(self._a)
if len(self._a) == 1:
return self._a[0].fluctuation_amplitude
return self._a[space].fluctuation_amplitude
return self._a[0].fluctuation_amplitude*self._azm
return self._a[space].fluctuation_amplitude*self._azm
@staticmethod
def offset_amplitude_realized(samples):
......
......@@ -42,12 +42,12 @@ def testAmplitudesConsistency(rseed, sspace, Astds, offset_std):
fsspace = ift.RGSpace((12,), (0.4,))
fa = ift.CorrelatedFieldMaker()
fa = ift.CorrelatedFieldMaker.make(offset_std, 1E-8, '')
fa.add_fluctuations(sspace, Astds[0], 1E-8, 1.1, 2., 2.1, .5,
-2, 1., 'spatial')
fa.add_fluctuations(fsspace, Astds[1], 1E-8, 3.1, 1., .5, .1,
-4, 1., 'freq')
op = fa.finalize(offset_std, 1E-8, '')
op = fa.finalize()
samples = [ift.from_random('normal',op.domain) for _ in range(nsam)]
tot_flm, _ = stats(fa.total_fluctuation,samples)
......@@ -73,14 +73,14 @@ def testAmplitudesConsistency(rseed, sspace, Astds, offset_std):
assert_allclose(slice_fluct_std0, sl_fluct_space, rtol=0.5)
assert_allclose(slice_fluct_std1, sl_fluct_freq, rtol=0.5)
fa = ift.CorrelatedFieldMaker()
fa = ift.CorrelatedFieldMaker.make(offset_std, .1, '')
fa.add_fluctuations(fsspace, Astds[1], 1., 3.1, 1., .5, .1,
-4, 1., 'freq')
m = 3.
x = fa.moment_slice_to_average(m)
fa.add_fluctuations(sspace, x, 1.5, 1.1, 2., 2.1, .5,
-2, 1., 'spatial', 0)
op = fa.finalize(offset_std, .1, '')
op = fa.finalize()
em, estd = stats(fa.slice_fluctuation(0),samples)
assert_allclose(m, em, rtol=0.5)
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment