3DPrinting/PrusaSlicer-NonPlainar
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity

Refactored mesh slicing code into a new TriangleMeshSlicer class

Browse source
This commit is contained in:
Alessandro Ranellucci 2014-01-15 20:31:38 +01:00
parent dfd9bc8958
commit 519ed91c68
4 changed files with 479 additions and 448 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

902
xs/src/TriangleMesh.cpp
View file

@ -164,450 +164,6 @@ void TriangleMesh::rotate(double angle, Point* center)
this->translate(+center->x, +center->y, 0);
}
void
TriangleMesh::slice(const std::vector<float> &z, std::vector<Polygons>* layers)
{
/*
This method gets called with a list of unscaled Z coordinates and outputs
a vector pointer having the same number of items as the original list.
Each item is a vector of polygons created by slicing our mesh at the
given heights.
This method should basically combine the behavior of the existing
Perl methods defined in lib/Slic3r/TriangleMesh.pm:
- analyze(): this creates the 'facets_edges' and the 'edges_facets'
tables (we don't need the 'edges' table)
- slice_facet(): this has to be done for each facet. It generates
intersection lines with each plane identified by the Z list.
The get_layer_range() binary search used to identify the Z range
of the facet is already ported to C++ (see Object.xsp)
- make_loops(): this has to be done for each layer. It creates polygons
from the lines generated by the previous step.
At the end, we free the tables generated by analyze() as we don't
need them anymore.
FUTURE: parallelize slice_facet() and make_loops()
NOTE: this method accepts a vector of floats because the mesh coordinate
type is float.
*/
// build a table to map a facet_idx to its three edge indices
this->require_shared_vertices();
typedef std::pair<int,int> t_edge;
typedef std::vector<t_edge> t_edges; // edge_idx => a_id,b_id
typedef std::map<t_edge,int> t_edges_map; // a_id,b_id => edge_idx
typedef std::vector< std::vector<int> > t_facets_edges;
t_facets_edges facets_edges;
facets_edges.resize(this->stl.stats.number_of_facets);
{
t_edges edges;
// reserve() instad of resize() because otherwise we couldn't read .size() below to assign edge_idx
edges.reserve(this->stl.stats.number_of_facets * 3); // number of edges = number of facets * 3
t_edges_map edges_map;
for (int facet_idx = 0; facet_idx < this->stl.stats.number_of_facets; facet_idx++) {
facets_edges[facet_idx].resize(3);
for (int i = 0; i <= 2; i++) {
int a_id = this->stl.v_indices[facet_idx].vertex[i];
int b_id = this->stl.v_indices[facet_idx].vertex[(i+1) % 3];
int edge_idx;
t_edges_map::const_iterator my_edge = edges_map.find(std::make_pair(b_id,a_id));
if (my_edge != edges_map.end()) {
edge_idx = my_edge->second;
} else {
/* admesh can assign the same edge ID to more than two facets (which is
still topologically correct), so we have to search for a duplicate of
this edge too in case it was already seen in this orientation */
my_edge = edges_map.find(std::make_pair(a_id,b_id));
if (my_edge != edges_map.end()) {
edge_idx = my_edge->second;
} else {
// edge isn't listed in table, so we insert it
edge_idx = edges.size();
edges.push_back(std::make_pair(a_id,b_id));
edges_map[ edges[edge_idx] ] = edge_idx;
}
}
facets_edges[facet_idx][i] = edge_idx;
#ifdef SLIC3R_DEBUG
printf(" [facet %d, edge %d] a_id = %d, b_id = %d --> edge %d\n", facet_idx, i, a_id, b_id, edge_idx);
#endif
}
}
}
std::vector<IntersectionLines> lines(z.size());
// clone shared vertices coordinates and scale them
stl_vertex* v_scaled_shared = (stl_vertex*)calloc(this->stl.stats.shared_vertices, sizeof(stl_vertex));
std::copy(this->stl.v_shared, this->stl.v_shared + this->stl.stats.shared_vertices, v_scaled_shared);
for (int i = 0; i < this->stl.stats.shared_vertices; i++) {
v_scaled_shared[i].x /= SCALING_FACTOR;
v_scaled_shared[i].y /= SCALING_FACTOR;
v_scaled_shared[i].z /= SCALING_FACTOR;
}
for (int facet_idx = 0; facet_idx < this->stl.stats.number_of_facets; facet_idx++) {
stl_facet* facet = &this->stl.facet_start[facet_idx];
// find facet extents
float min_z = fminf(facet->vertex[0].z, fminf(facet->vertex[1].z, facet->vertex[2].z));
float max_z = fmaxf(facet->vertex[0].z, fmaxf(facet->vertex[1].z, facet->vertex[2].z));
#ifdef SLIC3R_DEBUG
printf("\n==> FACET %d (%f,%f,%f - %f,%f,%f - %f,%f,%f):\n", facet_idx,
facet->vertex[0].x, facet->vertex[0].y, facet->vertex[0].z,
facet->vertex[1].x, facet->vertex[1].y, facet->vertex[1].z,
facet->vertex[2].x, facet->vertex[2].y, facet->vertex[2].z);
printf("z: min = %.2f, max = %.2f\n", min_z, max_z);
#endif
// find layer extents
std::vector<float>::const_iterator min_layer, max_layer;
min_layer = std::lower_bound(z.begin(), z.end(), min_z); // first layer whose slice_z is >= min_z
max_layer = std::upper_bound(z.begin() + (min_layer - z.begin()), z.end(), max_z) - 1; // last layer whose slice_z is <= max_z
#ifdef SLIC3R_DEBUG
printf("layers: min = %d, max = %d\n", (int)(min_layer - z.begin()), (int)(max_layer - z.begin()));
#endif
for (std::vector<float>::const_iterator it = min_layer; it != max_layer + 1; ++it) {
std::vector<float>::size_type layer_idx = it - z.begin();
float slice_z = *it / SCALING_FACTOR;
std::vector<IntersectionPoint> points;
std::vector< std::vector<IntersectionPoint>::size_type > points_on_layer;
bool found_horizontal_edge = false;
/* reorder vertices so that the first one is the one with lowest Z
this is needed to get all intersection lines in a consistent order
(external on the right of the line) */
int i = 0;
if (facet->vertex[1].z == min_z) {
// vertex 1 has lowest Z
i = 1;
} else if (facet->vertex[2].z == min_z) {
// vertex 2 has lowest Z
i = 2;
}
for (int j = i; (j-i) < 3; j++) { // loop through facet edges
int edge_id = facets_edges[facet_idx][j % 3];
int a_id = this->stl.v_indices[facet_idx].vertex[j % 3];
int b_id = this->stl.v_indices[facet_idx].vertex[(j+1) % 3];
stl_vertex* a = &v_scaled_shared[a_id];
stl_vertex* b = &v_scaled_shared[b_id];
if (a->z == b->z && a->z == slice_z) {
// edge is horizontal and belongs to the current layer
/* We assume that this method is never being called for horizontal
facets, so no other edge is going to be on this layer. */
stl_vertex* v0 = &v_scaled_shared[ this->stl.v_indices[facet_idx].vertex[0] ];
stl_vertex* v1 = &v_scaled_shared[ this->stl.v_indices[facet_idx].vertex[1] ];
stl_vertex* v2 = &v_scaled_shared[ this->stl.v_indices[facet_idx].vertex[2] ];
IntersectionLine line;
if (min_z == max_z) {
line.edge_type = feHorizontal;
} else if (v0->z < slice_z || v1->z < slice_z || v2->z < slice_z) {
line.edge_type = feTop;
std::swap(a, b);
std::swap(a_id, b_id);
} else {
line.edge_type = feBottom;
}
line.a.x = a->x;
line.a.y = a->y;
line.b.x = b->x;
line.b.y = b->y;
line.a_id = a_id;
line.b_id = b_id;
lines[layer_idx].push_back(line);
found_horizontal_edge = true;
// if this is a top or bottom edge, we can stop looping through edges
// because we won't find anything interesting
if (line.edge_type != feHorizontal) break;
} else if (a->z == slice_z) {
IntersectionPoint point;
point.x = a->x;
point.y = a->y;
point.point_id = a_id;
points.push_back(point);
points_on_layer.push_back(points.size()-1);
} else if (b->z == slice_z) {
IntersectionPoint point;
point.x = b->x;
point.y = b->y;
point.point_id = b_id;
points.push_back(point);
points_on_layer.push_back(points.size()-1);
} else if ((a->z < slice_z && b->z > slice_z) || (b->z < slice_z && a->z > slice_z)) {
// edge intersects the current layer; calculate intersection
IntersectionPoint point;
point.x = b->x + (a->x - b->x) * (slice_z - b->z) / (a->z - b->z);
point.y = b->y + (a->y - b->y) * (slice_z - b->z) / (a->z - b->z);
point.edge_id = edge_id;
points.push_back(point);
}
}
if (found_horizontal_edge) continue;
if (!points_on_layer.empty()) {
// we can't have only one point on layer because each vertex gets detected
// twice (once for each edge), and we can't have three points on layer because
// we assume this code is not getting called for horizontal facets
assert(points_on_layer.size() == 2);
assert( points[ points_on_layer[0] ].point_id == points[ points_on_layer[1] ].point_id );
if (points.size() < 3) continue; // no intersection point, this is a V-shaped facet tangent to plane
points.erase( points.begin() + points_on_layer[1] );
}
if (!points.empty()) {
assert(points.size() == 2); // facets must intersect each plane 0 or 2 times
IntersectionLine line;
line.a.x = points[1].x;
line.a.y = points[1].y;
line.b.x = points[0].x;
line.b.y = points[0].y;
line.a_id = points[1].point_id;
line.b_id = points[0].point_id;
line.edge_a_id = points[1].edge_id;
line.edge_b_id = points[0].edge_id;
lines[layer_idx].push_back(line);
}
}
}
free(v_scaled_shared);
// build loops
layers->resize(z.size());
for (std::vector<IntersectionLines>::iterator it = lines.begin(); it != lines.end(); ++it) {
int layer_idx = it - lines.begin();
#ifdef SLIC3R_DEBUG
printf("Layer %d:\n", layer_idx);
#endif
/*
SVG svg("lines.svg");
for (IntersectionLines::iterator line = it->begin(); line != it->end(); ++line)
svg.AddLine(*line);
svg.Close();
*/
// remove tangent edges
for (IntersectionLines::iterator line = it->begin(); line != it->end(); ++line) {
if (line->skip || line->edge_type == feNone) continue;
/* if the line is a facet edge, find another facet edge
having the same endpoints but in reverse order */
for (IntersectionLines::iterator line2 = line + 1; line2 != it->end(); ++line2) {
if (line2->skip || line2->edge_type == feNone) continue;
// are these facets adjacent? (sharing a common edge on this layer)
if (line->a_id == line2->a_id && line->b_id == line2->b_id) {
line2->skip = true;
/* if they are both oriented upwards or downwards (like a 'V')
then we can remove both edges from this layer since it won't
affect the sliced shape */
/* if one of them is oriented upwards and the other is oriented
downwards, let's only keep one of them (it doesn't matter which
one since all 'top' lines were reversed at slicing) */
if (line->edge_type == line2->edge_type) {
line->skip = true;
break;
}
} else if (line->a_id == line2->b_id && line->b_id == line2->a_id) {
/* if this edge joins two horizontal facets, remove both of them */
if (line->edge_type == feHorizontal && line2->edge_type == feHorizontal) {
line->skip = true;
line2->skip = true;
break;
}
}
}
}
// build a map of lines by edge_a_id and a_id
std::vector<IntersectionLinePtrs> by_edge_a_id, by_a_id;
by_edge_a_id.resize(this->stl.stats.number_of_facets * 3);
by_a_id.resize(this->stl.stats.shared_vertices);
for (IntersectionLines::iterator line = it->begin(); line != it->end(); ++line) {
if (line->skip) continue;
if (line->edge_a_id != -1) by_edge_a_id[line->edge_a_id].push_back(&(*line));
if (line->a_id != -1) by_a_id[line->a_id].push_back(&(*line));
}
CYCLE: while (1) {
// take first spare line and start a new loop
IntersectionLine* first_line = NULL;
for (IntersectionLines::iterator line = it->begin(); line != it->end(); ++line) {
if (line->skip) continue;
first_line = &(*line);
break;
}
if (first_line == NULL) break;
first_line->skip = true;
IntersectionLinePtrs loop;
loop.push_back(first_line);
/*
printf("first_line edge_a_id = %d, edge_b_id = %d, a_id = %d, b_id = %d, a = %d,%d, b = %d,%d\n",
first_line->edge_a_id, first_line->edge_b_id, first_line->a_id, first_line->b_id,
first_line->a.x, first_line->a.y, first_line->b.x, first_line->b.y);
*/
while (1) {
// find a line starting where last one finishes
IntersectionLine* next_line = NULL;
if (loop.back()->edge_b_id != -1) {
IntersectionLinePtrs* candidates = &(by_edge_a_id[loop.back()->edge_b_id]);
for (IntersectionLinePtrs::iterator lineptr = candidates->begin(); lineptr != candidates->end(); ++lineptr) {
if ((*lineptr)->skip) continue;
next_line = *lineptr;
break;
}
}
if (next_line == NULL && loop.back()->b_id != -1) {
IntersectionLinePtrs* candidates = &(by_a_id[loop.back()->b_id]);
for (IntersectionLinePtrs::iterator lineptr = candidates->begin(); lineptr != candidates->end(); ++lineptr) {
if ((*lineptr)->skip) continue;
next_line = *lineptr;
break;
}
}
if (next_line == NULL) {
// check whether we closed this loop
if ((loop.front()->edge_a_id != -1 && loop.front()->edge_a_id == loop.back()->edge_b_id)
|| (loop.front()->a_id != -1 && loop.front()->a_id == loop.back()->b_id)) {
// loop is complete
Polygon p;
p.points.reserve(loop.size());
for (IntersectionLinePtrs::iterator lineptr = loop.begin(); lineptr != loop.end(); ++lineptr) {
p.points.push_back((*lineptr)->a);
}
(*layers)[layer_idx].push_back(p);
#ifdef SLIC3R_DEBUG
printf(" Discovered %s polygon of %d points\n", (p.is_counter_clockwise() ? "ccw" : "cw"), (int)p.points.size());
#endif
goto CYCLE;
}
// we can't close this loop!
//// push @failed_loops, [@loop];
#ifdef SLIC3R_DEBUG
printf(" Unable to close this loop having %d points\n", (int)loop.size());
#endif
goto CYCLE;
}
/*
printf("next_line edge_a_id = %d, edge_b_id = %d, a_id = %d, b_id = %d, a = %d,%d, b = %d,%d\n",
next_line->edge_a_id, next_line->edge_b_id, next_line->a_id, next_line->b_id,
next_line->a.x, next_line->a.y, next_line->b.x, next_line->b.y);
*/
loop.push_back(next_line);
next_line->skip = true;
}
}
}
}
class _area_comp {
public:
_area_comp(std::vector<double>* _aa) : abs_area(_aa) {};
bool operator() (const size_t &a, const size_t &b) {
return (*this->abs_area)[a] > (*this->abs_area)[b];
}
private:
std::vector<double>* abs_area;
};
void
TriangleMesh::slice(const std::vector<float> &z, std::vector<ExPolygons>* layers)
{
std::vector<Polygons> layers_p;
this->slice(z, &layers_p);
/*
Input loops are not suitable for evenodd nor nonzero fill types, as we might get
two consecutive concentric loops having the same winding order - and we have to
respect such order. In that case, evenodd would create wrong inversions, and nonzero
would ignore holes inside two concentric contours.
So we're ordering loops and collapse consecutive concentric loops having the same
winding order.
TODO: find a faster algorithm for this, maybe with some sort of binary search.
If we computed a "nesting tree" we could also just remove the consecutive loops
having the same winding order, and remove the extra one(s) so that we could just
supply everything to offset_ex() instead of performing several union/diff calls.
we sort by area assuming that the outermost loops have larger area;
the previous sorting method, based on $b->contains_point($a->[0]), failed to nest
loops correctly in some edge cases when original model had overlapping facets
*/
layers->resize(z.size());
for (std::vector<Polygons>::const_iterator loops = layers_p.begin(); loops != layers_p.end(); ++loops) {
size_t layer_id = loops - layers_p.begin();
std::vector<double> area;
std::vector<double> abs_area;
std::vector<size_t> sorted_area; // vector of indices
for (Polygons::const_iterator loop = loops->begin(); loop != loops->end(); ++loop) {
double a = loop->area();
area.push_back(a);
abs_area.push_back(std::fabs(a));
sorted_area.push_back(loop - loops->begin());
}
std::sort(sorted_area.begin(), sorted_area.end(), _area_comp(&abs_area)); // outer first
// we don't perform a safety offset now because it might reverse cw loops
Polygons slices;
for (std::vector<size_t>::const_iterator loop_idx = sorted_area.begin(); loop_idx != sorted_area.end(); ++loop_idx) {
/* we rely on the already computed area to determine the winding order
of the loops, since the Orientation() function provided by Clipper
would do the same, thus repeating the calculation */
Polygons::const_iterator loop = loops->begin() + *loop_idx;
if (area[*loop_idx] >= 0) {
slices.push_back(*loop);
} else {
diff(slices, *loop, slices);
}
}
// perform a safety offset to merge very close facets (TODO: find test case for this)
double safety_offset = scale_(0.0499);
ExPolygons ex_slices;
offset2_ex(slices, ex_slices, +safety_offset, -safety_offset);
#ifdef SLIC3R_DEBUG
size_t holes_count = 0;
for (ExPolygons::const_iterator e = ex_slices.begin(); e != ex_slices.end(); ++e) {
holes_count += e->holes.size();
}
printf("Layer %zu (slice_z = %.2f): %zu surface(s) having %zu holes detected from %zu polylines\n",
layer_id, z[layer_id], ex_slices.size(), holes_count, loops->size());
#endif
ExPolygons* layer = &(*layers)[layer_id];
layer->insert(layer->end(), ex_slices.begin(), ex_slices.end());
}
}
TriangleMeshPtrs
TriangleMesh::split() const
{
@ -777,4 +333,462 @@ void TriangleMesh::ReadFromPerl(SV* vertices, SV* facets)
}
#endif
void
TriangleMeshSlicer::slice(const std::vector<float> &z, std::vector<Polygons>* layers)
{
/*
This method gets called with a list of unscaled Z coordinates and outputs
a vector pointer having the same number of items as the original list.
Each item is a vector of polygons created by slicing our mesh at the
given heights.
This method should basically combine the behavior of the existing
Perl methods defined in lib/Slic3r/TriangleMesh.pm:
- analyze(): this creates the 'facets_edges' and the 'edges_facets'
tables (we don't need the 'edges' table)
- slice_facet(): this has to be done for each facet. It generates
intersection lines with each plane identified by the Z list.
The get_layer_range() binary search used to identify the Z range
of the facet is already ported to C++ (see Object.xsp)
- make_loops(): this has to be done for each layer. It creates polygons
from the lines generated by the previous step.
At the end, we free the tables generated by analyze() as we don't
need them anymore.
FUTURE: parallelize slice_facet() and make_loops()
NOTE: this method accepts a vector of floats because the mesh coordinate
type is float.
*/
std::vector<IntersectionLines> lines(z.size());
for (int facet_idx = 0; facet_idx < this->mesh->stl.stats.number_of_facets; facet_idx++) {
stl_facet* facet = &this->mesh->stl.facet_start[facet_idx];
// find facet extents
float min_z = fminf(facet->vertex[0].z, fminf(facet->vertex[1].z, facet->vertex[2].z));
float max_z = fmaxf(facet->vertex[0].z, fmaxf(facet->vertex[1].z, facet->vertex[2].z));
#ifdef SLIC3R_DEBUG
printf("\n==> FACET %d (%f,%f,%f - %f,%f,%f - %f,%f,%f):\n", facet_idx,
facet->vertex[0].x, facet->vertex[0].y, facet->vertex[0].z,
facet->vertex[1].x, facet->vertex[1].y, facet->vertex[1].z,
facet->vertex[2].x, facet->vertex[2].y, facet->vertex[2].z);
printf("z: min = %.2f, max = %.2f\n", min_z, max_z);
#endif
// find layer extents
std::vector<float>::const_iterator min_layer, max_layer;
min_layer = std::lower_bound(z.begin(), z.end(), min_z); // first layer whose slice_z is >= min_z
max_layer = std::upper_bound(z.begin() + (min_layer - z.begin()), z.end(), max_z) - 1; // last layer whose slice_z is <= max_z
#ifdef SLIC3R_DEBUG
printf("layers: min = %d, max = %d\n", (int)(min_layer - z.begin()), (int)(max_layer - z.begin()));
#endif
for (std::vector<float>::const_iterator it = min_layer; it != max_layer + 1; ++it) {
std::vector<float>::size_type layer_idx = it - z.begin();
this->slice_facet(*it / SCALING_FACTOR, *facet, facet_idx, min_z, max_z, &lines[layer_idx]);
}
}
// v_scaled_shared could be freed here
// build loops
layers->resize(z.size());
for (std::vector<IntersectionLines>::iterator it = lines.begin(); it != lines.end(); ++it) {
int layer_idx = it - lines.begin();
#ifdef SLIC3R_DEBUG
printf("Layer %d:\n", layer_idx);
#endif
/*
SVG svg("lines.svg");
for (IntersectionLines::iterator line = it->begin(); line != it->end(); ++line) {
svg.AddLine(*line);
}
svg.Close();
*/
// remove tangent edges
for (IntersectionLines::iterator line = it->begin(); line != it->end(); ++line) {
if (line->skip || line->edge_type == feNone) continue;
/* if the line is a facet edge, find another facet edge
having the same endpoints but in reverse order */
for (IntersectionLines::iterator line2 = line + 1; line2 != it->end(); ++line2) {
if (line2->skip || line2->edge_type == feNone) continue;
// are these facets adjacent? (sharing a common edge on this layer)
if (line->a_id == line2->a_id && line->b_id == line2->b_id) {
line2->skip = true;
/* if they are both oriented upwards or downwards (like a 'V')
then we can remove both edges from this layer since it won't
affect the sliced shape */
/* if one of them is oriented upwards and the other is oriented
downwards, let's only keep one of them (it doesn't matter which
one since all 'top' lines were reversed at slicing) */
if (line->edge_type == line2->edge_type) {
line->skip = true;
break;
}
} else if (line->a_id == line2->b_id && line->b_id == line2->a_id) {
/* if this edge joins two horizontal facets, remove both of them */
if (line->edge_type == feHorizontal && line2->edge_type == feHorizontal) {
line->skip = true;
line2->skip = true;
break;
}
}
}
}
// build a map of lines by edge_a_id and a_id
std::vector<IntersectionLinePtrs> by_edge_a_id, by_a_id;
by_edge_a_id.resize(this->mesh->stl.stats.number_of_facets * 3);
by_a_id.resize(this->mesh->stl.stats.shared_vertices);
for (IntersectionLines::iterator line = it->begin(); line != it->end(); ++line) {
if (line->skip) continue;
if (line->edge_a_id != -1) by_edge_a_id[line->edge_a_id].push_back(&(*line));
if (line->a_id != -1) by_a_id[line->a_id].push_back(&(*line));
}
CYCLE: while (1) {
// take first spare line and start a new loop
IntersectionLine* first_line = NULL;
for (IntersectionLines::iterator line = it->begin(); line != it->end(); ++line) {
if (line->skip) continue;
first_line = &(*line);
break;
}
if (first_line == NULL) break;
first_line->skip = true;
IntersectionLinePtrs loop;
loop.push_back(first_line);
/*
printf("first_line edge_a_id = %d, edge_b_id = %d, a_id = %d, b_id = %d, a = %d,%d, b = %d,%d\n",
first_line->edge_a_id, first_line->edge_b_id, first_line->a_id, first_line->b_id,
first_line->a.x, first_line->a.y, first_line->b.x, first_line->b.y);
*/
while (1) {
// find a line starting where last one finishes
IntersectionLine* next_line = NULL;
if (loop.back()->edge_b_id != -1) {
IntersectionLinePtrs* candidates = &(by_edge_a_id[loop.back()->edge_b_id]);
for (IntersectionLinePtrs::iterator lineptr = candidates->begin(); lineptr != candidates->end(); ++lineptr) {
if ((*lineptr)->skip) continue;
next_line = *lineptr;
break;
}
}
if (next_line == NULL && loop.back()->b_id != -1) {
IntersectionLinePtrs* candidates = &(by_a_id[loop.back()->b_id]);
for (IntersectionLinePtrs::iterator lineptr = candidates->begin(); lineptr != candidates->end(); ++lineptr) {
if ((*lineptr)->skip) continue;
next_line = *lineptr;
break;
}
}
if (next_line == NULL) {
// check whether we closed this loop
if ((loop.front()->edge_a_id != -1 && loop.front()->edge_a_id == loop.back()->edge_b_id)
|| (loop.front()->a_id != -1 && loop.front()->a_id == loop.back()->b_id)) {
// loop is complete
Polygon p;
p.points.reserve(loop.size());
for (IntersectionLinePtrs::iterator lineptr = loop.begin(); lineptr != loop.end(); ++lineptr) {
p.points.push_back((*lineptr)->a);
}
(*layers)[layer_idx].push_back(p);
#ifdef SLIC3R_DEBUG
printf(" Discovered %s polygon of %d points\n", (p.is_counter_clockwise() ? "ccw" : "cw"), (int)p.points.size());
#endif
goto CYCLE;
}
// we can't close this loop!
//// push @failed_loops, [@loop];
//#ifdef SLIC3R_DEBUG
printf(" Unable to close this loop having %d points\n", (int)loop.size());
//#endif
goto CYCLE;
}
/*
printf("next_line edge_a_id = %d, edge_b_id = %d, a_id = %d, b_id = %d, a = %d,%d, b = %d,%d\n",
next_line->edge_a_id, next_line->edge_b_id, next_line->a_id, next_line->b_id,
next_line->a.x, next_line->a.y, next_line->b.x, next_line->b.y);
*/
loop.push_back(next_line);
next_line->skip = true;
}
}
}
}
void
TriangleMeshSlicer::slice_facet(float slice_z, const stl_facet &facet, const int &facet_idx, const float &min_z, const float &max_z, std::vector<IntersectionLine>* lines) const
{
std::vector<IntersectionPoint> points;
std::vector< std::vector<IntersectionPoint>::size_type > points_on_layer;
bool found_horizontal_edge = false;
/* reorder vertices so that the first one is the one with lowest Z
this is needed to get all intersection lines in a consistent order
(external on the right of the line) */
int i = 0;
if (facet.vertex[1].z == min_z) {
// vertex 1 has lowest Z
i = 1;
} else if (facet.vertex[2].z == min_z) {
// vertex 2 has lowest Z
i = 2;
}
for (int j = i; (j-i) < 3; j++) { // loop through facet edges
int edge_id = this->facets_edges[facet_idx][j % 3];
int a_id = this->mesh->stl.v_indices[facet_idx].vertex[j % 3];
int b_id = this->mesh->stl.v_indices[facet_idx].vertex[(j+1) % 3];
stl_vertex* a = &this->v_scaled_shared[a_id];
stl_vertex* b = &this->v_scaled_shared[b_id];
if (a->z == b->z && a->z == slice_z) {
// edge is horizontal and belongs to the current layer
/* We assume that this method is never being called for horizontal
facets, so no other edge is going to be on this layer. */
stl_vertex* v0 = &this->v_scaled_shared[ this->mesh->stl.v_indices[facet_idx].vertex[0] ];
stl_vertex* v1 = &this->v_scaled_shared[ this->mesh->stl.v_indices[facet_idx].vertex[1] ];
stl_vertex* v2 = &this->v_scaled_shared[ this->mesh->stl.v_indices[facet_idx].vertex[2] ];
IntersectionLine line;
if (min_z == max_z) {
line.edge_type = feHorizontal;
} else if (v0->z < slice_z || v1->z < slice_z || v2->z < slice_z) {
line.edge_type = feTop;
std::swap(a, b);
std::swap(a_id, b_id);
} else {
line.edge_type = feBottom;
}
line.a.x = a->x;
line.a.y = a->y;
line.b.x = b->x;
line.b.y = b->y;
line.a_id = a_id;
line.b_id = b_id;
lines->push_back(line);
found_horizontal_edge = true;
// if this is a top or bottom edge, we can stop looping through edges
// because we won't find anything interesting
if (line.edge_type != feHorizontal) return;
} else if (a->z == slice_z) {
IntersectionPoint point;
point.x = a->x;
point.y = a->y;
point.point_id = a_id;
points.push_back(point);
points_on_layer.push_back(points.size()-1);
} else if (b->z == slice_z) {
IntersectionPoint point;
point.x = b->x;
point.y = b->y;
point.point_id = b_id;
points.push_back(point);
points_on_layer.push_back(points.size()-1);
} else if ((a->z < slice_z && b->z > slice_z) || (b->z < slice_z && a->z > slice_z)) {
// edge intersects the current layer; calculate intersection
IntersectionPoint point;
point.x = b->x + (a->x - b->x) * (slice_z - b->z) / (a->z - b->z);
point.y = b->y + (a->y - b->y) * (slice_z - b->z) / (a->z - b->z);
point.edge_id = edge_id;
points.push_back(point);
}
}
if (found_horizontal_edge) return;
if (!points_on_layer.empty()) {
// we can't have only one point on layer because each vertex gets detected
// twice (once for each edge), and we can't have three points on layer because
// we assume this code is not getting called for horizontal facets
assert(points_on_layer.size() == 2);
assert( points[ points_on_layer[0] ].point_id == points[ points_on_layer[1] ].point_id );
if (points.size() < 3) return; // no intersection point, this is a V-shaped facet tangent to plane
points.erase( points.begin() + points_on_layer[1] );
}
if (!points.empty()) {
assert(points.size() == 2); // facets must intersect each plane 0 or 2 times
IntersectionLine line;
line.a.x = points[1].x;
line.a.y = points[1].y;
line.b.x = points[0].x;
line.b.y = points[0].y;
line.a_id = points[1].point_id;
line.b_id = points[0].point_id;
line.edge_a_id = points[1].edge_id;
line.edge_b_id = points[0].edge_id;
lines->push_back(line);
return;
}
}
class _area_comp {
public:
_area_comp(std::vector<double>* _aa) : abs_area(_aa) {};
bool operator() (const size_t &a, const size_t &b) {
return (*this->abs_area)[a] > (*this->abs_area)[b];
}
private:
std::vector<double>* abs_area;
};
void
TriangleMeshSlicer::slice(const std::vector<float> &z, std::vector<ExPolygons>* layers)
{
std::vector<Polygons> layers_p;
this->slice(z, &layers_p);
/*
Input loops are not suitable for evenodd nor nonzero fill types, as we might get
two consecutive concentric loops having the same winding order - and we have to
respect such order. In that case, evenodd would create wrong inversions, and nonzero
would ignore holes inside two concentric contours.
So we're ordering loops and collapse consecutive concentric loops having the same
winding order.
TODO: find a faster algorithm for this, maybe with some sort of binary search.
If we computed a "nesting tree" we could also just remove the consecutive loops
having the same winding order, and remove the extra one(s) so that we could just
supply everything to offset_ex() instead of performing several union/diff calls.
we sort by area assuming that the outermost loops have larger area;
the previous sorting method, based on $b->contains_point($a->[0]), failed to nest
loops correctly in some edge cases when original model had overlapping facets
*/
layers->resize(z.size());
for (std::vector<Polygons>::const_iterator loops = layers_p.begin(); loops != layers_p.end(); ++loops) {
size_t layer_id = loops - layers_p.begin();
std::vector<double> area;
std::vector<double> abs_area;
std::vector<size_t> sorted_area; // vector of indices
for (Polygons::const_iterator loop = loops->begin(); loop != loops->end(); ++loop) {
double a = loop->area();
area.push_back(a);
abs_area.push_back(std::fabs(a));
sorted_area.push_back(loop - loops->begin());
}
std::sort(sorted_area.begin(), sorted_area.end(), _area_comp(&abs_area)); // outer first
// we don't perform a safety offset now because it might reverse cw loops
Polygons slices;
for (std::vector<size_t>::const_iterator loop_idx = sorted_area.begin(); loop_idx != sorted_area.end(); ++loop_idx) {
/* we rely on the already computed area to determine the winding order
of the loops, since the Orientation() function provided by Clipper
would do the same, thus repeating the calculation */
Polygons::const_iterator loop = loops->begin() + *loop_idx;
if (area[*loop_idx] >= 0) {
slices.push_back(*loop);
} else {
diff(slices, *loop, slices);
}
}
// perform a safety offset to merge very close facets (TODO: find test case for this)
double safety_offset = scale_(0.0499);
ExPolygons ex_slices;
offset2_ex(slices, ex_slices, +safety_offset, -safety_offset);
#ifdef SLIC3R_DEBUG
size_t holes_count = 0;
for (ExPolygons::const_iterator e = ex_slices.begin(); e != ex_slices.end(); ++e) {
holes_count += e->holes.size();
}
printf("Layer %zu (slice_z = %.2f): %zu surface(s) having %zu holes detected from %zu polylines\n",
layer_id, z[layer_id], ex_slices.size(), holes_count, loops->size());
#endif
ExPolygons* layer = &(*layers)[layer_id];
layer->insert(layer->end(), ex_slices.begin(), ex_slices.end());
}
}
TriangleMeshSlicer::TriangleMeshSlicer(TriangleMesh* _mesh) : mesh(_mesh), v_scaled_shared(NULL)
{
// build a table to map a facet_idx to its three edge indices
this->mesh->require_shared_vertices();
typedef std::pair<int,int> t_edge;
typedef std::vector<t_edge> t_edges; // edge_idx => a_id,b_id
typedef std::map<t_edge,int> t_edges_map; // a_id,b_id => edge_idx
this->facets_edges.resize(this->mesh->stl.stats.number_of_facets);
{
t_edges edges;
// reserve() instad of resize() because otherwise we couldn't read .size() below to assign edge_idx
edges.reserve(this->mesh->stl.stats.number_of_facets * 3); // number of edges = number of facets * 3
t_edges_map edges_map;
for (int facet_idx = 0; facet_idx < this->mesh->stl.stats.number_of_facets; facet_idx++) {
this->facets_edges[facet_idx].resize(3);
for (int i = 0; i <= 2; i++) {
int a_id = this->mesh->stl.v_indices[facet_idx].vertex[i];
int b_id = this->mesh->stl.v_indices[facet_idx].vertex[(i+1) % 3];
int edge_idx;
t_edges_map::const_iterator my_edge = edges_map.find(std::make_pair(b_id,a_id));
if (my_edge != edges_map.end()) {
edge_idx = my_edge->second;
} else {
/* admesh can assign the same edge ID to more than two facets (which is
still topologically correct), so we have to search for a duplicate of
this edge too in case it was already seen in this orientation */
my_edge = edges_map.find(std::make_pair(a_id,b_id));
if (my_edge != edges_map.end()) {
edge_idx = my_edge->second;
} else {
// edge isn't listed in table, so we insert it
edge_idx = edges.size();
edges.push_back(std::make_pair(a_id,b_id));
edges_map[ edges[edge_idx] ] = edge_idx;
}
}
this->facets_edges[facet_idx][i] = edge_idx;
#ifdef SLIC3R_DEBUG
printf(" [facet %d, edge %d] a_id = %d, b_id = %d --> edge %d\n", facet_idx, i, a_id, b_id, edge_idx);
#endif
}
}
}
// clone shared vertices coordinates and scale them
this->v_scaled_shared = (stl_vertex*)calloc(this->mesh->stl.stats.shared_vertices, sizeof(stl_vertex));
std::copy(this->mesh->stl.v_shared, this->mesh->stl.v_shared + this->mesh->stl.stats.shared_vertices, this->v_scaled_shared);
for (int i = 0; i < this->mesh->stl.stats.shared_vertices; i++) {
this->v_scaled_shared[i].x /= SCALING_FACTOR;
this->v_scaled_shared[i].y /= SCALING_FACTOR;
this->v_scaled_shared[i].z /= SCALING_FACTOR;
}
}
TriangleMeshSlicer::~TriangleMeshSlicer()
{
if (this->v_scaled_shared != NULL) free(this->v_scaled_shared);
}
}

