The Graphics Pipeline & Rasterization

A 3D scene is described in multiple coordinate systems, and moving geometry between them is the first job of the pipeline.

Model space (also called object space): each mesh is defined relative to its own origin. A character's arm is at (0,0,0) in arm-space regardless of where the character stands.

World space: all objects are placed in a common scene coordinate system. The model matrix M (a 4×4 matrix) transforms from model to world space — it encodes the object's position, orientation, and scale.

View space (camera space): the camera is moved to the origin, looking down the −Z axis (by convention). The view matrix V transforms from world to camera space. Constructing V means inverting the camera's own model matrix, or equivalently, using the classic lookAt construction:

• forward = normalize(eye − target)
• right = normalize(up × forward)
• up' = forward × right

Then V is the 4×4 matrix that rotates and translates world coordinates into this camera frame.

Clip space: the projection matrix P transforms from camera space into a canonical volume (the view frustum gets remapped). After applying MVP = P·V·M, the resulting 4-vector (x, y, z, w) is in homogeneous clip coordinates. A vertex is inside the frustum when −w ≤ x, y, z ≤ w.

NDC (Normalized Device Coordinates): divide by w (the perspective divide) to get (x/w, y/w, z/w) ∈ [−1,1]³. This is the canonical cube.

Screen space: NDC coordinates are mapped to pixel coordinates by the viewport transform.

The key insight: all of these are just linear (or projective) transformations. Every vertex is multiplied by the same MVP matrix, making this trivially parallelizable — the GPU can transform millions of vertices simultaneously.

// GLSL vertex shader — runs once per vertex on the GPU
layout(location = 0) in vec3 aPosition;  // model space

uniform mat4 uModel;   // model → world
uniform mat4 uView;    // world → camera
uniform mat4 uProj;    // camera → clip

void main() {
    // gl_Position is the output clip-space position
    gl_Position = uProj * uView * uModel * vec4(aPosition, 1.0);
}

To represent perspective projection as a matrix multiplication (which GPUs can do efficiently), we use homogeneous coordinates: a point (x, y, z) in 3D is represented as the 4-vector (x, y, z, 1). Any scalar multiple (kx, ky, kz, k) for k≠0 represents the same point — you recover the 3D point by dividing by the last component.

This trick lets us encode the perspective divide (division by depth) inside a matrix. The perspective projection matrix for a frustum with field-of-view angle fovy, aspect ratio a, near plane n, and far plane f is:

P = | cot(fovy/2)/a   0          0              0       |
    | 0              cot(fovy/2)  0              0       |
    | 0               0       -(f+n)/(f-n)  -2fn/(f-n) |
    | 0               0          -1             0       |

After multiplying: the w component becomes −z_camera (the depth). After the perspective divide (x/w, y/w, z/w), points at depth z=n map to z_ndc=−1 and points at z=f map to z_ndc=+1 (OpenGL convention). This remapping is nonlinear in z — depth is compressed near the far plane and expanded near the near plane, which has crucial consequences for depth buffer precision.

Why do we do this? Because parallel lines in 3D should converge to a vanishing point in 2D (perspective foreshortening). A matrix multiplication followed by a divide achieves exactly this. An orthographic projection, by contrast, uses w=1 always (no perspective foreshortening) — useful for engineering drawings and 2D UI rendering.

The w-divide also enables homogeneous clipping: instead of clipping in 3D (complex, inconsistent), clip against the six planes −w ≤ x, y, z ≤ w in clip space (before the divide). This is why clipping happens before NDC.

import numpy as np

def perspective_matrix(fovy_deg, aspect, near, far):
    """Build an OpenGL-style perspective matrix."""
    f = 1.0 / np.tan(np.radians(fovy_deg) / 2.0)
    n, fa = near, far
    return np.array([
        [f/aspect, 0,              0,               0],
        [0,        f,              0,               0],
        [0,        0,  -(fa+n)/(fa-n),  -2*fa*n/(fa-n)],
        [0,        0,             -1,               0]
    ], dtype=np.float64)

# Transform a point at camera-space position (0, 0, -5)
p_cam = np.array([0, 0, -5, 1.0])
P = perspective_matrix(60, 16/9, 0.1, 100)
p_clip = P @ p_cam
print('clip:', p_clip)            # w = 5.0
print('ndc:', p_clip / p_clip[3]) # z_ndc between -1 and 1

