diff --git a/libsrc/core/python_ngcore.cpp b/libsrc/core/python_ngcore.cpp
index a0810ca6..f60e40d3 100644
--- a/libsrc/core/python_ngcore.cpp
+++ b/libsrc/core/python_ngcore.cpp
@@ -40,7 +40,7 @@ namespace ngcore
 
     if (py::isinstance<py::list>(value))
       {             
-        py::list vdl(value);
+        auto vdl = py::cast<py::list>(value);
         if (py::len(vdl) > 0)
           {
             if(py::isinstance<double>(vdl[0]))
@@ -59,7 +59,7 @@ namespace ngcore
 
     if (py::isinstance<py::tuple>(value))
       {
-        py::tuple vdt(value);
+        auto vdt = py::cast<py::tuple>(value);
         if (py::isinstance<py::float_>(value))
           flags.SetFlag(s, makeCArray<double>(vdt));
         if (py::isinstance<py::int_>(value))
@@ -69,29 +69,32 @@ namespace ngcore
       }
   }
 
-  Flags CreateFlagsFromKwArgs(py::object pyclass, const py::kwargs& kwargs, py::list info)
+  Flags CreateFlagsFromKwArgs(const py::kwargs& kwargs, py::object pyclass, py::list info)
   {
     static std::shared_ptr<Logger> logger = GetLogger("Flags");
-    auto flags_doc = pyclass.attr("__flags_doc__")();
     py::dict flags_dict;
 
     if (kwargs.contains("flags"))
       {
-        logger->warn("WARNING: using flags as kwarg is deprecated in {}, use the flag arguments as kwargs instead!",
-                     std::string(py::str(pyclass)));
+        logger->warn("WARNING: using flags as kwarg is deprecated{}, use the flag arguments as kwargs instead!",
+                     pyclass.is_none() ? "" : std::string(" in ") + std::string(py::str(pyclass)));
         auto addflags = py::cast<py::dict>(kwargs["flags"]);
         for (auto item : addflags)
           flags_dict[item.first.cast<string>().c_str()] = item.second;
       }
-    for (auto item : kwargs)
-      if (!flags_doc.contains(item.first.cast<string>().c_str()) &&
-          !(item.first.cast<string>() == "flags"))
-        logger->warn("WARNING: kwarg '{}' is an undocumented flags option for class {}, maybe there is a typo?",
-                     item.first.cast<string>(), std::string(py::str(pyclass)));
-
     py::dict special;
-    if(py::hasattr(pyclass,"__special_treated_flags__"))
-      special = pyclass.attr("__special_treated_flags__")();
+    if(!pyclass.is_none())
+      {
+        auto flags_doc = pyclass.attr("__flags_doc__")();
+        for (auto item : kwargs)
+          if (!flags_doc.contains(item.first.cast<string>().c_str()) &&
+              !(item.first.cast<string>() == "flags"))
+            logger->warn("WARNING: kwarg '{}' is an undocumented flags option for class {}, maybe there is a typo?",
+                         item.first.cast<string>(), std::string(py::str(pyclass)));
+      
+        if(py::hasattr(pyclass,"__special_treated_flags__"))
+          special = pyclass.attr("__special_treated_flags__")();
+      }
     for (auto item : kwargs)
       {
         auto name = item.first.cast<string>();
@@ -116,4 +119,47 @@ namespace ngcore
     return flags;
   }
 
+  py::dict CreateDictFromFlags(const Flags& flags)
+  {
+    py::dict d;
+    std::string key;
+    for(auto i : Range(flags.GetNFlagsFlags()))
+      {
+        auto& f = flags.GetFlagsFlag(i, key);
+        d[key.c_str()] = CreateDictFromFlags(f);
+      }
+    for(auto i : Range(flags.GetNStringListFlags()))
+      {
+        auto strlistflag = flags.GetStringListFlag(i, key);
+        py::list lst;
+        for(auto& val : *strlistflag)
+          lst.append(val);
+        d[key.c_str()] = lst;
+      }
+    for(auto i : Range(flags.GetNNumListFlags()))
+      {
+        auto numlistflag = flags.GetNumListFlag(i, key);
+        py::list lst;
+        for(auto& val : *numlistflag)
+          lst.append(val);
+        d[key.c_str()] = lst;
+      }
+    for(auto i : Range(flags.GetNStringFlags()))
+      {
+        auto val = flags.GetStringFlag(i, key);
+        d[key.c_str()] = val;
+      }
+    for(auto i : Range(flags.GetNNumFlags()))
+      {
+        auto val = flags.GetNumFlag(i, key);
+        d[key.c_str()] = val;
+      }
+    for(auto i : Range(flags.GetNDefineFlags()))
+      {
+        auto val = flags.GetDefineFlag(i, key);
+        d[key.c_str()] = val;
+      }
+    return d;
+  }
+
 } // namespace ngcore
diff --git a/libsrc/core/python_ngcore.hpp b/libsrc/core/python_ngcore.hpp
index 562467cc..7a2fe5b2 100644
--- a/libsrc/core/python_ngcore.hpp
+++ b/libsrc/core/python_ngcore.hpp
@@ -16,10 +16,9 @@ namespace ngcore
   Array<T> makeCArray(const py::object& obj)
   {
     Array<T> arr;
-    arr.SetAllocSize(py::len(obj));
     if(py::isinstance<py::list>(obj))
-      for(auto& val : py::cast<py::list>(obj))
-        arr.Append(py::cast<T>(val));
+        for(auto& val : py::cast<py::list>(obj))
+          arr.Append(py::cast<T>(val));
     else if(py::isinstance<py::tuple>(obj))
       for(auto& val : py::cast<py::tuple>(obj))
         arr.Append(py::cast<T>(val));
@@ -30,7 +29,10 @@ namespace ngcore
   
   void NGCORE_API SetFlag(Flags &flags, std::string s, py::object value);
   // Parse python kwargs to flags
-  Flags NGCORE_API CreateFlagsFromKwArgs(py::object pyclass, const py::kwargs& kwargs, py::list info = py::list());
+  Flags NGCORE_API CreateFlagsFromKwArgs(const py::kwargs& kwargs, py::object pyclass = py::none(),
+                                         py::list info = py::list());
+  // Create python dict from kwargs
+  py::dict NGCORE_API CreateDictFromFlags(const Flags& flags);
 
   // ***************  Archiving functionality  **************
   
diff --git a/libsrc/meshing/meshfunc.cpp b/libsrc/meshing/meshfunc.cpp
index b7d91a14..25a7e8b7 100644
--- a/libsrc/meshing/meshfunc.cpp
+++ b/libsrc/meshing/meshfunc.cpp
@@ -12,10 +12,11 @@ namespace netgen
 
 
   // extern double teterrpow; 
-  MESHING3_RESULT MeshVolume (MeshingParameters & mp, Mesh& mesh3d)
+  MESHING3_RESULT MeshVolume (const MeshingParameters & c_mp, Mesh& mesh3d)
   {
     static Timer t("MeshVolume"); RegionTimer reg(t);
-    
+
+     MeshingParameters mp = c_mp; // copy mp to change them here
      int oldne;
      int meshed;
 
@@ -637,7 +638,7 @@ namespace netgen
 
 
 
-  MESHING3_RESULT OptimizeVolume (MeshingParameters & mp, 
+  MESHING3_RESULT OptimizeVolume (const MeshingParameters & mp, 
 				  Mesh & mesh3d)
     //				  const CSGeometry * geometry)
   {
diff --git a/libsrc/meshing/meshfunc.hpp b/libsrc/meshing/meshfunc.hpp
index 7eaaa238..fdbdef4e 100644
--- a/libsrc/meshing/meshfunc.hpp
+++ b/libsrc/meshing/meshfunc.hpp
@@ -16,13 +16,13 @@ class Mesh;
 // class CSGeometry;
 
 /// Build tet-mesh
-DLL_HEADER MESHING3_RESULT MeshVolume (MeshingParameters & mp, Mesh& mesh3d);
+DLL_HEADER MESHING3_RESULT MeshVolume (const MeshingParameters & mp, Mesh& mesh3d);
 
 /// Build mixed-element mesh
 // MESHING3_RESULT MeshMixedVolume (MeshingParameters & mp, Mesh& mesh3d);
 
 /// Optimize tet-mesh
-DLL_HEADER MESHING3_RESULT OptimizeVolume (MeshingParameters & mp, Mesh& mesh3d);
+DLL_HEADER MESHING3_RESULT OptimizeVolume (const MeshingParameters & mp, Mesh& mesh3d);
 //			       const CSGeometry * geometry = NULL);
 
 DLL_HEADER void RemoveIllegalElements (Mesh & mesh3d);
diff --git a/libsrc/meshing/meshtype.hpp b/libsrc/meshing/meshtype.hpp
index 700727e7..aa6eb1d3 100644
--- a/libsrc/meshing/meshtype.hpp
+++ b/libsrc/meshing/meshtype.hpp
@@ -1280,6 +1280,8 @@ namespace netgen
     bool inverttrigs = false;
     ///
     bool autozrefine = false;
+
+    Flags geometrySpecificParameters;
     ///
     MeshingParameters ();
     ///
@@ -1307,7 +1309,7 @@ namespace netgen
     NgArray<MeshSizePoint> meshsize_points;
     
     void (*render_function)(bool) = NULL;
-    void Render(bool blocking = false)
+    void Render(bool blocking = false) const
     {
       if (render_function) 
         (*render_function)(blocking);
diff --git a/libsrc/meshing/python_mesh.cpp b/libsrc/meshing/python_mesh.cpp
index a58faa70..16f1b269 100644
--- a/libsrc/meshing/python_mesh.cpp
+++ b/libsrc/meshing/python_mesh.cpp
@@ -1019,17 +1019,14 @@ DLL_HEADER void ExportNetgenMeshing(py::module &m)
   typedef MeshingParameters MP;
   auto mp = py::class_<MP> (m, "MeshingParameters")
     .def(py::init<>())
-    .def(py::init([](py::kwargs kwargs)
+            .def(py::init([](MeshingParameters* other, py::kwargs kwargs)
                   {
                     MeshingParameters mp;
-                    CreateMPfromKwargs(mp, kwargs);
+                    if(other) mp = *other;
+                    CreateMPfromKwargs(mp, kwargs, false);
                     return mp;
-                  }), meshingparameter_description.c_str())
+                  }), py::arg("mp")=nullptr, meshingparameter_description.c_str())
     .def("__str__", &ToString<MP>)
-    .def_property("maxh", [](const MP & mp ) { return mp.maxh; },
-                  [](MP & mp, double maxh) { return mp.maxh = maxh; })
-    .def_property("grading", [](const MP & mp ) { return mp.grading; },
-                  [](MP & mp, double grading) { return mp.grading = grading; })
     .def("RestrictH", FunctionPointer
          ([](MP & mp, double x, double y, double z, double h)
           {
diff --git a/libsrc/meshing/python_mesh.hpp b/libsrc/meshing/python_mesh.hpp
index 8a35cd8b..62a85778 100644
--- a/libsrc/meshing/python_mesh.hpp
+++ b/libsrc/meshing/python_mesh.hpp
@@ -1,5 +1,9 @@
+#ifndef NETGEN_MESHING_PYTHON_MESH_HPP
+#define NETGEN_MESHING_PYTHON_MESH_HPP
 
 #include <pybind11/pybind11.h>
+
+#include <core/python_ngcore.hpp>
 #include "meshing.hpp"
 
 namespace netgen
@@ -80,89 +84,98 @@ elsizeweight: float = 0.2
 
 )delimiter";
 
-  inline void CreateMPfromKwargs(MeshingParameters& mp, py::kwargs kwargs)
+inline void CreateMPfromKwargs(MeshingParameters& mp, py::kwargs kwargs, bool throw_if_not_all_parsed=true)
   {
     if(kwargs.contains("optimize3d"))
-      mp.optimize3d = py::cast<string>(kwargs["optimize3d"]);
+      mp.optimize3d = py::cast<string>(kwargs.attr("pop")("optimize3d"));
     if(kwargs.contains("optsteps3d"))
-      mp.optsteps3d = py::cast<int>(kwargs["optsteps3d"]);
+      mp.optsteps3d = py::cast<int>(kwargs.attr("pop")("optsteps3d"));
     if(kwargs.contains("optimize2d"))
-      mp.optimize2d = py::cast<string>(kwargs["optimize2d"]);
+      mp.optimize2d = py::cast<string>(kwargs.attr("pop")("optimize2d"));
     if(kwargs.contains("optsteps2d"))
-      mp.optsteps2d = py::cast<int>(kwargs["optsteps2d"]);
+      mp.optsteps2d = py::cast<int>(kwargs.attr("pop")("optsteps2d"));
     if(kwargs.contains("opterrpow"))
-      mp.opterrpow = py::cast<double>(kwargs["opterrpow"]);
+      mp.opterrpow = py::cast<double>(kwargs.attr("pop")("opterrpow"));
     if(kwargs.contains("blockfill"))
-      mp.blockfill = py::cast<bool>(kwargs["blockfill"]);
+      mp.blockfill = py::cast<bool>(kwargs.attr("pop")("blockfill"));
     if(kwargs.contains("filldist"))
-      mp.filldist = py::cast<double>(kwargs["filldist"]);
+      mp.filldist = py::cast<double>(kwargs.attr("pop")("filldist"));
     if(kwargs.contains("safety"))
-      mp.safety = py::cast<double>(kwargs["safety"]);
+      mp.safety = py::cast<double>(kwargs.attr("pop")("safety"));
     if(kwargs.contains("relinnersafety"))
-      mp.relinnersafety = py::cast<double>(kwargs["relinnersafety"]);
+      mp.relinnersafety = py::cast<double>(kwargs.attr("pop")("relinnersafety"));
     if(kwargs.contains("uselocalh"))
-      mp.uselocalh = py::cast<bool>(kwargs["uselocalh"]);
+      mp.uselocalh = py::cast<bool>(kwargs.attr("pop")("uselocalh"));
     if(kwargs.contains("grading"))
-      mp.grading = py::cast<double>(kwargs["grading"]);
+      mp.grading = py::cast<double>(kwargs.attr("pop")("grading"));
     if(kwargs.contains("delaunay"))
-      mp.delaunay = py::cast<bool>(kwargs["delaunay"]);
+      mp.delaunay = py::cast<bool>(kwargs.attr("pop")("delaunay"));
     if(kwargs.contains("maxh"))
-      mp.maxh = py::cast<double>(kwargs["maxh"]);
+      mp.maxh = py::cast<double>(kwargs.attr("pop")("maxh"));
     if(kwargs.contains("minh"))
-      mp.minh = py::cast<double>(kwargs["minh"]);
+      mp.minh = py::cast<double>(kwargs.attr("pop")("minh"));
     if(kwargs.contains("meshsizefilename"))
-      mp.meshsizefilename = py::cast<string>(kwargs["meshsizefilename"]);
+      mp.meshsizefilename = py::cast<string>(kwargs.attr("pop")("meshsizefilename"));
     if(kwargs.contains("startinsurface"))
-      mp.startinsurface = py::cast<bool>(kwargs["startinsurface"]);
+      mp.startinsurface = py::cast<bool>(kwargs.attr("pop")("startinsurface"));
     if(kwargs.contains("checkoverlap"))
-      mp.checkoverlap = py::cast<bool>(kwargs["checkoverlap"]);
+      mp.checkoverlap = py::cast<bool>(kwargs.attr("pop")("checkoverlap"));
     if(kwargs.contains("checkoverlappingboundary"))
-      mp.checkoverlappingboundary = py::cast<bool>(kwargs["checkoverlappingboundary"]);
+      mp.checkoverlappingboundary = py::cast<bool>(kwargs.attr("pop")("checkoverlappingboundary"));
     if(kwargs.contains("checkchartboundary"))
-      mp.checkchartboundary = py::cast<bool>(kwargs["checkchartboundary"]);
+      mp.checkchartboundary = py::cast<bool>(kwargs.attr("pop")("checkchartboundary"));
     if(kwargs.contains("curvaturesafety"))
-      mp.curvaturesafety = py::cast<double>(kwargs["curvaturesafety"]);
+      mp.curvaturesafety = py::cast<double>(kwargs.attr("pop")("curvaturesafety"));
     if(kwargs.contains("segmentsperedge"))
-      mp.segmentsperedge = py::cast<double>(kwargs["segmentsperedge"]);
+      mp.segmentsperedge = py::cast<double>(kwargs.attr("pop")("segmentsperedge"));
     if(kwargs.contains("parthread"))
-      mp.parthread = py::cast<bool>(kwargs["parthread"]);
+      mp.parthread = py::cast<bool>(kwargs.attr("pop")("parthread"));
     if(kwargs.contains("elsizeweight"))
-      mp.elsizeweight = py::cast<double>(kwargs["elsizeweight"]);
+      mp.elsizeweight = py::cast<double>(kwargs.attr("pop")("elsizeweight"));
     if(kwargs.contains("perfstepsstart"))
-      mp.perfstepsstart = py::cast<int>(kwargs["perfstepsstart"]);
+      mp.perfstepsstart = py::cast<int>(kwargs.attr("pop")("perfstepsstart"));
     if(kwargs.contains("perfstepsend"))
-      mp.perfstepsend = py::cast<int>(kwargs["perfstepsend"]);
+      mp.perfstepsend = py::cast<int>(kwargs.attr("pop")("perfstepsend"));
     if(kwargs.contains("giveuptol2d"))
-      mp.giveuptol2d = py::cast<int>(kwargs["giveuptol2d"]);
+      mp.giveuptol2d = py::cast<int>(kwargs.attr("pop")("giveuptol2d"));
     if(kwargs.contains("giveuptol"))
-      mp.giveuptol = py::cast<int>(kwargs["giveuptol"]);
+      mp.giveuptol = py::cast<int>(kwargs.attr("pop")("giveuptol"));
     if(kwargs.contains("maxoutersteps"))
-      mp.maxoutersteps = py::cast<int>(kwargs["maxoutersteps"]);
+      mp.maxoutersteps = py::cast<int>(kwargs.attr("pop")("maxoutersteps"));
     if(kwargs.contains("starshapeclass"))
-      mp.starshapeclass = py::cast<int>(kwargs["starshapeclass"]);
+      mp.starshapeclass = py::cast<int>(kwargs.attr("pop")("starshapeclass"));
     if(kwargs.contains("baseelnp"))
-      mp.baseelnp = py::cast<int>(kwargs["baseelnp"]);
+      mp.baseelnp = py::cast<int>(kwargs.attr("pop")("baseelnp"));
     if(kwargs.contains("sloppy"))
-      mp.sloppy = py::cast<int>(kwargs["sloppy"]);
+      mp.sloppy = py::cast<int>(kwargs.attr("pop")("sloppy"));
     if(kwargs.contains("badellimit"))
-      mp.badellimit = py::cast<double>(kwargs["badellimit"]);
+      mp.badellimit = py::cast<double>(kwargs.attr("pop")("badellimit"));
     if(kwargs.contains("check_impossible"))
-      mp.check_impossible = py::cast<bool>(kwargs["check_impossible"]);
+      mp.check_impossible = py::cast<bool>(kwargs.attr("pop")("check_impossible"));
     if(kwargs.contains("only3D_domain_nr"))
-      mp.only3D_domain_nr = py::cast<int>(kwargs["only3D_domain_nr"]);
+      mp.only3D_domain_nr = py::cast<int>(kwargs.attr("pop")("only3D_domain_nr"));
     if(kwargs.contains("secondorder"))
-      mp.secondorder = py::cast<bool>(kwargs["secondorder"]);
+      mp.secondorder = py::cast<bool>(kwargs.attr("pop")("secondorder"));
     if(kwargs.contains("elementorder"))
-      mp.elementorder = py::cast<int>(kwargs["elementorder"]);
+      mp.elementorder = py::cast<int>(kwargs.attr("pop")("elementorder"));
     if(kwargs.contains("quad"))
-      mp.quad = py::cast<bool>(kwargs["quad"]);
+      mp.quad = py::cast<bool>(kwargs.attr("pop")("quad"));
     if(kwargs.contains("try_hexes"))
-      mp.try_hexes = py::cast<bool>(kwargs["try_hexes"]);
+      mp.try_hexes = py::cast<bool>(kwargs.attr("pop")("try_hexes"));
     if(kwargs.contains("inverttets"))
-      mp.inverttets = py::cast<bool>(kwargs["inverttets"]);
+      mp.inverttets = py::cast<bool>(kwargs.attr("pop")("inverttets"));
     if(kwargs.contains("inverttrigs"))
-      mp.inverttrigs = py::cast<bool>(kwargs["inverttrigs"]);
+      mp.inverttrigs = py::cast<bool>(kwargs.attr("pop")("inverttrigs"));
     if(kwargs.contains("autozrefine"))
-      mp.autozrefine = py::cast<bool>(kwargs["autozrefine"]);
+      mp.autozrefine = py::cast<bool>(kwargs.attr("pop")("autozrefine"));
+    if(kwargs.size())
+    {
+      if(throw_if_not_all_parsed)
+        throw Exception(string("Not all kwargs given to GenerateMesh could be parsed:") + string(py::str(kwargs)));
+      mp.geometrySpecificParameters = CreateFlagsFromKwArgs(kwargs);
+    }
   }
 } // namespace netgen
+
+#endif // NETGEN_MESHING_PYTHON_MESH_HPP
+ 
diff --git a/libsrc/stlgeom/meshstlsurface.cpp b/libsrc/stlgeom/meshstlsurface.cpp
index 4c1d940d..0abca632 100644
--- a/libsrc/stlgeom/meshstlsurface.cpp
+++ b/libsrc/stlgeom/meshstlsurface.cpp
@@ -13,15 +13,16 @@
 namespace netgen
 {
 
-static void STLFindEdges (STLGeometry & geom,
-			  class Mesh & mesh)
+static void STLFindEdges (STLGeometry & geom, Mesh & mesh,
+                          const MeshingParameters& mparam,
+                          const STLParameters& stlparam)
 {
   double h = mparam.maxh;
 
   // mark edge points:
   //int ngp = geom.GetNP();
 
-  geom.RestrictLocalH(mesh, h);
+  geom.RestrictLocalH(mesh, h, stlparam);
   
   PushStatusF("Mesh Lines");
 
@@ -229,18 +230,19 @@ static void STLFindEdges (STLGeometry & geom,
 
 
 
-void STLSurfaceMeshing1 (STLGeometry & geom, class Mesh & mesh,
-			 int retrynr);
+void STLSurfaceMeshing1 (STLGeometry & geom, class Mesh & mesh, const MeshingParameters& mparam,
+			 int retrynr, const STLParameters& stlparam);
 
 
-int STLSurfaceMeshing (STLGeometry & geom, class Mesh & mesh)
+int STLSurfaceMeshing (STLGeometry & geom, class Mesh & mesh, const MeshingParameters& mparam,
+                       const STLParameters& stlparam)
 {
   PrintFnStart("Do Surface Meshing");
 
   geom.PrepareSurfaceMeshing();
 
   if (mesh.GetNSeg() == 0)
-    STLFindEdges (geom, mesh);
+    STLFindEdges (geom, mesh, mparam, stlparam);
 
   int nopen;
   int outercnt = 20;
@@ -272,7 +274,7 @@ int STLSurfaceMeshing (STLGeometry & geom, class Mesh & mesh)
 	      if (multithread.terminate) { return MESHING3_TERMINATE; }
 
 	      trialcnt++;
-	      STLSurfaceMeshing1 (geom, mesh, trialcnt);
+	      STLSurfaceMeshing1 (geom, mesh, mparam, trialcnt, stlparam);
 
 	      mesh.FindOpenSegments();
 	      nopen = mesh.GetNOpenSegments();
@@ -526,8 +528,10 @@ int STLSurfaceMeshing (STLGeometry & geom, class Mesh & mesh)
 
 
 void STLSurfaceMeshing1 (STLGeometry & geom,
-			 class Mesh & mesh,
-			 int retrynr)
+			 Mesh & mesh,
+                         const MeshingParameters& mparam,
+			 int retrynr,
+                         const STLParameters& stlparam)
 {
   static int timer1 = NgProfiler::CreateTimer ("STL surface meshing1");
   static int timer1a = NgProfiler::CreateTimer ("STL surface meshing1a");
@@ -740,8 +744,8 @@ void STLSurfaceMeshing1 (STLGeometry & geom,
 
 
 void STLSurfaceOptimization (STLGeometry & geom,
-			     class Mesh & mesh,
-			     MeshingParameters & meshparam)
+			     Mesh & mesh,
+			     const MeshingParameters & mparam)
 {
   PrintFnStart("optimize STL Surface");
 
@@ -749,12 +753,12 @@ void STLSurfaceOptimization (STLGeometry & geom,
 
   optmesh.SetFaceIndex (0);
   optmesh.SetImproveEdges (0);
-  optmesh.SetMetricWeight (meshparam.elsizeweight);
+  optmesh.SetMetricWeight (mparam.elsizeweight);
 
-  PrintMessage(5,"optimize string = ", meshparam.optimize2d, " elsizew = ", meshparam.elsizeweight);
+  PrintMessage(5,"optimize string = ", mparam.optimize2d, " elsizew = ", mparam.elsizeweight);
 
-  for (int i = 1; i <= meshparam.optsteps2d; i++)
-    for (size_t j = 1; j <= meshparam.optimize2d.length(); j++)
+  for (int i = 1; i <= mparam.optsteps2d; i++)
+    for (size_t j = 1; j <= mparam.optimize2d.length(); j++)
       {
 	if (multithread.terminate)
 	  break;
@@ -762,7 +766,7 @@ void STLSurfaceOptimization (STLGeometry & geom,
 	//(*testout) << "optimize, before, step = " << meshparam.optimize2d[j-1] << mesh.Point (3679) << endl;
 
 	mesh.CalcSurfacesOfNode();
-	switch (meshparam.optimize2d[j-1])
+	switch (mparam.optimize2d[j-1])
 	  {
 	  case 's': 
 	    {
diff --git a/libsrc/stlgeom/python_stl.cpp b/libsrc/stlgeom/python_stl.cpp
index 39a52273..d98707ca 100644
--- a/libsrc/stlgeom/python_stl.cpp
+++ b/libsrc/stlgeom/python_stl.cpp
@@ -17,6 +17,118 @@ namespace netgen
   extern shared_ptr<NetgenGeometry> ng_geometry;
 }
 
+static string stlparameter_description = R"delimiter(
+STL Specific Meshing Parameters
+-------------------------------
+
+yangle: float =
+  Angle for edge detection
+
+contyangle: float =
+  Edges continue if angle > contyangle
+
+edgecornerangle: float =
+  Angle of geometry edge at which the mesher should set a point.
+
+)delimiter";
+
+void CreateSTLParametersFromKwargs(STLParameters& stlparam, py::dict kwargs)
+{
+  if(kwargs.contains("yangle"))
+    stlparam.yangle = py::cast<double>(kwargs.attr("pop")("yangle"));
+  if(kwargs.contains("contyangle"))
+    stlparam.contyangle = py::cast<double>(kwargs.attr("pop")("contyangle"));
+  if(kwargs.contains("edgecornerangle"))
+    stlparam.edgecornerangle = py::cast<double>(kwargs.attr("pop")("edgecornerangle"));
+  if(kwargs.contains("chartangle"))
+    stlparam.chartangle = py::cast<double>(kwargs.attr("pop")("chartangle"));
+  if(kwargs.contains("outerchartangle"))
+    stlparam.outerchartangle = py::cast<double>(kwargs.attr("pop")("outerchartangle"));
+  if(kwargs.contains("usesearchtree"))
+    stlparam.usesearchtree = py::cast<int>(kwargs.attr("pop")("usesearchtree"));
+  if(kwargs.contains("resthatlasfac"))
+  {
+    auto val = kwargs.attr("pop")("resthatlasfac");
+    if(val.is_none())
+      stlparam.resthatlasenable = false;
+    else
+    {
+      stlparam.resthatlasenable = true;
+      stlparam.resthatlasfac = py::cast<double>(val);
+    }
+  }
+  if(kwargs.contains("atlasminh"))
+    stlparam.atlasminh = py::cast<double>(kwargs.attr("pop")("atlasminh"));
+  if(kwargs.contains("resthsurfcurvfac"))
+  {
+    auto val = kwargs.attr("pop")("resthsurfcurvfac");
+    if(val.is_none())
+      stlparam.resthsurfcurvenable = false;
+    else
+    {
+      stlparam.resthsurfcurvenable = true;
+      stlparam.resthsurfcurvfac = py::cast<double>(val);
+    }
+  }
+  if(kwargs.contains("resthchartdistfac"))
+  {
+    auto val = kwargs.attr("pop")("resthchartdistfac");
+    if(val.is_none())
+      stlparam.resthchartdistenable = false;
+    else
+    {
+      stlparam.resthchartdistenable = true;
+      stlparam.resthchartdistfac = py::cast<double>(val);
+    }
+  }
+  if(kwargs.contains("resthcloseedgefac"))
+  {
+    auto val = kwargs.attr("pop")("resthcloseedgefac");
+    if(val.is_none())
+      stlparam.resthcloseedgeenable = false;
+    else
+    {
+      stlparam.resthcloseedgeenable = true;
+      stlparam.resthcloseedgefac = py::cast<double>(val);
+    }
+  }
+  if(kwargs.contains("resthedgeanglefac"))
+  {
+    auto val = kwargs.attr("pop")("resthedgeanglefac");
+    if(val.is_none())
+      stlparam.resthedgeangleenable = false;
+    else
+    {
+      stlparam.resthedgeangleenable = true;
+      stlparam.resthedgeanglefac = py::cast<double>(val);
+    }
+  }
+  if(kwargs.contains("resthsurfmeshcurvfac"))
+  {
+    auto val = kwargs.attr("pop")("resthsurfmeshcurvfac");
+    if(val.is_none())
+      stlparam.resthsurfmeshcurvenable = false;
+    else
+    {
+      stlparam.resthsurfmeshcurvenable = true;
+      stlparam.resthsurfmeshcurvfac = py::cast<double>(val);
+    }
+  }
+  if(kwargs.contains("resthlinelengthfac"))
+  {
+    auto val = kwargs.attr("pop")("resthlinelengthfac");
+    if(val.is_none())
+      stlparam.resthlinelengthenable = false;
+    else
+    {
+      stlparam.resthlinelengthenable = true;
+      stlparam.resthlinelengthfac = py::cast<double>(val);
+    }
+  }
+  if(kwargs.contains("recalc_h_opt"))
+    stlparam.recalc_h_opt = py::cast<bool>(kwargs.attr("pop")("recalc_h_opt"));
+}
+
 
 DLL_HEADER void ExportSTL(py::module & m)
 {
@@ -82,20 +194,27 @@ DLL_HEADER void ExportSTL(py::module & m)
                              MeshingParameters* pars, py::kwargs kwargs)
                          {
                            MeshingParameters mp;
-                           if(pars) mp = *pars;
-                           {
-                             py::gil_scoped_acquire aq;
-                             CreateMPfromKwargs(mp, kwargs);
+                           STLParameters stlparam;
+                           { py::gil_scoped_acquire aq;
+                             if(pars)
+                             {
+                               auto mp_flags = pars->geometrySpecificParameters;
+                               auto mp_kwargs = CreateDictFromFlags(mp_flags);
+                               CreateSTLParametersFromKwargs(stlparam, mp_kwargs);
+                               mp = *pars;
+                             }
+                             CreateSTLParametersFromKwargs(stlparam, kwargs);
+                             CreateMPfromKwargs(mp, kwargs); // this will throw if any kwargs are not passed
                            }
                            auto mesh = make_shared<Mesh>();
-                           SetGlobalMesh(mesh);
                            mesh->SetGeometry(geo);
                            ng_geometry = geo;
-                           geo->GenerateMesh(mesh,mp);
+                           SetGlobalMesh(mesh);
+                           STLMeshingDummy(geo.get(), mesh, mp, stlparam);
                            return mesh;
                          }, py::arg("mp") = nullptr,
       py::call_guard<py::gil_scoped_release>(),
-      meshingparameter_description.c_str())
+         (meshingparameter_description + stlparameter_description).c_str())
     ;
   m.def("LoadSTLGeometry", [] (const string & filename)
                            {
diff --git a/libsrc/stlgeom/stlgeom.cpp b/libsrc/stlgeom/stlgeom.cpp
index 2ccbf5f3..253ff02f 100644
--- a/libsrc/stlgeom/stlgeom.cpp
+++ b/libsrc/stlgeom/stlgeom.cpp
@@ -13,15 +13,17 @@ int usechartnormal = 1;
 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
 void STLMeshing (STLGeometry & geom,
-		 Mesh & mesh)
+		 Mesh & mesh,
+                 const MeshingParameters& mparam,
+                 const STLParameters& stlpar)
 {
   geom.Clear();
-  geom.BuildEdges();
-  geom.MakeAtlas(mesh);
+  geom.BuildEdges(stlpar);
+  geom.MakeAtlas(mesh, mparam, stlpar);
   if (multithread.terminate) { return; }
   geom.CalcFaceNums();
   geom.AddFaceEdges();
-  geom.LinkEdges();
+  geom.LinkEdges(stlpar);
 
   mesh.ClearFaceDescriptors();
   for (int i = 1; i <= geom.GetNOFaces(); i++)
@@ -93,9 +95,11 @@ void STLGeometry :: Save (string filename) const
 
 
 
+DLL_HEADER extern STLParameters stlparam;
 int STLGeometry :: GenerateMesh (shared_ptr<Mesh> & mesh, MeshingParameters & mparam)
 {
-  return STLMeshingDummy (this, mesh, mparam);
+  STLParameters stlpar = stlparam;
+  return STLMeshingDummy (this, mesh, mparam, stlpar);
 }
 
 
@@ -132,7 +136,7 @@ void STLGeometry :: STLInfo(double* data)
   data[7] = cons;
 }
 
-void STLGeometry :: MarkNonSmoothNormals()
+void STLGeometry :: MarkNonSmoothNormals(const STLParameters& stlparam)
 {
 
   PrintFnStart("Mark Non-Smooth Normals");
@@ -169,13 +173,13 @@ void STLGeometry :: MarkNonSmoothNormals()
 
 }
 
-void STLGeometry :: SmoothNormals()
+void STLGeometry :: SmoothNormals(const STLParameters& stlparam)
 {
   multithread.terminate = 0;
 
   //  UseExternalEdges();
 
-  BuildEdges();
+  BuildEdges(stlparam);
 
 
   DenseMatrix m(3), hm(3);
@@ -1240,13 +1244,13 @@ void STLGeometry :: ClearEdges()
 
 }
 
-void STLGeometry :: STLDoctorBuildEdges()
+void STLGeometry :: STLDoctorBuildEdges(const STLParameters& stlparam)
 {
   //  if (!trigsconverted) {return;}
   ClearEdges();
 
   meshlines.SetSize(0);
-  FindEdgesFromAngles();
+  FindEdgesFromAngles(stlparam);
 }
 
 void STLGeometry :: DeleteExternalEdgeAtSelected()
@@ -1737,7 +1741,7 @@ void STLGeometry :: InitMarkedTrigs()
     }
 }
 
-void STLGeometry :: MarkDirtyTrigs()
+void STLGeometry :: MarkDirtyTrigs(const STLParameters& stlparam)
 {
   PrintFnStart("mark dirty trigs");
   int i,j;
@@ -1813,12 +1817,12 @@ double STLGeometry :: CalcTrigBadness(int i)
 
 }
 
-void STLGeometry :: GeomSmoothRevertedTrigs()
+void STLGeometry :: GeomSmoothRevertedTrigs(const STLParameters& stlparam)
 {
   //double revertedangle = stldoctor.smoothangle/180.*M_PI;
   double fact = stldoctor.dirtytrigfact;
 
-  MarkRevertedTrigs();
+  MarkRevertedTrigs(stlparam);
 
   int i, j, k, l, p;
 
@@ -1860,13 +1864,13 @@ void STLGeometry :: GeomSmoothRevertedTrigs()
 	    }
 	}
     }
-  MarkRevertedTrigs();
+  MarkRevertedTrigs(stlparam);
 }
 
-void STLGeometry :: MarkRevertedTrigs()
+void STLGeometry :: MarkRevertedTrigs(const STLParameters& stlparam)
 {
   int i,j;
-  if (edgesperpoint.Size() != GetNP()) {BuildEdges();}
+  if (edgesperpoint.Size() != GetNP()) {BuildEdges(stlparam);}
 
   PrintFnStart("mark reverted trigs");
 
@@ -1906,11 +1910,11 @@ void STLGeometry :: MarkRevertedTrigs()
 
 }
 
-void STLGeometry :: SmoothDirtyTrigs()
+void STLGeometry :: SmoothDirtyTrigs(const STLParameters& stlparam)
 {
   PrintFnStart("smooth dirty trigs");
 
-  MarkDirtyTrigs();
+  MarkDirtyTrigs(stlparam);
 
   int i,j;
   int changed = 1;
@@ -1953,7 +1957,7 @@ void STLGeometry :: SmoothDirtyTrigs()
   calcedgedataanglesnew = 1;
 
 
-  MarkDirtyTrigs();
+  MarkDirtyTrigs(stlparam);
 
   int cnt = 0;
   for (i = 1; i <= GetNT(); i++)
@@ -2360,12 +2364,12 @@ int STLGeometry :: IsEdgeNum(int ap1, int ap2)
 }
 
 
-void STLGeometry :: BuildEdges()
+void STLGeometry :: BuildEdges(const STLParameters& stlparam)
 {
   //PrintFnStart("build edges");
   edges.SetSize(0);
   meshlines.SetSize(0);
-  FindEdgesFromAngles();
+  FindEdgesFromAngles(stlparam);
 }
 
 void STLGeometry :: UseExternalEdges()
@@ -2487,7 +2491,7 @@ void STLGeometry :: CalcEdgeDataAngles()
   PrintMessage (5,"calc edge data angles ... done");
 }
 
-void STLGeometry :: FindEdgesFromAngles()
+void STLGeometry :: FindEdgesFromAngles(const STLParameters& stlparam)
 {
   //  PrintFnStart("find edges from angles");
 
@@ -2714,7 +2718,7 @@ void STLGeometry :: AddFaceEdges()
   
 }
 
-void STLGeometry :: LinkEdges()
+void STLGeometry :: LinkEdges(const STLParameters& stlparam)
 {
   PushStatusF("Link Edges");
   PrintMessage(5,"have now ", GetNE(), " edges with yellow angle = ", stlparam.yangle, " degree");
@@ -3131,7 +3135,7 @@ int IsInArray(int n, const NgArray<int>& ia)
 }
 */
 
-void STLGeometry :: AddConeAndSpiralEdges()
+void STLGeometry :: AddConeAndSpiralEdges(const STLParameters& stlparam)
 {
   PrintMessage(5,"have now ", GetNE(), " edges with yellow angle = ", stlparam.yangle, " degree");
 
diff --git a/libsrc/stlgeom/stlgeom.hpp b/libsrc/stlgeom/stlgeom.hpp
index 5d7a5af1..7f31c151 100644
--- a/libsrc/stlgeom/stlgeom.hpp
+++ b/libsrc/stlgeom/stlgeom.hpp
@@ -40,7 +40,7 @@ namespace netgen
   }
   */
   
-  extern DLL_HEADER MeshingParameters mparam;
+// extern DLL_HEADER MeshingParameters mparam;
   
 
 
@@ -184,20 +184,20 @@ namespace netgen
     STLGeometry();
     virtual ~STLGeometry();
 
-    void DoArchive(Archive& ar)
+    void DoArchive(Archive& ar) override
     {
       STLTopology::DoArchive(ar);
     }
 
     void Clear();
 
-    virtual void Save (string filename) const;
+    virtual void Save (string filename) const override;
 
 
 	DLL_HEADER void STLInfo(double* data);
     //stldoctor:
-	DLL_HEADER void SmoothNormals();
-	DLL_HEADER void MarkNonSmoothNormals();
+	DLL_HEADER void SmoothNormals(const STLParameters& stlparam);
+	DLL_HEADER void MarkNonSmoothNormals(const STLParameters& stlparam);
 
 	DLL_HEADER void CalcEdgeData();
 	DLL_HEADER void CalcEdgeDataAngles();
@@ -251,7 +251,7 @@ namespace netgen
 	DLL_HEADER void AddClosedLinesToExternalEdges();
 	DLL_HEADER void AddLongLinesToExternalEdges();
 	DLL_HEADER void AddAllNotSingleLinesToExternalEdges();
-	DLL_HEADER void STLDoctorBuildEdges();
+	DLL_HEADER void STLDoctorBuildEdges(const STLParameters& stlparam);
 	DLL_HEADER void AddExternalEdgesFromGeomLine();
 	DLL_HEADER void DeleteDirtyExternalEdges();
 	DLL_HEADER void DeleteExternalEdgeAtSelected();
@@ -292,10 +292,10 @@ namespace netgen
 	DLL_HEADER int Vicinity(int trig) const;
 
 	DLL_HEADER void InitMarkedTrigs();
-	DLL_HEADER void MarkDirtyTrigs();
-	DLL_HEADER void SmoothDirtyTrigs();
-	DLL_HEADER void GeomSmoothRevertedTrigs();
-	DLL_HEADER void MarkRevertedTrigs();
+	DLL_HEADER void MarkDirtyTrigs(const STLParameters& stlparam);
+	DLL_HEADER void SmoothDirtyTrigs(const STLParameters& stlparam);
+	DLL_HEADER void GeomSmoothRevertedTrigs(const STLParameters& stlparam);
+	DLL_HEADER void MarkRevertedTrigs(const STLParameters& stlparam);
 	DLL_HEADER double CalcTrigBadness(int i);
 	DLL_HEADER int IsMarkedTrig(int trig) const;
 	DLL_HEADER void SetMarkedTrig(int trig, int num);
@@ -330,8 +330,8 @@ namespace netgen
     ///
 
     ///ReadTriangle->STLTriangle, initialise some important variables, always after load!!!
-    virtual void InitSTLGeometry (const NgArray<STLReadTriangle> & readtrigs);
-    virtual void TopologyChanged(); //do some things, if topology changed!
+    virtual void InitSTLGeometry (const NgArray<STLReadTriangle> & readtrigs) override;
+    virtual void TopologyChanged() override; //do some things, if topology changed!
     int CheckGeometryOverlapping();
 
     //get NO edges per point
@@ -353,18 +353,18 @@ namespace netgen
 
     ///Build EdgeSegments
     void ClearEdges();
-    void BuildEdges();
+    void BuildEdges(const STLParameters& stlparam);
     void BuildEdgesPerPoint();
     void UseExternalEdges();
 
 
-    void FindEdgesFromAngles();
+    void FindEdgesFromAngles(const STLParameters& stlparam);
     void CalcFaceNums();
     int GetNOBodys();
     int GetNOFaces() {return facecnt;}
-    void LinkEdges();
+    void LinkEdges(const STLParameters& stlparam);
 
-    void AddConeAndSpiralEdges();
+    void AddConeAndSpiralEdges(const STLParameters& stlparam);
     void AddFaceEdges(); //each face should have at least one starting edge (outherwise it won't be meshed)
 
     void GetDirtyChartTrigs(int chartnum, STLChart& chart, const NgArray<int>& outercharttrigs, 
@@ -382,7 +382,7 @@ namespace netgen
 
 
     //make charts with regions of a max. angle
-    void MakeAtlas(class Mesh & mesh);
+    void MakeAtlas(class Mesh & mesh, const MeshingParameters& mparam, const STLParameters& stlparam);
 
     //outerchartspertrig, sorted!
     int GetOCPTSize() const {return outerchartspertrig.Size();};
@@ -450,17 +450,16 @@ namespace netgen
     int LineEndPointsSet() const {return lineendpoints.Size() == GetNP();}
     void ClearLineEndPoints();
 
-	DLL_HEADER void RestrictLocalH(class Mesh & mesh, double gh);
-    void RestrictLocalHCurv(class Mesh & mesh, double gh);
+    DLL_HEADER void RestrictLocalH(class Mesh & mesh, double gh, const STLParameters& stlparam);
+    void RestrictLocalHCurv(class Mesh & mesh, double gh, const STLParameters& stlparam);
     void RestrictHChartDistOneChart(int chartnum, NgArray<int>& acttrigs, class Mesh & mesh, 
-				    double gh, double fact, double minh);
+				    double gh, double fact, double minh, const STLParameters& stlparam);
 
     friend class MeshingSTLSurface;
 
-
-    virtual int GenerateMesh (shared_ptr<Mesh> & mesh, MeshingParameters & mparam);
+    int GenerateMesh (shared_ptr<Mesh> & mesh, MeshingParameters & mparam) override;
     
-    virtual const Refinement & GetRefinement () const;
+    virtual const Refinement & GetRefinement () const override;
   };
  
 
@@ -468,7 +467,8 @@ namespace netgen
 
 
 
-  extern int STLMeshingDummy (STLGeometry* stlgeometry, shared_ptr<Mesh> & mesh, MeshingParameters & mparam);
+extern int STLMeshingDummy (STLGeometry* stlgeometry, shared_ptr<Mesh> & mesh, const MeshingParameters & mparam,
+                            const STLParameters& stlpar);
 
 
 }
diff --git a/libsrc/stlgeom/stlgeomchart.cpp b/libsrc/stlgeom/stlgeomchart.cpp
index 4186265a..b913ebde 100644
--- a/libsrc/stlgeom/stlgeomchart.cpp
+++ b/libsrc/stlgeom/stlgeomchart.cpp
@@ -17,7 +17,7 @@ int chartdebug = 0;
 
 
 
-void STLGeometry :: MakeAtlas(Mesh & mesh)
+void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, const STLParameters& stlparam)
 {
   // int timer1 = NgProfiler::CreateTimer ("makeatlas");
   /*
@@ -128,7 +128,7 @@ void STLGeometry :: MakeAtlas(Mesh & mesh)
 
       SetThreadPercent(100.0 * workedarea / atlasarea);
 
-      STLChart * chart = new STLChart(this);
+      STLChart * chart = new STLChart(this, stlparam);
       atlas.Append(chart);
 
       // *testout << "Chart " << atlas.Size() << endl;
@@ -520,7 +520,7 @@ void STLGeometry :: MakeAtlas(Mesh & mesh)
       // char key;
       // cin >> key;
       //calculate an estimate meshsize, not to produce too large outercharts, with factor 2 larger!
-      RestrictHChartDistOneChart(chartnum, chartdistacttrigs, mesh, h, 0.5, atlasminh);
+      RestrictHChartDistOneChart(chartnum, chartdistacttrigs, mesh, h, 0.5, atlasminh, stlparam);
       // NgProfiler::Print(stdout);
       // NgProfiler::StopTimer (timere2);      
       
@@ -572,7 +572,7 @@ void STLGeometry :: MakeAtlas(Mesh & mesh)
   mesh.SetMinimalH(mparam.minh);
   
   
-  AddConeAndSpiralEdges();
+  AddConeAndSpiralEdges(stlparam);
   
   PrintMessage(5,"Make Atlas finished");
 
diff --git a/libsrc/stlgeom/stlgeommesh.cpp b/libsrc/stlgeom/stlgeommesh.cpp
index a5e8cbcf..a5fcc0c2 100644
--- a/libsrc/stlgeom/stlgeommesh.cpp
+++ b/libsrc/stlgeom/stlgeommesh.cpp
@@ -705,7 +705,7 @@ int STLGeometry :: ProjectNearest(Point<3> & p3d) const
 
 	
 //Restrict local h due to curvature for make atlas
-void STLGeometry :: RestrictLocalHCurv(class Mesh & mesh, double gh)
+void STLGeometry :: RestrictLocalHCurv(class Mesh & mesh, double gh, const STLParameters& stlparam)
 {
   PushStatusF("Restrict H due to surface curvature");
 
@@ -810,7 +810,7 @@ void STLGeometry :: RestrictLocalHCurv(class Mesh & mesh, double gh)
 
 }
   //restrict local h due to near edges and due to outer chart distance
-void STLGeometry :: RestrictLocalH(class Mesh & mesh, double gh)
+void STLGeometry :: RestrictLocalH(class Mesh & mesh, double gh, const STLParameters& stlparam)
 {
   
   //bei jedem Dreieck alle Nachbardreiecke vergleichen, und, fallskein Kante dazwischen,
@@ -1077,7 +1077,7 @@ void STLGeometry :: RestrictLocalH(class Mesh & mesh, double gh)
 	  if (multithread.terminate)
 	    {PopStatus(); return;}
 
-	  RestrictHChartDistOneChart(i, acttrigs, mesh, gh, 1., 0.);
+	  RestrictHChartDistOneChart(i, acttrigs, mesh, gh, 1., 0., stlparam);
 	}
       
       PopStatus();
@@ -1117,7 +1117,8 @@ void STLGeometry :: RestrictLocalH(class Mesh & mesh, double gh)
 }
 
 void STLGeometry :: RestrictHChartDistOneChart(int chartnum, NgArray<int>& acttrigs, 
-					       class Mesh & mesh, double gh, double fact, double minh)
+					       class Mesh & mesh, double gh, double fact, double minh,
+                                               const STLParameters& stlparam)
 {
   static int timer1 = NgProfiler::CreateTimer ("restrictH OneChart 1");
   static int timer2 = NgProfiler::CreateTimer ("restrictH OneChart 2");
@@ -1345,7 +1346,8 @@ void STLGeometry :: RestrictHChartDistOneChart(int chartnum, NgArray<int>& acttr
 }
 
 
-int STLMeshingDummy (STLGeometry* stlgeometry, shared_ptr<Mesh> & mesh, MeshingParameters & mparam)
+int STLMeshingDummy (STLGeometry* stlgeometry, shared_ptr<Mesh> & mesh, const MeshingParameters & mparam,
+                     const STLParameters& stlparam)
 {
   if (mparam.perfstepsstart > mparam.perfstepsend) return 0;
 
@@ -1372,7 +1374,7 @@ int STLMeshingDummy (STLGeometry* stlgeometry, shared_ptr<Mesh> & mesh, MeshingP
   
       //mesh->DeleteMesh();
  
-      STLMeshing (*stlgeometry, *mesh);
+      STLMeshing (*stlgeometry, *mesh, mparam, stlparam);
 
       stlgeometry->edgesfound = 1;
       stlgeometry->surfacemeshed = 0;
@@ -1399,7 +1401,7 @@ int STLMeshingDummy (STLGeometry* stlgeometry, shared_ptr<Mesh> & mesh, MeshingP
 	}
 
       success = 0;
-      int retval = STLSurfaceMeshing (*stlgeometry, *mesh);
+      int retval = STLSurfaceMeshing (*stlgeometry, *mesh, mparam, stlparam);
       if (retval == MESHING3_OK)
 	{
 	  PrintMessage(3,"Success !!!!");
@@ -1471,13 +1473,14 @@ int STLMeshingDummy (STLGeometry* stlgeometry, shared_ptr<Mesh> & mesh, MeshingP
 	      mesh -> LoadLocalMeshSize (mparam.meshsizefilename);	      
 	      mesh -> CalcLocalHFromSurfaceCurvature (mparam.grading, 
 						      stlparam.resthsurfmeshcurvfac);
-	      mparam.optimize2d = "cmsmSm";
-	      STLSurfaceOptimization (*stlgeometry, *mesh, mparam);
+              MeshingParameters mpar = mparam;
+	      mpar.optimize2d = "cmsmSm";
+	      STLSurfaceOptimization (*stlgeometry, *mesh, mpar);
 #ifdef STAT_STREAM
 	      (*statout) << GetTime() << " & ";
 #endif
 
-	      mparam.Render();
+	      mpar.Render();
 	    }
 	  stlgeometry->surfaceoptimized = 1;
 	}
diff --git a/libsrc/stlgeom/stlpkg.cpp b/libsrc/stlgeom/stlpkg.cpp
index a016a71b..0858fa68 100644
--- a/libsrc/stlgeom/stlpkg.cpp
+++ b/libsrc/stlgeom/stlpkg.cpp
@@ -21,6 +21,8 @@ namespace netgen
 {
   DLL_HEADER extern shared_ptr<NetgenGeometry>  ng_geometry;
   DLL_HEADER extern shared_ptr<Mesh> mesh;
+  DLL_HEADER extern MeshingParameters mparam;
+  DLL_HEADER extern STLParameters stlparam;
 
   static VisualSceneSTLGeometry vsstlgeom;
   static VisualSceneSTLMeshing vsstlmeshing;
@@ -243,15 +245,15 @@ namespace netgen
 	  }
 	else if (strcmp (argv[1], "markdirtytrigs") == 0)
 	  {
-	    stlgeometry->MarkDirtyTrigs();
+	    stlgeometry->MarkDirtyTrigs(stlparam);
 	  }
 	else if (strcmp (argv[1], "smoothdirtytrigs") == 0)
 	  {
-	    stlgeometry->SmoothDirtyTrigs();
+	    stlgeometry->SmoothDirtyTrigs(stlparam);
 	  }
 	else if (strcmp (argv[1], "smoothrevertedtrigs") == 0)
 	  {
-	    stlgeometry->GeomSmoothRevertedTrigs();
+	    stlgeometry->GeomSmoothRevertedTrigs(stlparam);
 	  }
 	else if (strcmp (argv[1], "invertselectedtrig") == 0)
 	  {
@@ -306,11 +308,11 @@ namespace netgen
 	  }
 	else if (strcmp (argv[1], "smoothnormals") == 0)
 	  {
-	    stlgeometry->SmoothNormals();
+	    stlgeometry->SmoothNormals(stlparam);
 	  }
 	else if (strcmp (argv[1], "marknonsmoothnormals") == 0)
 	  {
-	    stlgeometry->MarkNonSmoothNormals();
+	    stlgeometry->MarkNonSmoothNormals(stlparam);
 	  }
 	else if (strcmp (argv[1], "addexternaledge") == 0)
 	  {
@@ -359,7 +361,7 @@ namespace netgen
 	  }
 	else if (strcmp (argv[1], "buildedges") == 0)
 	  {
-	    stlgeometry->STLDoctorBuildEdges();
+	    stlgeometry->STLDoctorBuildEdges(stlparam);
 	  }
 	else if (strcmp (argv[1], "confirmedge") == 0)
 	  {
@@ -526,7 +528,7 @@ namespace netgen
 	mesh -> SetLocalH (stlgeometry->GetBoundingBox().PMin() - Vec3d(10, 10, 10),
 			   stlgeometry->GetBoundingBox().PMax() + Vec3d(10, 10, 10),
 			   mparam.grading);
-	stlgeometry -> RestrictLocalH(*mesh, mparam.maxh);
+	stlgeometry -> RestrictLocalH(*mesh, mparam.maxh, stlparam);
 
 	if (stlparam.resthsurfmeshcurvenable)
 	  mesh -> CalcLocalHFromSurfaceCurvature (mparam.grading, 
diff --git a/libsrc/stlgeom/stltool.cpp b/libsrc/stlgeom/stltool.cpp
index 5a2eb61e..06ab66ac 100644
--- a/libsrc/stlgeom/stltool.cpp
+++ b/libsrc/stlgeom/stltool.cpp
@@ -607,7 +607,8 @@ STLTopEdge :: STLTopEdge (int p1, int p2, int trig1, int trig2)
 //+++++++++++++++++++   STL CHART   +++++++++++++++++++++++++++++++
 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
-STLChart :: STLChart(STLGeometry * ageometry)
+STLChart :: STLChart(STLGeometry * ageometry, const STLParameters& astlparam)
+    : geometry(ageometry), stlparam(astlparam)
 {
   // charttrigs = new NgArray<int> (0,0);
   // outertrigs = new NgArray<int> (0,0);
@@ -1462,7 +1463,6 @@ void STLParameters :: Print (ostream & ost) const
 }
 
 
+DLL_HEADER extern STLParameters stlparam;
 STLParameters stlparam;
-
-
 }
diff --git a/libsrc/stlgeom/stltool.hpp b/libsrc/stlgeom/stltool.hpp
index b602ad94..19ff12cc 100644
--- a/libsrc/stlgeom/stltool.hpp
+++ b/libsrc/stlgeom/stltool.hpp
@@ -38,6 +38,7 @@ extern void FIOWriteString(ostream& ios, char* str, int len);
 typedef NgArray <int> * ArrayINTPTR;
 
 class STLGeometry;
+class STLParameters;
 
 class STLChart
 {
@@ -49,11 +50,12 @@ private:
 
   NgArray<twoint> olimit; //outer limit of outer chart
   NgArray<twoint> ilimit; //outer limit of inner chart
+  const STLParameters& stlparam;
 
 
 public:
   
-  STLChart(STLGeometry * ageometry);
+  STLChart(STLGeometry * ageometry, const STLParameters& astlparam);
   ~STLChart();
   void AddChartTrig(int i);
   void AddOuterTrig(int i);
@@ -227,64 +229,67 @@ DLL_HEADER extern STLDoctorParams stldoctor;
 
 
 
-class STLParameters
+// TODO change enable flag to optional parameters
+class DLL_HEADER STLParameters
 {
 public:
   /// angle for edge detection
-  double yangle;
+  double yangle = 30.;
   double contyangle; //edges continued with contyangle
   /// angle of geometry edge at which the mesher should set a point
-  double edgecornerangle;
+  double edgecornerangle = 60.;
   /// angle inside on chart
-  double chartangle;
+  double chartangle = 15.;
   /// angle for overlapping parts of char
-  double outerchartangle;
+  double outerchartangle = 70.;
   /// 0 .. no, 1 .. local, (2 .. global)
   int usesearchtree;
   ///
   double resthatlasfac; 
-  int resthatlasenable;
+  bool resthatlasenable;
   double atlasminh;
 
-  double resthsurfcurvfac; 
-  int resthsurfcurvenable;
+  double resthsurfcurvfac = 1.; 
+  bool resthsurfcurvenable = false;
 
-  double resthchartdistfac;
-  int resthchartdistenable;
+  double resthchartdistfac = 1.5;
+  bool resthchartdistenable = true;
 
-  double resthcloseedgefac;
-  int resthcloseedgeenable;
+  double resthcloseedgefac = 2.;
+  bool resthcloseedgeenable = true;
   
-  double resthedgeanglefac;
-  int resthedgeangleenable;
+  double resthedgeanglefac = 1.;
+  bool resthedgeangleenable = false;
   
-  double resthsurfmeshcurvfac;
-  int resthsurfmeshcurvenable;
+  double resthsurfmeshcurvfac = 2.;
+  bool resthsurfmeshcurvenable = false;
   
-  double resthlinelengthfac;
-  int resthlinelengthenable;
+  double resthlinelengthfac = 0.5;
+  bool resthlinelengthenable = true;
 
   ///
-  int recalc_h_opt;
+  bool recalc_h_opt = true;
   ///
   STLParameters();
   ///
   void Print (ostream & ost) const;
 };
 
-DLL_HEADER extern STLParameters stlparam;
-
 
 void STLMeshing (STLGeometry & geom,
-		 class Mesh & mesh);
+		 Mesh & mesh,
+                 const MeshingParameters& mparam,
+                 const STLParameters& stlpar);
 
 
 int STLSurfaceMeshing (STLGeometry & geom,
-			class Mesh & mesh);
+                       Mesh & mesh,
+                       const MeshingParameters& mparam,
+                       const STLParameters& stlpar);
 
 void STLSurfaceOptimization (STLGeometry & geom,
-			     class Mesh & mesh,
-			     class MeshingParameters & mparam);
+			     Mesh & mesh,
+			     const MeshingParameters & mparam);
 
 
 
diff --git a/ng/netgenpy.cpp b/ng/netgenpy.cpp
index a4e7a8d1..d9d8a5c2 100644
--- a/ng/netgenpy.cpp
+++ b/ng/netgenpy.cpp
@@ -23,6 +23,7 @@ void DLL_HEADER ExportNgOCC(py::module &m);
 
 PYBIND11_MODULE(libngpy, ngpy)
 {
+  py::module::import("pyngcore");
     py::module meshing = ngpy.def_submodule("_meshing", "pybind meshing module");
     ExportNetgenMeshing(meshing);
     py::module csg = ngpy.def_submodule("_csg", "pybind csg module");
diff --git a/nglib/nglib.cpp b/nglib/nglib.cpp
index a177b9bd..98e37989 100644
--- a/nglib/nglib.cpp
+++ b/nglib/nglib.cpp
@@ -30,6 +30,8 @@ namespace netgen {
    extern void MeshFromSpline2D (SplineGeometry2d & geometry,
                                  shared_ptr<Mesh> & mesh, 
                                  MeshingParameters & mp);
+   extern MeshingParameters mparam;
+   DLL_HEADER extern STLParameters stlparam;
 }
 
 
@@ -667,7 +669,7 @@ namespace nglib
         }
       */
 
-      STLMeshing (*stlgeometry, *me);
+      STLMeshing (*stlgeometry, *me, mparam, stlparam);
 
       stlgeometry->edgesfound = 1;
       stlgeometry->surfacemeshed = 0;
@@ -709,7 +711,7 @@ namespace nglib
       stlgeometry->surfaceoptimized = 0;
       stlgeometry->volumemeshed = 0;
       */  
-      int retval = STLSurfaceMeshing (*stlgeometry, *me);
+      int retval = STLSurfaceMeshing (*stlgeometry, *me, mparam, stlparam);
       if (retval == MESHING3_OK)
       {
          (*mycout) << "Success !!!!" << endl;
diff --git a/python/meshing.py b/python/meshing.py
index 5fbebdf4..1af409a2 100644
--- a/python/meshing.py
+++ b/python/meshing.py
@@ -1,22 +1,66 @@
 from .libngpy._meshing import *
 
 class _MeshsizeObject:
-    pass
+    @property
+    def very_coarse(self):
+        return MeshingParameters(curvaturesafety=1,
+                                 segmentsperedge=0.3,
+                                 grading=0.7,
+                                 resthsurfcurvfac=0.25,
+                                 resthchartdistfac=0.8,
+                                 resthlinelengthfac=0.2,
+                                 resthcloseedgefac=0.5,
+                                 resthminedgelen=0.002,
+                                 resthedgeanglefac=0.25,
+                                 resthsurfmeshcurvfac=1.)
+    @property
+    def coarse(self):
+        return MeshingParameters(curvaturesafety=1.5,
+                                 segmentsperedge=0.5,
+                                 grading=0.5,
+                                 resthsurfcurvfac=0.5,
+                                 resthchartdistfac=1,
+                                 resthlinelengthfac=0.35,
+                                 resthcloseedgefac=1,
+                                 resthminedgelen=0.02,
+                                 resthedgeanglefac=0.5,
+                                 resthsurfmeshcurvfac=1.5)
+    @property
+    def moderate(self):
+        return MeshingParameters(curvaturesafety=2,
+                                 segmentsperedge=1,
+                                 grading=0.3,
+                                 resthsurfcurvfac=1.,
+                                 resthchartdistfac=1.5,
+                                 resthlinelengthfac=0.5,
+                                 resthcloseedgefac=2,
+                                 resthminedgelen=0.2,
+                                 resthedgeanglefac=1,
+                                 resthsurfmeshcurvfac=2.)
+    @property
+    def fine(self):
+        return MeshingParameters(curvaturesafety=3,
+                                 segmentsperedge=2,
+                                 grading=0.2,
+                                 resthsurfcurvfac=1.5,
+                                 resthchartdistfac=2,
+                                 resthlinelengthfac=1.5,
+                                 resthcloseedgefac=3.5,
+                                 resthminedgelen=1.,
+                                 resthedgeanglefac=1.5,
+                                 resthsurfmeshcurvfac=3.)
+
+    @property
+    def very_fine(self):
+        return  MeshingParameters(curvaturesafety=5,
+                                  segmentsperedge=3,
+                                  grading=0.1,
+                                  resthsurfcurvfac=3,
+                                  resthchartdistfac=5,
+                                  resthlinelengthfac=3,
+                                  resthcloseedgefac=5,
+                                  resthminedgelen=2.,
+                                  resthedgeanglefac=3.,
+                                  resthsurfmeshcurvfac=5.)
 
 meshsize = _MeshsizeObject()
-
-meshsize.very_coarse = MeshingParameters(curvaturesafety=1,
-                                         segmentsperedge=0.3,
-                                         grading=0.7)
-meshsize.coarse = MeshingParameters(curvaturesafety=1.5,
-                                    segmentsperedge=0.5,
-                                    grading=0.5)
-meshsize.moderate = MeshingParameters(curvaturesafety=2,
-                                      segmentsperedge=1,
-                                      grading=0.3)
-meshsize.fine = MeshingParameters(curvaturesafety=3,
-                                  segmentsperedge=2,
-                                  grading=0.2)
-meshsize.very_fine = MeshingParameters(curvaturesafety=5,
-                                       segmentsperedge=3,
-                                       grading=0.1)
diff --git a/tests/pytest/test_savemesh.py b/tests/pytest/test_savemesh.py
index 77ccd198..e98f2a79 100644
--- a/tests/pytest/test_savemesh.py
+++ b/tests/pytest/test_savemesh.py
@@ -28,7 +28,7 @@ def CreateGeo():
 def test_BBNDsave():
     mesh = CreateGeo().GenerateMesh(maxh=0.4,perfstepsend = meshing.MeshingStep.MESHSURFACE)
     for i in range(2):
-        mesh.GenerateVolumeMesh(only3D_domain=i+1,maxh=0.4)
+        mesh.GenerateVolumeMesh(only3D_domain_nr=i+1,maxh=0.4)
     mesh.SetGeometry(None)
     mesh.Save("test.vol")
     mesh2 = meshing.Mesh()