diff --git a/Core/Export/PyFmt.cpp b/Core/Export/PyFmt.cpp
index b008fd342cfe3ec18aa8b54f9a6952cb8294957c..2db857382385abd5bc9332db8ce3a31b702b077c 100644
--- a/Core/Export/PyFmt.cpp
+++ b/Core/Export/PyFmt.cpp
@@ -21,7 +21,7 @@
 
 namespace pyfmt {
 
-std::string preambled(const std::string& code)
+std::string printImportedSymbols(const std::string& code)
 {
     std::vector<std::string> to_declare;
     for (const std::string& key : {"angstrom", "deg", "nm", "nm2", "micrometer"})
@@ -30,10 +30,7 @@ std::string preambled(const std::string& code)
     for (const std::string& key : {"kvector_t"})
         if (code.find(key) != std::string::npos)
             to_declare.push_back(key);
-    return "import numpy, sys\n"
-           "import bornagain as ba\n"
-           "from bornagain import "
-           + StringUtils::join(to_declare, ", ") + "\n\n\n" + code;
+    return "from bornagain import " + StringUtils::join(to_declare, ", ") + "\n";
 }
 
 std::string printBool(double value)
diff --git a/Core/Export/PyFmt.h b/Core/Export/PyFmt.h
index bd174b9d4b91a2fe67d6dac41c4f5fc8ccd9442e..e3dfc741b1da73b1fe03d562d82be1f5d49985f2 100644
--- a/Core/Export/PyFmt.h
+++ b/Core/Export/PyFmt.h
@@ -27,7 +27,7 @@
 
 namespace pyfmt {
 
-std::string preambled(const std::string& code);
+std::string printImportedSymbols(const std::string& code);
 
 std::string printInt(int value);
 std::string printBool(double value);
diff --git a/Core/Export/SimulationToPython.cpp b/Core/Export/SimulationToPython.cpp
index ff0b0364f44e14e2b6e58af8cbbeabf44c3f422d..19483816bb9388da31ce1a7db91fea51a0873d62 100644
--- a/Core/Export/SimulationToPython.cpp
+++ b/Core/Export/SimulationToPython.cpp
@@ -464,7 +464,10 @@ std::string simulationCode(const ISimulation& simulation)
         throw std::runtime_error("Cannot export: Simulation has no sample");
     std::string code =
         SampleToPython().sampleCode(*simulation.sample()) + defineSimulate(&simulation);
-    return pyfmt::preambled(code);
+    return "import bornagain as ba\n"
+        + pyfmt::printImportedSymbols(code)
+        + "import ba_plot if  __name__ == '__main__'"
+        "\n\n" + code;
 }
 
 } // namespace
