diff --git a/src/libslic3r/Config.hpp b/src/libslic3r/Config.hpp index 151043a04..a8d0473b7 100644 --- a/src/libslic3r/Config.hpp +++ b/src/libslic3r/Config.hpp @@ -1786,6 +1786,17 @@ public: // Initialized by ConfigOptionEnum::get_enum_values() const t_config_enum_values *enum_keys_map = nullptr; + void set_enum_values(std::initializer_list> il) { + enum_values.clear(); + enum_values.reserve(il.size()); + enum_labels.clear(); + enum_labels.reserve(il.size()); + for (const std::pair p : il) { + enum_values.emplace_back(p.first); + enum_labels.emplace_back(p.second); + } + } + bool has_enum_value(const std::string &value) const { for (const std::string &v : enum_values) if (v == value) diff --git a/src/libslic3r/PrintConfig.cpp b/src/libslic3r/PrintConfig.cpp index 9a3bc3ce6..841bee5b9 100644 --- a/src/libslic3r/PrintConfig.cpp +++ b/src/libslic3r/PrintConfig.cpp @@ -138,7 +138,8 @@ CONFIG_OPTION_ENUM_DEFINE_STATIC_MAPS(SupportMaterialPattern) static const t_config_enum_values s_keys_map_SupportMaterialStyle { { "grid", smsGrid }, { "snug", smsSnug }, - { "tree", smsTree } + { "tree", smsTree }, + { "organic", smsOrganic } }; CONFIG_OPTION_ENUM_DEFINE_STATIC_MAPS(SupportMaterialStyle) @@ -2782,12 +2783,12 @@ void PrintConfigDef::init_fff_params() "will create more stable supports, while snug support towers will save material and reduce " "object scarring."); def->enum_keys_map = &ConfigOptionEnum::get_enum_values(); - def->enum_values.push_back("grid"); - def->enum_values.push_back("snug"); - def->enum_values.push_back("tree"); - def->enum_labels.push_back(L("Grid")); - def->enum_labels.push_back(L("Snug")); - def->enum_labels.push_back(L("Tree")); + def->set_enum_values({ + { "grid", L("Grid") }, + { "snug", L("Snug") }, + { "tree", L("Tree") }, + { "organic", L("Organic") } + }); def->mode = comAdvanced; def->set_default_value(new ConfigOptionEnum(smsGrid)); diff --git a/src/libslic3r/PrintConfig.hpp b/src/libslic3r/PrintConfig.hpp index 48493aa0b..0a6daf2d2 100644 --- a/src/libslic3r/PrintConfig.hpp +++ b/src/libslic3r/PrintConfig.hpp @@ -85,7 +85,7 @@ enum SupportMaterialPattern { }; enum SupportMaterialStyle { - smsGrid, smsSnug, smsTree, + smsGrid, smsSnug, smsTree, smsOrganic, }; enum SupportMaterialInterfacePattern { diff --git a/src/libslic3r/PrintObject.cpp b/src/libslic3r/PrintObject.cpp index 7ea55a3a4..9ffb57a72 100644 --- a/src/libslic3r/PrintObject.cpp +++ b/src/libslic3r/PrintObject.cpp @@ -2197,7 +2197,7 @@ void PrintObject::combine_infill() void PrintObject::_generate_support_material() { - if (m_config.support_material_style == smsTree) { + if (m_config.support_material_style == smsTree || m_config.support_material_style == smsOrganic) { fff_tree_support_generate(*this, std::function([this](){ this->throw_if_canceled(); })); } else { PrintObjectSupportMaterial support_material(this, m_slicing_params); diff --git a/src/libslic3r/SupportMaterial.cpp b/src/libslic3r/SupportMaterial.cpp index 491016887..b1b6600bc 100644 --- a/src/libslic3r/SupportMaterial.cpp +++ b/src/libslic3r/SupportMaterial.cpp @@ -800,6 +800,7 @@ public: { switch (m_style) { case smsTree: + case smsOrganic: assert(false); [[fallthrough]]; case smsGrid: diff --git a/src/libslic3r/TreeModelVolumes.cpp b/src/libslic3r/TreeModelVolumes.cpp index 41e78e45b..e098636fd 100644 --- a/src/libslic3r/TreeModelVolumes.cpp +++ b/src/libslic3r/TreeModelVolumes.cpp @@ -39,7 +39,7 @@ TreeSupportMeshGroupSettings::TreeSupportMeshGroupSettings(const PrintObject &pr // Support must be enabled and set to Tree style. assert(config.support_material); - assert(config.support_material_style == smsTree); + assert(config.support_material_style == smsTree || config.support_material_style == smsOrganic); // Calculate maximum external perimeter width over all printing regions, taking into account the default layer height. coordf_t external_perimeter_width = 0.; diff --git a/src/libslic3r/TreeSupport.cpp b/src/libslic3r/TreeSupport.cpp index c408b2f85..fecbef831 100644 --- a/src/libslic3r/TreeSupport.cpp +++ b/src/libslic3r/TreeSupport.cpp @@ -19,6 +19,9 @@ #include "Polyline.hpp" #include "MutablePolygon.hpp" #include "SupportMaterial.hpp" +#include "TriangleMeshSlicer.hpp" +#include "OpenVDBUtils.hpp" +#include #include #include @@ -26,6 +29,7 @@ #include #include #include +#include #ifdef _WIN32 #include //todo Remove! ONLY FOR PUBLIC BETA!! #endif // _WIN32 @@ -97,6 +101,20 @@ static inline void validate_range(const LineInformations &lines) validate_range(l); } +static inline void check_self_intersections(const Polygons &polygons, const std::string_view message) +{ +#ifdef _WIN32 + if (!intersecting_edges(polygons).empty()) + ::MessageBoxA(nullptr, (std::string("TreeSupport infill self intersections: ") + std::string(message)).c_str(), "Bug detected!", MB_OK | MB_SYSTEMMODAL | MB_SETFOREGROUND | MB_ICONWARNING); +#endif // _WIN32 +} +static inline void check_self_intersections(const ExPolygon &expoly, const std::string_view message) +{ +#ifdef _WIN32 + check_self_intersections(to_polygons(expoly), message); +#endif // _WIN32 +} + static inline void clip_for_diff(const Polygon &src, const BoundingBox &bbox, Polygon &out) { out.clear(); @@ -193,7 +211,7 @@ static std::vector>> group_me #endif // NDEBUG // Support must be enabled and set to Tree style. assert(object_config.support_material); - assert(object_config.support_material_style == smsTree); + assert(object_config.support_material_style == smsTree || object_config.support_material_style == smsOrganic); bool found_existing_group = false; TreeSupportSettings next_settings{ TreeSupportMeshGroupSettings{ print_object } }; @@ -312,17 +330,39 @@ void tree_supports_show_error(std::string message, bool critical) if (! (enforced_layer || blockers_layers.empty() || blockers_layers[layer_id].empty())) overhangs = diff(overhangs, blockers_layers[layer_id], ApplySafetyOffset::Yes); } - if (! enforcers_layers.empty() && ! enforcers_layers[layer_id].empty()) + //check_self_intersections(overhangs, "generate_overhangs1"); + if (! enforcers_layers.empty() && ! enforcers_layers[layer_id].empty()) { // Has some support enforcers at this layer, apply them to the overhangs, don't apply the support threshold angle. - if (Polygons enforced_overhangs = intersection(raw_overhangs_calculated ? raw_overhangs : diff(current_layer.lslices, lower_layer.lslices), enforcers_layers[layer_id]); + //enforcers_layers[layer_id] = union_(enforcers_layers[layer_id]); + //check_self_intersections(enforcers_layers[layer_id], "generate_overhangs - enforcers"); + //check_self_intersections(to_polygons(lower_layer.lslices), "generate_overhangs - lowerlayers"); + if (Polygons enforced_overhangs = intersection(raw_overhangs_calculated ? raw_overhangs : diff(current_layer.lslices, lower_layer.lslices), enforcers_layers[layer_id] /*, ApplySafetyOffset::Yes */); ! enforced_overhangs.empty()) { //FIXME this is a hack to make enforcers work on steep overhangs. - enforced_overhangs = diff(offset(enforced_overhangs, + //check_self_intersections(enforced_overhangs, "generate_overhangs - enforced overhangs1"); + //Polygons enforced_overhangs_prev = enforced_overhangs; + //check_self_intersections(to_polygons(union_ex(enforced_overhangs)), "generate_overhangs - enforced overhangs11"); + //check_self_intersections(offset(union_ex(enforced_overhangs), + //FIXME this is a fudge constant! + // scaled(0.4)), "generate_overhangs - enforced overhangs12"); + enforced_overhangs = diff(offset(union_ex(enforced_overhangs), //FIXME this is a fudge constant! scaled(0.4)), lower_layer.lslices); +#ifdef TREESUPPORT_DEBUG_SVG + if (! intersecting_edges(enforced_overhangs).empty()) { + static int irun = 0; + SVG::export_expolygons(debug_out_path("treesupport-self-intersections-%d.svg", ++irun), + { { { union_ex(enforced_overhangs_prev) }, { "prev", "yellow", 0.5f } }, + { { lower_layer.lslices }, { "lower_layer.lslices", "gray", 0.5f } }, + { { union_ex(enforced_overhangs) }, { "enforced_overhangs", "red", "black", "", scaled(0.1f), 0.5f } } }); + } +#endif // TREESUPPORT_DEBUG_SVG + //check_self_intersections(enforced_overhangs, "generate_overhangs - enforced overhangs2"); overhangs = overhangs.empty() ? std::move(enforced_overhangs) : union_(overhangs, enforced_overhangs); + //check_self_intersections(overhangs, "generate_overhangs - enforcers"); } + } out[layer_id] = std::move(overhangs); } }); @@ -675,7 +715,10 @@ static std::optional> polyline_sample_next_point_at_dis return lines; #else #ifdef _WIN32 - if (! BoundingBox(Point::new_scale(-170., -170.), Point::new_scale(170., 170.)).contains(get_extents(polygon))) + // Max dimensions for MK3 +// if (! BoundingBox(Point::new_scale(-170., -170.), Point::new_scale(170., 170.)).contains(get_extents(polygon))) + // Max dimensions for XL + if (! BoundingBox(Point::new_scale(-250., -250.), Point::new_scale(250., 250.)).contains(get_extents(polygon))) ::MessageBoxA(nullptr, "TreeSupport infill kravsky", "Bug detected!", MB_OK | MB_SYSTEMMODAL | MB_SETFOREGROUND | MB_ICONWARNING); #endif // _WIN32 @@ -702,10 +745,7 @@ static std::optional> polyline_sample_next_point_at_dis ::MessageBoxA(nullptr, "TreeSupport infill negative area", "Bug detected!", MB_OK | MB_SYSTEMMODAL | MB_SETFOREGROUND | MB_ICONWARNING); #endif // _WIN32 assert(intersecting_edges(to_polygons(expoly)).empty()); -#ifdef _WIN32 - if (! intersecting_edges(to_polygons(expoly)).empty()) - ::MessageBoxA(nullptr, "TreeSupport infill self intersections", "Bug detected!", MB_OK | MB_SYSTEMMODAL | MB_SETFOREGROUND | MB_ICONWARNING); -#endif // _WIN32 + check_self_intersections(expoly, "generate_support_infill_lines"); Surface surface(stInternal, std::move(expoly)); try { Polylines pl = filler->fill_surface(&surface, fill_params); @@ -831,6 +871,11 @@ static std::optional> polyline_sample_next_point_at_dis return union_(ret); } +static double layer_z(const SlicingParameters &slicing_params, const size_t layer_idx) +{ + return slicing_params.object_print_z_min + slicing_params.first_object_layer_height + layer_idx * slicing_params.layer_height; +} + static inline SupportGeneratorLayer& layer_initialize( SupportGeneratorLayer &layer_new, const SupporLayerType layer_type, @@ -838,7 +883,7 @@ static inline SupportGeneratorLayer& layer_initialize( const size_t layer_idx) { layer_new.layer_type = layer_type; - layer_new.print_z = slicing_params.object_print_z_min + slicing_params.first_object_layer_height + layer_idx * slicing_params.layer_height; + layer_new.print_z = layer_z(slicing_params, layer_idx); layer_new.height = layer_idx == 0 ? slicing_params.first_object_layer_height : slicing_params.layer_height; layer_new.bottom_z = layer_idx == 0 ? slicing_params.object_print_z_min : layer_new.print_z - layer_new.height; return layer_new; @@ -1082,6 +1127,8 @@ static void generate_initial_areas( overhang_regular = mesh_group_settings.support_offset == 0 ? overhang_raw : safe_offset_inc(overhang_raw, mesh_group_settings.support_offset, relevant_forbidden, mesh_config.min_radius * 1.75 + mesh_config.xy_min_distance, 0, 1); + //check_self_intersections(overhang_regular, "overhang_regular1"); + // offset ensures that areas that could be supported by a part of a support line, are not considered unsupported overhang Polygons remaining_overhang = intersection( diff(mesh_group_settings.support_offset == 0 ? @@ -1108,6 +1155,7 @@ static void generate_initial_areas( remaining_overhang = diff(remaining_overhang, safe_offset_inc(overhang_regular, 1.5 * extra_total_offset_acc, raw_collision, offset_step, 0, 1)); // Extending the overhangs by the inflated remaining overhangs. overhang_regular = union_(overhang_regular, diff(safe_offset_inc(remaining_overhang, extra_total_offset_acc, raw_collision, offset_step, 0, 1), relevant_forbidden)); + //check_self_intersections(overhang_regular, "overhang_regular2"); } // If the xy distance overrides the z distance, some support needs to be inserted further down. //=> Analyze which support points do not fit on this layer and check if they will fit a few layers down (while adding them an infinite amount of layers down would technically be closer the the setting description, it would not produce reasonable results. ) @@ -1159,6 +1207,7 @@ static void generate_initial_areas( if (mesh_group_settings.minimum_support_area > 0) remove_small(overhang_roofs, mesh_group_settings.minimum_roof_area); overhang_regular = diff(overhang_regular, overhang_roofs, ApplySafetyOffset::Yes); + //check_self_intersections(overhang_regular, "overhang_regular3"); for (ExPolygon &roof_part : union_ex(overhang_roofs)) overhang_processing.emplace_back(std::move(roof_part), true); } @@ -2369,6 +2418,8 @@ static void set_points_on_areas(const SupportElement &elem, SupportElements *lay next_elem.state.result_on_layer = move_inside_if_outside(next_elem.influence_area, elem.state.result_on_layer); // do not call recursive because then amount of layers would be restricted by the stack size } + // Mark the parent element as accessed from a valid child element. + next_elem.state.marked = true; } } @@ -2487,15 +2538,23 @@ static void create_nodes_from_area( { // Initialize points on layer 0, with a "random" point in the influence area. // Point is chosen based on an inaccurate estimate where the branches will split into two, but every point inside the influence area would produce a valid result. - for (SupportElement &init : move_bounds.front()) { - init.state.result_on_layer = move_inside_if_outside(init.influence_area, init.state.next_position); - // Also set the parent nodes, as these will be required for the first iteration of the loop below. - set_points_on_areas(init, move_bounds.size() > 1 ? &move_bounds[1] : nullptr); + { + SupportElements *layer_above = move_bounds.size() > 1 ? &move_bounds[1] : nullptr; + for (SupportElement &elem : *layer_above) + elem.state.marked = false; + for (SupportElement &init : move_bounds.front()) { + init.state.result_on_layer = move_inside_if_outside(init.influence_area, init.state.next_position); + // Also set the parent nodes, as these will be required for the first iteration of the loop below and mark the parent nodes. + set_points_on_areas(init, layer_above); + } } for (LayerIndex layer_idx = 1; layer_idx < LayerIndex(move_bounds.size()); ++ layer_idx) { auto &layer = move_bounds[layer_idx]; auto *layer_above = layer_idx + 1 < move_bounds.size() ? &move_bounds[layer_idx + 1] : nullptr; + if (layer_above) + for (SupportElement &elem : *layer_above) + elem.state.marked = false; for (SupportElement &elem : layer) { assert(! elem.state.deleted); assert(elem.state.layer_idx == layer_idx); @@ -2509,11 +2568,6 @@ static void create_nodes_from_area( } // we dont need to remove yet the parents as they will have a lower dtt and also no result_on_layer set elem.state.deleted = true; - for (int32_t parent_idx : elem.parents) - // When the roof was not able to generate downwards enough, the top elements may have not moved, and have result_on_layer already set. - // As this branch needs to be removed => all parents result_on_layer have to be invalidated. - (*layer_above)[parent_idx].state.result_on_layer_reset(); - continue; } else { // set the point where the branch will be placed on the model if (elem.state.to_model_gracious) @@ -2522,13 +2576,67 @@ static void create_nodes_from_area( set_to_model_contact_simple(elem); } } - if (! elem.state.deleted) - // element is valid now setting points in the layer above + if (! elem.state.deleted && ! elem.state.marked && elem.state.target_height == layer_idx) + // Just a tip surface with no supporting element. + elem.state.deleted = true; + if (elem.state.deleted) { + for (int32_t parent_idx : elem.parents) + // When the roof was not able to generate downwards enough, the top elements may have not moved, and have result_on_layer already set. + // As this branch needs to be removed => all parents result_on_layer have to be invalidated. + (*layer_above)[parent_idx].state.result_on_layer_reset(); + } + if (! elem.state.deleted) { + // Element is valid now setting points in the layer above and mark the parent nodes. set_points_on_areas(elem, layer_above); + } } } +#ifndef NDEBUG + // Verify the tree connectivity including the branch slopes. + for (LayerIndex layer_idx = 0; layer_idx + 1 < LayerIndex(move_bounds.size()); ++ layer_idx) { + auto &layer = move_bounds[layer_idx]; + auto &above = move_bounds[layer_idx + 1]; + for (SupportElement &elem : layer) + if (! elem.state.deleted) { + for (int32_t iparent : elem.parents) { + SupportElement &parent = above[iparent]; + assert(! parent.state.deleted); + assert(elem.state.result_on_layer_is_set() == parent.state.result_on_layer_is_set()); + if (elem.state.result_on_layer_is_set()) { + double radius_increase = config.getRadius(elem.state) - config.getRadius(parent.state); + assert(radius_increase >= 0); + double shift = (elem.state.result_on_layer - parent.state.result_on_layer).cast().norm(); + assert(shift < radius_increase + 2. * config.maximum_move_distance_slow); + } + } + } + } +#endif // NDEBUG + remove_deleted_elements(move_bounds); + +#ifndef NDEBUG + // Verify the tree connectivity including the branch slopes. + for (LayerIndex layer_idx = 0; layer_idx + 1 < LayerIndex(move_bounds.size()); ++ layer_idx) { + auto &layer = move_bounds[layer_idx]; + auto &above = move_bounds[layer_idx + 1]; + for (SupportElement &elem : layer) { + assert(! elem.state.deleted); + for (int32_t iparent : elem.parents) { + SupportElement &parent = above[iparent]; + assert(! parent.state.deleted); + assert(elem.state.result_on_layer_is_set() == parent.state.result_on_layer_is_set()); + if (elem.state.result_on_layer_is_set()) { + double radius_increase = config.getRadius(elem.state) - config.getRadius(parent.state); + assert(radius_increase >= 0); + double shift = (elem.state.result_on_layer - parent.state.result_on_layer).cast().norm(); + assert(shift < radius_increase + 2. * config.maximum_move_distance_slow); + } + } + } + } +#endif // NDEBUG } // For producing circular / elliptical areas from SupportElements (one DrawArea per one SupportElement) @@ -2649,7 +2757,8 @@ static void generate_branch_areas(const TreeModelVolumes &volumes, const TreeSup polygons_with_correct_center.emplace_back(std::move(part)); } // Increase the area again, to ensure the nozzle path when calculated later is very similar to the one assumed above. - polygons = diff_clipped(offset(polygons_with_correct_center, config.support_line_width / 2, jtMiter, 1.2), + assert(contains(polygons, draw_area.element->state.result_on_layer)); + polygons = diff_clipped(offset(polygons_with_correct_center, config.support_line_width / 2, jtMiter, 1.2), //FIXME Vojtech: Clipping may split the region into multiple pieces again, reversing the fixing effort. collision); } @@ -2695,10 +2804,12 @@ static void smooth_branch_areas( [&](const tbb::blocked_range &range) { for (size_t processing_idx = range.begin(); processing_idx < range.end(); ++ processing_idx) { DrawArea &draw_area = linear_data[processing_base + processing_idx]; + assert(draw_area.element->state.layer_idx == layer_idx); double max_outer_wall_distance = 0; bool do_something = false; for (int32_t parent_idx : draw_area.element->parents) { const SupportElement &parent = layer_above[parent_idx]; + assert(parent.state.layer_idx == layer_idx + 1); if (config.getRadius(parent.state) != config.getCollisionRadius(parent.state)) { do_something = true; max_outer_wall_distance = std::max(max_outer_wall_distance, (draw_area.element->state.result_on_layer - parent.state.result_on_layer).cast().norm() - (config.getRadius(*draw_area.element) - config.getRadius(parent))); @@ -2706,14 +2817,35 @@ static void smooth_branch_areas( } max_outer_wall_distance += max_radius_change_per_layer; // As this change is a bit larger than what usually appears, lost radius can be slowly reclaimed over the layers. if (do_something) { + assert(contains(draw_area.polygons, draw_area.element->state.result_on_layer)); Polygons max_allowed_area = offset(draw_area.polygons, float(max_outer_wall_distance), jtMiter, 1.2); for (int32_t parent_idx : draw_area.element->parents) { const SupportElement &parent = layer_above[parent_idx]; +#ifndef NDEBUG + assert(parent.state.layer_idx == layer_idx + 1); + assert(contains(linear_data[processing_base_above + parent_idx].polygons, parent.state.result_on_layer)); + double radius_increase = config.getRadius(draw_area.element->state) - config.getRadius(parent.state); + assert(radius_increase >= 0); + double shift = (draw_area.element->state.result_on_layer - parent.state.result_on_layer).cast().norm(); + assert(shift < radius_increase + 2. * config.maximum_move_distance_slow); +#endif // NDEBUG if (config.getRadius(parent.state) != config.getCollisionRadius(parent.state)) { // No other element on this layer than the current one may be connected to &parent, // thus it is safe to update parent's DrawArea directly. Polygons &dst = linear_data[processing_base_above + parent_idx].polygons; - dst = intersection(dst, max_allowed_area); +// Polygons orig = dst; + if (! dst.empty()) { + dst = intersection(dst, max_allowed_area); +#if 0 + if (dst.empty()) { + static int irun = 0; + SVG::export_expolygons(debug_out_path("treesupport-extrude_areas-smooth-error-%d.svg", irun ++), + { { { union_ex(max_allowed_area) }, { "max_allowed_area", "yellow", 0.5f } }, + { { union_ex(orig) }, { "orig", "red", "black", "", scaled(0.1f), 0.5f } } }); + ::MessageBoxA(nullptr, "TreeSupport smoothing bug", "Bug detected!", MB_OK | MB_SYSTEMMODAL | MB_SETFOREGROUND | MB_ICONWARNING); + } +#endif + } } } } @@ -2982,9 +3114,7 @@ static void draw_areas( // Only one link points to a node above from below. assert(! (++ it != map_downwards_old.end() && it->first == &elem)); } - if ((! child && elem.state.target_height == layer_idx) || (child && !child->state.result_on_layer_is_set())) - // We either come from nowhere at the final layer or we had invalid parents 2. should never happen but just to be sure - continue; + assert(child ? child->state.result_on_layer_is_set() : elem.state.target_height > layer_idx); } for (int32_t parent_idx : elem.parents) { SupportElement &parent = (*layer_above)[parent_idx]; @@ -2998,12 +3128,66 @@ static void draw_areas( linear_data_layers.emplace_back(linear_data.size()); } +#ifndef NDEBUG + for (size_t i = 0; i < move_bounds.size(); ++ i) { + size_t begin = linear_data_layers[i]; + size_t end = linear_data_layers[i + 1]; + for (size_t j = begin; j < end; ++ j) + assert(linear_data[j].element == &move_bounds[i][j - begin]); + } +#endif // NDEBUG + auto t_start = std::chrono::high_resolution_clock::now(); // Generate the circles that will be the branches. generate_branch_areas(volumes, config, move_bounds, linear_data); + +#if 0 + assert(linear_data_layers.size() == move_bounds.size() + 1); + for (const auto &draw_area : linear_data) + assert(contains(draw_area.polygons, draw_area.element->state.result_on_layer)); + for (size_t i = 0; i < move_bounds.size(); ++ i) { + size_t begin = linear_data_layers[i]; + size_t end = linear_data_layers[i + 1]; + for (size_t j = begin; j < end; ++ j) { + const auto &draw_area = linear_data[j]; + assert(draw_area.element == &move_bounds[i][j - begin]); + assert(contains(draw_area.polygons, draw_area.element->state.result_on_layer)); + } + } +#endif + +#if 0 + for (size_t area_layer_idx = 0; area_layer_idx + 1 < linear_data_layers.size(); ++ area_layer_idx) { + size_t begin = linear_data_layers[area_layer_idx]; + size_t end = linear_data_layers[area_layer_idx + 1]; + Polygons polygons; + for (size_t area_idx = begin; area_idx < end; ++ area_idx) { + DrawArea &area = linear_data[area_idx]; + append(polygons, area.polygons); + } + SVG::export_expolygons(debug_out_path("treesupport-extrude_areas-raw-%d.svg", area_layer_idx), + { { { union_ex(polygons) }, { "parent", "red", "black", "", scaled(0.1f), 0.5f } } }); + } +#endif + auto t_generate = std::chrono::high_resolution_clock::now(); // In some edgecases a branch may go though a hole, where the regular radius does not fit. This can result in an apparent jump in branch radius. As such this cases need to be caught and smoothed out. smooth_branch_areas(config, move_bounds, linear_data, linear_data_layers); + +#if 0 + for (size_t area_layer_idx = 0; area_layer_idx + 1 < linear_data_layers.size(); ++area_layer_idx) { + size_t begin = linear_data_layers[area_layer_idx]; + size_t end = linear_data_layers[area_layer_idx + 1]; + Polygons polygons; + for (size_t area_idx = begin; area_idx < end; ++area_idx) { + DrawArea& area = linear_data[area_idx]; + append(polygons, area.polygons); + } + SVG::export_expolygons(debug_out_path("treesupport-extrude_areas-smooth-%d.svg", area_layer_idx), + { { { union_ex(polygons) }, { "parent", "red", "black", "", scaled(0.1f), 0.5f } } }); + } +#endif + auto t_smooth = std::chrono::high_resolution_clock::now(); // drop down all trees that connect non gracefully with the model drop_non_gracious_areas(volumes, linear_data, support_layer_storage); @@ -3044,6 +3228,490 @@ static void draw_areas( "finalize_interface_and_support_areas " << dur_finalize << " ms"; } +#if 1 +// Test whether two circles, each on its own plane in 3D intersect. +// Circles are considered intersecting, if the lowest point on one circle is below the other circle's plane. +// Assumption: The two planes are oriented the same way. +static bool circles_intersect( + const Vec3d &p1, const Vec3d &n1, const double r1, + const Vec3d &p2, const Vec3d &n2, const double r2) +{ + assert(n1.dot(n2) >= 0); + + const Vec3d z = n1.cross(n2); + const Vec3d dir1 = z.cross(n1); + const Vec3d lowest_point1 = p1 + dir1 * (r1 / dir1.norm()); + assert(n2.dot(p1) >= n2.dot(lowest_point1)); + if (n2.dot(lowest_point1) <= 0) + return true; + const Vec3d dir2 = z.cross(n2); + const Vec3d lowest_point2 = p2 + dir2 * (r2 / dir2.norm()); + assert(n1.dot(p2) >= n1.dot(lowest_point2)); + return n1.dot(lowest_point2) <= 0; +} + +template +void triangulate_fan(indexed_triangle_set &its, int ifan, int ibegin, int iend) +{ + // at least 3 vertices, increasing order. + assert(ibegin + 3 <= iend); + assert(ibegin >= 0 && iend <= its.vertices.size()); + assert(ifan >= 0 && ifan < its.vertices.size()); + int num_faces = iend - ibegin; + its.indices.reserve(its.indices.size() + num_faces * 3); + for (int v = ibegin, u = iend - 1; v < iend; u = v ++) { + if (flip_normals) + its.indices.push_back({ ifan, u, v }); + else + its.indices.push_back({ ifan, v, u }); + } +} + +static void triangulate_strip(indexed_triangle_set &its, int ibegin1, int iend1, int ibegin2, int iend2) +{ + // at least 3 vertices, increasing order. + assert(ibegin1 + 3 <= iend1); + assert(ibegin1 >= 0 && iend1 <= its.vertices.size()); + assert(ibegin2 + 3 <= iend2); + assert(ibegin2 >= 0 && iend2 <= its.vertices.size()); + int n1 = iend1 - ibegin1; + int n2 = iend2 - ibegin2; + its.indices.reserve(its.indices.size() + (n1 + n2) * 3); + + // For the first vertex of 1st strip, find the closest vertex on the 2nd strip. + int istart2 = ibegin2; + { + const Vec3f &p1 = its.vertices[ibegin1]; + auto d2min = std::numeric_limits::max(); + for (int i = ibegin2; i < iend2; ++ i) { + const Vec3f &p2 = its.vertices[i]; + const float d2 = (p2 - p1).squaredNorm(); + if (d2 < d2min) { + d2min = d2; + istart2 = i; + } + } + } + + // Now triangulate the strip zig-zag fashion taking always the shortest connection if possible. + for (int u = ibegin1, v = istart2; n1 > 0 || n2 > 0;) { + bool take_first; + int u2, v2; + auto update_u2 = [&u2, u, ibegin1, iend1]() { + u2 = u; + if (++ u2 == iend1) + u2 = ibegin1; + }; + auto update_v2 = [&v2, v, ibegin2, iend2]() { + v2 = v; + if (++ v2 == iend2) + v2 = ibegin2; + }; + if (n1 == 0) { + take_first = false; + update_v2(); + } else if (n2 == 0) { + take_first = true; + update_u2(); + } else { + update_u2(); + update_v2(); + float l1 = (its.vertices[u2] - its.vertices[v]).squaredNorm(); + float l2 = (its.vertices[v2] - its.vertices[u]).squaredNorm(); + take_first = l1 < l2; + } + if (take_first) { + its.indices.push_back({ u, u2, v }); + -- n1; + u = u2; + } else { + its.indices.push_back({ u, v2, v }); + -- n2; + v = v2; + } + } +} + +// Discretize 3D circle, append to output vector, return ranges of indices of the points added. +static std::pair discretize_circle(const Vec3f ¢er, const Vec3f &normal, const float radius, const float eps, std::vector &pts) +{ + // Calculate discretization step and number of steps. + float angle_step = 2. * acos(1. - eps / radius); + auto nsteps = int(ceil(2 * M_PI / angle_step)); + angle_step = 2 * M_PI / nsteps; + + // Prepare coordinate system for the circle plane. + Vec3f x = normal.cross(Vec3f(0.f, -1.f, 0.f)).normalized(); + Vec3f y = normal.cross(x).normalized(); + assert(std::abs(x.cross(y).dot(normal) - 1.f) < EPSILON); + + // Discretize the circle. + int begin = int(pts.size()); + pts.reserve(pts.size() + nsteps); + float angle = 0; + x *= radius; + y *= radius; + for (int i = 0; i < nsteps; ++ i) { + pts.emplace_back(center + x * cos(angle) + y * sin(angle)); + angle += angle_step; + } + return { begin, int(pts.size()) }; +} + +static void extrude_branch( + const std::vector &path, + const TreeSupportSettings &config, + const SlicingParameters &slicing_params, + const std::vector &move_bounds, + indexed_triangle_set &result) +{ + Vec3d p1, p2, p3; + Vec3d v1, v2; + Vec3d nprev; + Vec3d ncurrent; + assert(path.size() >= 2); + static constexpr const float eps = 0.015f; + std::pair prev_strip; + +// char fname[2048]; +// static int irun = 0; + + for (size_t ipath = 1; ipath < path.size(); ++ ipath) { + const SupportElement &prev = *path[ipath - 1]; + const SupportElement ¤t = *path[ipath]; + assert(prev.state.layer_idx + 1 == current.state.layer_idx); + p1 = to_3d(unscaled(prev .state.result_on_layer), layer_z(slicing_params, prev .state.layer_idx)); + p2 = to_3d(unscaled(current.state.result_on_layer), layer_z(slicing_params, current.state.layer_idx)); + v1 = (p2 - p1).normalized(); + if (ipath == 1) { + nprev = v1; + // Extrude the bottom half sphere. + float radius = unscaled(config.getRadius(prev.state)); + float angle_step = 2. * acos(1. - eps / radius); + auto nsteps = int(ceil(M_PI / (2. * angle_step))); + angle_step = M_PI / (2. * nsteps); + int ifan = int(result.vertices.size()); + result.vertices.emplace_back((p1 - nprev * radius).cast()); + float angle = angle_step; + for (int i = 1; i < nsteps; ++ i, angle += angle_step) { + std::pair strip = discretize_circle((p1 - nprev * radius * cos(angle)).cast(), nprev.cast(), radius * sin(angle), eps, result.vertices); + if (i == 1) + triangulate_fan(result, ifan, strip.first, strip.second); + else + triangulate_strip(result, prev_strip.first, prev_strip.second, strip.first, strip.second); +// sprintf(fname, "d:\\temp\\meshes\\tree-partial-%d.obj", ++ irun); +// its_write_obj(result, fname); + prev_strip = strip; + } + } + if (ipath + 1 == path.size()) { + // End of the tube. + ncurrent = v1; + // Extrude the top half sphere. + float radius = unscaled(config.getRadius(current.state)); + float angle_step = 2. * acos(1. - eps / radius); + auto nsteps = int(ceil(M_PI / (2. * angle_step))); + angle_step = M_PI / (2. * nsteps); + auto angle = float(M_PI / 2.); + for (int i = 0; i < nsteps; ++ i, angle -= angle_step) { + std::pair strip = discretize_circle((p2 + ncurrent * radius * cos(angle)).cast(), ncurrent.cast(), radius * sin(angle), eps, result.vertices); + triangulate_strip(result, prev_strip.first, prev_strip.second, strip.first, strip.second); +// sprintf(fname, "d:\\temp\\meshes\\tree-partial-%d.obj", ++ irun); +// its_write_obj(result, fname); + prev_strip = strip; + } + int ifan = int(result.vertices.size()); + result.vertices.emplace_back((p2 + ncurrent * radius).cast()); + triangulate_fan(result, ifan, prev_strip.first, prev_strip.second); +// sprintf(fname, "d:\\temp\\meshes\\tree-partial-%d.obj", ++ irun); +// its_write_obj(result, fname); + } else { + const SupportElement &next = *path[ipath + 1]; + assert(current.state.layer_idx + 1 == next.state.layer_idx); + p3 = to_3d(unscaled(next.state.result_on_layer), layer_z(slicing_params, next.state.layer_idx)); + v2 = (p3 - p2).normalized(); + ncurrent = (v1 + v2).normalized(); + float radius = unscaled(config.getRadius(current.state)); + std::pair strip = discretize_circle(p2.cast(), ncurrent.cast(), radius, eps, result.vertices); + triangulate_strip(result, prev_strip.first, prev_strip.second, strip.first, strip.second); + prev_strip = strip; +// sprintf(fname, "d:\\temp\\meshes\\tree-partial-%d.obj", ++irun); +// its_write_obj(result, fname); + } +#if 0 + if (circles_intersect(p1, nprev, settings.getRadius(prev), p2, ncurrent, settings.getRadius(current))) { + // Cannot connect previous and current slice using a simple zig-zag triangulation, + // because the two circles intersect. + + } else { + // Continue with chaining. + + } +#endif + } +} +#endif + +static void draw_branches( + PrintObject &print_object, + const TreeModelVolumes &volumes, + const TreeSupportSettings &config, + const std::vector &overhangs, + std::vector &move_bounds, + + SupportGeneratorLayersPtr &bottom_contacts, + SupportGeneratorLayersPtr &top_contacts, + SupportGeneratorLayersPtr &intermediate_layers, + SupportGeneratorLayerStorage &layer_storage) +{ + static int irun = 0; + + const SlicingParameters& slicing_params = print_object.slicing_parameters(); + + // All SupportElements are put into a layer independent storage to improve parallelization. + std::vector> elements_with_link_down; + std::vector linear_data_layers; + { + std::vector> map_downwards_old; + std::vector> map_downwards_new; + linear_data_layers.emplace_back(0); + for (LayerIndex layer_idx = 0; layer_idx < LayerIndex(move_bounds.size()); ++ layer_idx) { + SupportElements *layer_above = layer_idx + 1 < move_bounds.size() ? &move_bounds[layer_idx + 1] : nullptr; + map_downwards_new.clear(); + std::sort(map_downwards_old.begin(), map_downwards_old.end(), [](auto& l, auto& r) { return l.first < r.first; }); + SupportElements &layer = move_bounds[layer_idx]; + for (size_t elem_idx = 0; elem_idx < layer.size(); ++ elem_idx) { + SupportElement &elem = layer[elem_idx]; + int child = -1; + if (layer_idx > 0) { + auto it = std::lower_bound(map_downwards_old.begin(), map_downwards_old.end(), &elem, [](auto& l, const SupportElement* r) { return l.first < r; }); + if (it != map_downwards_old.end() && it->first == &elem) { + child = it->second; + // Only one link points to a node above from below. + assert(!(++it != map_downwards_old.end() && it->first == &elem)); + } + const SupportElement *pchild = child == -1 ? nullptr : &move_bounds[layer_idx - 1][child]; + assert(pchild ? pchild->state.result_on_layer_is_set() : elem.state.target_height > layer_idx); + } + for (int32_t parent_idx : elem.parents) { + SupportElement &parent = (*layer_above)[parent_idx]; + if (parent.state.result_on_layer_is_set()) + map_downwards_new.emplace_back(&parent, elem_idx); + } + + elements_with_link_down.push_back({ &elem, int(child) }); + } + std::swap(map_downwards_old, map_downwards_new); + linear_data_layers.emplace_back(elements_with_link_down.size()); + } + } + + std::unique_ptr> closest_surface_point; + { + TriangleMesh mesh = print_object.model_object()->raw_mesh(); + mesh.transform(print_object.trafo_centered()); + double scale = 10.; + openvdb::FloatGrid::Ptr grid = mesh_to_grid(mesh.its, {}, scale, 0., 0.); + closest_surface_point = openvdb::tools::ClosestSurfacePoint::create(*grid); + std::vector pts, prev, projections; + std::vector distances; + for (const std::pair &element : elements_with_link_down) { + Vec3d pt = to_3d(unscaled(element.first->state.result_on_layer), layer_z(slicing_params, element.first->state.layer_idx)) * scale; + pts.push_back({ pt.x(), pt.y(), pt.z() }); + } + + const double collision_extra_gap = 1. * scale; + const double max_nudge_collision_avoidance = 2. * scale; + const double max_nudge_smoothing = 1. * scale; + + static constexpr const size_t num_iter = 100; // 1000; + for (size_t iter = 0; iter < num_iter; ++ iter) { + prev = pts; + projections = pts; + distances.assign(pts.size(), std::numeric_limits::max()); + closest_surface_point->searchAndReplace(projections, distances); + size_t num_moved = 0; + for (size_t i = 0; i < projections.size(); ++ i) { + const SupportElement &element = *elements_with_link_down[i].first; + const int below = elements_with_link_down[i].second; + if (pts[i] != projections[i]) { + // Nudge the circle center away from the collision. + Vec3d v{ projections[i].x() - pts[i].x(), projections[i].y() - pts[i].y(), projections[i].z() - pts[i].z() }; + double depth = v.norm(); + assert(std::abs(distances[i] - depth) < EPSILON); + double radius = unscaled(config.getRadius(element.state)) * scale; + if (depth < radius) { + // Collision detected to be removed. + ++ num_moved; + double dxy = sqrt(sqr(radius) - sqr(v.z())); + double nudge_dist_max = dxy - std::hypot(v.x(), v.y()) + //FIXME 1mm gap + + collision_extra_gap; + // Shift by maximum 2mm. + double nudge_dist = std::min(std::max(0., nudge_dist_max), max_nudge_collision_avoidance); + Vec2d nudge_v = to_2d(v).normalized() * (- nudge_dist); + pts[i].x() += nudge_v.x(); + pts[i].y() += nudge_v.y(); + } + } + // Laplacian smoothing + if (! element.parents.empty() && (below != -1 || element.state.layer_idx == 0)) { + Vec2d avg{ 0, 0 }; + const SupportElements &above = move_bounds[element.state.layer_idx + 1]; + const size_t offset_above = linear_data_layers[element.state.layer_idx + 1]; + double weight = 0.; + for (auto iparent : element.parents) { + double w = config.getRadius(above[iparent].state); + avg.x() += w * prev[offset_above + iparent].x(); + avg.y() += w * prev[offset_above + iparent].y(); + weight += w; + } + size_t cnt = element.parents.size(); + if (below != -1) { + const size_t offset_below = linear_data_layers[element.state.layer_idx - 1]; + const double w = weight; // config.getRadius(move_bounds[element.state.layer_idx - 1][below].state); + avg.x() += w * prev[offset_below + below].x(); + avg.y() += w * prev[offset_below + below].y(); + ++ cnt; + weight += w; + } + //avg /= double(cnt); + avg /= weight; + static constexpr const double smoothing_factor = 0.5; + Vec2d old_pos{ pts[i].x(), pts[i].y() }; + Vec2d new_pos = (1. - smoothing_factor) * old_pos + smoothing_factor * avg; + Vec2d shift = new_pos - old_pos; + double nudge_dist_max = shift.norm(); + // Shift by maximum 1mm, less than the collision avoidance factor. + double nudge_dist = std::min(std::max(0., nudge_dist_max), max_nudge_smoothing); + Vec2d nudge_v = shift.normalized() * nudge_dist; + pts[i].x() += nudge_v.x(); + pts[i].y() += nudge_v.y(); + } + } + printf("iteration: %d, moved: %d\n", int(iter), int(num_moved)); + if (num_moved == 0) + break; + } + +#if 1 + for (size_t i = 0; i < projections.size(); ++ i) { + elements_with_link_down[i].first->state.result_on_layer.x() = scaled(pts[i].x()) / scale; + elements_with_link_down[i].first->state.result_on_layer.y() = scaled(pts[i].y()) / scale; + } +#endif + } + + std::vector support_layer_storage(move_bounds.size()); + std::vector support_roof_storage(move_bounds.size()); + + // Unmark all nodes. + for (SupportElements &elements : move_bounds) + for (SupportElement &element : elements) + element.state.marked = false; + + // Traverse all nodes, generate tubes. + // Traversal stack with nodes and thier current parent + std::vector path; + indexed_triangle_set cummulative_mesh; + indexed_triangle_set partial_mesh; + indexed_triangle_set temp_mesh; + for (LayerIndex layer_idx = 0; layer_idx + 1 < LayerIndex(move_bounds.size()); ++ layer_idx) { + SupportElements &layer = move_bounds[layer_idx]; + SupportElements &layer_above = move_bounds[layer_idx + 1]; + + for (SupportElement &start_element : layer) + if (! start_element.state.marked && ! start_element.parents.empty()) { + // Collect elements up to a bifurcation above. + start_element.state.marked = true; + for (size_t parent_idx = 0; parent_idx < start_element.parents.size(); ++ parent_idx) { + path.clear(); + path.emplace_back(&start_element); + // Traverse each branch until it branches again. + SupportElement &first_parent = layer_above[start_element.parents[parent_idx]]; + assert(path.back()->state.layer_idx + 1 == first_parent.state.layer_idx); + path.emplace_back(&first_parent); + if (first_parent.parents.size() < 2) + first_parent.state.marked = true; + if (first_parent.parents.size() == 1) { + for (SupportElement *parent = &first_parent;;) { + SupportElement &next_parent = move_bounds[parent->state.layer_idx + 1][parent->parents.front()]; + assert(path.back()->state.layer_idx + 1 == next_parent.state.layer_idx); + path.emplace_back(&next_parent); + if (next_parent.parents.size() > 1) + break; + next_parent.state.marked = true; + if (next_parent.parents.size() == 0) + break; + parent = &next_parent; + } + } + // Triangulate the tube. + partial_mesh.clear(); + extrude_branch(path, config, slicing_params, move_bounds, partial_mesh); +#if 0 + char fname[2048]; + sprintf(fname, "d:\\temp\\meshes\\tree-raw-%d.obj", ++ irun); + its_write_obj(partial_mesh, fname); +#if 0 + temp_mesh.clear(); + cut_mesh(partial_mesh, layer_z(slicing_params, path.back()->state.layer_idx) + EPSILON, nullptr, &temp_mesh, false); + sprintf(fname, "d:\\temp\\meshes\\tree-trimmed1-%d.obj", irun); + its_write_obj(temp_mesh, fname); + partial_mesh.clear(); + cut_mesh(temp_mesh, layer_z(slicing_params, path.front()->state.layer_idx) - EPSILON, &partial_mesh, nullptr, false); + sprintf(fname, "d:\\temp\\meshes\\tree-trimmed2-%d.obj", irun); +#endif + its_write_obj(partial_mesh, fname); +#endif + its_merge(cummulative_mesh, partial_mesh); + } + } + } + + std::vector slice_z; + for (size_t layer_idx = 0; layer_idx < move_bounds.size(); ++ layer_idx) { + double print_z = slicing_params.object_print_z_min + slicing_params.first_object_layer_height + layer_idx * slicing_params.layer_height; + double layer_height = layer_idx == 0 ? slicing_params.first_object_layer_height : slicing_params.layer_height; + slice_z.emplace_back(float(print_z - layer_height * 0.5)); + } + // Remove the trailing slices. + while (! slice_z.empty()) + if (move_bounds[slice_z.size() - 1].empty()) + slice_z.pop_back(); + else + break; + +#if 0 + its_write_obj(cummulative_mesh, "d:\\temp\\meshes\\tree.obj"); +#endif + + MeshSlicingParamsEx params; + params.closing_radius = float(print_object.config().slice_closing_radius.value); + params.mode = MeshSlicingParams::SlicingMode::Positive; + std::vector slices = slice_mesh_ex(cummulative_mesh, slice_z, params); + for (size_t layer_idx = 0; layer_idx < slice_z.size(); ++ layer_idx) + if (! slices[layer_idx].empty()) { + SupportGeneratorLayer *&l = intermediate_layers[layer_idx]; + if (l == nullptr) + l = &layer_allocate(layer_storage, SupporLayerType::Base, slicing_params, layer_idx); + append(l->polygons, to_polygons(std::move(slices[layer_idx]))); + } + + // Trim the slices. + tbb::parallel_for(tbb::blocked_range(0, intermediate_layers.size()), + [&](const tbb::blocked_range &range) { + for (size_t layer_idx = range.begin(); layer_idx < range.end(); ++ layer_idx) + if (SupportGeneratorLayer *layer = intermediate_layers[layer_idx]; layer) { + Polygons &poly = intermediate_layers[layer_idx]->polygons; + poly = diff_clipped(poly, volumes.getCollision(0, layer_idx, true)); + } + }); + + finalize_interface_and_support_areas(print_object, volumes, config, overhangs, support_layer_storage, support_roof_storage, + bottom_contacts, top_contacts, intermediate_layers, layer_storage); +} + /*! * \brief Create the areas that need support. * @@ -3147,8 +3815,15 @@ static void generate_support_areas(Print &print, const BuildVolume &build_volume auto t_place = std::chrono::high_resolution_clock::now(); // ### draw these points as circles - draw_areas(*print.get_object(processing.second.front()), volumes, config, overhangs, move_bounds, - bottom_contacts, top_contacts, intermediate_layers, layer_storage); + + if (print_object.config().support_material_style == smsTree) + draw_areas(*print.get_object(processing.second.front()), volumes, config, overhangs, move_bounds, + bottom_contacts, top_contacts, intermediate_layers, layer_storage); + else { + assert(print_object.config().support_material_style == smsOrganic); + draw_branches(*print.get_object(processing.second.front()), volumes, config, overhangs, move_bounds, + bottom_contacts, top_contacts, intermediate_layers, layer_storage); + } auto t_draw = std::chrono::high_resolution_clock::now(); auto dur_pre_gen = 0.001 * std::chrono::duration_cast(t_precalc - t_start).count();