diff --git a/nifty5/__init__.py b/nifty5/__init__.py
index 0e1b17e390d3aef8234b72eba6195f8ee67158cc..c8ebede34edd9407a9c77f16ae371db8286c5710 100644
--- a/nifty5/__init__.py
+++ b/nifty5/__init__.py
@@ -1,31 +1,102 @@
 from .version import __version__
 
 from . import dobj
-from .domains import *
+
+from .domains.domain import Domain
+from .domains.structured_domain import StructuredDomain
+from .domains.unstructured_domain import UnstructuredDomain
+from .domains.rg_space import RGSpace
+from .domains.lm_space import LMSpace
+from .domains.gl_space import GLSpace
+from .domains.hp_space import HPSpace
+from .domains.power_space import PowerSpace
+from .domains.dof_space import DOFSpace
+from .domains.log_rg_space import LogRGSpace
+
 from .domain_tuple import DomainTuple
 from .field import Field
 
 from .nonlinearities import Exponential, Linear, PositiveTanh, Tanh
 
-from .models import *
-from .operators import *
+from .models.constant import Constant
+from .models.linear_model import LinearModel
+from .models.local_nonlinearity import (LocalModel, PointwiseExponential,
+                                 PointwisePositiveTanh, PointwiseTanh)
+from .models.model import Model
+from .models.multi_model import MultiModel
+from .models.variable import Variable
+
+from .operators.diagonal_operator import DiagonalOperator
+from .operators.dof_distributor import DOFDistributor
+from .operators.domain_distributor import DomainDistributor
+from .operators.endomorphic_operator import EndomorphicOperator
+from .operators.exp_transform import ExpTransform
+from .operators.fft_operator import FFTOperator
+from .operators.fft_smoothing_operator import FFTSmoothingOperator
+from .operators.geometry_remover import GeometryRemover
+from .operators.harmonic_transform_operator import HarmonicTransformOperator
+from .operators.inversion_enabler import InversionEnabler
+from .operators.laplace_operator import LaplaceOperator
+from .operators.linear_operator import LinearOperator
+from .operators.mask_operator import MaskOperator
+from .operators.multi_adaptor import MultiAdaptor
+from .operators.power_distributor import PowerDistributor
+from .operators.qht_operator import QHTOperator
+from .operators.sampling_enabler import SamplingEnabler
+from .operators.sandwich_operator import SandwichOperator
+from .operators.scaling_operator import ScalingOperator
+from .operators.selection_operator import SelectionOperator
+from .operators.slope_operator import SlopeOperator
+from .operators.smoothness_operator import SmoothnessOperator
+from .operators.symmetrizing_operator import SymmetrizingOperator
+
 from .probing.utils import probe_with_posterior_samples, probe_diagonal, \
     StatCalculator
-from .minimization import *
+
+from .minimization.line_search import LineSearch
+from .minimization.line_search_strong_wolfe import LineSearchStrongWolfe
+from .minimization.iteration_controller import IterationController
+from .minimization.gradient_norm_controller import GradientNormController
+from .minimization.minimizer import Minimizer
+from .minimization.conjugate_gradient import ConjugateGradient
+from .minimization.nonlinear_cg import NonlinearCG
+from .minimization.descent_minimizer import DescentMinimizer
+from .minimization.steepest_descent import SteepestDescent
+from .minimization.vl_bfgs import VL_BFGS
+from .minimization.l_bfgs import L_BFGS
+from .minimization.relaxed_newton import RelaxedNewton
+from .minimization.scipy_minimizer import ScipyMinimizer, NewtonCG, L_BFGS_B, ScipyCG
+from .minimization.energy import Energy
+from .minimization.quadratic_energy import QuadraticEnergy
+from .minimization.line_energy import LineEnergy
+from .minimization.energy_sum import EnergySum
 
 from .sugar import *
 from .plotting.plot import plot
 
 from .library import *
+from .library.amplitude_model import make_amplitude_model
+from .library.apply_data import ApplyData
+from .library.gaussian_energy import GaussianEnergy
+from .library.los_response import LOSResponse
+from .library.point_sources import PointSources
+from .library.poissonian_energy import PoissonianEnergy
+from .library.wiener_filter_curvature import WienerFilterCurvature
+from .library.wiener_filter_energy import WienerFilterEnergy
+from .library.correlated_fields import make_correlated_field, make_mf_correlated_field
+
 from . import extra
 
 from .utilities import memo
 
 from .logger import logger
 
