diff --git a/crates/paint/src/render.rs b/crates/paint/src/render.rs
index 9160072a..95f90f07 100644
--- a/crates/paint/src/render.rs
+++ b/crates/paint/src/render.rs
@@ -245,6 +245,29 @@ impl BrushContext {
         // position, and thin lines beaded.
         let small = mips.is_none() && !aliased && rm < SMALL_TIP_RADIUS;
         let small_reach2 = (rm + 0.75) * (rm + 0.75);
+        // Keep ordinary computed tips separate from the feature loop, so their dense rows
+        // can be vectorised without the per-pixel sampled-mask and stroke-depth branches.
+        if mips.is_none()
+            && !small
+            && !(aliased && b.square)
+            && !noise
+            && !tex_tip
+            && !wet
+            && depth.is_none()
+            && proj.is_none()
+            && d.angle == 0.0
+            && d.radius.is_finite()
+            && d.roundness.is_finite()
+            && rm.is_finite()
+            && cx.is_finite()
+            && cy.is_finite()
+            && hardness.is_finite()
+            && d.alpha.is_finite()
+            && w > 0
+        {
+            rasterize_round(d, rect, out, hardness, aliased);
+            return;
+        }
         for yy in 0..h {
             let y = rect.y0 + yy as i32;
             let dy = y as f32 + 0.5 - cy;
@@ -889,3 +912,48 @@ pub fn render_stroke(
     r.finish();
     r.composite(&pre, target, selection, lock_transparency, true)
 }
+
+// The zero-angle transform only changes signs; squaring gives the same radial distance.
+#[inline(never)]
+fn rasterize_round(d: &Dab, rect: Rect, out: &mut [f32], hardness: f32, aliased: bool) {
+    let (cx, cy) = if aliased { grid_center(d.center.x, d.center.y, 2.0 * d.radius) } else { (d.center.x as f32, d.center.y as f32) };
+    let ro = d.roundness.max(0.5 / d.radius).min(1.0);
+    let rm = (d.radius * ro).max(0.5);
+    let w = rect.width() as usize;
+    let reach2 = (rm + 0.5) * (rm + 0.5);
+    let inner = rm * hardness.clamp(0.0, 1.0);
+    for (yy, row) in out.chunks_exact_mut(w).enumerate() {
+        let dy = (rect.y0 + yy as i32) as f32 + 0.5 - cy;
+        if aliased {
+            for (xx, dst) in row.iter_mut().enumerate() {
+                let dx = (rect.x0 + xx as i32) as f32 + 0.5 - cx;
+                let d2 = (dx * ro).powi(2) + dy * dy;
+                *dst = if d2 < reach2 && d2.sqrt() <= rm { d.alpha } else { 0.0 };
+            }
+        } else if rm - inner < 1e-3 {
+            for (xx, dst) in row.iter_mut().enumerate() {
+                let dx = (rect.x0 + xx as i32) as f32 + 0.5 - cx;
+                let d2 = (dx * ro).powi(2) + dy * dy;
+                let edge = ((rm + 0.5 - d2.sqrt()).clamp(0.0, 1.0)).min(1.0);
+                *dst = if d2 >= reach2 || edge <= 0.0 { 0.0 } else { edge * d.alpha };
+            }
+        } else {
+            for (xx, dst) in row.iter_mut().enumerate() {
+                let dx = (rect.x0 + xx as i32) as f32 + 0.5 - cx;
+                let d2 = (dx * ro).powi(2) + dy * dy;
+                let dd = d2.sqrt();
+                let edge = ((rm + 0.5 - dd).clamp(0.0, 1.0)).min(1.0);
+                let t = ((dd - inner) / (rm - inner)).clamp(0.0, 1.0);
+                let s = 1.0 - t * t;
+                let s2 = s * s;
+                let falloff = if dd <= inner { 1.0 } else { s2 * s2 };
+                let val = edge * falloff;
+                *dst = if d2 >= reach2 || val <= 0.0 { 0.0 } else { val * d.alpha };
+            }
+        }
+    }
+}
+
+#[cfg(test)]
+#[path = "render_reference.rs"]
+mod reference;
diff --git a/crates/paint/src/render_reference.rs b/crates/paint/src/render_reference.rs
new file mode 100644
index 00000000..3811e7c9
--- /dev/null
+++ b/crates/paint/src/render_reference.rs
@@ -0,0 +1,249 @@
+//! Frozen scalar rasterizer used only as a pixel oracle for the computed-tip fast path.
+use super::*;
+use crate::tile::GrayTile;
+use photocraft_geom::Point;
+
+impl BrushContext {
+    fn rasterize_reference(&self, d: &Dab, dual: bool, rect: Rect, out: &mut Vec<f32>, mut depth: Option<&mut Vec<f32>>) {
+        let (w, h) = (rect.width() as usize, rect.height() as usize);
+        out.clear();
+        out.resize(w * h, 0.0);
+        if let Some(dp) = depth.as_deref_mut() {
+            dp.clear();
+            dp.resize(w * h, 0.0);
+        }
+        let b = &self.brush;
+        let (hardness, mips) = if dual { (b.dual_brush.hardness, self.dual_tip.as_deref()) } else { (b.hardness, self.tip.as_deref()) };
+        let aliased = b.aliased && !dual;
+        let wet = b.wet_edges && !dual && depth.is_none();
+        let noise = b.noise && !dual;
+        let tex_tip = !dual && b.texture.enabled && b.texture.each_tip && self.texture.is_some();
+        let (cx, cy) = if aliased { grid_center(d.center.x, d.center.y, 2.0 * d.radius) } else { (d.center.x as f32, d.center.y as f32) };
+        let (sn, cs) = d.angle.sin_cos();
+        let (fx, fy) = (if d.flip_x { -1.0 } else { 1.0 }, if d.flip_y { -1.0 } else { 1.0 });
+        // Brush Projection: stretch the sampling coordinates along the tilt direction, which
+        // foreshortens the tip there by `proj_scale`.
+        let proj = (!dual && d.proj_scale < 0.999).then(|| {
+            let (ps, pc) = d.proj_angle.sin_cos();
+            (pc, ps, 1.0 / d.proj_scale.max(0.05) - 1.0)
+        });
+        let r = d.radius;
+        let ro = d.roundness.max(0.5 / r).min(1.0);
+        let rm = (r * ro).max(0.5);
+        let reach2 = (rm + 0.5) * (rm + 0.5);
+        // Sampled-tip mapping.
+        let (level, inv_scale, tw, th) = match mips {
+            Some(m) => {
+                let (w0, h0, _) = &m.levels[0];
+                let big = *w0.max(h0) as f32;
+                (m.level_for(2.0 * r * ro.max(0.25)), big / (2.0 * r), *w0 as f32, *h0 as f32)
+            }
+            None => (0, 1.0, 1.0, 1.0),
+        };
+        // A pixel offset from the dab centre in tip space: to y-up, project, rotate by -angle, flip.
+        let to_tip = |dx: f32, dy: f32| {
+            let (mut ux, mut uy) = (dx, -dy);
+            if let Some((ax, ay, k)) = proj {
+                let t = (ux * ax + uy * ay) * k;
+                ux += t * ax;
+                uy += t * ay;
+            }
+            ((ux * cs + uy * sn) * fx, (-ux * sn + uy * cs) * fy)
+        };
+        // A tip a few pixels across takes each pixel's covered area (a 4×4 grid inside it):
+        // sampling only the pixel centre made a tiny dab's total coverage swing with its sub-pixel
+        // position, and thin lines beaded.
+        let small = mips.is_none() && !aliased && rm < SMALL_TIP_RADIUS;
+        let small_reach2 = (rm + 0.75) * (rm + 0.75);
+        for yy in 0..h {
+            let y = rect.y0 + yy as i32;
+            let dy = y as f32 + 0.5 - cy;
+            for xx in 0..w {
+                let x = rect.x0 + xx as i32;
+                let dx = x as f32 + 0.5 - cx;
+                let (u, v) = to_tip(dx, dy);
+                let (mut val, rn) = match mips {
+                    // Block: the whole-pixel square, hard (the Eraser's Block mode).
+                    _ if aliased && b.square => {
+                        if dx.abs() >= r || dy.abs() >= r {
+                            continue;
+                        }
+                        (1.0, 0.0)
+                    }
+                    None if small => {
+                        let d2 = (u * ro).powi(2) + v * v;
+                        if d2 >= small_reach2 {
+                            continue;
+                        }
+                        let mut sum = 0.0;
+                        for sy in SUBPIXEL {
+                            for sx in SUBPIXEL {
+                                let (su, sv) = to_tip(dx + sx, dy + sy);
+                                let sd = ((su * ro).powi(2) + sv * sv).sqrt();
+                                // Each sample anti-aliases over its own quarter pixel, so the
+                                // covered area changes smoothly with the dab's position.
+                                let inside = ((rm - sd) * 4.0 + 0.5).clamp(0.0, 1.0);
+                                if inside > 0.0 {
+                                    sum += inside * tip_falloff(sd.min(rm), rm, hardness);
+                                }
+                            }
+                        }
+                        (sum / (SUBPIXEL.len() * SUBPIXEL.len()) as f32, d2.sqrt() / rm)
+                    }
+                    None => {
+                        let d2 = (u * ro).powi(2) + v * v;
+                        if d2 >= reach2 {
+                            continue;
+                        }
+                        let dd = d2.sqrt();
+                        let val = if aliased { if dd <= rm { 1.0 } else { 0.0 } } else { dab_coverage(dd, rm, hardness) };
+                        (val, dd / rm)
+                    }
+                    Some(m) => {
+                        let (px, py) = (u * inv_scale + tw / 2.0, -(v / ro) * inv_scale + th / 2.0);
+                        if px < -1.0 || py < -1.0 || px > tw + 1.0 || py > th + 1.0 {
+                            continue;
+                        }
+                        let s = m.sample_clamped(level, px / tw, py / th);
+                        let val = if aliased { if s >= 0.5 { 1.0 } else { 0.0 } } else { s };
+                        (val, 1.0 - s)
+                    }
+                };
+                if val <= 0.0 {
+                    continue;
+                }
+                if noise && val < 1.0 {
+                    let n = hash2(x, y, NOISE_SALT);
+                    let hard = if n < val { 1.0 } else { 0.0 };
+                    val += (hard - val) * 0.7;
+                }
+                if tex_tip {
+                    val = mask_combine(b.texture.mode, val, self.texture_at(x, y), d.depth);
+                }
+                if wet {
+                    val *= wet_factor(rn);
+                }
+                if let Some(slot) = depth.as_deref_mut().and_then(|dp| dp.get_mut(yy * w + xx)) {
+                    *slot = (1.0 - rn).clamp(0.0, 1.0);
+                }
+                out[yy * w + xx] = val * d.alpha;
+            }
+        }
+    }
+}
+
+fn compare(ctx: &BrushContext, dab: &Dab, dual: bool, rect: Rect, with_depth: bool) {
+    // Reused, dirty buffers also check clearing and resizing.
+    let (mut actual, mut expected) = (vec![0.75; 19], vec![0.75; 19]);
+    let (mut actual_depth, mut expected_depth) = (vec![0.25; 23], vec![0.25; 23]);
+    ctx.rasterize_into(dab, dual, rect, &mut actual, with_depth.then_some(&mut actual_depth));
+    ctx.rasterize_reference(dab, dual, rect, &mut expected, with_depth.then_some(&mut expected_depth));
+    for (label, a, b) in [("coverage", actual, expected), ("depth", actual_depth, expected_depth)] {
+        assert_eq!(a.len(), b.len());
+        let mismatch = a.iter().zip(&b).enumerate().find(|(_, (a, b))| a.to_bits() != b.to_bits());
+        assert!(mismatch.is_none(), "{label} mismatch {mismatch:?}; {dab:?}; dual={dual}, rect={rect:?}");
+    }
+}
+
+#[test]
+fn computed_tips_match_scalar_bits() {
+    for aliased in [false, true] {
+        for hardness in [-0.5, 0.0, 0.1, 0.5, 0.9999, 1.0, 1.5] {
+            let ctx = BrushContext::new(&BrushSettings { aliased, hardness, ..Default::default() });
+            for radius in [0.5, 1.0, 2.99, 3.0, 3.01, 15.0, 31.5] {
+                for roundness in [0.0, 0.1, 0.4, 1.0] {
+                    for alpha in [-0.5, -0.0, 0.0, 0.37, 1.0] {
+                        let mut d = Dab::round(Point::new(-31.25, 63.75), radius, alpha);
+                        d.roundness = roundness;
+                        d.flip_x = true;
+                        d.flip_y = true;
+                        compare(&ctx, &d, false, ctx.dab_rect(&d, false), false);
+                    }
+                }
+            }
+        }
+    }
+}
+
+#[test]
+fn large_corner_dabs_keep_off_canvas_coverage() {
+    for (hardness, aliased) in [(0.0, false), (1.0, false), (1.0, true)] {
+        let ctx = BrushContext::new(&BrushSettings { hardness, aliased, ..Default::default() });
+        let d = Dab::round(Point::new(-1050.25, -1050.75), 1500.0, 0.73);
+        let rect = ctx.dab_rect(&d, false);
+        compare(&ctx, &d, false, rect, false);
+        let mut coverage = Vec::new();
+        ctx.rasterize(&d, false, rect, &mut coverage);
+        // The layer is an infinite plane: paint outside the visible document must survive.
+        assert!(coverage.iter().enumerate().any(|(i, v)| *v > 0.0 && rect.x0 + ((i % rect.width() as usize) as i32) < 0));
+        compare(&ctx, &d, false, Rect::new(0, 0, 32, 32), false);
+        compare(&ctx, &d, false, Rect::new(0, 0, 0, 0), false);
+    }
+}
+
+#[test]
+fn feature_and_exceptional_dabs_match_scalar_bits() {
+    let tip = GrayTile::from_fn(8, 5, |x, y| ((x + y) % 7) as f32 / 6.0);
+    for feature in 0..7 {
+        let mut b = BrushSettings { hardness: 0.25, ..Default::default() };
+        match feature {
+            0 => b.noise = true,
+            1 => b.wet_edges = true,
+            2 => {
+                b.aliased = true;
+                b.square = true;
+            }
+            3 => b.tip = TipShape::Sampled(tip.clone()),
+            4 => {
+                b.texture.enabled = true;
+                b.texture.each_tip = true;
+                b.texture.pattern = Pattern::Tile(tip.clone());
+            }
+            5 => b.dual_brush.tip = TipShape::Sampled(tip.clone()),
+            _ => {}
+        }
+        let ctx = BrushContext::new(&b);
+        for dual in [false, true] {
+            for with_depth in [false, true] {
+                for (angle, proj_scale) in [(0.0, 1.0), (0.7, 1.0), (0.0, 0.6)] {
+                    let mut d = Dab::round(Point::new(63.5, -0.25), 15.0, 0.47);
+                    d.angle = angle;
+                    d.roundness = 0.4;
+                    d.proj_scale = proj_scale;
+                    compare(&ctx, &d, dual, ctx.dab_rect(&d, dual), with_depth);
+                }
+            }
+        }
+    }
+    let ctx = BrushContext::new(&BrushSettings::default());
+    for value in [f32::NAN, f32::from_bits(0x7f80_0123), f32::INFINITY, f32::NEG_INFINITY] {
+        for field in 0..4 {
+            let mut d = Dab::round(Point::new(0.0, 0.0), 15.0, 1.0);
+            match field {
+                0 => d.alpha = value,
+                1 => d.angle = value,
+                2 => d.radius = value,
+                _ => d.roundness = value,
+            }
+            compare(&ctx, &d, false, Rect::new(-4, -4, 4, 4), false);
+        }
+        let ctx = BrushContext::new(&BrushSettings { hardness: value, ..Default::default() });
+        compare(&ctx, &Dab::round(Point::new(0.0, 0.0), 15.0, 1.0), false, Rect::new(-16, -16, 16, 16), false);
+    }
+}
+
+#[test]
+fn computed_tips_match_at_pixel_and_float_boundaries() {
+    for aliased in [false, true] {
+        for center in [Point::new(0.0, -0.0), Point::new(0.5, 0.5), Point::new(255.9999, -256.0001), Point::new(16_777_216.25, -16_777_216.75)] {
+            for radius in [3.0, 3.5, 15.25, 113.75] {
+                for hardness in [0.0, 0.6, 1.0] {
+                    let ctx = BrushContext::new(&BrushSettings { hardness, aliased, ..Default::default() });
+                    let mut d = Dab::round(center, radius, 0.137);
+                    d.roundness = 0.375;
+                    compare(&ctx, &d, false, ctx.dab_rect(&d, false), false);
+                }
+            }
+        }
+    }
+}
diff --git a/crates/raster/src/lib.rs b/crates/raster/src/lib.rs
index 64c74768..3315f0bf 100644
--- a/crates/raster/src/lib.rs
+++ b/crates/raster/src/lib.rs
@@ -524,12 +524,31 @@ impl Surface {
         for tc in r.tiles() {
             let tr = tc.rect().intersect(&r);
             let t = self.tile_mut(tc);
+            // Convert contiguous channel samples a row at a time. Dispatching on depth here
+            // removes the per-pixel channel loop and lets the compiler vectorise conversion.
+            let span = tr.width() as usize * n;
             for y in tr.y0..tr.y1 {
-                for x in tr.x0..tr.x1 {
-                    let o = (((y - r.y0) as usize) * w + (x - r.x0) as usize) * n;
-                    let base = (((y - tc.ty * TILE_SIZE) as usize) * TILE_SIZE as usize + (x - tc.tx * TILE_SIZE) as usize) * n;
-                    for i in 0..n {
-                        write_sample(&mut t.data, sample, base + i, data[o + i]);
+                let o = (((y - r.y0) as usize) * w + (tr.x0 - r.x0) as usize) * n;
+                let base = (((y - tc.ty * TILE_SIZE) as usize) * TILE_SIZE as usize + (tr.x0 - tc.tx * TILE_SIZE) as usize) * n;
+                let src = &data[o..o + span];
+                let bytes = sample.bytes();
+                let dst = &mut t.data[base * bytes..(base + span) * bytes];
+                match sample {
+                    SampleType::U8 => {
+                        for (d, v) in dst.iter_mut().zip(src) {
+                            *d = (v.clamp(0.0, 1.0) * 255.0 + 0.5) as u8;
+                        }
+                    }
+                    SampleType::U16 => {
+                        for (d, v) in dst.as_chunks_mut::<2>().0.iter_mut().zip(src) {
+                            let q = (v.clamp(0.0, 1.0) * 65535.0 + 0.5) as u16;
+                            d.copy_from_slice(&q.to_ne_bytes());
+                        }
+                    }
+                    SampleType::F32 => {
+                        for (d, v) in dst.as_chunks_mut::<4>().0.iter_mut().zip(src) {
+                            d.copy_from_slice(&v.to_ne_bytes());
+                        }
                     }
                 }
             }
@@ -871,6 +890,89 @@ mod tests {
     use photocraft_color::PixelFormat;
     use proptest::prelude::*;
 
+    fn write_region_reference(surface: &mut Surface, r: Rect, data: &[f32]) {
+        let n = surface.channels();
+        let w = r.width() as usize;
+        assert_eq!(data.len(), w * r.height() as usize * n, "region data length mismatch");
+        let sample = surface.format.sample;
+        for tc in r.tiles() {
+            let tr = tc.rect().intersect(&r);
+            let t = surface.tile_mut(tc);
+            for y in tr.y0..tr.y1 {
+                for x in tr.x0..tr.x1 {
+                    let o = (((y - r.y0) as usize) * w + (x - r.x0) as usize) * n;
+                    let base = (((y - tc.ty * TILE_SIZE) as usize) * TILE_SIZE as usize + (x - tc.tx * TILE_SIZE) as usize) * n;
+                    for i in 0..n {
+                        write_sample(&mut t.data, sample, base + i, data[o + i]);
+                    }
+                }
+            }
+        }
+    }
+
+    #[test]
+    fn region_writes_match_scalar_bytes_for_every_format() {
+        let special = [
+            0.0,
+            -0.0,
+            1.0,
+            -1.0,
+            2.0,
+            f32::INFINITY,
+            f32::NEG_INFINITY,
+            f32::NAN,
+            f32::from_bits(0x7fc12345),
+            f32::from_bits(0xffc12345),
+            f32::MIN_POSITIVE,
+            f32::from_bits(1),
+            0.5 / 255.0,
+            0.5 / 65535.0,
+        ];
+        for mode in [
+            ColorMode::Bitmap,
+            ColorMode::Grayscale,
+            ColorMode::Indexed,
+            ColorMode::Rgb,
+            ColorMode::Cmyk,
+            ColorMode::Lab,
+            ColorMode::Multichannel,
+            ColorMode::Duotone,
+        ] {
+            for sample in SampleType::ALL {
+                for alpha in [false, true] {
+                    let fmt = PixelFormat::new(mode, sample, alpha);
+                    let mut actual = Surface::with_default(fmt, &vec![0.25; fmt.channels()]);
+                    actual.write_pixel(900, 900, &vec![0.75; fmt.channels()]);
+                    let mut expected = actual.clone();
+                    for rect in
+                        [Rect::EMPTY, Rect::new(-270, -2, 270, 3), Rect::new(256, 256, 512, 257), Rect::new(-64, -64, 0, 0), Rect::new(255, 250, 258, 270)]
+                    {
+                        let data: Vec<_> = (0..rect.width() as usize * rect.height() as usize * fmt.channels())
+                            .map(|i| {
+                                if i % 3 == 0 {
+                                    special[i / 3 % special.len()]
+                                } else {
+                                    // Integer rounding boundaries and deterministic arbitrary f32 bit patterns.
+                                    let bits = (i as u32).wrapping_mul(2654435761);
+                                    if i % 3 == 1 { ((bits % 65536) as f32 + 0.5) / 65535.0 } else { f32::from_bits(bits) }
+                                }
+                            })
+                            .collect();
+                        let snapshot = actual.clone();
+                        let snapshot_bytes = snapshot.to_interleaved(rect.inflate(1));
+                        actual.write_region(rect, &data);
+                        write_region_reference(&mut expected, rect, &data);
+                        assert_eq!(actual, expected, "{fmt:?}, {rect:?}");
+                        assert_eq!(snapshot.to_interleaved(rect.inflate(1)), snapshot_bytes, "COW snapshot changed");
+                        // Repeat against already materialised, unshared tiles too.
+                        actual.write_region(rect, &data);
+                        assert_eq!(actual, expected);
+                    }
+                }
+            }
+        }
+    }
+
     #[test]
     fn from_rgba_into_matches_from_rgba() {
         for fmt in [PixelFormat::RGBA8, PixelFormat::GRAY8, PixelFormat::GRAYA8, PixelFormat::CMYKA8, PixelFormat::new(ColorMode::Lab, SampleType::F32, true)] {