Merge remote-tracking branch 'remotes/origin/lh_multi_material_segmentation' into vb_print_regions

39 files changed, 2690 insertions, 145 deletions

diff --git a/src/libslic3r/CMakeLists.txt b/src/libslic3r/CMakeLists.txt
index 6357011c1..47e9ff5ae 100644
--- a/src/libslic3r/CMakeLists.txt
+++ b/src/libslic3r/CMakeLists.txt
@@ -132,6 +132,8 @@ add_library(libslic3r STATIC
     Model.hpp
     ModelArrange.hpp
     ModelArrange.cpp
+    MultiMaterialSegmentation.cpp
+    MultiMaterialSegmentation.hpp
     CustomGCode.cpp
     CustomGCode.hpp
     Arrange.hpp
diff --git a/src/libslic3r/ClipperUtils.cpp b/src/libslic3r/ClipperUtils.cpp
index 724224604..a49665ae5 100644
--- a/src/libslic3r/ClipperUtils.cpp
+++ b/src/libslic3r/ClipperUtils.cpp
@@ -877,7 +877,7 @@ ClipperLib::Path mittered_offset_path_scaled(const Points &contour, const std::v
 				Vec2d nnext  = perp(ptnext - pt).normalized();
 
 				double delta  = deltas[i];
-				double sin_a  = clamp(-1., 1., cross2(nprev, nnext));
+				double sin_a  = std::clamp(cross2(nprev, nnext), -1., 1.);
 				double convex = sin_a * delta;
 				if (convex <= - sin_min_parallel) {
 					// Concave corner.
diff --git a/src/libslic3r/ElephantFootCompensation.cpp b/src/libslic3r/ElephantFootCompensation.cpp
index 0895e16d6..0adff1ba4 100644
--- a/src/libslic3r/ElephantFootCompensation.cpp
+++ b/src/libslic3r/ElephantFootCompensation.cpp
@@ -270,7 +270,7 @@ std::vector<float> contour_distance2(const EdgeGrid::Grid &grid, const size_t id
 				    const Vec2d   v  = (segment.second - segment.first).cast<double>();
 				    const Vec2d   va = (this->point - segment.first).cast<double>();
 				    const double  l2 = v.squaredNorm();  // avoid a sqrt
-				    const double  t  = (l2 == 0.0) ? 0. : clamp(0., 1., va.dot(v) / l2);
+				    const double  t  = (l2 == 0.0) ? 0. : std::clamp(va.dot(v) / l2, 0., 1.);
 				    // Closest point from this->point to the segment.
 				    const Vec2d   foot = segment.first.cast<double>() + t * v;
 					const Vec2d   bisector = foot - this->point.cast<double>();
diff --git a/src/libslic3r/ExPolygon.cpp b/src/libslic3r/ExPolygon.cpp
index 489023041..506ba8cb6 100644
--- a/src/libslic3r/ExPolygon.cpp
+++ b/src/libslic3r/ExPolygon.cpp
@@ -372,6 +372,27 @@ bool remove_sticks(ExPolygon &poly)
     return remove_sticks(poly.contour) || remove_sticks(poly.holes);
 }
 
+bool remove_small_and_small_holes(ExPolygons &expolygons, double min_area)
+{
+    bool   modified = false;
+    size_t free_idx = 0;
+    for (size_t expoly_idx = 0; expoly_idx < expolygons.size(); ++expoly_idx) {
+        if (std::abs(expolygons[expoly_idx].area()) >= min_area) {
+            // Expolygon is big enough, so also check all its holes
+            modified |= remove_small(expolygons[expoly_idx].holes, min_area);
+            if (free_idx < expoly_idx) {
+                std::swap(expolygons[expoly_idx].contour, expolygons[free_idx].contour);
+                std::swap(expolygons[expoly_idx].holes, expolygons[free_idx].holes);
+            }
+            ++free_idx;
+        } else
+            modified = true;
+    }
+    if (free_idx < expolygons.size())
+        expolygons.erase(expolygons.begin() + free_idx, expolygons.end());
+    return modified;
+}
+
 void keep_largest_contour_only(ExPolygons &polygons)
 {
 	if (polygons.size() > 1) {
diff --git a/src/libslic3r/ExPolygon.hpp b/src/libslic3r/ExPolygon.hpp
index 04b84f767..464310ac0 100644
--- a/src/libslic3r/ExPolygon.hpp
+++ b/src/libslic3r/ExPolygon.hpp
@@ -360,15 +360,11 @@ extern std::vector<BoundingBox> get_extents_vector(const ExPolygons &polygons);
 extern bool        remove_sticks(ExPolygon &poly);
 extern void 	   keep_largest_contour_only(ExPolygons &polygons);
 
-inline double area(const ExPolygon &poly) { return poly.area(); }
+inline double      area(const ExPolygon &poly) { return poly.area(); }
+inline double      area(const ExPolygons &polys) { double s = 0.; for (auto &p : polys) s += p.area(); return s; }
 
-inline double area(const ExPolygons &polys)
-{
-    double s = 0.;
-    for (auto &p : polys) s += p.area();
-
-    return s;
-}
+// Removes all expolygons smaller than min_area and also removes all holes smaller than min_area
+extern bool        remove_small_and_small_holes(ExPolygons &expolygons, double min_area);
 
 } // namespace Slic3r
 
diff --git a/src/libslic3r/ExtrusionSimulator.cpp b/src/libslic3r/ExtrusionSimulator.cpp
index 8c2ab037f..6b1f76abe 100644
--- a/src/libslic3r/ExtrusionSimulator.cpp
+++ b/src/libslic3r/ExtrusionSimulator.cpp
@@ -550,10 +550,10 @@ void gcode_paint_layer(
 		boost::geometry::expand(bboxLine, rect[2]);
 		boost::geometry::expand(bboxLine, rect[3]);
 		B2i bboxLinei(
-			V2i(clamp(0, nc-1, int(floor(bboxLine.min_corner().x()))),
-				clamp(0, nr-1, int(floor(bboxLine.min_corner().y())))),
-			V2i(clamp(0, nc-1, int(ceil (bboxLine.max_corner().x()))),
-				clamp(0, nr-1, int(ceil (bboxLine.max_corner().y())))));
+			V2i(std::clamp(int(floor(bboxLine.min_corner().x())), 0, nc-1),
+				std::clamp(int(floor(bboxLine.min_corner().y())), 0, nr-1)),
+			V2i(std::clamp(int(ceil(bboxLine.max_corner().x())), 0, nc-1),
+				std::clamp(int(ceil(bboxLine.max_corner().y())), 0, nr-1)));
 		// printf("bboxLinei %d,%d %d,%d\n", bboxLinei.min_corner().x(), bboxLinei.min_corner().y(), bboxLinei.max_corner().x(), bboxLinei.max_corner().y());
 #ifdef _DEBUG
 		float area = polyArea(rect, 4);
@@ -597,10 +597,10 @@ void gcode_paint_bitmap(
 		boost::geometry::expand(bboxLine, rect[2]);
 		boost::geometry::expand(bboxLine, rect[3]);
 		B2i bboxLinei(
-			V2i(clamp(0, nc-1, int(floor(bboxLine.min_corner().x()))),
-				clamp(0, nr-1, int(floor(bboxLine.min_corner().y())))),
-			V2i(clamp(0, nc-1, int(ceil (bboxLine.max_corner().x()))),
-				clamp(0, nr-1, int(ceil (bboxLine.max_corner().y())))));
+			V2i(std::clamp(int(floor(bboxLine.min_corner().x())), 0, nc-1),
+				std::clamp(int(floor(bboxLine.min_corner().y())), 0, nr-1)),
+			V2i(std::clamp(int(ceil(bboxLine.max_corner().x())), 0, nc-1),
+				std::clamp(int(ceil(bboxLine.max_corner().y())), 0, nr-1)));
 		// printf("bboxLinei %d,%d %d,%d\n", bboxLinei.min_corner().x(), bboxLinei.min_corner().y(), bboxLinei.max_corner().x(), bboxLinei.max_corner().y());
 		for (int j = bboxLinei.min_corner().y(); j + 1 < bboxLinei.max_corner().y(); ++ j) {
 			for (int i = bboxLinei.min_corner().x(); i + 1 < bboxLinei.max_corner().x(); ++i) {
@@ -664,10 +664,10 @@ void gcode_spread_points(
 		const float height_target = it->height;
 		B2f bbox(center - V2f(radius, radius), center + V2f(radius, radius));
 		B2i bboxi(
-			V2i(clamp(0, nc-1, int(floor(bbox.min_corner().x()))),
-				clamp(0, nr-1, int(floor(bbox.min_corner().y())))),
-			V2i(clamp(0, nc-1, int(ceil (bbox.max_corner().x()))),
-				clamp(0, nr-1, int(ceil (bbox.max_corner().y())))));
+			V2i(std::clamp(int(floor(bbox.min_corner().x())), 0, nc-1),
+				std::clamp(int(floor(bbox.min_corner().y())), 0, nr-1)),
+			V2i(std::clamp(int(ceil(bbox.max_corner().x())), 0, nc-1),
+				std::clamp(int(ceil(bbox.max_corner().y())), 0, nr-1)));
 		/*
 		// Fill in the spans, at which the circle intersects the rows.
 		int row_first = bboxi.min_corner().y();
@@ -758,7 +758,7 @@ void gcode_spread_points(
 				area_circle_total += area;
 				if (cell.area < area)
 					cell.area = area;
-				cell.fraction_covered = clamp(0.f, 1.f, (cell.area > 0) ? (area / cell.area) : 0);
+				cell.fraction_covered = std::clamp((cell.area > 0) ? (area / cell.area) : 0, 0.f, 1.f);
 				if (cell.fraction_covered == 0) {
 					-- n_cells;
 					continue;
@@ -1018,7 +1018,7 @@ void ExtrusionSimulator::evaluate_accumulator(ExtrusionSimulationType simulation
 			float p = mask[r][c];
 			#endif
 			int   idx = int(floor(p * float(pimpl->color_gradient.size()) + 0.5f));
-			V3uc  clr = pimpl->color_gradient[clamp(0, int(pimpl->color_gradient.size()-1), idx)];
+			V3uc  clr = pimpl->color_gradient[std::clamp(idx, 0, int(pimpl->color_gradient.size()-1))];
 			*ptr ++ = clr.get<0>();
 			*ptr ++ = clr.get<1>();
 			*ptr ++ = clr.get<2>();
diff --git a/src/libslic3r/Fill/Fill3DHoneycomb.cpp b/src/libslic3r/Fill/Fill3DHoneycomb.cpp
index 0dec8004b..e97a48680 100644
--- a/src/libslic3r/Fill/Fill3DHoneycomb.cpp
+++ b/src/libslic3r/Fill/Fill3DHoneycomb.cpp
@@ -55,8 +55,8 @@ static std::vector<coordf_t> perpendPoints(const coordf_t offset, const size_t b
 static inline void trim(Pointfs &pts, coordf_t minX, coordf_t minY, coordf_t maxX, coordf_t maxY)
 {
     for (Vec2d &pt : pts) {
-        pt(0) = clamp(minX, maxX, pt(0));
-        pt(1) = clamp(minY, maxY, pt(1));
+        pt.x() = std::clamp(pt.x(), minX, maxX);
+        pt.y() = std::clamp(pt.y(), minY, maxY);
     }
 }
 
diff --git a/src/libslic3r/Fill/FillGyroid.cpp b/src/libslic3r/Fill/FillGyroid.cpp
index c6510daa2..5c555df88 100644
--- a/src/libslic3r/Fill/FillGyroid.cpp
+++ b/src/libslic3r/Fill/FillGyroid.cpp
@@ -53,7 +53,7 @@ static inline Polyline make_wave(
     polyline.points.reserve(points.size());
     for (auto& point : points) {
         point(1) += offset;
-        point(1) = clamp(0., height, double(point(1)));
+        point(1) = std::clamp(double(point.y()), 0., height);
         if (vertical)
             std::swap(point(0), point(1));
         polyline.points.emplace_back((point * scaleFactor).cast<coord_t>());
diff --git a/src/libslic3r/Format/3mf.cpp b/src/libslic3r/Format/3mf.cpp
index b98cc5f04..1b59448fd 100644
--- a/src/libslic3r/Format/3mf.cpp
+++ b/src/libslic3r/Format/3mf.cpp
@@ -96,6 +96,7 @@ static constexpr const char* PRINTABLE_ATTR = "printable";
 static constexpr const char* INSTANCESCOUNT_ATTR = "instances_count";
 static constexpr const char* CUSTOM_SUPPORTS_ATTR = "slic3rpe:custom_supports";
 static constexpr const char* CUSTOM_SEAM_ATTR = "slic3rpe:custom_seam";
+static constexpr const char* MMU_SEGMENTATION_ATTR = "slic3rpe:mmu_segmentation";
 
 static constexpr const char* KEY_ATTR = "key";
 static constexpr const char* VALUE_ATTR = "value";
@@ -290,6 +291,7 @@ namespace Slic3r {
             std::vector<unsigned int> triangles;
             std::vector<std::string> custom_supports;
             std::vector<std::string> custom_seam;
+            std::vector<std::string> mmu_segmentation;
 
             bool empty() { return vertices.empty() || triangles.empty(); }
 
@@ -298,6 +300,7 @@ namespace Slic3r {
                 triangles.clear();
                 custom_supports.clear();
                 custom_seam.clear();
+                mmu_segmentation.clear();
             }
         };
 
@@ -1565,6 +1568,12 @@ namespace Slic3r {
 
         m_curr_object.geometry.custom_supports.push_back(get_attribute_value_string(attributes, num_attributes, CUSTOM_SUPPORTS_ATTR));
         m_curr_object.geometry.custom_seam.push_back(get_attribute_value_string(attributes, num_attributes, CUSTOM_SEAM_ATTR));
+//        m_curr_object.geometry.mmu_segmentation.push_back(get_attribute_value_string(attributes, num_attributes, MMU_SEGMENTATION_ATTR));
+        // FIXME Lukas H.: This is only for backward compatibility with older 3MF test files. Removes this when it is not necessary.
+        if(get_attribute_value_string(attributes, num_attributes, MMU_SEGMENTATION_ATTR) != "")
+            m_curr_object.geometry.mmu_segmentation.push_back(get_attribute_value_string(attributes, num_attributes, MMU_SEGMENTATION_ATTR));
+        else
+            m_curr_object.geometry.mmu_segmentation.push_back(get_attribute_value_string(attributes, num_attributes, "slic3rpe:mmu_painting"));
         return true;
     }
 
@@ -1889,15 +1898,18 @@ namespace Slic3r {
                 volume->source.transform = Slic3r::Geometry::Transformation(volume_matrix_to_object);
             volume->calculate_convex_hull();
 
-            // recreate custom supports and seam from previously loaded attribute
+            // recreate custom supports, seam and mmu segmentation from previously loaded attribute
             for (unsigned i=0; i<triangles_count; ++i) {
                 size_t index = src_start_id/3 + i;
                 assert(index < geometry.custom_supports.size());
                 assert(index < geometry.custom_seam.size());
+                assert(index < geometry.mmu_segmentation.size());
                 if (! geometry.custom_supports[index].empty())
                     volume->supported_facets.set_triangle_from_string(i, geometry.custom_supports[index]);
                 if (! geometry.custom_seam[index].empty())
                     volume->seam_facets.set_triangle_from_string(i, geometry.custom_seam[index]);
+                if (! geometry.mmu_segmentation[index].empty())
+                    volume->mmu_segmentation_facets.set_triangle_from_string(i, geometry.mmu_segmentation[index]);
             }
 
 
@@ -2532,6 +2544,15 @@ namespace Slic3r {
                     output_buffer += "\"";
                 }
 
+                std::string mmu_painting_data_string = volume->mmu_segmentation_facets.get_triangle_as_string(i);
+                if (! mmu_painting_data_string.empty()) {
+                    output_buffer += " ";
+                    output_buffer += MMU_SEGMENTATION_ATTR;
+                    output_buffer += "=\"";
+                    output_buffer += mmu_painting_data_string;
+                    output_buffer += "\"";
+                }
+
                 output_buffer += "/>\n";
 
                 if (! flush())
diff --git a/src/libslic3r/GCode/AvoidCrossingPerimeters.cpp b/src/libslic3r/GCode/AvoidCrossingPerimeters.cpp
index e00284fe1..49de854f2 100644
--- a/src/libslic3r/GCode/AvoidCrossingPerimeters.cpp
+++ b/src/libslic3r/GCode/AvoidCrossingPerimeters.cpp
@@ -620,7 +620,7 @@ static std::vector<float> contour_distance(const EdgeGrid::Grid     &grid,
                     const Vec2d  v        = (segment.second - segment.first).cast<double>();
                     const Vec2d  va       = (this->point - segment.first).cast<double>();
                     const double l2       = v.squaredNorm(); // avoid a sqrt
-                    const double t        = (l2 == 0.0) ? 0. : clamp(0., 1., va.dot(v) / l2);
+                    const double t        = (l2 == 0.0) ? 0. : std::clamp(va.dot(v) / l2, 0., 1.);
                     // Closest point from this->point to the segment.
                     const Vec2d  foot     = segment.first.cast<double>() + t * v;
                     const Vec2d  bisector = foot - this->point.cast<double>();
diff --git a/src/libslic3r/GCode/CoolingBuffer.cpp b/src/libslic3r/GCode/CoolingBuffer.cpp
index 6dd5e7e27..17f3e0b73 100644
--- a/src/libslic3r/GCode/CoolingBuffer.cpp
+++ b/src/libslic3r/GCode/CoolingBuffer.cpp
@@ -725,8 +725,8 @@ std::string CoolingBuffer::apply_layer_cooldown(
             if (int(layer_id) >= disable_fan_first_layers && int(layer_id) + 1 < full_fan_speed_layer) {
                 // Ramp up the fan speed from disable_fan_first_layers to full_fan_speed_layer.
                 float factor = float(int(layer_id + 1) - disable_fan_first_layers) / float(full_fan_speed_layer - disable_fan_first_layers);
-                fan_speed_new    = clamp(0, 255, int(float(fan_speed_new   ) * factor + 0.5f));
-                bridge_fan_speed = clamp(0, 255, int(float(bridge_fan_speed) * factor + 0.5f));
+                fan_speed_new    = std::clamp(int(float(fan_speed_new) * factor + 0.5f), 0, 255);
+                bridge_fan_speed = std::clamp(int(float(bridge_fan_speed) * factor + 0.5f), 0, 255);
             }
 #undef EXTRUDER_CONFIG
             bridge_fan_control = bridge_fan_speed > fan_speed_new;
diff --git a/src/libslic3r/Model.cpp b/src/libslic3r/Model.cpp
index 82d23021f..ce359e7bf 100644
--- a/src/libslic3r/Model.cpp
+++ b/src/libslic3r/Model.cpp
@@ -1064,6 +1064,7 @@ void ModelObject::convert_units(ModelObjectPtrs& new_objects, ConversionType con
 
             vol->supported_facets.assign(volume->supported_facets);
             vol->seam_facets.assign(volume->seam_facets);
+            vol->mmu_segmentation_facets.assign(volume->mmu_segmentation_facets);
 
             // Perform conversion only if the target "imperial" state is different from the current one.
             // This check supports conversion of "mixed" set of volumes, each with different "imperial" state.
@@ -1165,6 +1166,7 @@ ModelObjectPtrs ModelObject::cut(size_t instance, coordf_t z, bool keep_upper, b
 
         volume->supported_facets.clear();
         volume->seam_facets.clear();
+        volume->mmu_segmentation_facets.clear();
 
         if (! volume->is_model_part()) {
             // Modifiers are not cut, but we still need to add the instance transformation
@@ -1758,6 +1760,7 @@ void ModelVolume::assign_new_unique_ids_recursive()
     config.set_new_unique_id();
     supported_facets.set_new_unique_id();
     seam_facets.set_new_unique_id();
+    mmu_segmentation_facets.set_new_unique_id();
 }
 
 void ModelVolume::rotate(double angle, Axis axis)
@@ -2074,8 +2077,10 @@ bool model_volume_list_changed(const ModelObject &model_object_old, const ModelO
 
 bool model_custom_supports_data_changed(const ModelObject& mo, const ModelObject& mo_new) {
     assert(! model_volume_list_changed(mo, mo_new, ModelVolumeType::MODEL_PART));
-    assert(mo.volumes.size() == mo_new.volumes.size());
-    for (size_t i=0; i<mo.volumes.size(); ++i) {
+    // FIXME Lukas H.: Because of adding another mesh modifiers when slicing, then assert triggered and possible crash. It requires changing the integration of MMU segmentation.
+//    assert(mo.volumes.size() == mo_new.volumes.size());
+//    for (size_t i=0; i<mo.volumes.size(); ++i) {
+    for (size_t i=0; i<std::min(mo.volumes.size(), mo_new.volumes.size()); ++i) {
         if (! mo_new.volumes[i]->supported_facets.timestamp_matches(mo.volumes[i]->supported_facets))
             return true;
     }
@@ -2084,14 +2089,28 @@ bool model_custom_supports_data_changed(const ModelObject& mo, const ModelObject
 
 bool model_custom_seam_data_changed(const ModelObject& mo, const ModelObject& mo_new) {
     assert(! model_volume_list_changed(mo, mo_new, ModelVolumeType::MODEL_PART));
-    assert(mo.volumes.size() == mo_new.volumes.size());
-    for (size_t i=0; i<mo.volumes.size(); ++i) {
+    // FIXME Lukas H.: Because of adding another mesh modifiers when slicing, then assert triggered and possible crash. It requires changing the integration of MMU segmentation.
+//    assert(mo.volumes.size() == mo_new.volumes.size());
+//    for (size_t i=0; i<mo.volumes.size(); ++i) {
+    for (size_t i=0; i<std::min(mo.volumes.size(), mo_new.volumes.size()); ++i) {
         if (! mo_new.volumes[i]->seam_facets.timestamp_matches(mo.volumes[i]->seam_facets))
             return true;
     }
     return false;
 }
 
+bool model_mmu_segmentation_data_changed(const ModelObject& mo, const ModelObject& mo_new) {
+    assert(! model_volume_list_changed(mo, mo_new, ModelVolumeType::MODEL_PART));
+    // FIXME Lukas H.: Because of adding another mesh modifiers when slicing, then assert triggered and possible crash. It requires changing the integration of MMU segmentation.
+//    assert(mo.volumes.size() == mo_new.volumes.size());
+//    for (size_t i=0; i<mo.volumes.size(); ++i) {
+    for (size_t i=0; i<std::min(mo.volumes.size(), mo_new.volumes.size()); ++i) {
+        if (! mo_new.volumes[i]->mmu_segmentation_facets.timestamp_matches(mo.volumes[i]->mmu_segmentation_facets))
+            return true;
+    }
+    return false;
+}
+
 extern bool model_has_multi_part_objects(const Model &model)
 {
     for (const ModelObject *model_object : model.objects)
diff --git a/src/libslic3r/Model.hpp b/src/libslic3r/Model.hpp
index ea0043749..ed8030c9b 100644
--- a/src/libslic3r/Model.hpp
+++ b/src/libslic3r/Model.hpp
@@ -591,6 +591,9 @@ public:
     // List of seam enforcers/blockers.
     FacetsAnnotation    seam_facets;
 
+    // List of mesh facets painted for MMU segmentation.
+    FacetsAnnotation    mmu_segmentation_facets;
+
     // A parent object owning this modifier volume.
     ModelObject*        get_object() const { return this->object; }
     ModelVolumeType     type() const { return m_type; }
@@ -679,6 +682,7 @@ public:
         this->config.set_new_unique_id();
         this->supported_facets.set_new_unique_id();
         this->seam_facets.set_new_unique_id();
+        this->mmu_segmentation_facets.set_new_unique_id();
     }
 
 protected:
@@ -718,22 +722,26 @@ private:
 		assert(this->id().valid()); 
         assert(this->config.id().valid()); 
         assert(this->supported_facets.id().valid()); 
-        assert(this->seam_facets.id().valid()); 
+        assert(this->seam_facets.id().valid());
+        assert(this->mmu_segmentation_facets.id().valid());
         assert(this->id() != this->config.id());
         assert(this->id() != this->supported_facets.id());
         assert(this->id() != this->seam_facets.id());
+        assert(this->id() != this->mmu_segmentation_facets.id());
         if (mesh.stl.stats.number_of_facets > 1)
             calculate_convex_hull();
     }
     ModelVolume(ModelObject *object, TriangleMesh &&mesh, TriangleMesh &&convex_hull) :
 		m_mesh(new TriangleMesh(std::move(mesh))), m_convex_hull(new TriangleMesh(std::move(convex_hull))), m_type(ModelVolumeType::MODEL_PART), object(object) {
 		assert(this->id().valid()); 
-        assert(this->config.id().valid()); 
-        assert(this->supported_facets.id().valid()); 
-        assert(this->seam_facets.id().valid()); 
+        assert(this->config.id().valid());
+        assert(this->supported_facets.id().valid());
+        assert(this->seam_facets.id().valid());
+        assert(this->mmu_segmentation_facets.id().valid());
         assert(this->id() != this->config.id());
         assert(this->id() != this->supported_facets.id());
         assert(this->id() != this->seam_facets.id());
+        assert(this->id() != this->mmu_segmentation_facets.id());
 	}
 
     // Copying an existing volume, therefore this volume will get a copy of the ID assigned.
@@ -741,19 +749,22 @@ private:
         ObjectBase(other),
         name(other.name), source(other.source), m_mesh(other.m_mesh), m_convex_hull(other.m_convex_hull),
         config(other.config), m_type(other.m_type), object(object), m_transformation(other.m_transformation),
-        supported_facets(other.supported_facets), seam_facets(other.seam_facets)
+        supported_facets(other.supported_facets), seam_facets(other.seam_facets), mmu_segmentation_facets(other.mmu_segmentation_facets)
     {
 		assert(this->id().valid()); 
         assert(this->config.id().valid()); 
         assert(this->supported_facets.id().valid());
         assert(this->seam_facets.id().valid());
+        assert(this->mmu_segmentation_facets.id().valid());
         assert(this->id() != this->config.id());
         assert(this->id() != this->supported_facets.id());
         assert(this->id() != this->seam_facets.id());
+        assert(this->id() != this->mmu_segmentation_facets.id());
 		assert(this->id() == other.id());
         assert(this->config.id() == other.config.id());
         assert(this->supported_facets.id() == other.supported_facets.id());
         assert(this->seam_facets.id() == other.seam_facets.id());
+        assert(this->mmu_segmentation_facets.id() == other.mmu_segmentation_facets.id());
         this->set_material_id(other.material_id());
     }
     // Providing a new mesh, therefore this volume will get a new unique ID assigned.
@@ -764,9 +775,11 @@ private:
         assert(this->config.id().valid()); 
         assert(this->supported_facets.id().valid());
         assert(this->seam_facets.id().valid());
+        assert(this->mmu_segmentation_facets.id().valid());
         assert(this->id() != this->config.id());
         assert(this->id() != this->supported_facets.id());
         assert(this->id() != this->seam_facets.id());
+        assert(this->id() != this->mmu_segmentation_facets.id());
 		assert(this->id() != other.id());
         assert(this->config.id() == other.config.id());
         this->set_material_id(other.material_id());
@@ -777,9 +790,11 @@ private:
         assert(this->config.id() != other.config.id()); 
         assert(this->supported_facets.id() != other.supported_facets.id());
         assert(this->seam_facets.id() != other.seam_facets.id());
+        assert(this->mmu_segmentation_facets.id() != other.mmu_segmentation_facets.id());
         assert(this->id() != this->config.id());
         assert(this->supported_facets.empty());
         assert(this->seam_facets.empty());
+        assert(this->mmu_segmentation_facets.empty());
     }
 
     ModelVolume& operator=(ModelVolume &rhs) = delete;
@@ -787,17 +802,19 @@ private:
 	friend class cereal::access;
 	friend class UndoRedo::StackImpl;
 	// Used for deserialization, therefore no IDs are allocated.
-	ModelVolume() : ObjectBase(-1), config(-1), supported_facets(-1), seam_facets(-1), object(nullptr) {
+	ModelVolume() : ObjectBase(-1), config(-1), supported_facets(-1), seam_facets(-1), mmu_segmentation_facets(-1), object(nullptr) {
 		assert(this->id().invalid());
         assert(this->config.id().invalid());
         assert(this->supported_facets.id().invalid());
         assert(this->seam_facets.id().invalid());
+        assert(this->mmu_segmentation_facets.id().invalid());
 	}
 	template<class Archive> void load(Archive &ar) {
 		bool has_convex_hull;
         ar(name, source, m_mesh, m_type, m_material_id, m_transformation, m_is_splittable, has_convex_hull);
         cereal::load_by_value(ar, supported_facets);
         cereal::load_by_value(ar, seam_facets);
+        cereal::load_by_value(ar, mmu_segmentation_facets);
         cereal::load_by_value(ar, config);
 		assert(m_mesh);
 		if (has_convex_hull) {
@@ -813,6 +830,7 @@ private:
         ar(name, source, m_mesh, m_type, m_material_id, m_transformation, m_is_splittable, has_convex_hull);
         cereal::save_by_value(ar, supported_facets);
         cereal::save_by_value(ar, seam_facets);
+        cereal::save_by_value(ar, mmu_segmentation_facets);
         cereal::save_by_value(ar, config);
 		if (has_convex_hull)
 			cereal::save_optional(ar, m_convex_hull);
@@ -1094,6 +1112,10 @@ bool model_custom_supports_data_changed(const ModelObject& mo, const ModelObject
 // The function assumes that volumes list is synchronized.
 bool model_custom_seam_data_changed(const ModelObject& mo, const ModelObject& mo_new);
 
+// Test whether the now ModelObject has newer MMU segmentation data than the old one.
+// The function assumes that volumes list is synchronized.
+extern bool model_mmu_segmentation_data_changed(const ModelObject& mo, const ModelObject& mo_new);
+
 // If the model has multi-part objects, then it is currently not supported by the SLA mode.
 // Either the model cannot be loaded, or a SLA printer has to be activated.
 bool model_has_multi_part_objects(const Model &model);
diff --git a/src/libslic3r/MultiMaterialSegmentation.cpp b/src/libslic3r/MultiMaterialSegmentation.cpp
new file mode 100644
index 000000000..778a64018
--- /dev/null
+++ b/src/libslic3r/MultiMaterialSegmentation.cpp
@@ -0,0 +1,1523 @@
+#include "BoundingBox.hpp"
+#include "ClipperUtils.hpp"
+#include "EdgeGrid.hpp"
+#include "Layer.hpp"
+#include "Print.hpp"
+#include "VoronoiVisualUtils.hpp"
+
+#include <utility>
+#include <cfloat>
+#include <unordered_set>
+
+#include <boost/log/trivial.hpp>
+
+#include <tbb/parallel_for.h>
+
+#include <boost/geometry.hpp>
+#include <boost/geometry/geometries/point.hpp>
+#include <boost/geometry/geometries/segment.hpp>
+#include <boost/geometry/index/rtree.hpp>
+
+namespace Slic3r {
+struct ColoredLine {
+    Line line;
+    int color;
+    int poly_idx = -1;
+    int local_line_idx = -1;
+};
+}
+
+#include <boost/polygon/polygon.hpp>
+namespace boost::polygon {
+template <>
+struct geometry_concept<Slic3r::ColoredLine> { typedef segment_concept type; };
+
+template <>
+struct segment_traits<Slic3r::ColoredLine> {
+    typedef coord_t coordinate_type;
+    typedef Slic3r::Point point_type;
+
+    static inline point_type get(const Slic3r::ColoredLine& line, const direction_1d& dir) {
+        return dir.to_int() ? line.line.b : line.line.a;
+    }
+};
+}
+
+namespace Slic3r {
+
+// Assumes that is at most same projected_l length or below than projection_l
+static bool project_line_on_line(const Line &projection_l, const Line &projected_l, Line *new_projected)
+{
+    const Vec2d  v1 = (projection_l.b - projection_l.a).cast<double>();
+    const Vec2d  va = (projected_l.a - projection_l.a).cast<double>();
+    const Vec2d  vb = (projected_l.b - projection_l.a).cast<double>();
+    const double l2 = v1.squaredNorm(); // avoid a sqrt
+    if (l2 == 0.0)
+        return false;
+    double t1 = va.dot(v1) / l2;
+    double t2 = vb.dot(v1) / l2;
+    t1        = std::clamp(t1, 0., 1.);
+    t2        = std::clamp(t2, 0., 1.);
+    assert(t1 >= 0.);
+    assert(t2 >= 0.);
+    assert(t1 <= 1.);
+    assert(t2 <= 1.);
+
+    Point p1       = projection_l.a + (t1 * v1).cast<coord_t>();
+    Point p2       = projection_l.a + (t2 * v1).cast<coord_t>();
+    *new_projected = Line(p1, p2);
+    return true;
+}
+
+struct PaintedLine
+{
+    size_t contour_idx;
+    size_t line_idx;
+    Line   projected_line;
+    int    color;
+};
+
+struct PaintedLineVisitor
+{
+    PaintedLineVisitor(const EdgeGrid::Grid &grid, std::vector<PaintedLine> &painted_lines) : grid(grid), painted_lines(painted_lines)
+    {
+        painted_lines_set.reserve(painted_lines.capacity());
+    }
+
+    void reset() { painted_lines_set.clear(); }
+
+    bool operator()(coord_t iy, coord_t ix)
+    {
+        // Called with a row and column of the grid cell, which is intersected by a line.
+        auto         cell_data_range = grid.cell_data_range(iy, ix);
+        const Vec2d  v1              = line_to_test.vector().cast<double>();
+        for (auto it_contour_and_segment = cell_data_range.first; it_contour_and_segment != cell_data_range.second; ++it_contour_and_segment) {
+            Line        grid_line = grid.line(*it_contour_and_segment);
+            const Vec2d v2        = grid_line.vector().cast<double>();
+            // When lines have too different length, it is necessary to normalize them
+            if (Slic3r::sqr(v1.dot(v2)) > cos_threshold2 * v1.squaredNorm() * v2.squaredNorm()) {
+                // The two vectors are nearly collinear (their mutual angle is lower than 30 degrees)
+                if (painted_lines_set.find(*it_contour_and_segment) == painted_lines_set.end()) {
+                    double dist_1     = grid_line.distance_to(line_to_test.a);
+                    double dist_2     = grid_line.distance_to(line_to_test.b);
+                    double dist_3     = line_to_test.distance_to(grid_line.a);
+                    double dist_4     = line_to_test.distance_to(grid_line.b);
+                    double total_dist = std::min(std::min(dist_1, dist_2), std::min(dist_3, dist_4));
+
+                    if (total_dist < 50 * SCALED_EPSILON) {
+                        Line line_to_test_projected;
+                        project_line_on_line(grid_line, line_to_test, &line_to_test_projected);
+
+                        if (Line(grid_line.a, line_to_test_projected.a).length() > Line(grid_line.a, line_to_test_projected.b).length()) {
+                            line_to_test_projected.reverse();
+                        }
+                        painted_lines.push_back({it_contour_and_segment->first, it_contour_and_segment->second, line_to_test_projected, this->color});
+                        painted_lines_set.insert(*it_contour_and_segment);
+                    }
+                }
+            }
+        }
+        // Continue traversing the grid along the edge.
+        return true;
+    }
+
+    const EdgeGrid::Grid                                                                 &grid;
+    std::vector<PaintedLine>                                                             &painted_lines;
+    Line                                                                                  line_to_test;
+    std::unordered_set<std::pair<size_t, size_t>, boost::hash<std::pair<size_t, size_t>>> painted_lines_set;
+    int                                                                                   color = -1;
+
+    static inline const double                                                            cos_threshold2 = Slic3r::sqr(cos(M_PI * 30. / 180.));
+};
+
+static std::vector<ColoredLine> to_colored_lines(const Polygon &polygon, int color)
+{
+    std::vector<ColoredLine> lines;
+    if (polygon.points.size() > 2) {
+        lines.reserve(polygon.points.size());
+        for (auto it = polygon.points.begin(); it != polygon.points.end() - 1; ++it)
+            lines.push_back({Line(*it, *(it + 1)), color});
+        lines.push_back({Line(polygon.points.back(), polygon.points.front()), color});
+    }
+    return lines;
+}
+
+static Polygon colored_points_to_polygon(const std::vector<ColoredLine> &lines)
+{
+    Polygon out;
+    out.points.reserve(lines.size());
+    for (const ColoredLine &l : lines)
+        out.points.emplace_back(l.line.a);
+    return out;
+}
+
+static Polygons colored_points_to_polygon(const std::vector<std::vector<ColoredLine>> &lines)
+{
+    Polygons out;
+    for (const std::vector<ColoredLine> &l : lines)
+        out.emplace_back(colored_points_to_polygon(l));
+    return out;
+}
+
+// Flatten the vector of vectors into a vector.
+static inline std::vector<ColoredLine> to_lines(const std::vector<std::vector<ColoredLine>> &c_lines)
+{
+    size_t n_lines = 0;
+    for (const auto &c_line : c_lines)
+        n_lines += c_line.size();
+    std::vector<ColoredLine> lines;
+    lines.reserve(n_lines);
+    for (const auto &c_line : c_lines)
+        lines.insert(lines.end(), c_line.begin(), c_line.end());
+    return lines;
+}
+
+// Double vertex equal to a coord_t point after conversion to double.
+static bool vertex_equal_to_point(const Voronoi::VD::vertex_type &vertex, const Point &ipt)
+{
+    // Convert ipt to doubles, force the 80bit FPU temporary to 64bit and then compare.
+    // This should work with any settings of math compiler switches and the C++ compiler
+    // shall understand the memcpies as type punning and it shall optimize them out.
+    using ulp_cmp_type = boost::polygon::detail::ulp_comparison<double>;
+    ulp_cmp_type ulp_cmp;
+    static constexpr int ULPS = boost::polygon::voronoi_diagram_traits<double>::vertex_equality_predicate_type::ULPS;
+    return ulp_cmp(vertex.x(), double(ipt.x()), ULPS) == ulp_cmp_type::EQUAL &&
+           ulp_cmp(vertex.y(), double(ipt.y()), ULPS) == ulp_cmp_type::EQUAL;
+}
+
+static inline bool vertex_equal_to_point(const Voronoi::VD::vertex_type *vertex, const Point &ipt) {
+    return vertex_equal_to_point(*vertex, ipt);
+}
+
+static std::vector<std::pair<size_t, size_t>> get_segments(const std::vector<ColoredLine> &polygon)
+{
+    std::vector<std::pair<size_t, size_t>> segments;
+
+    size_t segment_end = 0;
+    while (segment_end + 1 < polygon.size() && polygon[segment_end].color == polygon[segment_end + 1].color)
+        segment_end++;
+
+    if (segment_end == polygon.size() - 1)
+        return {std::make_pair(0, polygon.size() - 1)};
+
+    size_t first_different_color = (segment_end + 1) % polygon.size();
+    for (size_t line_offset_idx = 0; line_offset_idx < polygon.size(); ++line_offset_idx) {
+        size_t start_s = (first_different_color + line_offset_idx) % polygon.size();
+        size_t end_s   = start_s;
+
+        while (line_offset_idx + 1 < polygon.size() && polygon[start_s].color == polygon[(first_different_color + line_offset_idx + 1) % polygon.size()].color) {
+            end_s = (first_different_color + line_offset_idx + 1) % polygon.size();
+            line_offset_idx++;
+        }
+        segments.emplace_back(start_s, end_s);
+    }
+    return segments;
+}
+
+static std::vector<std::vector<std::pair<size_t, size_t>>> get_all_segments(const std::vector<std::vector<ColoredLine>> &color_poly)
+{
+    std::vector<std::vector<std::pair<size_t, size_t>>> all_segments(color_poly.size());
+    for (size_t poly_idx = 0; poly_idx < color_poly.size(); ++poly_idx) {
+        const std::vector<ColoredLine> &c_polygon = color_poly[poly_idx];
+        all_segments[poly_idx]                    = get_segments(c_polygon);
+    }
+    return all_segments;
+}
+
+static std::vector<ColoredLine> colorize_line(const Line &              line_to_process,
+                                              const size_t              start_idx,
+                                              const size_t              end_idx,
+                                              std::vector<PaintedLine> &painted_lines)
+{
+    std::vector<PaintedLine> internal_painted;
+    for (size_t line_idx = start_idx; line_idx <= end_idx; ++line_idx) { internal_painted.emplace_back(painted_lines[line_idx]); }
+    const int                filter_eps_value = scale_(0.1f);
+    std::vector<PaintedLine> filtered_lines;
+    filtered_lines.emplace_back(internal_painted.front());
+    for (size_t line_idx = 1; line_idx < internal_painted.size(); ++line_idx) {
+        PaintedLine &prev = filtered_lines.back();
+        PaintedLine &curr = internal_painted[line_idx];
+
+        double prev_length        = prev.projected_line.length();
+        double curr_dist_start    = (curr.projected_line.a - prev.projected_line.a).cast<double>().norm();
+        double dist_between_lines = curr_dist_start - prev_length;
+
+        if (dist_between_lines >= 0) {
+            if (prev.color == curr.color) {
+                if (dist_between_lines <= filter_eps_value) {
+                    prev.projected_line.b = curr.projected_line.b;
+                } else {
+                    filtered_lines.emplace_back(curr);
+                }
+            } else {
+                filtered_lines.emplace_back(curr);
+            }
+        } else {
+            double curr_dist_end = (curr.projected_line.b - prev.projected_line.a).cast<double>().norm();
+            if (curr_dist_end <= prev_length) {
+            } else {
+                if (prev.color == curr.color) {
+                    prev.projected_line.b = curr.projected_line.b;
+                } else {
+                    curr.projected_line.a = prev.projected_line.b;
+                    filtered_lines.emplace_back(curr);
+                }
+            }
+        }
+    }
+
+    std::vector<ColoredLine> final_lines;
+    double                   dist_to_start = (filtered_lines.front().projected_line.a - line_to_process.a).cast<double>().norm();
+    if (dist_to_start <= filter_eps_value) {
+        filtered_lines.front().projected_line.a = line_to_process.a;
+        final_lines.push_back({filtered_lines.front().projected_line, filtered_lines.front().color});
+    } else {
+        final_lines.push_back({Line(line_to_process.a, filtered_lines.front().projected_line.a), 0});
+        final_lines.push_back({filtered_lines.front().projected_line, filtered_lines.front().color});
+    }
+
+    for (size_t line_idx = 1; line_idx < filtered_lines.size(); ++line_idx) {
+        ColoredLine &prev = final_lines.back();
+        PaintedLine &curr = filtered_lines[line_idx];
+
+        double line_dist = (curr.projected_line.a - prev.line.b).cast<double>().norm();
+        if (line_dist <= filter_eps_value) {
+            if (prev.color == curr.color) {
+                prev.line.b = curr.projected_line.b;
+            } else {
+                prev.line.b = curr.projected_line.a;
+                final_lines.push_back({curr.projected_line, curr.color});
+            }
+        } else {
+            final_lines.push_back({Line(prev.line.b, curr.projected_line.a), 0});
+            final_lines.push_back({curr.projected_line, curr.color});
+        }
+    }
+
+    double dist_to_end = (final_lines.back().line.b - line_to_process.b).cast<double>().norm();
+    if (dist_to_end <= filter_eps_value)
+        final_lines.back().line.b = line_to_process.b;
+    else
+        final_lines.push_back({Line(final_lines.back().line.b, line_to_process.b), 0});
+
+    for (size_t line_idx = 1; line_idx < final_lines.size(); ++line_idx)
+        assert(final_lines[line_idx - 1].line.b == final_lines[line_idx].line.a);
+
+    for (size_t line_idx = 2; line_idx < final_lines.size(); ++line_idx) {
+        const ColoredLine &line_0 = final_lines[line_idx - 2];
+        ColoredLine &      line_1 = final_lines[line_idx - 1];
+        const ColoredLine &line_2 = final_lines[line_idx - 0];
+
+        if (line_0.color == line_2.color && line_0.color != line_1.color)
+            if (line_1.line.length() <= scale_(0.2)) line_1.color = line_0.color;
+    }
+
+    std::vector<ColoredLine> colored_lines_simpl;
+    colored_lines_simpl.emplace_back(final_lines.front());
+    for (size_t line_idx = 1; line_idx < final_lines.size(); ++line_idx) {
+        const ColoredLine &line_0 = final_lines[line_idx];
+
+        if (colored_lines_simpl.back().color == line_0.color)
+            colored_lines_simpl.back().line.b = line_0.line.b;
+        else
+            colored_lines_simpl.emplace_back(line_0);
+    }
+
+    final_lines = colored_lines_simpl;
+
+    if (final_lines.size() > 1) {
+        if (final_lines.front().color != final_lines[1].color && final_lines.front().line.length() <= scale_(0.2)) {
+            final_lines[1].line.a = final_lines.front().line.a;
+            final_lines.erase(final_lines.begin());
+        }
+    }
+
+    if (final_lines.size() > 1) {
+        if (final_lines.back().color != final_lines[final_lines.size() - 2].color && final_lines.back().line.length() <= scale_(0.2)) {
+            final_lines[final_lines.size() - 2].line.b = final_lines.back().line.b;
+            final_lines.pop_back();
+        }
+    }
+
+    return final_lines;
+}
+
+static std::vector<ColoredLine> colorize_polygon(const Polygon &poly, const size_t start_idx, const size_t end_idx, std::vector<PaintedLine> &painted_lines)
+{
+    std::vector<ColoredLine> new_lines;
+    Lines                    lines = poly.lines();
+
+    for (size_t idx = 0; idx < painted_lines[start_idx].line_idx; ++idx)
+        new_lines.emplace_back(ColoredLine{lines[idx], 0});
+
+    for (size_t first_idx = start_idx; first_idx <= end_idx; ++first_idx) {
+        size_t second_idx = first_idx;
+        while (second_idx <= end_idx && painted_lines[first_idx].line_idx == painted_lines[second_idx].line_idx) ++second_idx;
+        --second_idx;
+
+        assert(painted_lines[first_idx].line_idx == painted_lines[second_idx].line_idx);
+        std::vector<ColoredLine> lines_c_line = colorize_line(lines[painted_lines[first_idx].line_idx], first_idx, second_idx, painted_lines);
+        new_lines.insert(new_lines.end(), lines_c_line.begin(), lines_c_line.end());
+
+        if (second_idx + 1 <= end_idx)
+            for (size_t idx = painted_lines[second_idx].line_idx + 1; idx < painted_lines[second_idx + 1].line_idx; ++idx)
+                new_lines.emplace_back(ColoredLine{lines[idx], 0});
+
+        first_idx = second_idx;
+    }
+
+    for (size_t idx = painted_lines[end_idx].line_idx + 1; idx < poly.size(); ++idx)
+        new_lines.emplace_back(ColoredLine{lines[idx], 0});
+
+    for (size_t line_idx = 2; line_idx < new_lines.size(); ++line_idx) {
+        const ColoredLine &line_0 = new_lines[line_idx - 2];
+        ColoredLine &      line_1 = new_lines[line_idx - 1];
+        const ColoredLine &line_2 = new_lines[line_idx - 0];
+
+        if (line_0.color == line_2.color && line_0.color != line_1.color && line_0.color >= 1) {
+            if (line_1.line.length() <= scale_(0.5)) line_1.color = line_0.color;
+        }
+    }
+
+    for (size_t line_idx = 3; line_idx < new_lines.size(); ++line_idx) {
+        const ColoredLine &line_0 = new_lines[line_idx - 3];
+        ColoredLine &      line_1 = new_lines[line_idx - 2];
+        ColoredLine &      line_2 = new_lines[line_idx - 1];
+        const ColoredLine &line_3 = new_lines[line_idx - 0];
+
+        if (line_0.color == line_3.color && (line_0.color != line_1.color || line_0.color != line_2.color) && line_0.color >= 1 && line_3.color >= 1) {
+            if ((line_1.line.length() + line_2.line.length()) <= scale_(0.5)) {
+                line_1.color = line_0.color;
+                line_2.color = line_0.color;
+            }
+        }
+    }
+
+    std::vector<std::pair<size_t, size_t>> segments       = get_segments(new_lines);
+    auto                                   segment_length = [&new_lines](const std::pair<size_t, size_t> &segment) {
+        double total_length = 0;
+        for (size_t seg_start_idx = segment.first; seg_start_idx != segment.second; seg_start_idx = (seg_start_idx + 1 < new_lines.size()) ? seg_start_idx + 1 : 0)
+            total_length += new_lines[seg_start_idx].line.length();
+        total_length += new_lines[segment.second].line.length();
+        return total_length;
+    };
+
+    for (size_t pair_idx = 1; pair_idx < segments.size(); ++pair_idx) {
+        int color0 = new_lines[segments[pair_idx - 1].first].color;
+        int color1 = new_lines[segments[pair_idx - 0].first].color;
+
+        double seg0l = segment_length(segments[pair_idx - 1]);
+        double seg1l = segment_length(segments[pair_idx - 0]);
+
+        if (color0 != color1 && seg0l >= scale_(0.1) && seg1l <= scale_(0.2)) {
+            for (size_t seg_start_idx = segments[pair_idx].first; seg_start_idx != segments[pair_idx].second; seg_start_idx = (seg_start_idx + 1 < new_lines.size()) ? seg_start_idx + 1 : 0)
+                new_lines[seg_start_idx].color = color0;
+            new_lines[segments[pair_idx].second].color = color0;
+        }
+    }
+
+    segments = get_segments(new_lines);
+    for (size_t pair_idx = 1; pair_idx < segments.size(); ++pair_idx) {
+        int    color0 = new_lines[segments[pair_idx - 1].first].color;
+        int    color1 = new_lines[segments[pair_idx - 0].first].color;
+        double seg1l  = segment_length(segments[pair_idx - 0]);
+
+        if (color0 >= 1 && color0 != color1 && seg1l <= scale_(0.2)) {
+            for (size_t seg_start_idx = segments[pair_idx].first; seg_start_idx != segments[pair_idx].second; seg_start_idx = (seg_start_idx + 1 < new_lines.size()) ? seg_start_idx + 1 : 0)
+                new_lines[seg_start_idx].color = color0;
+            new_lines[segments[pair_idx].second].color = color0;
+        }
+    }
+
+    for (size_t pair_idx = 2; pair_idx < segments.size(); ++pair_idx) {
+        int color0 = new_lines[segments[pair_idx - 2].first].color;
+        int color1 = new_lines[segments[pair_idx - 1].first].color;
+        int color2 = new_lines[segments[pair_idx - 0].first].color;
+
+        if (color0 > 0 && color0 == color2 && color0 != color1 && segment_length(segments[pair_idx - 1]) <= scale_(0.5)) {
+            for (size_t seg_start_idx = segments[pair_idx].first; seg_start_idx != segments[pair_idx].second; seg_start_idx = (seg_start_idx + 1 < new_lines.size()) ? seg_start_idx + 1 : 0)
+                new_lines[seg_start_idx].color = color0;
+            new_lines[segments[pair_idx].second].color = color0;
+        }
+    }
+
+    return new_lines;
+}
+
+static std::vector<std::vector<ColoredLine>> colorize_polygons(const Polygons &polygons, std::vector<PaintedLine> &painted_lines)
+{
+    const size_t start_idx = 0;
+    const size_t end_idx   = painted_lines.size() - 1;
+
+    std::vector<std::vector<ColoredLine>> new_polygons;
+
+    for (size_t idx = 0; idx < painted_lines[start_idx].contour_idx; ++idx)
+        new_polygons.emplace_back(to_colored_lines(polygons[idx], 0));
+
+    for (size_t first_idx = start_idx; first_idx <= end_idx; ++first_idx) {
+        size_t second_idx = first_idx;
+        while (second_idx <= end_idx && painted_lines[first_idx].contour_idx == painted_lines[second_idx].contour_idx)
+            ++second_idx;
+        --second_idx;
+
+        assert(painted_lines[first_idx].contour_idx == painted_lines[second_idx].contour_idx);
+        std::vector<ColoredLine> polygon_c = colorize_polygon(polygons[painted_lines[first_idx].contour_idx], first_idx, second_idx, painted_lines);
+        new_polygons.emplace_back(polygon_c);
+
+        if (second_idx + 1 <= end_idx)
+            for (size_t idx = painted_lines[second_idx].contour_idx + 1; idx < painted_lines[second_idx + 1].contour_idx; ++idx)
+                new_polygons.emplace_back(to_colored_lines(polygons[idx], 0));
+
+        first_idx = second_idx;
+    }
+
+    for (size_t idx = painted_lines[end_idx].contour_idx + 1; idx < polygons.size(); ++idx)
+        new_polygons.emplace_back(to_colored_lines(polygons[idx], 0));
+
+    return new_polygons;
+}
+
+using boost::polygon::voronoi_diagram;
+
+struct MMU_Graph
+{
+    enum class ARC_TYPE { BORDER, NON_BORDER };
+
+    struct Arc
+    {
+        size_t   from_idx;
+        size_t   to_idx;
+        int      color;
+        ARC_TYPE type;
+        bool     used{false};
+
+        bool operator==(const Arc &rhs) const { return (from_idx == rhs.from_idx) && (to_idx == rhs.to_idx) && (color == rhs.color) && (type == rhs.type); }
+        bool operator!=(const Arc &rhs) const { return !operator==(rhs); }
+    };
+
+    struct Node
+    {
+        Point                     point;
+        std::list<MMU_Graph::Arc> neighbours;
+
+        void remove_edge(const size_t to_idx)
+        {
+            for (auto arc_it = this->neighbours.begin(); arc_it != this->neighbours.end(); ++arc_it) {
+                if (arc_it->to_idx == to_idx) {
+                    assert(arc_it->type != ARC_TYPE::BORDER);
+                    this->neighbours.erase(arc_it);
+                    break;
+                }
+            }
+        }
+    };
+
+    std::vector<MMU_Graph::Node> nodes;
+    std::vector<MMU_Graph::Arc>  arcs;
+    size_t                       all_border_points{};
+
+    std::vector<size_t> polygon_idx_offset;
+    std::vector<size_t> polygon_sizes;
+
+    void remove_edge(const size_t from_idx, const size_t to_idx)
+    {
+        nodes[from_idx].remove_edge(to_idx);
+        nodes[to_idx].remove_edge(from_idx);
+    }
+
+    size_t get_global_index(const size_t poly_idx, const size_t point_idx) const { return polygon_idx_offset[poly_idx] + point_idx; }
+
+    void append_edge(const size_t &from_idx, const size_t &to_idx, int color = -1, ARC_TYPE type = ARC_TYPE::NON_BORDER)
+    {
+        // Don't append duplicate edges between the same nodes.
+        for (const MMU_Graph::Arc &arc : this->nodes[from_idx].neighbours)
+            if (arc.to_idx == to_idx)
+                return;
+        for (const MMU_Graph::Arc &arc : this->nodes[to_idx].neighbours)
+            if (arc.to_idx == to_idx)
+                return;
+
+        this->nodes[from_idx].neighbours.push_back({from_idx, to_idx, color, type});
+        this->nodes[to_idx].neighbours.push_back({to_idx, from_idx, color, type});
+        this->arcs.push_back({from_idx, to_idx, color, type});
+        this->arcs.push_back({to_idx, from_idx, color, type});
+    }
+
+    // Ignoring arcs in the opposite direction
+    MMU_Graph::Arc get_arc(size_t idx) { return this->arcs[idx * 2]; }
+
+    size_t nodes_count() const { return this->nodes.size(); }
+
+    void remove_nodes_with_one_arc()
+    {
+        std::queue<size_t> update_queue;
+        for (const MMU_Graph::Node &node : this->nodes)
+            if (node.neighbours.size() == 1) update_queue.emplace(&node - &this->nodes.front());
+
+        while (!update_queue.empty()) {
+            size_t           node_from_idx = update_queue.front();
+            MMU_Graph::Node &node_from     = this->nodes[update_queue.front()];
+            update_queue.pop();
+            if (node_from.neighbours.empty())
+                continue;
+
+            assert(node_from.neighbours.size() == 1);
+            size_t           node_to_idx = node_from.neighbours.front().to_idx;
+            MMU_Graph::Node &node_to     = this->nodes[node_to_idx];
+            this->remove_edge(node_from_idx, node_to_idx);
+            if (node_to.neighbours.size() == 1)
+                update_queue.emplace(node_to_idx);
+        }
+    }
+
+    void add_contours(const std::vector<std::vector<ColoredLine>> &color_poly)
+    {
+        this->all_border_points = nodes.size();
+        this->polygon_sizes     = std::vector<size_t>(color_poly.size());
+        for (size_t polygon_idx = 0; polygon_idx < color_poly.size(); ++polygon_idx)
+            this->polygon_sizes[polygon_idx] = color_poly[polygon_idx].size();
+        this->polygon_idx_offset    = std::vector<size_t>(color_poly.size());
+        this->polygon_idx_offset[0] = 0;
+        for (size_t polygon_idx = 1; polygon_idx < color_poly.size(); ++polygon_idx) {
+            this->polygon_idx_offset[polygon_idx] = this->polygon_idx_offset[polygon_idx - 1] + color_poly[polygon_idx - 1].size();
+        }
+
+        size_t poly_idx = 0;
+        for (const std::vector<ColoredLine> &color_lines : color_poly) {
+            size_t line_idx = 0;
+            for (const ColoredLine &color_line : color_lines) {
+                size_t from_idx = this->get_global_index(poly_idx, line_idx);
+                size_t to_idx   = this->get_global_index(poly_idx, (line_idx + 1) % color_lines.size());
+                this->append_edge(from_idx, to_idx, color_line.color, ARC_TYPE::BORDER);
+                ++line_idx;
+            }
+            ++poly_idx;
+        }
+    }
+
+    // Nodes 0..all_border_points are only one with are on countour. Other vertexis are consider as not on coouter. So we check if base on attach index
+    inline bool is_vertex_on_contour(const Voronoi::VD::vertex_type *vertex) const
+    {
+        assert(vertex != nullptr);
+        return vertex->color() < this->all_border_points;
+    }
+
+    inline bool is_edge_attach_to_contour(const voronoi_diagram<double>::const_edge_iterator &edge_iterator) const
+    {
+        return this->is_vertex_on_contour(edge_iterator->vertex0()) || this->is_vertex_on_contour(edge_iterator->vertex1());
+    }
+
+    inline bool is_edge_connecting_two_contour_vertices(const voronoi_diagram<double>::const_edge_iterator &edge_iterator) const
+    {
+        return this->is_vertex_on_contour(edge_iterator->vertex0()) && this->is_vertex_on_contour(edge_iterator->vertex1());
+    }
+};
+
+namespace bg  = boost::geometry;
+namespace bgm = boost::geometry::model;
+namespace bgi = boost::geometry::index;
+
+// float is needed because for coord_t bgi::intersects throws "bad numeric conversion: positive overflow"
+using rtree_point_t = bgm::point<float, 2, boost::geometry::cs::cartesian>;
+using rtree_t       = bgi::rtree<std::pair<rtree_point_t, size_t>, bgi::rstar<16, 4>>;
+
+static inline rtree_point_t mk_rtree_point(const Point &pt) { return rtree_point_t(float(pt.x()), float(pt.y())); }
+
+static inline Point mk_point(const Voronoi::VD::vertex_type *point) { return Point(coord_t(point->x()), coord_t(point->y())); }
+
+static inline Point mk_point(const Voronoi::Internal::point_type &point) { return Point(coord_t(point.x()), coord_t(point.y())); }
+
+static inline Point mk_point(const voronoi_diagram<double>::vertex_type &point) { return Point(coord_t(point.x()), coord_t(point.y())); }
+
+static inline void mark_processed(const voronoi_diagram<double>::const_edge_iterator &edge_iterator)
+{
+    edge_iterator->color(true);
+    edge_iterator->twin()->color(true);
+}
+
+// Return true, if "p" is closer to line.a, then line.b
+static inline bool is_point_closer_to_beginning_of_line(const Line &line, const Point &p)
+{
+    return (p - line.a).cast<double>().squaredNorm() < (p - line.b).cast<double>().squaredNorm();
+}
+
+static inline bool has_same_color(const ColoredLine &cl1, const ColoredLine &cl2) { return cl1.color == cl2.color; }
+
+// Determines if the line points from the point between two contour lines is pointing inside polygon or outside.
+static inline bool points_inside(const Line &contour_first, const Line &contour_second, const Point &new_point)
+{
+    // Used in points_inside for decision if line leading thought the common point of two lines is pointing inside polygon or outside
+    auto three_points_inward_normal = [](const Point &left, const Point &middle, const Point &right) -> Vec2d {
+        assert(left != middle);
+        assert(middle != right);
+        return (perp(Point(middle - left)).cast<double>().normalized() + perp(Point(right - middle)).cast<double>().normalized()).normalized();
+    };
+
+    assert(contour_first.b == contour_second.a);
+    Vec2d  inward_normal = three_points_inward_normal(contour_first.a, contour_first.b, contour_second.b);
+    Vec2d  edge_norm     = (new_point - contour_first.b).cast<double>().normalized();
+    double side          = inward_normal.dot(edge_norm);
+    //    assert(side != 0.);
+    return side > 0.;
+}
+
+static inline bool line_intersection_with_epsilon(const Line &line_to_extend, const Line &other, Point *intersection)
+{
+    Line extended_line = line_to_extend;
+    extended_line.extend(15 * SCALED_EPSILON);
+    return extended_line.intersection(other, intersection);
+}
+
+// For every ColoredLine in lines_colored_out, assign the index of the polygon to which belongs and also the index of this line inside of the polygon.
+static inline void init_polygon_indices(const MMU_Graph                             &graph,
+                                        const std::vector<std::vector<ColoredLine>> &color_poly,
+                                        std::vector<ColoredLine>                    &lines_colored_out)
+{
+    size_t poly_idx = 0;
+    for (const std::vector<ColoredLine> &color_lines : color_poly) {
+        size_t line_idx = 0;
+        for (size_t color_line_idx = 0; color_line_idx < color_lines.size(); ++color_line_idx) {
+            size_t from_idx                            = graph.get_global_index(poly_idx, line_idx);
+            lines_colored_out[from_idx].poly_idx       = int(poly_idx);
+            lines_colored_out[from_idx].local_line_idx = int(line_idx);
+            ++line_idx;
+        }
+        ++poly_idx;
+    }
+}
+
+static MMU_Graph build_graph(size_t layer_idx, const std::vector<std::vector<ColoredLine>> &color_poly)
+{
+    Geometry::VoronoiDiagram vd;
+    std::vector<ColoredLine> lines_colored  = to_lines(color_poly);
+    Polygons                 color_poly_tmp = colored_points_to_polygon(color_poly);
+    const Points             points         = to_points(color_poly_tmp);
+    const Lines              lines          = to_lines(color_poly_tmp);
+
+    // The algorithm adds edges to the graph that are between two different colors.
+    // If a polygon is colored entirely with one color, we need to add at least one edge from that polygon artificially.
+    // Adding this edge is necessary for cases where the expolygon has an outer contour colored whole with one color
+    // and a hole colored with a different color. If an edge wasn't added to the graph,
+    // the entire expolygon would be colored with single random color instead of two different.
+    std::vector<bool>        force_edge_adding(color_poly.size());
+
+    // For each polygon, check if it is all colored with the same color. If it is, we need to force adding one edge to it.
+    for (const std::vector<ColoredLine> &c_poly : color_poly) {
+        bool force_edge = true;
+        for (const ColoredLine &c_line : c_poly)
+            if (c_line.color != c_poly.front().color) {
+                force_edge = false;
+                break;
+            }
+        force_edge_adding[&c_poly - &color_poly.front()] = force_edge;
+    }
+
+    boost::polygon::construct_voronoi(lines_colored.begin(), lines_colored.end(), &vd);
+    MMU_Graph graph;
+    for (const Point &point : points)
+        graph.nodes.push_back({point});
+
+    graph.add_contours(color_poly);
+    init_polygon_indices(graph, color_poly, lines_colored);
+
+    assert(graph.nodes.size() == lines_colored.size());
+
+    // All Voronoi vertices are post-processes to merge very close vertices to single. Witch Eliminates issues with intersection edges.
+    // Also, Voronoi vertices outside of the bounding of input polygons are throw away by marking them.
+    auto append_voronoi_vertices_to_graph = [&graph, &color_poly_tmp, &vd]() -> void {
+        auto is_equal_points = [](const Point &p1, const Point &p2) { return p1 == p2 || (p1 - p2).cast<double>().norm() <= 3 * SCALED_EPSILON; };
+
+        BoundingBox bbox = get_extents(color_poly_tmp);
+        bbox.offset(SCALED_EPSILON);
+        // EdgeGrid is used for vertices near to contour and rtree for other vertices
+        // FIXME Lukas H.: Get rid of EdgeGrid and rtree. Use only one structure for both cases.
+        EdgeGrid::Grid grid;
+        grid.set_bbox(bbox);
+        grid.create(color_poly_tmp, coord_t(scale_(10.)));
+        rtree_t rtree;
+        for (const voronoi_diagram<double>::vertex_type &vertex : vd.vertices()) {
+            vertex.color(-1);
+            Point vertex_point = mk_point(vertex);
+
+            const Point &first_point  = graph.nodes[graph.get_arc(vertex.incident_edge()->cell()->source_index()).from_idx].point;
+            const Point &second_point = graph.nodes[graph.get_arc(vertex.incident_edge()->twin()->cell()->source_index()).from_idx].point;
+
+            if (vertex_equal_to_point(&vertex, first_point)) {
+                assert(vertex.color() != vertex.incident_edge()->cell()->source_index());
+                assert(vertex.color() != vertex.incident_edge()->twin()->cell()->source_index());
+                vertex.color(graph.get_arc(vertex.incident_edge()->cell()->source_index()).from_idx);
+            } else if (vertex_equal_to_point(&vertex, second_point)) {
+                assert(vertex.color() != vertex.incident_edge()->cell()->source_index());
+                assert(vertex.color() != vertex.incident_edge()->twin()->cell()->source_index());
+                vertex.color(graph.get_arc(vertex.incident_edge()->twin()->cell()->source_index()).from_idx);
+            } else if (bbox.contains(vertex_point)) {
+                EdgeGrid::Grid::ClosestPointResult cp = grid.closest_point_signed_distance(vertex_point, coord_t(3 * SCALED_EPSILON));
+                if (cp.valid()) {
+                    size_t global_idx      = graph.get_global_index(cp.contour_idx, cp.start_point_idx);
+                    size_t global_idx_next = graph.get_global_index(cp.contour_idx, (cp.start_point_idx + 1) % color_poly_tmp[cp.contour_idx].points.size());
+                    vertex.color(is_equal_points(vertex_point, graph.nodes[global_idx].point) ? global_idx : global_idx_next);
+                } else {
+                    if (rtree.empty()) {
+                        rtree.insert(std::make_pair(mk_rtree_point(vertex_point), graph.nodes_count()));
+                        vertex.color(graph.nodes_count());
+                        graph.nodes.push_back({vertex_point});
+                    } else {
+                        std::vector<std::pair<rtree_point_t, size_t>> closest;
+                        rtree.query(bgi::nearest(mk_rtree_point(vertex_point), 1), std::back_inserter(closest));
+                        assert(!closest.empty());
+                        rtree_point_t r_point = closest.front().first;
+                        Point         closest_p(bg::get<0>(r_point), bg::get<1>(r_point));
+                        if (Line(vertex_point, closest_p).length() > 3 * SCALED_EPSILON) {
+                            rtree.insert(std::make_pair(mk_rtree_point(vertex_point), graph.nodes_count()));
+                            vertex.color(graph.nodes_count());
+                            graph.nodes.push_back({vertex_point});
+                        } else {
+                            vertex.color(closest.front().second);
+                        }
+                    }
+                }
+            }
+        }
+    };
+
+    append_voronoi_vertices_to_graph();
+
+    auto get_prev_contour_line = [&lines_colored, &color_poly, &graph](const voronoi_diagram<double>::const_edge_iterator &edge_it) -> ColoredLine {
+        size_t contour_line_local_idx = lines_colored[edge_it->cell()->source_index()].local_line_idx;
+        size_t contour_line_size      = color_poly[lines_colored[edge_it->cell()->source_index()].poly_idx].size();
+        size_t contour_prev_idx       = graph.get_global_index(lines_colored[edge_it->cell()->source_index()].poly_idx,
+                                                         (contour_line_local_idx > 0) ? contour_line_local_idx - 1 : contour_line_size - 1);
+        return lines_colored[contour_prev_idx];
+    };
+
+    auto get_next_contour_line = [&lines_colored, &color_poly, &graph](const voronoi_diagram<double>::const_edge_iterator &edge_it) -> ColoredLine {
+        size_t contour_line_local_idx = lines_colored[edge_it->cell()->source_index()].local_line_idx;
+        size_t contour_line_size      = color_poly[lines_colored[edge_it->cell()->source_index()].poly_idx].size();
+        size_t contour_next_idx       = graph.get_global_index(lines_colored[edge_it->cell()->source_index()].poly_idx,
+                                                         (contour_line_local_idx + 1) % contour_line_size);
+        return lines_colored[contour_next_idx];
+    };
+
+    BoundingBox bbox = get_extents(color_poly_tmp);
+    bbox.offset(scale_(10.));
+    const double bbox_dim_max = double(std::max(bbox.size().x(), bbox.size().y()));
+
+    // Make a copy of the input segments with the double type.
+    std::vector<Voronoi::Internal::segment_type> segments;
+    for (const Line &line : lines)
+        segments.emplace_back(Voronoi::Internal::point_type(double(line.a(0)), double(line.a(1))),
+                              Voronoi::Internal::point_type(double(line.b(0)), double(line.b(1))));
+
+    for (auto edge_it = vd.edges().begin(); edge_it != vd.edges().end(); ++edge_it) {
+        // Skip second half-edge
+        if (edge_it->cell()->source_index() > edge_it->twin()->cell()->source_index() || edge_it->color())
+            continue;
+
+        if (edge_it->is_infinite() && (edge_it->vertex0() != nullptr || edge_it->vertex1() != nullptr)) {
+            // Infinite edge is leading through a point on the counter, but there are no Voronoi vertices.
+            // So we could fix this case by computing the intersection between the contour line and infinity edge.
+            std::vector<Voronoi::Internal::point_type> samples;
+            Voronoi::Internal::clip_infinite_edge(points, segments, *edge_it, bbox_dim_max, &samples);
+            if (samples.empty())
+                continue;
+
+            const Line         edge_line(mk_point(samples[0]), mk_point(samples[1]));
+            const ColoredLine &contour_line = lines_colored[edge_it->cell()->source_index()];
+            Point              contour_intersection;
+
+            if (line_intersection_with_epsilon(contour_line.line, edge_line, &contour_intersection)) {
+                const MMU_Graph::Arc &graph_arc = graph.get_arc(edge_it->cell()->source_index());
+                const size_t          from_idx  = (edge_it->vertex1() != nullptr) ? edge_it->vertex1()->color() : edge_it->vertex0()->color();
+                size_t                to_idx    = ((contour_line.line.a - contour_intersection).cast<double>().squaredNorm() <
+                                 (contour_line.line.b - contour_intersection).cast<double>().squaredNorm()) ?
+                                                      graph_arc.from_idx :
+                                                      graph_arc.to_idx;
+                if (from_idx != to_idx && from_idx < graph.nodes_count() && to_idx < graph.nodes_count()) {
+                    graph.append_edge(from_idx, to_idx);
+                    mark_processed(edge_it);
+                }
+            }
+        } else if (edge_it->is_finite()) {
+            const Point  v0       = mk_point(edge_it->vertex0());
+            const Point  v1       = mk_point(edge_it->vertex1());
+            const size_t from_idx = edge_it->vertex0()->color();
+            const size_t to_idx   = edge_it->vertex1()->color();
+
+            // Both points are on contour, so skip them. In cases of duplicate Voronoi vertices, skip edges between the same two points.
+            if (graph.is_edge_connecting_two_contour_vertices(edge_it) || (edge_it->vertex0()->color() == edge_it->vertex1()->color())) continue;
+
+            const Line        edge_line(v0, v1);
+            const Line        contour_line      = lines_colored[edge_it->cell()->source_index()].line;
+            const ColoredLine colored_line      = lines_colored[edge_it->cell()->source_index()];
+            const ColoredLine contour_line_prev = get_prev_contour_line(edge_it);
+            const ColoredLine contour_line_next = get_next_contour_line(edge_it);
+
+            Point intersection;
+            if (edge_it->vertex0()->color() >= graph.nodes_count() || edge_it->vertex1()->color() >= graph.nodes_count()) {
+//                if(edge_it->vertex0()->color() < graph.nodes_count() && !graph.is_vertex_on_contour(edge_it->vertex0())) {
+//
+//                }
+                if (edge_it->vertex1()->color() < graph.nodes_count() && !graph.is_vertex_on_contour(edge_it->vertex1())) {
+                    Line contour_line_twin = lines_colored[edge_it->twin()->cell()->source_index()].line;
+                    if (line_intersection_with_epsilon(contour_line_twin, edge_line, &intersection)) {
+                        const MMU_Graph::Arc &graph_arc = graph.get_arc(edge_it->twin()->cell()->source_index());
+                        const size_t          to_idx_l  = is_point_closer_to_beginning_of_line(contour_line_twin, intersection) ? graph_arc.from_idx :
+                                                                                                                                  graph_arc.to_idx;
+                        graph.append_edge(edge_it->vertex1()->color(), to_idx_l);
+                    } else if (line_intersection_with_epsilon(contour_line, edge_line, &intersection)) {
+                        const MMU_Graph::Arc &graph_arc = graph.get_arc(edge_it->cell()->source_index());
+                        const size_t to_idx_l = is_point_closer_to_beginning_of_line(contour_line, intersection) ? graph_arc.from_idx : graph_arc.to_idx;
+                        graph.append_edge(edge_it->vertex1()->color(), to_idx_l);
+                    }
+                    mark_processed(edge_it);
+                }
+            } else if (graph.is_edge_attach_to_contour(edge_it)) {
+                mark_processed(edge_it);
+                // Skip edges witch connection two points on a contour
+                if (graph.is_edge_connecting_two_contour_vertices(edge_it))
+                    continue;
+
+                if (graph.is_vertex_on_contour(edge_it->vertex0())) {
+                    if (is_point_closer_to_beginning_of_line(contour_line, v0)) {
+                        if ((!has_same_color(contour_line_prev, colored_line) || force_edge_adding[colored_line.poly_idx]) && points_inside(contour_line_prev.line, contour_line, v1)) {
+                            graph.append_edge(from_idx, to_idx);
+                            force_edge_adding[colored_line.poly_idx] = false;
+                        }
+                    } else {
+                        if ((!has_same_color(contour_line_next, colored_line) || force_edge_adding[colored_line.poly_idx]) && points_inside(contour_line, contour_line_next.line, v1)) {
+                            graph.append_edge(from_idx, to_idx);
+                            force_edge_adding[colored_line.poly_idx] = false;
+                        }
+                    }
+                } else {
+                    assert(graph.is_vertex_on_contour(edge_it->vertex1()));
+                    if (is_point_closer_to_beginning_of_line(contour_line, v1)) {
+                        if ((!has_same_color(contour_line_prev, colored_line) || force_edge_adding[colored_line.poly_idx]) && points_inside(contour_line_prev.line, contour_line, v0)) {
+                            graph.append_edge(from_idx, to_idx);
+                            force_edge_adding[colored_line.poly_idx] = false;
+                        }
+                    } else {
+                        if ((!has_same_color(contour_line_next, colored_line) || force_edge_adding[colored_line.poly_idx]) && points_inside(contour_line, contour_line_next.line, v0)) {
+                            graph.append_edge(from_idx, to_idx);
+                            force_edge_adding[colored_line.poly_idx] = false;
+                        }
+                    }
+                }
+            } else if (line_intersection_with_epsilon(contour_line, edge_line, &intersection)) {
+                mark_processed(edge_it);
+                Point real_v0 = graph.nodes[edge_it->vertex0()->color()].point;
+                Point real_v1 = graph.nodes[edge_it->vertex1()->color()].point;
+
+                if (is_point_closer_to_beginning_of_line(contour_line, intersection)) {
+                    Line first_part(intersection, real_v0);
+                    Line second_part(intersection, real_v1);
+
+                    if (!has_same_color(contour_line_prev, colored_line)) {
+                        if (points_inside(contour_line_prev.line, contour_line, first_part.b)) {
+                            graph.append_edge(edge_it->vertex0()->color(), graph.get_arc(edge_it->cell()->source_index()).from_idx);
+                        }
+                        if (points_inside(contour_line_prev.line, contour_line, second_part.b)) {
+                            graph.append_edge(edge_it->vertex1()->color(), graph.get_arc(edge_it->cell()->source_index()).from_idx);
+                        }
+                    }
+                } else {
+                    const size_t int_point_idx = graph.get_arc(edge_it->cell()->source_index()).to_idx;
+                    const Point  int_point     = graph.nodes[int_point_idx].point;
+
+                    const Line first_part(int_point, real_v0);
+                    const Line second_part(int_point, real_v1);
+
+                    if (!has_same_color(contour_line_next, colored_line)) {
+                        if (points_inside(contour_line, contour_line_next.line, first_part.b)) {
+                            graph.append_edge(edge_it->vertex0()->color(), int_point_idx);
+                        }
+                        if (points_inside(contour_line, contour_line_next.line, second_part.b)) {
+                            graph.append_edge(edge_it->vertex1()->color(), int_point_idx);
+                        }
+                    }
+                }
+            }
+        }
+    }
+
+    for (auto edge_it = vd.edges().begin(); edge_it != vd.edges().end(); ++edge_it) {
+        // Skip second half-edge and processed edges
+        if (edge_it->cell()->source_index() > edge_it->twin()->cell()->source_index() || edge_it->color())
+            continue;
+
+        if (edge_it->is_finite() && !bool(edge_it->color()) && edge_it->vertex0()->color() < graph.nodes_count() &&
+            edge_it->vertex1()->color() < graph.nodes_count()) {
+            // Skip cases, when the edge is between two same vertices, which is in cases two near vertices were merged together.
+            if (edge_it->vertex0()->color() == edge_it->vertex1()->color())
+                continue;
+
+            size_t from_idx = edge_it->vertex0()->color();
+            size_t to_idx   = edge_it->vertex1()->color();
+            graph.append_edge(from_idx, to_idx);
+        }
+        mark_processed(edge_it);
+    }
+
+    graph.remove_nodes_with_one_arc();
+    return graph;
+}
+
+static inline Polygon to_polygon(const Lines &lines)
+{
+    Polygon poly_out;
+    poly_out.points.reserve(lines.size());
+    for (const Line &line : lines)
+        poly_out.points.emplace_back(line.a);
+    return poly_out;
+}
+
+// Returns list of polygons and assigned colors.
+// It iterates through all nodes on the border between two different colors, and from this point,
+// start selection always left most edges for every node to construct CCW polygons.
+// Assumes that graph is planar (without self-intersection edges)
+static std::vector<std::pair<Polygon, size_t>> extract_colored_segments(MMU_Graph &graph)
+{
+    // When there is no next arc, then is returned original_arc or edge with is marked as used
+    auto get_next = [&graph](const Line &process_line, MMU_Graph::Arc &original_arc) -> MMU_Graph::Arc & {
+        std::vector<std::pair<MMU_Graph::Arc *, double>> sorted_arcs;
+        for (MMU_Graph::Arc &arc : graph.nodes[original_arc.to_idx].neighbours) {
+            if (graph.nodes[arc.to_idx].point == process_line.a || arc.used)
+                continue;
+
+            assert(original_arc.to_idx == arc.from_idx);
+            Vec2d process_line_vec_n   = (process_line.a - process_line.b).cast<double>().normalized();
+            Vec2d neighbour_line_vec_n = (graph.nodes[arc.to_idx].point - graph.nodes[arc.from_idx].point).cast<double>().normalized();
+
+            double angle = ::acos(std::clamp(neighbour_line_vec_n.dot(process_line_vec_n), -1.0, 1.0));
+            if (Slic3r::cross2(neighbour_line_vec_n, process_line_vec_n) < 0.0)
+                angle = 2.0 * (double) PI - angle;
+
+            sorted_arcs.emplace_back(&arc, angle);
+        }
+
+        std::sort(sorted_arcs.begin(), sorted_arcs.end(),
+                  [](std::pair<MMU_Graph::Arc *, double> &l, std::pair<MMU_Graph::Arc *, double> &r) -> bool { return l.second < r.second; });
+
+        // Try to return left most edge witch is unused
+        for (auto &sorted_arc : sorted_arcs)
+            if (!sorted_arc.first->used)
+                return *sorted_arc.first;
+
+        if (sorted_arcs.empty())
+            return original_arc;
+
+        return *(sorted_arcs.front().first);
+    };
+
+    std::vector<std::pair<Polygon, size_t>> polygons_segments;
+    for (MMU_Graph::Node &node : graph.nodes)
+        for (MMU_Graph::Arc &arc : node.neighbours)
+            arc.used = false;
+
+    for (size_t node_idx = 0; node_idx < graph.all_border_points; ++node_idx) {
+        MMU_Graph::Node &node = graph.nodes[node_idx];
+
+        for (MMU_Graph::Arc &arc : node.neighbours) {
+            if (arc.type == MMU_Graph::ARC_TYPE::NON_BORDER || arc.used) continue;
+
+            Line process_line(node.point, graph.nodes[arc.to_idx].point);
+            arc.used = true;
+
+            Lines face_lines;
+            face_lines.emplace_back(process_line);
+            Point start_p = process_line.a;
+
+            Line            p_vec = process_line;
+            MMU_Graph::Arc *p_arc = &arc;
+            do {
+                MMU_Graph::Arc &next = get_next(p_vec, *p_arc);
+                face_lines.emplace_back(Line(graph.nodes[next.from_idx].point, graph.nodes[next.to_idx].point));
+                if (next.used) break;
+
+                next.used = true;
+                p_vec     = Line(graph.nodes[next.from_idx].point, graph.nodes[next.to_idx].point);
+                p_arc     = &next;
+            } while (graph.nodes[p_arc->to_idx].point != start_p);
+
+            Polygon poly = to_polygon(face_lines);
+            if (poly.is_counter_clockwise() && poly.is_valid())
+                polygons_segments.emplace_back(poly, arc.color);
+        }
+    }
+    return polygons_segments;
+}
+
+// Used in remove_multiple_edges_in_vertices()
+// Returns length of edge with is connected to contour. To this length is include other edges with follows it if they are almost straight (with the
+// tolerance of 15) And also if node between two subsequent edges is connected only to these two edges.
+static inline double compute_edge_length(MMU_Graph &graph, size_t start_idx, MMU_Graph::Arc &start_edge)
+{
+    for (MMU_Graph::Node &node : graph.nodes)
+        for (MMU_Graph::Arc &arc : node.neighbours)
+            arc.used = false;
+
+    start_edge.used                   = true;
+    MMU_Graph::Arc *arc               = &start_edge;
+    size_t          idx               = start_idx;
+    double          line_total_length = Line(graph.nodes[idx].point, graph.nodes[arc->to_idx].point).length();
+    while (graph.nodes[arc->to_idx].neighbours.size() == 2) {
+        bool found = false;
+        for (MMU_Graph::Arc &arc_n : graph.nodes[arc->to_idx].neighbours) {
+            if (arc_n.type == MMU_Graph::ARC_TYPE::NON_BORDER && !arc_n.used && arc_n.to_idx != idx) {
+                Line first_line(graph.nodes[idx].point, graph.nodes[arc->to_idx].point);
+                Line second_line(graph.nodes[arc->to_idx].point, graph.nodes[arc_n.to_idx].point);
+
+                Vec2d  first_line_vec    = (first_line.a - first_line.b).cast<double>();
+                Vec2d  second_line_vec   = (second_line.b - second_line.a).cast<double>();
+                Vec2d  first_line_vec_n  = first_line_vec.normalized();
+                Vec2d  second_line_vec_n = second_line_vec.normalized();
+                double angle             = ::acos(std::clamp(first_line_vec_n.dot(second_line_vec_n), -1.0, 1.0));
+                if (Slic3r::cross2(first_line_vec_n, second_line_vec_n) < 0.0)
+                    angle = 2.0 * (double) PI - angle;
+
+                if (std::abs(angle - PI) >= (PI / 12))
+                    continue;
+
+                idx = arc->to_idx;
+                arc = &arc_n;
+
+                line_total_length += Line(graph.nodes[idx].point, graph.nodes[arc->to_idx].point).length();
+                arc_n.used = true;
+                found      = true;
+                break;
+            }
+        }
+        if (!found)
+            break;
+    }
+
+    return line_total_length;
+}
+
+// Used for fixing double Voronoi edges for concave parts of the polygon.
+static void remove_multiple_edges_in_vertices(MMU_Graph &graph, const std::vector<std::vector<ColoredLine>> &color_poly)
+{
+    std::vector<std::vector<std::pair<size_t, size_t>>> colored_segments = get_all_segments(color_poly);
+    for (const std::vector<std::pair<size_t, size_t>> &colored_segment_p : colored_segments) {
+        size_t poly_idx = &colored_segment_p - &colored_segments.front();
+        for (const std::pair<size_t, size_t> &colored_segment : colored_segment_p) {
+            size_t first_idx  = graph.get_global_index(poly_idx, colored_segment.first);
+            size_t second_idx = graph.get_global_index(poly_idx, (colored_segment.second + 1) % graph.polygon_sizes[poly_idx]);
+            Line   seg_line(graph.nodes[first_idx].point, graph.nodes[second_idx].point);
+
+            if (graph.nodes[first_idx].neighbours.size() >= 3) {
+                std::vector<std::pair<MMU_Graph::Arc *, double>> arc_to_check;
+                for (MMU_Graph::Arc &n_arc : graph.nodes[first_idx].neighbours) {
+                    if (n_arc.type == MMU_Graph::ARC_TYPE::NON_BORDER) {
+                        double total_len = compute_edge_length(graph, first_idx, n_arc);
+                        arc_to_check.emplace_back(&n_arc, total_len);
+                    }
+                }
+                std::sort(arc_to_check.begin(), arc_to_check.end(),
+                          [](std::pair<MMU_Graph::Arc *, double> &l, std::pair<MMU_Graph::Arc *, double> &r) -> bool { return l.second > r.second; });
+
+                while (arc_to_check.size() > 1) {
+                    graph.remove_edge(first_idx, arc_to_check.back().first->to_idx);
+                    arc_to_check.pop_back();
+                }
+            }
+        }
+    }
+}
+
+static void cut_segmented_layers(const std::vector<ExPolygons> &input_expolygons, std::vector<std::vector<std::pair<ExPolygon, size_t>>> &segmented_regions, const float cut_width) {
+    tbb::parallel_for(tbb::blocked_range<size_t>(0, segmented_regions.size()),[&](const tbb::blocked_range<size_t>& range) {
+        for (size_t layer_idx = range.begin(); layer_idx < range.end(); ++layer_idx) {
+            std::vector<std::pair<ExPolygon, size_t>> segmented_regions_cuts;
+            for (const std::pair<ExPolygon, size_t> &colored_expoly : segmented_regions[layer_idx]) {
+                ExPolygons cut_colored_expoly = diff_ex(colored_expoly.first, offset_ex(input_expolygons[layer_idx], cut_width));
+                for (const ExPolygon &expoly : cut_colored_expoly) {
+                    segmented_regions_cuts.emplace_back(expoly, colored_expoly.second);
+                }
+            }
+            segmented_regions[layer_idx] = segmented_regions_cuts;
+        }
+    }); // end of parallel_for
+}
+
+// Returns MMU segmentation of top and bottom layers based on painting in MMU segmentation gizmo
+static inline std::vector<std::vector<ExPolygons>> mmu_segmentation_top_and_bottom_layers(const PrintObject &print_object, const std::vector<ExPolygons> &input_expolygons)
+{
+    const size_t num_extruders = print_object.print()->config().nozzle_diameter.size();
+    const ConstLayerPtrsAdaptor layers = print_object.layers();
+    std::vector<std::vector<ExPolygons>> triangles_by_color(num_extruders);
+    triangles_by_color.assign(num_extruders, std::vector<ExPolygons>(layers.size()));
+    for (const ModelVolume *mv : print_object.model_object()->volumes) {
+        for (size_t extruder_idx = 0; extruder_idx < num_extruders; ++extruder_idx) {
+            const indexed_triangle_set custom_facets = mv->mmu_segmentation_facets.get_facets(*mv, EnforcerBlockerType(extruder_idx));
+            if (!mv->is_model_part() || custom_facets.indices.empty())
+                continue;
+
+            const Transform3f tr = print_object.trafo().cast<float>() * mv->get_matrix().cast<float>();
+            for (size_t facet_idx = 0; facet_idx < custom_facets.indices.size(); ++facet_idx) {
+                float min_z = std::numeric_limits<float>::max();
+                float max_z = std::numeric_limits<float>::lowest();
+
+                std::array<Vec3f, 3> facet;
+                Points               projected_facet(3);
+                for (int p_idx = 0; p_idx < 3; ++p_idx) {
+                    facet[p_idx] = tr * custom_facets.vertices[custom_facets.indices[facet_idx](p_idx)];
+                    max_z        = std::max(max_z, facet[p_idx].z());
+                    min_z        = std::min(min_z, facet[p_idx].z());
+                }
+
+                // Sort the vertices by z-axis for simplification of projected_facet on slices
+                std::sort(facet.begin(), facet.end(), [](const Vec3f &p1, const Vec3f &p2) { return p1.z() < p2.z(); });
+
+                for (int p_idx = 0; p_idx < 3; ++p_idx) {
+                    projected_facet[p_idx] = Point(scale_(facet[p_idx].x()), scale_(facet[p_idx].y()));
+                    projected_facet[p_idx] = projected_facet[p_idx] - print_object.center_offset();
+                }
+
+                ExPolygon triangle = ExPolygon(projected_facet);
+
+                // Find lowest slice not below the triangle.
+                auto first_layer = std::upper_bound(layers.begin(), layers.end(), float(min_z - EPSILON),
+                                                    [](float z, const Layer *l1) { return z < l1->slice_z + l1->height * 0.5; });
+                auto last_layer  = std::upper_bound(layers.begin(), layers.end(), float(max_z - EPSILON),
+                                                   [](float z, const Layer *l1) { return z < l1->slice_z + l1->height * 0.5; });
+
+                if (last_layer == layers.end())
+                    --last_layer;
+
+                if (first_layer == layers.end() || (first_layer != layers.begin() && facet[0].z() < (*first_layer)->print_z - EPSILON))
+                    --first_layer;
+
+                for (auto layer_it = first_layer; (layer_it != (last_layer + 1) && layer_it != layers.end()); ++layer_it) {
+                    size_t layer_idx = layer_it - layers.begin();
+                    triangles_by_color[extruder_idx][layer_idx].emplace_back(triangle);
+                }
+            }
+        }
+    }
+
+    auto get_extrusion_width = [&layers = std::as_const(layers)](const size_t layer_idx) -> float {
+        auto extrusion_width_it = std::max_element(layers[layer_idx]->regions().begin(), layers[layer_idx]->regions().end(),
+                                                   [](const LayerRegion *l1, const LayerRegion *l2) {
+                                                       return l1->region().config().perimeter_extrusion_width <
+                                                              l2->region().config().perimeter_extrusion_width;
+                                                   });
+        assert(extrusion_width_it != layers[layer_idx]->regions().end());
+        return float((*extrusion_width_it)->region().config().perimeter_extrusion_width);
+    };
+
+    auto get_top_solid_layers = [&layers = std::as_const(layers)](const size_t layer_idx) -> int {
+        auto top_solid_layer_it = std::max_element(layers[layer_idx]->regions().begin(), layers[layer_idx]->regions().end(),
+                                                   [](const LayerRegion *l1, const LayerRegion *l2) {
+                                                       return l1->region().config().top_solid_layers < l2->region().config().top_solid_layers;
+                                                   });
+        assert(top_solid_layer_it != layers[layer_idx]->regions().end());
+        return (*top_solid_layer_it)->region().config().top_solid_layers;
+    };
+
+    auto get_bottom_solid_layers = [&layers = std::as_const(layers)](const size_t layer_idx) -> int {
+        auto top_bottom_layer_it = std::max_element(layers[layer_idx]->regions().begin(), layers[layer_idx]->regions().end(),
+                                                    [](const LayerRegion *l1, const LayerRegion *l2) {
+                                                        return l1->region().config().bottom_solid_layers < l2->region().config().bottom_solid_layers;
+                                                    });
+        assert(top_bottom_layer_it != layers[layer_idx]->regions().end());
+        return (*top_bottom_layer_it)->region().config().bottom_solid_layers;
+    };
+
+    std::vector<ExPolygons> top_layers(input_expolygons.size());
+    top_layers.back() = input_expolygons.back();
+    tbb::parallel_for(tbb::blocked_range<size_t>(1, input_expolygons.size()), [&](const tbb::blocked_range<size_t> &range) {
+        for (size_t layer_idx = range.begin(); layer_idx < range.end(); ++layer_idx) {
+            float extrusion_width     = 0.1f * float(scale_(get_extrusion_width(layer_idx)));
+            top_layers[layer_idx - 1] = diff_ex(input_expolygons[layer_idx - 1], offset_ex(input_expolygons[layer_idx], extrusion_width));
+        }
+    }); // end of parallel_for
+
+    std::vector<ExPolygons> bottom_layers(input_expolygons.size());
+    bottom_layers.front() = input_expolygons.front();
+    tbb::parallel_for(tbb::blocked_range<size_t>(0, input_expolygons.size() - 1), [&](const tbb::blocked_range<size_t> &range) {
+        for (size_t layer_idx = range.begin(); layer_idx < range.end(); ++layer_idx) {
+            float extrusion_width        = 0.1f * float(scale_(get_extrusion_width(layer_idx)));
+            bottom_layers[layer_idx + 1] = diff_ex(input_expolygons[layer_idx + 1], offset_ex(input_expolygons[layer_idx], extrusion_width));
+        }
+    }); // end of parallel_for
+
+    tbb::parallel_for(tbb::blocked_range<size_t>(0, input_expolygons.size()), [&](const tbb::blocked_range<size_t> &range) {
+        for (size_t layer_idx = range.begin(); layer_idx < range.end(); ++layer_idx) {
+            float extrusion_width = 0.1f * float(scale_(get_extrusion_width(layer_idx)));
+            for (std::vector<ExPolygons> &triangles : triangles_by_color) {
+                if (!triangles[layer_idx].empty() && (!top_layers[layer_idx].empty() || !bottom_layers[layer_idx].empty())) {
+                    ExPolygons connected = union_ex(offset_ex(triangles[layer_idx], float(10 * SCALED_EPSILON)));
+                    triangles[layer_idx] = union_ex(offset_ex(offset_ex(connected, -extrusion_width / 1), extrusion_width / 1));
+                } else {
+                    triangles[layer_idx].clear();
+                }
+            }
+        }
+    }); // end of parallel_for
+
+    std::vector<std::vector<ExPolygons>> triangles_by_color_bottom(num_extruders);
+    std::vector<std::vector<ExPolygons>> triangles_by_color_top(num_extruders);
+    triangles_by_color_bottom.assign(num_extruders, std::vector<ExPolygons>(input_expolygons.size()));
+    triangles_by_color_top.assign(num_extruders, std::vector<ExPolygons>(input_expolygons.size()));
+
+    for (size_t layer_idx = 0; layer_idx < input_expolygons.size(); ++layer_idx) {
+        BOOST_LOG_TRIVIAL(debug) << "MMU segmentation of top layer: " << layer_idx;
+        float      extrusion_width  = scale_(get_extrusion_width(layer_idx));
+        int        top_solid_layers = get_top_solid_layers(layer_idx);
+        ExPolygons top_expolygon    = top_layers[layer_idx];
+        if (top_expolygon.empty())
+            continue;
+
+        for (size_t color_idx = 0; color_idx < triangles_by_color.size(); ++color_idx) {
+            if (triangles_by_color[color_idx][layer_idx].empty())
+                continue;
+            ExPolygons intersection_poly = intersection_ex(triangles_by_color[color_idx][layer_idx], top_expolygon);
+            if (!intersection_poly.empty()) {
+                triangles_by_color_top[color_idx][layer_idx].insert(triangles_by_color_top[color_idx][layer_idx].end(), intersection_poly.begin(),
+                                                                    intersection_poly.end());
+                for (int last_idx = int(layer_idx) - 1; last_idx >= std::max(int(layer_idx - top_solid_layers), int(0)); --last_idx) {
+                    float offset_value = float(layer_idx - last_idx) * (-1.0f) * extrusion_width;
+                    if (offset_ex(top_expolygon, offset_value).empty())
+                        continue;
+                    ExPolygons layer_slices_trimmed = input_expolygons[last_idx];
+
+                    for (int last_idx_1 = last_idx; last_idx_1 < int(layer_idx); ++last_idx_1) {
+                        layer_slices_trimmed = intersection_ex(layer_slices_trimmed, input_expolygons[last_idx_1 + 1]);
+                    }
+
+                    ExPolygons offset_e            = offset_ex(layer_slices_trimmed, offset_value);
+                    ExPolygons intersection_poly_2 = intersection_ex(triangles_by_color_top[color_idx][layer_idx], offset_e);
+                    triangles_by_color_top[color_idx][last_idx].insert(triangles_by_color_top[color_idx][last_idx].end(), intersection_poly_2.begin(),
+                                                                       intersection_poly_2.end());
+                }
+            }
+        }
+    }
+
+    for (size_t layer_idx = 0; layer_idx < input_expolygons.size(); ++layer_idx) {
+        BOOST_LOG_TRIVIAL(debug) << "MMU segmentation of bottom layer: " << layer_idx;
+        float             extrusion_width     = scale_(get_extrusion_width(layer_idx));
+        int               bottom_solid_layers = get_bottom_solid_layers(layer_idx);
+        const ExPolygons &bottom_expolygon    = bottom_layers[layer_idx];
+        if (bottom_expolygon.empty())
+            continue;
+
+        for (size_t color_idx = 0; color_idx < triangles_by_color.size(); ++color_idx) {
+            if (triangles_by_color[color_idx][layer_idx].empty())
+                continue;
+
+            ExPolygons intersection_poly = intersection_ex(triangles_by_color[color_idx][layer_idx], bottom_expolygon);
+            if (!intersection_poly.empty()) {
+                triangles_by_color_bottom[color_idx][layer_idx].insert(triangles_by_color_bottom[color_idx][layer_idx].end(), intersection_poly.begin(),
+                                                                       intersection_poly.end());
+                for (size_t last_idx = layer_idx + 1; last_idx < std::min(layer_idx + bottom_solid_layers, input_expolygons.size()); ++last_idx) {
+                    float offset_value = float(last_idx - layer_idx) * (-1.0f) * extrusion_width;
+                    if (offset_ex(bottom_expolygon, offset_value).empty())
+                        continue;
+                    ExPolygons layer_slices_trimmed = input_expolygons[last_idx];
+                    for (int last_idx_1 = int(last_idx); last_idx_1 > int(layer_idx); --last_idx_1) {
+                        layer_slices_trimmed = intersection_ex(layer_slices_trimmed, offset_ex(input_expolygons[last_idx_1 - 1], offset_value));
+                    }
+
+                    ExPolygons offset_e            = offset_ex(layer_slices_trimmed, offset_value);
+                    ExPolygons intersection_poly_2 = intersection_ex(triangles_by_color_bottom[color_idx][layer_idx], offset_e);
+                    append(triangles_by_color_bottom[color_idx][last_idx], std::move(intersection_poly_2));
+                }
+            }
+        }
+    }
+
+    std::vector<std::vector<ExPolygons>> triangles_by_color_merged(num_extruders);
+    triangles_by_color_merged.assign(num_extruders, std::vector<ExPolygons>(input_expolygons.size()));
+    for (size_t layer_idx = 0; layer_idx < input_expolygons.size(); ++layer_idx) {
+        for (size_t color_idx = 0; color_idx < triangles_by_color_merged.size(); ++color_idx) {
+            auto &self = triangles_by_color_merged[color_idx][layer_idx];
+            append(self, std::move(triangles_by_color_bottom[color_idx][layer_idx]));
+            append(self, std::move(triangles_by_color_top[color_idx][layer_idx]));
+            self = union_ex(self);
+        }
+
+        // Cut all colors for cases when two colors are overlapping
+        for (size_t color_idx = 1; color_idx < triangles_by_color_merged.size(); ++color_idx) {
+            triangles_by_color_merged[color_idx][layer_idx] = diff_ex(triangles_by_color_merged[color_idx][layer_idx],
+                                                                      triangles_by_color_merged[color_idx - 1][layer_idx]);
+        }
+    }
+
+    return triangles_by_color_merged;
+}
+
+static std::vector<std::vector<std::pair<ExPolygon, size_t>>> merge_segmented_layers(
+    const std::vector<std::vector<std::pair<ExPolygon, size_t>>> &segmented_regions, std::vector<std::vector<ExPolygons>> &&top_and_bottom_layers)
+{
+    std::vector<std::vector<std::pair<ExPolygon, size_t>>> segmented_regions_merged(segmented_regions.size());
+
+    tbb::parallel_for(tbb::blocked_range<size_t>(0, segmented_regions.size()), [&](const tbb::blocked_range<size_t> &range) {
+        for (size_t layer_idx = range.begin(); layer_idx < range.end(); ++layer_idx) {
+            BOOST_LOG_TRIVIAL(debug) << "MMU segmentation - merging region: " << layer_idx;
+            for (const std::pair<ExPolygon, size_t> &colored_expoly : segmented_regions[layer_idx]) {
+                ExPolygons cut_colored_expoly = {colored_expoly.first};
+                for (const std::vector<ExPolygons> &top_and_bottom_layer : top_and_bottom_layers)
+                    cut_colored_expoly = diff_ex(cut_colored_expoly, top_and_bottom_layer[layer_idx]);
+                for (ExPolygon &ex_poly : cut_colored_expoly)
+                    segmented_regions_merged[layer_idx].emplace_back(std::move(ex_poly), colored_expoly.second);
+            }
+
+            for (size_t color_idx = 0; color_idx < top_and_bottom_layers.size(); ++color_idx)
+                for (ExPolygon &expoly : top_and_bottom_layers[color_idx][layer_idx])
+                    segmented_regions_merged[layer_idx].emplace_back(std::move(expoly), color_idx);
+        }
+    }); // end of parallel_for
+
+    return segmented_regions_merged;
+}
+
+std::vector<std::vector<std::pair<ExPolygon, size_t>>> multi_material_segmentation_by_painting(const PrintObject &print_object)
+{
+    std::vector<std::vector<std::pair<ExPolygon, size_t>>> segmented_regions(print_object.layers().size());
+    std::vector<std::vector<PaintedLine>>                  painted_lines(print_object.layers().size());
+    std::vector<EdgeGrid::Grid>                            edge_grids(print_object.layers().size());
+    const ConstLayerPtrsAdaptor                            layers = print_object.layers();
+    std::vector<ExPolygons>                                input_expolygons(layers.size());
+    std::vector<Polygons>                                  input_polygons(layers.size());
+
+    // Merge all regions and remove small holes
+    tbb::parallel_for(tbb::blocked_range<size_t>(0, layers.size()), [&](const tbb::blocked_range<size_t> &range) {
+        for (size_t layer_idx = range.begin(); layer_idx < range.end(); ++layer_idx) {
+            ExPolygons ex_polygons;
+            for (LayerRegion *region : layers[layer_idx]->regions())
+                for (const Surface &surface : region->slices.surfaces)
+                    Slic3r::append(ex_polygons, offset_ex(surface.expolygon, float(SCALED_EPSILON)));
+            // All expolygons are expanded by SCALED_EPSILON, merged, and then shrunk again by SCALED_EPSILON
+            // to ensure that very close polygons will be merged.
+            ex_polygons = union_ex(ex_polygons);
+            // Remove all expolygons and holes with an area less than 0.01mm^2
+            remove_small_and_small_holes(ex_polygons, Slic3r::sqr(scale_(0.1f)));
+            // Occasionally, some input polygons contained self-intersections that caused problems with Voronoi diagrams
+            // and consequently with the extraction of colored segments by function extract_colored_segments.
+            // Calling simplify_polygons removes these self-intersections.
+            // Also, occasionally input polygons contained several points very close together (distance between points is 1 or so).
+            // Such close points sometimes caused that the Voronoi diagram has self-intersecting edges around these vertices.
+            // This consequently leads to issues with the extraction of colored segments by function extract_colored_segments.
+            // Calling expolygons_simplify fixed these issues.
+            input_expolygons[layer_idx] = simplify_polygons_ex(to_polygons(expolygons_simplify(offset_ex(ex_polygons, float(-SCALED_EPSILON)), SCALED_EPSILON)));
+            input_polygons[layer_idx] = to_polygons(input_expolygons[layer_idx]);
+        }
+    }); // end of parallel_for
+
+    for (size_t layer_idx = 0; layer_idx < layers.size(); ++layer_idx) {
+        BoundingBox  bbox(get_extents(input_expolygons[layer_idx]));
+        bbox.offset(SCALED_EPSILON);
+        edge_grids[layer_idx].set_bbox(bbox);
+        edge_grids[layer_idx].create(input_expolygons[layer_idx], coord_t(scale_(10.)));
+    }
+
+    for (const ModelVolume *mv : print_object.model_object()->volumes) {
+        const size_t num_extruders = print_object.print()->config().nozzle_diameter.size();
+        for (size_t extruder_idx = 1; extruder_idx < num_extruders; ++extruder_idx) {
+            const indexed_triangle_set custom_facets = mv->mmu_segmentation_facets.get_facets(*mv, EnforcerBlockerType(extruder_idx));
+            if (!mv->is_model_part() || custom_facets.indices.empty())
+                continue;
+
+            const Transform3f tr = print_object.trafo().cast<float>() * mv->get_matrix().cast<float>();
+            for (size_t facet_idx = 0; facet_idx < custom_facets.indices.size(); ++facet_idx) {
+                float min_z = std::numeric_limits<float>::max();
+                float max_z = std::numeric_limits<float>::lowest();
+
+                std::array<Vec3f, 3> facet;
+                for (int p_idx = 0; p_idx < 3; ++p_idx) {
+                    facet[p_idx] = tr * custom_facets.vertices[custom_facets.indices[facet_idx](p_idx)];
+                    max_z        = std::max(max_z, facet[p_idx].z());
+                    min_z        = std::min(min_z, facet[p_idx].z());
+                }
+
+                // Sort the vertices by z-axis for simplification of projected_facet on slices
+                std::sort(facet.begin(), facet.end(), [](const Vec3f &p1, const Vec3f &p2) { return p1.z() < p2.z(); });
+
+                // Find lowest slice not below the triangle.
+                auto first_layer = std::upper_bound(print_object.layers().begin(), print_object.layers().end(), float(min_z - EPSILON),
+                                                    [](float z, const Layer *l1) { return z < l1->slice_z; });
+                auto last_layer  = std::upper_bound(print_object.layers().begin(), print_object.layers().end(), float(max_z + EPSILON),
+                                                   [](float z, const Layer *l1) { return z < l1->slice_z; });
+                --last_layer;
+
+                for (auto layer_it = first_layer; layer_it != (last_layer + 1); ++layer_it) {
+                    const Layer *layer     = *layer_it;
+                    size_t       layer_idx = layer_it - print_object.layers().begin();
+                    if (facet[0].z() > layer->slice_z || layer->slice_z > facet[2].z())
+                        continue;
+
+                    // https://kandepet.com/3d-printing-slicing-3d-objects/
+                    float t            = (float(layer->slice_z) - facet[0].z()) / (facet[2].z() - facet[0].z());
+                    Vec3f line_start_f = facet[0] + t * (facet[2] - facet[0]);
+                    Vec3f line_end_f;
+
+                    if (facet[1].z() > layer->slice_z) {
+                        // [P0, P2] a [P0, P1]
+                        float t1   = (float(layer->slice_z) - facet[0].z()) / (facet[1].z() - facet[0].z());
+                        line_end_f = facet[0] + t1 * (facet[1] - facet[0]);
+                    } else if (facet[1].z() <= layer->slice_z) {
+                        // [P0, P2] a [P1, P2]
+                        float t2   = (float(layer->slice_z) - facet[1].z()) / (facet[2].z() - facet[1].z());
+                        line_end_f = facet[1] + t2 * (facet[2] - facet[1]);
+                    }
+
+                    Point line_start(scale_(line_start_f.x()), scale_(line_start_f.y()));
+                    Point line_end(scale_(line_end_f.x()), scale_(line_end_f.y()));
+                    line_start -= print_object.center_offset();
+                    line_end   -= print_object.center_offset();
+
+                    std::vector<PaintedLine> painted_line_tmp;
+                    PaintedLineVisitor       visitor(edge_grids[layer_idx], painted_line_tmp);
+                    visitor.reset();
+                    visitor.line_to_test.a = line_start;
+                    visitor.line_to_test.b = line_end;
+                    visitor.color          = int(extruder_idx);
+                    edge_grids[layer_idx].visit_cells_intersecting_line(line_start, line_end, visitor);
+
+                    append(painted_lines[layer_idx], std::move(painted_line_tmp));
+                }
+            }
+        }
+    }
+
+    tbb::parallel_for(tbb::blocked_range<size_t>(0, print_object.layers().size()), [&](const tbb::blocked_range<size_t> &range) {
+        for (size_t layer_idx = range.begin(); layer_idx < range.end(); ++layer_idx) {
+            //    for(size_t layer_idx = 0; layer_idx < print_object.layers().size(); ++layer_idx) {
+            BOOST_LOG_TRIVIAL(debug) << "MMU segmentation of layer: " << layer_idx;
+            auto comp = [&edge_grids, layer_idx](const PaintedLine &first, const PaintedLine &second) {
+                Point first_start_p = *(edge_grids[layer_idx].contours()[first.contour_idx].begin() + first.line_idx);
+
+                return first.contour_idx < second.contour_idx ||
+                       (first.contour_idx == second.contour_idx &&
+                        (first.line_idx < second.line_idx ||
+                         (first.line_idx == second.line_idx &&
+                          Line(first_start_p, first.projected_line.a).length() < Line(first_start_p, second.projected_line.a).length())));
+            };
+
+            std::sort(painted_lines[layer_idx].begin(), painted_lines[layer_idx].end(), comp);
+            std::vector<PaintedLine> &painted_lines_single = painted_lines[layer_idx];
+
+            if (!painted_lines_single.empty()) {
+                std::vector<std::vector<ColoredLine>> color_poly = colorize_polygons(input_polygons[layer_idx], painted_lines_single);
+                MMU_Graph                             graph      = build_graph(layer_idx, color_poly);
+                remove_multiple_edges_in_vertices(graph, color_poly);
+                graph.remove_nodes_with_one_arc();
+                std::vector<std::pair<Polygon, size_t>> segmentation = extract_colored_segments(graph);
+                for (const std::pair<Polygon, size_t> &region : segmentation)
+                    segmented_regions[layer_idx].emplace_back(region);
+            }
+        }
+    }); // end of parallel_for
+
+    if (auto w = print_object.print()->config().mmu_segmented_region_max_width; w > 0.f)
+        cut_segmented_layers(input_expolygons, segmented_regions, float(-scale_(w)));
+
+//    return segmented_regions;
+    std::vector<std::vector<ExPolygons>>                   top_and_bottom_layers    = mmu_segmentation_top_and_bottom_layers(print_object, input_expolygons);
+    std::vector<std::vector<std::pair<ExPolygon, size_t>>> segmented_regions_merged = merge_segmented_layers(segmented_regions, std::move(top_and_bottom_layers));
+    return segmented_regions_merged;
+}
+
+} // namespace Slic3r
diff --git a/src/libslic3r/MultiMaterialSegmentation.hpp b/src/libslic3r/MultiMaterialSegmentation.hpp
new file mode 100644
index 000000000..213744f15
--- /dev/null
+++ b/src/libslic3r/MultiMaterialSegmentation.hpp
@@ -0,0 +1,18 @@
+#ifndef slic3r_MultiMaterialSegmentation_hpp_
+#define slic3r_MultiMaterialSegmentation_hpp_
+
+#include <utility>
+#include <vector>
+
+namespace Slic3r {
+
+
+class PrintObject;
+class ExPolygon;
+
+// Returns MMU segmentation based on painting in MMU segmentation gizmo
+std::vector<std::vector<std::pair<ExPolygon, size_t>>> multi_material_segmentation_by_painting(const PrintObject &print_object);
+
+} // namespace Slic3r
+
+#endif // slic3r_MultiMaterialSegmentation_hpp_
diff --git a/src/libslic3r/Preset.cpp b/src/libslic3r/Preset.cpp
index 7a258182e..3be867d53 100644
--- a/src/libslic3r/Preset.cpp
+++ b/src/libslic3r/Preset.cpp
@@ -445,7 +445,7 @@ const std::vector<std::string>& Preset::print_options()
         "perimeter_extrusion_width", "external_perimeter_extrusion_width", "infill_extrusion_width", "solid_infill_extrusion_width",
         "top_infill_extrusion_width", "support_material_extrusion_width", "infill_overlap", "infill_anchor", "infill_anchor_max", "bridge_flow_ratio", "clip_multipart_objects",
         "elefant_foot_compensation", "xy_size_compensation", "threads", "resolution", "wipe_tower", "wipe_tower_x", "wipe_tower_y",
-        "wipe_tower_width", "wipe_tower_rotation_angle", "wipe_tower_brim_width", "wipe_tower_bridging", "single_extruder_multi_material_priming",
+        "wipe_tower_width", "wipe_tower_rotation_angle", "wipe_tower_brim_width", "wipe_tower_bridging", "single_extruder_multi_material_priming", "mmu_segmented_region_max_width",
         "wipe_tower_no_sparse_layers", "compatible_printers", "compatible_printers_condition", "inherits"
     };
     return s_opts;
diff --git a/src/libslic3r/Print.cpp b/src/libslic3r/Print.cpp
index 99119f566..8288ddb19 100644
--- a/src/libslic3r/Print.cpp
+++ b/src/libslic3r/Print.cpp
@@ -216,6 +216,8 @@ bool Print::invalidate_state_by_config_options(const ConfigOptionResolver & /* n
             osteps.emplace_back(posSupportMaterial);
             steps.emplace_back(psSkirt);
             steps.emplace_back(psBrim);
+        } else if (opt_key == "mmu_segmented_region_max_width") {
+            invalidated |= this->invalidate_all_steps();
         } else {
             // for legacy, if we can't handle this option let's invalidate all steps
             //FIXME invalidate all steps of all objects as well?
diff --git a/src/libslic3r/Print.hpp b/src/libslic3r/Print.hpp
index 05d417084..073779ae6 100644
--- a/src/libslic3r/Print.hpp
+++ b/src/libslic3r/Print.hpp
@@ -546,7 +546,6 @@ public:
             [object_id](const PrintObject *obj) { return obj->id() == object_id; });
         return (it == m_objects.end()) ? nullptr : *it;
     }
-//    ConstPrintRegionPtrsAdaptor regions() const { return ConstPrintRegionPtrsAdaptor(&m_regions); }
     // How many of PrintObject::copies() over all print objects are there?
     // If zero, then the print is empty and the print shall not be executed.
     unsigned int                num_object_instances() const;
diff --git a/src/libslic3r/PrintConfig.cpp b/src/libslic3r/PrintConfig.cpp
index 72722e6fc..3fffe7f7d 100644
--- a/src/libslic3r/PrintConfig.cpp
+++ b/src/libslic3r/PrintConfig.cpp
@@ -1433,6 +1433,14 @@ void PrintConfigDef::init_fff_params()
     def->mode = comExpert;
     def->set_default_value(new ConfigOptionBool(false));
 
+    def = this->add("mmu_segmented_region_max_width", coFloat);
+    def->label = L("Maximum width of a segmented region");
+    def->tooltip = L("Maximum width of a segmented region. Zero disables this feature.");
+    def->sidetext = L("mm (zero to disable)");
+    def->min = 0;
+    def->mode = comExpert;
+    def->set_default_value(new ConfigOptionFloat(0.f));
+
     def = this->add("ironing", coBool);
     def->label = L("Enable ironing");
     def->tooltip = L("Enable ironing of the top layers with the hot print head for smooth surface");
diff --git a/src/libslic3r/PrintConfig.hpp b/src/libslic3r/PrintConfig.hpp
index 0d94513d6..e7bc1ded5 100644
--- a/src/libslic3r/PrintConfig.hpp
+++ b/src/libslic3r/PrintConfig.hpp
@@ -699,6 +699,7 @@ PRINT_CONFIG_CLASS_DERIVED_DEFINE(
     ((ConfigOptionFloat,              max_print_height))
     ((ConfigOptionFloats,             min_print_speed))
     ((ConfigOptionFloat,              min_skirt_length))
+    ((ConfigOptionFloat,              mmu_segmented_region_max_width))
     ((ConfigOptionString,             notes))
     ((ConfigOptionFloats,             nozzle_diameter))
     ((ConfigOptionBool,               only_retract_when_crossing_perimeters))
diff --git a/src/libslic3r/PrintObject.cpp b/src/libslic3r/PrintObject.cpp
index 5422908ba..d5808a1f3 100644
--- a/src/libslic3r/PrintObject.cpp
+++ b/src/libslic3r/PrintObject.cpp
@@ -6,6 +6,7 @@
 #include "Geometry.hpp"
 #include "I18N.hpp"
 #include "Layer.hpp"
+#include "MultiMaterialSegmentation.hpp"
 #include "SupportMaterial.hpp"
 #include "Surface.hpp"
 #include "Slicing.hpp"
diff --git a/src/libslic3r/Slicing.cpp b/src/libslic3r/Slicing.cpp
index 82b3cf1b6..c39f07a1f 100644
--- a/src/libslic3r/Slicing.cpp
+++ b/src/libslic3r/Slicing.cpp
@@ -298,7 +298,7 @@ std::vector<double> layer_height_profile_adaptive(const SlicingParameters& slici
     if (z_gap > 0.0)
     {
         layer_height_profile.push_back(slicing_params.object_print_z_height());
-        layer_height_profile.push_back(clamp(slicing_params.min_layer_height, slicing_params.max_layer_height, z_gap));
+        layer_height_profile.push_back(std::clamp(z_gap, slicing_params.min_layer_height, slicing_params.max_layer_height));
     }
 
     return layer_height_profile;
@@ -376,7 +376,7 @@ std::vector<double> smooth_height_profile(const std::vector<double>& profile, co
                 }
             }
 
-            height = clamp(slicing_params.min_layer_height, slicing_params.max_layer_height, (weight_total != 0.0) ? height /= weight_total : hi);
+            height = std::clamp((weight_total != 0.0) ? height /= weight_total : hi, slicing_params.min_layer_height, slicing_params.max_layer_height);
             if (smoothing_params.keep_min)
                 height = std::min(height, hi);
         }
@@ -502,7 +502,7 @@ void adjust_layer_height_profile(
                 assert(false);
                 break;
         }
-        height = clamp(slicing_params.min_layer_height, slicing_params.max_layer_height, height);
+        height = std::clamp(height, slicing_params.min_layer_height, slicing_params.max_layer_height);
         if (zz == z_span_variable.second) {
             // This is the last point of the profile.
             if (profile_new[profile_new.size() - 2] + EPSILON > zz) {
@@ -670,11 +670,11 @@ int generate_layer_height_texture(
 		assert(mid <= slicing_params.object_print_z_height());
 		coordf_t h = hi - lo;
 		hi = std::min(hi, slicing_params.object_print_z_height());
-        int cell_first = clamp(0, ncells-1, int(ceil(lo * z_to_cell)));
-        int cell_last  = clamp(0, ncells-1, int(floor(hi * z_to_cell)));
+        int cell_first = std::clamp(int(ceil(lo * z_to_cell)), 0, ncells-1);
+        int cell_last  = std::clamp(int(floor(hi * z_to_cell)), 0, ncells-1);
         for (int cell = cell_first; cell <= cell_last; ++ cell) {
             coordf_t idxf = (0.5 * hscale + (h - slicing_params.layer_height)) * coordf_t(palette_raw.size()-1) / hscale;
-            int idx1 = clamp(0, int(palette_raw.size() - 1), int(floor(idxf)));
+            int idx1 = std::clamp(int(floor(idxf)), 0, int(palette_raw.size() - 1));
             int idx2 = std::min(int(palette_raw.size() - 1), idx1 + 1);
 			coordf_t t = idxf - coordf_t(idx1);
             const Vec3crd &color1 = palette_raw[idx1];
@@ -693,9 +693,9 @@ int generate_layer_height_texture(
 			assert(row >= 0 && row < rows);
 			assert(col >= 0 && col < cols);
             unsigned char *ptr = (unsigned char*)data + (row * cols + col) * 4;
-            ptr[0] = (unsigned char)clamp<int>(0, 255, int(floor(color(0) + 0.5)));
-            ptr[1] = (unsigned char)clamp<int>(0, 255, int(floor(color(1) + 0.5)));
-            ptr[2] = (unsigned char)clamp<int>(0, 255, int(floor(color(2) + 0.5)));
+            ptr[0] = (unsigned char)std::clamp(int(floor(color(0) + 0.5)), 0, 255);
+            ptr[1] = (unsigned char)std::clamp(int(floor(color(1) + 0.5)), 0, 255);
+            ptr[2] = (unsigned char)std::clamp(int(floor(color(2) + 0.5)), 0, 255);
             ptr[3] = 255;
             if (col == 0 && row > 0) {
                 // Duplicate the first value in a row as a last value of the preceding row.
@@ -706,11 +706,11 @@ int generate_layer_height_texture(
             }
         }
         if (level_of_detail_2nd_level) {
-            cell_first = clamp(0, ncells1-1, int(ceil(lo * z_to_cell1))); 
-            cell_last  = clamp(0, ncells1-1, int(floor(hi * z_to_cell1)));
+            cell_first = std::clamp(int(ceil(lo * z_to_cell1)), 0, ncells1-1);
+            cell_last  = std::clamp(int(floor(hi * z_to_cell1)), 0, ncells1-1);
             for (int cell = cell_first; cell <= cell_last; ++ cell) {
                 coordf_t idxf = (0.5 * hscale + (h - slicing_params.layer_height)) * coordf_t(palette_raw.size()-1) / hscale;
-                int idx1 = clamp(0, int(palette_raw.size() - 1), int(floor(idxf)));
+                int idx1 = std::clamp(int(floor(idxf)), 0, int(palette_raw.size() - 1));
                 int idx2 = std::min(int(palette_raw.size() - 1), idx1 + 1);
     			coordf_t t = idxf - coordf_t(idx1);
                 const Vec3crd &color1 = palette_raw[idx1];
@@ -725,9 +725,9 @@ int generate_layer_height_texture(
     			assert(row >= 0 && row < rows/2);
     			assert(col >= 0 && col < cols/2);
                 unsigned char *ptr = data1 + (row * cols1 + col) * 4;
-                ptr[0] = (unsigned char)clamp<int>(0, 255, int(floor(color(0) + 0.5)));
-                ptr[1] = (unsigned char)clamp<int>(0, 255, int(floor(color(1) + 0.5)));
-                ptr[2] = (unsigned char)clamp<int>(0, 255, int(floor(color(2) + 0.5)));
+                ptr[0] = (unsigned char)std::clamp(int(floor(color(0) + 0.5)), 0, 255);
+                ptr[1] = (unsigned char)std::clamp(int(floor(color(1) + 0.5)), 0, 255);
+                ptr[2] = (unsigned char)std::clamp(int(floor(color(2) + 0.5)), 0, 255);
                 ptr[3] = 255;
                 if (col == 0 && row > 0) {
                     // Duplicate the first value in a row as a last value of the preceding row.
diff --git a/src/libslic3r/TriangleSelector.cpp b/src/libslic3r/TriangleSelector.cpp
index 7570075e5..95f2984e3 100644
--- a/src/libslic3r/TriangleSelector.cpp
+++ b/src/libslic3r/TriangleSelector.cpp
@@ -36,7 +36,7 @@ void TriangleSelector::Triangle::set_division(int sides_to_split, int special_si
 void TriangleSelector::select_patch(const Vec3f& hit, int facet_start,
                                     const Vec3f& source, float radius,
                                     CursorType cursor_type, EnforcerBlockerType new_state,
-                                    const Transform3d& trafo)
+                                    const Transform3d& trafo, bool triangle_splitting)
 {
     assert(facet_start < m_orig_size_indices);
 
@@ -59,7 +59,7 @@ void TriangleSelector::select_patch(const Vec3f& hit, int facet_start,
     while (facet_idx < int(facets_to_check.size())) {
         int facet = facets_to_check[facet_idx];
         if (! visited[facet]) {
-            if (select_triangle(facet, new_state)) {
+            if (select_triangle(facet, new_state, false, triangle_splitting)) {
                 // add neighboring facets to list to be proccessed later
                 for (int n=0; n<3; ++n) {
                     int neighbor_idx = m_mesh->stl.neighbors_start[facet].neighbor[n];
@@ -73,13 +73,56 @@ void TriangleSelector::select_patch(const Vec3f& hit, int facet_start,
     }
 }
 
+void TriangleSelector::seed_fill_select_triangles(const Vec3f& hit, int facet_start, float seed_fill_angle)
+{
+    this->seed_fill_unselect_all_triangles();
+
+    std::vector<bool> visited(m_triangles.size(), false);
+    std::queue<size_t> facet_queue;
+    facet_queue.push(facet_start);
+
+    // Check if neighbour_facet_idx is satisfies angle in seed_fill_angle and append it to facet_queue if it do.
+    auto check_angle_and_append = [this, &facet_queue](const size_t facet_idx, const size_t neighbour_facet_idx, const float seed_fill_angle) -> void {
+        double dot_product       = m_triangles[neighbour_facet_idx].normal.dot(m_triangles[facet_idx].normal);
+        dot_product              = std::clamp(dot_product, 0., 1.);
+        double facet_angle_limit = cos(Geometry::deg2rad(seed_fill_angle));
+        if ((dot_product + EPSILON) >= facet_angle_limit)
+            facet_queue.push(neighbour_facet_idx);
+    };
+
+    while(!facet_queue.empty()) {
+        size_t current_facet = facet_queue.front();
+        facet_queue.pop();
 
+        if (!visited[current_facet]) {
+            if (!m_triangles[current_facet].is_split())
+                m_triangles[current_facet].select_by_seed_fill();
+
+            if (m_triangles[current_facet].is_split())
+                for (int split_triangle_idx = 0; split_triangle_idx <= m_triangles[current_facet].number_of_split_sides(); ++split_triangle_idx) {
+                    assert(split_triangle_idx < int(m_triangles[current_facet].children.size()));
+                    assert(m_triangles[current_facet].children[split_triangle_idx] < int(m_triangles.size()));
+
+                    if (!visited[m_triangles[current_facet].children[split_triangle_idx]])
+                        check_angle_and_append(current_facet, m_triangles[current_facet].children[split_triangle_idx], seed_fill_angle);
+                }
+
+            if (int(current_facet) < m_orig_size_indices)
+                for (int neighbor_idx : m_mesh->stl.neighbors_start[current_facet].neighbor) {
+                    assert(neighbor_idx >= 0);
+                    if (neighbor_idx >= 0 && !visited[neighbor_idx])
+                        check_angle_and_append(current_facet, neighbor_idx, seed_fill_angle);
+                }
+        }
+        visited[current_facet] = true;
+    }
+}
 
 // Selects either the whole triangle (discarding any children it had), or divides
 // the triangle recursively, selecting just subtriangles truly inside the circle.
 // This is done by an actual recursive call. Returns false if the triangle is
 // outside the cursor.
-bool TriangleSelector::select_triangle(int facet_idx, EnforcerBlockerType type, bool recursive_call)
+bool TriangleSelector::select_triangle(int facet_idx, EnforcerBlockerType type, bool recursive_call, bool triangle_splitting)
 {
     assert(facet_idx < int(m_triangles.size()));
 
@@ -108,7 +151,10 @@ bool TriangleSelector::select_triangle(int facet_idx, EnforcerBlockerType type,
             return true;
         }
 
-        split_triangle(facet_idx);
+        if(triangle_splitting)
+            split_triangle(facet_idx);
+        else if(!m_triangles[facet_idx].is_split())
+            m_triangles[facet_idx].set_state(type);
         tr = &m_triangles[facet_idx]; // might have been invalidated
 
 
@@ -118,7 +164,7 @@ bool TriangleSelector::select_triangle(int facet_idx, EnforcerBlockerType type,
                 assert(i < int(tr->children.size()));
                 assert(tr->children[i] < int(m_triangles.size()));
 
-                select_triangle(tr->children[i], type, true);
+                select_triangle(tr->children[i], type, true, triangle_splitting);
                 tr = &m_triangles[facet_idx]; // might have been invalidated
             }
         }
@@ -417,7 +463,7 @@ TriangleSelector::TriangleSelector(const TriangleMesh& mesh)
 }
 
 
-void TriangleSelector::reset()
+void TriangleSelector::reset(const EnforcerBlockerType reset_state)
 {
     if (m_orig_size_indices != 0) // unless this is run from constructor
         garbage_collect();
@@ -428,7 +474,7 @@ void TriangleSelector::reset()
     for (size_t i=0; i<m_mesh->its.indices.size(); ++i) {
         const stl_triangle_vertex_indices& ind = m_mesh->its.indices[i];
         const Vec3f& normal = m_mesh->stl.facet_start[i].normal;
-        push_triangle(ind[0], ind[1], ind[2], normal);
+        push_triangle(ind[0], ind[1], ind[2], normal, reset_state);
     }
     m_orig_size_vertices = m_vertices.size();
     m_orig_size_indices = m_triangles.size();
@@ -454,13 +500,13 @@ void TriangleSelector::set_edge_limit(float edge_limit)
 
 
 
-void TriangleSelector::push_triangle(int a, int b, int c, const Vec3f& normal)
+void TriangleSelector::push_triangle(int a, int b, int c, const Vec3f& normal, const EnforcerBlockerType state)
 {
     for (int i : {a, b, c}) {
         assert(i >= 0 && i < int(m_vertices.size()));
         ++m_vertices[i].ref_cnt;
     }
-    m_triangles.emplace_back(a, b, c, normal);
+    m_triangles.emplace_back(a, b, c, normal, state);
 }
 
 
@@ -550,8 +596,9 @@ indexed_triangle_set TriangleSelector::get_facets(EnforcerBlockerType state) con
 std::map<int, std::vector<bool>> TriangleSelector::serialize() const
 {
     // Each original triangle of the mesh is assigned a number encoding its state
-    // or how it is split. Each triangle is encoded by 4 bits (xxyy):
-    // leaf triangle: xx = EnforcerBlockerType, yy = 0
+    // or how it is split. Each triangle is encoded by 4 bits (xxyy) or 8 bits (zzzzxxyy):
+    // leaf triangle: xx = EnforcerBlockerType (Only values 0, 1, and 2. Value 3 is used as an indicator for additional 4 bits.), yy = 0
+    // leaf triangle: xx = 0b11, yy = 0b00, zzzz = EnforcerBlockerType (subtracted by 3)
     // non-leaf:      xx = special side, yy = number of split sides
     // These are bitwise appended and formed into one 64-bit integer.
 
@@ -594,9 +641,17 @@ std::map<int, std::vector<bool>> TriangleSelector::serialize() const
                     serialize_recursive(tr.children[child_idx]);
             } else {
                 // In case this is leaf, we better save information about its state.
-                assert(int(tr.get_state()) <= 3);
-                data.push_back(int(tr.get_state()) & 0b01);
-                data.push_back(int(tr.get_state()) & 0b10);
+                assert(int(tr.get_state()) <= 15);
+                if (3 <= int(tr.get_state()) && int(tr.get_state()) <= 15) {
+                    data.insert(data.end(), {true, true});
+                    for (size_t bit_idx = 0; bit_idx < 4; ++bit_idx) {
+                        size_t bit_mask = uint64_t(0b0001) << bit_idx;
+                        data.push_back((int(tr.get_state()) - 3) & bit_mask);
+                    }
+                } else {
+                    data.push_back(int(tr.get_state()) & 0b01);
+                    data.push_back(int(tr.get_state()) & 0b10);
+                }
                 ++stored_triangles;
             }
         };
@@ -608,13 +663,13 @@ std::map<int, std::vector<bool>> TriangleSelector::serialize() const
     return out;
 }
 
-void TriangleSelector::deserialize(const std::map<int, std::vector<bool>> data)
+void TriangleSelector::deserialize(const std::map<int, std::vector<bool>> data, const EnforcerBlockerType init_state)
 {
-    reset(); // dump any current state
+    reset(init_state); // dump any current state
     for (const auto& [triangle_id, code] : data) {
         assert(triangle_id < int(m_triangles.size()));
         assert(! code.empty());
-        int processed_triangles = 0;
+        int processed_nibbles = 0;
         struct ProcessingInfo {
             int facet_id = 0;
             int processed_children = 0;
@@ -626,18 +681,26 @@ void TriangleSelector::deserialize(const std::map<int, std::vector<bool>> data)
 
         while (true) {
             // Read next triangle info.
-            int next_code = 0;
-            for (int i=3; i>=0; --i) {
-                next_code = next_code << 1;
-                next_code |= int(code[4 * processed_triangles + i]);
+            std::array<int, 2> next_code{};
+            for(size_t nibble_idx = 0; nibble_idx < 2; ++nibble_idx) {
+                assert(nibble_idx < 2);
+                if(nibble_idx >= 1 && (next_code[0] >> 2) != 0b11)
+                    break;
+
+                for (int i = 3; i >= 0; --i) {
+                    next_code[nibble_idx] = next_code[nibble_idx] << 1;
+                    next_code[nibble_idx] |= int(code[4 * processed_nibbles + i]);
+                }
+
+                ++processed_nibbles;
             }
-            ++processed_triangles;
 
-            int num_of_split_sides = (next_code & 0b11);
+            int num_of_split_sides = (next_code[0] & 0b11);
             int num_of_children = num_of_split_sides != 0 ? num_of_split_sides + 1 : 0;
             bool is_split = num_of_children != 0;
-            EnforcerBlockerType state = EnforcerBlockerType(next_code >> 2);
-            int special_side = (next_code >> 2);
+            // Value of the second nibble was subtracted by 3, so it is added back.
+            auto state = EnforcerBlockerType(next_code[0] >> 2 == 0b11 ? next_code[1] + 3 : next_code[0] >> 2);
+            int special_side = (next_code[0] >> 2);
 
             // Take care of the first iteration separately, so handling of the others is simpler.
             if (parents.empty()) {
@@ -693,6 +756,18 @@ void TriangleSelector::deserialize(const std::map<int, std::vector<bool>> data)
     }
 }
 
+void TriangleSelector::seed_fill_unselect_all_triangles() {
+    for (Triangle &triangle : m_triangles)
+        if (!triangle.is_split())
+            triangle.unselect_by_seed_fill();
+}
+
+void TriangleSelector::seed_fill_apply_on_triangles(EnforcerBlockerType new_state)
+{
+    for (Triangle &triangle : m_triangles)
+        if (!triangle.is_split() && triangle.is_selected_by_seed_fill())
+            triangle.set_state(new_state);
+}
 
 TriangleSelector::Cursor::Cursor(
         const Vec3f& center_, const Vec3f& source_, float radius_world,
diff --git a/src/libslic3r/TriangleSelector.hpp b/src/libslic3r/TriangleSelector.hpp
index 9d1590070..383aab480 100644
--- a/src/libslic3r/TriangleSelector.hpp
+++ b/src/libslic3r/TriangleSelector.hpp
@@ -29,13 +29,18 @@ public:
     explicit TriangleSelector(const TriangleMesh& mesh);
 
     // Select all triangles fully inside the circle, subdivide where needed.
-    void select_patch(const Vec3f& hit,    // point where to start
-                      int facet_start,     // facet that point belongs to
-                      const Vec3f& source, // camera position (mesh coords)
-                      float radius,        // radius of the cursor
-                      CursorType type,     // current type of cursor
+    void select_patch(const Vec3f        &hit,         // point where to start
+                      int                 facet_start, // facet that point belongs to
+                      const Vec3f        &source,      // camera position (mesh coords)
+                      float               radius,      // radius of the cursor
+                      CursorType          type,        // current type of cursor
                       EnforcerBlockerType new_state,   // enforcer or blocker?
-                      const Transform3d& trafo); // matrix to get from mesh to world
+                      const Transform3d  &trafo,       // matrix to get from mesh to world
+                      bool                triangle_splitting); // If triangles will be split base on the cursor or not
+
+    void seed_fill_select_triangles(const Vec3f &hit,               // point where to start
+                                    int          facet_start,       // facet that point belongs to
+                                    float        seed_fill_angle);  // the maximal angle between two facets to be painted by the same color
 
     // Get facets currently in the given state.
     indexed_triangle_set get_facets(EnforcerBlockerType state) const;
@@ -44,7 +49,7 @@ public:
     void set_facet(int facet_idx, EnforcerBlockerType state);
 
     // Clear everything and make the tree empty.
-    void reset();
+    void reset(const EnforcerBlockerType reset_state = EnforcerBlockerType{0});
 
     // Remove all unnecessary data.
     void garbage_collect();
@@ -54,8 +59,13 @@ public:
     std::map<int, std::vector<bool>> serialize() const;
 
     // Load serialized data. Assumes that correct mesh is loaded.
-    void deserialize(const std::map<int, std::vector<bool>> data);
+    void deserialize(const std::map<int, std::vector<bool>> data, const EnforcerBlockerType init_state = EnforcerBlockerType{0});
+
+    // For all triangles, remove the flag indicating that the triangle was selected by seed fill.
+    void seed_fill_unselect_all_triangles();
 
+    // For all triangles selected by seed fill, set new EnforcerBlockerType and remove flag indicating that triangle was selected by seed fill.
+    void seed_fill_apply_on_triangles(EnforcerBlockerType new_state);
 
 protected:
     // Triangle and info about how it's split.
@@ -63,10 +73,10 @@ protected:
     public:
         // Use TriangleSelector::push_triangle to create a new triangle.
         // It increments/decrements reference counter on vertices.
-        Triangle(int a, int b, int c, const Vec3f& normal_)
+        Triangle(int a, int b, int c, const Vec3f& normal_, const EnforcerBlockerType init_state)
             : verts_idxs{a, b, c},
               normal{normal_},
-              state{EnforcerBlockerType(0)},
+              state{init_state},
               number_of_splits{0},
               special_side_idx{0},
               old_number_of_splits{0}
@@ -90,6 +100,12 @@ protected:
         void set_state(EnforcerBlockerType type) { assert(! is_split()); state = type; }
         EnforcerBlockerType get_state() const { assert(! is_split()); return state; }
 
+        // Set if the triangle has been selected or unselected by seed fill.
+        void select_by_seed_fill() { assert(! is_split()); m_selected_by_seed_fill = true; }
+        void unselect_by_seed_fill() { assert(! is_split()); m_selected_by_seed_fill = false; }
+        // Get if the triangle has been selected or not by seed fill.
+        bool is_selected_by_seed_fill() const { assert(! is_split()); return m_selected_by_seed_fill; }
+
         // Get info on how it's split.
         bool is_split() const { return number_of_split_sides() != 0; }
         int number_of_split_sides() const { return number_of_splits; }
@@ -101,6 +117,7 @@ protected:
         int number_of_splits;
         int special_side_idx;
         EnforcerBlockerType state;
+        bool m_selected_by_seed_fill = false;
 
         // How many children were spawned during last split?
         // Is not reset on remerging the triangle.
@@ -153,8 +170,7 @@ protected:
     float m_old_cursor_radius_sqr;
 
     // Private functions:
-    bool select_triangle(int facet_idx, EnforcerBlockerType type,
-                         bool recursive_call = false);
+    bool select_triangle(int facet_idx, EnforcerBlockerType type, bool recursive_call = false, bool triangle_splitting = true);
     int  vertices_inside(int facet_idx) const;
     bool faces_camera(int facet) const;
     void undivide_triangle(int facet_idx);
@@ -162,7 +178,7 @@ protected:
     void remove_useless_children(int facet_idx); // No hidden meaning. Triangles are meant.
     bool is_pointer_in_triangle(int facet_idx) const;
     bool is_edge_inside_cursor(int facet_idx) const;
-    void push_triangle(int a, int b, int c, const Vec3f& normal);
+    void push_triangle(int a, int b, int c, const Vec3f &normal, const EnforcerBlockerType state = EnforcerBlockerType{0});
     void perform_split(int facet_idx, EnforcerBlockerType old_state);
 };
 
diff --git a/src/libslic3r/VoronoiOffset.cpp b/src/libslic3r/VoronoiOffset.cpp
index 2108388f5..e8d13a6ad 100644
--- a/src/libslic3r/VoronoiOffset.cpp
+++ b/src/libslic3r/VoronoiOffset.cpp
@@ -41,7 +41,7 @@ namespace detail {
 			// Degenerate to a single closest point.
 			t = - b / (2. * a);
 			assert(t >= - EPSILON && t <= 1. + EPSILON);
-			return Slic3r::clamp(0., 1., t);
+			return std::clamp(t, 0., 1.);
 		} else {
 			u = sqrt(u);
 			out.first = 2;
@@ -1142,7 +1142,7 @@ std::vector<Vec2d> edge_offset_contour_intersections(
 #endif // NDEBUG
                 if (! bisector || (dmin != dmax && offset_distance >= dmin)) {
                     double t = (offset_distance - dmin) / (dmax - dmin);
-                    t = clamp(0., 1., t);
+                    t = std::clamp(t, 0., 1.);
                     if (d1 < d0) {
                         out[edge_idx2] = Slic3r::lerp(vertex_point(v1), vertex_point(v0), t);
                         // mark visited
diff --git a/src/libslic3r/libslic3r.h b/src/libslic3r/libslic3r.h
index 2b0c94161..d2c1d0fc1 100644
--- a/src/libslic3r/libslic3r.h
+++ b/src/libslic3r/libslic3r.h
@@ -239,26 +239,20 @@ template<typename T> inline bool one_of(const T& v, const std::initializer_list<
     { return contains(il, v); }
 
 template<typename T>
-static inline T sqr(T x)
+constexpr inline T sqr(T x)
 {
     return x * x;
 }
 
-template <typename T>
-static inline T clamp(const T low, const T high, const T value)
-{
-    return std::max(low, std::min(high, value));
-}
-
 template <typename T, typename Number>
-static inline T lerp(const T& a, const T& b, Number t)
+constexpr inline T lerp(const T& a, const T& b, Number t)
 {
     assert((t >= Number(-EPSILON)) && (t <= Number(1) + Number(EPSILON)));
     return (Number(1) - t) * a + t * b;
 }
 
 template <typename Number>
-static inline bool is_approx(Number value, Number test_value)
+constexpr inline bool is_approx(Number value, Number test_value)
 {
     return std::fabs(double(value) - double(test_value)) < double(EPSILON);
 }
diff --git a/src/slic3r/CMakeLists.txt b/src/slic3r/CMakeLists.txt
index aa3bf01c4..333f15113 100644
--- a/src/slic3r/CMakeLists.txt
+++ b/src/slic3r/CMakeLists.txt
@@ -57,6 +57,8 @@ set(SLIC3R_GUI_SOURCES
     GUI/Gizmos/GLGizmoPainterBase.hpp
     GUI/Gizmos/GLGizmoSeam.cpp
     GUI/Gizmos/GLGizmoSeam.hpp
+    GUI/Gizmos/GLGizmoMmuSegmentation.cpp
+    GUI/Gizmos/GLGizmoMmuSegmentation.hpp
     GUI/GLSelectionRectangle.cpp
     GUI/GLSelectionRectangle.hpp
     GUI/GLModel.hpp
diff --git a/src/slic3r/GUI/GLCanvas3D.cpp b/src/slic3r/GUI/GLCanvas3D.cpp
index 455eb8a78..6f4650e98 100644
--- a/src/slic3r/GUI/GLCanvas3D.cpp
+++ b/src/slic3r/GUI/GLCanvas3D.cpp
@@ -259,7 +259,7 @@ void GLCanvas3D::LayersEditing::render_overlay(const GLCanvas3D& canvas) const
     ImGui::SameLine();
     float widget_align = ImGui::GetCursorPosX();
     ImGui::PushItemWidth(imgui.get_style_scaling() * 120.0f);
-    m_adaptive_quality = clamp(0.0f, 1.f, m_adaptive_quality);
+    m_adaptive_quality = std::clamp(m_adaptive_quality, 0.0f, 1.f);
     ImGui::SliderFloat("", &m_adaptive_quality, 0.0f, 1.f, "%.2f");
 
     ImGui::Separator();
@@ -1062,7 +1062,8 @@ void GLCanvas3D::toggle_model_objects_visibility(bool visible, const ModelObject
                     const GLGizmosManager& gm = get_gizmos_manager();
                     auto gizmo_type = gm.get_current_type();
                     if (    (gizmo_type == GLGizmosManager::FdmSupports
-                          || gizmo_type == GLGizmosManager::Seam)
+                          || gizmo_type == GLGizmosManager::Seam
+                          || gizmo_type == GLGizmosManager::MmuSegmentation)
                         && ! vol->is_modifier)
                         vol->force_neutral_color = true;
                     else
@@ -2926,7 +2927,8 @@ void GLCanvas3D::on_mouse(wxMouseEvent& evt)
         else if (evt.LeftDown() && (evt.ShiftDown() || evt.AltDown()) && m_picking_enabled) {
             if (m_gizmos.get_current_type() != GLGizmosManager::SlaSupports
              && m_gizmos.get_current_type() != GLGizmosManager::FdmSupports
-             && m_gizmos.get_current_type() != GLGizmosManager::Seam) {
+             && m_gizmos.get_current_type() != GLGizmosManager::Seam
+             && m_gizmos.get_current_type() != GLGizmosManager::MmuSegmentation) {
                 m_rectangle_selection.start_dragging(m_mouse.position, evt.ShiftDown() ? GLSelectionRectangle::Select : GLSelectionRectangle::Deselect);
                 m_dirty = true;
             }
@@ -4808,7 +4810,8 @@ void GLCanvas3D::_render_bed(bool bottom, bool show_axes) const
             (m_gizmos.get_current_type() != GLGizmosManager::FdmSupports
           && m_gizmos.get_current_type() != GLGizmosManager::SlaSupports
           && m_gizmos.get_current_type() != GLGizmosManager::Hollow
-          && m_gizmos.get_current_type() != GLGizmosManager::Seam);
+          && m_gizmos.get_current_type() != GLGizmosManager::Seam
+          && m_gizmos.get_current_type() != GLGizmosManager::MmuSegmentation);
 
     wxGetApp().plater()->get_bed().render(const_cast<GLCanvas3D&>(*this), bottom, scale_factor, show_axes, show_texture);
 }
@@ -4868,7 +4871,8 @@ void GLCanvas3D::_render_objects() const
             const GLGizmosManager& gm = get_gizmos_manager();
             GLGizmosManager::EType type = gm.get_current_type();
             if (type == GLGizmosManager::FdmSupports
-             || type == GLGizmosManager::Seam) {
+                || type == GLGizmosManager::Seam
+                || type == GLGizmosManager::MmuSegmentation) {
                 shader->stop_using();
                 gm.render_painter_gizmo();
                 shader->start_using();
diff --git a/src/slic3r/GUI/Gizmos/GLGizmoMmuSegmentation.cpp b/src/slic3r/GUI/Gizmos/GLGizmoMmuSegmentation.cpp
new file mode 100644
index 000000000..e5acda8c6
--- /dev/null
+++ b/src/slic3r/GUI/Gizmos/GLGizmoMmuSegmentation.cpp
@@ -0,0 +1,671 @@
+#include "GLGizmoMmuSegmentation.hpp"
+
+#include "slic3r/GUI/GLCanvas3D.hpp"
+#include "slic3r/GUI/GUI_App.hpp"
+#include "slic3r/GUI/ImGuiWrapper.hpp"
+#include "slic3r/GUI/Camera.hpp"
+#include "slic3r/GUI/Plater.hpp"
+#include "slic3r/GUI/BitmapCache.hpp"
+#include "libslic3r/PresetBundle.hpp"
+#include "libslic3r/Model.hpp"
+
+
+#include <GL/glew.h>
+
+namespace Slic3r::GUI {
+
+void GLGizmoMmuSegmentation::on_shutdown()
+{
+    m_parent.use_slope(false);
+}
+
+std::string GLGizmoMmuSegmentation::on_get_name() const
+{
+    // FIXME Lukas H.: Discuss and change shortcut
+    return (_L("MMU painting") + " [N]").ToUTF8().data();
+}
+
+bool GLGizmoMmuSegmentation::on_is_selectable() const
+{
+    return (wxGetApp().preset_bundle->printers.get_edited_preset().printer_technology() == ptFFF
+            && wxGetApp().get_mode() != comSimple && wxGetApp().extruders_edited_cnt() > 1);
+}
+
+static std::vector<std::array<uint8_t, 3>> get_extruders_colors()
+{
+    unsigned char                       rgb_color[3] = {};
+    std::vector<std::string>            colors       = Slic3r::GUI::wxGetApp().plater()->get_extruder_colors_from_plater_config();
+    std::vector<std::array<uint8_t, 3>> colors_out(colors.size());
+    for (const std::string &color : colors) {
+        Slic3r::GUI::BitmapCache::parse_color(color, rgb_color);
+        size_t color_idx      = &color - &colors.front();
+        colors_out[color_idx] = {rgb_color[0], rgb_color[1], rgb_color[2]};
+    }
+
+    return colors_out;
+}
+
+static std::vector<std::string> get_extruders_names()
+{
+    size_t                   extruders_count = wxGetApp().extruders_edited_cnt();
+    std::vector<std::string> extruders_out;
+    extruders_out.reserve(extruders_count);
+    for (size_t extruder_idx = 1; extruder_idx <= extruders_count; ++extruder_idx)
+        extruders_out.emplace_back("Extruder " + std::to_string(extruder_idx));
+
+    return extruders_out;
+}
+
+void GLGizmoMmuSegmentation::init_extruders_data()
+{
+    m_original_extruders_names     = get_extruders_names();
+    m_original_extruders_colors    = get_extruders_colors();
+    m_modified_extruders_colors    = m_original_extruders_colors;
+    m_first_selected_extruder_idx  = 0;
+    m_second_selected_extruder_idx = 1;
+}
+
+bool GLGizmoMmuSegmentation::on_init()
+{
+    // FIXME Lukas H.: Discuss and change shortcut
+    m_shortcut_key = WXK_CONTROL_N;
+
+    m_desc["reset_direction"]      = _L("Reset direction");
+    m_desc["clipping_of_view"]     = _L("Clipping of view") + ": ";
+    m_desc["cursor_size"]          = _L("Brush size") + ": ";
+    m_desc["cursor_type"]          = _L("Brush shape") + ": ";
+    m_desc["first_color_caption"]  = _L("Left mouse button") + ": ";
+    m_desc["first_color"]          = _L("First color");
+    m_desc["second_color_caption"] = _L("Right mouse button") + ": ";
+    m_desc["second_color"]         = _L("Second color");
+    m_desc["remove_caption"]       = _L("Shift + Left mouse button") + ": ";
+    m_desc["remove"]               = _L("Remove painted color");
+    m_desc["remove_all"]           = _L("Remove all painted colors");
+    m_desc["circle"]               = _L("Circle");
+    m_desc["sphere"]               = _L("Sphere");
+    m_desc["seed_fill_angle"]      = _L("Seed fill angle");
+
+    init_extruders_data();
+
+    return true;
+}
+
+void GLGizmoMmuSegmentation::render_painter_gizmo() const
+{
+    const Selection& selection = m_parent.get_selection();
+
+    glsafe(::glEnable(GL_BLEND));
+    glsafe(::glEnable(GL_DEPTH_TEST));
+
+    render_triangles(selection);
+
+    m_c->object_clipper()->render_cut();
+    render_cursor();
+
+    glsafe(::glDisable(GL_BLEND));
+}
+
+bool GLGizmoMmuSegmentation::gizmo_event(SLAGizmoEventType action, const Vec2d& mouse_position, bool shift_down, bool alt_down, bool control_down)
+{
+    if (action == SLAGizmoEventType::MouseWheelUp
+        || action == SLAGizmoEventType::MouseWheelDown) {
+        if (control_down) {
+            double pos = m_c->object_clipper()->get_position();
+            pos = action == SLAGizmoEventType::MouseWheelDown
+                  ? std::max(0., pos - 0.01)
+                  : std::min(1., pos + 0.01);
+            m_c->object_clipper()->set_position(pos, true);
+            return true;
+        }
+        else if (alt_down) {
+            m_cursor_radius = action == SLAGizmoEventType::MouseWheelDown
+                              ? std::max(m_cursor_radius - CursorRadiusStep, CursorRadiusMin)
+                              : std::min(m_cursor_radius + CursorRadiusStep, CursorRadiusMax);
+            m_parent.set_as_dirty();
+            return true;
+        }
+    }
+
+    if (action == SLAGizmoEventType::ResetClippingPlane) {
+        m_c->object_clipper()->set_position(-1., false);
+        return true;
+    }
+
+    if (action == SLAGizmoEventType::LeftDown
+        || action == SLAGizmoEventType::RightDown
+        || (action == SLAGizmoEventType::Dragging && m_button_down != Button::None)) {
+
+        if (m_triangle_selectors.empty())
+            return false;
+
+        EnforcerBlockerType new_state = EnforcerBlockerType::NONE;
+        if (! shift_down) {
+            if (action == SLAGizmoEventType::Dragging)
+                new_state = m_button_down == Button::Left
+                            ? EnforcerBlockerType(m_first_selected_extruder_idx)
+                            : EnforcerBlockerType(m_second_selected_extruder_idx);
+            else
+                new_state = action == SLAGizmoEventType::LeftDown
+                            ? EnforcerBlockerType(m_first_selected_extruder_idx)
+                            : EnforcerBlockerType(m_second_selected_extruder_idx);
+        }
+
+        const Camera& camera = wxGetApp().plater()->get_camera();
+        const Selection& selection = m_parent.get_selection();
+        const ModelObject* mo = m_c->selection_info()->model_object();
+        const ModelInstance* mi = mo->instances[selection.get_instance_idx()];
+        const Transform3d& instance_trafo = mi->get_transformation().get_matrix();
+
+        // List of mouse positions that will be used as seeds for painting.
+        std::vector<Vec2d> mouse_positions{mouse_position};
+
+        // In case current mouse position is far from the last one,
+        // add several positions from between into the list, so there
+        // are no gaps in the painted region.
+        {
+            if (m_last_mouse_click == Vec2d::Zero())
+                m_last_mouse_click = mouse_position;
+            // resolution describes minimal distance limit using circle radius
+            // as a unit (e.g., 2 would mean the patches will be touching).
+            double resolution = 0.7;
+            double diameter_px =  resolution  * m_cursor_radius * camera.get_zoom();
+            int patches_in_between = int(((mouse_position - m_last_mouse_click).norm() - diameter_px) / diameter_px);
+            if (patches_in_between > 0) {
+                Vec2d diff = (mouse_position - m_last_mouse_click)/(patches_in_between+1);
+                for (int i=1; i<=patches_in_between; ++i)
+                    mouse_positions.emplace_back(m_last_mouse_click + i*diff);
+            }
+        }
+        m_last_mouse_click = Vec2d::Zero(); // only actual hits should be saved
+
+        // Precalculate transformations of individual meshes.
+        std::vector<Transform3d> trafo_matrices;
+        for (const ModelVolume* mv : mo->volumes) {
+            if (mv->is_model_part())
+                trafo_matrices.emplace_back(instance_trafo * mv->get_matrix());
+        }
+
+        // Now "click" into all the prepared points and spill paint around them.
+        for (const Vec2d& mp : mouse_positions) {
+            update_raycast_cache(mp, camera, trafo_matrices);
+
+            bool dragging_while_painting = (action == SLAGizmoEventType::Dragging && m_button_down != Button::None);
+
+            // The mouse button click detection is enabled when there is a valid hit.
+            // Missing the object entirely
+            // shall not capture the mouse.
+            if (m_rr.mesh_id != -1) {
+                if (m_button_down == Button::None)
+                    m_button_down = ((action == SLAGizmoEventType::LeftDown) ? Button::Left : Button::Right);
+            }
+
+            if (m_rr.mesh_id == -1) {
+                // In case we have no valid hit, we can return. The event will be stopped when
+                // dragging while painting (to prevent scene rotations and moving the object)
+                return dragging_while_painting;
+            }
+
+            const Transform3d& trafo_matrix = trafo_matrices[m_rr.mesh_id];
+
+            // Calculate direction from camera to the hit (in mesh coords):
+            Vec3f camera_pos = (trafo_matrix.inverse() * camera.get_position()).cast<float>();
+
+            assert(m_rr.mesh_id < int(m_triangle_selectors.size()));
+            if (m_seed_fill_enabled)
+                m_triangle_selectors[m_rr.mesh_id]->seed_fill_apply_on_triangles(new_state);
+            else
+                m_triangle_selectors[m_rr.mesh_id]->select_patch(m_rr.hit, int(m_rr.facet), camera_pos, m_cursor_radius, m_cursor_type,
+                                                                 new_state, trafo_matrix, m_triangle_splitting_enabled);
+            m_last_mouse_click = mouse_position;
+        }
+
+        return true;
+    }
+
+    if (action == SLAGizmoEventType::Moving && m_seed_fill_enabled) {
+        if (m_triangle_selectors.empty())
+            return false;
+
+        const Camera &       camera         = wxGetApp().plater()->get_camera();
+        const Selection &    selection      = m_parent.get_selection();
+        const ModelObject *  mo             = m_c->selection_info()->model_object();
+        const ModelInstance *mi             = mo->instances[selection.get_instance_idx()];
+        const Transform3d &  instance_trafo = mi->get_transformation().get_matrix();
+
+        // Precalculate transformations of individual meshes.
+        std::vector<Transform3d> trafo_matrices;
+        for (const ModelVolume *mv : mo->volumes)
+            if (mv->is_model_part())
+                trafo_matrices.emplace_back(instance_trafo * mv->get_matrix());
+
+        // Now "click" into all the prepared points and spill paint around them.
+        update_raycast_cache(mouse_position, camera, trafo_matrices);
+
+        if (m_rr.mesh_id == -1) {
+            // Clean selected by seed fill for all triangles
+            for (auto &triangle_selector : m_triangle_selectors)
+                triangle_selector->seed_fill_unselect_all_triangles();
+
+            // In case we have no valid hit, we can return.
+            return false;
+        }
+
+        assert(m_rr.mesh_id < int(m_triangle_selectors.size()));
+        m_triangle_selectors[m_rr.mesh_id]->seed_fill_select_triangles(m_rr.hit, int(m_rr.facet), m_seed_fill_angle);
+        return true;
+    }
+
+    if ((action == SLAGizmoEventType::LeftUp || action == SLAGizmoEventType::RightUp)
+        && m_button_down != Button::None) {
+        // Take snapshot and update ModelVolume data.
+        wxString action_name;
+        if (get_painter_type() == PainterGizmoType::FDM_SUPPORTS) {
+            if (shift_down)
+                action_name = _L("Remove selection");
+            else {
+                if (m_button_down == Button::Left)
+                    action_name = _L("Add supports");
+                else
+                    action_name = _L("Block supports");
+            }
+        }
+        if (get_painter_type() == PainterGizmoType::SEAM) {
+            if (shift_down)
+                action_name = _L("Remove selection");
+            else {
+                if (m_button_down == Button::Left)
+                    action_name = _L("Enforce seam");
+                else
+                    action_name = _L("Block seam");
+            }
+        }
+
+        activate_internal_undo_redo_stack(true);
+        Plater::TakeSnapshot snapshot(wxGetApp().plater(), action_name);
+        update_model_object();
+
+        m_button_down = Button::None;
+        m_last_mouse_click = Vec2d::Zero();
+        return true;
+    }
+
+    return false;
+}
+
+void GLGizmoMmuSegmentation::set_painter_gizmo_data(const Selection &selection)
+{
+    GLGizmoPainterBase::set_painter_gizmo_data(selection);
+
+    if (m_state != On)
+        return;
+
+    int prev_extruders_count = m_original_extruders_colors.size();
+    if (prev_extruders_count != wxGetApp().extruders_edited_cnt() || get_extruders_colors() != m_original_extruders_colors) {
+        this->init_extruders_data();
+        // Reinitialize triangle selectors because of change of extruder count need also change the size of GLIndexedVertexArray
+        if (prev_extruders_count != wxGetApp().extruders_edited_cnt())
+            this->init_model_triangle_selectors();
+    }
+}
+
+static void render_extruders_combo(const std::string                         &label,
+                                   const std::vector<std::string>            &extruders,
+                                   const std::vector<std::array<uint8_t, 3>> &extruders_colors,
+                                   size_t                                    &selection_idx)
+{
+    assert(!extruders_colors.empty());
+    assert(extruders_colors.size() == extruders_colors.size());
+
+    size_t selection_out = selection_idx;
+
+    // It is necessary to use BeginGroup(). Otherwise, when using SameLine() is called, then other items will be drawn inside the combobox.
+    ImGui::BeginGroup();
+    ImVec2 combo_pos = ImGui::GetCursorScreenPos();
+    if (ImGui::BeginCombo(label.c_str(), "")) {
+        for (size_t extruder_idx = 0; extruder_idx < extruders.size(); ++extruder_idx) {
+            ImGui::PushID(int(extruder_idx));
+            ImVec2 start_position = ImGui::GetCursorScreenPos();
+
+            if (ImGui::Selectable("", extruder_idx == selection_idx))
+                selection_out = extruder_idx;
+
+            ImGui::SameLine();
+            ImGuiStyle &style  = ImGui::GetStyle();
+            float       height = ImGui::GetTextLineHeight();
+            ImGui::GetWindowDrawList()->AddRectFilled(start_position, ImVec2(start_position.x + height + height / 2, start_position.y + height),
+                                                      IM_COL32(extruders_colors[extruder_idx][0], extruders_colors[extruder_idx][1], extruders_colors[extruder_idx][2], 255));
+            ImGui::GetWindowDrawList()->AddRect(start_position, ImVec2(start_position.x + height + height / 2, start_position.y + height), IM_COL32_BLACK);
+
+            ImGui::SetCursorScreenPos(ImVec2(start_position.x + height + height / 2 + style.FramePadding.x, start_position.y));
+            ImGui::Text("%s", extruders[extruder_idx].c_str());
+            ImGui::PopID();
+        }
+
+        ImGui::EndCombo();
+    }
+
+    ImVec2      backup_pos = ImGui::GetCursorScreenPos();
+    ImGuiStyle &style      = ImGui::GetStyle();
+
+    ImGui::SetCursorScreenPos(ImVec2(combo_pos.x + style.FramePadding.x, combo_pos.y + style.FramePadding.y));
+    ImVec2 p      = ImGui::GetCursorScreenPos();
+    float  height = ImGui::GetTextLineHeight();
+
+    ImGui::GetWindowDrawList()->AddRectFilled(p, ImVec2(p.x + height + height / 2, p.y + height),
+                                              IM_COL32(extruders_colors[selection_idx][0], extruders_colors[selection_idx][1],
+                                                       extruders_colors[selection_idx][2], 255));
+    ImGui::GetWindowDrawList()->AddRect(p, ImVec2(p.x + height + height / 2, p.y + height), IM_COL32_BLACK);
+
+    ImGui::SetCursorScreenPos(ImVec2(p.x + height + height / 2 + style.FramePadding.x, p.y));
+    ImGui::Text("%s", extruders[selection_out].c_str());
+    ImGui::SetCursorScreenPos(backup_pos);
+    ImGui::EndGroup();
+
+    selection_idx = selection_out;
+}
+
+void GLGizmoMmuSegmentation::on_render_input_window(float x, float y, float bottom_limit)
+{
+    if (!m_c->selection_info()->model_object())
+        return;
+
+    const float approx_height = m_imgui->scaled(23.0f);
+                            y = std::min(y, bottom_limit - approx_height);
+    m_imgui->set_next_window_pos(x, y, ImGuiCond_Always);
+
+    m_imgui->begin(on_get_name(), ImGuiWindowFlags_NoMove | ImGuiWindowFlags_AlwaysAutoResize | ImGuiWindowFlags_NoCollapse);
+
+    // First calculate width of all the texts that are could possibly be shown. We will decide set the dialog width based on that:
+    const float clipping_slider_left = std::max(m_imgui->calc_text_size(m_desc.at("clipping_of_view")).x,
+                                                m_imgui->calc_text_size(m_desc.at("reset_direction")).x) + m_imgui->scaled(1.5f);
+    const float cursor_slider_left       = m_imgui->calc_text_size(m_desc.at("cursor_size")).x + m_imgui->scaled(1.f);
+    const float autoset_slider_left      = m_imgui->calc_text_size(m_desc.at("seed_fill_angle")).x + m_imgui->scaled(1.f);
+    const float cursor_type_radio_left   = m_imgui->calc_text_size(m_desc.at("cursor_type")).x + m_imgui->scaled(1.f);
+    const float cursor_type_radio_width1 = m_imgui->calc_text_size(m_desc["circle"]).x + m_imgui->scaled(2.5f);
+    const float cursor_type_radio_width2 = m_imgui->calc_text_size(m_desc["sphere"]).x + m_imgui->scaled(2.5f);
+    const float button_width             = m_imgui->calc_text_size(m_desc.at("remove_all")).x + m_imgui->scaled(1.f);
+    const float buttons_width            = m_imgui->scaled(0.5f);
+    const float minimal_slider_width     = m_imgui->scaled(4.f);
+    const float color_button_width       = m_imgui->calc_text_size("").x + m_imgui->scaled(1.75f);
+    const float combo_label_width        = std::max(m_imgui->calc_text_size(m_desc.at("first_color")).x,
+                                                    m_imgui->calc_text_size(m_desc.at("second_color")).x) + m_imgui->scaled(1.f);
+
+    float caption_max    = 0.f;
+    float total_text_max = 0.;
+    for (const std::string &t : {"first_color", "second_color", "remove"}) {
+        caption_max    = std::max(caption_max, m_imgui->calc_text_size(m_desc.at(t + "_caption")).x);
+        total_text_max = std::max(total_text_max, caption_max + m_imgui->calc_text_size(m_desc.at(t)).x);
+    }
+    caption_max += m_imgui->scaled(1.f);
+    total_text_max += m_imgui->scaled(1.f);
+
+    float window_width = minimal_slider_width + std::max(autoset_slider_left, std::max(cursor_slider_left, clipping_slider_left));
+    window_width       = std::max(window_width, total_text_max);
+    window_width       = std::max(window_width, button_width);
+    window_width       = std::max(window_width, cursor_type_radio_left + cursor_type_radio_width1 + cursor_type_radio_width2);
+    window_width       = std::max(window_width, 2.f * buttons_width + m_imgui->scaled(1.f));
+
+    auto draw_text_with_caption = [this, &caption_max](const wxString &caption, const wxString &text) {
+        m_imgui->text_colored(ImGuiWrapper::COL_ORANGE_LIGHT, caption);
+        ImGui::SameLine(caption_max);
+        m_imgui->text(text);
+    };
+
+    for (const std::string &t : {"first_color", "second_color", "remove"})
+        draw_text_with_caption(m_desc.at(t + "_caption"), m_desc.at(t));
+
+    m_imgui->text("");
+    ImGui::Separator();
+
+    const std::array<uint8_t, 3> &select_first_color  = m_modified_extruders_colors[m_first_selected_extruder_idx];
+    const std::array<uint8_t, 3> &select_second_color = m_modified_extruders_colors[m_second_selected_extruder_idx];
+
+    m_imgui->text(m_desc.at("first_color"));
+    ImGui::SameLine(combo_label_width);
+    ImGui::PushItemWidth(window_width - combo_label_width - color_button_width);
+    render_extruders_combo("##first_color_combo", m_original_extruders_names, m_original_extruders_colors, m_first_selected_extruder_idx);
+    ImGui::SameLine();
+
+    ImVec4 first_color = ImVec4(float(select_first_color[0]) / 255.0f, float(select_first_color[1]) / 255.0f, float(select_first_color[2]) / 255.0f, 1.0f);
+    ImVec4 second_color = ImVec4(float(select_second_color[0]) / 255.0f, float(select_second_color[1]) / 255.0f, float(select_second_color[2]) / 255.0f, 1.0f);
+    if(ImGui::ColorEdit4("First color##color_picker", (float*)&first_color, ImGuiColorEditFlags_NoAlpha | ImGuiColorEditFlags_NoInputs | ImGuiColorEditFlags_NoLabel))
+        m_modified_extruders_colors[m_first_selected_extruder_idx] = {uint8_t(first_color.x * 255.0f), uint8_t(first_color.y * 255.0f), uint8_t(first_color.z * 255.0f)};
+
+    m_imgui->text(m_desc.at("second_color"));
+    ImGui::SameLine(combo_label_width);
+    ImGui::PushItemWidth(window_width - combo_label_width - color_button_width);
+    render_extruders_combo("##second_color_combo", m_original_extruders_names, m_original_extruders_colors, m_second_selected_extruder_idx);
+    ImGui::SameLine();
+    if(ImGui::ColorEdit4("Second color##color_picker", (float*)&second_color, ImGuiColorEditFlags_NoAlpha | ImGuiColorEditFlags_NoInputs | ImGuiColorEditFlags_NoLabel))
+        m_modified_extruders_colors[m_second_selected_extruder_idx] = {uint8_t(second_color.x * 255.0f), uint8_t(second_color.y * 255.0f), uint8_t(second_color.z * 255.0f)};
+
+    ImGui::Separator();
+
+    if (m_imgui->checkbox(_L("Seed fill"), m_seed_fill_enabled))
+        if (!m_seed_fill_enabled)
+            for (auto &triangle_selector : m_triangle_selectors)
+                triangle_selector->seed_fill_unselect_all_triangles();
+
+    m_imgui->text(m_desc["seed_fill_angle"] + ":");
+    ImGui::AlignTextToFramePadding();
+    std::string format_str = std::string("%.f") + I18N::translate_utf8("°", "Degree sign to use in the respective slider in FDM supports gizmo,"
+                                                                            "placed after the number with no whitespace in between.");
+    ImGui::SameLine(autoset_slider_left);
+    ImGui::PushItemWidth(window_width - autoset_slider_left);
+    m_imgui->disabled_begin(!m_seed_fill_enabled);
+    m_imgui->slider_float("##seed_fill_angle", &m_seed_fill_angle, 0.f, 90.f, format_str.data());
+    m_imgui->disabled_end();
+
+    ImGui::Separator();
+
+    if (m_imgui->button(m_desc.at("remove_all"))) {
+        Plater::TakeSnapshot snapshot(wxGetApp().plater(), wxString(_L("Reset selection")));
+        ModelObject *mo  = m_c->selection_info()->model_object();
+        int          idx = -1;
+        for (ModelVolume *mv : mo->volumes) {
+            if (mv->is_model_part()) {
+                ++idx;
+                size_t extruder_id = (mv->extruder_id() > 0) ? mv->extruder_id() - 1 : 0;
+                m_triangle_selectors[idx]->reset(EnforcerBlockerType(extruder_id));
+            }
+        }
+
+        update_model_object();
+        m_parent.set_as_dirty();
+    }
+
+    const float max_tooltip_width = ImGui::GetFontSize() * 20.0f;
+
+    ImGui::AlignTextToFramePadding();
+    m_imgui->text(m_desc.at("cursor_size"));
+    ImGui::SameLine(cursor_slider_left);
+    ImGui::PushItemWidth(window_width - cursor_slider_left);
+    ImGui::SliderFloat(" ", &m_cursor_radius, CursorRadiusMin, CursorRadiusMax, "%.2f");
+    if (ImGui::IsItemHovered()) {
+        ImGui::BeginTooltip();
+        ImGui::PushTextWrapPos(max_tooltip_width);
+        ImGui::TextUnformatted(_L("Alt + Mouse wheel").ToUTF8().data());
+        ImGui::PopTextWrapPos();
+        ImGui::EndTooltip();
+    }
+
+    ImGui::AlignTextToFramePadding();
+    m_imgui->text(m_desc.at("cursor_type"));
+    ImGui::SameLine(cursor_type_radio_left + m_imgui->scaled(0.f));
+    ImGui::PushItemWidth(cursor_type_radio_width1);
+
+    bool sphere_sel = m_cursor_type == TriangleSelector::CursorType::SPHERE;
+    if (m_imgui->radio_button(m_desc["sphere"], sphere_sel))
+        sphere_sel = true;
+
+    if (ImGui::IsItemHovered()) {
+        ImGui::BeginTooltip();
+        ImGui::PushTextWrapPos(max_tooltip_width);
+        ImGui::TextUnformatted(_L("Paints all facets inside, regardless of their orientation.").ToUTF8().data());
+        ImGui::PopTextWrapPos();
+        ImGui::EndTooltip();
+    }
+
+    ImGui::SameLine(cursor_type_radio_left + cursor_type_radio_width2 + m_imgui->scaled(0.f));
+    ImGui::PushItemWidth(cursor_type_radio_width2);
+
+    if (m_imgui->radio_button(m_desc["circle"], !sphere_sel))
+        sphere_sel = false;
+
+    if (ImGui::IsItemHovered()) {
+        ImGui::BeginTooltip();
+        ImGui::PushTextWrapPos(max_tooltip_width);
+        ImGui::TextUnformatted(_L("Ignores facets facing away from the camera.").ToUTF8().data());
+        ImGui::PopTextWrapPos();
+        ImGui::EndTooltip();
+    }
+
+    m_cursor_type = sphere_sel ? TriangleSelector::CursorType::SPHERE : TriangleSelector::CursorType::CIRCLE;
+
+    m_imgui->checkbox(_L("Split triangles"), m_triangle_splitting_enabled);
+
+    ImGui::Separator();
+    if (m_c->object_clipper()->get_position() == 0.f) {
+        ImGui::AlignTextToFramePadding();
+        m_imgui->text(m_desc.at("clipping_of_view"));
+    } else {
+        if (m_imgui->button(m_desc.at("reset_direction"))) {
+            wxGetApp().CallAfter([this]() { m_c->object_clipper()->set_position(-1., false); });
+        }
+    }
+
+    ImGui::SameLine(clipping_slider_left);
+    ImGui::PushItemWidth(window_width - clipping_slider_left);
+    auto clp_dist = float(m_c->object_clipper()->get_position());
+    if (ImGui::SliderFloat("  ", &clp_dist, 0.f, 1.f, "%.2f"))
+        m_c->object_clipper()->set_position(clp_dist, true);
+    if (ImGui::IsItemHovered()) {
+        ImGui::BeginTooltip();
+        ImGui::PushTextWrapPos(max_tooltip_width);
+        ImGui::TextUnformatted(_L("Ctrl + Mouse wheel").ToUTF8().data());
+        ImGui::PopTextWrapPos();
+        ImGui::EndTooltip();
+    }
+    m_imgui->end();
+}
+
+void GLGizmoMmuSegmentation::update_model_object() const
+{
+    bool updated = false;
+    ModelObject* mo = m_c->selection_info()->model_object();
+    int idx = -1;
+    for (ModelVolume* mv : mo->volumes) {
+        if (! mv->is_model_part())
+            continue;
+        ++idx;
+        updated |= mv->mmu_segmentation_facets.set(*m_triangle_selectors[idx].get());
+    }
+
+    if (updated)
+        m_parent.post_event(SimpleEvent(EVT_GLCANVAS_SCHEDULE_BACKGROUND_PROCESS));
+}
+
+void GLGizmoMmuSegmentation::init_model_triangle_selectors()
+{
+    const ModelObject *mo = m_c->selection_info()->model_object();
+    m_triangle_selectors.clear();
+
+    for (const ModelVolume *mv : mo->volumes) {
+        if (!mv->is_model_part())
+            continue;
+
+        // This mesh does not account for the possible Z up SLA offset.
+        const TriangleMesh *mesh = &mv->mesh();
+
+        size_t extruder_id = (mv->extruder_id() > 0) ? mv->extruder_id() - 1 : 0;
+        m_triangle_selectors.emplace_back(std::make_unique<TriangleSelectorMmuGui>(*mesh, m_modified_extruders_colors));
+        m_triangle_selectors.back()->deserialize(mv->mmu_segmentation_facets.get_data(), EnforcerBlockerType(extruder_id));
+    }
+}
+
+void GLGizmoMmuSegmentation::update_from_model_object()
+{
+    wxBusyCursor wait;
+    this->init_model_triangle_selectors();
+}
+
+PainterGizmoType GLGizmoMmuSegmentation::get_painter_type() const
+{
+    return PainterGizmoType::MMU_SEGMENTATION;
+}
+
+void TriangleSelectorMmuGui::render(ImGuiWrapper *imgui)
+{
+    std::vector<int> color_cnt(m_iva_colors.size());
+    int              seed_fill_cnt = 0;
+    for (auto &iva_color : m_iva_colors)
+        iva_color.release_geometry();
+    m_iva_seed_fill.release_geometry();
+
+    for (size_t color_idx = 0; color_idx < m_iva_colors.size(); ++color_idx) {
+        for (const Triangle &tr : m_triangles) {
+            if (!tr.valid || tr.is_split() || tr.is_selected_by_seed_fill() || tr.get_state() != EnforcerBlockerType(color_idx))
+                continue;
+
+            for (int i = 0; i < 3; ++i)
+                m_iva_colors[color_idx].push_geometry(double(m_vertices[tr.verts_idxs[i]].v[0]),
+                                                      double(m_vertices[tr.verts_idxs[i]].v[1]),
+                                                      double(m_vertices[tr.verts_idxs[i]].v[2]),
+                                                      double(tr.normal[0]),
+                                                      double(tr.normal[1]),
+                                                      double(tr.normal[2]));
+            m_iva_colors[color_idx].push_triangle(color_cnt[color_idx], color_cnt[color_idx] + 1, color_cnt[color_idx] + 2);
+            color_cnt[color_idx] += 3;
+        }
+    }
+
+    for (const Triangle &tr : m_triangles) {
+        if (!tr.valid || tr.is_split() || !tr.is_selected_by_seed_fill()) continue;
+
+        for (int i = 0; i < 3; ++i)
+            m_iva_seed_fill.push_geometry(double(m_vertices[tr.verts_idxs[i]].v[0]),
+                                          double(m_vertices[tr.verts_idxs[i]].v[1]),
+                                          double(m_vertices[tr.verts_idxs[i]].v[2]),
+                                          double(tr.normal[0]),
+                                          double(tr.normal[1]),
+                                          double(tr.normal[2]));
+        m_iva_seed_fill.push_triangle(seed_fill_cnt, seed_fill_cnt + 1, seed_fill_cnt + 2);
+        seed_fill_cnt += 3;
+    }
+
+    for (auto &iva_color : m_iva_colors)
+        iva_color.finalize_geometry(true);
+    m_iva_seed_fill.finalize_geometry(true);
+
+    std::vector<bool> render_colors(m_iva_colors.size());
+    for (size_t color_idx = 0; color_idx < m_iva_colors.size(); ++color_idx)
+        render_colors[color_idx] = m_iva_colors[color_idx].has_VBOs();
+    bool render_seed_fill = m_iva_seed_fill.has_VBOs();
+
+    auto *shader = wxGetApp().get_shader("gouraud");
+    if (!shader) return;
+
+    shader->start_using();
+    ScopeGuard guard([shader]() {
+        if (shader)
+            shader->stop_using();
+    });
+    shader->set_uniform("slope.actived", false);
+
+    for (size_t color_idx = 0; color_idx < m_iva_colors.size(); ++color_idx) {
+        if (render_colors[color_idx]) {
+            std::array<float, 4> color = {float(m_colors[color_idx][0]) / 255.0f, float(m_colors[color_idx][1]) / 255.0f,
+                                          float(m_colors[color_idx][2]) / 255.0f, 1.f};
+            shader->set_uniform("uniform_color", color);
+            m_iva_colors[color_idx].render();
+        }
+    }
+
+    if (render_seed_fill) {
+        std::array<float, 4> color = {0.f, 1.f, 0.44f, 1.f};
+        shader->set_uniform("uniform_color", color);
+        m_iva_seed_fill.render();
+    }
+}
+
+} // namespace Slic3r
diff --git a/src/slic3r/GUI/Gizmos/GLGizmoMmuSegmentation.hpp b/src/slic3r/GUI/Gizmos/GLGizmoMmuSegmentation.hpp
new file mode 100644
index 000000000..35ade52fb
--- /dev/null
+++ b/src/slic3r/GUI/Gizmos/GLGizmoMmuSegmentation.hpp
@@ -0,0 +1,69 @@
+#ifndef slic3r_GLGizmoMmuSegmentation_hpp_
+#define slic3r_GLGizmoMmuSegmentation_hpp_
+
+#include "GLGizmoPainterBase.hpp"
+
+namespace Slic3r::GUI {
+
+class TriangleSelectorMmuGui : public TriangleSelectorGUI {
+public:
+    explicit TriangleSelectorMmuGui(const TriangleMesh& mesh, const std::vector<std::array<uint8_t, 3>> &colors)
+        : TriangleSelectorGUI(mesh), m_colors(colors) {
+        m_iva_colors = std::vector<GLIndexedVertexArray>(colors.size());
+    }
+
+    // Render current selection. Transformation matrices are supposed
+    // to be already set.
+    void render(ImGuiWrapper* imgui) override;
+
+private:
+    const std::vector<std::array<uint8_t, 3>> &m_colors;
+    std::vector<GLIndexedVertexArray> m_iva_colors;
+};
+
+class GLGizmoMmuSegmentation : public GLGizmoPainterBase
+{
+public:
+    GLGizmoMmuSegmentation(GLCanvas3D& parent, const std::string& icon_filename, unsigned int sprite_id)
+        : GLGizmoPainterBase(parent, icon_filename, sprite_id) {}
+
+    void render_painter_gizmo() const override;
+
+    bool gizmo_event(SLAGizmoEventType action, const Vec2d& mouse_position, bool shift_down, bool alt_down, bool control_down) override;
+
+    virtual void set_painter_gizmo_data(const Selection& selection) override;
+
+protected:
+    void on_render_input_window(float x, float y, float bottom_limit) override;
+    std::string on_get_name() const override;
+
+    bool on_is_selectable() const override;
+
+    size_t                              m_first_selected_extruder_idx  = 0;
+    size_t                              m_second_selected_extruder_idx = 1;
+    std::vector<std::string>            m_original_extruders_names;
+    std::vector<std::array<uint8_t, 3>> m_original_extruders_colors;
+    std::vector<std::array<uint8_t, 3>> m_modified_extruders_colors;
+
+private:
+    bool on_init() override;
+
+    void update_model_object() const override;
+    void update_from_model_object() override;
+
+    void on_opening() override {}
+    void on_shutdown() override;
+    PainterGizmoType get_painter_type() const override;
+
+    void init_model_triangle_selectors();
+    void init_extruders_data();
+
+    // This map holds all translated description texts, so they can be easily referenced during layout calculations
+    // etc. When language changes, GUI is recreated and this class constructed again, so the change takes effect.
+    std::map<std::string, wxString> m_desc;
+};
+
+} // namespace Slic3r
+
+
+#endif // slic3r_GLGizmoMmuSegmentation_hpp_
diff --git a/src/slic3r/GUI/Gizmos/GLGizmoPainterBase.cpp b/src/slic3r/GUI/Gizmos/GLGizmoPainterBase.cpp
index 3a932c598..adb95b045 100644
--- a/src/slic3r/GUI/Gizmos/GLGizmoPainterBase.cpp
+++ b/src/slic3r/GUI/Gizmos/GLGizmoPainterBase.cpp
@@ -73,7 +73,7 @@ void GLGizmoPainterBase::activate_internal_undo_redo_stack(bool activate)
     if (activate && ! m_internal_stack_active) {
         wxString str = get_painter_type() == PainterGizmoType::FDM_SUPPORTS
                            ? _L("Entering Paint-on supports")
-                           : _L("Entering Seam painting");
+                           : (get_painter_type() == PainterGizmoType::MMU_SEGMENTATION ? _L("Entering MMU segmentation") : _L("Entering Seam painting"));
         Plater::TakeSnapshot(wxGetApp().plater(), str);
         wxGetApp().plater()->enter_gizmos_stack();
         m_internal_stack_active = true;
@@ -81,7 +81,7 @@ void GLGizmoPainterBase::activate_internal_undo_redo_stack(bool activate)
     if (! activate && m_internal_stack_active) {
         wxString str = get_painter_type() == PainterGizmoType::SEAM
                            ? _L("Leaving Seam painting")
-                           : _L("Leaving Paint-on supports");
+                           : (get_painter_type() == PainterGizmoType::MMU_SEGMENTATION ? _L("Leaving MMU segmentation") : _L("Leaving Paint-on supports"));
         wxGetApp().plater()->leave_gizmos_stack();
         Plater::TakeSnapshot(wxGetApp().plater(), str);
         m_internal_stack_active = false;
@@ -180,11 +180,12 @@ void GLGizmoPainterBase::render_cursor() const
     if (m_rr.mesh_id == -1)
         return;
 
-
-    if (m_cursor_type == TriangleSelector::SPHERE)
-        render_cursor_sphere(trafo_matrices[m_rr.mesh_id]);
-    else
-        render_cursor_circle();
+    if (!m_seed_fill_enabled) {
+        if (m_cursor_type == TriangleSelector::SPHERE)
+            render_cursor_sphere(trafo_matrices[m_rr.mesh_id]);
+        else
+            render_cursor_circle();
+    }
 }
 
 
@@ -388,14 +389,50 @@ bool GLGizmoPainterBase::gizmo_event(SLAGizmoEventType action, const Vec2d& mous
             Vec3f camera_pos = (trafo_matrix.inverse() * camera.get_position()).cast<float>();
 
             assert(m_rr.mesh_id < int(m_triangle_selectors.size()));
-            m_triangle_selectors[m_rr.mesh_id]->select_patch(m_rr.hit, m_rr.facet, camera_pos,
-                                       m_cursor_radius, m_cursor_type, new_state, trafo_matrix);
+            if (m_seed_fill_enabled)
+                m_triangle_selectors[m_rr.mesh_id]->seed_fill_apply_on_triangles(new_state);
+            else
+                m_triangle_selectors[m_rr.mesh_id]->select_patch(m_rr.hit, int(m_rr.facet), camera_pos, m_cursor_radius, m_cursor_type,
+                                                                 new_state, trafo_matrix, m_triangle_splitting_enabled);
             m_last_mouse_click = mouse_position;
         }
 
         return true;
     }
 
+    if (action == SLAGizmoEventType::Moving && m_seed_fill_enabled) {
+        if (m_triangle_selectors.empty())
+            return false;
+
+        const Camera &       camera         = wxGetApp().plater()->get_camera();
+        const Selection &    selection      = m_parent.get_selection();
+        const ModelObject *  mo             = m_c->selection_info()->model_object();
+        const ModelInstance *mi             = mo->instances[selection.get_instance_idx()];
+        const Transform3d &  instance_trafo = mi->get_transformation().get_matrix();
+
+        // Precalculate transformations of individual meshes.
+        std::vector<Transform3d> trafo_matrices;
+        for (const ModelVolume *mv : mo->volumes)
+            if (mv->is_model_part())
+                trafo_matrices.emplace_back(instance_trafo * mv->get_matrix());
+
+        // Now "click" into all the prepared points and spill paint around them.
+        update_raycast_cache(mouse_position, camera, trafo_matrices);
+
+        if (m_rr.mesh_id == -1) {
+            // Clean selected by seed fill for all triangles
+            for (auto &triangle_selector : m_triangle_selectors)
+                triangle_selector->seed_fill_unselect_all_triangles();
+
+            // In case we have no valid hit, we can return.
+            return false;
+        }
+
+        assert(m_rr.mesh_id < int(m_triangle_selectors.size()));
+        m_triangle_selectors[m_rr.mesh_id]->seed_fill_select_triangles(m_rr.hit, int(m_rr.facet), m_seed_fill_angle);
+        return true;
+    }
+
     if ((action == SLAGizmoEventType::LeftUp || action == SLAGizmoEventType::RightUp)
       && m_button_down != Button::None) {
         // Take snapshot and update ModelVolume data.
@@ -558,12 +595,14 @@ void TriangleSelectorGUI::render(ImGuiWrapper* imgui)
 {
     int enf_cnt = 0;
     int blc_cnt = 0;
+    int seed_fill_cnt = 0;
 
     m_iva_enforcers.release_geometry();
     m_iva_blockers.release_geometry();
+    m_iva_seed_fill.release_geometry();
 
     for (const Triangle& tr : m_triangles) {
-        if (! tr.valid || tr.is_split() || tr.get_state() == EnforcerBlockerType::NONE)
+        if (!tr.valid || tr.is_split() || tr.get_state() == EnforcerBlockerType::NONE || tr.is_selected_by_seed_fill())
             continue;
 
         GLIndexedVertexArray& va = tr.get_state() == EnforcerBlockerType::ENFORCER
@@ -580,17 +619,32 @@ void TriangleSelectorGUI::render(ImGuiWrapper* imgui)
                              double(tr.normal[0]),
                              double(tr.normal[1]),
                              double(tr.normal[2]));
-        va.push_triangle(cnt,
-                         cnt+1,
-                         cnt+2);
+        va.push_triangle(cnt, cnt + 1, cnt + 2);
         cnt += 3;
     }
 
+    for (const Triangle &tr : m_triangles) {
+        if (!tr.valid || tr.is_split() || !tr.is_selected_by_seed_fill())
+            continue;
+
+        for (int i = 0; i < 3; ++i)
+            m_iva_seed_fill.push_geometry(double(m_vertices[tr.verts_idxs[i]].v[0]),
+                                          double(m_vertices[tr.verts_idxs[i]].v[1]),
+                                          double(m_vertices[tr.verts_idxs[i]].v[2]),
+                                          double(tr.normal[0]),
+                                          double(tr.normal[1]),
+                                          double(tr.normal[2]));
+        m_iva_seed_fill.push_triangle(seed_fill_cnt, seed_fill_cnt + 1, seed_fill_cnt + 2);
+        seed_fill_cnt += 3;
+    }
+
     m_iva_enforcers.finalize_geometry(true);
     m_iva_blockers.finalize_geometry(true);
+    m_iva_seed_fill.finalize_geometry(true);
 
     bool render_enf = m_iva_enforcers.has_VBOs();
     bool render_blc = m_iva_blockers.has_VBOs();
+    bool render_seed_fill = m_iva_seed_fill.has_VBOs();
 
     auto* shader = wxGetApp().get_shader("gouraud");
     if (! shader)
@@ -612,6 +666,12 @@ void TriangleSelectorGUI::render(ImGuiWrapper* imgui)
         m_iva_blockers.render();
     }
 
+    if (render_seed_fill) {
+        std::array<float, 4> color = { 0.f, 1.00f, 0.44f, 1.f };
+        shader->set_uniform("uniform_color", color);
+        m_iva_seed_fill.render();
+    }
+
 
 #ifdef PRUSASLICER_TRIANGLE_SELECTOR_DEBUG
     if (imgui)
diff --git a/src/slic3r/GUI/Gizmos/GLGizmoPainterBase.hpp b/src/slic3r/GUI/Gizmos/GLGizmoPainterBase.hpp
index da415ce09..55e4df865 100644
--- a/src/slic3r/GUI/Gizmos/GLGizmoPainterBase.hpp
+++ b/src/slic3r/GUI/Gizmos/GLGizmoPainterBase.hpp
@@ -12,7 +12,6 @@
 
 
 
-
 namespace Slic3r {
 
 enum class EnforcerBlockerType : int8_t;
@@ -22,10 +21,12 @@ namespace GUI {
 enum class SLAGizmoEventType : unsigned char;
 class ClippingPlane;
 struct Camera;
+class GLGizmoMmuSegmentation;
 
 enum class PainterGizmoType {
     FDM_SUPPORTS,
-    SEAM
+    SEAM,
+    MMU_SEGMENTATION
 };
 
 
@@ -36,7 +37,7 @@ public:
 
     // Render current selection. Transformation matrices are supposed
     // to be already set.
-    void render(ImGuiWrapper* imgui = nullptr);
+    virtual void render(ImGuiWrapper* imgui = nullptr);
 
 #ifdef PRUSASLICER_TRIANGLE_SELECTOR_DEBUG
     void render_debug(ImGuiWrapper* imgui);
@@ -48,6 +49,8 @@ private:
     GLIndexedVertexArray m_iva_enforcers;
     GLIndexedVertexArray m_iva_blockers;
     std::array<GLIndexedVertexArray, 3> m_varrays;
+protected:
+    GLIndexedVertexArray m_iva_seed_fill;
 };
 
 
@@ -65,8 +68,8 @@ private:
 public:
     GLGizmoPainterBase(GLCanvas3D& parent, const std::string& icon_filename, unsigned int sprite_id);
     ~GLGizmoPainterBase() override {}
-    void set_painter_gizmo_data(const Selection& selection);
-    bool gizmo_event(SLAGizmoEventType action, const Vec2d& mouse_position, bool shift_down, bool alt_down, bool control_down);
+    virtual void set_painter_gizmo_data(const Selection& selection);
+    virtual bool gizmo_event(SLAGizmoEventType action, const Vec2d& mouse_position, bool shift_down, bool alt_down, bool control_down);
 
     // Following function renders the triangles and cursor. Having this separated
     // from usual on_render method allows to render them before transparent objects,
@@ -94,6 +97,9 @@ protected:
 
     TriangleSelector::CursorType m_cursor_type = TriangleSelector::SPHERE;
 
+    bool  m_triangle_splitting_enabled = true;
+    bool  m_seed_fill_enabled          = false;
+    float m_seed_fill_angle            = 0.f;
 
 private:
     bool is_mesh_point_clipped(const Vec3d& point, const Transform3d& trafo) const;
@@ -141,6 +147,8 @@ protected:
     void on_load(cereal::BinaryInputArchive& ar) override;
     void on_save(cereal::BinaryOutputArchive& ar) const override {}
     CommonGizmosDataID on_get_requirements() const override;
+
+    friend class ::Slic3r::GUI::GLGizmoMmuSegmentation;
 };
 
 
diff --git a/src/slic3r/GUI/Gizmos/GLGizmoRotate.cpp b/src/slic3r/GUI/Gizmos/GLGizmoRotate.cpp
index 6b6905e4d..63cb051aa 100644
--- a/src/slic3r/GUI/Gizmos/GLGizmoRotate.cpp
+++ b/src/slic3r/GUI/Gizmos/GLGizmoRotate.cpp
@@ -110,7 +110,7 @@ void GLGizmoRotate::on_update(const UpdateData& data)
     Vec2d orig_dir = Vec2d::UnitX();
     Vec2d new_dir = mouse_pos.normalized();
 
-    double theta = ::acos(clamp(-1.0, 1.0, new_dir.dot(orig_dir)));
+    double theta = ::acos(std::clamp(new_dir.dot(orig_dir), -1.0, 1.0));
     if (cross2(orig_dir, new_dir) < 0.0)
         theta = 2.0 * (double)PI - theta;
 
diff --git a/src/slic3r/GUI/Gizmos/GLGizmos.hpp b/src/slic3r/GUI/Gizmos/GLGizmos.hpp
index e8e73959c..0252dee15 100644
--- a/src/slic3r/GUI/Gizmos/GLGizmos.hpp
+++ b/src/slic3r/GUI/Gizmos/GLGizmos.hpp
@@ -32,6 +32,7 @@ enum class SLAGizmoEventType : unsigned char {
 #include "slic3r/GUI/Gizmos/GLGizmoFlatten.hpp"
 #include "slic3r/GUI/Gizmos/GLGizmoSlaSupports.hpp"
 #include "slic3r/GUI/Gizmos/GLGizmoFdmSupports.hpp"
+#include "slic3r/GUI/Gizmos/GLGizmoMmuSegmentation.hpp"
 #include "slic3r/GUI/Gizmos/GLGizmoCut.hpp"
 #include "slic3r/GUI/Gizmos/GLGizmoHollow.hpp"
 
diff --git a/src/slic3r/GUI/Gizmos/GLGizmosCommon.hpp b/src/slic3r/GUI/Gizmos/GLGizmosCommon.hpp
index 28be1b97f..6f59dc95e 100644
--- a/src/slic3r/GUI/Gizmos/GLGizmosCommon.hpp
+++ b/src/slic3r/GUI/Gizmos/GLGizmosCommon.hpp
@@ -32,7 +32,8 @@ enum class SLAGizmoEventType : unsigned char {
     ManualEditing,
     MouseWheelUp,
     MouseWheelDown,
-    ResetClippingPlane
+    ResetClippingPlane,
+    Moving
 };
 
 
diff --git a/src/slic3r/GUI/Gizmos/GLGizmosManager.cpp b/src/slic3r/GUI/Gizmos/GLGizmosManager.cpp
index bd02cbbbb..351c12901 100644
--- a/src/slic3r/GUI/Gizmos/GLGizmosManager.cpp
+++ b/src/slic3r/GUI/Gizmos/GLGizmosManager.cpp
@@ -19,6 +19,7 @@
 #include "slic3r/GUI/Gizmos/GLGizmoCut.hpp"
 #include "slic3r/GUI/Gizmos/GLGizmoHollow.hpp"
 #include "slic3r/GUI/Gizmos/GLGizmoSeam.hpp"
+#include "slic3r/GUI/Gizmos/GLGizmoMmuSegmentation.hpp"
 
 #include "libslic3r/Model.hpp"
 #include "libslic3r/PresetBundle.hpp"
@@ -108,6 +109,7 @@ bool GLGizmosManager::init()
     m_gizmos.emplace_back(new GLGizmoSlaSupports(m_parent, "sla_supports.svg", 6));
     m_gizmos.emplace_back(new GLGizmoFdmSupports(m_parent, "fdm_supports.svg", 7));
     m_gizmos.emplace_back(new GLGizmoSeam(m_parent, "seam.svg", 8));
+    m_gizmos.emplace_back(new GLGizmoMmuSegmentation(m_parent, "fdm_supports.svg", 9));
 
     m_common_gizmos_data.reset(new CommonGizmosDataPool(&m_parent));
 
@@ -401,6 +403,7 @@ void GLGizmosManager::set_painter_gizmo_data()
 
     dynamic_cast<GLGizmoFdmSupports*>(m_gizmos[FdmSupports].get())->set_painter_gizmo_data(m_parent.get_selection());
     dynamic_cast<GLGizmoSeam*>(m_gizmos[Seam].get())->set_painter_gizmo_data(m_parent.get_selection());
+    dynamic_cast<GLGizmoMmuSegmentation*>(m_gizmos[MmuSegmentation].get())->set_painter_gizmo_data(m_parent.get_selection());
 }
 
 // Returns true if the gizmo used the event to do something, false otherwise.
@@ -417,6 +420,8 @@ bool GLGizmosManager::gizmo_event(SLAGizmoEventType action, const Vec2d& mouse_p
         return dynamic_cast<GLGizmoFdmSupports*>(m_gizmos[FdmSupports].get())->gizmo_event(action, mouse_position, shift_down, alt_down, control_down);
     else if (m_current == Seam)
         return dynamic_cast<GLGizmoSeam*>(m_gizmos[Seam].get())->gizmo_event(action, mouse_position, shift_down, alt_down, control_down);
+    else if (m_current == MmuSegmentation)
+        return dynamic_cast<GLGizmoMmuSegmentation*>(m_gizmos[MmuSegmentation].get())->gizmo_event(action, mouse_position, shift_down, alt_down, control_down);
     else
         return false;
 }
@@ -493,7 +498,7 @@ bool GLGizmosManager::on_mouse_wheel(wxMouseEvent& evt)
 {
     bool processed = false;
 
-    if (m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam) {
+    if (m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam || m_current == MmuSegmentation) {
         float rot = (float)evt.GetWheelRotation() / (float)evt.GetWheelDelta();
         if (gizmo_event((rot > 0.f ? SLAGizmoEventType::MouseWheelUp : SLAGizmoEventType::MouseWheelDown), Vec2d::Zero(), evt.ShiftDown(), evt.AltDown(), evt.ControlDown()))
             processed = true;
@@ -518,9 +523,11 @@ bool GLGizmosManager::on_mouse(wxMouseEvent& evt)
     bool control_down = evt.CmdDown();
 
     // mouse anywhere
-    if (evt.Moving())
+    if (evt.Moving()) {
         m_tooltip = update_hover_state(mouse_pos);
-    else if (evt.LeftUp()) {
+        if (m_current == MmuSegmentation)
+            gizmo_event(SLAGizmoEventType::Moving, mouse_pos, evt.ShiftDown(), evt.AltDown());
+    } else if (evt.LeftUp()) {
         if (m_mouse_capture.left) {
             processed = true;
             m_mouse_capture.left = false;
@@ -626,7 +633,7 @@ bool GLGizmosManager::on_mouse(wxMouseEvent& evt)
         m_tooltip = "";
 
         if (evt.LeftDown() && (!control_down || grabber_contains_mouse())) {
-            if ((m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam)
+            if ((m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam || m_current == MmuSegmentation)
                 && gizmo_event(SLAGizmoEventType::LeftDown, mouse_pos, evt.ShiftDown(), evt.AltDown()))
                 // the gizmo got the event and took some action, there is no need to do anything more
                 processed = true;
@@ -652,16 +659,16 @@ bool GLGizmosManager::on_mouse(wxMouseEvent& evt)
             // event was taken care of by the SlaSupports gizmo
             processed = true;
         }
-        else if (evt.RightDown() && !control_down && selected_object_idx != -1 && (m_current == FdmSupports || m_current == Seam)
+        else if (evt.RightDown() && !control_down && selected_object_idx != -1 && (m_current == FdmSupports || m_current == Seam || m_current == MmuSegmentation)
             && gizmo_event(SLAGizmoEventType::RightDown, mouse_pos)) {
-            // event was taken care of by the FdmSupports / Seam gizmo
+            // event was taken care of by the FdmSupports / Seam / MMUPainting gizmo
             processed = true;
         }
         else if (evt.Dragging() && m_parent.get_move_volume_id() != -1
-             && (m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam))
+            && (m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam || m_current == MmuSegmentation))
             // don't allow dragging objects with the Sla gizmo on
             processed = true;
-        else if (evt.Dragging() && !control_down && (m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam)
+        else if (evt.Dragging() && !control_down && (m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam  || m_current == MmuSegmentation)
             && gizmo_event(SLAGizmoEventType::Dragging, mouse_pos, evt.ShiftDown(), evt.AltDown())) {
             // the gizmo got the event and took some action, no need to do anything more here
             m_parent.set_as_dirty();
@@ -674,7 +681,7 @@ bool GLGizmosManager::on_mouse(wxMouseEvent& evt)
             else if (evt.RightIsDown())
                 gizmo_event(SLAGizmoEventType::RightUp, mouse_pos, evt.ShiftDown(), evt.AltDown(), true);
         }
-        else if (evt.LeftUp() && (m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam) && !m_parent.is_mouse_dragging()) {
+        else if (evt.LeftUp() && (m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam || m_current == MmuSegmentation) && !m_parent.is_mouse_dragging()) {
             // in case SLA/FDM gizmo is selected, we just pass the LeftUp event and stop processing - neither
             // object moving or selecting is suppressed in that case
             gizmo_event(SLAGizmoEventType::LeftUp, mouse_pos, evt.ShiftDown(), evt.AltDown(), control_down);
@@ -684,7 +691,7 @@ bool GLGizmosManager::on_mouse(wxMouseEvent& evt)
             // to avoid to loose the selection when user clicks an the white faces of a different object while the Flatten gizmo is active
             processed = true;
         }
-        else if (evt.RightUp() && (m_current == FdmSupports || m_current == Seam) && !m_parent.is_mouse_dragging()) {
+        else if (evt.RightUp() && (m_current == FdmSupports || m_current == Seam || m_current == MmuSegmentation) && !m_parent.is_mouse_dragging()) {
             gizmo_event(SLAGizmoEventType::RightUp, mouse_pos, evt.ShiftDown(), evt.AltDown(), control_down);
             processed = true;
         }
@@ -768,7 +775,7 @@ bool GLGizmosManager::on_char(wxKeyEvent& evt)
         case 'r' :
         case 'R' :
         {
-            if ((m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam) && gizmo_event(SLAGizmoEventType::ResetClippingPlane))
+            if ((m_current == SlaSupports || m_current == Hollow || m_current == FdmSupports || m_current == Seam || m_current == MmuSegmentation) && gizmo_event(SLAGizmoEventType::ResetClippingPlane))
                 processed = true;
 
             break;
diff --git a/src/slic3r/GUI/Gizmos/GLGizmosManager.hpp b/src/slic3r/GUI/Gizmos/GLGizmosManager.hpp
index c8951966e..01d7ea85c 100644
--- a/src/slic3r/GUI/Gizmos/GLGizmosManager.hpp
+++ b/src/slic3r/GUI/Gizmos/GLGizmosManager.hpp
@@ -68,6 +68,7 @@ public:
         SlaSupports,
         FdmSupports,
         Seam,
+        MmuSegmentation,
         Undefined
     };
 
diff --git a/src/slic3r/GUI/Plater.cpp b/src/slic3r/GUI/Plater.cpp
index 3aff82fb4..66e1c9ed6 100644
--- a/src/slic3r/GUI/Plater.cpp
+++ b/src/slic3r/GUI/Plater.cpp
@@ -3273,6 +3273,7 @@ void Plater::priv::reload_from_disk()
                         new_volume->convert_from_meters();
                     new_volume->supported_facets.assign(old_volume->supported_facets);
                     new_volume->seam_facets.assign(old_volume->seam_facets);
+                    new_volume->mmu_segmentation_facets.assign(old_volume->mmu_segmentation_facets);
                     std::swap(old_model_object->volumes[sel_v.volume_idx], old_model_object->volumes.back());
                     old_model_object->delete_volume(old_model_object->volumes.size() - 1);
 #if ENABLE_ALLOW_NEGATIVE_Z
@@ -3345,13 +3346,14 @@ void Plater::priv::fix_through_netfabb(const int obj_idx, const int vol_idx/* =
 
     ModelObject* mo = model.objects[obj_idx];
 
-    // If there are custom supports/seams, remove them. Fixed mesh
+    // If there are custom supports/seams/mmu segmentation, remove them. Fixed mesh
     // may be different and they would make no sense.
     bool paint_removed = false;
     for (ModelVolume* mv : mo->volumes) {
-        paint_removed |= ! mv->supported_facets.empty() || ! mv->seam_facets.empty();
+        paint_removed |= ! mv->supported_facets.empty() || ! mv->seam_facets.empty() || ! mv->mmu_segmentation_facets.empty();
         mv->supported_facets.clear();
         mv->seam_facets.clear();
+        mv->mmu_segmentation_facets.clear();
     }
     if (paint_removed) {
         // snapshot_time is captured by copy so the lambda knows where to undo/redo to.
diff --git a/src/slic3r/GUI/Tab.cpp b/src/slic3r/GUI/Tab.cpp
index 5920ab5e3..ff2a70a94 100644
--- a/src/slic3r/GUI/Tab.cpp
+++ b/src/slic3r/GUI/Tab.cpp
@@ -1589,6 +1589,7 @@ void TabPrint::build()
 
         optgroup = page->new_optgroup(L("Advanced"));
         optgroup->append_single_option_line("interface_shells");
+        optgroup->append_single_option_line("mmu_segmented_region_max_width");
 
     page = add_options_page(L("Advanced"), "wrench");
         optgroup = page->new_optgroup(L("Extrusion width"));