merge master; currently broken

.gitignore

@@ -10,6 +10,8 @@ setup.cfg
 .document
 .svn/
 *.csv
+.pytest_cache/
+*.png
 
 # from https://github.com/github/gitignore/blob/master/Python.gitignore
 

.gitlab-ci.yml

 image: $CONTAINER_TEST_IMAGE
 
 variables:
-  CONTAINER_TEST_IMAGE: gitlab-registry.mpcdf.mpg.de/ift/nifty:$CI_BUILD_REF_NAME
+  CONTAINER_TEST_IMAGE: gitlab-registry.mpcdf.mpg.de/$CI_PROJECT_PATH:$CI_BUILD_REF_NAME
   OMP_NUM_THREADS: 1
 
 stages:
@@ -15,6 +15,8 @@ build_docker_from_scratch:
     - schedules
   image: docker:stable
   stage: build_docker
+  before_script:
+    - ls
   script:
     - docker login -u gitlab-ci-token -p $CI_BUILD_TOKEN gitlab-registry.mpcdf.mpg.de
     - docker build -t $CONTAINER_TEST_IMAGE --no-cache .
@@ -25,33 +27,32 @@ build_docker_from_cache:
     - schedules
   image: docker:stable
   stage: build_docker
+  before_script:
+    - ls
   script:
     - docker login -u gitlab-ci-token -p $CI_BUILD_TOKEN gitlab-registry.mpcdf.mpg.de
     - docker build -t $CONTAINER_TEST_IMAGE .
     - docker push $CONTAINER_TEST_IMAGE
 
-test_python2_with_coverage:
+test_serial:
   stage: test
   script:
-    - python setup.py install --user -f
-    - mpiexec -n 2 --bind-to none nosetests -q 2> /dev/null
-    - nosetests -q --with-coverage --cover-package=nifty5 --cover-erase
+    - pytest-3 -q --cov=nifty5 test
     - >
-      coverage report --omit "*plotting*,*distributed_do*"
+      python3 -m coverage report --omit "*plot*,*distributed_do*" | tee coverage.txt
     - >
-      coverage report --omit "*plotting*,*distributed_do*" | grep TOTAL | awk '{ print "TOTAL: "$4; }'
+      grep TOTAL coverage.txt | awk '{ print "TOTAL: "$4; }'
 
-test_python3:
+test_mpi:
   stage: test
+  variables:
+    OMPI_MCA_btl_vader_single_copy_mechanism: none
   script:
-    - python3 setup.py install --user -f
-    - mpiexec -n 2 --bind-to none nosetests3 -q 2> /dev/null
-    - nosetests3 -q
+    - mpiexec -n 2 --bind-to none pytest-3 -q test
 
 pages:
   stage: release
   script:
-    - python setup.py install --user -f
     - sh docs/generate.sh
     - mv docs/build/ public/
   artifacts:
@@ -60,142 +61,53 @@ pages:
   only:
   - NIFTy_5
 
-before_script:
-  - export MPLBACKEND="agg"
-
-run_critical_filtering:
-  stage: demo_runs
-  script:
-    - ls
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - python demos/critical_filtering.py
-    - python3 demos/critical_filtering.py
-  artifacts:
-    paths:
-      - '*.png'
-
-run_nonlinear_critical_filter:
-  stage: demo_runs
-  script:
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - python demos/nonlinear_critical_filter.py
-    - python3 demos/nonlinear_critical_filter.py
-  artifacts:
-    paths:
-      - '*.png'
 
-run_nonlinear_wiener_filter:
-  stage: demo_runs
-  script:
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - python demos/nonlinear_wiener_filter.py
-    - python3 demos/nonlinear_wiener_filter.py
-  only:
-    - run_demos
-  artifacts:
-    paths:
-      - '*.png'
-
-# FIXME: disable for now. Fixing it is part of issue #244.
-#run_poisson_demo:
-#  stage: demo_runs
-#  script:
-#    - python setup.py install --user -f
-#    - python3 setup.py install --user -f
-#    - python demos/poisson_demo.py
-#    - python3 demos/poisson_demo.py
-#  artifacts:
-#    paths:
-#      - '*.png'
-
-run_probing:
-  stage: demo_runs
-  script:
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - python demos/probing.py
-    - python3 demos/probing.py
-  artifacts:
-    paths:
-      - '*.png'
-
-run_sampling:
-  stage: demo_runs
-  script:
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - python demos/sampling.py
-    - python3 demos/sampling.py
-  artifacts:
-    paths:
-      - '*.png'
+before_script:
+  - python3 setup.py install --user -f
 