20
xs/src/TriangleMesh.hpp
View file

@ -12,6 +12,7 @@
namespace Slic3r {
class TriangleMesh;
class TriangleMeshSlicer;
typedef std::vector<TriangleMesh*> TriangleMeshPtrs;
class TriangleMesh
@ -30,8 +31,6 @@ class TriangleMesh
void translate(float x, float y, float z);
void align_to_origin();
void rotate(double angle, Point* center);
void slice(const std::vector<float> &z, std::vector<Polygons>* layers);
void slice(const std::vector<float> &z, std::vector<ExPolygons>* layers);
TriangleMeshPtrs split() const;
void merge(const TriangleMesh* mesh);
void horizontal_projection(ExPolygons &retval) const;
@ -47,6 +46,7 @@ class TriangleMesh
private:
void require_shared_vertices();
friend class TriangleMeshSlicer;
};
enum FacetEdgeType { feNone, feTop, feBottom, feHorizontal };
@ -75,6 +75,22 @@ class IntersectionLine
typedef std::vector<IntersectionLine> IntersectionLines;
typedef std::vector<IntersectionLine*> IntersectionLinePtrs;
class TriangleMeshSlicer
{
public:
TriangleMesh* mesh;
TriangleMeshSlicer(TriangleMesh* _mesh);
~TriangleMeshSlicer();
void slice(const std::vector<float> &z, std::vector<Polygons>* layers);
void slice(const std::vector<float> &z, std::vector<ExPolygons>* layers);
void slice_facet(float slice_z, const stl_facet &facet, const int &facet_idx, const float &min_z, const float &max_z, std::vector<IntersectionLine>* lines) const;
private:
typedef std::vector< std::vector<int> > t_facets_edges;
t_facets_edges facets_edges;
stl_vertex* v_scaled_shared;
};
}
#endif

2
xs/t/01_trianglemesh.t
View file

@ -83,7 +83,7 @@ my $cube = {
my $result = $m->slice(\@z);
my $SCALING_FACTOR = 0.000001;
for my $i (0..$#z) {
is scalar(@{$result->[$i]}), 1, 'number of returned polygons per layer';
is scalar(@{$result->[$i]}), 1, "number of returned polygons per layer (z = " . $z[$i] . ")";
is $result->[$i][0]->area, 20*20/($SCALING_FACTOR**2), 'size of returned polygon';
}
}

3
xs/xsp/TriangleMesh.xsp
View file

@ -142,7 +142,8 @@ TriangleMesh::slice(z)
delete z;
std::vector<ExPolygons> layers;
THIS->slice(z_f, &layers);
TriangleMeshSlicer mslicer(THIS);
mslicer.slice(z_f, &layers);
AV* layers_av = newAV();
av_extend(layers_av, layers.size()-1);