-from .multi import *
+from .multi.multi_domain import MultiDomain
+from .multi.multi_field import MultiField
+from .multi.block_diagonal_operator import BlockDiagonalOperator
 
-from .energies import *
+from .energies.kl import SampledKullbachLeiblerDivergence
+from .energies.hamiltonian import Hamiltonian
 
 # We deliberately don't set __all__ here, because we don't want people to do a
 # "from nifty5 import *"; that would swamp the global namespace.
diff --git a/nifty5/domains/__init__.py b/nifty5/domains/__init__.py
index 549be30cbb0357ef9dbe5e7307482bc64b12b164..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 100644
--- a/nifty5/domains/__init__.py
+++ b/nifty5/domains/__init__.py
@@ -1,14 +0,0 @@
-from .domain import Domain
-from .unstructured_domain import UnstructuredDomain
-from .structured_domain import StructuredDomain
-from .rg_space import RGSpace
-from .lm_space import LMSpace
-from .hp_space import HPSpace
-from .gl_space import GLSpace
-from .dof_space import DOFSpace
-from .power_space import PowerSpace
-from .log_rg_space import LogRGSpace
-
-__all__ = ["Domain", "UnstructuredDomain", "StructuredDomain", "RGSpace",
-           "LMSpace", "HPSpace", "GLSpace", "DOFSpace", "PowerSpace",
-           "LogRGSpace"]
diff --git a/nifty5/energies/__init__.py b/nifty5/energies/__init__.py
index 49f40904af27a9e3e787b351de63fe68302dfbdd..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 100644
--- a/nifty5/energies/__init__.py
+++ b/nifty5/energies/__init__.py
@@ -1,2 +0,0 @@
-from .kl import SampledKullbachLeiblerDivergence
-from .hamiltonian import Hamiltonian
diff --git a/nifty5/library/__init__.py b/nifty5/library/__init__.py
index 2678799c9b90ce4c4879a418a54102cfb5f504ae..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 100644
--- a/nifty5/library/__init__.py
+++ b/nifty5/library/__init__.py
@@ -1,9 +0,0 @@
-from .amplitude_model import make_amplitude_model
-from .apply_data import ApplyData
-from .gaussian_energy import GaussianEnergy
-from .los_response import LOSResponse
-from .point_sources import PointSources
-from .poissonian_energy import PoissonianEnergy
-from .wiener_filter_curvature import WienerFilterCurvature
-from .wiener_filter_energy import WienerFilterEnergy
-from .correlated_fields import make_correlated_field, make_mf_correlated_field
diff --git a/nifty5/library/point_sources.py b/nifty5/library/point_sources.py
index 2f5d8f1ab2f19795402e553ae3cc147a078f5617..1bb69c3ee1b2415631f8a4ff8e050c8384f5becc 100644
--- a/nifty5/library/point_sources.py
+++ b/nifty5/library/point_sources.py
@@ -2,10 +2,10 @@ import numpy as np
 from scipy.stats import invgamma, norm
 from ..field import Field
 from ..sugar import makeOp
-from ..multi import MultiField
-from ..models import Model
+from ..multi.multi_field import MultiField
+from ..models.model import Model
 
-from ..operators import SelectionOperator
+from ..operators.selection_operator import SelectionOperator
 from ..utilities import memo
 
 
