Merge branch 'master' into plotting

README.md

@@ -101,7 +101,7 @@ Starting with a fresh Ubuntu installation move to a folder like
 
         git clone https://gitlab.mpcdf.mpg.de/ift/pyHealpix.git
         cd pyHealpix
-        autoreconf -i && ./configure && sudo make install
+        autoreconf -i && ./configure && make -j4 && sudo make install
         cd ..
 
 -   Finally, NIFTy:
@@ -135,7 +135,7 @@ For pyHealpix, use:
 
     git clone https://gitlab.mpcdf.mpg.de/ift/pyHealpix.git
     cd pyHealpix
-    autoreconf -i && ./configure --prefix=$HOME/.local && make install
+    autoreconf -i && ./configure --prefix=$HOME/.local && make -j4 && make install
     cd ..
 
 ### Installation on OS X 10.11
@@ -156,7 +156,7 @@ may cause trouble.
 
         git clone https://gitlab.mpcdf.mpg.de/ift/pyHealpix.git
         cd pyHealpix
-        autoreconf -i && ./configure && sudo make install
+        autoreconf -i && ./configure --prefix=`python-config --prefix` && make -j4 && sudo make install
         cd ..
 
 -   Install NIFTy:

demos/wiener_filter_hamiltonian.py

@@ -107,7 +107,7 @@ if __name__ == "__main__":
 #                        callback=distance_measure,
 #                        max_history_length=3)
 
-    m0 = Field(s_space, val=1)
+    m0 = Field(s_space, val=1.)
 
     energy = WienerFilterEnergy(position=m0, D=D, j=j)
 

nifty/domain_object.py

@@ -17,13 +17,14 @@
 # along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
 import abc
+from nifty.nifty_meta import NiftyMeta
 
 from keepers import Loggable,\
                     Versionable
 
 
 class DomainObject(Versionable, Loggable, object):
-    __metaclass__ = abc.ABCMeta
+    __metaclass__ = NiftyMeta
 
     def __init__(self):
         # _global_id is used in the Versioning module from keepers
@@ -40,6 +41,20 @@ class DomainObject(Versionable, Loggable, object):
         return result_hash
 
     def __eq__(self, x):
+        """ Checks if two domain_objects are equal.
+
+        Parameters
+        ----------
+        x: domain_object
+            The domain_object `self` is compared to.
+
+        Returns
+        -------
+        bool
+            True if `self` and x describe the same manifold.
+
+        """
+
         if isinstance(x, type(self)):
             for key in vars(self).keys():
                 item1 = vars(self)[key]
@@ -57,23 +72,137 @@ class DomainObject(Versionable, Loggable, object):
 
     @abc.abstractproperty
     def shape(self):
+        """ Returns the shape of the underlying array-like object.
+
+        Returns
+        -------
+        tuple of ints
+            The shape of the underlying array-like object.
+
+        Raises
+        ------
+        NotImplementedError
+            If called for this abstract class.
+
+        """
+
         raise NotImplementedError(
             "There is no generic shape for DomainObject.")
 
     @abc.abstractproperty
     def dim(self):
+        """ Returns the number of pixel-dimensions the object has.
+
+        Returns
+        -------
+        int
+            An Integer representing the number of pixels the discretized
+            manifold has.
+
+        Raises
+        ------
+        NotImplementedError
+            If called for this abstract class.
+
+        """
+
         raise NotImplementedError(
             "There is no generic dim for DomainObject.")
 
     @abc.abstractmethod
     def weight(self, x, power=1, axes=None, inplace=False):
+        """ Weights the field on this domain with the space's volume-weights.
+
+        Weights hereby refer to integration weights, as they appear in
+        discretized integrals. Per default, this function mutliplies each bin
+        of the field x by its volume, which lets it behave like a density
+        (top form). However, different powers of the volume can be applied
+        with the power parameter. The axes parameter specifies which of the
+        field array's indices correspond to this domain.
+
+        Parameters
+        ----------
+        x : distributed_data_object
+            The fields data array.
+        power : int, *optional*
+            The power to which the volume-weight is raised (default: 1).
+        axes : {int, tuple}, *optional*
+            Specifies the axes of x which represent this domain
+            (default: None).
+            If axes==None:
+                weighting is applied with respect to all axes
+        inplace : bool, *optional*
+            If this is True, the weighting is done on the values of x,
+            if it is False, x is not modified and this method returns a
+            weighted copy of x (default: False).
+
+        Returns
+        -------
+        distributed_data_object
+            A weighted version of x, with volume-weights raised to the
+            given power.
+
+        Raises
+        ------
+        NotImplementedError
+            If called for this abstract class.
+
+        """
         raise NotImplementedError(
             "There is no generic weight-method for DomainObject.")
 
