Changed precision from ~25b (2 147.48mm -> 0.1um) to ~40b (17km -> 0.1um)

it's scaled by 10^6 (~20b), with 10^2 (~7b) under the epsilon
also, CLIPPER_OFFSET_SCALE scale the number by 2^17, and clipper disallow the use of the 2 higher bits, so it can't go higher than 45b
This commit is contained in:
supermerill 2019-04-08 19:47:26 +02:00
parent 782d3706e1
commit 18a57ee137
5 changed files with 48 additions and 37 deletions

View File

@ -1258,9 +1258,9 @@ Polygons EdgeGrid::Grid::contours_simplified(coord_t offset, bool fill_holes) co
// 1) Collect the lines. // 1) Collect the lines.
std::vector<Line> lines; std::vector<Line> lines;
EndPointMapType start_point_to_line_idx; EndPointMapType start_point_to_line_idx;
for (int r = 0; r <= int(m_rows); ++ r) { for (coord_t r = 0; r <= coord_t(m_rows); ++ r) {
for (int c = 0; c <= int(m_cols); ++ c) { for (coord_t c = 0; c <= coord_t(m_cols); ++ c) {
int addr = (r + 1) * cell_cols + c + 1; size_t addr = (r + 1) * cell_cols + c + 1;
bool left = cell_inside[addr - 1]; bool left = cell_inside[addr - 1];
bool top = cell_inside[addr - cell_cols]; bool top = cell_inside[addr - cell_cols];
bool current = cell_inside[addr]; bool current = cell_inside[addr];
@ -1364,9 +1364,9 @@ inline int segments_could_intersect(
const Slic3r::Point &ip1, const Slic3r::Point &ip2, const Slic3r::Point &ip1, const Slic3r::Point &ip2,
const Slic3r::Point &jp1, const Slic3r::Point &jp2) const Slic3r::Point &jp1, const Slic3r::Point &jp2)
{ {
Vec2i64 iv = (ip2 - ip1).cast<int64_t>(); Vec2crd iv = (ip2 - ip1);
Vec2i64 vij1 = (jp1 - ip1).cast<int64_t>(); Vec2crd vij1 = (jp1 - ip1);
Vec2i64 vij2 = (jp2 - ip1).cast<int64_t>(); Vec2crd vij2 = (jp2 - ip1);
int64_t tij1 = cross2(iv, vij1); int64_t tij1 = cross2(iv, vij1);
int64_t tij2 = cross2(iv, vij2); int64_t tij2 = cross2(iv, vij2);
int sij1 = (tij1 > 0) ? 1 : ((tij1 < 0) ? -1 : 0); // signum int sij1 = (tij1 > 0) ? 1 : ((tij1 < 0) ? -1 : 0); // signum

View File

@ -276,13 +276,13 @@ int SegmentIntersection::ordering_along_line(const SegmentIntersection &other) c
// other.iSegment succeeds this->iSegment // other.iSegment succeeds this->iSegment
assert(seg_end_a == seg_start_b); assert(seg_end_a == seg_start_b);
// Avoid calling the 128bit x 128bit multiplication below if this->line intersects the common point. // Avoid calling the 128bit x 128bit multiplication below if this->line intersects the common point.
if (cross2(Vec2i64(this->line->dir.cast<int64_t>()), (seg_end_b - this->line->pos).cast<int64_t>()) == 0) if (cross2(Vec2crd(this->line->dir), (seg_end_b - this->line->pos)) == 0)
return 0; return 0;
} else if ((other.iSegment + 1) % poly_a.points.size() == this->iSegment) { } else if ((other.iSegment + 1) % poly_a.points.size() == this->iSegment) {
// this->iSegment succeeds other.iSegment // this->iSegment succeeds other.iSegment
assert(seg_start_a == seg_end_b); assert(seg_start_a == seg_end_b);
// Avoid calling the 128bit x 128bit multiplication below if this->line intersects the common point. // Avoid calling the 128bit x 128bit multiplication below if this->line intersects the common point.
if (cross2(Vec2i64(this->line->dir.cast<int64_t>()), (seg_start_a - this->line->pos).cast<int64_t>()) == 0) if (cross2(Vec2crd(this->line->dir), (seg_start_a - this->line->pos)) == 0)
return 0; return 0;
} else { } else {
// General case. // General case.
@ -290,33 +290,33 @@ int SegmentIntersection::ordering_along_line(const SegmentIntersection &other) c
} }
// First test, whether both points of one segment are completely in one half-plane of the other line. // First test, whether both points of one segment are completely in one half-plane of the other line.
const Vec2i64 vec_b = (seg_end_b - seg_start_b).cast<int64_t>(); const Vec2crd vec_b = (seg_end_b - seg_start_b);
int side_start = signum(cross2(vec_b, (seg_start_a - seg_start_b).cast<int64_t>())); int side_start = signum(cross2(vec_b, (seg_start_a - seg_start_b)));
int side_end = signum(cross2(vec_b, (seg_end_a - seg_start_b).cast<int64_t>())); int side_end = signum(cross2(vec_b, (seg_end_a - seg_start_b)));
int side = side_start * side_end; int side = side_start * side_end;
if (side > 0) if (side > 0)
// This segment is completely inside one half-plane of the other line, therefore the ordering is trivial. // This segment is completely inside one half-plane of the other line, therefore the ordering is trivial.
return signum(cross2(vec_b, this->line->dir.cast<int64_t>())) * side_start; return signum(cross2(vec_b, this->line->dir)) * side_start;
const Vec2i64 vec_a = (seg_end_a - seg_start_a).cast<int64_t>(); const Vec2crd vec_a = (seg_end_a - seg_start_a);
int side_start2 = signum(cross2(vec_a, (seg_start_b - seg_start_a).cast<int64_t>())); int side_start2 = signum(cross2(vec_a, (seg_start_b - seg_start_a)));
int side_end2 = signum(cross2(vec_a, (seg_end_b - seg_start_a).cast<int64_t>())); int side_end2 = signum(cross2(vec_a, (seg_end_b - seg_start_a)));
int side2 = side_start2 * side_end2; int side2 = side_start2 * side_end2;
//if (side == 0 && side2 == 0) //if (side == 0 && side2 == 0)
// The segments share one of their end points. // The segments share one of their end points.
if (side2 > 0) if (side2 > 0)
// This segment is completely inside one half-plane of the other line, therefore the ordering is trivial. // This segment is completely inside one half-plane of the other line, therefore the ordering is trivial.
return signum(cross2(this->line->dir.cast<int64_t>(), vec_a)) * side_start2; return signum(cross2(this->line->dir, vec_a)) * side_start2;
// The two segments intersect and they are not sucessive segments of the same contour. // The two segments intersect and they are not sucessive segments of the same contour.
// Ordering of the points depends on the position of the segment intersection (left / right from this->line), // Ordering of the points depends on the position of the segment intersection (left / right from this->line),
// therefore a simple test over the input segment end points is not sufficient. // therefore a simple test over the input segment end points is not sufficient.
// Find the parameters of intersection of the two segmetns with this->line. // Find the parameters of intersection of the two segmetns with this->line.
int64_t denom1 = cross2(this->line->dir.cast<int64_t>(), vec_a); int64_t denom1 = cross2(this->line->dir, vec_a);
int64_t denom2 = cross2(this->line->dir.cast<int64_t>(), vec_b); int64_t denom2 = cross2(this->line->dir, vec_b);
Vec2i64 vx_a = (seg_start_a - this->line->pos).cast<int64_t>(); Vec2crd vx_a = (seg_start_a - this->line->pos);
Vec2i64 vx_b = (seg_start_b - this->line->pos).cast<int64_t>(); Vec2crd vx_b = (seg_start_b - this->line->pos);
int64_t t1_times_denom1 = vx_a(0) * vec_a(1) - vx_a(1) * vec_a(0); int64_t t1_times_denom1 = vx_a(0) * vec_a(1) - vx_a(1) * vec_a(0);
int64_t t2_times_denom2 = vx_b(0) * vec_b(1) - vx_b(1) * vec_b(0); int64_t t2_times_denom2 = vx_b(0) * vec_b(1) - vx_b(1) * vec_b(0);
assert(denom1 != 0); assert(denom1 != 0);
@ -330,7 +330,7 @@ bool SegmentIntersection::operator<(const SegmentIntersection &other) const
#ifdef _DEBUG #ifdef _DEBUG
Point p1 = this->pos(); Point p1 = this->pos();
Point p2 = other.pos(); Point p2 = other.pos();
int64_t d = this->line->dir.cast<int64_t>().dot((p2 - p1).cast<int64_t>()); int64_t d = this->line->dir.dot((p2 - p1));
#endif /* _DEBUG */ #endif /* _DEBUG */
int ordering = this->ordering_along_line(other); int ordering = this->ordering_along_line(other);
#ifdef _DEBUG #ifdef _DEBUG
@ -510,7 +510,7 @@ static bool prepare_infill_hatching_segments(
for (size_t i = 1; i < sil.intersections.size(); ++ i) { for (size_t i = 1; i < sil.intersections.size(); ++ i) {
Point p1 = sil.intersections[i - 1].pos(); Point p1 = sil.intersections[i - 1].pos();
Point p2 = sil.intersections[i].pos(); Point p2 = sil.intersections[i].pos();
int64_t d = sil.dir.cast<int64_t>().dot((p2 - p1).cast<int64_t>()); int64_t d = sil.dir.dot((p2 - p1));
assert(d >= - int64_t(SCALED_EPSILON)); assert(d >= - int64_t(SCALED_EPSILON));
} }
#endif /* _DEBUG */ #endif /* _DEBUG */

View File

@ -1913,11 +1913,12 @@ static Points::iterator project_point_to_polygon_and_insert(Polygon &polygon, co
const Point &p2 = polygon.points[j]; const Point &p2 = polygon.points[j];
const Slic3r::Point v_seg = p2 - p1; const Slic3r::Point v_seg = p2 - p1;
const Slic3r::Point v_pt = pt - p1; const Slic3r::Point v_pt = pt - p1;
const int64_t l2_seg = int64_t(v_seg(0)) * int64_t(v_seg(0)) + int64_t(v_seg(1)) * int64_t(v_seg(1)); //going into double because coord*coord can overflow.
int64_t t_pt = int64_t(v_seg(0)) * int64_t(v_pt(0)) + int64_t(v_seg(1)) * int64_t(v_pt(1)); const double l2_seg = double(v_seg(0)) * double(v_seg(0)) + double(v_seg(1)) * double(v_seg(1));
double t_pt = double(v_seg(0)) * double(v_pt(0)) + double(v_seg(1)) * double(v_pt(1));
if (t_pt < 0) { if (t_pt < 0) {
// Closest to p1. // Closest to p1.
double dabs = sqrt(int64_t(v_pt(0)) * int64_t(v_pt(0)) + int64_t(v_pt(1)) * int64_t(v_pt(1))); double dabs = sqrt(double(v_pt(0)) * double(v_pt(0)) + double(v_pt(1)) * double(v_pt(1)));
if (dabs < d_min) { if (dabs < d_min) {
d_min = dabs; d_min = dabs;
i_min = i; i_min = i;
@ -1930,7 +1931,7 @@ static Points::iterator project_point_to_polygon_and_insert(Polygon &polygon, co
} else { } else {
// Closest to the segment. // Closest to the segment.
assert(t_pt >= 0 && t_pt <= l2_seg); assert(t_pt >= 0 && t_pt <= l2_seg);
int64_t d_seg = int64_t(v_seg(1)) * int64_t(v_pt(0)) - int64_t(v_seg(0)) * int64_t(v_pt(1)); double d_seg = double(v_seg(1)) * double(v_pt(0)) - double(v_seg(0)) * double(v_pt(1));
double d = double(d_seg) / sqrt(double(l2_seg)); double d = double(d_seg) / sqrt(double(l2_seg));
double dabs = std::abs(d); double dabs = std::abs(d);
if (dabs < d_min) { if (dabs < d_min) {
@ -2010,9 +2011,9 @@ std::vector<float> polygon_angles_at_vertices(const Polygon &polygon, const std:
const Point &p2 = polygon.points[idx_next]; const Point &p2 = polygon.points[idx_next];
const Point v1 = p1 - p0; const Point v1 = p1 - p0;
const Point v2 = p2 - p1; const Point v2 = p2 - p1;
int64_t dot = int64_t(v1(0))*int64_t(v2(0)) + int64_t(v1(1))*int64_t(v2(1)); double dot = double(v1(0))*double(v2(0)) + double(v1(1))*double(v2(1));
int64_t cross = int64_t(v1(0))*int64_t(v2(1)) - int64_t(v1(1))*int64_t(v2(0)); double cross = double(v1(0))*double(v2(1)) - double(v1(1))*double(v2(0));
float angle = float(atan2(double(cross), double(dot))); float angle = float(atan2(cross, dot));
angles[idx_curr] = angle; angles[idx_curr] = angle;
} }

View File

@ -25,14 +25,24 @@ enum Orientation
// As the points are limited to 30 bits + signum, // As the points are limited to 30 bits + signum,
// the temporaries u, v, w are limited to 61 bits + signum, // the temporaries u, v, w are limited to 61 bits + signum,
// and d is limited to 63 bits + signum and we are good. // and d is limited to 63 bits + signum and we are good.
static inline Orientation orient(const Point &a, const Point &b, const Point &c) //note: now coord_t is int64_t, so the algorithm is now adjusted to fallback to double is too big.
{ static inline Orientation orient(const Point &a, const Point &b, const Point &c) {
// BOOST_STATIC_ASSERT(sizeof(coord_t) * 2 == sizeof(int64_t)); // BOOST_STATIC_ASSERT(sizeof(coord_t) * 2 == sizeof(int64_t));
int64_t u = int64_t(b(0)) * int64_t(c(1)) - int64_t(b(1)) * int64_t(c(0)); // BOOST_STATIC_ASSERT(sizeof(coord_t) == sizeof(int64_t));
int64_t v = int64_t(a(0)) * int64_t(c(1)) - int64_t(a(1)) * int64_t(c(0)); if (a.x() <= 0xffffffff && b.x() <= 0xffffffff && c.x() <= 0xffffffff &&
int64_t w = int64_t(a(0)) * int64_t(b(1)) - int64_t(a(1)) * int64_t(b(0)); a.y() <= 0xffffffff && b.y() <= 0xffffffff && c.y() <= 0xffffffff) {
int64_t d = u - v + w; int64_t u = int64_t(b(0)) * int64_t(c(1)) - int64_t(b(1)) * int64_t(c(0));
return (d > 0) ? ORIENTATION_CCW : ((d == 0) ? ORIENTATION_COLINEAR : ORIENTATION_CW); int64_t v = int64_t(a(0)) * int64_t(c(1)) - int64_t(a(1)) * int64_t(c(0));
int64_t w = int64_t(a(0)) * int64_t(b(1)) - int64_t(a(1)) * int64_t(b(0));
int64_t d = u - v + w;
return (d > 0) ? ORIENTATION_CCW : ((d == 0) ? ORIENTATION_COLINEAR : ORIENTATION_CW);
} else {
double u = double(b(0)) * double(c(1)) - double(b(1)) * double(c(0));
double v = double(a(0)) * double(c(1)) - double(a(1)) * double(c(0));
double w = double(a(0)) * double(b(1)) - double(a(1)) * double(b(0));
double d = u - v + w;
return (d > 0) ? ORIENTATION_CCW : ((d == 0) ? ORIENTATION_COLINEAR : ORIENTATION_CW);
}
} }
// Return orientation of the polygon by checking orientation of the left bottom corner of the polygon // Return orientation of the polygon by checking orientation of the left bottom corner of the polygon

View File

@ -20,7 +20,7 @@
#include "Technologies.hpp" #include "Technologies.hpp"
typedef int32_t coord_t; typedef int64_t coord_t;
typedef double coordf_t; typedef double coordf_t;
//FIXME This epsilon value is used for many non-related purposes: //FIXME This epsilon value is used for many non-related purposes: