ffmpeg/libswscale/csputils.c

/*
 * Copyright (C) 2024 Niklas Haas
 *
 * This file is part of FFmpeg.
 *
 * FFmpeg is free software; you can redistribute it and/or
 * modify it under the terms of the GNU Lesser General Public
 * License as published by the Free Software Foundation; either
 * version 2.1 of the License, or (at your option) any later version.
 *
 * FFmpeg is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * Lesser General Public License for more details.
 *
 * You should have received a copy of the GNU Lesser General Public
 * License along with FFmpeg; if not, write to the Free Software
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
 */

#include "libavutil/csp.h"

#include "csputils.h"
#include "utils.h"

void ff_sws_matrix3x3_mul(SwsMatrix3x3 *a, const SwsMatrix3x3 *b)
{
    float a00 = a->m[0][0], a01 = a->m[0][1], a02 = a->m[0][2],
          a10 = a->m[1][0], a11 = a->m[1][1], a12 = a->m[1][2],
          a20 = a->m[2][0], a21 = a->m[2][1], a22 = a->m[2][2];

    for (int i = 0; i < 3; i++) {
        a->m[0][i] = a00 * b->m[0][i] + a01 * b->m[1][i] + a02 * b->m[2][i];
        a->m[1][i] = a10 * b->m[0][i] + a11 * b->m[1][i] + a12 * b->m[2][i];
        a->m[2][i] = a20 * b->m[0][i] + a21 * b->m[1][i] + a22 * b->m[2][i];
    }
}

void ff_sws_matrix3x3_rmul(const SwsMatrix3x3 *a, SwsMatrix3x3 *b)
{
    float b00 = b->m[0][0], b01 = b->m[0][1], b02 = b->m[0][2],
          b10 = b->m[1][0], b11 = b->m[1][1], b12 = b->m[1][2],
          b20 = b->m[2][0], b21 = b->m[2][1], b22 = b->m[2][2];

    for (int i = 0; i < 3; i++) {
        b->m[i][0] = a->m[i][0] * b00 + a->m[i][1] * b10 + a->m[i][2] * b20;
        b->m[i][1] = a->m[i][0] * b01 + a->m[i][1] * b11 + a->m[i][2] * b21;
        b->m[i][2] = a->m[i][0] * b02 + a->m[i][1] * b12 + a->m[i][2] * b22;
    }
}

void ff_sws_matrix3x3_invert(SwsMatrix3x3 *mat)
{
    double m00 = mat->m[0][0], m01 = mat->m[0][1], m02 = mat->m[0][2],
           m10 = mat->m[1][0], m11 = mat->m[1][1], m12 = mat->m[1][2],
           m20 = mat->m[2][0], m21 = mat->m[2][1], m22 = mat->m[2][2];

    // calculate the adjoint
    double a00 =  (m11 * m22 - m21 * m12);
    double a01 = -(m01 * m22 - m21 * m02);
    double a02 =  (m01 * m12 - m11 * m02);
    double a10 = -(m10 * m22 - m20 * m12);
    double a11 =  (m00 * m22 - m20 * m02);
    double a12 = -(m00 * m12 - m10 * m02);
    double a20 =  (m10 * m21 - m20 * m11);
    double a21 = -(m00 * m21 - m20 * m01);
    double a22 =  (m00 * m11 - m10 * m01);

    // calculate the determinant (as inverse == 1/det * adjoint,
    // adjoint * m == identity * det, so this calculates the det)
    double det = m00 * a00 + m10 * a01 + m20 * a02;
    det = 1.0 / det;

    mat->m[0][0] = det * a00;
    mat->m[0][1] = det * a01;
    mat->m[0][2] = det * a02;
    mat->m[1][0] = det * a10;
    mat->m[1][1] = det * a11;
    mat->m[1][2] = det * a12;
    mat->m[2][0] = det * a20;
    mat->m[2][1] = det * a21;
    mat->m[2][2] = det * a22;
}

void ff_sws_matrix3x3_apply(const SwsMatrix3x3 *mat, float vec[3])
{
    float x = vec[0], y = vec[1], z = vec[2];

    for (int i = 0; i < 3; i++)
        vec[i] = mat->m[i][0] * x + mat->m[i][1] * y + mat->m[i][2] * z;
}

/* Recovers (X / Y) from a CIE xy value. */
static inline AVRational cie_X(AVCIExy xy)
{
    return av_div_q(xy.x, xy.y);
}

/* Recovers (Z / Y) from a CIE xy value. */
static inline AVRational cie_Z(AVCIExy xy)
{
    return av_div_q(av_sub_q(av_sub_q(av_make_q(1, 1), xy.x), xy.y), xy.y);
}

SwsMatrix3x3 ff_sws_rgb2xyz(const AVColorPrimariesDesc *desc)
{
    SwsMatrix3x3 out = {{{0}}};
    float S[3], X[3], Z[3], Xw, Zw;

    X[0] = av_q2d(cie_X(desc->prim.r));
    X[1] = av_q2d(cie_X(desc->prim.g));
    X[2] = av_q2d(cie_X(desc->prim.b));

    Z[0] = av_q2d(cie_Z(desc->prim.r));
    Z[1] = av_q2d(cie_Z(desc->prim.g));
    Z[2] = av_q2d(cie_Z(desc->prim.b));

    Xw = av_q2d(cie_X(desc->wp));
    Zw = av_q2d(cie_Z(desc->wp));

    /* S = XYZ^-1 * W */
    for (int i = 0; i < 3; i++) {
        out.m[0][i] = X[i];
        out.m[1][i] = 1.0f;
        out.m[2][i] = Z[i];
    }

    ff_sws_matrix3x3_invert(&out);

    for (int i = 0; i < 3; i++)
        S[i] = out.m[i][0] * Xw + out.m[i][1] + out.m[i][2] * Zw;

    /* M = [Sc * XYZc] */
    for (int i = 0; i < 3; i++) {
        out.m[0][i] = S[i] * X[i];
        out.m[1][i] = S[i];
        out.m[2][i] = S[i] * Z[i];
    }

    return out;
}

SwsMatrix3x3 ff_sws_xyz2rgb(const AVColorPrimariesDesc *prim)
{
    SwsMatrix3x3 out = ff_sws_rgb2xyz(prim);
    ff_sws_matrix3x3_invert(&out);
    return out;
}

/* Matrix used in CAT16, a revised one-step linear transform method */
static const SwsMatrix3x3 m_cat16 = {{
    {  0.401288f, 0.650173f, -0.051461f },
    { -0.250268f, 1.204414f,  0.045854f },
    { -0.002079f, 0.048952f,  0.953127f },
}};

static const SwsMatrix3x3 m_cat16_inv = {{
    {  1.862068f, -1.011255f,  0.149187f },
    {  0.387527f,  0.621447f, -0.008974f },
    { -0.015841f, -0.034123f,  1.049964f },
}};

// M := M * XYZd<-XYZs
static void apply_chromatic_adaptation(AVWhitepointCoefficients src,
                                       AVWhitepointCoefficients dst,
                                       SwsMatrix3x3 *mat)
{
    SwsMatrix3x3 tmp = {0};
    float C[3][2];
    const float srcX = av_q2d(cie_X(src)),
                srcZ = av_q2d(cie_Z(src)),
                dstX = av_q2d(cie_X(dst)),
                dstZ = av_q2d(cie_Z(dst));

    if (ff_cie_xy_equal(src, dst))
        return;

    // Linear "von Kries" method, adapted from CIECAM16
    // http://www.brucelindbloom.com/index.html?Eqn_ChromAdapt.html

    for (int i = 0; i < 3; i++) {
        // source cone
        C[i][0] = m_cat16.m[i][0] * srcX
                + m_cat16.m[i][1] * 1
                + m_cat16.m[i][2] * srcZ,

        // dest cone
        C[i][1] = m_cat16.m[i][0] * dstX
                + m_cat16.m[i][1] * 1
                + m_cat16.m[i][2] * dstZ;
    }

    // tmp := I * [Cd/Cs] * Ma
    for (int i = 0; i < 3; i++)
        tmp.m[i][i] = C[i][1] / C[i][0];

    ff_sws_matrix3x3_mul(&tmp, &m_cat16);

    // M := M * Ma^-1 * tmp
    ff_sws_matrix3x3_mul(mat, &m_cat16_inv);
    ff_sws_matrix3x3_mul(mat, &tmp);
}

SwsMatrix3x3 ff_sws_get_adaptation(const AVPrimaryCoefficients *prim,
                                   AVWhitepointCoefficients from,
                                   AVWhitepointCoefficients to)
{
    SwsMatrix3x3 rgb2xyz, xyz2rgb;
    const AVColorPrimariesDesc csp = {
        .prim = *prim,
        .wp = from,
    };

    rgb2xyz = ff_sws_rgb2xyz(&csp);
    xyz2rgb = ff_sws_xyz2rgb(&csp);
    apply_chromatic_adaptation(from, to, &xyz2rgb);
    ff_sws_matrix3x3_mul(&xyz2rgb, &rgb2xyz);
    return xyz2rgb;
}

static const AVWhitepointCoefficients d65 = {
    .x = {3127, 10000},
    .y = {3290, 10000}
};

static const SwsMatrix3x3 hpe = {{ /* HPE XYZ->LMS (D65) method */
    {  0.40024f, 0.70760f, -0.08081f },
    { -0.22630f, 1.16532f,  0.04570f },
    {  0.00000f, 0.00000f,  0.91822f },
}};

SwsMatrix3x3 ff_sws_ipt_rgb2lms(const AVColorPrimariesDesc *prim)
{
    const float c = 0.04f; // 4% crosstalk
    SwsMatrix3x3 rgb2xyz;
    SwsMatrix3x3 m = {{
        { 1 - 2*c,       c,       c },
        {       c, 1 - 2*c,       c },
        {       c,       c, 1 - 2*c },
    }};

    ff_sws_matrix3x3_mul(&m, &hpe);

    // Apply chromatic adaptation to D65 if the input white point differs
    apply_chromatic_adaptation(prim->wp, d65, &m);

    rgb2xyz = ff_sws_rgb2xyz(prim);
    ff_sws_matrix3x3_mul(&m, &rgb2xyz);
    return m;
}

SwsMatrix3x3 ff_sws_ipt_lms2rgb(const AVColorPrimariesDesc *prim)
{
    SwsMatrix3x3 rgb2lms = ff_sws_ipt_rgb2lms(prim);
    ff_sws_matrix3x3_invert(&rgb2lms);
    return rgb2lms;
}

/* Test the sign of 'p' relative to the line 'ab' (barycentric coordinates) */
static int test_point_line(AVCIExy p, AVCIExy a, AVCIExy b)
{
    return av_cmp_q(av_mul_q(av_sub_q(p.x, b.x), av_sub_q(a.y, b.y)),
                    av_mul_q(av_sub_q(a.x, b.x), av_sub_q(p.y, b.y)));
}


/* Test if a point is entirely inside a gamut */
static float test_point_gamut(const AVCIExy point,
                              const AVPrimaryCoefficients *prim)
{
    int d1 = test_point_line(point, prim->r, prim->g),
        d2 = test_point_line(point, prim->g, prim->b),
        d3 = test_point_line(point, prim->b, prim->r);

    bool has_neg = d1 < 0 || d2 < 0 || d3 < 0,
         has_pos = d1 > 0 || d2 > 0 || d3 > 0;

    return !(has_neg && has_pos);
}

bool ff_prim_superset(const AVPrimaryCoefficients *a, const AVPrimaryCoefficients *b)
{
    return test_point_gamut(b->r, a) &&
           test_point_gamut(b->g, a) &&
           test_point_gamut(b->b, a);
}