-run_tomography:
+run_ipynb:
   stage: demo_runs
   script:
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - python demos/tomography.py
-    - python3 demos/tomography.py
-  artifacts:
-    paths:
-      - '*.png'
+    - jupyter nbconvert --execute --ExecutePreprocessor.timeout=None demos/Wiener_Filter.ipynb
 
-run_wiener_filter_data_space_noiseless:
+run_getting_started_1:
   stage: demo_runs
   script:
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - python demos/wiener_filter_data_space_noiseless.py
-    - python3 demos/wiener_filter_data_space_noiseless.py
+    - python3 demos/getting_started_1.py
+    - mpiexec -n 2 --bind-to none python3 demos/getting_started_1.py 2> /dev/null
   artifacts:
     paths:
       - '*.png'
 
-run_wiener_filter_easy.py:
+run_getting_started_2:
   stage: demo_runs
   script:
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - python demos/wiener_filter_easy.py
-    - python3 demos/wiener_filter_easy.py
+    - python3 demos/getting_started_2.py
+    - mpiexec -n 2 --bind-to none python3 demos/getting_started_2.py 2> /dev/null
   artifacts:
     paths:
       - '*.png'
 
-run_wiener_filter_via_curvature.py:
+run_getting_started_3:
   stage: demo_runs
   script:
-    - pip install --user numericalunits
-    - pip3 install --user numericalunits
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - python demos/wiener_filter_via_curvature.py
-    - python3 demos/wiener_filter_via_curvature.py
+    - python3 demos/getting_started_3.py
   artifacts:
     paths:
       - '*.png'
 
-run_wiener_filter_via_hamiltonian.py:
+run_bernoulli:
   stage: demo_runs
   script:
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - python demos/wiener_filter_via_hamiltonian.py
-    - python3 demos/wiener_filter_via_hamiltonian.py
+    - python3 demos/bernoulli_demo.py
   artifacts:
     paths:
       - '*.png'
 
-run_ipynb:
+run_curve_fitting:
   stage: demo_runs
   script:
-    - python setup.py install --user -f
-    - python3 setup.py install --user -f
-    - jupyter nbconvert --execute --ExecutePreprocessor.timeout=None demos/Wiener_Filter.ipynb
+    - python3 demos/polynomial_fit.py
   artifacts:
     paths:
       - '*.png'

Dockerfile

 FROM debian:testing-slim
 
 RUN apt-get update && apt-get install -y \
-    # Needed for gitlab tests
-    git \
+    # Needed for setup
+    git python3-pip \
     # Packages needed for NIFTy
-    libfftw3-dev \
-    python python-pip python-dev python-future python-scipy \
-    python3 python3-pip python3-dev python3-future python3-scipy \
+    python3-scipy \
     # Documentation build dependencies
-    python-sphinx python-sphinx-rtd-theme python-numpydoc \
+    python3-sphinx-rtd-theme dvipng texlive-latex-base texlive-latex-extra \
     # Testing dependencies
-    python-nose python-parameterized \
-    python3-nose python3-parameterized \
+    python3-pytest-cov jupyter \
     # Optional NIFTy dependencies
-    openmpi-bin libopenmpi-dev python-mpi4py python3-mpi4py \
-    # Packages needed for NIFTy
-  && pip install pyfftw \
+    libfftw3-dev python3-mpi4py python3-matplotlib \
+  # more optional NIFTy dependencies
   && pip3 install pyfftw \
-  # Optional NIFTy dependencies
-  && pip install git+https://gitlab.mpcdf.mpg.de/ift/pyHealpix.git \
   && pip3 install git+https://gitlab.mpcdf.mpg.de/ift/pyHealpix.git \
