LibGL: Use integer math in rasterizer coverage calculations

This makes the software rasterizer use integers for triangle coverage
calculations. The previously used floating point algorithm was not
precise enough in certain situations and showed gaps between triangles.

This is not yet subpixel accurate.

Userland/Libraries/LibGL/SoftwareRasterizer.cpp

@@ -7,23 +7,23 @@
 #include "SoftwareRasterizer.h"
 #include <AK/Function.h>
 #include <LibGfx/Painter.h>
+#include <LibGfx/Vector2.h>
+#include <LibGfx/Vector3.h>
 
 namespace GL {
 
-static constexpr size_t RASTERIZER_BLOCK_SIZE = 16;
+using IntVector2 = Gfx::Vector2<int>;
+using IntVector3 = Gfx::Vector3<int>;
 
-struct FloatVector2 {
-    float x;
-    float y;
-};
+static constexpr int RASTERIZER_BLOCK_SIZE = 16;
 
-constexpr static float triangle_area(const FloatVector2& a, const FloatVector2& b, const FloatVector2& c)
+constexpr static int edge_function(const IntVector2& a, const IntVector2& b, const IntVector2& c)
 {
-    return ((c.x - a.x) * (b.y - a.y) - (c.y - a.y) * (b.x - a.x)) / 2;
+    return ((c.x() - a.x()) * (b.y() - a.y()) - (c.y() - a.y()) * (b.x() - a.x()));
 }
 
 template<typename T>
-constexpr static T interpolate(const T& v0, const T& v1, const T& v2, const FloatVector4& barycentric_coords)
+constexpr static T interpolate(const T& v0, const T& v1, const T& v2, const FloatVector3& barycentric_coords)
 {
     return v0 * barycentric_coords.x() + v1 * barycentric_coords.y() + v2 * barycentric_coords.z();
 }
@@ -47,141 +47,131 @@ static void rasterize_triangle(const RasterizerOptions& options, Gfx::Bitmap& re
     VERIFY((render_target.height() % RASTERIZER_BLOCK_SIZE) == 0);
 
     // Calculate area of the triangle for later tests
-    FloatVector2 v0 = { triangle.vertices[0].x, triangle.vertices[0].y };
-    FloatVector2 v1 = { triangle.vertices[1].x, triangle.vertices[1].y };
-    FloatVector2 v2 = { triangle.vertices[2].x, triangle.vertices[2].y };
+    IntVector2 v0 { (int)triangle.vertices[0].x, (int)triangle.vertices[0].y };
+    IntVector2 v1 { (int)triangle.vertices[1].x, (int)triangle.vertices[1].y };
+    IntVector2 v2 { (int)triangle.vertices[2].x, (int)triangle.vertices[2].y };
 
-    float area = triangle_area(v0, v1, v2);
+    int area = edge_function(v0, v1, v2);
     if (area == 0)
         return;
 
-    float one_over_area = 1 / area;
+    float one_over_area = 1.0f / area;
 
     // Obey top-left rule:
     // This sets up "zero" for later pixel coverage tests.
     // Depending on where on the triangle the edge is located
-    // it is either tested against 0 or float epsilon, effectively
+    // it is either tested against 0 or 1, effectively
     // turning "< 0" into "<= 0"
-    float constexpr epsilon = AK::NumericLimits<float>::epsilon();
-    FloatVector4 zero { epsilon, epsilon, epsilon, 0.0f };
-    if (v1.y > v0.y || (v1.y == v0.y && v1.x < v0.x))
+    IntVector3 zero { 1, 1, 1 };
+    if (v1.y() > v0.y() || (v1.y() == v0.y() && v1.x() < v0.x()))
         zero.set_z(0);
-    if (v2.y > v1.y || (v2.y == v1.y && v2.x < v1.x))
+    if (v2.y() > v1.y() || (v2.y() == v1.y() && v2.x() < v1.x()))
         zero.set_x(0);
-    if (v0.y > v2.y || (v0.y == v2.y && v0.x < v2.x))
+    if (v0.y() > v2.y() || (v0.y() == v2.y() && v0.x() < v2.x()))
         zero.set_y(0);
 
-    // This function calculates the barycentric coordinates for the pixel relative to the triangle.
-    auto barycentric_coordinates = [v0, v1, v2, one_over_area](float x, float y) -> FloatVector4 {
-        FloatVector2 p { x, y };
+    // This function calculates the 3 edge values for the pixel relative to the triangle.
+    auto calculate_edge_values = [v0, v1, v2](const IntVector2& p) -> IntVector3 {
         return {
-            triangle_area(v1, v2, p) * one_over_area,
-            triangle_area(v2, v0, p) * one_over_area,
-            triangle_area(v0, v1, p) * one_over_area,
-            0.0f
+            edge_function(v1, v2, p),
+            edge_function(v2, v0, p),
+            edge_function(v0, v1, p),
         };
     };
 
-    // This function tests whether a point lies within the triangle
-    auto test_point = [zero](const FloatVector4& point) -> bool {
-        return point.x() >= zero.x()
-            && point.y() >= zero.y()
-            && point.z() >= zero.z();
+    // This function tests whether a point as identified by its 3 edge values lies within the triangle
+    auto test_point = [zero](const IntVector3& edges) -> bool {
+        return edges.x() >= zero.x()
+            && edges.y() >= zero.y()
+            && edges.z() >= zero.z();
     };
 
-    // Calculate bounds
-    FloatVector2 min { AK::min(v0.x, AK::min(v1.x, v2.x)), AK::min(v0.y, AK::min(v1.y, v2.y)) };
-    FloatVector2 max { AK::max(v0.x, AK::max(v1.x, v2.x)), AK::max(v0.y, AK::max(v1.y, v2.y)) };
-
     // Calculate block-based bounds
-    int iminx = floorf(min.x);
-    int iminy = floorf(min.y);
-    int imaxx = ceilf(max.x);
-    int imaxy = ceilf(max.y);
-
-    iminx = clamp(iminx, 0, render_target.width() - 1);
-    imaxx = clamp(imaxx, 0, render_target.width() - 1);
-    iminy = clamp(iminy, 0, render_target.height() - 1);
-    imaxy = clamp(imaxy, 0, render_target.height() - 1);
-
-    int bx0 = iminx / RASTERIZER_BLOCK_SIZE;
-    int bx1 = imaxx / RASTERIZER_BLOCK_SIZE + 1;
-    int by0 = iminy / RASTERIZER_BLOCK_SIZE;
-    int by1 = imaxy / RASTERIZER_BLOCK_SIZE + 1;
+    // clang-format off
+    const int bx0 = max(0,                      min(min(v0.x(), v1.x()), v2.x())                            ) / RASTERIZER_BLOCK_SIZE;
+    const int bx1 = min(render_target.width(),  max(max(v0.x(), v1.x()), v2.x()) + RASTERIZER_BLOCK_SIZE - 1) / RASTERIZER_BLOCK_SIZE;
+    const int by0 = max(0,                      min(min(v0.y(), v1.y()), v2.y())                            ) / RASTERIZER_BLOCK_SIZE;
+    const int by1 = min(render_target.height(), max(max(v0.y(), v1.y()), v2.y()) + RASTERIZER_BLOCK_SIZE - 1) / RASTERIZER_BLOCK_SIZE;
+    // clang-format on
 
     // Iterate over all blocks within the bounds of the triangle
     for (int by = by0; by < by1; by++) {
         for (int bx = bx0; bx < bx1; bx++) {
 
-            // The 4 block corners
-            int x0 = bx * RASTERIZER_BLOCK_SIZE;
-            int y0 = by * RASTERIZER_BLOCK_SIZE;
-            int x1 = bx * RASTERIZER_BLOCK_SIZE + RASTERIZER_BLOCK_SIZE;
-            int y1 = by * RASTERIZER_BLOCK_SIZE + RASTERIZER_BLOCK_SIZE;
-
-            // Barycentric coordinates of the 4 block corners
-            auto a = barycentric_coordinates(x0, y0);
-            auto b = barycentric_coordinates(x1, y0);
-            auto c = barycentric_coordinates(x0, y1);
-            auto d = barycentric_coordinates(x1, y1);
+            // Edge values of the 4 block corners
+            // clang-format off
+            auto b0 = calculate_edge_values({ bx * RASTERIZER_BLOCK_SIZE,                         by * RASTERIZER_BLOCK_SIZE });
+            auto b1 = calculate_edge_values({ bx * RASTERIZER_BLOCK_SIZE + RASTERIZER_BLOCK_SIZE, by * RASTERIZER_BLOCK_SIZE });
+            auto b2 = calculate_edge_values({ bx * RASTERIZER_BLOCK_SIZE,                         by * RASTERIZER_BLOCK_SIZE + RASTERIZER_BLOCK_SIZE });
+            auto b3 = calculate_edge_values({ bx * RASTERIZER_BLOCK_SIZE + RASTERIZER_BLOCK_SIZE, by * RASTERIZER_BLOCK_SIZE + RASTERIZER_BLOCK_SIZE });
+            // clang-format on
 
             // If the whole block is outside any of the triangle edges we can discard it completely
-            if ((a.x() < zero.x() && b.x() < zero.x() && c.x() < zero.x() && d.x() < zero.x())
-                || (a.y() < zero.y() && b.y() < zero.y() && c.y() < zero.y() && d.y() < zero.y())
-                || (a.z() < zero.z() && b.z() < zero.z() && c.z() < zero.z() && d.z() < zero.z()))
+            // We test this by and'ing the relevant edge function values together for all block corners
+            // and checking if the negative sign bit is set for all of them
+            if ((b0.x() & b1.x() & b2.x() & b3.x()) & 0x80000000)
+                continue;
+
+            if ((b0.y() & b1.y() & b2.y() & b3.y()) & 0x80000000)
+                continue;
+
+            if ((b0.z() & b1.z() & b2.z() & b3.z()) & 0x80000000)
                 continue;
 
-            // barycentric coordinate derivatives
-            auto dcdx = (b - a) / RASTERIZER_BLOCK_SIZE;
-            auto dcdy = (c - a) / RASTERIZER_BLOCK_SIZE;
+            // edge value derivatives
+            auto dbdx = (b1 - b0) / RASTERIZER_BLOCK_SIZE;
+            auto dbdy = (b2 - b0) / RASTERIZER_BLOCK_SIZE;
+            // step edge value after each horizontal span: 1 down, BLOCK_SIZE left
+            auto step_y = dbdy - dbdx * RASTERIZER_BLOCK_SIZE;
 
-            if (test_point(a) && test_point(b) && test_point(c) && test_point(d)) {
+            int x0 = bx * RASTERIZER_BLOCK_SIZE;
+            int y0 = by * RASTERIZER_BLOCK_SIZE;
+            int x1 = x0 + RASTERIZER_BLOCK_SIZE;
+            int y1 = y0 + RASTERIZER_BLOCK_SIZE;
+
+            if (test_point(b0) && test_point(b1) && test_point(b2) && test_point(b3)) {
                 // The block is fully contained within the triangle
                 // Fill the block without further coverage tests
-                for (int y = y0; y < y1; y++) {
-                    auto coords = a;
+                auto coords = b0;
+                for (int y = y0; y < y1; y++, coords += step_y) {
                     auto* pixel = &render_target.scanline(y)[x0];
                     auto* depth = &depth_buffer.scanline(y)[x0];
-                    for (int x = x0; x < x1; x++) {
+                    for (int x = x0; x < x1; x++, coords += dbdx, pixel++, depth++) {
+                        auto barycentric = FloatVector3(coords.x(), coords.y(), coords.z()) * one_over_area;
                         if (options.enable_depth_test) {
-                            float z = interpolate(triangle.vertices[0].z, triangle.vertices[1].z, triangle.vertices[2].z, coords);
+                            float z = interpolate(triangle.vertices[0].z, triangle.vertices[1].z, triangle.vertices[2].z, barycentric);
                             if (z < *depth) {
-                                *pixel = to_rgba32(pixel_shader(coords, triangle));
+                                *pixel = to_rgba32(pixel_shader(barycentric, triangle));
                                 *depth = z;
                             }
                         } else {
-                            *pixel = to_rgba32(pixel_shader(coords, triangle));
+                            *pixel = to_rgba32(pixel_shader(barycentric, triangle));
                         }
-                        pixel++;
-                        depth++;
-                        coords = coords + dcdx;
                     }
-                    a = a + dcdy;
                 }
             } else {
                 // The block overlaps at least one triangle edge
                 // We need to test coverage of every pixel within the block
-                for (int y = y0; y < y1; y++) {
-                    auto coords = a;
+                auto coords = b0;
+                for (int y = y0; y < y1; y++, coords += step_y) {
                     auto* pixel = &render_target.scanline(y)[x0];
                     auto* depth = &depth_buffer.scanline(y)[x0];
-                    for (int x = x0; x < x1; x++) {
-                        if (test_point(coords)) {
-                            if (options.enable_depth_test) {
-                                float z = interpolate(triangle.vertices[0].z, triangle.vertices[1].z, triangle.vertices[2].z, coords);
-                                if (z < *depth) {
-                                    *pixel = to_rgba32(pixel_shader(coords, triangle));
-                                    *depth = z;
-                                }
-                            } else {
-                                *pixel = to_rgba32(pixel_shader(coords, triangle));
+                    for (int x = x0; x < x1; x++, coords += dbdx, pixel++, depth++) {
+                        
+                        if (!test_point(coords))
+                            continue;
+
+                        auto barycentric = FloatVector3(coords.x(), coords.y(), coords.z()) * one_over_area;
+                        if (options.enable_depth_test) {
+                            float z = interpolate(triangle.vertices[0].z, triangle.vertices[1].z, triangle.vertices[2].z, barycentric);
+                            if (z < *depth) {
+                                *pixel = to_rgba32(pixel_shader(barycentric, triangle));
+                                *depth = z;
                             }
+                        } else {
+                            *pixel = to_rgba32(pixel_shader(barycentric, triangle));
                         }
-                        pixel++;
-                        depth++;
-                        coords = coords + dcdx;
                     }
-                    a = a + dcdy;
                 }
             }
         }
@@ -204,7 +194,7 @@ SoftwareRasterizer::SoftwareRasterizer(const Gfx::IntSize& min_size)
 void SoftwareRasterizer::submit_triangle(const GLTriangle& triangle)
 {
     if (m_options.shade_smooth) {
-        rasterize_triangle(m_options, *m_render_target, *m_depth_buffer, triangle, [](const FloatVector4& v, const GLTriangle& t) -> FloatVector4 {
+        rasterize_triangle(m_options, *m_render_target, *m_depth_buffer, triangle, [](const FloatVector3& v, const GLTriangle& t) -> FloatVector4 {
             const float r = t.vertices[0].r * v.x() + t.vertices[1].r * v.y() + t.vertices[2].r * v.z();
             const float g = t.vertices[0].g * v.x() + t.vertices[1].g * v.y() + t.vertices[2].g * v.z();
             const float b = t.vertices[0].b * v.x() + t.vertices[1].b * v.y() + t.vertices[2].b * v.z();
@@ -212,7 +202,7 @@ void SoftwareRasterizer::submit_triangle(const GLTriangle& triangle)
             return { r, g, b, a };
         });
     } else {
-        rasterize_triangle(m_options, *m_render_target, *m_depth_buffer, triangle, [](const FloatVector4&, const GLTriangle& t) -> FloatVector4 {
+        rasterize_triangle(m_options, *m_render_target, *m_depth_buffer, triangle, [](const FloatVector3&, const GLTriangle& t) -> FloatVector4 {
             return { t.vertices[0].r, t.vertices[0].g, t.vertices[0].b, t.vertices[0].a };
         });
     }