diff --git a/nifty5/minimization/__init__.py b/nifty5/minimization/__init__.py
index 21616060311b8a2d3afc760e7334ab8803dcc5d5..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 100644
--- a/nifty5/minimization/__init__.py
+++ b/nifty5/minimization/__init__.py
@@ -1,24 +0,0 @@
-from .line_search import LineSearch
-from .line_search_strong_wolfe import LineSearchStrongWolfe
-from .iteration_controller import IterationController
-from .gradient_norm_controller import GradientNormController
-from .minimizer import Minimizer
-from .conjugate_gradient import ConjugateGradient
-from .nonlinear_cg import NonlinearCG
-from .descent_minimizer import DescentMinimizer
-from .steepest_descent import SteepestDescent
-from .vl_bfgs import VL_BFGS
-from .l_bfgs import L_BFGS
-from .relaxed_newton import RelaxedNewton
-from .scipy_minimizer import ScipyMinimizer, NewtonCG, L_BFGS_B, ScipyCG
-from .energy import Energy
-from .quadratic_energy import QuadraticEnergy
-from .line_energy import LineEnergy
-from .energy_sum import EnergySum
-
-__all__ = ["LineSearch", "LineSearchStrongWolfe",
-           "IterationController", "GradientNormController",
-           "Minimizer", "ConjugateGradient", "NonlinearCG", "DescentMinimizer",
-           "SteepestDescent", "VL_BFGS", "RelaxedNewton", "ScipyMinimizer",
-           "NewtonCG", "L_BFGS_B", "ScipyCG", "Energy", "QuadraticEnergy",
-           "LineEnergy", "L_BFGS", "EnergySum"]
diff --git a/nifty5/minimization/energy.py b/nifty5/minimization/energy.py
index 0104bbbf69439b6932d38b03b93cc0dfc06755da..aaf9ccf4a357873e0739dd3e3c0bfaebbcc73721 100644
--- a/nifty5/minimization/energy.py
+++ b/nifty5/minimization/energy.py
@@ -18,7 +18,6 @@
 
 import numpy as np
 from ..field import Field
-from ..multi import MultiField
 from ..utilities import NiftyMetaBase, memo
 
 
diff --git a/nifty5/models/__init__.py b/nifty5/models/__init__.py
index 6a3dd17e75cc605e01f25f22a95be48b2dbe7c8f..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 100644
--- a/nifty5/models/__init__.py
+++ b/nifty5/models/__init__.py
@@ -1,11 +0,0 @@
-from .constant import Constant
-from .linear_model import LinearModel
-from .local_nonlinearity import (LocalModel, PointwiseExponential,
-                                 PointwisePositiveTanh, PointwiseTanh)
-from .model import Model
-from .multi_model import MultiModel
-from .variable import Variable
-
-__all__ = ['Model', 'Constant', 'LocalModel', 'Variable',
-           'LinearModel', 'PointwiseTanh', 'PointwisePositiveTanh',
-           'PointwiseExponential', 'MultiModel']
diff --git a/nifty5/multi/__init__.py b/nifty5/multi/__init__.py
index f627ea25774fcf14d2d8e12bcbddbaa4fe21ac0e..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 100644
--- a/nifty5/multi/__init__.py
+++ b/nifty5/multi/__init__.py
@@ -1,5 +0,0 @@
-from .multi_domain import MultiDomain
-from .multi_field import MultiField
-from .block_diagonal_operator import BlockDiagonalOperator
-
-__all__ = ["MultiDomain", "MultiField", "BlockDiagonalOperator"]
diff --git a/nifty5/operators/__init__.py b/nifty5/operators/__init__.py
index 4cc51fbe602b941dbf07f2e681a007ae29dd1190..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 100644
--- a/nifty5/operators/__init__.py
+++ b/nifty5/operators/__init__.py
@@ -1,32 +0,0 @@
-from .diagonal_operator import DiagonalOperator
-from .dof_distributor import DOFDistributor
-from .domain_distributor import DomainDistributor
-from .endomorphic_operator import EndomorphicOperator
-from .exp_transform import ExpTransform
-from .fft_operator import FFTOperator
-from .fft_smoothing_operator import FFTSmoothingOperator
-from .geometry_remover import GeometryRemover
-from .harmonic_transform_operator import HarmonicTransformOperator
-from .inversion_enabler import InversionEnabler
-from .laplace_operator import LaplaceOperator
-from .linear_operator import LinearOperator
-from .mask_operator import MaskOperator
-from .multi_adaptor import MultiAdaptor
-from .power_distributor import PowerDistributor
-from .qht_operator import QHTOperator
-from .sampling_enabler import SamplingEnabler
-from .sandwich_operator import SandwichOperator
-from .scaling_operator import ScalingOperator
-from .selection_operator import SelectionOperator
-from .slope_operator import SlopeOperator
-from .smoothness_operator import SmoothnessOperator
-from .symmetrizing_operator import SymmetrizingOperator
-
-__all__ = ["LinearOperator", "EndomorphicOperator", "ScalingOperator",
-           "DiagonalOperator", "HarmonicTransformOperator", "FFTOperator",
-           "FFTSmoothingOperator", "GeometryRemover", "MaskOperator",
-           "LaplaceOperator", "SmoothnessOperator", "PowerDistributor",
-           "InversionEnabler", "SandwichOperator", "SamplingEnabler",
-           "DOFDistributor", "SelectionOperator", "MultiAdaptor",
-           "ExpTransform", "SymmetrizingOperator", "QHTOperator",
-           "SlopeOperator", "DomainDistributor"]
diff --git a/nifty5/operators/exp_transform.py b/nifty5/operators/exp_transform.py
index 4f56a07d1815ae24f8066661a2ac1147631e9fb3..cdf988c4ca54040dc1ea747bd3a8b911594d0718 100644
--- a/nifty5/operators/exp_transform.py
+++ b/nifty5/operators/exp_transform.py
@@ -1,7 +1,8 @@
 import numpy as np
 
 from ..domain_tuple import DomainTuple