-    def pre_cast(self, x, axes=None):
+    def pre_cast(self, x, axes):
+        """ Casts input for Field.val before Field performs the cast.
+
+        Parameters
+        ----------
+        x : {array-like, castable}
+            an array-like object or anything that can be cast to arrays.
+        axes : tuple of ints
+            Specifies the axes of x which correspond to this domain.
+
+        Returns
+        -------
+        {array-like, castable}
+            Processed input where casting that needs Space-specific knowledge
+            (for example location of pixels on the manifold) was performed.
+
+
+        See Also
+        --------
+        post_cast
+
+        Notes
+        -----
+            Usually returns x, except if a power spectrum is given to a
+            PowerSpace, where this spectrum is evaluated at the power indices.
+
+        """
+
         return x
 
-    def post_cast(self, x, axes=None):
+    def post_cast(self, x, axes):
+        """ Performs casting operations that are done after Field's cast.
+
+        Parameters
+        ----------
+        x : {array-like, castable}
+            an array-like object or anything that can be cast to arrays.
+        axes : tuple of ints
+            Specifies the axes of x which correspond to this domain.
+
+        See Also
+        --------
+        pre_cast
+
+        Returns
+        -------
+        distributed_data_object
+            Processed input where casting that needs Space-specific knowledge
+            (for example location of pixels on the manifold) was performed.
+
+        """
+
         return x
 
     # ---Serialization---

nifty/energies/energy.py

@@ -16,10 +16,14 @@
 # You should have received a copy of the GNU General Public License
 # along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
+from nifty.nifty_meta import NiftyMeta
+
 from keepers import Loggable
 
 
 class Energy(Loggable, object):
+    __metaclass__ = NiftyMeta
+
     def __init__(self, position):
         self._cache = {}
         try:

nifty/field.py

@@ -167,8 +167,8 @@ class Field(Loggable, Versionable, object):
 
     # ---Powerspectral methods---
 
-    def power_analyze(self, spaces=None, log=False, nbin=None, binbounds=None,
-                      real_signal=True):
+    def power_analyze(self, spaces=None, logarithmic=False, nbin=None, binbounds=None,
+                      decompose_power=False):
         # check if all spaces in `self.domain` are either harmonic or
         # power_space instances
         for sp in self.domain:
@@ -210,18 +210,18 @@ class Field(Loggable, Versionable, object):
             self.val.get_axes_local_distribution_strategy(
                 self.domain_axes[space_index])
 
-        harmonic_domain = self.domain[space_index]
-        power_domain = PowerSpace(harmonic_domain=harmonic_domain,
+        harmonic_partner = self.domain[space_index]
+        power_domain = PowerSpace(harmonic_partner=harmonic_partner,
                                   distribution_strategy=distribution_strategy,
-                                  log=log, nbin=nbin, binbounds=binbounds)
+                                  logarithmic=logarithmic, nbin=nbin, binbounds=binbounds)
 
         # extract pindex and rho from power_domain
         pindex = power_domain.pindex
         rho = power_domain.rho
 
-        if real_signal:
+        if decompose_power:
             hermitian_part, anti_hermitian_part = \
-                harmonic_domain.hermitian_decomposition(
+                harmonic_partner.hermitian_decomposition(
                                             self.val,
                                             axes=self.domain_axes[space_index])
 
@@ -245,7 +245,7 @@ class Field(Loggable, Versionable, object):
         result_domain = list(self.domain)
         result_domain[space_index] = power_domain
 
-        if real_signal:
+        if decompose_power:
             result_dtype = np.complex
         else:
             result_dtype = np.float
@@ -303,8 +303,8 @@ class Field(Loggable, Versionable, object):
 
         return result_obj
 
-    def power_synthesize(self, spaces=None, real_power=True, real_signal=True,
-                         mean=None, std=None):
+    def power_synthesize(self, spaces=None, real_power=True,
+                         real_signal=False, mean=None, std=None):
 
         # check if the `spaces` input is valid
         spaces = utilities.cast_axis_to_tuple(spaces, len(self.domain))
@@ -322,8 +322,8 @@ class Field(Loggable, Versionable, object):
         result_domain = list(self.domain)
         for power_space_index in spaces:
             power_space = self.domain[power_space_index]
-            harmonic_domain = power_space.harmonic_domain
-            result_domain[power_space_index] = harmonic_domain
+            harmonic_partner = power_space.harmonic_partner
+            result_domain[power_space_index] = harmonic_partner
 
         # create random samples: one or two, depending on whether the
         # power spectrum is real or complex
@@ -365,8 +365,8 @@ class Field(Loggable, Versionable, object):
 
         if real_signal:
             for power_space_index in spaces:
-                harmonic_domain = result_domain[power_space_index]
-                result_val_list = [harmonic_domain.hermitian_decomposition(
+                harmonic_partner = result_domain[power_space_index]
+                result_val_list = [harmonic_partner.hermitian_decomposition(
                                     result_val,
                                     axes=result.domain_axes[power_space_index],
                                     preserve_gaussian_variance=True)[0]

nifty/minimization/quasi_newton_minimizer.py

@@ -17,6 +17,7 @@
 # along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
 import abc
+from nifty.nifty_meta import NiftyMeta
 
 import numpy as np
 
@@ -26,7 +27,7 @@ from .line_searching import LineSearchStrongWolfe
 
 
 class QuasiNewtonMinimizer(Loggable, object):
-    __metaclass__ = abc.ABCMeta
+    __metaclass__ = NiftyMeta
 
     def __init__(self, line_searcher=LineSearchStrongWolfe(), callback=None,
                  convergence_tolerance=1E-4, convergence_level=3,

nifty/nifty_meta.py

+# -*- coding: utf-8 -*-
+
+import abc
+
+
+class DocStringInheritor(type):
+    """
+    A variation on
+    http://groups.google.com/group/comp.lang.python/msg/26f7b4fcb4d66c95
+    by Paul McGuire
+    """
+    def __new__(meta, name, bases, clsdict):
+        if not('__doc__' in clsdict and clsdict['__doc__']):
+            for mro_cls in (mro_cls for base in bases
+                            for mro_cls in base.mro()):
+                doc = mro_cls.__doc__
+                if doc:
+                    clsdict['__doc__'] = doc
+                    break
+        for attr, attribute in clsdict.items():
+            if not attribute.__doc__:
+                for mro_cls in (mro_cls for base in bases
+                                for mro_cls in base.mro()
+                                if hasattr(mro_cls, attr)):
+                    doc = getattr(getattr(mro_cls, attr), '__doc__')
+                    if doc:
+                        if isinstance(attribute, property):
+                            clsdict[attr] = property(attribute.fget,
+                                                     attribute.fset,
+                                                     attribute.fdel,
+                                                     doc)
+                        else:
+                            attribute.__doc__ = doc
+                        break
+        return super(DocStringInheritor, meta).__new__(meta, name,
+                                                       bases, clsdict)
+
+
+class NiftyMeta(DocStringInheritor, abc.ABCMeta):
+    pass

nifty/nifty_utilities.py

@@ -19,6 +19,7 @@
 import numpy as np
 from itertools import product
 
+
 def get_slice_list(shape, axes):
     """
     Helper function which generates slice list(s) to traverse over all
@@ -65,8 +66,7 @@ def get_slice_list(shape, axes):
         return
 
 
-
-def cast_axis_to_tuple(axis, length):
+def cast_axis_to_tuple(axis, length=None):
     if axis is None:
         return None
     try:
@@ -78,16 +78,17 @@ def cast_axis_to_tuple(axis, length):
             raise TypeError(
                 "Could not convert axis-input to tuple of ints")
 
-    # shift negative indices to positive ones
-    axis = tuple(item if (item >= 0) else (item + length) for item in axis)
+    if length is not None:
+        # shift negative indices to positive ones
+        axis = tuple(item if (item >= 0) else (item + length) for item in axis)
 
-    # Deactivated this, in order to allow for the ComposedOperator
-    # remove duplicate entries
-    # axis = tuple(set(axis))
+        # Deactivated this, in order to allow for the ComposedOperator
+        # remove duplicate entries
+        # axis = tuple(set(axis))
 
-    # assert that all entries are elements in [0, length]
-    for elem in axis:
-        assert (0 <= elem < length)
+        # assert that all entries are elements in [0, length]
+        for elem in axis:
+            assert (0 <= elem < length)
 
     return axis
 

nifty/operators/composed_operator/composed_operator.py

@@ -21,7 +21,9 @@ from nifty.operators.linear_operator import LinearOperator
 
 class ComposedOperator(LinearOperator):
     # ---Overwritten properties and methods---
-    def __init__(self, operators):
+    def __init__(self, operators, default_spaces=None):
+        super(ComposedOperator, self).__init__(default_spaces)
+
         self._operator_store = ()
         for op in operators:
             if not isinstance(op, LinearOperator):
@@ -54,10 +56,6 @@ class ComposedOperator(LinearOperator):
                 self._target += op.target
         return self._target
 
-    @property
-    def implemented(self):
-        return True
-
     @property
     def unitary(self):
         return False

nifty/operators/diagonal_operator/diagonal_operator.py

@@ -30,12 +30,11 @@ class DiagonalOperator(EndomorphicOperator):
 
     # ---Overwritten properties and methods---
 
-    def __init__(self, domain=(), implemented=True,
-                 diagonal=None, bare=False, copy=True,
-                 distribution_strategy=None):
-        self._domain = self._parse_domain(domain)
+    def __init__(self, domain=(), diagonal=None, bare=False, copy=True,
+                 distribution_strategy=None, default_spaces=None):
+        super(DiagonalOperator, self).__init__(default_spaces)
 
-        self._implemented = bool(implemented)
+        self._domain = self._parse_domain(domain)
 
         if distribution_strategy is None:
             if isinstance(diagonal, distributed_data_object):
@@ -100,10 +99,6 @@ class DiagonalOperator(EndomorphicOperator):
     def domain(self):
         return self._domain
 
-    @property
-    def implemented(self):
-        return self._implemented
-
     @property
     def self_adjoint(self):
         if self._self_adjoint is None:
@@ -144,16 +139,12 @@ class DiagonalOperator(EndomorphicOperator):
                   distribution_strategy=self.distribution_strategy,
                   copy=copy)
 
-        # weight if the given values were `bare` and `implemented` is True
+        # weight if the given values were `bare` is True
         # do inverse weightening if the other way around
-        if bare and self.implemented:
+        if bare:
             # If `copy` is True, we won't change external data by weightening
             # Otherwise, inplace weightening would change the external field
             f.weight(inplace=copy)
-        elif not bare and not self.implemented:
-            # If `copy` is True, we won't change external data by weightening
-            # Otherwise, inplace weightening would change the external field
-            f.weight(inplace=copy, power=-1)
 
         # Reset the self_adjoint property:
         self._self_adjoint = None

nifty/operators/fft_operator/fft_operator.py

@@ -34,13 +34,72 @@ from transformations import RGRGTransformation,\
 
 
 class FFTOperator(LinearOperator):
+    """Transforms between a pair of position and harmonic domains.
+
+    Built-in domain pairs are
+      - a harmonic and a non-harmonic RGSpace (with matching distances)
+      - a HPSpace and a LMSpace
+      - a GLSpace and a LMSpace
+    Within a domain pair, both orderings are possible.
+
+    The operator provides a "times" and an "adjoint_times" operation.
+    For a pair of RGSpaces, the "adjoint_times" operation is equivalent to
+    "inverse_times"; for the sphere-related domains this is not the case, since
+    the operator matrix is not square.
+
+    Parameters
+    ----------
+    domain: Space or single-element tuple of Spaces
+        The domain of the data that is input by "times" and output by
+        "adjoint_times".
+    target: Space or single-element tuple of Spaces (optional)
+        The domain of the data that is output by "times" and input by
+        "adjoint_times".
+        If omitted, a co-domain will be chosen automatically.
+        Whenever "domain" is an RGSpace, the codomain (and its parameters) are
+        uniquely determined (except for "zerocenter").
+        For GLSpace, HPSpace, and LMSpace, a sensible (but not unique)
+        co-domain is chosen that should work satisfactorily in most situations,
+        but for full control, the user should explicitly specify a codomain.
+    module: String (optional)
+        Software module employed for carrying out the transform operations.
+        For RGSpace pairs this can be "numpy" or "fftw", where "numpy" is
+        always available, but "fftw" offers higher performance and
+        parallelization. For sphere-related domains, only "pyHealpix" is
+        available. If omitted, "fftw" is selected for RGSpaces if available,
+        else "numpy"; on the sphere the default is "pyHealpix".
+    domain_dtype: data type (optional)
+        Data type of the fields that go into "times" and come out of
+        "adjoint_times". Default is "numpy.complex".
+    target_dtype: data type (optional)
+        Data type of the fields that go into "adjoint_times" and come out of
+        "times". Default is "numpy.complex".
+
+    Attributes
+    ----------
+    domain: Tuple of Spaces (with one entry)
+        The domain of the data that is input by "times" and output by
+        "adjoint_times".
+    target: Tuple of Spaces (with one entry)
+        The domain of the data that is output by "times" and input by
+        "adjoint_times".
+    unitary: bool
+        Returns True if the operator is unitary (currently only the case if
+        the domain and codomain are RGSpaces), else False.
+
+    Raises
+    ------
+    ValueError:
+        if "domain" or "target" are not of the proper type.
+
+    """
 
     # ---Class attributes---
 
     default_codomain_dictionary = {RGSpace: RGSpace,
                                    HPSpace: LMSpace,
                                    GLSpace: LMSpace,
-                                   LMSpace: HPSpace,
+                                   LMSpace: GLSpace,
                                    }
 
     transformation_dictionary = {(RGSpace, RGSpace): RGRGTransformation,
@@ -52,8 +111,9 @@ class FFTOperator(LinearOperator):
 
     # ---Overwritten properties and methods---
 
-    def __init__(self, domain=(), target=None, module=None,
-                 domain_dtype=None, target_dtype=None):
+    def __init__(self, domain, target=None, module=None,
+                 domain_dtype=None, target_dtype=None, default_spaces=None):
+        super(FFTOperator, self).__init__(default_spaces)
 
         # Initialize domain and target
 
@@ -83,8 +143,8 @@ class FFTOperator(LinearOperator):
 
         # Store the dtype information
         if domain_dtype is None:
-            self.logger.info("Setting domain_dtype to np.float.")
-            self.domain_dtype = np.float
+            self.logger.info("Setting domain_dtype to np.complex.")
+            self.domain_dtype = np.complex
         else:
             self.domain_dtype = np.dtype(domain_dtype)
 
@@ -114,11 +174,11 @@ class FFTOperator(LinearOperator):
 
         result_field = x.copy_empty(domain=result_domain,
                                     dtype=self.target_dtype)
-        result_field.set_val(new_val=new_val, copy=False)
+        result_field.set_val(new_val=new_val, copy=True)
 
         return result_field
 
-    def _inverse_times(self, x, spaces):
+    def _adjoint_times(self, x, spaces):
         spaces = utilities.cast_axis_to_tuple(spaces, len(x.domain))
         if spaces is None:
             # this case means that x lives on only one space, which is
@@ -138,7 +198,7 @@ class FFTOperator(LinearOperator):
 
         result_field = x.copy_empty(domain=result_domain,
                                     dtype=self.domain_dtype)
-        result_field.set_val(new_val=new_val, copy=False)
+        result_field.set_val(new_val=new_val, copy=True)
 
         return result_field
 
@@ -152,18 +212,40 @@ class FFTOperator(LinearOperator):
     def target(self):
         return self._target
 
-    @property
-    def implemented(self):
-        return True
-
     @property
     def unitary(self):
-        return True
+        return (self._forward_transformation.unitary and
+                self._backward_transformation.unitary)
 
     # ---Added properties and methods---
 
     @classmethod
     def get_default_codomain(cls, domain):
+        """Returns a codomain to the given domain.
+
+        Parameters
+        ----------
+        domain: Space
+            An instance of RGSpace, HPSpace, GLSpace or LMSpace.
+
+        Returns
+        -------
+        target: Space
+            A (more or less perfect) counterpart to "domain" with respect
+            to a FFT operation.