easyinterface.py 1.87 KB
Newer Older
Philipp Arras's avatar
Philipp Arras committed
1
import nifty5 as ift
Philipp Arras's avatar
Philipp Arras committed
2
3
import numpy as np
np.random.seed(42)
Philipp Arras's avatar
Philipp Arras committed
4
5
6
7
8

sspace = ift.RGSpace((128, 128), (0.2, 0.2))
hspace = sspace.get_default_codomain()
target = ift.PowerSpace(hspace)

Philipp Arras's avatar
Work    
Philipp Arras committed
9
10
11
12
13
14
15
fa = ift.FinalAmplitude()
fa.add_fluctuations(target, 1, 0.1, 1, 0.1, 1, 0.1, -2, 0.1, 'fst')
fa.add_fluctuations(target, 1, 0.1, 1, 0.1, 1, 0.1, -2, 0.1, 'snd')
fa.finalize(1, 0.1, '', offset=-10)

exit()

Philipp Arras's avatar
Philipp Arras committed
16
A = ift.NormalizedAmplitude(target, 16, 1, 1, -3, 1, 0, 1, 0, 1)
Philipp Arras's avatar
Philipp Arras committed
17
18
A = ift.WPAmplitude(target, [0, -2, 0], [1E-5, 1, 1], 1, 0.99,
                    ['rest', 'smooth', 'wienersigma'])
Philipp Arras's avatar
Philipp Arras committed
19
20
21
22
23
24
25
26
27

avgA = ift.full(A.target, 0.)
n = 1000
for _ in range(n):
    avgA = avgA + A(ift.from_random('normal', A.domain))
avgA = avgA/n
corfldfixA = ift.CorrelatedField(sspace, avgA)
corfld = ift.CorrelatedField(sspace, A)

Philipp Arras's avatar
Philipp Arras committed
28
29
30
31
32
corfld = ift.CorrelatedFieldNormAmplitude(sspace, A, 0, 1)
corfldfixA = ift.CorrelatedFieldNormAmplitude(sspace, avgA, 0, 1)

cstpos = ift.from_random('normal', corfld.domain)
p, p1, p2 = [ift.Plot() for _ in range(3)]
Philipp Arras's avatar
Philipp Arras committed
33
lst, lst1 = [avgA**2], []
Philipp Arras's avatar
Philipp Arras committed
34
skys, skys1, skys2 = [], [], []
Philipp Arras's avatar
Philipp Arras committed
35
36
37
38
for _ in range(8):
    pos = ift.from_random('normal', corfld.domain)

    skyfixA = corfldfixA.force(pos)
Philipp Arras's avatar
Philipp Arras committed
39
    skys.append(skyfixA)
Philipp Arras's avatar
Philipp Arras committed
40
    ft = ift.HartleyOperator(hspace, sspace).scale(hspace.scalar_dvol**-0.5)
Philipp Arras's avatar
Philipp Arras committed
41
42
    lst.append(ift.power_analyze(ft.inverse(skyfixA)))

Philipp Arras's avatar
Philipp Arras committed
43
    foo = ift.MultiField.union([cstpos, pos.extract(A.domain)])
Philipp Arras's avatar
Philipp Arras committed
44
    skys2.append(corfld(foo))
Philipp Arras's avatar
Philipp Arras committed
45

Philipp Arras's avatar
Philipp Arras committed
46
    sky = corfld(pos)
Philipp Arras's avatar
Philipp Arras committed
47
    skys1.append(sky)
Philipp Arras's avatar
Philipp Arras committed
48
    lst1.append(A.force(pos))
Philipp Arras's avatar
Philipp Arras committed
49
50
51
52
53
54
55
56
57
58
59

for pp, ll in [(p, skys), (p1, skys1), (p2, skys2)]:
    mi, ma = None, None
    if True:
        mi, ma = np.inf, -np.inf
        for ss in ll:
            mi = min([mi, np.amin(ss.val)])
            ma = max([ma, np.amax(ss.val)])
    for ss in ll:
        pp.add(ss, zmin=mi, zmax=ma)

Philipp Arras's avatar
Philipp Arras committed
60
61
p.add(lst)
p1.add(lst1)
Philipp Arras's avatar
Philipp Arras committed
62
63
64
65
p2.add(lst1)
p1.output(name='full.png')
p.output(name='A_fixed.png')
p2.output(name='xi_fixed.png')