Commit 5b3640bb authored by Martin Reinecke's avatar Martin Reinecke
Browse files

stage 3

parent f98653ce
Pipeline #64968 failed with stages
in 7 minutes and 27 seconds
......@@ -123,15 +123,15 @@ class LMSpace(StructuredDomain):
lm0 = gl.get_default_codomain()
theta = pyHealpix.GL_thetas(gl.nlat)
# evaluate the kernel function at the required thetas
kernel_sphere = Field.from_global_data(gl, func(theta))
kernel_sphere = Field.from_arr(gl, func(theta))
# normalize the kernel such that the integral over the sphere is 4pi
kernel_sphere = kernel_sphere * (4 * np.pi / kernel_sphere.integrate())
# compute the spherical harmonic coefficients of the kernel
op = HarmonicTransformOperator(lm0, gl)
kernel_lm = op.adjoint_times(kernel_sphere.weight(1)).to_global_data()
kernel_lm = op.adjoint_times(kernel_sphere.weight(1)).val
# evaluate the k lengths of the harmonic space
k_lengths = self.get_k_length_array().to_global_data().astype(np.int)
return Field.from_global_data(self, kernel_lm[k_lengths])
k_lengths = self.get_k_length_array().val.astype(np.int)
return Field.from_arr(self, kernel_lm[k_lengths])
@property
def lmax(self):
......
......@@ -103,7 +103,7 @@ def _stats(op, samples):
sc = StatCalculator()
for s in samples:
sc.add(op(s.extract(op.domain)))
return sc.mean.to_global_data(), sc.var.sqrt().to_global_data()
return sc.mean.val, sc.var.sqrt().val
class _LognormalMomentMatching(Operator):
......@@ -140,7 +140,7 @@ class _SlopeRemover(EndomorphicOperator):
def apply(self, x, mode):
self._check_input(x, mode)
x = x.to_global_data()
x = x.val
if mode == self.TIMES:
res = x - x[self._last]*self._sc[self._extender]
else:
......@@ -175,14 +175,14 @@ class _TwoLogIntegrations(LinearOperator):
reverse = sl + (slice(None, None, -1),)
if mode == self.TIMES:
x = x.to_global_data()
x = x.val
res = np.empty(self._target.shape)
res[first] = res[second] = 0
res[from_third] = np.cumsum(x[second], axis=axis)
res[from_third] = (res[from_third] + res[no_border])/2*self._log_vol[extender_sl] + x[first]
res[from_third] = np.cumsum(res[from_third], axis=axis)
else:
x = x.to_global_data_rw()
x = x.val.copy()
res = np.zeros(self._domain.shape)
x[from_third] = np.cumsum(x[from_third][reverse], axis=axis)[reverse]
res[first] += x[from_third]
......@@ -200,7 +200,7 @@ class _Normalization(Operator):
hspace[space] = hspace[space].harmonic_partner
hspace = makeDomain(hspace)
pd = PowerDistributor(hspace, power_space=self._domain[space], space=space)
mode_multiplicity = pd.adjoint(full(pd.target, 1.)).to_global_data_rw()
mode_multiplicity = pd.adjoint(full(pd.target, 1.)).val.copy()
zero_mode = (slice(None),)*self._domain.axes[space][0] + (0,)
mode_multiplicity[zero_mode] = 0
self._mode_multiplicity = from_global_data(self._domain, mode_multiplicity)
......@@ -237,7 +237,7 @@ class _Distributor(LinearOperator):
def apply(self, x, mode):
self._check_input(x, mode)
x = x.to_global_data()
x = x.val
if mode == self.TIMES:
res = x[self._dofdex]
else:
......@@ -504,7 +504,7 @@ class CorrelatedFieldMaker:
scm = 1.
for a in self._a:
op = a.fluctuation_amplitude*self._azm.one_over()
res = np.array([op(from_random('normal', op.domain)).to_global_data()
res = np.array([op(from_random('normal', op.domain)).val
for _ in range(nsamples)])
scm *= res**2 + 1.
return fluctuations_slice_mean/np.mean(np.sqrt(scm))
......
......@@ -121,11 +121,11 @@ def _make_dynamic_operator(target,
elif len(sigc) != len(m.target.shape) - 1:
raise (ValueError, "Shape mismatch!")
c = FieldAdapter(UnstructuredDomain(len(sigc)), keys[1])
c = makeOp(Field.from_global_data(c.target, np.array(sigc)))(c)
c = makeOp(Field(c.target, np.array(sigc)))(c)
lightspeed = ScalingOperator(-0.5, c.target)(c).exp()
scaling = np.array(m.target[0].distances[1:])/m.target[0].distances[0]
scaling = DiagonalOperator(Field.from_global_data(c.target, scaling))
scaling = DiagonalOperator(Field(c.target, scaling))
ops['lightspeed'] = scaling(lightspeed)
c = LightConeOperator(c.target, m.target, quant)(c.exp())
......
......@@ -80,7 +80,7 @@ class _RestOperator(LinearOperator):
def apply(self, x, mode):
self._check_input(x, mode)
res = x.to_global_data()
res = x.val
if mode == self.TIMES:
res = self._gconf.grid2dirty(res)
else:
......@@ -102,10 +102,10 @@ class RadioGridder(LinearOperator):
import nifty_gridder
self._check_input(x, mode)
if mode == self.TIMES:
x = self._bl.ms2vis(x.to_global_data().reshape((-1, 1)), self._idx)
x = self._bl.ms2vis(x.val.reshape((-1, 1)), self._idx)
res = nifty_gridder.vis2grid(self._bl, self._gconf, self._idx, x)
else:
res = nifty_gridder.grid2vis(self._bl, self._gconf, self._idx,
x.to_global_data())
x.val)
res = self._bl.vis2ms(res, self._idx).reshape((-1,))
return from_global_data(self._tgt(mode), res)
......@@ -27,7 +27,7 @@ from ..operators.operator import Operator
def _field_from_function(domain, func, absolute=False):
domain = DomainTuple.make(domain)
k_array = _make_coords(domain, absolute=absolute)
return Field.from_global_data(domain, func(k_array))
return Field(domain, func(k_array))
def _make_coords(domain, absolute=False):
......@@ -57,14 +57,14 @@ class _LightConeDerivative(LinearOperator):
def apply(self, x, mode):
self._check_input(x, mode)
x = x.to_global_data()
x = x.val
res = np.zeros(self._tgt(mode).shape, dtype=self._derivatives.dtype)
for i in range(self.domain.shape[0]):
if mode == self.TIMES:
res += self._derivatives[i]*x[i]
else:
res[i] = np.sum(self._derivatives[i]*x.real)
return Field.from_global_data(self._tgt(mode), res)
return Field(self._tgt(mode), res)
def _cone_arrays(c, domain, sigx, want_gradient):
......@@ -133,9 +133,9 @@ class LightConeOperator(Operator):
def apply(self, x):
islin = isinstance(x, Linearization)
val = x.val.to_global_data() if islin else x.to_global_data()
val = x.val.val if islin else x.val
a, derivs = _cone_arrays(val, self.target, self._sigx, islin)
res = Field.from_global_data(self.target, a)
res = Field(self.target, a)
if not islin:
return res
jac = _LightConeDerivative(x.jac.target, self.target, derivs)(x.jac)
......
......@@ -178,7 +178,7 @@ class LOSResponse(LinearOperator):
# get the shape of the local data slice
w_i = _comp_traverse(localized_pixel_starts,
localized_pixel_ends,
self._local_shape,
self.domain[0].shape,
np.array(self.domain[0].distances),
1./(1./difflen+truncation*sigmas),
difflen,
......@@ -188,7 +188,7 @@ class LOSResponse(LinearOperator):
boxsz = 16
nlos = len(w_i)
npix = np.prod(self._local_shape)
npix = np.prod(self.domain[0].shape)
ntot = 0
for i in w_i:
ntot += len(i[1])
......@@ -203,12 +203,12 @@ class LOSResponse(LinearOperator):
ilos[ofs:ofs+nval] = cnt
iarr[ofs:ofs+nval] = i[0]
xwgt[ofs:ofs+nval] = i[1]
fullidx = np.unravel_index(i[0], self._local_shape)
fullidx = np.unravel_index(i[0], self.domain[0].shape)
tmp = np.zeros(nval, dtype=np.float64)
fct = 1.
for j in range(ndim):
tmp += (fullidx[j]//boxsz)*fct
fct *= self._local_shape[j]
fct *= self.domain[0].shape[j]
tmp += cnt/float(nlos)
tmp += iarr[ofs:ofs+nval]/(float(nlos)*float(npix))
pri[ofs:ofs+nval] = tmp
......@@ -220,7 +220,7 @@ class LOSResponse(LinearOperator):
xwgt = xwgt[xtmp]
self._smat = aslinearoperator(
coo_matrix((xwgt, (ilos, iarr)),
shape=(nlos, np.prod(self._local_shape))))
shape=(nlos, np.prod(self.domain[0].shape))))
self._target = DomainTuple.make(UnstructuredDomain(nlos))
......@@ -228,8 +228,7 @@ class LOSResponse(LinearOperator):
self._check_input(x, mode)
if mode == self.TIMES:
result_arr = self._smat.matvec(x.val.reshape(-1))
return Field.from_global_data(self._target, result_arr,
sum_up=True)
local_input_data = x.to_global_data().reshape(-1)
res = self._smat.rmatvec(local_input_data).reshape(self._local_shape)
return Field(self._target, result_arr)
local_input_data = x.val.reshape(-1)
res = self._smat.rmatvec(local_input_data).reshape(self.domain[0].shape)
return Field(self._domain, res)
......@@ -132,14 +132,10 @@ class MultiField(object):
return MultiField(domain, tuple(Field(dom, val)
for dom in domain._domains))
def to_global_data(self):
return {key: val.to_global_data()
for key, val in zip(self._domain.keys(), self._val)}
@staticmethod
def from_global_data(domain, arr, sum_up=False):
def from_global_data(domain, arr):
return MultiField(
domain, tuple(Field.from_global_data(domain[key], arr[key], sum_up)
domain, tuple(Field(domain[key], arr[key])
for key in domain.keys()))
def norm(self, ord=2):
......
......@@ -137,7 +137,7 @@ def operator_spectrum(A, k, hermitian, which='LM', tol=0):
Ar = SandwichOperator.make(_DomRemover(A.domain).adjoint, A)
M = ssl.LinearOperator(
shape=2*(size,),
matvec=lambda x: Ar(from_global_data(Ar.domain, x)).to_global_data())
matvec=lambda x: Ar(from_global_data(Ar.domain, x)).val)
f = ssl.eigsh if hermitian else ssl.eigs
eigs = f(M, k=k, tol=tol, return_eigenvectors=False, which=which)
return np.flip(np.sort(eigs), axis=0)
......@@ -74,7 +74,7 @@ class DOFDistributor(LinearOperator):
wgt = np.bincount(dofdex.val.ravel(), minlength=nbin)
wgt = wgt*partner.scalar_dvol
else:
dvol = Field.from_global_data(partner, partner.dvol).val
dvol = Field.from_arr(partner, partner.dvol).val
wgt = np.bincount(dofdex.val.ravel(),
minlength=nbin, weights=dvol)
# The explicit conversion to float64 is necessary because bincount
......
......@@ -58,11 +58,9 @@ class DomainTupleFieldInserter(LinearOperator):
def apply(self, x, mode):
self._check_input(x, mode)
# FIXME Make fully MPI compatible without global_data
if mode == self.TIMES:
res = np.zeros(self.target.shape, dtype=x.dtype)
res[self._slc] = x.to_global_data()
return Field.from_global_data(self.target, res)
res[self._slc] = x.val
return Field(self.target, res)
else:
return Field.from_global_data(self.domain,
x.to_global_data()[self._slc])
return Field(self.domain, x.val[self._slc])
......@@ -92,7 +92,7 @@ class FieldZeroPadder(LinearOperator):
# i1 = idx+(-Nyquist,)
# xnew[i1] *= 0.5
else:
xnew[idx + (slice(0, v.shape[d]),)] = x
xnew[idx + (slice(0, v.shape[d]),)] = v
else: # ADJOINT_TIMES
if self._central:
shp = list(x.shape)
......@@ -110,4 +110,5 @@ class FieldZeroPadder(LinearOperator):
xnew = v[idx + (slice(0, tgtshp[d]),)]
curshp[d] = xnew.shape[d]
return Field(self._tgt(mode), xnew)
v = xnew.copy()
return Field(self._tgt(mode), v)
......@@ -95,9 +95,9 @@ class LinearInterpolator(LinearOperator):
def apply(self, x, mode):
self._check_input(x, mode)
x_val = x.to_global_data()
x_val = x.val
if mode == self.TIMES:
res = self._mat.matvec(x_val.reshape(-1))
else:
res = self._mat.rmatvec(x_val).reshape(self.domain.shape)
return Field.from_global_data(self._tgt(mode), res)
return Field(self._tgt(mode), res)
......@@ -23,9 +23,6 @@ from ..field import Field
from .linear_operator import LinearOperator
# MR FIXME: this needs a redesign to avoid most _global_data() calls
# Possible approach: keep everything defined on `domain` distributed and only
# collect the unstructured Fields.
class MaskOperator(LinearOperator):
"""Implementation of a mask response
......@@ -41,17 +38,17 @@ class MaskOperator(LinearOperator):
if not isinstance(flags, Field):
raise TypeError
self._domain = DomainTuple.make(flags.domain)
self._flags = np.logical_not(flags.to_global_data())
self._flags = np.logical_not(flags.val)
self._target = DomainTuple.make(UnstructuredDomain(self._flags.sum()))
self._capability = self.TIMES | self.ADJOINT_TIMES
def apply(self, x, mode):
self._check_input(x, mode)
x = x.to_global_data()
x = x.val
if mode == self.TIMES:
res = x[self._flags]
return Field.from_global_data(self.target, res)
return Field(self.target, res)
res = np.empty(self.domain.shape, x.dtype)
res[self._flags] = x
res[~self._flags] = 0
return Field.from_global_data(self.domain, res)
return Field(self.domain, res)
......@@ -91,4 +91,5 @@ class RegriddingOperator(LinearOperator):
xnew += v[idx + (self._bindex[d-d0]+1,)] * wgt
curshp[d] = xnew.shape[d]
v = xnew.copy()
return Field(self._tgt(mode), xnew)
......@@ -374,7 +374,7 @@ class MatrixProductOperator(EndomorphicOperator):
def apply(self, x, mode):
self._check_input(x, mode)
res = x.to_global_data()
res = x.val
f = self._mat.dot if mode == self.TIMES else self._mat_tr.dot
res = f(res)
return Field.from_global_data(self._domain, res)
return Field(self._domain, res)
......@@ -298,7 +298,7 @@ def _plot1D(f, ax, **kwargs):
dist = dom.distances[0]
xcoord = np.arange(npoints, dtype=np.float64)*dist
for i, fld in enumerate(f):
ycoord = fld.to_global_data()
ycoord = fld.val
plt.plot(xcoord, ycoord, label=label[i],
linewidth=linewidth[i], alpha=alpha[i])
_limit_xy(**kwargs)
......@@ -310,7 +310,7 @@ def _plot1D(f, ax, **kwargs):
plt.yscale(kwargs.pop("yscale", "log"))
xcoord = dom.k_lengths
for i, fld in enumerate(f):
ycoord = fld.to_global_data_rw()
ycoord = fld.val.copy()
ycoord[0] = ycoord[1]
plt.plot(xcoord, ycoord, label=label[i],
linewidth=linewidth[i], alpha=alpha[i])
......@@ -336,9 +336,9 @@ def _plot2D(f, ax, **kwargs):
raise TypeError("need 1D RGSpace as second domain")
if dom[1].shape[0] == 1:
from .sugar import from_global_data
f = from_global_data(f.domain[0], f.to_global_data()[..., 0])
f = from_global_data(f.domain[0], f.val[..., 0])
else:
rgb = _rgb_data(f.to_global_data())
rgb = _rgb_data(f.val)
have_rgb = True
foo = kwargs.pop("norm", None)
......@@ -363,7 +363,7 @@ def _plot2D(f, ax, **kwargs):
**aspect)
else:
im = ax.imshow(
f.to_global_data().T, extent=[0, nx*dx, 0, ny*dy],
f.val.T, extent=[0, nx*dx, 0, ny*dy],
vmin=kwargs.get("zmin"), vmax=kwargs.get("zmax"),
cmap=cmap, origin="lower", **norm, **aspect)
plt.colorbar(im)
......@@ -385,7 +385,7 @@ def _plot2D(f, ax, **kwargs):
if have_rgb:
res[mask] = rgb[base.ang2pix(ptg)]
else:
res[mask] = f.to_global_data()[base.ang2pix(ptg)]
res[mask] = f.val[base.ang2pix(ptg)]
else:
ra = np.linspace(0, 2*np.pi, dom.nlon+1)
dec = pyHealpix.GL_thetas(dom.nlat)
......@@ -395,7 +395,7 @@ def _plot2D(f, ax, **kwargs):
if have_rgb:
res[mask] = rgb[ilat*dom[0].nlon + ilon]
else:
res[mask] = f.to_global_data()[ilat*dom.nlon + ilon]
res[mask] = f.val[ilat*dom.nlon + ilon]
plt.axis('off')
if have_rgb:
plt.imshow(res, origin="lower")
......
......@@ -35,7 +35,7 @@ from .plot import Plot
__all__ = ['PS_field', 'power_analyze', 'create_power_operator',
'create_harmonic_smoothing_operator', 'from_random',
'full', 'from_global_data', 'from_local_data',
'full', 'from_global_data',
'makeDomain', 'sqrt', 'exp', 'log', 'tanh', 'sigmoid',
'sin', 'cos', 'tan', 'sinh', 'cosh', 'log10',
'absolute', 'one_over', 'clip', 'sinc', "log1p", "expm1",
......@@ -297,28 +297,7 @@ def from_global_data(domain, arr):
The newly created random field
"""
if isinstance(domain, (dict, MultiDomain)):
return MultiField.from_global_data(domain, arr, sum_up)
return Field.from_arr(domain, arr)
def from_local_data(domain, arr):
"""Convenience function creating Fields/MultiFields from Numpy arrays or
dicts of Numpy arrays.
Parameters
----------
domain : Domainoid
the intended domain of the output field
arr : Numpy array if `domain` corresponds to a `DomainTuple`,
dictionary of Numpy arrays if `domain` corresponds to a `MultiDomain`
Returns
-------
Field or MultiField
The newly created field
"""
if isinstance(domain, (dict, MultiDomain)):
return MultiField.from_local_data(domain, arr)
return MultiField.from_global_data(domain, arr)
return Field.from_arr(domain, arr)
......@@ -512,7 +491,7 @@ def calculate_position(operator, output):
raise TypeError
if output.domain != operator.target:
raise TypeError
cov = 1e-3*output.to_global_data().max()**2
cov = 1e-3*output.val.max()**2
invcov = ScalingOperator(cov, output.domain).inverse
d = output + invcov.draw_sample(from_inverse=True)
lh = GaussianEnergy(d, invcov)(operator)
......
......@@ -175,7 +175,6 @@ def test_dataconv():
s1 = ift.RGSpace((10,))
ld = np.arange(s1.shape[0])
gd = np.arange(s1.shape[0])
assert_equal(ld, ift.from_local_data(s1, ld).val)
assert_equal(gd, ift.from_global_data(s1, gd).val)
......
......@@ -57,10 +57,10 @@ def test_kl(constants, point_estimates, mirror_samples):
# Test gradient
for kk in h.domain.keys():
res0 = klpure.gradient.to_global_data()[kk]
res0 = klpure.gradient[kk].val
if kk in constants:
res0 = 0*res0
res1 = kl.gradient.to_global_data()[kk]
res1 = kl.gradient[kk].val
assert_allclose(res0, res1)
# Test number of samples
......@@ -70,13 +70,13 @@ def test_kl(constants, point_estimates, mirror_samples):
# Test point_estimates (after drawing samples)
for kk in point_estimates:
for ss in kl.samples:
ss = ss.to_global_data()[kk]
ss = ss[kk].val
assert_allclose(ss, 0*ss)
# Test constants (after some minimization)
cg = ift.GradientNormController(iteration_limit=5)
minimizer = ift.NewtonCG(cg)
kl, _ = minimizer(kl)
diff = (mean0 - kl.position).to_global_data()
diff = (mean0 - kl.position).to_dict()
for kk in constants:
assert_allclose(diff[kk], 0*diff[kk])
assert_allclose(diff[kk].val, 0*diff[kk].val)
......@@ -25,7 +25,7 @@ pmp = pytest.mark.parametrize
def _lin2grad(lin):
return lin.jac(ift.full(lin.domain, 1.)).to_global_data()
return lin.jac(ift.full(lin.domain, 1.)).val
def jt(lin, check):
......@@ -36,7 +36,7 @@ def test_special_gradients():
dom = ift.UnstructuredDomain((1,))
f = ift.full(dom, 2.4)
var = ift.Linearization.make_var(f)
s = f.to_global_data()
s = f.val
jt(var.clip(0, 10), np.ones_like(s))
jt(var.clip(-1, 0), np.zeros_like(s))
......@@ -62,8 +62,8 @@ def test_actual_gradients(f):
eps = 1e-8
var0 = ift.Linearization.make_var(fld)
var1 = ift.Linearization.make_var(fld + eps)
f0 = getattr(var0, f)().val.to_global_data()
f1 = getattr(var1, f)().val.to_global_data()
f0 = getattr(var0, f)().val.val
f1 = getattr(var1, f)().val.val
df0 = (f1 - f0)/eps
df1 = _lin2grad(getattr(var0, f)())
assert_allclose(df0, df1, rtol=100*eps)
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