smooth_sky.py 2.53 KB
Newer Older
Jakob Knollmueller's avatar
Jakob Knollmueller committed
1
def make_correlated_field(s_space, amplitude_model):
2
    '''
Martin Reinecke's avatar
Martin Reinecke committed
3
    Method for construction of correlated fields
4
5
6

    Parameters
    ----------
Martin Reinecke's avatar
Martin Reinecke committed
7
    s_space : Field domain
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24

    amplitude_model : model for correlation structure
    '''
    from .. import (FFTOperator, Field, MultiField, PointwiseExponential,
                    PowerDistributor, Variable)
    h_space = s_space.get_default_codomain()
    ht = FFTOperator(h_space, s_space)
    p_space = amplitude_model.value.domain[0]
    power_distributor = PowerDistributor(h_space, p_space)
    position = {}
    position['xi'] = Field.from_random('normal', h_space)
    position['tau'] = amplitude_model.position['tau']
    position['phi'] = amplitude_model.position['phi']
    position = MultiField(position)

    xi = Variable(position)['xi']
    A = power_distributor(amplitude_model)
Jakob Knollmueller's avatar
Jakob Knollmueller committed
25
26
27
    correlated_field_h = A * xi
    correlated_field = ht(correlated_field_h)
    internals = {'correlated_field_h': correlated_field_h,
28
29
                 'power_distributor': power_distributor,
                 'ht': ht}
Jakob Knollmueller's avatar
Jakob Knollmueller committed
30
    return correlated_field, internals
31
32


Martin Reinecke's avatar
Martin Reinecke committed
33
def make_mf_correlated_field(s_space_spatial, s_space_energy,
34
35
                             amplitude_model_spatial, amplitude_model_energy):
    '''
Martin Reinecke's avatar
Martin Reinecke committed
36
    Method for construction of correlated multi-frequency fields
37
    '''
Philipp Arras's avatar
Fixups    
Philipp Arras committed
38
39
40
41
    from .. import (DomainTuple, Field, MultiField,
                    PointwiseExponential, Variable)
    from ..operators import (DomainDistributor, PowerDistributor,
                             HarmonicTransformOperator)
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
    h_space_spatial = s_space_spatial.get_default_codomain()
    h_space_energy = s_space_energy.get_default_codomain()
    h_space = DomainTuple.make((h_space_spatial, h_space_energy))
    ht1 = HarmonicTransformOperator(h_space, space=0)
    ht2 = HarmonicTransformOperator(ht1.target, space=1)
    ht = ht2*ht1

    p_space_spatial = amplitude_model_spatial.value.domain[0]
    p_space_energy = amplitude_model_energy.value.domain[0]

    pd_spatial = PowerDistributor(h_space, p_space_spatial, 0)
    pd_energy = PowerDistributor(pd_spatial.domain, p_space_energy, 1)
    pd = pd_spatial*pd_energy

    dom_distr_0 = DomainDistributor(pd.domain, 0)
    dom_distr_1 = DomainDistributor(pd.domain, 1)
    a_spatial = dom_distr_1(amplitude_model_spatial)
    a_energy = dom_distr_0(amplitude_model_energy)
    a = a_spatial*a_energy
    A = pd(a)

    position = MultiField({'xi': Field.from_random('normal', h_space)})
    xi = Variable(position)['xi']
    logsky_h = A*xi
    logsky = ht(logsky_h)
    return PointwiseExponential(logsky)