smooth_sky.py 2.61 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
def make_smooth_sky_model(s_space, amplitude_model):
    '''
    Method for construction of correlated sky model

    Parameters
    ----------
    s_space : domain of sky model

    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)
    logsky_h = A * xi
    logsky = ht(logsky_h)
    internals = {'logsky_h': logsky_h,
                 'power_distributor': power_distributor,
                 'ht': ht}
    return PointwiseExponential(logsky), internals


def make_smooth_mf_sky_model(s_space_spatial, s_space_energy,
                             amplitude_model_spatial, amplitude_model_energy):
    '''
    Method for construction of correlated sky model

    Parameters
    ----------
    s_space : domain of sky model

    amplitude_model : model for correlation structure
    '''
Philipp Arras's avatar
Fixups    
Philipp Arras committed
44
45
46
47
    from .. import (DomainTuple, Field, MultiField,
                    PointwiseExponential, Variable)
    from ..operators import (DomainDistributor, PowerDistributor,
                             HarmonicTransformOperator)
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
    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)