/* See LICENSE file for copyright and license details. */
#include "common.h"
#ifndef TEST
/*
* The implementation in this fill is not supposed to be fast,
* it is supposed to be absolutely correct. The implementation
* shall be as clear and obvious as possible. It is important
* that everything written to cells is as precise as possible
* and that the function behavious correctly when outside the
* raster. This implementation shall be usable to generate an
* as perfect as possible rasterisation and shall be usable to
* test the correctness and accuracy of the fast implementation.
*/
static void
draw_cell(RASTER *restrict raster, double x1, double y1, double x2, double y2, double dx, double dy,
double cxmin, double cxmax, double rxmin, double rxmax, double rymin, double rymax, ssize_t x, size_t cellY)
{
double cy1, cy2, ydir = (double)SIGNUM(dy);
double rightedge, opposite, midadjacent;
size_t cell;
rightedge = cxmax;
cxmin = fmax(cxmin, rxmin);
cxmax = fmin(cxmax, rxmax);
cy1 = dx ? y1 + (cxmin - x1) * dy / dx : y1;
cy2 = dx ? y2 + (cxmax - x2) * dy / dx : y2;
cy1 = fmin(fmax(cy1, rymin), rymax);
cy2 = fmin(fmax(cy2, rymin), rymax);
opposite = ydir * fabs(cy2 - cy1);
midadjacent = (cxmin + cxmax) / 2.0;
if (x < 0) {
raster->cells[cellY].opposite_coverage += opposite;
raster->cells[cellY].cell_coverage += opposite;
} else if ((size_t)x >= raster->width) {
cell = cellY + (raster->width - 1);
raster->cells[cell].opposite_coverage += opposite;
} else {
cell = cellY + (size_t)x;
raster->cells[cell].opposite_coverage += opposite;
raster->cells[cell].cell_coverage += opposite * fabs(rightedge - midadjacent);
}
}
void
draw_linear_bezier_reference(RASTER *restrict raster, double x1, double y1, double x2, double y2)
{
ssize_t iwidth = (ssize_t)raster->width;
double width = (double)raster->width;
double height = (double)raster->height;
double dx = x2 - x1;
double dy = y2 - y1;
double fymin = fmin(fmax(fmin(y1, y2), 0.0), height);
double fymax = fmin(fmax(fmax(y1, y2), 0.0), height);
double floor_fymin = floor(fymin);
double ceil_fymax = ceil(fymax);
size_t ymin = (size_t)floor_fymin;
size_t ymax = (size_t)ceil_fymax;
size_t y, cellY;
double rowx1, rowx2;
double rxmin, rxmax, floor_rxmin;
double rymax, rymin, ceil_rxmax;
ssize_t x, xmin, xmax, xmin_cropped, xmax_cropped;
double cxmin, cxmax;
ymax += (size_t)(ymax == ymin);
for (y = ymin; y < ymax; y++) {
rymin = (double)y;
rymax = (double)y + 1.0;
rymin = fmax(rymin, fymin);
rymax = fmin(rymax, fymax);
rowx1 = dy ? x1 + (rymin - y1) * dx / dy : x1;
rowx2 = dy ? x2 + (rymax - y2) * dx / dy : x2;
rxmin = fmin(rowx1, rowx2);
rxmax = fmax(rowx1, rowx2);
cellY = y * raster->width;
floor_rxmin = floor(rxmin);
ceil_rxmax = ceil(rxmax);
xmin = (ssize_t)floor_rxmin;
xmax = (ssize_t)ceil_rxmax;
xmax += (ssize_t)(xmax == xmin);
xmin_cropped = xmin;
xmax_cropped = xmax;
if (xmin < 0.0) {
cxmin = (double)xmin;
cxmax = xmax < 0 ? (double)(xmax + 1) : 0.0;
draw_cell(raster, x1, y1, x2, y2, dx, dy, cxmin, cxmax, rxmin, rxmax, rymin, rymax, xmin, cellY);
xmin_cropped = 0;
}
if (xmax > iwidth) {
cxmin = fmax(floor_rxmin, width);
cxmax = (double)xmax;
draw_cell(raster, x1, y1, x2, y2, dx, dy, cxmin, cxmax, rxmin, rxmax, rymin, rymax, (ssize_t)cxmin, cellY);
xmax_cropped = iwidth;
}
for (x = xmin_cropped; x < xmax_cropped; x++) {
cxmin = (double)x;
cxmax = (double)x + 1.0;
draw_cell(raster, x1, y1, x2, y2, dx, dy, cxmin, cxmax, rxmin, rxmax, rymin, rymax, x, cellY);
}
}
}
#else
#define CHECK_ZEROED()\
do {\
size_t i__;\
for (i__ = 0; i__ < 100; i__++) {\
ASSERT(!raster->cells[i__].cell_coverage);\
ASSERT(!raster->cells[i__].opposite_coverage);\
}\
} while (0)
#define CHECK_CELL(X, Y, AREA, SHADOW)\
do {\
size_t i__ = (size_t)(Y) * raster->width + (size_t)(X);\
ASSERT(raster->cells[i__].cell_coverage == (AREA));\
ASSERT(raster->cells[i__].opposite_coverage == (SHADOW));\
raster->cells[i__].cell_coverage = 0;\
raster->cells[i__].opposite_coverage = 0;\
} while (0)
#define SUB_CELL(X, Y, AREA, SHADOW)\
do {\
size_t i__ = (size_t)(Y) * raster->width + (size_t)(X);\
raster->cells[i__].cell_coverage -= (AREA);\
raster->cells[i__].opposite_coverage -= (SHADOW);\
} while (0)
static RASTER *raster;
static void
draw(double x1, double y1, double x2, double y2)
{
alarm(5);
draw_linear_bezier_reference(raster, x1, y1, x2, y2);
alarm(0);
}
static void
check_horizontal(void)
{
ssize_t y, x1, x2;
double fy, fx1, fx2;
for (y = -2; y < 12; y++) {
for (x1 = -2; x1 < 12; x1++) {
for (x2 = -2; x2 < 12; x2++) {
fy = (double)y;
fx1 = (double)x1;
fx2 = (double)x2;
draw(fx1, fy, fx2, fy);
CHECK_ZEROED();
draw(fx1 + 0.25, fy, fx2 + 0.25, fy);
CHECK_ZEROED();
draw(fx1, fy + 0.25, fx2, fy + 0.25);
CHECK_ZEROED();
draw(fx1 + 0.25, fy + 0.25, fx2 + 0.25, fy + 0.25);
CHECK_ZEROED();
}
}
}
}
static void
check_vertical(void)
{
ssize_t x, y1, y2, i, ydir;
double yhit;
double fx, fy1, fy2;
for (x = -2; x < 12; x++) {
for (y1 = -2; y1 < 12; y1++) {
for (y2 = -2; y2 < 12; y2++) {
ydir = SIGNUM(y2 - y1);
fx = (double)x;
fy1 = (double)y1;
fy2 = (double)y2;
draw(fx, fy1, fx, fy2);
for (i = 0; i < 10; i++) {
if (i < y1 && i < y2)
continue;
if (i >= y1 && i >= y2)
continue;
CHECK_CELL(x < 0 ? 0 : x < 10 ? x : 9,
i,
x < 0 ? (double)ydir : x < 10 ? (double)ydir : 0,
(double)ydir);
}
CHECK_ZEROED();
draw(fx + 0.25, fy1, fx + 0.25, fy2);
for (i = 0; i < 10; i++) {
if (i < y1 && i < y2)
continue;
if (i >= y1 && i >= y2)
continue;
CHECK_CELL(x < 0 ? 0 : x < 10 ? x : 9,
i,
x < 0 ? (double)ydir : x < 10 ? (double)ydir * 0.75 : 0,
(double)ydir);
}
CHECK_ZEROED();
draw(fx, fy1 + 0.25, fx, fy2 + 0.25);
for (i = 0; i < 10; i++) {
if (i < y1 + (y1 > y2) && i < y2 + (y2 > y1))
continue;
if (i >= y1 + (y1 > y2) && i >= y2 + (y2 > y1))
continue;
yhit = i == (y1 < y2 ? y1 : y2) ? 0.75 : i == (y2 > y1 ? y2 : y1) ? 0.25 : 1.0;
CHECK_CELL(x < 0 ? 0 : x < 10 ? x : 9,
i,
(x < 0 ? (double)ydir : x < 10 ? (double)ydir : 0) * yhit,
(double)ydir * yhit);
}
CHECK_ZEROED();
draw(fx + 0.25, fy1 + 0.25, fx + 0.25, fy2 + 0.25);
for (i = 0; i < 10; i++) {
if (i < y1 + (y1 > y2) && i < y2 + (y2 > y1))
continue;
if (i >= y1 + (y1 > y2) && i >= y2 + (y2 > y1))
continue;
yhit = i == (y1 < y2 ? y1 : y2) ? 0.75 : i == (y2 > y1 ? y2 : y1) ? 0.25 : 1.0;
CHECK_CELL(x < 0 ? 0 : x < 10 ? x : 9,
i,
(x < 0 ? (double)ydir : x < 10 ? (double)ydir * 0.75 : 0) * yhit,
(double)ydir * yhit);
}
CHECK_ZEROED();
}
}
}
}
static void
check_11_diagonal(ssize_t xshift, ssize_t yshift)
{
ssize_t i1, i2, i, imin, imax, x, y;
double fi1, fi2, idir;
double fxshift = (double)xshift;
double fyshift = (double)yshift;
for (i1 = -2; i1 < 12; i1++) {
for (i2 = -2; i2 < 12; i2++) {
imin = i1 < i2 ? i1 : i2;
imax = i1 > i2 ? i1 : i2;
fi1 = (double)i1;
fi2 = (double)i2;
idir = SIGNUM(fi2 - fi1);
draw(fi1 + fxshift, fi1 + fyshift, fi2 + fxshift, fi2 + fyshift);
for (i = imin; i < imax; i++) {
x = i + xshift;
y = i + yshift;
if (y < 0 || y >= 10)
continue;
if (x >= 10)
CHECK_CELL(9, y, 0, idir);
else if (x >= 0)
CHECK_CELL(x, y, 0.5 * idir, idir);
else
CHECK_CELL(0, y, idir, idir);
}
CHECK_ZEROED();
draw(fi1 + fxshift, fi1 + fyshift + 0.25, fi2 + fxshift, fi2 + fyshift + 0.25);
for (i = imin; i < imax; i++) {
x = i + xshift;
y = i + yshift;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, idir * 0.75);
else if (x >= 0)
SUB_CELL(x, y, idir * (0.75 - 0.75 * 0.75 / 2), idir * 0.75);
else
SUB_CELL(0, y, idir * 0.75, idir * 0.75);
}
x = i + xshift;
y = i + yshift + 1;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, idir * 0.25);
else if (x >= 0)
SUB_CELL(x, y, idir * 0.25 * 0.25 / 2, idir * 0.25);
else
SUB_CELL(0, y, idir * 0.25, idir * 0.25);
}
}
CHECK_ZEROED();
draw(fi1 + fxshift + 0.25, fi1 + fyshift, fi2 + fxshift + 0.25, fi2 + fyshift);
for (i = imin; i < imax; i++) {
x = i + xshift;
y = i + yshift;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, idir * 0.75);
else if (x >= 0)
SUB_CELL(x, y, idir * 0.75 * 0.75 / 2, idir * 0.75);
else
SUB_CELL(0, y, idir * 0.75, idir * 0.75);
}
x = i + xshift + 1;
y = i + yshift;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, idir * 0.25);
else if (x >= 0)
SUB_CELL(x, y, idir * (0.25 - 0.25 * 0.25 / 2), idir * 0.25);
else
SUB_CELL(0, y, idir * 0.25, idir * 0.25);
}
}
CHECK_ZEROED();
draw(fi1 + fxshift + 0.25, fi1 + fyshift + 0.25, fi2 + fxshift + 0.25, fi2 + fyshift + 0.25);
for (i = imin; i < imax; i++) {
x = i + xshift;
y = i + yshift;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, idir * 0.75);
else if (x >= 0)
SUB_CELL(x, y, idir * 0.75 * 0.75 / 2, idir * 0.75);
else
SUB_CELL(0, y, idir * 0.75, idir * 0.75);
}
x = i + xshift + 1;
y = i + yshift + 1;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, idir * 0.25);
else if (x >= 0)
SUB_CELL(x, y, idir * (0.25 - 0.25 * 0.25 / 2), idir * 0.25);
else
SUB_CELL(0, y, idir * 0.25, idir * 0.25);
}
}
CHECK_ZEROED();
}
}
}
static void
check_11_antidiagonal(ssize_t xshift, ssize_t yshift)
{
ssize_t i1, i2, i, imin, imax, x, y;
double fi1, fi2, ydir;
double fxshift = (double)xshift;
double fyshift = (double)yshift;
for (i1 = -2; i1 < 12; i1++) {
for (i2 = -2; i2 < 12; i2++) {
imin = i1 < i2 ? i1 : i2;
imax = i1 > i2 ? i1 : i2;
fi1 = (double)i1;
fi2 = (double)i2;
ydir = -SIGNUM(fi2 - fi1);
draw(fi1 + fxshift, fi2 + fyshift, fi2 + fxshift, fi1 + fyshift);
for (i = imin; i < imax; i++) {
x = i + xshift;
y = imax - i + imin + yshift - 1;
if (y < 0 || y >= 10)
continue;
if (x >= 10)
CHECK_CELL(9, y, 0, ydir);
else if (x >= 0)
CHECK_CELL(x, y, 0.5 * ydir, ydir);
else
CHECK_CELL(0, y, ydir, ydir);
}
CHECK_ZEROED();
draw(fi1 + fxshift, fi2 + fyshift + 0.25, fi2 + fxshift, fi1 + fyshift + 0.25);
for (i = imin; i < imax; i++) {
x = i + xshift;
y = imax - i + imin + yshift - 1;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, ydir * 0.75);
else if (x >= 0)
SUB_CELL(x, y, ydir * 0.75 * 0.75 / 2, ydir * 0.75);
else
SUB_CELL(0, y, ydir * 0.75, ydir * 0.75);
}
x = i + xshift;
y = imax - i + imin + yshift - 1 + 1;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, ydir * 0.25);
else if (x >= 0)
SUB_CELL(x, y, ydir * (0.25 - 0.25 * 0.25 / 2), ydir * 0.25);
else
SUB_CELL(0, y, ydir * 0.25, ydir * 0.25);
}
}
CHECK_ZEROED();
draw(fi1 + fxshift + 0.25, fi2 + fyshift, fi2 + fxshift + 0.25, fi1 + fyshift);
for (i = imin; i < imax; i++) {
x = i + xshift;
y = imax - i + imin + yshift - 1;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, ydir * 0.75);
else if (x >= 0)
SUB_CELL(x, y, ydir * 0.75 * 0.75 / 2, ydir * 0.75);
else
SUB_CELL(0, y, ydir * 0.75, ydir * 0.75);
}
x = i + xshift + 1;
y = imax - i + imin + yshift - 1;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, ydir * 0.25);
else if (x >= 0)
SUB_CELL(x, y, ydir * (0.25 - 0.25 * 0.25 / 2), ydir * 0.25);
else
SUB_CELL(0, y, ydir * 0.25, ydir * 0.25);
}
}
CHECK_ZEROED();
draw(fi1 + fxshift + 0.25, fi2 + fyshift + 0.25, fi2 + fxshift + 0.25, fi1 + fyshift + 0.25);
for (i = imin; i < imax; i++) {
x = i + xshift;
y = imax - i + imin + yshift - 1;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, ydir * 0.5);
else if (x >= 0)
SUB_CELL(x, y, ydir / 8.0, ydir * 0.5);
else
SUB_CELL(0, y, ydir * 0.5, ydir * 0.5);
}
x = i + xshift;
y = imax - i + imin + yshift - 1 + 1;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, ydir * 0.25);
else if (x >= 0)
SUB_CELL(x, y, ydir * (2.5 / 16.0), ydir * 0.25);
else
SUB_CELL(0, y, ydir * 0.25, ydir * 0.25);
}
x = i + xshift + 1;
y = imax - i + imin + yshift - 1;
if (0 <= y && y < 10) {
if (x >= 10)
SUB_CELL(9, y, 0, ydir * 0.25);
else if (x >= 0)
SUB_CELL(x, y, ydir * (3.5 / 16.0), ydir * 0.25);
else
SUB_CELL(0, y, ydir * 0.25, ydir * 0.25);
}
}
CHECK_ZEROED();
}
}
}
static double
frand(double a, double b)
{
int i = rand();
double f = (double)i / (double)RAND_MAX;
double max = fmax(a, b);
double min = fmin(a, b);
return fma(f, max - min, min);
}
static void
try_random_line(double xmin, double xmax, double ymin, double ymax)
{
double x1 = frand(xmin, xmax);
double y1 = frand(ymin, ymax);
double x2 = frand(xmin, xmax);
double y2 = frand(ymin, ymax);
draw(x1, y1, x2, y2);
}
static void
check_random_line(double rand_xmin, double rand_xmax, double rand_ymin, double rand_ymax)
{
double x1 = frand(rand_xmin, rand_xmax);
double y1 = frand(rand_ymin, rand_ymax);
double x2 = frand(rand_xmin, rand_xmax);
double y2 = frand(rand_ymin, rand_ymax);
double dx = x2 - x1;
double dy = y2 - y1;
double ymin = fmin(y1, y2);
double ymax = fmax(y1, y2);
double xmin = fmin(x1, x2);
double xmax = fmax(x1, x2);
double ydir = SIGNUM(y2 - y1);
double inklevel, area, expected_area, inkdir, dist;
size_t y, x, jumps;
int hit;
draw(x1, y1, x2, y2);
area = 0;
for (y = 0; y < 10; y++) {
inklevel = 0;
jumps = 0;
inkdir = 0;
for (x = 0; x < 10; x++) {
area += raster->cells[y * 10 + x].cell_coverage + inklevel;
inklevel += raster->cells[y * 10 + x].opposite_coverage;
if (inkdir == 0) {
inkdir = raster->cells[y * 10 + x].opposite_coverage;
} else if (inkdir < 0 && raster->cells[y * 10 + x].opposite_coverage > 0) {
ASSERT(!jumps);
jumps += 1;
inkdir = raster->cells[y * 10 + x].opposite_coverage;
} else if (inkdir > 0 && raster->cells[y * 10 + x].opposite_coverage < 0) {
ASSERT(!jumps);
jumps += 1;
inkdir = raster->cells[y * 10 + x].opposite_coverage;
}
if (!dy)
continue;
hit = (double)x >= floor(xmin) && (double)x <= ceil(xmax);
hit |= !x && xmin < 0;
hit &= (double)y >= floor(ymin) && (double)y <= ceil(ymax);
ASSERT(hit || !raster->cells[y * 10 + x].cell_coverage);
hit = (double)x >= floor(xmin) && (double)x <= ceil(xmax);
hit &= (double)y >= floor(ymin) && (double)y <= ceil(ymax);
if (!hit || !dx || !x)
continue;
dist = ((double)x+0.5)*dy/dx - ((double)y+0.5) + (y1 - x1*dy/dx);
dist /= hypot(dy/dx, 1);
dist = fabs(dist);
if (dist > 0.708)
ASSERT(!raster->cells[y * 10 + x].cell_coverage);
}
if ((double)y < floor(ymin) || (double)y >= ceil(ymax))
ASSERT(inklevel == 0.0);
else if ((double)y == floor(ymin))
ASSERT(eq(inklevel, ydir * (fmin(ymax, (double)y + 1.0) - ymin)));
else if ((double)y == floor(ymax))
ASSERT(eq(inklevel, ydir * fmod(ymax, 1.0)));
else
ASSERT(eq(inklevel, ydir));
}
if (ymax <= 0 || ymin >= 10 || !dy || fmin(x1, x2) >= 10) {
expected_area = 0;
} else if (fmax(x1, x2) < 0) {
expected_area = ydir * 10 * (fmin(10, ymax) - fmax(ymin, 0));
} else if (!dx) {
expected_area = ydir * (10 - fmin(fmax(0, x1), 10)) * (fmin(10, ymax) - fmax(ymin, 0));
} else {
#define y(x) (y1 + ((x) - x1) * k)
#define x(y) (x1 + ((y) - y1) / k)
double k = dy / dx;
double a = xmin;
double b = xmax;
double mid;
a = y(a) > 10 ? x(10) : y(a) < 0 ? x(0) : a;
b = y(b) > 10 ? x(10) : y(b) < 0 ? x(0) : b;
b = fmin(b, 10);
a = fmin(a, 10);
mid = (a + b) / 2;
if (mid < 0) {
expected_area = (10 - 0) * fabs(y(b) - y(a));
if (b > 0)
expected_area -= fabs(y(b) - y(0)) * (b - 0) / 2;
} else {
expected_area = (10 - mid) * fabs(y(b) - y(a));
if (a < 0)
expected_area -= fabs(y(0) - y(a)) * (0 - a) / 2;
}
expected_area *= ydir;
#undef y
#undef x
}
ASSERT(eq(area, expected_area));
rtgrpblib_reset_raster(raster, 10, 10);
}
int
main(void)
{
size_t i;
raster = rtgrpblib_create_raster(10, 10);
ASSERT(raster);
srand((unsigned)(uintptr_t)raster);
check_horizontal();
check_vertical();
check_11_diagonal(0, 0);
check_11_diagonal(0, 2);
check_11_diagonal(0, -2);
check_11_diagonal(2, 0);
check_11_diagonal(2, 2);
check_11_diagonal(2, -2);
check_11_diagonal(-2, 0);
check_11_diagonal(-2, 2);
check_11_diagonal(-2, -2);
check_11_antidiagonal(0, 0);
check_11_antidiagonal(0, 2);
check_11_antidiagonal(0, -2);
check_11_antidiagonal(2, 0);
check_11_antidiagonal(2, 2);
check_11_antidiagonal(2, -2);
check_11_antidiagonal(-2, 0);
check_11_antidiagonal(-2, 2);
check_11_antidiagonal(-2, -2);
for (i = 0; i < 100000UL; i++) {
try_random_line(0., 1., 0., 1.);
try_random_line(-5., 15., -5., 15.);
}
rtgrpblib_reset_raster(raster, 10, 10);
for (i = 0; i < 10000UL; i++) {
check_random_line(0., 1., 0., 1.);
check_random_line(-5., 15., -5., 15.);
}
free(raster);
return 0;
}
#endif