@@ -477,7 +480,8 @@ std::string SimulationToPython::simulationPlotCode(const ISimulation& simulation
 {
     return simulationCode(simulation)
            + "if __name__ == '__main__':\n"
-             "    ba.run_and_plot(get_simulation(get_sample()))\n";
+             "    import ba_plot\n"
+             "    ba_plot.run_and_plot(get_simulation(get_sample()))\n";
 }
 
 std::string SimulationToPython::simulationSaveCode(const ISimulation& simulation,
@@ -485,6 +489,7 @@ std::string SimulationToPython::simulationSaveCode(const ISimulation& simulation
 {
     return simulationCode(simulation)
            + "if __name__ == '__main__':\n"
-             "    ba.run_and_save(get_simulation(get_sample()), \""
+             "    import ba_plot\n"
+             "    baplot.run_and_save(get_simulation(get_sample()), \""
            + fname + "\")\n";
 }
diff --git a/Examples/scatter2d/ApproximationDA.py b/Examples/scatter2d/ApproximationDA.py
index f6d60aa1180586eeeab8b939f12b8ceb6d32bb24..e3d7392c9ebb38544f54ebe673b1f32056110bc7 100644
--- a/Examples/scatter2d/ApproximationDA.py
+++ b/Examples/scatter2d/ApproximationDA.py
@@ -59,4 +59,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/ApproximationLMA.py b/Examples/scatter2d/ApproximationLMA.py
index 1d69886e5bec82e76e1fd209f4b692dbd72abe92..db46e4af07013235692c5052412ec9f0baea6986 100644
--- a/Examples/scatter2d/ApproximationLMA.py
+++ b/Examples/scatter2d/ApproximationLMA.py
@@ -67,4 +67,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/ApproximationSSCA.py b/Examples/scatter2d/ApproximationSSCA.py
index 287dc47c68c67264dc589b8afcd11b9d6c380857..663e194a90f6c6dc48a2e5a6c5f4618ee8e7358c 100644
--- a/Examples/scatter2d/ApproximationSSCA.py
+++ b/Examples/scatter2d/ApproximationSSCA.py
@@ -60,4 +60,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/BeamDivergence.py b/Examples/scatter2d/BeamDivergence.py
index e8d4a04fc849a9c83b6db60e7a0d6dae763f73f2..07e6f1b285891c5684805c2f568ce6564e01fbf3 100644
--- a/Examples/scatter2d/BeamDivergence.py
+++ b/Examples/scatter2d/BeamDivergence.py
@@ -57,4 +57,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/BiMaterialCylinders.py b/Examples/scatter2d/BiMaterialCylinders.py
index ac991ce7e7984f89c7c0947be71ff589998533cd..0d8e8d351af7ab0daed48a3d017513ab11e63d11 100644
--- a/Examples/scatter2d/BiMaterialCylinders.py
+++ b/Examples/scatter2d/BiMaterialCylinders.py
@@ -65,4 +65,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/BoxesWithSpecularPeak.py b/Examples/scatter2d/BoxesWithSpecularPeak.py
index ca87675492c734fa9eba542365b52daf66401d9f..8c0184df47922c040e47eae2230bf1535b812f0a 100644
--- a/Examples/scatter2d/BoxesWithSpecularPeak.py
+++ b/Examples/scatter2d/BoxesWithSpecularPeak.py
@@ -60,4 +60,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/BuriedParticles.py b/Examples/scatter2d/BuriedParticles.py
index 6fb9cd0efa19b48c9610a586817e22e4e5bdec8d..cbe1b459e5b3d370c24842ec6906d34b81b97f33 100644
--- a/Examples/scatter2d/BuriedParticles.py
+++ b/Examples/scatter2d/BuriedParticles.py
@@ -56,4 +56,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/ConstantBackground.py b/Examples/scatter2d/ConstantBackground.py
index 5244d0ceb44ecb4422243b28553e06a4547a048e..7b1f7ed27581d0f54f6eb8f8ca7b8c246203899a 100644
--- a/Examples/scatter2d/ConstantBackground.py
+++ b/Examples/scatter2d/ConstantBackground.py
@@ -51,4 +51,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/CoreShellNanoparticles.py b/Examples/scatter2d/CoreShellNanoparticles.py
index 19162de038f9f1adab3e4f0dc3cf408886ad6a95..8c2148dc7e5958b51c51a41e49e395ab3cc4f830 100644
--- a/Examples/scatter2d/CoreShellNanoparticles.py
+++ b/Examples/scatter2d/CoreShellNanoparticles.py
@@ -56,4 +56,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/CorrelatedRoughness.py b/Examples/scatter2d/CorrelatedRoughness.py
index 916215b5f3d6495c0c40acd3e3af733d49bd805c..1ce98d904cfc5cd4c91ee550253914b95e176609 100644
--- a/Examples/scatter2d/CorrelatedRoughness.py
+++ b/Examples/scatter2d/CorrelatedRoughness.py
@@ -51,4 +51,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/CosineRipplesAtRectLattice.py b/Examples/scatter2d/CosineRipplesAtRectLattice.py
index 694f055b18d42239bd4c27d9755f64ea177dffd2..b34cdb3fdefbc8f747883f76082880f0194a4f92 100644
--- a/Examples/scatter2d/CosineRipplesAtRectLattice.py
+++ b/Examples/scatter2d/CosineRipplesAtRectLattice.py
@@ -59,4 +59,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/CylindersAndPrisms.py b/Examples/scatter2d/CylindersAndPrisms.py
index 216d9afa21fc7e527523825747332f02a575f617..cf8db3228dd171cb9476614094014899da7a7b4c 100644
--- a/Examples/scatter2d/CylindersAndPrisms.py
+++ b/Examples/scatter2d/CylindersAndPrisms.py
@@ -56,4 +56,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/CylindersInAverageLayer.py b/Examples/scatter2d/CylindersInAverageLayer.py
index e13db805276b53b4f7ef36268695407b952b2719..d2b8925c3a2a24bd42fe961a6c84ecae12c043f5 100644
--- a/Examples/scatter2d/CylindersInAverageLayer.py
+++ b/Examples/scatter2d/CylindersInAverageLayer.py
@@ -51,4 +51,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/CylindersInBA.py b/Examples/scatter2d/CylindersInBA.py
index ec8d50a4d19c6726c6712150220df56ddfc9f6f4..17fc834da2e7e81049cfe6c5a479514603677277 100644
--- a/Examples/scatter2d/CylindersInBA.py
+++ b/Examples/scatter2d/CylindersInBA.py
@@ -47,4 +47,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/CylindersInDWBA.py b/Examples/scatter2d/CylindersInDWBA.py
index 076dd5960ea6f8c30f582d65bec6fcfe29418433..b8420ad4a09ef06d2322c4fba4b82f42aec3a537 100644
--- a/Examples/scatter2d/CylindersInDWBA.py
+++ b/Examples/scatter2d/CylindersInDWBA.py
@@ -49,4 +49,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/CylindersWithSizeDistribution.py b/Examples/scatter2d/CylindersWithSizeDistribution.py
index 65439c35962bcc4d4dc761d58ec248073f05f2ec..4cfce432bfefc1981cbd3094dd29680f82e5356f 100644
--- a/Examples/scatter2d/CylindersWithSizeDistribution.py
+++ b/Examples/scatter2d/CylindersWithSizeDistribution.py
@@ -53,4 +53,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/DetectorResolutionFunction.py b/Examples/scatter2d/DetectorResolutionFunction.py
index d6bde346d70437cfb77b3464cce1d7995823d38c..c9d2a61f956cf16c6d4a3aa6a24c37415352c676 100644
--- a/Examples/scatter2d/DetectorResolutionFunction.py
+++ b/Examples/scatter2d/DetectorResolutionFunction.py
@@ -51,4 +51,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/HalfSpheresInAverageTopLayer.py b/Examples/scatter2d/HalfSpheresInAverageTopLayer.py
index d8933487c15ee857e3fe4d076279a20aa13c5d1e..96fdcad873b8545886d087281dd8a8e7d57c869e 100644
--- a/Examples/scatter2d/HalfSpheresInAverageTopLayer.py
+++ b/Examples/scatter2d/HalfSpheresInAverageTopLayer.py
@@ -61,4 +61,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/HexagonalLatticesWithBasis.py b/Examples/scatter2d/HexagonalLatticesWithBasis.py
index dd99d3f2c12d94698c014dbb3744acd052252230..62c91b2c42adaf263bf26772d8deb4479a299321 100644
--- a/Examples/scatter2d/HexagonalLatticesWithBasis.py
+++ b/Examples/scatter2d/HexagonalLatticesWithBasis.py
@@ -68,4 +68,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/Interference1DRadialParaCrystal.py b/Examples/scatter2d/Interference1DRadialParaCrystal.py
index b17da17637cf6aeaee65c27d0425fb0ad0e890be..616ccfa7aac0b581cd0fda64d3bfbc69e0cb2491 100644
--- a/Examples/scatter2d/Interference1DRadialParaCrystal.py
+++ b/Examples/scatter2d/Interference1DRadialParaCrystal.py
@@ -55,4 +55,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/Interference2DCenteredSquareLattice.py b/Examples/scatter2d/Interference2DCenteredSquareLattice.py
index d56deadc3ca34f85a50c50437d04b6f0d6af1517..cbd0cf87f17784032859a44f2bb14435798857f7 100644
--- a/Examples/scatter2d/Interference2DCenteredSquareLattice.py
+++ b/Examples/scatter2d/Interference2DCenteredSquareLattice.py
@@ -68,4 +68,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/Interference2DLatticeSumOfRotated.py b/Examples/scatter2d/Interference2DLatticeSumOfRotated.py
index ac2a432a9da7e192364851d3d6305380eaaab203..32c1d01c575228c608889517d8216f1ab5e3f818 100644
--- a/Examples/scatter2d/Interference2DLatticeSumOfRotated.py
+++ b/Examples/scatter2d/Interference2DLatticeSumOfRotated.py
@@ -46,4 +46,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/Interference2DParaCrystal.py b/Examples/scatter2d/Interference2DParaCrystal.py
index 013a65d4b458132bea58253b2009e6eea19ecc72..7570ea6a668ebbfa8745028e77aa175f22f8b720 100644
--- a/Examples/scatter2d/Interference2DParaCrystal.py
+++ b/Examples/scatter2d/Interference2DParaCrystal.py
@@ -61,4 +61,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/Interference2DRotatedSquareLattice.py b/Examples/scatter2d/Interference2DRotatedSquareLattice.py
index f428b4409637b43848f0d11f24738b227de4494d..6979000034164d4ed2b4bffb7e8387630dc6a474 100644
--- a/Examples/scatter2d/Interference2DRotatedSquareLattice.py
+++ b/Examples/scatter2d/Interference2DRotatedSquareLattice.py
@@ -59,4 +59,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/Interference2DSquareFiniteLattice.py b/Examples/scatter2d/Interference2DSquareFiniteLattice.py
index 5fc42c81f4853502b19feb56df8f809956b47f50..ac3bdaf8e1fc635c81c6d8cf3bfaecdb7a2c835f 100644
--- a/Examples/scatter2d/Interference2DSquareFiniteLattice.py
+++ b/Examples/scatter2d/Interference2DSquareFiniteLattice.py
@@ -57,4 +57,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/Interference2DSquareLattice.py b/Examples/scatter2d/Interference2DSquareLattice.py
index cb64a4e5415e4e78e4b21bd2883147c6f3a03859..1476f817786d05cc736c26f8a998c9430a39723c 100644
--- a/Examples/scatter2d/Interference2DSquareLattice.py
+++ b/Examples/scatter2d/Interference2DSquareLattice.py
@@ -58,4 +58,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/MagneticSpheres.py b/Examples/scatter2d/MagneticSpheres.py
index ef10a3bf24d124e387b2659f3c98bea79e3937e5..ef60d2b05d92ef3db5952d4fb8d30e9cad488bcd 100644
--- a/Examples/scatter2d/MagneticSpheres.py
+++ b/Examples/scatter2d/MagneticSpheres.py
@@ -57,4 +57,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/MesoCrystal.py b/Examples/scatter2d/MesoCrystal.py
index c238716ca40d635712e945e67c0a2f217c124cc6..5d7f473319659179c10fb71877a604e928a3e605 100644
--- a/Examples/scatter2d/MesoCrystal.py
+++ b/Examples/scatter2d/MesoCrystal.py
@@ -61,4 +61,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/ParticlesCrossingInterface.py b/Examples/scatter2d/ParticlesCrossingInterface.py
index a3be3b18f528851b620ad51d834f92e1a768bf2e..dedbeb3a6505ebfd3626de9d9aca808b8b495ec1 100644
--- a/Examples/scatter2d/ParticlesCrossingInterface.py
+++ b/Examples/scatter2d/ParticlesCrossingInterface.py
@@ -75,4 +75,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/RotatedPyramids.py b/Examples/scatter2d/RotatedPyramids.py
index a4e73820ca7c2b52e4464c1112b5aaf77f4ca06c..100965a8a4e42ab110ef0e1232844eabc612214e 100644
--- a/Examples/scatter2d/RotatedPyramids.py
+++ b/Examples/scatter2d/RotatedPyramids.py
@@ -51,4 +51,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/SpheresAtHexLattice.py b/Examples/scatter2d/SpheresAtHexLattice.py
index 961213c595bfb5553f4b08f3e02a0a7378265d44..d2bffe000f467af00a9b4a109e6a97057f94eb50 100644
--- a/Examples/scatter2d/SpheresAtHexLattice.py
+++ b/Examples/scatter2d/SpheresAtHexLattice.py
@@ -59,4 +59,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/TriangularRipple.py b/Examples/scatter2d/TriangularRipple.py
index 628e36846f17f7cb7155402f7e996c4cf1556311..6074955fca6b2fb81aef44f43b975e4a6b1e285c 100644
--- a/Examples/scatter2d/TriangularRipple.py
+++ b/Examples/scatter2d/TriangularRipple.py
@@ -60,4 +60,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Examples/scatter2d/TwoTypesOfCylindersWithSizeDistribution.py b/Examples/scatter2d/TwoTypesOfCylindersWithSizeDistribution.py
index 7f4450d0fd4f96d33513ae5479dcc635b1bd1db1..de8116d82376ff6cee3ef0aede2a321ea95f258a 100644
--- a/Examples/scatter2d/TwoTypesOfCylindersWithSizeDistribution.py
+++ b/Examples/scatter2d/TwoTypesOfCylindersWithSizeDistribution.py
@@ -62,4 +62,5 @@ def get_simulation(sample):
 
 
 if __name__ == '__main__':
-    ba.run_and_plot(get_simulation(get_sample()))
+    import ba_plot
+    ba_plot.run_and_plot(get_simulation(get_sample()))
diff --git a/Tests/Functional/PyEmbedded/Tests.cpp b/Tests/Functional/PyEmbedded/Tests.cpp
index d54021a9b14ac3905c02353f9d76606068f0e346..27c9e31d4e7c437bf1117d4c9ae85e24541db836 100644
--- a/Tests/Functional/PyEmbedded/Tests.cpp
+++ b/Tests/Functional/PyEmbedded/Tests.cpp
@@ -372,7 +372,9 @@ TEST_F(PyEmbedded, ExportToPythonAndBack)
     std::unique_ptr<MultiLayer> sample(factory.createSampleByName("CylindersAndPrismsBuilder"));
 
     const std::string code = ExportToPython::sampleCode(*sample);
-    const std::string snippet = pyfmt::preambled(code);
+    const std::string snippet = "import bornagain as ba\n"
+        + pyfmt::printImportedSymbols(code)
+        + "\n\n" + code;
 
     const auto multilayer =
         PyImport::createFromPython(snippet, "get_sample", BABuild::buildLibDir());
diff --git a/auto/Wrap/doxygenCore.i b/auto/Wrap/doxygenCore.i
index 479b47904f7e200e58aa704a47d114b82870fcbf..2a2f3c967dafe774cffb3416ba844630a2b7f5e2 100644
--- a/auto/Wrap/doxygenCore.i
+++ b/auto/Wrap/doxygenCore.i
@@ -2406,7 +2406,7 @@ Returns default metric name.
 
 
 // File: namespacepyfmt.xml
-%feature("docstring")  pyfmt::preambled "std::string pyfmt::preambled(const std::string &code)
+%feature("docstring")  pyfmt::printImportedSymbols "std::string pyfmt::printImportedSymbols(const std::string &code)
 ";
 
 %feature("docstring")  pyfmt::printBool "std::string pyfmt::printBool(double value)