diff --git a/demos/critical_filtering.py b/demos/critical_filtering.py
index 4e96544086c60a1de8efa6ba0459ef62e060bb31..d37984f401d61768efcf3bed5e2b354601e4455c 100644
--- a/demos/critical_filtering.py
+++ b/demos/critical_filtering.py
@@ -51,7 +51,7 @@ if __name__ == "__main__":
     data_power = ift.log(ift.power_analyze(fft.adjoint_times(d),
                                            binbounds=p_space.binbounds))
     d_data = d.val
-    ift.plotting.plot(d, name="data.png")
+    ift.plot(d, name="data.png")
 
     IC1 = ift.GradientNormController(verbose=True, iteration_limit=100,
                                      tol_abs_gradnorm=0.1)
@@ -90,5 +90,5 @@ if __name__ == "__main__":
 
         # Plot current estimate
         ift.dobj.mprint(i)
-        if i % 5 == 0:
-            ift.plotting.plot(fft(m0), name='map.png')
+        if i % 50 == 0:
+            ift.plot(fft(m0), name='map.png')
diff --git a/demos/log_normal_wiener_filter.py b/demos/log_normal_wiener_filter.py
index 94a79a52ef0cf9f67a6a591d8251f0eb05f60b93..87563debbfae8a130727a7f2aaa9673ee30fb759 100644
--- a/demos/log_normal_wiener_filter.py
+++ b/demos/log_normal_wiener_filter.py
@@ -68,14 +68,13 @@ if __name__ == "__main__":
     # Plotting
     plotdict = {"xlabel": "Pixel index", "ylabel": "Pixel index",
                 "colormap": "Planck-like"}
-    ift.plotting.plot(mock_signal, name="mock_signal.png", **plotdict)
+    ift.plot(mock_signal, name="mock_signal.png", **plotdict)
     logdata = np.log(ift.dobj.to_global_data(data.val)).reshape(signal_space.shape)
-    ift.plotting.plot(ift.Field(signal_space,
-                                val=ift.dobj.from_global_data(logdata)),
-                      name="log_of_data.png", **plotdict)
-    # ift.plotting.plot(m1,name='m_LBFGS.png', **plotdict)
-    ift.plotting.plot(m2, name='m_Newton.png', **plotdict)
-    # ift.plotting.plot(m3, name='m_SteepestDescent.png', **plotdict)
+    ift.plot(ift.Field(signal_space, val=ift.dobj.from_global_data(logdata)),
+             name="log_of_data.png", **plotdict)
+    # ift.plot(m1,name='m_LBFGS.png', **plotdict)
+    ift.plot(m2, name='m_Newton.png', **plotdict)
+    # ift.plot(m3, name='m_SteepestDescent.png', **plotdict)
 
     # Probing the variance
     class Proby(ift.DiagonalProberMixin, ift.Prober):
@@ -85,4 +84,4 @@ if __name__ == "__main__":
 
     sm = ift.FFTSmoothingOperator(signal_space, sigma=0.02)
     variance = sm(proby.diagonal.weight(-1))
-    ift.plotting.plot(variance, name='variance.png', **plotdict)
+    ift.plot(variance, name='variance.png', **plotdict)
diff --git a/demos/nonlinear_critical_filter.py b/demos/nonlinear_critical_filter.py
index 2d9e092f7dcf12cd9edb54b3b53c4e7a5696e2ae..7feef7420d90a66f0c4903e5aea20e3fb0e65031 100644
--- a/demos/nonlinear_critical_filter.py
+++ b/demos/nonlinear_critical_filter.py
@@ -110,7 +110,7 @@ if __name__ == "__main__":
         # excitation monopole to 1
         m0, t0 = adjust_zero_mode(m0, t0)
 
