From 6f2a83e99057227f0e4409af238fdec3b3ab6add Mon Sep 17 00:00:00 2001
From: Raph Levien <raph@google.com>
Date: Mon, 30 Oct 2023 11:58:01 -0700
Subject: [PATCH] Move numerical robustness into tiling

Change representation to tile-relative coordinates for segment endpoints.
---
 crates/encoding/src/path.rs   |  4 +-
 shader/fine.wgsl              | 84 +++++++++++++----------------------
 shader/path_tiling.wgsl       | 41 +++++++++++++----
 shader/shared/segment.wgsl    |  5 ++-
 src/cpu_shader/fine.rs        | 20 +++++----
 src/cpu_shader/path_tiling.rs | 54 +++++++++++++++-------
 6 files changed, 118 insertions(+), 90 deletions(-)

diff --git a/crates/encoding/src/path.rs b/crates/encoding/src/path.rs
index 83e537042..c1c521d30 100644
--- a/crates/encoding/src/path.rs
+++ b/crates/encoding/src/path.rs
@@ -207,8 +207,8 @@ pub struct SegmentCount {
 #[derive(Clone, Copy, Debug, Zeroable, Pod, Default)]
 #[repr(C)]
 pub struct PathSegment {
-    pub origin: [f32; 2],
-    pub delta: [f32; 2],
+    pub point0: [f32; 2],
+    pub point1: [f32; 2],
     pub y_edge: f32,
     pub _padding: u32,
 }
diff --git a/shader/fine.wgsl b/shader/fine.wgsl
index cc5c1c8e3..016095327 100644
--- a/shader/fine.wgsl
+++ b/shader/fine.wgsl
@@ -162,29 +162,18 @@ fn fill_path_ms(fill: CmdFill, local_id: vec2<u32>, result: ptr<function, array<
         // TODO: might save a register rewriting this in terms of limit
         if th_ix < slice_size {
             let segment = segments[seg_off];
-            let xy0 = segment.origin;
-            let xy1 = xy0 + segment.delta;
+            let xy0 = segment.point0;
+            let xy1 = segment.point1;
             var y_edge_f = f32(TILE_HEIGHT);
             var delta = select(-1, 1, xy1.x <= xy0.x);
-            if xy0.x == 0.0 && xy1.x == 0.0 {
-                if xy0.y == 0.0 {
-                    y_edge_f = 0.0;
-                } else if xy1.y == 0.0 {
-                    y_edge_f = 0.0;
-                    delta = -delta;
-                }
-            } else {
-                if xy0.x == 0.0 {
-                    if xy0.y != 0.0 {
-                        y_edge_f = xy0.y;
-                    }
-                } else if xy1.x == 0.0 && xy1.y != 0.0 {
-                    y_edge_f = xy1.y;
-                }
-                // discard horizontal lines aligned to pixel grid
-                if !(xy0.y == xy1.y && xy0.y == floor(xy0.y)) {
-                    count = span(xy0.x, xy1.x) + span(xy0.y, xy1.y) - 1u;
-                }
+            if xy0.x == 0.0 {
+                y_edge_f = xy0.y;
+            } else if xy1.x == 0.0 {
+                y_edge_f = xy1.y;
+            }
+            // discard horizontal lines aligned to pixel grid
+            if !(xy0.y == xy1.y && xy0.y == floor(xy0.y)) {
+                count = span(xy0.x, xy1.x) + span(xy0.y, xy1.y) - 1u;
             }
             let y_edge = u32(ceil(y_edge_f));
             if y_edge < TILE_HEIGHT {
@@ -222,8 +211,8 @@ fn fill_path_ms(fill: CmdFill, local_id: vec2<u32>, result: ptr<function, array<
             let seg_off = fill.seg_data + batch * WG_SIZE + el_ix;
             let segment = segments[seg_off];
             // Coordinates are relative to tile origin
-            let xy0_in = segment.origin;
-            let xy1_in = xy0_in + segment.delta;
+            let xy0_in = segment.point0;
+            let xy1_in = segment.point1;
             let is_down = xy1_in.y >= xy0_in.y;
             let xy0 = select(xy1_in, xy0_in, is_down);
             let xy1 = select(xy0_in, xy1_in, is_down);
@@ -514,27 +503,17 @@ fn fill_path_ms_evenodd(fill: CmdFill, local_id: vec2<u32>, result: ptr<function
         if th_ix < slice_size {
             let segment = segments[seg_off];
             // Coordinates are relative to tile origin
-            let xy0 = segment.origin;
-            let xy1 = xy0 + segment.delta;
+            let xy0 = segment.point0;
+            let xy1 = segment.point1;
             var y_edge_f = f32(TILE_HEIGHT);
-            if xy0.x == 0.0 && xy1.x == 0.0 {
-                if xy0.y == 0.0 {
-                    y_edge_f = 0.0;
-                } else if xy1.y == 0.0 {
-                    y_edge_f = 0.0;
-                }
-            } else {
-                if xy0.x == 0.0 {
-                    if xy0.y != 0.0 {
-                        y_edge_f = xy0.y;
-                    }
-                } else if xy1.x == 0.0 && xy1.y != 0.0 {
-                    y_edge_f = xy1.y;
-                }
-                // discard horizontal lines aligned to pixel grid
-                if !(xy0.y == xy1.y && xy0.y == floor(xy0.y)) {
-                    count = span(xy0.x, xy1.x) + span(xy0.y, xy1.y) - 1u;
-                }
+            if xy0.x == 0.0 {
+                y_edge_f = xy0.y;
+            } else if xy1.x == 0.0 {
+                y_edge_f = xy1.y;
+            }
+            // discard horizontal lines aligned to pixel grid
+            if !(xy0.y == xy1.y && xy0.y == floor(xy0.y)) {
+                count = span(xy0.x, xy1.x) + span(xy0.y, xy1.y) - 1u;
             }
             let y_edge = u32(ceil(y_edge_f));
             if y_edge < TILE_HEIGHT {
@@ -571,8 +550,8 @@ fn fill_path_ms_evenodd(fill: CmdFill, local_id: vec2<u32>, result: ptr<function
             let sub_ix = i - select(0u, sh_count[el_ix - 1u], el_ix > 0u);
             let seg_off = fill.seg_data + batch * WG_SIZE + el_ix;
             let segment = segments[seg_off];
-            let xy0_in = segment.origin;
-            let xy1_in = xy0_in + segment.delta;
+            let xy0_in = segment.point0;
+            let xy1_in = segment.point1;
             let is_down = xy1_in.y >= xy0_in.y;
             let xy0 = select(xy1_in, xy0_in, is_down);
             let xy1 = select(xy0_in, xy1_in, is_down);
@@ -807,17 +786,18 @@ fn fill_path(fill: CmdFill, xy: vec2<f32>, result: ptr<function, array<f32, PIXE
     for (var i = 0u; i < n_segs; i++) {
         let seg_off = fill.seg_data + i;
         let segment = segments[seg_off];
-        let y = segment.origin.y - xy.y;
+        let y = segment.point0.y - xy.y;
+        let delta = segment.point1 - segment.point0;
         let y0 = clamp(y, 0.0, 1.0);
-        let y1 = clamp(y + segment.delta.y, 0.0, 1.0);
+        let y1 = clamp(y + delta.y, 0.0, 1.0);
         let dy = y0 - y1;
         if dy != 0.0 {
-            let vec_y_recip = 1.0 / segment.delta.y;
+            let vec_y_recip = 1.0 / delta.y;
             let t0 = (y0 - y) * vec_y_recip;
             let t1 = (y1 - y) * vec_y_recip;
-            let startx = segment.origin.x - xy.x;
-            let x0 = startx + t0 * segment.delta.x;
-            let x1 = startx + t1 * segment.delta.x;
+            let startx = segment.point0.x - xy.x;
+            let x0 = startx + t0 * delta.x;
+            let x1 = startx + t1 * delta.x;
             let xmin0 = min(x0, x1);
             let xmax0 = max(x0, x1);
             for (var i = 0u; i < PIXELS_PER_THREAD; i += 1u) {
@@ -831,7 +811,7 @@ fn fill_path(fill: CmdFill, xy: vec2<f32>, result: ptr<function, array<f32, PIXE
                 area[i] += a * dy;
             }
         }
-        let y_edge = sign(segment.delta.x) * clamp(xy.y - segment.y_edge + 1.0, 0.0, 1.0);
+        let y_edge = sign(delta.x) * clamp(xy.y - segment.y_edge + 1.0, 0.0, 1.0);
         for (var i = 0u; i < PIXELS_PER_THREAD; i += 1u) {
             area[i] += y_edge;
         }
diff --git a/shader/path_tiling.wgsl b/shader/path_tiling.wgsl
index 04a640d02..8cd260ef4 100644
--- a/shader/path_tiling.wgsl
+++ b/shader/path_tiling.wgsl
@@ -120,19 +120,42 @@ fn main(
                 xy1 = vec2(x_clip, yt);
             }
         }
-        // See comments in CPU version of shader
         var y_edge = 1e9;
-        if xy0.x == tile_xy.x && xy1.x != tile_xy.x && xy0.y != tile_xy.y {
-            y_edge = xy0.y;
-        } else if xy1.x == tile_xy.x && xy1.y != tile_xy.y {
-            y_edge = xy1.y;
+        // Apply numerical robustness logic
+        var p0 = xy0 - tile_xy;
+        var p1 = xy1 - tile_xy;
+        // When we move to f16, this will be f16::MIN_POSITIVE
+        let EPSILON = 1e-6;
+        if p0.x == 0.0 {
+            if p1.x == 0.0 {
+                p0.x = EPSILON;
+                if p0.y == 0.0 {
+                    // Entire tile
+                    p1.x = EPSILON;
+                    p1.y = f32(TILE_HEIGHT);
+                } else {
+                    // Make segment disappear
+                    p1.x = 2.0 * EPSILON;
+                    p1.y = p0.y;
+                }
+            } else if p0.y == 0.0 {
+                p0.x = EPSILON;
+            } else {
+                y_edge = p0.y;
+            }
+        } else if p1.x == 0.0 {
+            if p1.y == 0.0 {
+                p1.x = EPSILON;
+            } else {
+                y_edge = p1.y;
+            }
         }
         if !is_down {
-            let tmp = xy0;
-            xy0 = xy1;
-            xy1 = tmp;
+            let tmp = p0;
+            p0 = p1;
+            p1 = tmp;
         }
-        let segment = Segment(xy0 - tile_xy, xy1 - xy0, y_edge - tile_xy.y);
+        let segment = Segment(p0, p1, y_edge);
         segments[seg_start + seg_within_slice] = segment;
     }
 }
diff --git a/shader/shared/segment.wgsl b/shader/shared/segment.wgsl
index 9fe134478..03a68e4c9 100644
--- a/shader/shared/segment.wgsl
+++ b/shader/shared/segment.wgsl
@@ -2,8 +2,9 @@
 
 // Segments laid out for contiguous storage
 struct Segment {
-    origin: vec2<f32>,
-    delta: vec2<f32>,
+    // Points are relative to tile origin
+    point0: vec2<f32>,
+    point1: vec2<f32>,
     y_edge: f32,
 }
 
diff --git a/src/cpu_shader/fine.rs b/src/cpu_shader/fine.rs
index c64c87627..dc3203959 100644
--- a/src/cpu_shader/fine.rs
+++ b/src/cpu_shader/fine.rs
@@ -58,20 +58,24 @@ fn fill_path(area: &mut [f32], segments: &[PathSegment], fill: &CmdFill, x_tile:
         *a = backdrop_f;
     }
     for segment in &segments[fill.seg_data as usize..][..n_segs as usize] {
+        let delta = [
+            segment.point1[0] - segment.point0[0],
+            segment.point1[1] - segment.point0[1],
+        ];
         for yi in 0..TILE_HEIGHT {
-            let y = segment.origin[1] - (y_tile + yi as f32);
+            let y = segment.point0[1] - (y_tile + yi as f32);
             let y0 = y.clamp(0.0, 1.0);
-            let y1 = (y + segment.delta[1]).clamp(0.0, 1.0);
+            let y1 = (y + delta[1]).clamp(0.0, 1.0);
             let dy = y0 - y1;
-            let y_edge = segment.delta[0].signum()
-                * (y_tile + yi as f32 - segment.y_edge + 1.0).clamp(0.0, 1.0);
+            let y_edge =
+                delta[0].signum() * (y_tile + yi as f32 - segment.y_edge + 1.0).clamp(0.0, 1.0);
             if dy != 0.0 {
-                let vec_y_recip = segment.delta[1].recip();
+                let vec_y_recip = delta[1].recip();
                 let t0 = (y0 - y) * vec_y_recip;
                 let t1 = (y1 - y) * vec_y_recip;
-                let startx = segment.origin[0] - x_tile;
-                let x0 = startx + t0 * segment.delta[0];
-                let x1 = startx + t1 * segment.delta[0];
+                let startx = segment.point0[0] - x_tile;
+                let x0 = startx + t0 * delta[0];
+                let x1 = startx + t1 * delta[0];
                 let xmin0 = x0.min(x1);
                 let xmax0 = x0.max(x1);
                 for i in 0..TILE_WIDTH {
diff --git a/src/cpu_shader/path_tiling.rs b/src/cpu_shader/path_tiling.rs
index bf2975233..849474cf9 100644
--- a/src/cpu_shader/path_tiling.rs
+++ b/src/cpu_shader/path_tiling.rs
@@ -115,28 +115,48 @@ fn path_tiling_main(
                 xy1 = Vec2::new(x_clip, yt);
             }
         }
+        let mut y_edge = 1e9;
+        // Apply numerical robustness logic
+        let mut p0 = xy0 - tile_xy;
+        let mut p1 = xy1 - tile_xy;
+        const EPSILON: f32 = 1e-6;
+        if p0.x == 0.0 {
+            if p1.x == 0.0 {
+                p0.x = EPSILON;
+                if p0.y == 0.0 {
+                    // Entire tile
+                    p1.x = EPSILON;
+                    p1.y = TILE_HEIGHT as f32;
+                } else {
+                    // Make segment disappear
+                    p1.x = 2.0 * EPSILON;
+                    p1.y = p0.y;
+                }
+            } else if p0.y == 0.0 {
+                p0.x = EPSILON;
+            } else {
+                y_edge = p0.y;
+            }
+        } else if p1.x == 0.0 {
+            if p1.y == 0.0 {
+                p1.x = EPSILON;
+            } else {
+                y_edge = p1.y;
+            }
+        }
         if !is_down {
-            (xy0, xy1) = (xy1, xy0);
+            (p0, p1) = (p1, p0);
         }
-        // TODO (part of move to 8 byte encoding for segments): don't store y_edge at all,
-        // resolve this in fine.
-        let y_edge = if xy0.x == tile_xy.x && xy1.x != tile_xy.x && xy0.y != tile_xy.y {
-            xy0.y
-        } else if xy1.x == tile_xy.x && xy1.y != tile_xy.y {
-            xy1.y
-        } else {
-            1e9
-        };
         let segment = PathSegment {
-            origin: (xy0 - tile_xy).to_array(),
-            delta: (xy1 - xy0).to_array(),
-            y_edge: y_edge - tile_xy.y,
+            point0: p0.to_array(),
+            point1: p1.to_array(),
+            y_edge,
             _padding: Default::default(),
         };
-        assert!(xy0.x >= tile_xy.x && xy0.x <= tile_xy1.x);
-        assert!(xy0.y >= tile_xy.y && xy0.y <= tile_xy1.y);
-        assert!(xy1.x >= tile_xy.x && xy1.x <= tile_xy1.x);
-        assert!(xy1.y >= tile_xy.y && xy1.y <= tile_xy1.y);
+        assert!(p0.x >= 0.0 && p0.x <= TILE_WIDTH as f32);
+        assert!(p0.y >= 0.0 && p0.y <= TILE_HEIGHT as f32);
+        assert!(p1.x >= 0.0 && p1.x <= TILE_WIDTH as f32);
+        assert!(p1.y >= 0.0 && p1.y <= TILE_HEIGHT as f32);
         segments[(seg_start + seg_within_slice) as usize] = segment;
     }
 }