Before we can rasterize a triangle, we must ensure it lies within the view volume. Triangles that straddle the frustum boundary must be clipped: cut so that only the visible portion remains.

Why clip? If a vertex has w=0 or w<0 (behind the camera), the perspective divide is undefined or flips the geometry inside-out. We must remove those parts before dividing.

The Sutherland-Hodgman algorithm clips a polygon against a half-space one at a time. For each clipping plane:

• If both endpoints of an edge are inside → keep the second vertex.
• If going outside → compute and output the intersection point.
• If going inside → compute the intersection point and output both it and the second vertex.
• If both outside → output nothing.

Applying this for all six frustum planes (left, right, bottom, top, near, far) handles the complete clip.

The near plane is the most critical. If a triangle vertex is behind the camera (z > 0 in camera space, i.e., w < 0 after projection), rasterizing it without clipping would project it to the wrong side of the screen. The near clip plane prevents this.

After clipping, clipped triangles may have 4, 5, or 6 vertices (a clipped triangle can become a pentagon if two planes clip it). These are trivially re-triangulated (fan from the first vertex) before rasterization.

Guard-band clipping: modern GPUs use a wider clip rectangle (2× or 4× screen size) for X/Y clipping, deferring the scissor test to the rasterizer. Only the near/far and trivial reject cases require exact Sutherland-Hodgman. This reduces clipping overhead significantly for large triangles that just barely extend off-screen.

# Sutherland-Hodgman: clip polygon against one half-space
# inside(v) returns True if v is on the inside of the plane
# intersect(a, b) returns the clip-plane intersection point

def clip_polygon_one_plane(polygon, inside, intersect):
    output = []
    if not polygon:
        return output
    s = polygon[-1]  # start with last vertex
    for e in polygon:
        if inside(e):
            if not inside(s):
                output.append(intersect(s, e))  # edge crosses in
            output.append(e)
        elif inside(s):
            output.append(intersect(s, e))       # edge crosses out
        s = e
    return output

After the perspective divide we have NDC coordinates: (x, y, z) ∈ [−1, 1]³. The viewport transform maps these to pixel coordinates:

x_screen = (x_ndc + 1) / 2 * width  + x_offset
y_screen = (1 - y_ndc) / 2 * height + y_offset   // Y flipped: NDC +Y is up, screen +Y is down
z_screen = (z_ndc + 1) / 2           // depth mapped to [0, 1]

Note the Y flip: in NDC, +Y is up (mathematical convention). On screen, +Y is down (raster convention). The depth value is remapped to [0, 1] for the depth buffer.

After the viewport transform, the three vertices of a triangle are now in screen space: floating-point (x, y) pixel positions and a depth z ∈ [0, 1]. The rasterizer's job is to determine which pixels lie inside this triangle.

Back-face culling happens before rasterization: compute the signed area of the screen-space triangle using the cross product of two edges. If the area is negative (or positive, depending on winding convention), the triangle faces away from the camera and is discarded. This eliminates ~half of all triangles with zero cost.

Signed area = ((v1 - v0) × (v2 - v0)).z / 2

If your convention is counter-clockwise = front-facing, discard when signed area ≤ 0.

def signed_area_2d(v0, v1, v2):
    """Signed area of a 2D triangle (positive = CCW)."""
    return 0.5 * ((v1[0]-v0[0])*(v2[1]-v0[1])
                 -(v2[0]-v0[0])*(v1[1]-v0[1]))

def viewport_transform(x_ndc, y_ndc, z_ndc, width, height):
    xs = (x_ndc + 1) / 2 * width
    ys = (1 - y_ndc) / 2 * height  # flip Y
    zs = (z_ndc + 1) / 2           # [0, 1] depth
    return xs, ys, zs

The core rasterization question is: for each pixel, is it inside the triangle? And if so, what are the interpolated attributes (color, texture coordinate, normal) at that pixel?

Barycentric coordinates are the key tool. Any point P inside triangle (V0, V1, V2) can be written as:

P = λ₀V₀ + λ₁V₁ + λ₂V₂, with λ₀ + λ₁ + λ₂ = 1

The λ values are the barycentric coordinates. Their geometric meaning: λ₀ = (area of triangle PV1V2) / (area of triangle V0V1V2). A point is inside the triangle iff all three λᵢ ≥ 0.

Computing λ from 2D screen-space positions uses ratios of signed triangle areas (or equivalently, edge functions).