-from ..domains import PowerSpace, RGSpace
+from ..domains.power_space import PowerSpace
+from ..domains.rg_space import RGSpace
 from ..field import Field
 from .linear_operator import LinearOperator
 from .. import dobj
diff --git a/nifty5/operators/linear_operator.py b/nifty5/operators/linear_operator.py
index 6cf5924b90abab3365c5d72833276f94b2d6b72e..988cfbddbbd7172202ca9e7787726ff88ea2a33b 100644
--- a/nifty5/operators/linear_operator.py
+++ b/nifty5/operators/linear_operator.py
@@ -198,7 +198,8 @@ class LinearOperator(NiftyMetaBase()):
 
     def __call__(self, x):
         """Same as :meth:`times`"""
-        from ..models import LinearModel, Model
+        from ..models.model import Model
+        from ..models.linear_model import LinearModel
         if isinstance(x, Model):
             return LinearModel(x, self)
         return self.apply(x, self.TIMES)
diff --git a/nifty5/operators/multi_adaptor.py b/nifty5/operators/multi_adaptor.py
index 7bc8a644228640f5bc9e8a2255863b3a2b0139dd..b54cc4cd1c23068dafe1bd9ebba02fbed5c22eee 100644
--- a/nifty5/operators/multi_adaptor.py
+++ b/nifty5/operators/multi_adaptor.py
@@ -16,12 +16,13 @@
 # NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik
 # and financially supported by the Studienstiftung des deutschen Volkes.
 
-from ..multi import MultiDomain, MultiField
+from ..multi.multi_domain import MultiDomain
+from ..multi.multi_field import MultiField
 from .linear_operator import LinearOperator
 
 
 class MultiAdaptor(LinearOperator):
-    """Transforms a Field into a MultiField and vise versa when 
+    """Transforms a Field into a MultiField and vise versa when
     using adjoint_times.
 
     Parameters
diff --git a/nifty5/operators/selection_operator.py b/nifty5/operators/selection_operator.py
index e2aa8f90885b2fab47f6ac5d710bf7da2be34631..4d88e442903aedd4eaa23063574050c9e61f1512 100644
--- a/nifty5/operators/selection_operator.py
+++ b/nifty5/operators/selection_operator.py
@@ -16,7 +16,7 @@
 # NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik
 # and financially supported by the Studienstiftung des deutschen Volkes.
 
-from ..operators import LinearOperator
+from .linear_operator import LinearOperator
 
 
 class SelectionOperator(LinearOperator):
@@ -31,7 +31,7 @@ class SelectionOperator(LinearOperator):
         String identifier of the wanted subdomain
     """
     def __init__(self, domain, key):
-        from ..multi import MultiDomain
+        from ..multi.multi_domain import MultiDomain
         if not isinstance(domain, MultiDomain):
             raise TypeError("Domain must be a MultiDomain")
         self._target = domain[key]
@@ -56,5 +56,5 @@ class SelectionOperator(LinearOperator):
         if mode == self.TIMES:
             return x[self._key].copy()
         else:
-            from ..multi import MultiField
+            from ..multi.multi_field import MultiField
             return MultiField({self._key: x.copy()})
diff --git a/test/__init__.py b/test/__init__.py
index e4193cf0582fc150f03646dbec76d0930622bdc6..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 100644
--- a/test/__init__.py
+++ b/test/__init__.py
@@ -1 +0,0 @@
-from . import common