-    ift.plotting.plot(true_sky)
-    ift.plotting.plot(nonlinearity(FFT.adjoint_times(power0*m0)),
-                      title='reconstructed_sky')
-    ift.plotting.plot(MeasurementOperator.adjoint_times(d))
+    ift.plot(true_sky)
+    ift.plot(nonlinearity(FFT.adjoint_times(power0*m0)),
+             title='reconstructed_sky')
+    ift.plot(MeasurementOperator.adjoint_times(d))
diff --git a/demos/paper_demos/cartesian_wiener_filter.py b/demos/paper_demos/cartesian_wiener_filter.py
index e2b7eec5906536665780f6d3c60e4ce62deaa955..4fd696d093a656d10289b7347925e6c376fb49c5 100644
--- a/demos/paper_demos/cartesian_wiener_filter.py
+++ b/demos/paper_demos/cartesian_wiener_filter.py
@@ -107,13 +107,12 @@ if __name__ == "__main__":
     sm = ift.FFTSmoothingOperator(plot_space, sigma=0.03)
     plotdict = {"xlabel": "Pixel index", "ylabel": "Pixel index",
                 "colormap": "Planck-like"}
-    ift.plotting.plot(
+    ift.plot(
         ift.log(ift.sqrt(sm(ift.Field(plot_space, val=variance.val.real)))),
         name='uncertainty.png', zmin=0., zmax=3., title="Uncertainty map",
         **plotdict)
-    ift.plotting.plot(ift.Field(plot_space, val=mock_signal.val.real),
-                      name='mock_signal.png', **plotdict)
-    ift.plotting.plot(ift.Field(plot_space, val=data.val.real),
-                      name='data.png', **plotdict)
-    ift.plotting.plot(ift.Field(plot_space, val=m.val.real),
-                      name='map.png', **plotdict)
+    ift.plot(ift.Field(plot_space, val=mock_signal.val.real),
+             name='mock_signal.png', **plotdict)
+    ift.plot(ift.Field(plot_space, val=data.val.real),
+             name='data.png', **plotdict)
+    ift.plot(ift.Field(plot_space, val=m.val.real), name='map.png', **plotdict)
diff --git a/demos/paper_demos/wiener_filter.py b/demos/paper_demos/wiener_filter.py
index 04b131b634d7fbffb5ef0ce60f3e7b18e9b15186..ddf5a63a74468b782716bb5b21fc60d327046309 100644
--- a/demos/paper_demos/wiener_filter.py
+++ b/demos/paper_demos/wiener_filter.py
@@ -69,8 +69,8 @@ if __name__ == "__main__":
     # Plotting
     plotdict = {"xlabel": "Pixel index", "ylabel": "Pixel index",
                 "colormap": "Planck-like"}
-    ift.plotting.plot(variance, name="uncertainty.png", **plotdict)
-    ift.plotting.plot(mock_signal, name="mock_signal.png", **plotdict)
-    ift.plotting.plot(ift.Field(signal_space, val=data.val),
-                      name="data.png", **plotdict)
-    ift.plotting.plot(m, name="map.png", **plotdict)
+    ift.plot(variance, name="uncertainty.png", **plotdict)
+    ift.plot(mock_signal, name="mock_signal.png", **plotdict)
+    ift.plot(ift.Field(signal_space, val=data.val),
+             name="data.png", **plotdict)
+    ift.plot(m, name="map.png", **plotdict)
diff --git a/demos/wiener_filter_via_curvature.py b/demos/wiener_filter_via_curvature.py
index eb38204a0bee6eb23a6f91d02511082134cd3348..b32f05632704956adc804de26d5ac431e12933aa 100644
--- a/demos/wiener_filter_via_curvature.py
+++ b/demos/wiener_filter_via_curvature.py
@@ -77,9 +77,8 @@ if __name__ == "__main__":
 
     sspace2 = ift.RGSpace(shape, distances=L/N_pixels/nu.m)
 
-    ift.plotting.plot(ift.Field(sspace2, mock_signal.val)/nu.K,
-                      name="mock_signal.png")
+    ift.plot(ift.Field(sspace2, mock_signal.val)/nu.K, name="mock_signal.png")
     data = ift.dobj.to_global_data(data.val).reshape(sspace2.shape)/nu.K
     data = ift.Field(sspace2, val=ift.dobj.from_global_data(data))/nu.K
-    ift.plotting.plot(ift.Field(sspace2, val=data), name="data.png")
-    ift.plotting.plot(ift.Field(sspace2, m_s.val)/nu.K, name="map.png")
+    ift.plot(ift.Field(sspace2, val=data), name="data.png")
+    ift.plot(ift.Field(sspace2, m_s.val)/nu.K, name="map.png")
diff --git a/demos/wiener_filter_via_hamiltonian.py b/demos/wiener_filter_via_hamiltonian.py
index d45a2ee4ce54a5df6cf7070ca8fb9f887a01fae9..b1d62dbbe9a4b7756e1e4e4c0809e85600d682dd 100644
--- a/demos/wiener_filter_via_hamiltonian.py
+++ b/demos/wiener_filter_via_hamiltonian.py
@@ -62,8 +62,8 @@ if __name__ == "__main__":
     energy, convergence = minimizer(energy)
     m = energy.position
     D = energy.curvature
-    ift.plotting.plot(ss, name="signal.png", colormap="Planck-like")
-    ift.plotting.plot(fft(m), name="m.png", colormap="Planck-like")
+    ift.plot(ss, name="signal.png", colormap="Planck-like")
+    ift.plot(fft(m), name="m.png", colormap="Planck-like")
 
     # sampling the uncertainty map
     sample_variance = ift.Field.zeros(s_space)
@@ -77,4 +77,4 @@ if __name__ == "__main__":
     sample_mean /= n_samples
     sample_variance /= n_samples
     variance = sample_variance - sample_mean**2
-    ift.plotting.plot(variance, name="variance.png", colormap="Planck-like")
+    ift.plot(variance, name="variance.png", colormap="Planck-like")
diff --git a/nifty4/__init__.py b/nifty4/__init__.py
index 8345aae2b7d1c06f2b32eca5a6d7f83aac2b64bb..199905dcf63b7cf7bd0ae5ce27a27e1bd39e7129 100644
--- a/nifty4/__init__.py
+++ b/nifty4/__init__.py
@@ -9,6 +9,6 @@ from .spaces import *
 from .operators import *
 from .probing import *
 from .sugar import *
-from . import plotting
+from .plotting import plot
 from . import library
 from . import dobj
diff --git a/nifty4/domain_object.py b/nifty4/domain_object.py
index 85e107fe5ceec02cd1d51a2689b23e06e54309f5..1e57592b0218f4d26250492682dc5cd0a92edb01 100644
--- a/nifty4/domain_object.py
+++ b/nifty4/domain_object.py
@@ -114,7 +114,4 @@ class DomainObject(with_metaclass(
     @property
     def total_volume(self):
         tmp = self.dvol()
-        if np.isscalar(tmp):
-            return self.dim * tmp
-        else:
-            return np.sum(tmp)
+        return self.dim * tmp if np.isscalar(tmp) else np.sum(tmp)
diff --git a/nifty4/field.py b/nifty4/field.py
index b57d80a2a50550533ad30ee5e3c26f3372ebb4de..ee6ad6d797b2c74906b2a60d2b349d49f6de0f7b 100644
--- a/nifty4/field.py
+++ b/nifty4/field.py
@@ -541,7 +541,7 @@ class Field(object):
         minmax = [self.min(), self.max()]
         mean = self.mean()
         return "nifty4.Field instance\n- domain      = " + \
-               repr(self._domain) + \
+               self._domain.__str__() + \
                "\n- val         = " + repr(self.val) + \
                "\n  - min.,max. = " + str(minmax) + \
                "\n  - mean = " + str(mean)