const float ff_pq_eotf_lut[PQ_LUT_SIZE+1] = {
    0.0000000e+00f, 4.0422718e-05f, 1.3111372e-04f, 2.6236826e-04f, 4.3151495e-04f, 6.3746885e-04f, 8.7982383e-04f, 1.1585362e-03f,
    1.4737819e-03f, 1.8258818e-03f, 2.2152586e-03f, 2.6424098e-03f, 3.1078907e-03f, 3.6123021e-03f, 4.1562821e-03f, 4.7405001e-03f,
    5.3656521e-03f, 6.0324583e-03f, 6.7416568e-03f, 7.4940095e-03f, 8.2902897e-03f, 9.1312924e-03f, 1.0017822e-02f, 1.0950702e-02f,
    1.1930764e-02f, 1.2958861e-02f, 1.4035847e-02f, 1.5162600e-02f, 1.6340000e-02f, 1.7568948e-02f, 1.8850346e-02f, 2.0185119e-02f,
    2.1574192e-02f, 2.3018509e-02f, 2.4519029e-02f, 2.6076704e-02f, 2.7692516e-02f, 2.9367449e-02f, 3.1102509e-02f, 3.2898690e-02f,
    3.4757019e-02f, 3.6678526e-02f, 3.8664261e-02f, 4.0715262e-02f, 4.2832601e-02f, 4.5017354e-02f, 4.7270617e-02f, 4.9593473e-02f,
    5.1987040e-02f, 5.4452441e-02f, 5.6990819e-02f, 5.9603301e-02f, 6.2291055e-02f, 6.5055251e-02f, 6.7897080e-02f, 7.0817717e-02f,
    7.3818379e-02f, 7.6900283e-02f, 8.0064675e-02f, 8.3312774e-02f, 8.6645849e-02f, 9.0065169e-02f, 9.3572031e-02f, 9.7167704e-02f,
    1.0085351e-01f, 1.0463077e-01f, 1.0850082e-01f, 1.1246501e-01f, 1.1652473e-01f, 1.2068130e-01f, 1.2493614e-01f, 1.2929066e-01f,
    1.3374626e-01f, 1.3830439e-01f, 1.4296648e-01f, 1.4773401e-01f, 1.5260848e-01f, 1.5759132e-01f, 1.6268405e-01f, 1.6788821e-01f,
    1.7320534e-01f, 1.7863697e-01f, 1.8418467e-01f, 1.8985004e-01f, 1.9563470e-01f, 2.0154019e-01f, 2.0756818e-01f, 2.1372031e-01f,
    2.1999824e-01f, 2.2640365e-01f, 2.3293824e-01f, 2.3960372e-01f, 2.4640186e-01f, 2.5333431e-01f, 2.6040288e-01f, 2.6760935e-01f,
    2.7495552e-01f, 2.8244319e-01f, 2.9007421e-01f, 2.9785041e-01f, 3.0577373e-01f, 3.1384594e-01f, 3.2206899e-01f, 3.3044481e-01f,
    3.3897533e-01f, 3.4766253e-01f, 3.5650838e-01f, 3.6551487e-01f, 3.7468409e-01f, 3.8401794e-01f, 3.9351855e-01f, 4.0318799e-01f,
    4.1302836e-01f, 4.2304177e-01f, 4.3323036e-01f, 4.4359629e-01f, 4.5414181e-01f, 4.6486897e-01f, 4.7578006e-01f, 4.8687732e-01f,
    4.9816302e-01f, 5.0963944e-01f, 5.2130889e-01f, 5.3317369e-01f, 5.4523628e-01f, 5.5749886e-01f, 5.6996391e-01f, 5.8263384e-01f,
    5.9551111e-01f, 6.0859816e-01f, 6.2189750e-01f, 6.3541162e-01f, 6.4914307e-01f, 6.6309439e-01f, 6.7726819e-01f, 6.9166705e-01f,
    7.0629384e-01f, 7.2115077e-01f, 7.3624074e-01f, 7.5156646e-01f, 7.6713065e-01f, 7.8293608e-01f, 7.9898553e-01f, 8.1528181e-01f,
    8.3182776e-01f, 8.4862623e-01f, 8.6568012e-01f, 8.8299235e-01f, 9.0056585e-01f, 9.1840360e-01f, 9.3650860e-01f, 9.5488388e-01f,
    9.7353277e-01f, 9.9245779e-01f, 1.0116623e+00f, 1.0311496e+00f, 1.0509226e+00f, 1.0709847e+00f, 1.0913391e+00f, 1.1119889e+00f,
    1.1329376e+00f, 1.1541885e+00f, 1.1757448e+00f, 1.1976100e+00f, 1.2197875e+00f, 1.2422807e+00f, 1.2650931e+00f, 1.2882282e+00f,
    1.3116900e+00f, 1.3354812e+00f, 1.3596059e+00f, 1.3840676e+00f, 1.4088701e+00f, 1.4340170e+00f, 1.4595121e+00f, 1.4853593e+00f,
    1.5115622e+00f, 1.5381247e+00f, 1.5650507e+00f, 1.5923442e+00f, 1.6200090e+00f, 1.6480492e+00f, 1.6764688e+00f, 1.7052718e+00f,
    1.7344629e+00f, 1.7640451e+00f, 1.7940233e+00f, 1.8244015e+00f, 1.8551840e+00f, 1.8863752e+00f, 1.9179792e+00f, 1.9500006e+00f,
    1.9824437e+00f, 2.0153130e+00f, 2.0486129e+00f, 2.0823479e+00f, 2.1165227e+00f, 2.1511419e+00f, 2.1862101e+00f, 2.2217319e+00f,
    2.2577128e+00f, 2.2941563e+00f, 2.3310679e+00f, 2.3684523e+00f, 2.4063146e+00f, 2.4446597e+00f, 2.4834925e+00f, 2.5228182e+00f,
    2.5626417e+00f, 2.6029683e+00f, 2.6438031e+00f, 2.6851514e+00f, 2.7270184e+00f, 2.7694094e+00f, 2.8123299e+00f, 2.8557852e+00f,
    2.8997815e+00f, 2.9443230e+00f, 2.9894159e+00f, 3.0350657e+00f, 3.0812783e+00f, 3.1280593e+00f, 3.1754144e+00f, 3.2233495e+00f,
    3.2718705e+00f, 3.3209833e+00f, 3.3706938e+00f, 3.4210082e+00f, 3.4719324e+00f, 3.5234727e+00f, 3.5756351e+00f, 3.6284261e+00f,
    3.6818526e+00f, 3.7359195e+00f, 3.7906340e+00f, 3.8460024e+00f, 3.9020315e+00f, 3.9587277e+00f, 4.0160977e+00f, 4.0741483e+00f,
    4.1328861e+00f, 4.1923181e+00f, 4.2524511e+00f, 4.3132921e+00f, 4.3748480e+00f, 4.4371260e+00f, 4.5001332e+00f, 4.5638768e+00f,
    4.6283650e+00f, 4.6936032e+00f, 4.7595999e+00f, 4.8263624e+00f, 4.8938982e+00f, 4.9622151e+00f, 5.0313205e+00f, 5.1012223e+00f,
    5.1719283e+00f, 5.2434463e+00f, 5.3157843e+00f, 5.3889502e+00f, 5.4629521e+00f, 5.5377982e+00f, 5.6134968e+00f, 5.6900560e+00f,
    5.7674843e+00f, 5.8457900e+00f, 5.9249818e+00f, 6.0050682e+00f, 6.0860578e+00f, 6.1679595e+00f, 6.2507819e+00f, 6.3345341e+00f,
    6.4192275e+00f, 6.5048661e+00f, 6.5914616e+00f, 6.6790231e+00f, 6.7675600e+00f, 6.8570816e+00f, 6.9475975e+00f, 7.0391171e+00f,
    7.1316500e+00f, 7.2252060e+00f, 7.3197948e+00f, 7.4154264e+00f, 7.5121107e+00f, 7.6098577e+00f, 7.7086777e+00f, 7.8085807e+00f,
    7.9095772e+00f, 8.0116775e+00f, 8.1148922e+00f, 8.2192318e+00f, 8.3247071e+00f, 8.4313287e+00f, 8.5391076e+00f, 8.6480548e+00f,
    8.7581812e+00f, 8.8694982e+00f, 8.9820168e+00f, 9.0957485e+00f, 9.2107048e+00f, 9.3268971e+00f, 9.4443372e+00f, 9.5630368e+00f,
    9.6830115e+00f, 9.8042658e+00f, 9.9268155e+00f, 1.0050673e+01f, 1.0175850e+01f, 1.0302359e+01f, 1.0430213e+01f, 1.0559425e+01f,
    1.0690006e+01f, 1.0821970e+01f, 1.0955331e+01f, 1.1090100e+01f, 1.1226290e+01f, 1.1363917e+01f, 1.1502992e+01f, 1.1643529e+01f,
    1.1785542e+01f, 1.1929044e+01f, 1.2074050e+01f, 1.2220573e+01f, 1.2368628e+01f, 1.2518229e+01f, 1.2669390e+01f, 1.2822125e+01f,
    1.2976449e+01f, 1.3132377e+01f, 1.3289925e+01f, 1.3449105e+01f, 1.3609935e+01f, 1.3772429e+01f, 1.3936602e+01f, 1.4102470e+01f,
    1.4270054e+01f, 1.4439360e+01f, 1.4610407e+01f, 1.4783214e+01f, 1.4957794e+01f, 1.5134166e+01f, 1.5312345e+01f, 1.5492348e+01f,
    1.5674192e+01f, 1.5857894e+01f, 1.6043471e+01f, 1.6230939e+01f, 1.6420317e+01f, 1.6611622e+01f, 1.6804871e+01f, 1.7000083e+01f,
    1.7197275e+01f, 1.7396465e+01f, 1.7597672e+01f, 1.7800914e+01f, 1.8006210e+01f, 1.8213578e+01f, 1.8423038e+01f, 1.8634608e+01f,
    1.8848308e+01f, 1.9064157e+01f, 1.9282175e+01f, 1.9502381e+01f, 1.9724796e+01f, 1.9949439e+01f, 2.0176331e+01f, 2.0405492e+01f,
    2.0636950e+01f, 2.0870711e+01f, 2.1106805e+01f, 2.1345250e+01f, 2.1586071e+01f, 2.1829286e+01f, 2.2074919e+01f, 2.2322992e+01f,
    2.2573525e+01f, 2.2826542e+01f, 2.3082066e+01f, 2.3340118e+01f, 2.3600721e+01f, 2.3863900e+01f, 2.4129676e+01f, 2.4398074e+01f,
    2.4669117e+01f, 2.4942828e+01f, 2.5219233e+01f, 2.5498355e+01f, 2.5780219e+01f, 2.6064849e+01f, 2.6352271e+01f, 2.6642509e+01f,
    2.6935589e+01f, 2.7231536e+01f, 2.7530377e+01f, 2.7832137e+01f, 2.8136843e+01f, 2.8444520e+01f, 2.8755196e+01f, 2.9068898e+01f,
    2.9385662e+01f, 2.9705496e+01f, 3.0028439e+01f, 3.0354517e+01f, 3.0683758e+01f, 3.1016192e+01f, 3.1351846e+01f, 3.1690750e+01f,
    3.2032932e+01f, 3.2378422e+01f, 3.2727250e+01f, 3.3079445e+01f, 3.3435038e+01f, 3.3794058e+01f, 3.4156537e+01f, 3.4522505e+01f,
    3.4891993e+01f, 3.5265034e+01f, 3.5641658e+01f, 3.6021897e+01f, 3.6405785e+01f, 3.6793353e+01f, 3.7184634e+01f, 3.7579661e+01f,
    3.7978468e+01f, 3.8381088e+01f, 3.8787555e+01f, 3.9197904e+01f, 3.9612169e+01f, 4.0030385e+01f, 4.0452587e+01f, 4.0878810e+01f,
    4.1309104e+01f, 4.1743478e+01f, 4.2181981e+01f, 4.2624651e+01f, 4.3071525e+01f, 4.3522639e+01f, 4.3978031e+01f, 4.4437739e+01f,
    4.4901803e+01f, 4.5370259e+01f, 4.5843148e+01f, 4.6320508e+01f, 4.6802379e+01f, 4.7288801e+01f, 4.7779815e+01f, 4.8275461e+01f,
    4.8775780e+01f, 4.9280813e+01f, 4.9790603e+01f, 5.0305191e+01f, 5.0824620e+01f, 5.1348933e+01f, 5.1878172e+01f, 5.2412382e+01f,
    5.2951607e+01f, 5.3495890e+01f, 5.4045276e+01f, 5.4599811e+01f, 5.5159540e+01f, 5.5724510e+01f, 5.6294765e+01f, 5.6870353e+01f,
    5.7451339e+01f, 5.8037735e+01f, 5.8629606e+01f, 5.9227001e+01f, 5.9829968e+01f, 6.0438557e+01f, 6.1052818e+01f, 6.1672799e+01f,
    6.2298552e+01f, 6.2930128e+01f, 6.3567578e+01f, 6.4210953e+01f, 6.4860306e+01f, 6.5515690e+01f, 6.6177157e+01f, 6.6844762e+01f,
    6.7518558e+01f, 6.8198599e+01f, 6.8884942e+01f, 6.9577641e+01f, 7.0276752e+01f, 7.0982332e+01f, 7.1694438e+01f, 7.2413127e+01f,
    7.3138457e+01f, 7.3870486e+01f, 7.4609273e+01f, 7.5354878e+01f, 7.6107361e+01f, 7.6866782e+01f, 7.7633203e+01f, 7.8406684e+01f,
    7.9187312e+01f, 7.9975101e+01f, 8.0770139e+01f, 8.1572490e+01f, 8.2382216e+01f, 8.3199385e+01f, 8.4024059e+01f, 8.4856307e+01f,
    8.5696193e+01f, 8.6543786e+01f, 8.7399153e+01f, 8.8262362e+01f, 8.9133482e+01f, 9.0012582e+01f, 9.0899733e+01f, 9.1795005e+01f,
    9.2698470e+01f, 9.3610199e+01f, 9.4530265e+01f, 9.5458741e+01f, 9.6395701e+01f, 9.7341219e+01f, 9.8295370e+01f, 9.9258231e+01f,
    1.0022988e+02f, 1.0121039e+02f, 1.0219984e+02f, 1.0319830e+02f, 1.0420587e+02f, 1.0522261e+02f, 1.0624862e+02f, 1.0728396e+02f,
    1.0832872e+02f, 1.0938299e+02f, 1.1044684e+02f, 1.1152036e+02f, 1.1260365e+02f, 1.1369677e+02f, 1.1479982e+02f, 1.1591288e+02f,
    1.1703605e+02f, 1.1816941e+02f, 1.1931305e+02f, 1.2046706e+02f, 1.2163153e+02f, 1.2280656e+02f, 1.2399223e+02f, 1.2518864e+02f,
    1.2639596e+02f, 1.2761413e+02f, 1.2884333e+02f, 1.3008365e+02f, 1.3133519e+02f, 1.3259804e+02f, 1.3387231e+02f, 1.3515809e+02f,
    1.3645549e+02f, 1.3776461e+02f, 1.3908555e+02f, 1.4041841e+02f, 1.4176331e+02f, 1.4312034e+02f, 1.4448961e+02f, 1.4587123e+02f,
    1.4726530e+02f, 1.4867194e+02f, 1.5009126e+02f, 1.5152336e+02f, 1.5296837e+02f, 1.5442638e+02f, 1.5589753e+02f, 1.5738191e+02f,
    1.5887965e+02f, 1.6039087e+02f, 1.6191567e+02f, 1.6345419e+02f, 1.6500655e+02f, 1.6657285e+02f, 1.6815323e+02f, 1.6974781e+02f,
    1.7135672e+02f, 1.7298007e+02f, 1.7461800e+02f, 1.7627063e+02f, 1.7793810e+02f, 1.7962053e+02f, 1.8131805e+02f, 1.8303080e+02f,
    1.8475891e+02f, 1.8650252e+02f, 1.8826176e+02f, 1.9003676e+02f, 1.9182767e+02f, 1.9363463e+02f, 1.9545777e+02f, 1.9729724e+02f,
    1.9915319e+02f, 2.0102575e+02f, 2.0291507e+02f, 2.0482131e+02f, 2.0674460e+02f, 2.0868510e+02f, 2.1064296e+02f, 2.1261833e+02f,
    2.1461136e+02f, 2.1662222e+02f, 2.1865105e+02f, 2.2069802e+02f, 2.2276328e+02f, 2.2484700e+02f, 2.2694934e+02f, 2.2907045e+02f,
    2.3121064e+02f, 2.3336982e+02f, 2.3554827e+02f, 2.3774618e+02f, 2.3996370e+02f, 2.4220102e+02f, 2.4445831e+02f, 2.4673574e+02f,
    2.4903349e+02f, 2.5135174e+02f, 2.5369067e+02f, 2.5605046e+02f, 2.5843129e+02f, 2.6083336e+02f, 2.6325684e+02f, 2.6570192e+02f,
    2.6816880e+02f, 2.7065767e+02f, 2.7316872e+02f, 2.7570215e+02f, 2.7825815e+02f, 2.8083692e+02f, 2.8343867e+02f, 2.8606359e+02f,
    2.8871189e+02f, 2.9138378e+02f, 2.9407946e+02f, 2.9679914e+02f, 2.9954304e+02f, 3.0231137e+02f, 3.0510434e+02f, 3.0792217e+02f,
    3.1076508e+02f, 3.1363330e+02f, 3.1652704e+02f, 3.1944653e+02f, 3.2239199e+02f, 3.2536367e+02f, 3.2836178e+02f, 3.3138657e+02f,
    3.3443826e+02f, 3.3751710e+02f, 3.4062333e+02f, 3.4375718e+02f, 3.4691890e+02f, 3.5010874e+02f, 3.5332694e+02f, 3.5657377e+02f,
    3.5984946e+02f, 3.6315428e+02f, 3.6648848e+02f, 3.6985233e+02f, 3.7324608e+02f, 3.7667000e+02f, 3.8012436e+02f, 3.8360942e+02f,
    3.8712547e+02f, 3.9067276e+02f, 3.9425159e+02f, 3.9786223e+02f, 4.0150496e+02f, 4.0518006e+02f, 4.0888783e+02f, 4.1262855e+02f,
    4.1640274e+02f, 4.2021025e+02f, 4.2405159e+02f, 4.2792707e+02f, 4.3183699e+02f, 4.3578166e+02f, 4.3976138e+02f, 4.4377647e+02f,
    4.4782724e+02f, 4.5191401e+02f, 4.5603709e+02f, 4.6019681e+02f, 4.6439350e+02f, 4.6862749e+02f, 4.7289910e+02f, 4.7720867e+02f,
    4.8155654e+02f, 4.8594305e+02f, 4.9036854e+02f, 4.9483336e+02f, 4.9933787e+02f, 5.0388240e+02f, 5.0846733e+02f, 5.1309301e+02f,
    5.1775981e+02f, 5.2246808e+02f, 5.2721821e+02f, 5.3201056e+02f, 5.3684551e+02f, 5.4172344e+02f, 5.4664473e+02f, 5.5160978e+02f,
    5.5661897e+02f, 5.6167269e+02f, 5.6677135e+02f, 5.7191535e+02f, 5.7710508e+02f, 5.8234097e+02f, 5.8762342e+02f, 5.9295285e+02f,
    5.9832968e+02f, 6.0375433e+02f, 6.0922723e+02f, 6.1474882e+02f, 6.2031952e+02f, 6.2593979e+02f, 6.3161006e+02f, 6.3733078e+02f,
    6.4310241e+02f, 6.4892540e+02f, 6.5480021e+02f, 6.6072730e+02f, 6.6670715e+02f, 6.7274023e+02f, 6.7882702e+02f, 6.8496800e+02f,
    6.9116365e+02f, 6.9741447e+02f, 7.0372096e+02f, 7.1008361e+02f, 7.1650293e+02f, 7.2297942e+02f, 7.2951361e+02f, 7.3610602e+02f,
    7.4275756e+02f, 7.4946797e+02f, 7.5623818e+02f, 7.6306873e+02f, 7.6996016e+02f, 7.7691302e+02f, 7.8392787e+02f, 7.9100526e+02f,
    7.9814576e+02f, 8.0534993e+02f, 8.1261837e+02f, 8.1995163e+02f, 8.2735032e+02f, 8.3481501e+02f, 8.4234632e+02f, 8.4994483e+02f,
    8.5761116e+02f, 8.6534592e+02f, 8.7314974e+02f, 8.8102323e+02f, 8.8896702e+02f, 8.9698176e+02f, 9.0506809e+02f, 9.1322665e+02f,
    9.2145810e+02f, 9.2976310e+02f, 9.3814232e+02f, 9.4659643e+02f, 9.5512612e+02f, 9.6373206e+02f, 9.7241496e+02f, 9.8117550e+02f,
    9.9001441e+02f, 9.9893238e+02f, 1.0079301e+03f, 1.0170084e+03f, 1.0261679e+03f, 1.0354094e+03f, 1.0447337e+03f, 1.0541414e+03f,
    1.0636334e+03f, 1.0732104e+03f, 1.0828731e+03f, 1.0926225e+03f, 1.1024592e+03f, 1.1123841e+03f, 1.1223979e+03f, 1.1325016e+03f,
    1.1426958e+03f, 1.1529814e+03f, 1.1633594e+03f, 1.1738304e+03f, 1.1843954e+03f, 1.1950552e+03f, 1.2058107e+03f, 1.2166627e+03f,
    1.2276122e+03f, 1.2386601e+03f, 1.2498072e+03f, 1.2610544e+03f, 1.2724027e+03f, 1.2838531e+03f, 1.2954063e+03f, 1.3070635e+03f,
    1.3188262e+03f, 1.3306940e+03f, 1.3426686e+03f, 1.3547509e+03f, 1.3669420e+03f, 1.3792428e+03f, 1.3916544e+03f, 1.4041778e+03f,
    1.4168140e+03f, 1.4295640e+03f, 1.4424289e+03f, 1.4554098e+03f, 1.4685078e+03f, 1.4817238e+03f, 1.4950591e+03f, 1.5085147e+03f,
    1.5220916e+03f, 1.5357912e+03f, 1.5496144e+03f, 1.5635624e+03f, 1.5776364e+03f, 1.5918375e+03f, 1.6061670e+03f, 1.6206260e+03f,
    1.6352156e+03f, 1.6499372e+03f, 1.6647920e+03f, 1.6797811e+03f, 1.6949059e+03f, 1.7101676e+03f, 1.7255674e+03f, 1.7411067e+03f,
    1.7567867e+03f, 1.7726087e+03f, 1.7885742e+03f, 1.8046844e+03f, 1.8209406e+03f, 1.8373443e+03f, 1.8538967e+03f, 1.8705994e+03f,
    1.8874536e+03f, 1.9044608e+03f, 1.9216225e+03f, 1.9389401e+03f, 1.9564150e+03f, 1.9740486e+03f, 1.9918426e+03f, 2.0097984e+03f,
    2.0279175e+03f, 2.0462014e+03f, 2.0646517e+03f, 2.0832699e+03f, 2.1020577e+03f, 2.1210165e+03f, 2.1401481e+03f, 2.1594540e+03f,
    2.1789359e+03f, 2.1985954e+03f, 2.2184342e+03f, 2.2384540e+03f, 2.2586565e+03f, 2.2790434e+03f, 2.2996165e+03f, 2.3203774e+03f,
    2.3413293e+03f, 2.3624714e+03f, 2.3838068e+03f, 2.4053372e+03f, 2.4270646e+03f, 2.4489908e+03f, 2.4711177e+03f, 2.4934471e+03f,
    2.5159811e+03f, 2.5387214e+03f, 2.5616702e+03f, 2.5848293e+03f, 2.6082007e+03f, 2.6317866e+03f, 2.6555888e+03f, 2.6796095e+03f,
    2.7038507e+03f, 2.7283145e+03f, 2.7530031e+03f, 2.7779186e+03f, 2.8030631e+03f, 2.8284388e+03f, 2.8540479e+03f, 2.8798927e+03f,
    2.9059754e+03f, 2.9322983e+03f, 2.9588635e+03f, 2.9856736e+03f, 3.0127308e+03f, 3.0400374e+03f, 3.0675959e+03f, 3.0954086e+03f,
    3.1234780e+03f, 3.1518066e+03f, 3.1803969e+03f, 3.2092512e+03f, 3.2383723e+03f, 3.2677625e+03f, 3.2974246e+03f, 3.3273611e+03f,
    3.3575747e+03f, 3.3880680e+03f, 3.4188437e+03f, 3.4499045e+03f, 3.4812533e+03f, 3.5128926e+03f, 3.5448255e+03f, 3.5770546e+03f,
    3.6095828e+03f, 3.6424131e+03f, 3.6755483e+03f, 3.7089914e+03f, 3.7427454e+03f, 3.7768132e+03f, 3.8111979e+03f, 3.8459027e+03f,
    3.8809304e+03f, 3.9162844e+03f, 3.9519678e+03f, 3.9879837e+03f, 4.0243354e+03f, 4.0610261e+03f, 4.0980592e+03f, 4.1354380e+03f,
    4.1731681e+03f, 4.2112483e+03f, 4.2496844e+03f, 4.2884798e+03f, 4.3276381e+03f, 4.3671627e+03f, 4.4070572e+03f, 4.4473253e+03f,
    4.4879706e+03f, 4.5289968e+03f, 4.5704076e+03f, 4.6122068e+03f, 4.6543981e+03f, 4.6969854e+03f, 4.7399727e+03f, 4.7833637e+03f,
    4.8271625e+03f, 4.8713731e+03f, 4.9159995e+03f, 4.9610458e+03f, 5.0065162e+03f, 5.0524147e+03f, 5.0987457e+03f, 5.1455133e+03f,
    5.1927219e+03f, 5.2403759e+03f, 5.2884795e+03f, 5.3370373e+03f, 5.3860537e+03f, 5.4355333e+03f, 5.4854807e+03f, 5.5359004e+03f,
    5.5867972e+03f, 5.6381757e+03f, 5.6900408e+03f, 5.7423972e+03f, 5.7952499e+03f, 5.8486037e+03f, 5.9024637e+03f, 5.9568349e+03f,
    6.0117223e+03f, 6.0671311e+03f, 6.1230664e+03f, 6.1795336e+03f, 6.2365379e+03f, 6.2940847e+03f, 6.3521793e+03f, 6.4108273e+03f,
    6.4700342e+03f, 6.5298056e+03f, 6.5901471e+03f, 6.6510643e+03f, 6.7125632e+03f, 6.7746495e+03f, 6.8373290e+03f, 6.9006078e+03f,
    6.9644918e+03f, 7.0289872e+03f, 7.0941001e+03f, 7.1598366e+03f, 7.2262031e+03f, 7.2932059e+03f, 7.3608513e+03f, 7.4291460e+03f,
    7.4981006e+03f, 7.5677134e+03f, 7.6379952e+03f, 7.7089527e+03f, 7.7805929e+03f, 7.8529226e+03f, 7.9259489e+03f, 7.9996786e+03f,
    8.0741191e+03f, 8.1492774e+03f, 8.2251609e+03f, 8.3017769e+03f, 8.3791329e+03f, 8.4572364e+03f, 8.5360950e+03f, 8.6157163e+03f,
    8.6961082e+03f, 8.7772786e+03f, 8.8592352e+03f, 8.9419862e+03f, 9.0255397e+03f, 9.1099038e+03f, 9.1950869e+03f, 9.2810973e+03f,
    9.3679435e+03f, 9.4556340e+03f, 9.5441776e+03f, 9.6335829e+03f, 9.7238588e+03f, 9.8150143e+03f, 9.9070583e+03f, 1.0000000e+04f,
    1e4f, /* extra padding to avoid out of bounds access */
};