draw_linear_bezier_reference: check for random lines

Signed-off-by: Mattias Andrée <maandree@kth.se>

1 files changed, 131 insertions, 10 deletions

diff --git a/draw_linear_bezier_reference.c b/draw_linear_bezier_reference.c
index a395c18..271e0a8 100644
--- a/draw_linear_bezier_reference.c
+++ b/draw_linear_bezier_reference.c
@@ -101,6 +101,7 @@ draw_linear_bezier_reference(RASTER *restrict raster, double x1, double y1, doub
 
 #else
 
+
 #define CHECK_ZEROED()\
 	do {\
 		size_t i__;\
@@ -130,7 +131,7 @@ draw_linear_bezier_reference(RASTER *restrict raster, double x1, double y1, doub
 static RASTER *raster;
 
 
-#if 0
+#if 1
 static void
 print_raster(void)
 {
@@ -497,21 +498,134 @@ check_11_antidiagonal(ssize_t xshift, ssize_t yshift)
 
 
 static double
-frand(void)
+frand(double a, double b)
 {
 	int i = rand();
-	return (double)i / (double)RAND_MAX;
+	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
-try_random_line(void)
+check_random_line(double rand_xmin, double rand_xmax, double rand_ymin, double rand_ymax)
 {
-	double x1 = frand();
-	double y1 = frand();
-	double x2 = frand();
-	double y2 = frand();
+	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  = x >= floor(xmin) && x <= ceil(xmax);
+			hit |= !x && xmin < 0;
+			hit &= y >= floor(ymin) && y <= ceil(ymax);
+			ASSERT(hit || !raster->cells[y * 10 + x].cell_coverage);
+
+			hit  = x >= floor(xmin) && x <= ceil(xmax);
+			hit &= y >= floor(ymin) && 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);
 }
 
 
@@ -547,10 +661,17 @@ main(void)
 	check_11_antidiagonal(-2, 2);
 	check_11_antidiagonal(-2, -2);
 
-	for (i = 0; i < 100000UL; i++)
-		try_random_line();
+	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;
 }