Commit d6894d23 authored by Martin Reinecke's avatar Martin Reinecke

stage1

parent 766aa216
......@@ -24,16 +24,11 @@ from .version import __version__
from keepers import MPILogger
logger = MPILogger()
# it is important to import config before d2o such that NIFTy is able to
# pre-create d2o's configuration object with the corrected path
from .config import dependency_injector,\
nifty_configuration,\
d2o_configuration
nifty_configuration
logger.logger.setLevel(nifty_configuration['loglevel'])
from d2o import distributed_data_object, d2o_librarian
from .field import Field
from .random import Random
......
......@@ -18,7 +18,6 @@
from __future__ import division
import numpy as np
from d2o import distributed_data_object
from .field import Field
......@@ -32,8 +31,6 @@ def _math_helper(x, function):
result_val = x.val.apply_scalar_function(function)
result = x.copy_empty(dtype=result_val.dtype)
result.val = result_val
elif isinstance(x, distributed_data_object):
result = x.apply_scalar_function(function, inplace=False)
else:
result = function(np.asarray(x))
......
......@@ -142,7 +142,7 @@ class DomainObject(with_metaclass(
Parameters
----------
x : distributed_data_object
x : numpy.ndarray
The fields data array.
power : int, *optional*
The power to which the volume-weight is raised (default: 1).
......@@ -158,7 +158,7 @@ class DomainObject(with_metaclass(
Returns
-------
distributed_data_object
numpy.ndarray
A weighted version of x, with volume-weights raised to the
given power.
......@@ -217,7 +217,7 @@ class DomainObject(with_metaclass(
Returns
-------
distributed_data_object
numpy.ndarray
Processed input where casting that needs Space-specific knowledge
(for example location of pixels on the manifold) was performed.
......
This diff is collapsed.
......@@ -20,9 +20,6 @@ from __future__ import division
from builtins import range
import numpy as np
from d2o import distributed_data_object,\
STRATEGIES as DISTRIBUTION_STRATEGIES
from ...config import nifty_configuration as gc
from ...field import Field
from ..endomorphic_operator import EndomorphicOperator
......@@ -47,10 +44,6 @@ class DiagonalOperator(EndomorphicOperator):
(default: False).
copy : boolean
Internal copy of the diagonal (default: True)
distribution_strategy : string
setting the prober distribution_strategy of the
diagonal (default : None). In case diagonal is d2o-object or Field,
their distribution_strategy is used as a fallback.
default_spaces : tuple of ints *optional*
Defines on which space(s) of a given field the Operator acts by
default (default: None)
......@@ -66,9 +59,6 @@ class DiagonalOperator(EndomorphicOperator):
Indicates whether the Operator is unitary or not.
self_adjoint : boolean
Indicates whether the operator is self_adjoint or not.
distribution_strategy : string
Defines the distribution_strategy of the distributed_data_object
in which the diagonal entries are stored in.
Raises
------
......@@ -83,16 +73,6 @@ class DiagonalOperator(EndomorphicOperator):
deals with the bare entries that allow for correct interpretation
of the matrix entries; e.g., as variance in case of an covariance operator.
Examples
--------
>>> x_space = RGSpace(5)
>>> D = DiagonalOperator(x_space, diagonal=[1., 3., 2., 4., 6.])
>>> f = Field(x_space, val=2.)
>>> res = D.times(f)
>>> res.val
<distributed_data_object>
array([ 2., 6., 4., 8., 12.])
See Also
--------
EndomorphicOperator
......@@ -102,28 +82,17 @@ class DiagonalOperator(EndomorphicOperator):
# ---Overwritten properties and methods---
def __init__(self, domain=(), diagonal=None, bare=False, copy=True,
distribution_strategy=None, default_spaces=None):
default_spaces=None):
super(DiagonalOperator, self).__init__(default_spaces)
self._domain = self._parse_domain(domain)
if distribution_strategy is None:
if isinstance(diagonal, distributed_data_object):
distribution_strategy = diagonal.distribution_strategy
elif isinstance(diagonal, Field):
distribution_strategy = diagonal.distribution_strategy
self._distribution_strategy = self._parse_distribution_strategy(
distribution_strategy=distribution_strategy,
val=diagonal)
self._self_adjoint = None
self._unitary = None
self.set_diagonal(diagonal=diagonal, bare=bare, copy=copy)
def _add_attributes_to_copy(self, copy, **kwargs):
copy._domain = self._domain
copy._distribution_strategy = self._distribution_strategy
copy.set_diagonal(diagonal=self.diagonal(bare=True), bare=True)
copy._self_adjoint = self._self_adjoint
copy._unitary = self._unitary
......@@ -207,35 +176,6 @@ class DiagonalOperator(EndomorphicOperator):
# ---Added properties and methods---
@property
def distribution_strategy(self):
"""
distribution_strategy : string
Defines the way how the diagonal operator is distributed
among the nodes. Available distribution_strategies are:
'fftw', 'equal' and 'not'.
Notes :
https://arxiv.org/abs/1606.05385
"""
return self._distribution_strategy
def _parse_distribution_strategy(self, distribution_strategy, val):
if distribution_strategy is None:
if isinstance(val, distributed_data_object):
distribution_strategy = val.distribution_strategy
elif isinstance(val, Field):
distribution_strategy = val.distribution_strategy
else:
self.logger.info("Datamodel set to default!")
distribution_strategy = gc['default_distribution_strategy']
elif distribution_strategy not in DISTRIBUTION_STRATEGIES['all']:
raise ValueError(
"Invalid distribution_strategy!")
return distribution_strategy
def set_diagonal(self, diagonal, bare=False, copy=True):
""" Sets the diagonal of the Operator.
......@@ -254,7 +194,6 @@ class DiagonalOperator(EndomorphicOperator):
# use the casting functionality from Field to process `diagonal`
f = Field(domain=self.domain,
val=diagonal,
distribution_strategy=self.distribution_strategy,
copy=copy)
# weight if the given values were `bare` is True
......@@ -280,8 +219,6 @@ class DiagonalOperator(EndomorphicOperator):
# here the actual multiplication takes place
return operation(self.diagonal(copy=False))(x)
# if the distribution_strategy of self is sub-slice compatible to
# the one of x, reshape the local data of self and apply it directly
active_axes = []
if spaces is None:
active_axes = list(range(len(x.shape)))
......@@ -289,17 +226,7 @@ class DiagonalOperator(EndomorphicOperator):
for space_index in spaces:
active_axes += x.domain_axes[space_index]
axes_local_distribution_strategy = \
x.val.get_axes_local_distribution_strategy(active_axes)
if axes_local_distribution_strategy == self.distribution_strategy:
local_diagonal = self._diagonal.val.get_local_data(copy=False)
else:
# create an array that is sub-slice compatible
self.logger.warn("The input field is not sub-slice compatible to "
"the distribution strategy of the operator.")
redistr_diagonal_val = self._diagonal.val.copy(
distribution_strategy=axes_local_distribution_strategy)
local_diagonal = redistr_diagonal_val.get_local_data(copy=False)
local_diagonal = self._diagonal.val.get_local_data(copy=False)
reshaper = [x.val.data.shape[i] if i in active_axes else 1
for i in range(len(x.shape))]
......
......@@ -21,7 +21,6 @@ from builtins import object
import warnings
import numpy as np
from d2o import distributed_data_object, STRATEGIES
from ....config import dependency_injector as gdi
from ....config import nifty_configuration as gc
from .... import nifty_utilities as utilities
......@@ -331,15 +330,13 @@ class MPIFFT(Transform):
return return_val
def _repack_to_fftw_and_transform(self, val, axes, **kwargs):
temp_val = val.copy_empty(distribution_strategy='fftw')
self.logger.info("Repacking d2o to fftw distribution strategy")
temp_val = val.copy_empty()
temp_val.set_full_data(val, copy=False)
# Recursive call to transform
result = self.transform(temp_val, axes, **kwargs)
return_val = result.copy_empty(
distribution_strategy=val.distribution_strategy)
return_val = result.copy_empty()
return_val.set_full_data(data=result, copy=False)
return return_val
......@@ -433,32 +430,10 @@ class MPIFFT(Transform):
raise ValueError("Provided axes does not match array shape")
# If the input is a numpy array we transform it locally
if not isinstance(val, distributed_data_object):
# Cast to a np.ndarray
temp_val = np.asarray(val)
temp_val = np.asarray(val)
# local transform doesn't apply transforms inplace
return_val = self._local_transform(temp_val, axes)
else:
if val.distribution_strategy in STRATEGIES['slicing']:
if axes is None or 0 in axes:
if val.distribution_strategy != 'fftw':
return_val = \
self._repack_to_fftw_and_transform(
val, axes, **kwargs
)
else:
return_val = self._mpi_transform(
val, axes, **kwargs
)
else:
return_val = self._local_transform(
val, axes, **kwargs
)
else:
return_val = self._repack_to_fftw_and_transform(
val, axes, **kwargs
)
# local transform doesn't apply transforms inplace
return_val = self._local_transform(temp_val, axes)
return return_val
......@@ -590,35 +565,14 @@ class SerialFFT(Transform):
not all(axis in range(len(val.shape)) for axis in axes):
raise ValueError("Provided axes does not match array shape")
return_val = val.copy_empty(global_shape=val.shape,
dtype=np.complex)
if (axes is None) or (0 in axes) or \
(val.distribution_strategy not in STRATEGIES['slicing']):
return_val = np.empty(val.shape, dtype=np.complex)
if val.distribution_strategy == 'not':
local_val = val.get_local_data(copy=False)
else:
local_val = val.get_full_data()
local_val = val
result_data = self._atomic_transform(local_val=local_val,
axes=axes,
local_offset_Q=False)
return_val.set_full_data(result_data, copy=False)
else:
local_offset_list = np.cumsum(
np.concatenate([[0, ],
val.distributor.all_local_slices[:, 2]]))
local_offset_Q = \
bool(local_offset_list[val.distributor.comm.rank] % 2)
local_val = val.get_local_data()
result_data = self._atomic_transform(local_val=local_val,
axes=axes,
local_offset_Q=local_offset_Q)
return_val.set_local_data(result_data, copy=False)
result_data = self._atomic_transform(local_val=local_val,
axes=axes,
local_offset_Q=False)
return_val=result_data
return return_val
......
......@@ -156,21 +156,21 @@ class RGRGTransformation(Transformation):
Tval_imag = self._transform.transform(val.imag, axes,
**kwargs)
if self.codomain.harmonic:
Tval_real.data.real += Tval_real.data.imag
Tval_real.data.imag = \
Tval_imag.data.real + Tval_imag.data.imag
Tval_real.real += Tval_real.imag
Tval_real.imag = \
Tval_imag.real + Tval_imag.imag
else:
Tval_real.data.real -= Tval_real.data.imag
Tval_real.data.imag = \
Tval_imag.data.real - Tval_imag.data.imag
Tval_real.real -= Tval_real.imag
Tval_real.imag = \
Tval_imag.real - Tval_imag.imag
Tval = Tval_real
else:
Tval = self._transform.transform(val, axes, **kwargs)
if self.codomain.harmonic:
Tval.data.real += Tval.data.imag
Tval.real += Tval.imag
else:
Tval.data.real -= Tval.data.imag
Tval.real -= Tval.imag
Tval = Tval.real
if not self._transform.codomain.harmonic:
......
......@@ -35,14 +35,12 @@ class SlicingTransformation(Transformation):
for slice_list in utilities.get_slice_list(val.shape, axes):
if return_val is None:
return_val = val.copy_empty(global_shape=return_shape)
data = val.get_data(slice_list, copy=False)
data = data.get_full_data()
return_val = np.empty(return_shape,dtype=val.dtype)
data = val[slice_list]
data = self._transformation_of_slice(data, **kwargs)
return_val.set_data(data=data, to_key=slice_list, copy=False)
return_val[slice_list] = data
return return_val
......
......@@ -44,27 +44,6 @@ class ResponseOperator(LinearOperator):
raised if:
* len of sigma-list and exposure-list are not equal
Notes
-----
Examples
--------
>>> x1 = RGSpace(5)
>>> x2 = RGSpace(10)
>>> R = ResponseOperator(domain=(x1,x2), sigma=[.5, .25],
exposure=[2.,3.])
>>> f = Field((x1,x2), val=4.)
>>> R.times(f)
<distributed_data_object>
array([[ 24., 24., 24., 24., 24., 24., 24., 24., 24., 24.],
[ 24., 24., 24., 24., 24., 24., 24., 24., 24., 24.],
[ 24., 24., 24., 24., 24., 24., 24., 24., 24., 24.],
[ 24., 24., 24., 24., 24., 24., 24., 24., 24., 24.],
[ 24., 24., 24., 24., 24., 24., 24., 24., 24., 24.]])
See Also
--------
"""
def __init__(self, domain, sigma=[1.], exposure=[1.],
......
......@@ -22,8 +22,6 @@ import numpy as np
from ..space import Space
from d2o import arange, distributed_data_object
class LMSpace(Space):
"""
......@@ -126,33 +124,18 @@ class LMSpace(Space):
else:
return x.copy()
def get_distance_array(self, distribution_strategy):
if distribution_strategy == 'not': # short cut
lmax = self.lmax
ldist = np.empty((self.dim,), dtype=np.float64)
ldist[0:lmax+1] = np.arange(lmax+1, dtype=np.float64)
tmp = np.empty((2*lmax+2), dtype=np.float64)
tmp[0::2] = np.arange(lmax+1)
tmp[1::2] = np.arange(lmax+1)
idx = lmax+1
for l in range(1, lmax+1):
ldist[idx:idx+2*(lmax+1-l)] = tmp[2*l:]
idx += 2*(lmax+1-l)
dists = distributed_data_object(
global_shape=self.shape,
dtype=np.float,
distribution_strategy=distribution_strategy)
dists.set_local_data(ldist)
return dists
dists = arange(start=0, stop=self.shape[0],
distribution_strategy=distribution_strategy)
dists = dists.apply_scalar_function(
lambda x: self._distance_array_helper(x, self.lmax),
dtype=np.float64)
return dists
def get_distance_array(self):
lmax = self.lmax
ldist = np.empty((self.dim,), dtype=np.float64)
ldist[0:lmax+1] = np.arange(lmax+1, dtype=np.float64)
tmp = np.empty((2*lmax+2), dtype=np.float64)
tmp[0::2] = np.arange(lmax+1)
tmp[1::2] = np.arange(lmax+1)
idx = lmax+1
for l in range(1, lmax+1):
ldist[idx:idx+2*(lmax+1-l)] = tmp[2*l:]
idx += 2*(lmax+1-l)
return ldist
def get_unique_distances(self):
return np.arange(self.lmax+1, dtype=np.float64)
......
......@@ -20,9 +20,6 @@
import ast
import numpy as np
from d2o import distributed_data_object,\
STRATEGIES as DISTRIBUTION_STRATEGIES
from ...spaces.space import Space
from functools import reduce
from ...config import nifty_configuration as gc
......@@ -35,10 +32,6 @@ class PowerSpace(Space):
----------
harmonic_partner : Space
The harmonic Space of which this is the power space.
distribution_strategy : str *optional*
The distribution strategy used for the distributed_data_objects
derived from this PowerSpace, e.g. the pindex.
(default : 'not')
logarithmic : bool *optional*
True if logarithmic binning should be used (default : None).
nbin : {int, None} *optional*
......@@ -61,7 +54,7 @@ class PowerSpace(Space):
Attributes
----------
pindex : distributed_data_object
pindex : numpy.ndarray
This holds the information which pixel of the harmonic partner gets
mapped to which power bin
kindex : numpy.ndarray
......@@ -91,18 +84,11 @@ class PowerSpace(Space):
# ---Overwritten properties and methods---
def __init__(self, harmonic_partner, distribution_strategy=None,
def __init__(self, harmonic_partner,
logarithmic=None, nbin=None, binbounds=None):
super(PowerSpace, self).__init__()
self._ignore_for_hash += ['_pindex', '_kindex', '_rho']
if distribution_strategy is None:
distribution_strategy = gc['default_distribution_strategy']
elif distribution_strategy not in DISTRIBUTION_STRATEGIES['global']:
raise ValueError(
"distribution_strategy must be a global-type "
"strategy.")
if not (isinstance(harmonic_partner, Space) and
harmonic_partner.harmonic):
raise ValueError("harmonic_partner must be a harmonic space.")
......@@ -116,26 +102,25 @@ class PowerSpace(Space):
if binbounds is not None:
binbounds = tuple(binbounds)
key = (harmonic_partner, distribution_strategy, logarithmic, nbin,
binbounds)
key = (harmonic_partner, logarithmic, nbin, binbounds)
if self._powerIndexCache.get(key) is None:
distance_array = \
self.harmonic_partner.get_distance_array(distribution_strategy)
self.harmonic_partner.get_distance_array()
temp_binbounds = self._compute_binbounds(
harmonic_partner=self.harmonic_partner,
distribution_strategy=distribution_strategy,
logarithmic=logarithmic,
nbin=nbin,
binbounds=binbounds)
temp_pindex = self._compute_pindex(
harmonic_partner=self.harmonic_partner,
distance_array=distance_array,
binbounds=temp_binbounds,
distribution_strategy=distribution_strategy)
temp_rho = temp_pindex.bincount().get_full_data()
temp_kindex = \
(temp_pindex.bincount(weights=distance_array).get_full_data() /
temp_rho)
binbounds=temp_binbounds)
#temp_rho = temp_pindex.bincount().get_full_data()
temp_rho = np.bincount(temp_pindex.flatten())
#temp_kindex = \
# (temp_pindex.bincount(weights=distance_array).get_full_data() /
# temp_rho)
temp_kindex = np.bincount(temp_pindex.flatten(),weights=distance_array.flatten())/temp_rho
self._powerIndexCache[key] = (temp_binbounds,
temp_pindex,
temp_kindex,
......@@ -145,7 +130,7 @@ class PowerSpace(Space):
self._powerIndexCache[key]
@staticmethod
def _compute_binbounds(harmonic_partner, distribution_strategy,
def _compute_binbounds(harmonic_partner,
logarithmic, nbin, binbounds):
if logarithmic is None and nbin is None and binbounds is None:
......@@ -182,19 +167,10 @@ class PowerSpace(Space):
return result
@staticmethod
def _compute_pindex(harmonic_partner, distance_array, binbounds,
distribution_strategy):
# Compute pindex, kindex and rho according to bb
pindex = distributed_data_object(
global_shape=distance_array.shape,
dtype=np.int,
distribution_strategy=distribution_strategy)
def _compute_pindex(harmonic_partner, distance_array, binbounds):
if binbounds is None:
binbounds = harmonic_partner.get_natural_binbounds()
pindex.set_local_data(
np.searchsorted(binbounds, distance_array.get_local_data()))
return pindex
return np.searchsorted(binbounds, distance_array)
def pre_cast(self, x, axes):
""" Casts power spectrum functions to discretized power spectra.
......@@ -224,10 +200,8 @@ class PowerSpace(Space):
# ---Mandatory properties and methods---
def __repr__(self):
return ("PowerSpace(harmonic_partner=%r, distribution_strategy=%r, "
"binbounds=%r)"
% (self.harmonic_partner, self.pindex.distribution_strategy,
self._binbounds))
return ("PowerSpace(harmonic_partner=%r, binbounds=%r)"
% (self.harmonic_partner, self._binbounds))
@property
def harmonic(self):
......@@ -247,9 +221,7 @@ class PowerSpace(Space):
return float(reduce(lambda x, y: x*y, self.pindex.shape))
def copy(self):
distribution_strategy = self.pindex.distribution_strategy
return self.__class__(harmonic_partner=self.harmonic_partner,
distribution_strategy=distribution_strategy,
binbounds=self._binbounds)
def weight(self, x, power, axes, inplace=False):
......@@ -269,10 +241,8 @@ class PowerSpace(Space):
return result_x
def get_distance_array(self, distribution_strategy):
return distributed_data_object(
self.kindex, dtype=np.float64,
distribution_strategy=distribution_strategy)
def get_distance_array(self):
return self.kindex.copy()
def get_fft_smoothing_kernel_function(self, sigma):
raise NotImplementedError(
......@@ -292,7 +262,7 @@ class PowerSpace(Space):
@property
def pindex(self):
""" A distributed_data_object having the shape of the harmonic partner
""" A numpy.ndarray having the shape of the harmonic partner
space containing the indices of the power bin a pixel belongs to.
"""
return self._pindex
......@@ -313,8 +283,6 @@ class PowerSpace(Space):
def _to_hdf5(self, hdf5_group):
hdf5_group.attrs['binbounds'] = str(self._binbounds)
hdf5_group.attrs['distribution_strategy'] = \
self._pindex.distribution_strategy
return {
'harmonic_partner': self.harmonic_partner,
......@@ -324,5 +292,4 @@ class PowerSpace(Space):
def _from_hdf5(cls, hdf5_group, repository):
hp = repository.get('harmonic_partner', hdf5_group)
bb = ast.literal_eval(hdf5_group.attrs['binbounds'])
ds = hdf5_group.attrs['distribution_strategy']
return PowerSpace(hp, ds, binbounds=bb)
return PowerSpace(hp, binbounds=bb)
......@@ -34,9 +34,6 @@ from functools import reduce
import numpy as np
from d2o import distributed_data_object,\
STRATEGIES