Skip to content
Snippets Groups Projects
Commit 255f0ce5 authored by Philipp Frank's avatar Philipp Frank
Browse files

add intro scrips

parent b4590ea7
Branches
No related tags found
1 merge request!3Update to latest nifty and add intro files
import nifty8 as ift
import numpy as np
# Domains / spaces
sp = ift.RGSpace(shape=(2, 3),
distances=(1., 1.5))
print(sp.shape)
print(sp.distances)
print(sp.total_volume)
print(sp.scalar_dvol)
print()
sp1 = ift.UnstructuredDomain(4)
print(sp1)
print()
# Fields
d_arr = np.array([12.1, 4.3, 21., 1110.])
d = ift.makeField(sp1, d_arr)
print(d)
s = ift.from_random(sp) # draws std normal samples
print(s)
print(s.integrate())
print(s.sum())
print()
print(s.domain)
print(s.val)
#print(s.val_rw())
# Product domains
dom = ift.makeDomain((sp, sp1))
print(dom)
print(dom.shape)
print(dom.size)
# MultiDomain: dictionaries of Domains
dom = ift.makeDomain({"key0": sp, "key1": sp1})
print(dom)
print(dom.size)
mfld = ift.from_random(dom)
mfld = mfld.sin()
print(mfld.domain)
print(mfld.val)
print()
print()
print()
print()
# Operators
sp = ift.RGSpace((2, 2), (1, 1))
mask = ift.makeField(sp,
np.array([[False, False],
[False, True]])
)
op = ift.MaskOperator(mask)
print(op)
print(op.domain)
print(op.target)
inp = ift.from_random(op.domain)
print("Input")
print(inp)
print("Output")
print(op(inp))
import nifty8 as ift
from matplotlib import pyplot as plt
import numpy as np
plt.ion()
position_space = ift.RGSpace([128,128])
harmonic_space = position_space.get_default_codomain()
HT = ift.HarmonicTransformOperator(harmonic_space, position_space)
def sqrtpspec(k):
return 1./(20. + k**2)
p_space = ift.PowerSpace(harmonic_space)
pd = ift.PowerDistributor(harmonic_space, p_space)
a = ift.PS_field(p_space, sqrtpspec)
A_field = pd(a)
A = ift.makeOp(A_field)
xi = ift.FieldAdapter(harmonic_space,'xi')
args = {
'offset_mean': 0,
'offset_std': (1e-3, 1e-6),
# Amplitude of field fluctuations
'fluctuations': (1., 0.8), # 1.0, 1e-2
# Exponent of power law power spectrum component
'loglogavgslope': (-3., 1), # -6.0, 1
# Amplitude of integrated Wiener process power spectrum component
'flexibility': (2, 1.), # 1.0, 0.5
# How ragged the integrated Wiener process component is
'asperity': (0.5, 0.4) # 0.1, 0.5
}
correlated_field = ift.SimpleCorrelatedField(position_space, **args)
# GP = HT @ A @ xi
xi_true = ift.from_random(correlated_field.domain)
s = correlated_field(xi_true)
sky = ift.exp(correlated_field)
sky_true = sky(xi_true)
data_raw = np.random.poisson(sky_true.val)
data = ift.makeField(sky.target, data_raw)
log_likelihood = ift.PoissonianEnergy(data) @ sky
ic_sampling = ift.AbsDeltaEnergyController(name="Sampling (linear)",
deltaE=0.05, iteration_limit=100)
ic_sampling_nl = ift.AbsDeltaEnergyController(name="Sampling (non-linear)",
deltaE=0.05, iteration_limit=10)
ic_newton = ift.AbsDeltaEnergyController(name='Newton', deltaE=0.5,
convergence_level=2, iteration_limit=5)
minimizer = ift.NewtonCG(ic_newton)
minimizer_sampling = (lambda iiter: None if iiter < 3 else ift.NewtonCG(ic_sampling_nl))
n_iterations = 5
n_samples = 5
samples = ift.optimize_kl(log_likelihood, n_iterations, n_samples, minimizer,
ic_sampling, minimizer_sampling)
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment