From 9cd9d0f3de06099dff988c4f9abdcd5eb61633c3 Mon Sep 17 00:00:00 2001 From: Michael Ludwig Date: Tue, 3 Aug 2021 15:05:26 -0400 Subject: [PATCH] Fix nearly flat edge detection in GrTriangulator Adds a GM reproducing triangulator failure. The stroker generated geometry shown in the GM had two edges intersecting on the inner contour. One of the edges was practically vertical, such that its top and bottom vertices X coordinates differed by 1ulp (approx 4.8e-7 when x = 7.5). The intersection was computed correctly in double but then clamped to exactly 7.5 when stored back in float at (7.5, ~11), which is still a reasonable intersection point along the edge represented by (7.5,0) to (7.5000048,~68). However, with a sorting direction of horizontal for the vertices, vertices are sorted on increasing x and then *decreasing* y. The computed intersection is then technically above the edge's top vertex. In the original code, the edge was not considered flat, so it assumes that the intersection vertex is below the top vertex and above the bottom vertex. The issue was with the nearly_flat function only checking the difference against machine epsilon without taking the magnitude of the vertices into account. However, just changing nearly_flat to use ulp instead of epsilon led to infinite recursion and asserts on a number of tests due to inaccuracies in edge intersection tests. To address these issues, I reworked clamp() to work on X and Y independently, which avoids the need to check nearly_flat entirely. Additionally, edge intersection detects when it might be inaccurate due to really large vertex coords and recursively splits to reject false positives. Bug: skia:12244 Change-Id: I9a92a163e8c53af799332f1f9369fccd3a9fbf31 Reviewed-on: https://skia-review.googlesource.com/c/skia/+/432196 Commit-Queue: Michael Ludwig Reviewed-by: Robert Phillips --- gm/strokes.cpp | 26 ++++++ src/gpu/GrTriangulator.cpp | 167 ++++++++++++++++++++++++++++--------- src/gpu/GrTriangulator.h | 12 ++- 3 files changed, 162 insertions(+), 43 deletions(-) diff --git a/gm/strokes.cpp b/gm/strokes.cpp index 8673ac2596..455e30a247 100644 --- a/gm/strokes.cpp +++ b/gm/strokes.cpp @@ -618,3 +618,29 @@ DEF_SIMPLE_GM(inner_join_geometry, canvas, 1000, 700) { } } } + +DEF_SIMPLE_GM(skbug12244, canvas, 150, 150) { + // Should look like a stroked triangle; these vertices are the results of the SkStroker + // but we draw as a filled path in order to highlight that it's the GPU triangulating path + // renderer that's the source of the problem, and not the stroking operation. The original + // path was a simple: + // m(0,0), l(100, 40), l(0, 80), l(0,0) with a stroke width of 15px + SkPath path; + path.moveTo(2.7854299545288085938, -6.9635753631591796875); + path.lineTo( 120.194366455078125, 40); + path.lineTo(-7.5000004768371582031, 91.07775115966796875); + path.lineTo(-7.5000004768371582031, -11.077748298645019531); + path.lineTo(2.7854299545288085938, -6.9635753631591796875); + path.moveTo(-2.7854299545288085938, 6.9635753631591796875); + path.lineTo( 0, 0); + path.lineTo( 7.5, 0); + path.lineTo(7.5000004768371582031, 68.92224884033203125); + path.lineTo( 79.805633544921875, 40); + path.lineTo(-2.7854299545288085938, 6.9635753631591796875); + + SkPaint p; + p.setColor(SK_ColorGREEN); + + canvas->translate(20.f, 20.f); + canvas->drawPath(path, p); +} diff --git a/src/gpu/GrTriangulator.cpp b/src/gpu/GrTriangulator.cpp index 35f2131682..95342cd043 100644 --- a/src/gpu/GrTriangulator.cpp +++ b/src/gpu/GrTriangulator.cpp @@ -145,34 +145,127 @@ bool GrTriangulator::Line::intersect(const Line& other, SkPoint* point) const { double scale = 1.0 / denom; point->fX = double_to_clamped_scalar((fB * other.fC - other.fB * fC) * scale); point->fY = double_to_clamped_scalar((other.fA * fC - fA * other.fC) * scale); - round(point); + round(point); return point->isFinite(); } -bool GrTriangulator::Edge::intersect(const Edge& other, SkPoint* p, uint8_t* alpha) const { - TESS_LOG("intersecting %g -> %g with %g -> %g\n", - fTop->fID, fBottom->fID, other.fTop->fID, other.fBottom->fID); - if (fTop == other.fTop || fBottom == other.fBottom) { +// If the edge's vertices differ by many orders of magnitude, the computed line equation can have +// significant error in its distance and intersection tests. To avoid this, we recursively subdivide +// long edges and effectively perform a binary search to perform a more accurate intersection test. +static bool edge_line_needs_recursion(const SkPoint& p0, const SkPoint& p1) { + // ilogbf(0) returns an implementation-defined constant, but we are choosing to saturate + // negative exponents to 0 for comparisons sake. We're only trying to recurse on lines with + // very large coordinates. + int expDiffX = std::abs((std::abs(p0.fX) < 1.f ? 0 : std::ilogbf(p0.fX)) - + (std::abs(p1.fX) < 1.f ? 0 : std::ilogbf(p1.fX))); + int expDiffY = std::abs((std::abs(p0.fY) < 1.f ? 0 : std::ilogbf(p0.fY)) - + (std::abs(p1.fY) < 1.f ? 0 : std::ilogbf(p1.fY))); + // Differ by more than 2^20, or roughly a factor of one million. + return expDiffX > 20 || expDiffY > 20; +} + +static bool recursive_edge_intersect(const Line& u, SkPoint u0, SkPoint u1, + const Line& v, SkPoint v0, SkPoint v1, + SkPoint* p, double* s, double* t) { + // First check if the bounding boxes of [u0,u1] intersects [v0,v1]. If they do not, then the + // two line segments cannot intersect in their domain (even if the lines themselves might). + // - don't use SkRect::intersect since the vertices aren't sorted and horiz/vertical lines + // appear as empty rects, which then never "intersect" according to SkRect. + if (std::min(u0.fX, u1.fX) > std::max(v0.fX, v1.fX) || + std::max(u0.fX, u1.fX) < std::min(v0.fX, v1.fX) || + std::min(u0.fY, u1.fY) > std::max(v0.fY, v1.fY) || + std::max(u0.fY, u1.fY) < std::min(v0.fY, v1.fY)) { return false; } - double denom = fLine.fA * other.fLine.fB - fLine.fB * other.fLine.fA; + + // Compute intersection based on current segment vertices; if an intersection is found but the + // vertices differ too much in magnitude, we recurse using the midpoint of the segment to + // reject false positives. We don't currently try to avoid false negatives (e.g. large magnitude + // line reports no intersection but there is one). + double denom = u.fA * v.fB - u.fB * v.fA; if (denom == 0.0) { return false; } - double dx = static_cast(other.fTop->fPoint.fX) - fTop->fPoint.fX; - double dy = static_cast(other.fTop->fPoint.fY) - fTop->fPoint.fY; - double sNumer = dy * other.fLine.fB + dx * other.fLine.fA; - double tNumer = dy * fLine.fB + dx * fLine.fA; + double dx = static_cast(v0.fX) - u0.fX; + double dy = static_cast(v0.fY) - u0.fY; + double sNumer = dy * v.fB + dx * v.fA; + double tNumer = dy * u.fB + dx * u.fA; // If (sNumer / denom) or (tNumer / denom) is not in [0..1], exit early. // This saves us doing the divide below unless absolutely necessary. if (denom > 0.0 ? (sNumer < 0.0 || sNumer > denom || tNumer < 0.0 || tNumer > denom) : (sNumer > 0.0 || sNumer < denom || tNumer > 0.0 || tNumer < denom)) { return false; } - double s = sNumer / denom; - SkASSERT(s >= 0.0 && s <= 1.0); - p->fX = double_to_clamped_scalar(fTop->fPoint.fX - s * fLine.fB); - p->fY = double_to_clamped_scalar(fTop->fPoint.fY + s * fLine.fA); + + *s = sNumer / denom; + *t = tNumer / denom; + SkASSERT(*s >= 0.0 && *s <= 1.0 && *t >= 0.0 && *t <= 1.0); + + const bool uNeedsSplit = edge_line_needs_recursion(u0, u1); + const bool vNeedsSplit = edge_line_needs_recursion(v0, v1); + if (!uNeedsSplit && !vNeedsSplit) { + p->fX = double_to_clamped_scalar(u0.fX - (*s) * u.fB); + p->fY = double_to_clamped_scalar(u0.fY + (*s) * u.fA); + return true; + } else { + double sScale = 1.0, sShift = 0.0; + double tScale = 1.0, tShift = 0.0; + + if (uNeedsSplit) { + SkPoint uM = {(float) (0.5 * u0.fX + 0.5 * u1.fX), + (float) (0.5 * u0.fY + 0.5 * u1.fY)}; + sScale = 0.5; + if (*s >= 0.5) { + u1 = uM; + sShift = 0.5; + } else { + u0 = uM; + } + } + if (vNeedsSplit) { + SkPoint vM = {(float) (0.5 * v0.fX + 0.5 * v1.fX), + (float) (0.5 * v0.fY + 0.5 * v1.fY)}; + tScale = 0.5; + if (*t >= 0.5) { + v1 = vM; + tShift = 0.5; + } else { + v0 = vM; + } + } + + // Just recompute both lines, even if only one was split; we're already in a slow path. + if (recursive_edge_intersect(Line(u0, u1), u0, u1, Line(v0, v1), v0, v1, p, s, t)) { + // Adjust s and t back to full range + *s = sScale * (*s) + sShift; + *t = tScale * (*t) + tShift; + return true; + } else { + // False positive + return false; + } + } +} + +bool GrTriangulator::Edge::intersect(const Edge& other, SkPoint* p, uint8_t* alpha) const { + TESS_LOG("intersecting %g -> %g with %g -> %g\n", + fTop->fID, fBottom->fID, other.fTop->fID, other.fBottom->fID); + if (fTop == other.fTop || fBottom == other.fBottom || + fTop == other.fBottom || fBottom == other.fTop) { + // If the two edges share a vertex by construction, they have already been split and + // shouldn't be considered "intersecting" anymore. + return false; + } + + double s, t; // needed to interpolate vertex alpha + const bool intersects = recursive_edge_intersect( + fLine, fTop->fPoint, fBottom->fPoint, + other.fLine, other.fTop->fPoint, other.fBottom->fPoint, + p, &s, &t); + if (!intersects) { + return false; + } + if (alpha) { if (fType == EdgeType::kInner || other.fType == EdgeType::kInner) { // If the intersection is on any interior edge, it needs to stay fully opaque or later @@ -186,7 +279,6 @@ bool GrTriangulator::Edge::intersect(const Edge& other, SkPoint* p, uint8_t* alp // Could be two connectors crossing, or a connector crossing an outer edge. // Take the max interpolated alpha SkASSERT(fType == EdgeType::kConnector || other.fType == EdgeType::kConnector); - double t = tNumer / denom; *alpha = std::max((1.0 - s) * fTop->fAlpha + s * fBottom->fAlpha, (1.0 - t) * other.fTop->fAlpha + t * other.fBottom->fAlpha); } @@ -912,22 +1004,21 @@ Vertex* GrTriangulator::makeSortedVertex(const SkPoint& p, uint8_t alpha, Vertex return v; } -// If an edge's top and bottom points differ only by 1/2 machine epsilon in the primary -// sort criterion, it may not be possible to split correctly, since there is no point which is -// below the top and above the bottom. This function detects that case. -static bool nearly_flat(const Comparator& c, Edge* edge) { - SkPoint diff = edge->fBottom->fPoint - edge->fTop->fPoint; - float primaryDiff = c.fDirection == Comparator::Direction::kHorizontal ? diff.fX : diff.fY; - return fabs(primaryDiff) <= std::numeric_limits::epsilon() && primaryDiff != 0.0f; -} - +// Clamps x and y coordinates independently, so the returned point will lie within the bounding +// box formed by the corners of 'min' and 'max' (although min/max here refer to the ordering +// imposed by 'c'). static SkPoint clamp(SkPoint p, SkPoint min, SkPoint max, const Comparator& c) { - if (c.sweep_lt(p, min)) { - return min; - } else if (c.sweep_lt(max, p)) { - return max; + if (c.fDirection == Comparator::Direction::kHorizontal) { + // With horizontal sorting, we know min.x <= max.x, but there's no relation between + // Y components unless min.x == max.x. + return {SkTPin(p.fX, min.fX, max.fX), + min.fY < max.fY ? SkTPin(p.fY, min.fY, max.fY) + : SkTPin(p.fY, max.fY, min.fY)}; } else { - return p; + // And with vertical sorting, we know Y's relation but not necessarily X's. + return {min.fX < max.fX ? SkTPin(p.fX, min.fX, max.fX) + : SkTPin(p.fX, max.fX, min.fX), + SkTPin(p.fY, min.fY, max.fY)}; } } @@ -968,17 +1059,13 @@ bool GrTriangulator::checkForIntersection(Edge* left, Edge* right, EdgeList* act while (top && c.sweep_lt(p, top->fPoint)) { top = top->fPrev; } - bool leftFlat = nearly_flat(c, left); - bool rightFlat = nearly_flat(c, right); - if (leftFlat && rightFlat) { - return false; - } - if (!leftFlat) { - p = clamp(p, left->fTop->fPoint, left->fBottom->fPoint, c); - } - if (!rightFlat) { - p = clamp(p, right->fTop->fPoint, right->fBottom->fPoint, c); - } + + // Always clamp the intersection to lie between the vertices of each segment, since + // in theory that's where the intersection is, but in reality, floating point error may + // have computed an intersection beyond a vertex's component(s). + p = clamp(p, left->fTop->fPoint, left->fBottom->fPoint, c); + p = clamp(p, right->fTop->fPoint, right->fBottom->fPoint, c); + if (coincident(p, left->fTop->fPoint)) { v = left->fTop; } else if (coincident(p, left->fBottom->fPoint)) { diff --git a/src/gpu/GrTriangulator.h b/src/gpu/GrTriangulator.h index 97f97e59cf..f1603aaf08 100644 --- a/src/gpu/GrTriangulator.h +++ b/src/gpu/GrTriangulator.h @@ -413,9 +413,15 @@ struct GrTriangulator::Edge { bool fUsedInLeftPoly; bool fUsedInRightPoly; Line fLine; - double dist(const SkPoint& p) const { return fLine.dist(p); } - bool isRightOf(Vertex* v) const { return fLine.dist(v->fPoint) < 0.0; } - bool isLeftOf(Vertex* v) const { return fLine.dist(v->fPoint) > 0.0; } + + double dist(const SkPoint& p) const { + // Coerce points coincident with the vertices to have dist = 0, since converting from + // a double intersection point back to float storage might construct a point that's no + // longer on the ideal line. + return (p == fTop->fPoint || p == fBottom->fPoint) ? 0.0 : fLine.dist(p); + } + bool isRightOf(Vertex* v) const { return this->dist(v->fPoint) < 0.0; } + bool isLeftOf(Vertex* v) const { return this->dist(v->fPoint) > 0.0; } void recompute() { fLine = Line(fTop, fBottom); } void insertAbove(Vertex*, const Comparator&); void insertBelow(Vertex*, const Comparator&);