Merge branch 'master' into 'docstrings_energy_minimization'

Master

See merge request !131

.gitlab-ci.yml

@@ -13,7 +13,7 @@ before_script:
   - apt-get update
   - >
     apt-get install -y build-essential python python-pip python-dev git
-    autoconf gsl-bin libgsl-dev wget python-numpy
+    autoconf gsl-bin libgsl-dev wget python-numpy cython
   - pip install --upgrade -r ci/requirements_base.txt
   - chmod +x ci/*.sh
 

README.md

@@ -38,25 +38,25 @@ certain grids, **fields** that are defined on spaces, and **operators**
 that apply to fields.
 
 -   [Spaces](http://www.mpa-garching.mpg.de/ift/nifty/space.html)
-    -   `rg_space` - *n*-dimensional regular Euclidean grid
-    -   `lm_space` - spherical harmonics
-    -   `gl_space` - Gauss-Legendre grid on the 2-sphere
-    -   `hp_space` - [HEALPix](http://sourceforge.net/projects/healpix/)
+    -   `RGSpace` - *n*-dimensional regular Euclidean grid
+    -   `LMSpace` - spherical harmonics
+    -   `GLSpace` - Gauss-Legendre grid on the 2-sphere
+    -   `HPSpace` - [HEALPix](http://sourceforge.net/projects/healpix/)
         grid on the 2-sphere
 -   [Fields](http://www.mpa-garching.mpg.de/ift/nifty/field.html)
-    -   `field` - generic class for (discretized) fields
+    -   `Field` - generic class for (discretized) fields
 
 <!-- -->
 
-    field.conjugate     field.dim          field.norm
-    field.dot           field.set_val      field.weight
+    Field.conjugate     Field.dim          Field.norm
+    Field.dot           Field.set_val      Field.weight
 
 -   [Operators](http://www.mpa-garching.mpg.de/ift/nifty/operator.html)
-    -   `diagonal_operator` - purely diagonal matrices in a specified
+    -   `DiagonalOperator` - purely diagonal matrices in a specified
         basis
-    -   `projection_operator` - projections onto subsets of a specified
+    -   `ProjectionOperator` - projections onto subsets of a specified
         basis
-    -   `propagator_operator` - information propagator in Wiener filter
+    -   `PropagatorOperator` - information propagator in Wiener filter
         theory
     -   (and more)
 -   (and more)
@@ -74,9 +74,7 @@ Installation
 
 ### Download
 
-The latest release is tagged **v1.0.7** and is available as a source
-package at [](https://gitlab.mpcdf.mpg.de/ift/NIFTy/tags). The current
-version can be obtained by cloning the repository:
+The current version of Nifty3 can be obtained by cloning the repository:
 
     git clone https://gitlab.mpcdf.mpg.de/ift/NIFTy.git
 
@@ -89,53 +87,30 @@ Starting with a fresh Ubuntu installation move to a folder like
 
 -   Install basic packages like python, python-dev, gsl and others:
 
-        sudo apt-get install curl git autoconf
-        sudo apt-get install python-dev python-pip gsl-bin libgsl0-dev libfreetype6-dev libpng-dev libatlas-base-dev
-
--   Using pip install numpy etc...:
-
-        sudo pip install numpy
+        sudo apt-get install curl git autoconf python-dev python-pip python-numpy
 
 -   Install pyHealpix:
 
         git clone https://gitlab.mpcdf.mpg.de/ift/pyHealpix.git
         cd pyHealpix
-        autoreconf -i && ./configure && make -j4 && sudo make install
+        autoreconf -i && ./configure --prefix=$HOME/.local --enable-openmp --enable-native-optimizations && make -j4 && make install
         cd ..
 
 -   Finally, NIFTy:
 
         git clone https://gitlab.mpcdf.mpg.de/ift/NIFTy.git
-        cd nifty
-        sudo python setup.py install
+        cd NIFTy
+        python setup.py install --user
         cd ..
 
-### Installation on a Linux cluster
-
-This is for you if you want to install NIFTy on a HPC machine or cluster
-that is hosted by your university or institute. Most of the dependencies
-will most likely already be there, but you won't have superuser
-privileges. In this case, instead of:
-
-    sudo python setup.py install
-
-use:
-
-    python setup.py install --user
-
-or:
-
-    python setup.py install --install-lib=/SOMEWHERE
+### Installation on Linux systems in general
 
-in the instruction above. This will install the python packages into
-your local user directory.
-
-For pyHealpix, use:
-
-    git clone https://gitlab.mpcdf.mpg.de/ift/pyHealpix.git
-    cd pyHealpix
-    autoreconf -i && ./configure --prefix=$HOME/.local && make -j4 && make install
-    cd ..
+Since all the "unconventional" packages (i.e. pyHealpix and NIFTy) listed in the
+section above are installed
+within the home directory of the user, the installation instructions for these
+should also work on any Linux machine where you do not have root access.
+In this case you have to ensure with your system administrators that the
+"standard" dependencies (python, numpy, etc.) are installed system-wide.
 
 ### Installation on OS X 10.11
 
@@ -154,27 +129,20 @@ may cause trouble.
 
         git clone https://gitlab.mpcdf.mpg.de/ift/pyHealpix.git
         cd pyHealpix
-        autoreconf -i && ./configure --prefix=`python-config --prefix` && make -j4 && sudo make install
+        autoreconf -i && ./configure --prefix=`python-config --prefix` --enable-openmp --enable-native-optimizations && make -j4 && sudo make install
         cd ..
 
+    (The third command installs the package system-wide. User-specific
+    installation would be preferrable, but we haven't found a simple recipe yet
+    how to determine the installation prefix ...)
+
 -   Install NIFTy:
 
         git clone https://gitlab.mpcdf.mpg.de/ift/NIFTy.git
-        cd nifty
-        sudo python setup.py install
+        cd NIFTy
+        python setup.py install --user
         cd ..
 
-### Installation using pypi
-
-NIFTY can be installed using [PyPI](https://pypi.python.org/pypi) and
-**pip** by running the following command:
-
-    pip install ift_nifty
-
-Alternatively, a private or user specific installation can be done by:
-
-    pip install --user ift_nifty
-
 ### Running the tests
 
 In oder to run the tests one needs two additional packages:

ci/install_h5py.sh

 #!/bin/bash
 
-wget https://api.github.com/repos/h5py/h5py/tags -O - | grep tarball_url | grep -v rc | head -n 1 | cut -d '"' -f 4 | wget -i - -O h5py.tar.gz
-tar xzf h5py.tar.gz
-cd h5py-h5py*
 export CC=mpicc
 export HDF5_DIR=/usr/lib/x86_64-linux-gnu/hdf5/openmpi
-python setup.py configure --mpi
-python setup.py build
-python setup.py install
-cd ..
-rm -r h5py-h5py*
-rm h5py.tar.gz
+export HDF5_MPI="ON"
+pip install --no-binary=h5py h5py

ci/install_pyfftw.sh

 #!/bin/bash
 
-git clone -b mpi https://github.com/ultimanet/pyFFTW.git
+git clone -b mpi https://github.com/fredros/pyFFTW.git
 cd pyFFTW/
 CC=mpicc python setup.py build_ext install
 cd ..

ci/requirements.txt

 numpy
+cython
 mpi4py
 matplotlib
 plotly

ci/requirements_base.txt

 numpy
+cython
 nose
 parameterized
 coverage

demos/wiener_filter_harmonic - Kopie.py

-from nifty import *
-from mpi4py import MPI
-import plotly.offline as py
-import plotly.graph_objs as go
-
-
-comm = MPI.COMM_WORLD
-rank = comm.rank
-
-
-def plot_maps(x, name):
-
-    trace = [None]*len(x)
-
-    keys = x.keys()
-    field = x[keys[0]]
-    domain = field.domain[0]
-    shape = len(domain.shape)
-    max_n = domain.shape[0]*domain.distances[0]
-    step = domain.distances[0]
-    x_axis = np.arange(0, max_n, step)
-
-    if shape == 1:
-        for ii in xrange(len(x)):
-            trace[ii] = go.Scatter(x= x_axis, y=x[keys[ii]].val.get_full_data(), name=keys[ii])
-        fig = go.Figure(data=trace)
-
-        py.plot(fig, filename=name)
-
-    elif shape == 2:
-        for ii in xrange(len(x)):
-            py.plot([go.Heatmap(z=x[keys[ii]].val.get_full_data().real)], filename=keys[ii])
-    else:
-        raise TypeError("Only 1D and 2D field plots are supported")
-
-def plot_power(x, name):
-
-    layout = go.Layout(
-        xaxis=dict(
-            type='log',
-            autorange=True
-        ),
-        yaxis=dict(
-            type='log',
-            autorange=True
-        )
-    )
-
-    trace = [None]*len(x)
-
-    keys = x.keys()
-    field = x[keys[0]]
-    domain = field.domain[0]
-    x_axis = domain.kindex
-
-    for ii in xrange(len(x)):
-        trace[ii] = go.Scatter(x= x_axis, y=x[keys[ii]].val.get_full_data(), name=keys[ii])
-
-    fig = go.Figure(data=trace, layout=layout)
-    py.plot(fig, filename=name)
-
-np.random.seed(42)
-
-
-
-if __name__ == "__main__":
-
-    distribution_strategy = 'not'
-
-    # setting spaces
-    npix = np.array([500])  # number of pixels
-    total_volume = 1.  # total length
-
-    # setting signal parameters
-    lambda_s = .05  # signal correlation length
-    sigma_s = 10.  # signal variance
-
-
-    #setting response operator parameters
-    length_convolution = .025
-    exposure = 1.
-
-    # calculating parameters
-    k_0 = 4. / (2 * np.pi * lambda_s)
-    a_s = sigma_s ** 2. * lambda_s * total_volume
-
-    # creation of spaces
-#    x1 = RGSpace([npix,npix], distances=total_volume / npix,
-#                 zerocenter=False)
-#    k1 = RGRGTransformation.get_codomain(x1)
-    x1 = HPSpace(32)
-    k1 = HPLMTransformation.get_codomain(x1)
-
-    p1 = PowerSpace(harmonic_partner=k1, logarithmic=False)
-
-
-    # creating Power Operator with given spectrum
-    spec = (lambda k: a_s / (1 + (k / k_0) ** 2) ** 2)
-    p_field = Field(p1, val=spec)
-    S_op = create_power_operator(k1, spec)
-
-    # creating FFT-Operator and Response-Operator with Gaussian convolution
-    Fft_op = FFTOperator(domain=x1, target=k1,
-                        domain_dtype=np.float64,
-                        target_dtype=np.complex128)
-    R_op = ResponseOperator(x1, sigma=[length_convolution],
-                            exposure=[exposure])
-
-    # drawing a random field
-    sk = p_field.power_synthesize(real_signal=True, mean=0.)
-    s = Fft_op.adjoint_times(sk)
-
-    signal_to_noise = 1
-    N_op = DiagonalOperator(R_op.target, diagonal=s.var()/signal_to_noise, bare=True)
-    n = Field.from_random(domain=R_op.target,
-                          random_type='normal',
-                          std=s.std()/np.sqrt(signal_to_noise),
-                          mean=0.)
-
-    d = R_op(s) + n
-
-    # Wiener filter
-    j = Fft_op.times(R_op.adjoint_times(N_op.inverse_times(d)))
-    D = HarmonicPropagatorOperator(S=S_op, N=N_op, R=R_op)
-
-    mk = D(j)
-
-    m = Fft_op.adjoint_times(mk)
-
-#    z={}
-#    z["signal"] = s
-#    z["reconstructed_map"] = m
-#    z["data"] = d
-#    z["lambda"] = R_op(s)
-#    z["j"] = j
-#
-#    plot_maps(z, "Wiener_filter.html")
-
-

nifty/domain_object.py

@@ -87,7 +87,7 @@ class DomainObject(Versionable, Loggable, object):
 
     @abc.abstractproperty
     def shape(self):
-        """ Returns the shape of the underlying array-like object.
+        """ The domain-object's shape contribution to the underlying array.
 
         Returns
         -------

nifty/energies/energy.py

@@ -71,7 +71,7 @@ class Energy(Loggable, object):
             position = position.copy()
         except AttributeError:
             pass
-        self.position = position
+        self._position = position
 
     def at(self, position):
         """ Initializes and returns new Energy object at new position.

nifty/field_types/field_type.py

@@ -27,16 +27,3 @@ class FieldType(DomainObject):
         else:
             result = x.copy()
         return result
-
-    def process(self, method_name, array, inplace=True, **kwargs):
-        try:
-            result_array = self.__getattr__(method_name)(array,
-                                                         inplace,
-                                                         **kwargs)
-        except AttributeError:
-            if inplace:
-                result_array = array
-            else:
-                result_array = array.copy()
-
-        return result_array

nifty/spaces/lm_space/lm_space.py

@@ -94,9 +94,12 @@ class LMSpace(Space):
         hermitian_part = x.copy_empty()
         anti_hermitian_part = x.copy_empty()
         hermitian_part[:] = x.real
-        anti_hermitian_part[:] = x.imag
+        anti_hermitian_part[:] = x.imag * 1j
         return (hermitian_part, anti_hermitian_part)
 
+    def hermitian_fixed_points(self):
+        return None
+
     # ---Mandatory properties and methods---
 
     @property

nifty/spaces/power_space/power_space.py

@@ -114,6 +114,11 @@ class PowerSpace(Space):
         self._pundex = power_index['pundex']
         self._k_array = power_index['k_array']
 
+        if self.config['nbin'] is not None:
+            if self.config['nbin'] > len(self.kindex):
+                self.logger.warn("nbin was set to a value being larger than "
+                                 "the length of kindex!")
+
     def pre_cast(self, x, axes):
         """ Casts power spectrum functions to discretized power spectra.
 
@@ -249,12 +254,12 @@ class PowerSpace(Space):
 
     def _to_hdf5(self, hdf5_group):
         hdf5_group['kindex'] = self.kindex
-        hdf5_group['rho'] = self.config["rho"]
-        hdf5_group['pundex'] = self.config["pundex"]
+        hdf5_group['rho'] = self.rho
+        hdf5_group['pundex'] = self.pundex
         hdf5_group['logarithmic'] = self.config["logarithmic"]
         # Store nbin as string, since it can be None
-        hdf5_group.attrs['nbin'] = str(self.nbin)
-        hdf5_group.attrs['binbounds'] = str(self.binbounds)
+        hdf5_group.attrs['nbin'] = str(self.config["nbin"])
+        hdf5_group.attrs['binbounds'] = str(self.config["binbounds"])
 
         return {
             'harmonic_partner': self.harmonic_partner,
@@ -274,10 +279,10 @@ class PowerSpace(Space):
         new_ps._harmonic_partner = repository.get('harmonic_partner',
                                                   hdf5_group)
 
-        new_ps.config = {}
-        new_ps.config['logarithmic'] = hdf5_group['logarithmic'][()]
-        exec("new_ps.config['nbin'] = " + hdf5_group.attrs['nbin'])
-        exec("new_ps.config['binbounds'] = " + hdf5_group.attrs['binbounds'])
+        new_ps._config = {}
+        new_ps._config['logarithmic'] = hdf5_group['logarithmic'][()]
+        exec("new_ps._config['nbin'] = " + hdf5_group.attrs['nbin'])
+        exec("new_ps._config['binbounds'] = " + hdf5_group.attrs['binbounds'])
 
         new_ps._pindex = repository.get('pindex', hdf5_group)
         new_ps._kindex = hdf5_group['kindex'][:]

nifty/spaces/rg_space/rg_space.py

@@ -110,7 +110,7 @@ class RGSpace(Space):
         hermitian_part /= 2.
 
         # use subtraction since it is faster than flipping another time
-        anti_hermitian_part = (x-hermitian_part)/1j
+        anti_hermitian_part = (x-hermitian_part)
 
         if preserve_gaussian_variance:
             hermitian_part, anti_hermitian_part = \
@@ -120,6 +120,18 @@ class RGSpace(Space):
 
         return (hermitian_part, anti_hermitian_part)
 
+    def hermitian_fixed_points(self):
+        shape = self.shape
+        mid_index = np.array(shape)//2
+        ndlist = [2 if (shape[i] % 2 == 0) else 1 for i in xrange(len(shape))]
+        ndlist = tuple(ndlist)
+        odd_axes_list = np.array([1 if (shape[i] % 2 == 1) else 0
+                                  for i in xrange(len(shape))])
+        fixed_points = []
+        for i in np.ndindex(ndlist):
+            fixed_points += [tuple((i+odd_axes_list) * mid_index)]
+        return fixed_points
+
     def _hermitianize_correct_variance(self, hermitian_part,
                                        anti_hermitian_part, axes):
         # Correct the variance by multiplying sqrt(2)
@@ -145,8 +157,9 @@ class RGSpace(Space):
         return hermitian_part, anti_hermitian_part
 
     def _hermitianize_inverter(self, x, axes):
+        shape = x.shape
         # calculate the number of dimensions the input array has
-        dimensions = len(x.shape)
+        dimensions = len(shape)
         # prepare the slicing object which will be used for mirroring
         slice_primitive = [slice(None), ] * dimensions
         # copy the input data
@@ -158,11 +171,17 @@ class RGSpace(Space):
         # flip in the desired directions
         for i in axes:
             slice_picker = slice_primitive[:]
-            slice_picker[i] = slice(1, None, None)
+            if shape[i] % 2 == 0:
+                slice_picker[i] = slice(1, None, None)
+            else:
+                slice_picker[i] = slice(None)
             slice_picker = tuple(slice_picker)
 
             slice_inverter = slice_primitive[:]
-            slice_inverter[i] = slice(None, 0, -1)
+            if shape[i] % 2 == 0:
+                slice_inverter[i] = slice(None, 0, -1)
+            else:
+                slice_inverter[i] = slice(None, None, -1)
             slice_inverter = tuple(slice_inverter)
 
             try:

nifty/spaces/space/space.py

@@ -58,7 +58,7 @@ class Space(DomainObject):
     """
 
     def __init__(self):
-        
+
         super(Space, self).__init__()
 
     @abc.abstractproperty

test/test_serialization.py

@@ -21,31 +21,18 @@ import unittest
 from numpy.testing import assert_equal
 from keepers import Repository
 from test.common import expand, generate_spaces
+from nose.plugins.skip import SkipTest
 import os
 
-try:
-    import h5py
-except ImportError:
-    h5py_available = False
-else:
-    h5py_available = True
-
-if h5py_available:
-    class SpaceSerializationTests(unittest.TestCase):
-        # variable to hold the repository
-        _repo = None
-
-        @classmethod
-        def setUpClass(cls):
-            cls._repo = Repository('test.h5')
-
-        @expand([[space] for space in generate_spaces()])
-        def test_serialization(self, space):
-            self._repo.add(space, 'space')
-            self._repo.commit()
-
-            assert_equal(space, self._repo.get('space'))
-
-        @classmethod
-        def tearDownClass(cls):
-            os.remove('test.h5')
+class SpaceSerializationTests(unittest.TestCase):
+    @expand([[space] for space in generate_spaces()])
+    def test_serialization(self, space):
+        try:
+            import h5py
+        except ImportError:
+            raise SkipTest
+        repo = Repository('test.h5')
+        repo.add(space, 'space')
+        repo.commit()
+        assert_equal(space, repo.get('space'))
+        os.remove('test.h5')

test/test_spaces/test_lm_space.py

@@ -116,7 +116,7 @@ class LMSpaceFunctionalityTests(unittest.TestCase):
             real)
         assert_almost_equal(
             l.hermitian_decomposition(distributed_data_object(x))[1],
-            imag)
+            imag*1j)
 
     @expand(get_weight_configs())
     def test_weight(self, x, power, axes, inplace, expected):

test/test_spaces/test_rg_space.py

@@ -162,7 +162,7 @@ def get_hermitian_configs():
         [0.37928780+0.j, 0.33013206+0.26511434j, 0.48235714+0.j,
             0.33013206-0.26511434j],
         [0.47773437+0.17342852j, 0.31059374+0.02379381j, 0.64003551-0.23838932j,
-            0.27824437-0.0197826j]])
+            0.27824437-0.0197826j]])*1j
 
     h_1_x = np.array([
         [[0.23987021+0.41617749j, 0.34605012+0.55462234j, 0.07947035+0.73360723j,
@@ -205,7 +205,7 @@ def get_hermitian_configs():
              0.64466122-0.10318736j],
          [0.19076379+0.j, 0.49389371+0.03664104j, 0.85091112+0.j,
              0.49389371-0.03664104j]]
-        ])
+        ])*1j
     return [
         [h_0_x, None, h_0_res_real, h_0_res_imag],
         [h_1_x, (2,), h_1_res_real, h_1_res_imag]