[tessellator] Direct comparison of result in edges_compare_x_for_y.

We need to compare the x-coordinate of a line at a for a particular y,
without loss of precision.

The x-coordinate along an edge for a given y is:
  X = A_x + (Y - A_y) * A_dx / A_dy

So the inequality we wish to test is:
  A_x + (Y - A_y) * A_dx / A_dy -?- B_x + (Y - B_y) * B_dx / B_dy,
where -?- is our inequality operator.

By construction we know that A_dy and B_dy (and (Y - A_y), (Y - B_y)) are
all positive, so we can rearrange it thus without causing a sign
change:
  A_dy * B_dy * (A_x - B_x) -?- (Y - B_y) * B_dx * A_dy
                                - (Y - A_y) * A_dx * B_dy

Given the assumption that all the deltas fit within 32 bits, we can compute
this comparison directly using 128 bit arithmetic.
This commit is contained in:
Chris Wilson 2008-10-01 12:18:16 +01:00
parent 7db03ac68c
commit ab23c29953
2 changed files with 55 additions and 41 deletions

View file

@ -219,42 +219,66 @@ _slope_compare (cairo_bo_edge_t *a,
}
}
static cairo_quorem64_t
edge_x_for_y (cairo_bo_edge_t *edge,
int32_t y)
/*
* We need to compare the x-coordinate of a line for a particular y,
* without loss of precision.
*
* The x-coordinate along an edge for a given y is:
* X = A_x + (Y - A_y) * A_dx / A_dy
*
* So the inequality we wish to test is:
* A_x + (Y - A_y) * A_dx / A_dy -?- B_x + (Y - B_y) * B_dx / B_dy,
* where -?- is our inequality operator.
*
* By construction, we know that A_dy and B_dy (and (Y - A_y), (Y - B_y)) are
* all positive, so we can rearrange it thus without causing a sign change:
* A_dy * B_dy * (A_x - B_x) -?- (Y - B_y) * B_dx * A_dy
* - (Y - A_y) * A_dx * B_dy
*
* Given the assumption that all the deltas fit within 32 bits, we can compute
* this comparison directly using 128 bit arithmetic.
*
* (And put the burden of the work on developing fast 128 bit ops, which are
* required throughout the tessellator.)
*
* See the similar discussion for _slope_compare().
*/
static int
edges_compare_x_for_y (const cairo_bo_edge_t *a,
const cairo_bo_edge_t *b,
int32_t y)
{
/* XXX: We're assuming here that dx and dy will still fit in 32
* bits. That's not true in general as there could be overflow. We
* should prevent that before the tessellation algorithm
* begins.
*/
int32_t dx = edge->bottom.x - edge->top.x;
int32_t dy = edge->bottom.y - edge->top.y;
int64_t numerator;
cairo_quorem64_t quorem;
int32_t adx, ady;
int32_t bdx, bdy;
cairo_int128_t L, R;
if (edge->middle.y == y) {
quorem.quo = edge->middle.x;
quorem.rem = 0;
return quorem;
}
if (edge->bottom.y == y) {
quorem.quo = edge->bottom.x;
quorem.rem = 0;
return quorem;
}
if (dy == 0) {
quorem.quo = _cairo_int32_to_int64 (edge->top.x);
quorem.rem = 0;
return quorem;
}
adx = a->bottom.x - a->top.x;
ady = a->bottom.y - a->top.y;
/* edge->top.x + (y - edge->top.y) * dx / dy */
numerator = _cairo_int32x32_64_mul ((y - edge->top.y), dx);
quorem = _cairo_int64_divrem (numerator, dy);
quorem.quo += edge->top.x;
bdx = b->bottom.x - b->top.x;
bdy = b->bottom.y - b->top.y;
return quorem;
L = _cairo_int64x32_128_mul (_cairo_int32x32_64_mul (ady, bdy),
a->top.x - b->top.x);
R = _cairo_int128_sub (_cairo_int64x32_128_mul (_cairo_int32x32_64_mul (bdx,
ady),
y - b->top.y),
_cairo_int64x32_128_mul (_cairo_int32x32_64_mul (adx,
bdy),
y - a->top.y));
/* return _cairo_int128_cmp (L, R); */
if (_cairo_int128_lt (L, R))
return -1;
if (_cairo_int128_gt (L, R))
return 1;
return 0;
}
static int
@ -262,8 +286,6 @@ _cairo_bo_sweep_line_compare_edges (cairo_bo_sweep_line_t *sweep_line,
cairo_bo_edge_t *a,
cairo_bo_edge_t *b)
{
cairo_quorem64_t ax;
cairo_quorem64_t bx;
int cmp;
if (a == b)
@ -292,18 +314,9 @@ _cairo_bo_sweep_line_compare_edges (cairo_bo_sweep_line_t *sweep_line,
if (amin > bmax) return +1;
}
ax = edge_x_for_y (a, sweep_line->current_y);
bx = edge_x_for_y (b, sweep_line->current_y);
if (ax.quo > bx.quo)
return 1;
else if (ax.quo < bx.quo)
return -1;
/* Quotients are identical, test remainder. */
if (ax.rem > bx.rem)
return 1;
else if (ax.rem < bx.rem)
return -1;
cmp = edges_compare_x_for_y (a, b, sweep_line->current_y);
if (cmp)
return cmp;
/* The two edges intersect exactly at y, so fall back on slope
* comparison. We know that this compare_edges function will be

View file

@ -181,6 +181,7 @@ cairo_int128_t I _cairo_int64_to_int128 (cairo_int64_t i);
#define _cairo_int128_sub(a,b) _cairo_uint128_sub(a,b)
#define _cairo_int128_mul(a,b) _cairo_uint128_mul(a,b)
cairo_int128_t I _cairo_int64x64_128_mul (cairo_int64_t a, cairo_int64_t b);
#define _cairo_int64x32_128_mul(a, b) _cairo_int64x64_128_mul(a, _cairo_int32_to_int64(b))
#define _cairo_int128_lsl(a,b) _cairo_uint128_lsl(a,b)
#define _cairo_int128_rsl(a,b) _cairo_uint128_rsl(a,b)
#define _cairo_int128_rsa(a,b) _cairo_uint128_rsa(a,b)