The edge function test: define edge functions e₀₁(P) = (P − V0) × (V1 − V0) for each edge. A point is inside iff all three edge functions have the same sign as the triangle's signed area. This is equivalent to checking barycentric coordinates and is what GPU rasterizers compute.

The incremental algorithm: for a bounding box scan, edge functions change by a constant amount per pixel step (Δx adds a constant, Δy adds a constant). This means we can compute e(P) for the first pixel and then use additions only for subsequent pixels — no multiplications per pixel.

Attribute interpolation: once we have (λ₀, λ₁, λ₂), any per-vertex attribute A interpolates as A(P) = λ₀A₀ + λ₁A₁ + λ₂A₂. This gives smooth shading, texture coordinates, and normals across the triangle surface.

Filling rules: to avoid double-counting pixels on shared edges between adjacent triangles, GPUs use the top-left rule: a pixel center is included if the edge is a top edge (horizontal, going right) or a left edge (going up). This ensures each shared edge is owned by exactly one triangle.

import numpy as np

def rasterize_triangle(v0, v1, v2, attr0, attr1, attr2, framebuffer, zbuffer):
    """
    v0,v1,v2: (x, y, z) in screen space
    attr0..2: per-vertex color (r,g,b)
    """
    h, w = framebuffer.shape[:2]
    # Bounding box
    x0 = max(0, int(min(v0[0], v1[0], v2[0])))
    x1 = min(w-1, int(max(v0[0], v1[0], v2[0])))
    y0 = max(0, int(min(v0[1], v1[1], v2[1])))
    y1 = min(h-1, int(max(v0[1], v1[1], v2[1])))

    total_area = signed_area_2d(v0, v1, v2)
    if abs(total_area) < 1e-8:
        return  # degenerate

    for y in range(y0, y1+1):
        for x in range(x0, x1+1):
            p = (x + 0.5, y + 0.5)  # pixel center
            lam0 = signed_area_2d(p, v1, v2) / total_area
            lam1 = signed_area_2d(v0, p,  v2) / total_area
            lam2 = 1 - lam0 - lam1
            if lam0 >= 0 and lam1 >= 0 and lam2 >= 0:
                z = lam0*v0[2] + lam1*v1[2] + lam2*v2[2]
                if z < zbuffer[y, x]:          # depth test
                    zbuffer[y, x] = z
                    color = lam0*attr0 + lam1*attr1 + lam2*attr2
                    framebuffer[y, x] = np.clip(color, 0, 255)

There is a subtle but critical bug in naive barycentric interpolation: the screen-space barycentric coordinates λᵢ are computed after the perspective divide, but the attributes (texture coordinates, colors, normals) live in 3D world space. Linear interpolation in screen space does not correspond to linear interpolation in 3D — it produces texture swimming and incorrect shading.

The problem: after perspective projection, equal steps in screen space correspond to larger steps in 3D near the camera and smaller steps far away. So linearly interpolating a texture coordinate in screen space produces distorted results (think of a tiled floor in perspective — tiles near the camera look stretched if you interpolate naively).

The fix: for perspective-correct interpolation, you must interpolate attributes divided by the clip-space w (before the perspective divide), then divide by the interpolated 1/w:

If the per-vertex value is tᵢ and the per-vertex clip w is wᵢ, then at a screen-space point with barycentric coordinates (λ₀, λ₁, λ₂):

t_correct = (λ₀·t₀/w₀ + λ₁·t₁/w₁ + λ₂·t₂/w₂)
            / (λ₀/w₀   + λ₁/w₁   + λ₂/w₂)

This is perspective-correct interpolation. The denominator gives a perspective-correct 1/w at the pixel, and the numerator gives the attribute weighted by 1/w.

Modern GPUs implement this automatically when you declare in variables in a fragment shader — they do the /w division before interpolating and multiply back before delivering the value to your fragment shader. The keyword noperspective disables this and gives raw screen-space linear interpolation (useful for 2D UI effects).

The mathematical derivation uses the fact that 1/z is linear in screen space (a classic result from projective geometry), while z itself is not.

# Perspective-correct interpolation of texture coord (u, v)
# v0, v1, v2 are (x_screen, y_screen, z_screen, w_clip)
# uv0, uv1, uv2 are (u, v) at each vertex

