"""Unlit orthographic raster evidence; bounding-box alignment is not registration."""
import numpy as np


def rasterize(parts, bounds, shape, front):
    height, width = shape
    front = np.asarray(front, dtype=float)
    right = -np.cross(front, [0, 0, 1])
    all_vertices = np.concatenate([p['vertices'] for p in parts])
    lo = np.array([(all_vertices @ right).min(), all_vertices[:, 2].min()])
    span = np.array([np.ptp(all_vertices @ right), np.ptp(all_vertices[:, 2])])
    if np.any(span <= 0):
        raise ValueError('Degenerate projection bounds')
    x0, y0, x1, y1 = bounds
    rgb = np.full((height, width, 3), 220, dtype=np.uint8)
    depth = np.full((height, width), -np.inf)
    front_facing = np.zeros((height, width), dtype=bool)
    covered = np.zeros((height, width), dtype=bool)
    occluded = np.zeros((height, width), dtype=bool)
    for part in parts:
        vertices, faces, uvs, texture = [part[k] for k in ('vertices', 'faces', 'uvs', 'texture')]
        mode = part.get('alpha_mode', 'OPAQUE')
        if mode not in {'OPAQUE', 'MASK'}:
            raise ValueError('Unsupported diagnostic alpha mode: ' + mode)
        factor = np.asarray(part.get('base_color_factor', [1, 1, 1, 1]), dtype=float)
        screen = np.c_[(vertices @ right - lo[0]) / span[0] * (x1-x0) + x0,
                       y1 - (vertices[:, 2] - lo[1]) / span[1] * (y1-y0)]
        for face, uv in zip(faces, uvs, strict=True):
            p = screen[face]
            xmin, ymin = np.maximum(np.floor(p.min(0)).astype(int), 0)
            xmax, ymax = np.minimum(np.ceil(p.max(0)).astype(int), [width-1, height-1])
            if xmin > xmax or ymin > ymax:
                continue
            xx, yy = np.meshgrid(np.arange(xmin, xmax+1), np.arange(ymin, ymax+1))
            d = (p[1,1]-p[2,1])*(p[0,0]-p[2,0])+(p[2,0]-p[1,0])*(p[0,1]-p[2,1])
            if abs(d) < 1e-12:
                continue
            a = ((p[1,1]-p[2,1])*(xx-p[2,0])+(p[2,0]-p[1,0])*(yy-p[2,1]))/d
            b = ((p[2,1]-p[0,1])*(xx-p[2,0])+(p[0,0]-p[2,0])*(yy-p[2,1]))/d
            bary = np.stack([a,b,1-a-b], axis=-1)
            inside = (bary >= -1e-7).all(-1)
            xx, yy, bary = xx[inside], yy[inside], bary[inside]
            z = bary @ (vertices[face] @ front)
            coords = bary @ uv
            tx = np.clip(np.rint(coords[:,0]*(texture.shape[1]-1)).astype(int),0,texture.shape[1]-1)
            ty = np.clip(np.rint((1-coords[:,1])*(texture.shape[0]-1)).astype(int),0,texture.shape[0]-1)
            if mode == 'MASK':
                visible = texture[ty, tx, 3] / 255 * factor[3] >= part.get('alpha_cutoff', .5)
                xx, yy, z, tx, ty = xx[visible], yy[visible], z[visible], tx[visible], ty[visible]
            previous = depth[yy, xx]
            occluded[yy, xx] |= np.isfinite(previous) & (np.abs(previous-z)>1e-6)
            keep = z > previous
            xx, yy, z, tx, ty = xx[keep], yy[keep], z[keep], tx[keep], ty[keep]
            # Convert only visible sampled texels, never the full potentially 16K atlas.
            srgb = texture[ty, tx, :3].astype(float) / 255
            linear = np.where(srgb <= .04045, srgb / 12.92, ((srgb + .055) / 1.055) ** 2.4)
            linear *= factor[:3]
            corrected = np.where(linear <= .0031308, 12.92 * linear,
                                 1.055 * np.maximum(linear, 0) ** (1/2.4) - .055)
            rgb[yy,xx] = np.clip(np.rint(corrected * 255), 0, 255).astype(np.uint8)
            depth[yy,xx] = z
            normal = np.cross(vertices[face[1]]-vertices[face[0]],vertices[face[2]]-vertices[face[0]])
            front_facing[yy,xx] = normal @ front > 0
            covered[yy,xx] = True
    return dict(rgb=rgb, coverage=covered, front_facing=front_facing, occluded=occluded)