-  # Testing dependencies
-  && pip install coverage \
+  && pip3 install jupyter \
   && rm -rf /var/lib/apt/lists/*
 
-# Needed for demos to be running
-RUN apt-get update && apt-get install -y python-matplotlib python3-matplotlib \
-  && python3 -m pip install --upgrade pip && python3 -m pip install jupyter && python -m pip install --upgrade pip && python -m pip install jupyter \
-  && rm -rf /var/lib/apt/lists/*
+# Set matplotlib backend
+ENV MPLBACKEND agg
 
 # Create user (openmpi does not like to be run as root)
 RUN useradd -ms /bin/bash testinguser

README.md

@@ -4,7 +4,7 @@ NIFTy - Numerical Information Field Theory
 [![coverage report](https://gitlab.mpcdf.mpg.de/ift/NIFTy/badges/NIFTy_5/coverage.svg)](https://gitlab.mpcdf.mpg.de/ift/NIFTy/commits/NIFTy_5)
 
 **NIFTy** project homepage:
-[http://ift.pages.mpcdf.de/NIFTy](http://ift.pages.mpcdf.de/NIFTy)
+[http://ift.pages.mpcdf.de/nifty](http://ift.pages.mpcdf.de/nifty)
 
 Summary
 -------
@@ -13,23 +13,31 @@ Summary
 
 **NIFTy**, "**N**umerical **I**nformation **F**ield **T**heor<strong>y</strong>", is
 a versatile library designed to enable the development of signal
-inference algorithms that operate regardless of the underlying spatial
-grid and its resolution. Its object-oriented framework is written in
-Python, although it accesses libraries written in C++ and C for
-efficiency.
+inference algorithms that operate regardless of the underlying grids
+(spatial, spectral, temporal, …) and their resolutions.
+Its object-oriented framework is written in Python, although it accesses
+libraries written in C++ and C for efficiency.
 
 NIFTy offers a toolkit that abstracts discretized representations of
 continuous spaces, fields in these spaces, and operators acting on
-fields into classes. The correct normalization of operations on
-fields is taken care of automatically without concerning the user. This
-allows for an abstract formulation and programming of inference
+these fields into classes.
+This allows for an abstract formulation and programming of inference
 algorithms, including those derived within information field theory.
-Thus, NIFTy permits its user to rapidly prototype algorithms in 1D, and
-then apply the developed code in higher-dimensional settings of real
-world problems. The set of spaces on which NIFTy operates comprises
-point sets, *n*-dimensional regular grids, spherical spaces, their
-harmonic counterparts, and product spaces constructed as combinations of
-those.
+NIFTy's interface is designed to resemble IFT formulae in the sense
+that the user implements algorithms in NIFTy independent of the topology
+of the underlying spaces and the discretization scheme.
+Thus, the user can develop algorithms on subsets of problems and on
+spaces where the detailed performance of the algorithm can be properly
+evaluated and then easily generalize them to other, more complex spaces
+and the full problem, respectively.
+
+The set of spaces on which NIFTy operates comprises point sets,
+*n*-dimensional regular grids, spherical spaces, their harmonic
+counterparts, and product spaces constructed as combinations of those.
+NIFTy takes care of numerical subtleties like the normalization of
+operations on fields and the numerical representation of model
+components, allowing the user to focus on formulating the abstract
+inference procedures and process-specific model properties.
 
 
 Installation
@@ -37,11 +45,11 @@ Installation
 
 ### Requirements
 
-- [Python](https://www.python.org/) (v2.7.x or 3.5.x)
+- [Python 3](https://www.python.org/) (3.5.x or later)
 - [SciPy](https://www.scipy.org/)
-- [pyFFTW](https://pypi.python.org/pypi/pyFFTW)
 
 Optional dependencies:
+- [pyFFTW](https://pypi.python.org/pypi/pyFFTW) for faster Fourier transforms
 - [pyHealpix](https://gitlab.mpcdf.mpg.de/ift/pyHealpix) (for harmonic
     transforms involving domains on the sphere)
 - [mpi4py](https://mpi4py.scipy.org) (for MPI-parallel execution)
@@ -61,57 +69,60 @@ distributions, the "apt" lines will need slight changes.
 
 NIFTy5 and its mandatory dependencies can be installed via:
 
-    sudo apt-get install git libfftw3-dev python python-pip python-dev
-    pip install --user git+https://gitlab.mpcdf.mpg.de/ift/NIFTy.git@NIFTy_5
+    sudo apt-get install git python3 python3-pip python3-dev
+    pip3 install --user git+https://gitlab.mpcdf.mpg.de/ift/NIFTy.git@NIFTy_5
+
+Plotting support is added via:
+
+    sudo apt-get install python3-matplotlib
+
+NIFTy uses Numpy's FFT implementation by default. For large problems FFTW may be
+used because of its higher performance. It can be installed via:
+
+    sudo apt-get install libfftw3-dev
+    pip3 install --user pyfftw
+
+To enable FFTW usage in NIFTy, call
+
+    nifty5.fft.enable_fftw()
+
+at the beginning of your code.
 
 (Note: If you encounter problems related to `pyFFTW`, make sure that you are
 using a pip-installed `pyFFTW` package. Unfortunately, some distributions are
 shipping an incorrectly configured `pyFFTW` package, which does not cooperate
 with the installed `FFTW3` libraries.)
 
-Plotting support is added via:
-
-    pip install --user matplotlib
-
 Support for spherical harmonic transforms is added via:
 
-    pip install --user git+https://gitlab.mpcdf.mpg.de/ift/pyHealpix.git
+    pip3 install --user git+https://gitlab.mpcdf.mpg.de/ift/pyHealpix.git
 
 MPI support is added via:
 
-    sudo apt-get install openmpi-bin libopenmpi-dev
-    pip install --user mpi4py
-
-### Installation for Python 3
-
-If you want to run NIFTy with Python 3, you need to make the following changes
-to the instructions above:
-
-- in all `apt-get` commands, replace `python-*` by `python3-*`
-- in all `pip` commands, replace `pip` by `pip3`
+    sudo apt-get install python3-mpi4py
 
 ### Running the tests
 
-In oder to run the tests one needs two additional packages:
+To run the tests, additional packages are required:
 
-    pip install --user nose parameterized coverage
+    sudo apt-get install python3-pytest-cov
 
 Afterwards the tests (including a coverage report) can be run using the
 following command in the repository root:
 
-    nosetests -x --with-coverage --cover-html --cover-package=nifty5
+    pytest-3 --cov=nifty5 test
 
 
 ### First Steps
 
 For a quick start, you can browse through the [informal
-introduction](http://ift.pages.mpcdf.de/NIFTy/code.html) or
+introduction](http://ift.pages.mpcdf.de/nifty/code.html) or
 dive into NIFTy by running one of the demonstrations, e.g.:
 
-    python demos/wiener_filter_via_curvature.py
+    python3 demos/getting_started_1.py
 
 
-### Acknowledgement
+### Acknowledgements
 
 Please acknowledge the use of NIFTy in your publication(s) by using a
 phrase such as the following:
@@ -119,10 +130,10 @@ phrase such as the following:
 > "Some of the results in this publication have been derived using the
 > NIFTy package [(https://gitlab.mpcdf.mpg.de/ift/NIFTy)](https://gitlab.mpcdf.mpg.de/ift/NIFTy)"
 
-and a citation to one of the [publications](http://ift.pages.mpcdf.de/NIFTy/citations.html).
+and a citation to one of the [publications](http://ift.pages.mpcdf.de/nifty/citations.html).
 
 
-### Release Notes
+### Licensing terms
 
 The NIFTy package is licensed under the terms of the
 [GPLv3](https://www.gnu.org/licenses/gpl.html) and is distributed

demos/Wiener_Filter.ipynb

@@ -58,7 +58,7 @@
     "### Posterior\n",
     "The Posterior is given by:\n",
     "\n",
-    "$$\\mathcal P (s|d) \\propto P(s,d) = \\mathcal G(d-Rs,N) \\,\\mathcal G(s,S) \\propto \\mathcal G (m,D) $$\n",
+    "$$\\mathcal P (s|d) \\propto P(s,d) = \\mathcal G(d-Rs,N) \\,\\mathcal G(s,S) \\propto \\mathcal G (s-m,D) $$\n",
     "\n",
     "where\n",
     "$$\\begin{align}\n",
@@ -172,7 +172,7 @@
     "                                    tol_abs_gradnorm=0.1)\n",
     "    # WienerFilterCurvature is (R.adjoint*N.inverse*R + Sh.inverse) plus some handy\n",
     "    # helper methods.\n",
-    "    return ift.library.WienerFilterCurvature(R,N,Sh,iteration_controller=IC,iteration_controller_sampling=IC)"
+    "    return ift.WienerFilterCurvature(R,N,Sh,iteration_controller=IC,iteration_controller_sampling=IC)"
    ]
   },
   {
@@ -429,8 +429,8 @@
     "mask[l:h] = 0\n",
     "mask = ift.Field.from_global_data(s_space, mask)\n",
     "\n",
-    "R = ift.DiagonalOperator(mask)*HT\n",
-    "n = n.to_global_data()\n",
+    "R = ift.DiagonalOperator(mask)(HT)\n",
+    "n = n.to_global_data_rw()\n",
     "n[l:h] = 0\n",
     "n = ift.Field.from_global_data(s_space, n)\n",
     "\n",
@@ -501,7 +501,7 @@
     "m_data = HT(m).to_global_data()\n",
     "m_var_data = m_var.to_global_data()\n",
     "uncertainty = np.sqrt(m_var_data)\n",
-    "d_data = d.to_global_data()\n",
+    "d_data = d.to_global_data_rw()\n",
     "\n",
     "# Set lost data to NaN for proper plotting\n",
     "d_data[d_data == 0] = np.nan"
@@ -585,8 +585,8 @@
     "mask[l:h,l:h] = 0.\n",
     "mask = ift.Field.from_global_data(s_space, mask)\n",
     "\n",
-    "R = ift.DiagonalOperator(mask)*HT\n",
-    "n = n.to_global_data()\n",
+    "R = ift.DiagonalOperator(mask)(HT)\n",
+    "n = n.to_global_data_rw()\n",
     "n[l:h, l:h] = 0\n",
     "n = ift.Field.from_global_data(s_space, n)\n",
     "curv = Curvature(R=R, N=N, Sh=Sh)\n",
@@ -717,21 +717,21 @@
  "metadata": {
   "celltoolbar": "Slideshow",
   "kernelspec": {
-   "display_name": "Python 2",
+   "display_name": "Python 3",
    "language": "python",
-   "name": "python2"
+   "name": "python3"
   },
   "language_info": {
    "codemirror_mode": {
     "name": "ipython",
-    "version": 2
+    "version": 3
    },
    "file_extension": ".py",
    "mimetype": "text/x-python",
    "name": "python",
    "nbconvert_exporter": "python",
-   "pygments_lexer": "ipython2",
-   "version": "2.7.15"
+   "pygments_lexer": "ipython3",
+   "version": "3.6.6"
   }
  },
  "nbformat": 4,

demos/bernoulli_demo.py

+# This program is free software: you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation, either version 3 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+# GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License
+# along with this program.  If not, see <http://www.gnu.org/licenses/>.
+#
+# Copyright(C) 2013-2019 Max-Planck-Society
+#
+# NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik.
+
+#####################################################################
+# Bernoulli reconstruction
+# Reconstruct an event probability field with values between 0 and 1
+# from the observed events
+# 1D (set mode=0), 2D (mode=1), or on the sphere (mode=2)
+#####################################################################
+
+import numpy as np
+
+import nifty5 as ift
+
+if __name__ == '__main__':
+    np.random.seed(41)
+
+    # Set up the position space of the signal
+    mode = 2
+    if mode == 0:
+        # One-dimensional regular grid
+        position_space = ift.RGSpace(1024)
+    elif mode == 1:
+        # Two-dimensional regular grid
+        position_space = ift.RGSpace([512, 512])
+    else:
+        # Sphere
+        position_space = ift.HPSpace(128)
+
+    # Define harmonic space and transform
+    harmonic_space = position_space.get_default_codomain()
+    HT = ift.HarmonicTransformOperator(harmonic_space, position_space)
+
+    position = ift.from_random('normal', harmonic_space)
+
+    # Define power spectrum and amplitudes
+    def sqrtpspec(k):
+        return 1./(20. + k**2)
+
+    A = ift.create_power_operator(harmonic_space, sqrtpspec)
+
+    # Set up a sky operator and instrumental response
+    sky = ift.sigmoid(HT(A))
+    GR = ift.GeometryRemover(position_space)
+    R = GR
+
+    # Generate mock data
+    p = R(sky)
+    mock_position = ift.from_random('normal', harmonic_space)
+    tmp = p(mock_position).to_global_data().astype(np.float64)
+    data = np.random.binomial(1, tmp)
+    data = ift.Field.from_global_data(R.target, data)
+
+    # Compute likelihood and Hamiltonian
+    position = ift.from_random('normal', harmonic_space)
+    likelihood = ift.BernoulliEnergy(data)(p)
+    ic_newton = ift.DeltaEnergyController(
+        name='Newton', iteration_limit=100, tol_rel_deltaE=1e-8)
+    minimizer = ift.NewtonCG(ic_newton)
+    ic_sampling = ift.GradientNormController(iteration_limit=100)
+
+    # Minimize the Hamiltonian
+    H = ift.StandardHamiltonian(likelihood, ic_sampling)
+    H = ift.EnergyAdapter(position, H, want_metric=True)
+    # minimizer = ift.L_BFGS(ic_newton)
+    H, convergence = minimizer(H)
+
+    reconstruction = sky(H.position)
+
+    plot = ift.Plot()
+    plot.add(reconstruction, title='reconstruction')
+    plot.add(GR.adjoint_times(data), title='data')
+    plot.add(sky(mock_position), title='truth')
+    plot.output(nx=3, xsize=16, ysize=9, title="results", name="bernoulli.png")

demos/critical_filtering.py

-import nifty5 as ift
-from nifty5.library.nonlinearities import Linear, Tanh, Exponential
-import numpy as np
-np.random.seed(42)
-
-
-def adjust_zero_mode(m0, t0):
-    mtmp = m0.to_global_data().copy()
-    zero_position = len(m0.shape)*(0,)
-    zero_mode = mtmp[zero_position]
-    mtmp[zero_position] = zero_mode / abs(zero_mode)
-    ttmp = t0.to_global_data().copy()
-    ttmp[0] += 2 * np.log(abs(zero_mode))
-    return (ift.Field.from_global_data(m0.domain, mtmp),
-            ift.Field.from_global_data(t0.domain, ttmp))
-
-
-if __name__ == "__main__":
-    noise_level = 1.
-    p_spec = (lambda k: (.5 / (k + 1) ** 3))
-
-    nonlinearity = Linear()
-    # Set up position space
-    s_space = ift.RGSpace((128, 128))
-    h_space = s_space.get_default_codomain()
-
-    # Define harmonic transformation and associated harmonic space
-    HT = ift.HarmonicTransformOperator(h_space, target=s_space)
-
-    # Setting up power space
-    p_space = ift.PowerSpace(h_space,
-                             binbounds=ift.PowerSpace.useful_binbounds(
-                                 h_space, logarithmic=True))
-
-    # Choosing the prior correlation structure and defining
-    # correlation operator
-    p = ift.PS_field(p_space, p_spec)
-    log_p = ift.log(p)
-    S = ift.create_power_operator(h_space, power_spectrum=lambda k: 1e-5)
-
-    # Drawing a sample sh from the prior distribution in harmonic space
-    sh = S.draw_sample()
-
-    # Choosing the measurement instrument
-    # Instrument = SmoothingOperator(s_space, sigma=0.01)
-    mask = np.ones(s_space.shape)
-    # mask[6000:8000] = 0.
-    mask[30:70, 30:70] = 0.
-    mask = ift.Field.from_global_data(s_space, mask)
-
-    MaskOperator = ift.DiagonalOperator(mask)
-    R = ift.GeometryRemover(s_space)
-    R = R*MaskOperator
-    # R = R*HT
-    # R = R * ift.create_harmonic_smoothing_operator((harmonic_space,), 0,
-    #                                                response_sigma)
-    MeasurementOperator = R
-
-    d_space = MeasurementOperator.target
-
-    Distributor = ift.PowerDistributor(target=h_space, power_space=p_space)
-    power = Distributor(ift.exp(0.5*log_p))
-    # Creating the mock data
-    true_sky = nonlinearity(HT(power*sh))
-    noiseless_data = MeasurementOperator(true_sky)
-    noise_amplitude = noiseless_data.val.std()*noise_level
-    N = ift.ScalingOperator(noise_amplitude**2, d_space)
-    n = N.draw_sample()
-    # Creating the mock data
-    d = noiseless_data + n
-
-    m0 = ift.full(h_space, 1e-7)
-    t0 = ift.full(p_space, -4.)
-    power0 = Distributor.times(ift.exp(0.5 * t0))
-
-    plotdict = {"colormap": "Planck-like"}
-    zmin = true_sky.min()
-    zmax = true_sky.max()
-    ift.plot(true_sky, title="True sky", name="true_sky.png", **plotdict)
-    ift.plot(MeasurementOperator.adjoint_times(d), title="Data",
-             name="data.png", **plotdict)
-
-    IC1 = ift.GradientNormController(name="IC1", iteration_limit=100,
-                                     tol_abs_gradnorm=1e-3)
-    LS = ift.LineSearchStrongWolfe(c2=0.02)
-    minimizer = ift.RelaxedNewton(IC1,