builtin-programs/gpu/gpu-fns.folk

Wish the GPU compiles function "rotate" {{vec2 v float a} vec2 {
    float s = sin(a);
    float c = cos(a);
    mat2 m = mat2(c, s, -s, c);
    return m * v;
}}

Wish the GPU compiles function "wedge2d" {{vec2 v vec2 w} float {
    return v.x * w.y - v.y * w.x;
}}
# See https://www.reedbeta.com/blog/quadrilateral-interpolation-part-2/
Wish the GPU compiles function "invBilinear" {
      {vec2 pos vec2 p2 vec2 p3 vec2 p1 vec2 p0 fn wedge2d} vec2 {

    vec2 q = pos - p0;
    vec2 b1 = p1 - p0;
    vec2 b2 = p2 - p0;
    vec2 b3 = p0 - p1 - p2 + p3;

    // Set up quadratic formula
    float A = wedge2d(b2, b3);
    float B = wedge2d(b3, q) - wedge2d(b1, b2);
    float C = wedge2d(b1, q);

    // Solve for v
    vec2 uv;
    if (abs(A) < 0.001) {
        // Linear form
        uv.y = -C / B;
    } else {
        // Quadratic form. Take positive root for CCW winding with V-up
        float discrim = B * B - 4.0 * A * C;
        // Fix from https://www.reedbeta.com/blog/quadrilateral-interpolation-part-2/#comment-3629608573
        // (makes images not break with certain vertex inputs)
        uv.y = 2.0 * C / (-B - sqrt(discrim));
        // uv.y = 0.5 * (-B + sqrt(discrim)) / A;
    }
    
    // Solve for u, using largest-magnitude component
    vec2 denom = b1 + uv.y * b3;
    if (abs(denom.x) > abs(denom.y))
        uv.x = (q.x - b2.x * uv.y) / denom.x;
    else
        uv.x = (q.y - b2.y * uv.y) / denom.y;

    uv.y = 1.0 - uv.y;
    return uv;
}}