def persp_correct_interp(lam, uvs, ws):
    """
    lam: (lam0, lam1, lam2) screen-space barycentric coords
    uvs: [(u0,v0), (u1,v1), (u2,v2)]
    ws:  [w0, w1, w2]  clip-space w values
    """
    inv_ws = [1.0/w for w in ws]
    denom = sum(lam[i] * inv_ws[i] for i in range(3))
    u = sum(lam[i] * uvs[i][0] * inv_ws[i] for i in range(3)) / denom
    v = sum(lam[i] * uvs[i][1] * inv_ws[i] for i in range(3)) / denom
    return u, v

When multiple triangles overlap in screen space, we need to show only the closest one — this is the visibility problem. The naïve solution (painter's algorithm: draw back-to-front) fails for cyclic overlaps and is expensive to sort.

The z-buffer algorithm (Catmull, 1974) solves this elegantly: maintain a 2D array of depth values (the z-buffer or depth buffer), initialized to the maximum depth (1.0 or +∞). When rasterizing a fragment:

1. Compute the interpolated depth z ∈ [0, 1] at this pixel.
2. If z < zbuffer[y][x]: this fragment is closer than whatever was drawn before. Update zbuffer[y][x] = z and write the color to the framebuffer.
3. Otherwise: discard this fragment — something closer was already drawn.

The algorithm is O(n) in the number of fragments, not O(n log n) like sorting. It handles all cases correctly — even cyclic overlaps — because it operates per-pixel, not per-object.

Depth precision issues: the perspective matrix maps depth non-linearly to [0, 1]. Most precision is concentrated near the near plane, and very little near the far plane. With a 24-bit depth buffer and a near/far ratio of 1000:1, you get significant z-fighting (flickering between two close surfaces far from the camera). Solutions:

• Use a 32-bit float depth buffer.
• Keep near/far ratio below 1000:1.
• Use reversed-Z: map near to 1 and far to 0. Combined with floating-point depth buffers, this distributes precision more uniformly because floats have more precision near 0.

Early-Z / depth prepass: modern GPUs execute the depth test before the fragment shader when possible (hardware early-Z). If the fragment will be discarded by the depth test, we skip the (expensive) fragment shader. A depth prepass renders all geometry with no fragment shader first to fill the z-buffer, then renders again with shading — this converts most fragment shader invocations to early-Z discards.

# Minimal software z-buffer rasterizer core
import numpy as np

W, H = 800, 600
framebuf = np.zeros((H, W, 3), dtype=np.uint8)
zbuf     = np.full((H, W), np.inf)  # +inf = nothing drawn yet

def draw_fragment(x, y, z, color):
    if 0 <= x < W and 0 <= y < H:
        if z < zbuf[y, x]:          # depth test: closer wins
            zbuf[y, x]       = z
            framebuf[y, x]   = color

After rasterization produces a fragment (a candidate pixel with interpolated attributes), the fragment shader (or pixel shader in Direct3D terminology) runs to compute the final color. This is where lighting, texturing, shadows, and all visual effects are computed.

A fragment shader receives:

• Interpolated data from the vertex shader: texture coordinates, normals, tangents, view direction
• Uniform data: textures (samplers), material properties, light positions
• Built-in values: gl_FragCoord (screen position + depth), gl_FrontFacing

And outputs:

• One or more colors written to render targets
• Optionally, a modified depth value

The Phong lighting model — the classic real-time lighting approximation:

I = kₐ·Iₐ + kd·(L·N)·Id + ks·(R·V)ⁿ·Is

Where:

• kₐ, kd, ks = ambient, diffuse, specular material coefficients
• L = unit vector toward light, N = surface normal
• R = reflect(−L, N) = reflection direction
• V = unit vector toward camera
• n = shininess exponent
• Iₐ, Id, Is = ambient, diffuse, specular light colors

This runs per fragment — every pixel computes this equation. With millions of pixels per frame, fragment shaders dominate GPU workload in complex scenes.

Post-fragment operations (in fixed hardware after the shader):

• Scissor test: discard fragments outside a screen rectangle
• Alpha test / discard: fragment shader can call discard to kill a fragment (used for cutouts like leaves)
• Depth test: the z-buffer test described above
• Blending: combine fragment color with existing framebuffer color using configurable blend equations (used for transparency, particle effects, UI overlays)

// GLSL fragment shader — Phong shading
precision highp float;

in vec3 vNormal;      // interpolated world-space normal
in vec3 vFragPos;     // interpolated world-space position
in vec2 vTexCoord;    // interpolated texture coordinate

uniform sampler2D uDiffuseMap;
uniform vec3 uLightPos;
uniform vec3 uViewPos;
uniform vec3 uLightColor;

out vec4 fragColor;

void main() {
    vec3 N = normalize(vNormal);
    vec3 L = normalize(uLightPos - vFragPos);
    vec3 V = normalize(uViewPos  - vFragPos);
    vec3 R = reflect(-L, N);

    vec3 albedo  = texture(uDiffuseMap, vTexCoord).rgb;
    float diff   = max(dot(L, N), 0.0);
    float spec   = pow(max(dot(R, V), 0.0), 32.0);

    vec3 ambient  = 0.1  * albedo;
    vec3 diffuse  = diff * albedo * uLightColor;
    vec3 specular = spec * uLightColor;

    fragColor = vec4(ambient + diffuse + specular, 1.0);
}

The graphics pipeline is not just an abstract sequence of operations — it maps directly onto the physical architecture of a GPU. Understanding this explains why GPUs are so fast and what their limitations are.

The GPU is a SIMD machine with thousands of cores: NVIDIA's CUDA architecture groups threads into warps of 32. All 32 threads in a warp execute the same instruction at the same time (SIMD). If different threads take different branches (if/else), both branches execute and results are masked — this is called warp divergence and halves throughput.

The vertex shader runs once per vertex. Vertices are independent so any warp can process any subset of vertices. Vertex shader occupancy is typically high.

Rasterization is implemented in fixed-function hardware (not programmable) on modern GPUs. It uses a tile-based approach: the screen is divided into small tiles (16×16 or 32×32 pixels), and each tile can be tested against a triangle's bounding box in parallel. Tiles that overlap a triangle are handed off to a rasterization unit that generates fragments.

The fragment shader runs once per fragment. Fragments from the same triangle primitive are scheduled together for cache efficiency (their texture coordinates are spatially close, so texture cache hits are likely). Fragments are also launched in 2×2 quads even at triangle edges — the extra fragments are used to compute texture coordinate derivatives (dFdx/dFdy) for mipmap selection, then discarded.

Texture caches and bandwidth: texture fetching (sampling) is the most memory-bandwidth-intensive operation. GPUs have L1 and L2 texture caches with specialized addressing hardware (bilinear filtering, mipmapping) implemented in hardware. A texture cache miss stalls the warp; the GPU hides this latency by switching to another ready warp — this is latency hiding through massive concurrency, the core architectural principle.

Modern additions: mesh shaders (replace vertex/geometry shaders with a more flexible compute-like stage), ray tracing cores (BVH traversal in hardware), tensor cores (matrix multiplication for AI workloads). But the rasterization pipeline described here remains the core rendering path.

// Mental model: the GPU's inner loop (vastly simplified)
// For each draw call:

// 1. VERTEX STAGE (massively parallel, one thread per vertex)
for each vertex v in parallel:
    v.clip_pos = uProj * uView * uModel * v.position
    v.normal   = normalMatrix * v.normal
    // ... other per-vertex computations

// 2. PRIMITIVE ASSEMBLY + CLIPPING (fixed function)
for each triangle (v0, v1, v2):
    clipped_tris = sutherland_hodgman_clip(v0, v1, v2)
    for each clipped_tri: emit to rasterizer

// 3. RASTERIZATION (fixed function, tile-parallel)
for each triangle in parallel (tile-based):
    for each fragment (x,y) in triangle bounding box:
        if inside_triangle(x, y, triangle):
            emit fragment with interpolated attributes

// 4. FRAGMENT STAGE (massively parallel, one thread per fragment)
for each fragment f in parallel:
    if early_z_pass(f):         // early depth test
        f.color = fragment_shader(f)  // programmable
        depth_test_and_write(f)

Triangles in 3D are mapped to 2D images (textures) using UV coordinates — 2D coordinates in [0,1]² that are stored per-vertex and interpolated (perspective-correctly!) to each fragment. The fragment shader then samples the texture at the interpolated UV to get the surface color.

Texture filtering: a texture has discrete texels, but fragments map to non-integer UV positions. Nearest-neighbor filtering just takes the closest texel — fast but aliased (blocky). Bilinear filtering blends the four nearest texels weighted by distance — smooth and very hardware-efficient.

Mipmapping solves a different problem: aliasing from minification. When a large texture is viewed from far away, many texels map to the same pixel — averaging all of them is expensive, and sampling only one produces aliasing (shimmering noise). The solution is to precompute a mip chain: a sequence of textures at half resolution (1×, 1/2×, 1/4×, ...) down to 1×1. The hardware selects the mip level whose texels are approximately pixel-sized, then samples from that level.

The mip level is selected from the texture coordinate gradient — how quickly the UV changes across neighboring pixels. The GPU computes this for a 2×2 quad of fragments simultaneously using the dFdx(uv) and dFdy(uv) derivative instructions. LOD = log2(max(|dFdx(uv)|, |dFdy(uv)|) * texture_resolution). Trilinear filtering blends between adjacent mip levels to avoid popping.

Anisotropic filtering (AF) handles the case where the surface is viewed at a glancing angle: the pixel footprint is elongated in texture space, so we need to sample multiple times along the elongated axis. 16× AF takes up to 16 samples — it's the most expensive standard filter mode but gives dramatically sharper textures on floors and walls viewed at low angles.

This is why you see options like "Texture Quality: Bilinear / Trilinear / Anisotropic 4× / 16×" in game settings — they trade sample count for texture sharpness.

// GLSL: explicit LOD computation (what hardware does automatically)
vec2 uv = vTexCoord * textureSize(uTex, 0); // UV in texels
vec2 ddx = dFdx(uv);
vec2 ddy = dFdy(uv);
float lod = 0.5 * log2(max(dot(ddx, ddx), dot(ddy, ddy)));

// Sample with explicit LOD (for custom filtering):
vec4 color = textureLod(uTex, vTexCoord, lod);

// Or just let the hardware do it automatically (trilinear):
vec4 color2 = texture(uTex, vTexCoord); // hardware selects LOD

The fundamental problem: triangles have sharp edges, but the framebuffer is a discrete grid of pixels. A pixel is either inside or outside the triangle — this binary decision at pixel centers causes jagged edges (jaggies). This is spatial aliasing, the graphics manifestation of the Nyquist sampling theorem.

MSAA (Multisample Anti-Aliasing): the most common hardware solution. Instead of one sample per pixel, take N samples at different subpixel positions (4×, 8×, 16× are common). The triangle coverage test runs at the full sample rate, but the fragment shader runs once per pixel (not once per sample) — the depth and coverage mask are stored per sample, and the final color is the average over covered samples. This gives smooth edges at the cost of N× more depth/stencil buffer storage and the coverage test overhead, but NOT N× more fragment shader invocations.

SSAA (Supersample Anti-Aliasing): render at 2× or 4× resolution then downsample. Highest quality but fully N× more expensive (fragment shaders included). Too expensive for real-time.

TAA (Temporal Anti-Aliasing): jitter the camera sub-pixel each frame and blend with previous frames using a reprojection buffer. Essentially samples the same pixel at different subpixel positions across time. Very effective (smooth, essentially free), but introduces ghosting on fast-moving objects and requires careful velocity buffer handling.

FXAA / SMAA: post-process filters that detect and blur edges in the final image based on local contrast. Cheap and effective on static images but blurs textures and misses geometric edges not visible in the color buffer.

FSR / DLSS: upscaling techniques (not strictly AA). Render at lower resolution and upscale using spatial filters (FSR) or neural networks trained on high-resolution references (DLSS). Dramatically improve performance while maintaining reasonable quality.

The choice of AA technique encapsulates tradeoffs between quality, performance, memory bandwidth, and temporal stability — a recurring theme in real-time rendering.

// MSAA conceptual model (4x MSAA)
// Each pixel has 4 sample slots at subpixel offsets
// Rasterizer runs coverage test at all 4 samples
// Fragment shader runs ONCE if any sample is covered
// Final output = average of covered samples

struct Pixel {
    vec3 color;           // result of ONE fragment shader call
    float depth[4];       // depth at each of 4 MSAA samples
    bool  covered[4];     // coverage mask: which samples hit?
};

// Resolve pass: average covered sample colors
vec3 resolve(Pixel p) {
    int count = 0;
    for (int i = 0; i < 4; i++) if (p.covered[i]) count++;
    // Full coverage (all 4) → pure color; partial → blend with background
    return (p.color * count + background * (4 - count)) / 4;
}