#ifndef FASTLED_INTERNAL
#define FASTLED_INTERNAL
#endif

#include "fl/system/fastled.h"

#include "fl/math/ease.h"
#include "lib8tion.h" // This is the problematic header that's hard to include

#include "fl/math/math.h"
#include "fl/math/intmap.h"
#include "fl/math/sin32.h"
#include "fl/stl/int.h"
#include "fl/math/math.h"
#include "fl/stl/weak_ptr.h"
#include "fl/stl/flat_map.h"
#include "fl/stl/align.h"
#include "fl/math/fixed_point.h"

namespace fl {

// Static PROGMEM lookup table for gamma 2.8 correction (8-bit -> 16-bit).
// Generated by: round(pow(i / 255.0, 2.8) * 65535.0) for i in 0..255.
// Uses flash storage on AVR/ESP, zero heap allocation on all platforms.
FL_ALIGN_PROGMEM(64) const u16 GAMMA_2_8_LUT[256] FL_PROGMEM = {
        0,     0,     0,     0,     1,     1,     2,     3,     4,     6,
        8,    10,    13,    16,    19,    24,    28,    33,    39,    46,
       53,    60,    69,    78,    88,    98,   110,   122,   135,   149,
      164,   179,   196,   214,   232,   252,   273,   295,   317,   341,
      366,   393,   420,   449,   478,   510,   542,   575,   610,   647,
      684,   723,   764,   806,   849,   894,   940,   988,  1037,  1088,
     1140,  1194,  1250,  1307,  1366,  1427,  1489,  1553,  1619,  1686,
     1756,  1827,  1900,  1975,  2051,  2130,  2210,  2293,  2377,  2463,
     2552,  2642,  2734,  2829,  2925,  3024,  3124,  3227,  3332,  3439,
     3548,  3660,  3774,  3890,  4008,  4128,  4251,  4376,  4504,  4634,
     4766,  4901,  5038,  5177,  5319,  5464,  5611,  5760,  5912,  6067,
     6224,  6384,  6546,  6711,  6879,  7049,  7222,  7397,  7576,  7757,
     7941,  8128,  8317,  8509,  8704,  8902,  9103,  9307,  9514,  9723,
     9936, 10151, 10370, 10591, 10816, 11043, 11274, 11507, 11744, 11984,
    12227, 12473, 12722, 12975, 13230, 13489, 13751, 14017, 14285, 14557,
    14833, 15111, 15393, 15678, 15967, 16259, 16554, 16853, 17155, 17461,
    17770, 18083, 18399, 18719, 19042, 19369, 19700, 20034, 20372, 20713,
    21058, 21407, 21759, 22115, 22475, 22838, 23206, 23577, 23952, 24330,
    24713, 25099, 25489, 25884, 26282, 26683, 27089, 27499, 27913, 28330,
    28752, 29178, 29608, 30041, 30479, 30921, 31367, 31818, 32272, 32730,
    33193, 33660, 34131, 34606, 35085, 35569, 36057, 36549, 37046, 37547,
    38052, 38561, 39075, 39593, 40116, 40643, 41175, 41711, 42251, 42796,
    43346, 43899, 44458, 45021, 45588, 46161, 46737, 47319, 47905, 48495,
    49091, 49691, 50295, 50905, 51519, 52138, 52761, 53390, 54023, 54661,
    55303, 55951, 56604, 57261, 57923, 58590, 59262, 59939, 60621, 61308,
    62000, 62697, 63399, 64106, 64818, 65535};

u16 gamma_2_8(u8 value) {
    return FL_PGM_READ_WORD_ALIGNED(&GAMMA_2_8_LUT[value]);
}

// 8-bit easing functions
u8 easeInQuad8(u8 i) {
    // Simple quadratic ease-in: i^2 scaled to 8-bit range
    // Using scale8(i, i) which computes (i * i) / 255
    return scale8(i, i);
}

u8 easeInOutQuad8(u8 i) {
    constexpr u16 MAX = 0xFF;            // 255
    constexpr u16 HALF = (MAX + 1) >> 1; // 128
    constexpr u16 DENOM = MAX;           // divisor for scaling
    constexpr u16 ROUND = DENOM >> 1;    // for rounding

    if (i < HALF) {
        // first half: y = 2·(i/MAX)² → y_i = 2·i² / MAX
        u32 t = i;
        u32 num = 2 * t * t + ROUND; // 2*i², +half for rounding
        return u8(num / DENOM);
    } else {
        // second half: y = 1 − 2·(1−i/MAX)²
        // → y_i = MAX − (2·(MAX−i)² / MAX)
        u32 d = MAX - i;
        u32 num = 2 * d * d + ROUND; // 2*(MAX−i)², +half for rounding
        return u8(MAX - (num / DENOM));
    }
}

u8 easeInOutCubic8(u8 i) {
    constexpr u16 MAX = 0xFF;                  // 255
    constexpr u16 HALF = (MAX + 1) >> 1;       // 128
    constexpr u32 DENOM = (u32)MAX * MAX; // 255*255 = 65025
    constexpr u32 ROUND = DENOM >> 1;          // for rounding

    if (i < HALF) {
        // first half: y = 4·(i/MAX)³ → y_i = 4·i³ / MAX²
        u32 ii = i;
        u32 cube = ii * ii * ii;    // i³
        u32 num = 4 * cube + ROUND; // 4·i³, +half denom for rounding
        return u8(num / DENOM);
    } else {
        // second half: y = 1 − ((−2·t+2)³)/2
        // where t = i/MAX; equivalently:
        // y_i = MAX − (4·(MAX−i)³ / MAX²)
        u32 d = MAX - i;
        u32 cube = d * d * d; // (MAX−i)³
        u32 num = 4 * cube + ROUND;
        return u8(MAX - (num / DENOM));
    }
}

u8 easeOutQuad8(u8 i) {
    // ease-out is the inverse of ease-in: 1 - (1-t)²
    // For 8-bit: y = MAX - (MAX-i)² / MAX
    constexpr u16 MAX = 0xFF;
    u32 d = MAX - i;              // (MAX - i)
    u32 num = d * d + (MAX >> 1); // (MAX-i)² + rounding
    return u8(MAX - (num / MAX));
}

u8 easeInCubic8(u8 i) {
    // Simple cubic ease-in: i³ scaled to 8-bit range
    // y = i³ / MAX²
    constexpr u16 MAX = 0xFF;
    constexpr u32 DENOM = (u32)MAX * MAX;
    constexpr u32 ROUND = DENOM >> 1;

    u32 ii = i;
    u32 cube = ii * ii * ii; // i³
    u32 num = cube + ROUND;
    return u8(num / DENOM);
}

u8 easeOutCubic8(u8 i) {
    // ease-out cubic: 1 - (1-t)³
    // For 8-bit: y = MAX - (MAX-i)³ / MAX²
    constexpr u16 MAX = 0xFF;
    constexpr u32 DENOM = (u32)MAX * MAX;
    constexpr u32 ROUND = DENOM >> 1;

    u32 d = MAX - i;      // (MAX - i)
    u32 cube = d * d * d; // (MAX-i)³
    u32 num = cube + ROUND;
    return u8(MAX - (num / DENOM));
}

u8 easeInSine8(u8 i) {

    static const u8 easeInSineTable[256] FL_PROGMEM = {
        0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   0,   1,   1,   1,
        1,   1,   1,   1,   2,   2,   2,   2,   2,   3,   3,   3,   3,   4,
        4,   4,   4,   5,   5,   5,   6,   6,   6,   7,   7,   7,   8,   8,
        8,   9,   9,   10,  10,  11,  11,  12,  12,  12,  13,  13,  14,  14,
        15,  16,  16,  17,  17,  18,  18,  19,  20,  20,  21,  21,  22,  23,
        23,  24,  25,  25,  26,  27,  27,  28,  29,  30,  30,  31,  32,  33,
        33,  34,  35,  36,  37,  37,  38,  39,  40,  41,  42,  42,  43,  44,
        45,  46,  47,  48,  49,  50,  51,  52,  52,  53,  54,  55,  56,  57,
        58,  59,  60,  61,  62,  63,  64,  65,  67,  68,  69,  70,  71,  72,
        73,  74,  75,  76,  77,  79,  80,  81,  82,  83,  84,  86,  87,  88,
        89,  90,  91,  93,  94,  95,  96,  98,  99,  100, 101, 103, 104, 105,
        106, 108, 109, 110, 112, 113, 114, 115, 117, 118, 119, 121, 122, 123,
        125, 126, 127, 129, 130, 132, 133, 134, 136, 137, 139, 140, 141, 143,
        144, 146, 147, 148, 150, 151, 153, 154, 156, 157, 159, 160, 161, 163,
        164, 166, 167, 169, 170, 172, 173, 175, 176, 178, 179, 181, 182, 184,
        185, 187, 188, 190, 191, 193, 194, 196, 197, 199, 200, 202, 204, 205,
        207, 208, 210, 211, 213, 214, 216, 217, 219, 221, 222, 224, 225, 227,
        228, 230, 231, 233, 235, 236, 238, 239, 241, 242, 244, 246, 247, 249,
        250, 252, 253, 255};

    // ease-in sine: 1 - cos(t * π/2)
    // Handle boundary conditions explicitly
    return FL_PGM_READ_BYTE_NEAR(&easeInSineTable[i]);
}

u8 easeOutSine8(u8 i) {
    // ease-out sine: sin(t * π/2)
    // Delegate to 16-bit version for consistency and accuracy
    // Scale 8-bit input to 16-bit range, call 16-bit function, scale result back
    u16 input16 = map8_to_16(i);
    u16 result16 = easeOutSine16(input16);
    return map16_to_8(result16);
}

u8 easeInOutSine8(u8 i) {
    // ease-in-out sine: -(cos(π*t) - 1) / 2
    // Delegate to 16-bit version for consistency and accuracy
    // Scale 8-bit input to 16-bit range, call 16-bit function, scale result back
    u16 input16 = map8_to_16(i);
    u16 result16 = easeInOutSine16(input16);
    return map16_to_8(result16);
}

// 16-bit easing functions
u16 easeInQuad16(u16 i) {
    // Simple quadratic ease-in: i^2 scaled to 16-bit range
    // Using scale16(i, i) which computes (i * i) / 65535
    return scale16(i, i);
}

u16 easeInOutQuad16(u16 x) {
    // 16-bit quadratic ease-in / ease-out function
    constexpr u32 MAX = 0xFFFF;          // 65535
    constexpr u32 HALF = (MAX + 1) >> 1; // 32768
    constexpr u32 DENOM = MAX;           // divisor
    constexpr u32 ROUND = DENOM >> 1;    // for rounding

    if (x < HALF) {
        // first half: y = 2·(x/MAX)² → y_i = 2·x² / MAX
        fl::u64 xi = x;
        fl::u64 num = 2 * xi * xi + ROUND; // 2*x², +half for rounding
        return u16(num / DENOM);
    } else {
        // second half: y = 1 − 2·(1−x/MAX)² → y_i = MAX − (2·(MAX−x)² / MAX)
        fl::u64 d = MAX - x;
        fl::u64 num = 2 * d * d + ROUND; // 2*(MAX−x)², +half for rounding
        return u16(MAX - (num / DENOM));
    }
}

u16 easeInOutCubic16(u16 x) {
    const u32 MAX = 0xFFFF;             // 65535
    const u32 HALF = (MAX + 1) >> 1;    // 32768
    const fl::u64 M2 = (fl::u64)MAX * MAX; // 65535² = 4 294 836 225

    if (x < HALF) {
        // first half:  y = 4·(x/MAX)³  →  y_i = 4·x³ / MAX²
        fl::u64 xi = x;
        fl::u64 cube = xi * xi * xi; // x³
        // add M2/2 for rounding
        fl::u64 num = 4 * cube + (M2 >> 1);
        return (u16)(num / M2);
    } else {
        // second half: y = 1 − ((2·(1−x/MAX))³)/2
        // → y_i = MAX − (4·(MAX−x)³ / MAX²)
        fl::u64 d = MAX - x;
        fl::u64 cube = d * d * d; // (MAX−x)³
        fl::u64 num = 4 * cube + (M2 >> 1);
        return (u16)(MAX - (num / M2));
    }
}

u16 easeOutQuad16(u16 i) {
    // ease-out quadratic: 1 - (1-t)²
    // For 16-bit: y = MAX - (MAX-i)² / MAX
    constexpr u32 MAX = 0xFFFF;     // 65535
    constexpr u32 ROUND = MAX >> 1; // for rounding

    fl::u64 d = MAX - i;         // (MAX - i)
    fl::u64 num = d * d + ROUND; // (MAX-i)² + rounding
    return u16(MAX - (num / MAX));
}

u16 easeInCubic16(u16 i) {
    // Simple cubic ease-in: i³ scaled to 16-bit range
    // y = i³ / MAX²
    constexpr u32 MAX = 0xFFFF;                // 65535
    constexpr fl::u64 DENOM = (fl::u64)MAX * MAX; // 65535²
    constexpr fl::u64 ROUND = DENOM >> 1;          // for rounding

    fl::u64 ii = i;
    fl::u64 cube = ii * ii * ii; // i³
    fl::u64 num = cube + ROUND;
    return u16(num / DENOM);
}

u16 easeOutCubic16(u16 i) {
    // ease-out cubic: 1 - (1-t)³
    // For 16-bit: y = MAX - (MAX-i)³ / MAX²
    constexpr u32 MAX = 0xFFFF;                // 65535
    constexpr fl::u64 DENOM = (fl::u64)MAX * MAX; // 65535²
    constexpr fl::u64 ROUND = DENOM >> 1;          // for rounding

    fl::u64 d = MAX - i;      // (MAX - i)
    fl::u64 cube = d * d * d; // (MAX-i)³
    fl::u64 num = cube + ROUND;
    return u16(MAX - (num / DENOM));
}

u16 easeInSine16(u16 i) {
    // ease-in sine: 1 - cos(t * π/2)
    // Handle boundary conditions explicitly
    if (i == 0)
        return 0;
    // Remove the hard-coded boundary for 65535 and let math handle it
        
    // For 16-bit: use cos32 for efficiency and accuracy
    // Map i from [0,65535] to [0,4194304] in cos32 space (zero to quarter wave)
    // Formula: 1 - cos(t * π/2) where t goes from 0 to 1
    // sin32/cos32 quarter cycle is 16777216/4 = 4194304
    u32 angle = ((fl::u64)i * 4194304ULL) / 65535ULL;
    i32 cos_result = fl::cos32(angle);

    // Convert cos32 output and apply easing formula: 1 - cos(t * π/2)
    // cos32 output range is [-2147418112, 2147418112]
    // At t=0: cos(0) = 2147418112, result should be 0
    // At t=1: cos(π/2) = 0, result should be 65535
    
    const fl::i64 MAX_COS32 = 2147418112LL;
    
    // Calculate: (MAX_COS32 - cos_result) and scale to [0, 65535]
    fl::i64 adjusted = MAX_COS32 - (fl::i64)cos_result;
    
    // Scale from [0, 2147418112] to [0, 65535]
    fl::u64 result = (fl::u64)adjusted * 65535ULL + (MAX_COS32 >> 1); // Add half for rounding
    u16 final_result = (u16)(result / (fl::u64)MAX_COS32);
    
    return final_result;
}

u16 easeOutSine16(u16 i) {
    // ease-out sine: sin(t * π/2)
    // Handle boundary conditions explicitly
    if (i == 0)
        return 0;
    if (i == 65535)
        return 65535;

    // For 16-bit: use sin32 for efficiency and accuracy
    // Map i from [0,65535] to [0,4194304] in sin32 space (zero to quarter wave)
    // Formula: sin(t * π/2) where t goes from 0 to 1
    // sin32 quarter cycle is 16777216/4 = 4194304
    u32 angle = ((fl::u64)i * 4194304ULL) / 65535ULL;
    i32 sin_result = fl::sin32(angle);
    
    // Convert sin32 output range [-2147418112, 2147418112] to [0, 65535]
    // sin32 output is in range -32767*65536 to +32767*65536
    // For ease-out sine, we only use positive portion [0, 2147418112] -> [0, 65535]
    return (u16)((fl::u64)sin_result * 65535ULL / 2147418112ULL);
}

u16 easeInOutSine16(u16 i) {
    // ease-in-out sine: -(cos(π*t) - 1) / 2
    // Handle boundary conditions explicitly
    if (i == 0)
        return 0;
    if (i == 65535)
        return 65535;

    // For 16-bit: use cos32 for efficiency and accuracy
    // Map i from [0,65535] to [0,8388608] in cos32 space (0 to half wave)
    // Formula: (1 - cos(π*t)) / 2 where t goes from 0 to 1
    // sin32/cos32 half cycle is 16777216/2 = 8388608
    u32 angle = ((fl::u64)i * 8388608ULL) / 65535ULL;
    i32 cos_result = fl::cos32(angle);

    // Convert cos32 output and apply easing formula: (1 - cos(π*t)) / 2
    // cos32 output range is [-2147418112, 2147418112]
    // We want: (2147418112 - cos_result) / 2, then scale to [0, 65535]
    fl::i64 adjusted = (2147418112LL - (fl::i64)cos_result) / 2;
    return (u16)((fl::u64)adjusted * 65535ULL / 2147418112ULL);
}

// --- Gamma8 implementation ---

// Fixed-point key for gamma cache: unsigned 4.12 (range [0, 15.999], 1/4096 resolution).
using GammaKey = fl::ufixed_point<4, 12>;

class Gamma8Impl : public Gamma8 {
public:
    explicit Gamma8Impl(float gamma) {
        // i=0 is mathematically exact regardless of gamma: pow(0, any) = 0.
        // The s8x24::pow short-circuit covers exact 0 input, but we set it
        // here directly anyway because (a) the loop below avoids the
        // round-trip math and (b) `mLut` is otherwise uninitialized.
        // i=255 used to need a workaround for the log2(1+t) endpoint
        // residual; that snap is now handled inside s8x24::pow itself
        // (see #2969).
        mLut[0] = 0;
        // Compute the 256-entry u16 gamma LUT in fixed-point so we don't
        // pull `__ieee754_pow` (libm, ~2.7 KB) into release builds — the
        // double-precision pow chain dominates the top-9 bytes attributed
        // to libm in the post-#2908 ESP32-S3 NEOPIXEL Blink audit
        // (see #2886 / #2910).
        //
        // `s8x24` is 8-integer + 24-fractional bits (same 32-bit storage as
        // s16x16, but 256× the sub-LSB resolution). Both log2_fp/exp2_fp use
        // the same 4-term minimax polynomial, but s8x24 carries the full
        // 24-bit intermediate precision end-to-end instead of truncating
        // back to 16 frac bits at each stage — bringing the runtime output
        // within ~1 LSB of true float pow() at the u16 output. Combined
        // with the special-case for gamma=2.8 in Gamma8::getOrCreate(),
        // this closes the divergence with the precomputed GAMMA_2_8_LUT
        // (see #2963 audit + ucs7604 "default gamma 2.8" subcase).
        //
        // Bit budget check: max intermediate is exp*log2_fp(1/255) =
        // 16 * -7.994 ≈ -127.9, fits in s8x24's signed [-128, 128) range.
        // (GammaKey caps user gamma at 16.)
        //
        // libm-free: log2_fp/exp2_fp are pure integer-polynomial impls.
        // The only float kept is the one-shot s8x24(gamma) constructor
        // call (pulls __mulsf3 / __fixsfsi helpers, both << 100 B).
        const fl::s8x24 gamma_fp(gamma);
        constexpr fl::s8x24 inv_255_fp(1.0f / 255.0f);
        for (int i = 1; i < 256; ++i) {
            const fl::s8x24 x = static_cast<i32>(i) * inv_255_fp;  // (0, 1]
            const fl::s8x24 r = fl::s8x24::pow(x, gamma_fp);       // (0, 1]
            // r.raw() is the s8x24 raw with FRAC_BITS=24, range [0, 2^24].
            // Scale to u16 [0, 65535] with round-half-up. 24+16 bit
            // multiplication needs a u64 intermediate:
            //   result = ((u64)raw * 65535 + (1<<23)) >> 24
            const fl::u64 scaled =
                (static_cast<fl::u64>(static_cast<fl::u32>(r.raw())) * 65535ull
                 + (1ull << 23)) >> 24;
            mLut[i] = static_cast<u16>(scaled > 65535u ? 65535u : scaled);
        }
    }

    // Construct a Gamma8Impl by copying a precomputed PROGMEM LUT directly
    // into mLut. Used by Gamma8::getOrCreate(2.8f) to alias GAMMA_2_8_LUT,
    // guaranteeing bit-identical output with fl::gamma_2_8() and skipping
    // the runtime pow chain entirely. The `from_progmem_lut_tag` tag
    // disambiguates from the float ctor; we can't add a ctor taking just
    // `const u16*` because that would silently bind to integer-literal
    // gammas (`Gamma8::getOrCreate(2)` would be a footgun).
    struct from_progmem_lut_tag {};
    Gamma8Impl(from_progmem_lut_tag, const u16* progmem_lut) {
        for (int i = 0; i < 256; ++i) {
            mLut[i] = FL_PGM_READ_WORD_ALIGNED(&progmem_lut[i]);
        }
    }

    void convert(fl::span<const u8> input,
                  fl::span<u16> output) const override {
        const int n =
            input.size() < output.size() ? input.size() : output.size();
        for (int i = 0; i < n; ++i) {
            output[i] = mLut[input[i]];
        }
    }

    void convert(fl::span<const fl::ufixed_point<8, 8>> input,
                  fl::span<u16> output) const override {
        const int n =
            input.size() < output.size() ? input.size() : output.size();
        for (int i = 0; i < n; ++i) {
            output[i] = lerpLut(input[i]);
        }
    }

    void convert(fl::span<const fl::ufixed_point<8, 8>> input,
                  fl::span<fl::ufixed_point<8, 8>> output) const override {
        using FP = fl::ufixed_point<8, 8>;
        const int n =
            input.size() < output.size() ? input.size() : output.size();
        for (int i = 0; i < n; ++i) {
            output[i] = FP::from_raw(lerpLut(input[i]));
        }
    }

private:
    u16 lerpLut(const fl::ufixed_point<8, 8>& fp) const {
        u16 raw = fp.raw();
        u8 idx = static_cast<u8>(raw >> 8);    // integer part (LUT index)
        u8 frac = static_cast<u8>(raw & 0xFF); // fractional part (0-255)
        u16 a = mLut[idx];
        u16 b = (idx < 255) ? mLut[idx + 1] : mLut[idx];
        i32 diff = static_cast<i32>(b) - static_cast<i32>(a);
        return static_cast<u16>(a + ((diff * frac) >> 8));
    }

    FL_ALIGNAS(64) u16 mLut[256];
};

fl::shared_ptr<const Gamma8> Gamma8::getOrCreate(float gamma) {
    // Clamp gamma to the cache-key domain BEFORE deriving the key and
    // before feeding it to s8x24. `GammaKey` is `ufixed_point<4, 12>`,
    // so its representable range is [0, ~16). A negative input would
    // wrap unsigned in the GammaKey ctor; an input ≥ 16 would either
    // saturate or wrap. Clamping up-front also keeps the s8x24 bit
    // budget safe: `exp * log2_fp(1/255) ≈ -7.994 * gamma` must stay
    // within s8x24's signed [-128, 128) integer range — gamma ≤ 16
    // gives -127.9, right at the edge.
    constexpr float kGammaMax = 15.99975f;  // just under 4.12 ufixed max
    const float gamma_clamped = (gamma < 0.0f)         ? 0.0f
                                : (gamma > kGammaMax)  ? kGammaMax
                                                       : gamma;
    GammaKey key(gamma_clamped);

    // Fast path: gamma 2.8 is the canonical default (every legacy
    // `addLeds<NEOPIXEL>` flow, every WS281x chipset). We have a
    // precomputed `GAMMA_2_8_LUT` already in .rodata (used by the
    // free function `fl::gamma_2_8`). Construct a Gamma8 instance
    // that copies that table directly instead of recomputing 256
    // fixed-point pow operations (which would diverge from the
    // precomputed values by ~30-50 LSB). Cached via weak_ptr to
    // match the documented `getOrCreate` lifetime contract — the
    // instance disappears once all callers drop their shared_ptr,
    // and the next call rebuilds it.
    constexpr GammaKey k28(2.8f);
    if (key == k28) {
        static fl::weak_ptr<const Gamma8> s28_weak;
        if (auto existing = s28_weak.lock()) {
            return existing;
        }
        auto ptr = fl::make_shared<Gamma8Impl>(
            Gamma8Impl::from_progmem_lut_tag{}, GAMMA_2_8_LUT);
        s28_weak = ptr;
        return ptr;
    }

    // Slow path: arbitrary gamma. Single-entry weak_ptr cache. Users
    // who switch gamma values at runtime pay one extra `Gamma8Impl`
    // reconstruction per switch (256 fixed-point pow operations +
    // ~512 B alloc). See #2928 / #2886.
    static GammaKey sCachedKey{};
    static fl::weak_ptr<const Gamma8> sCachedPtr;
    static bool sCacheValid = false;

    if (sCacheValid && key == sCachedKey) {
        if (auto existing = sCachedPtr.lock()) {
            return existing;
        }
    }

    // Use the clamped gamma so the LUT we build matches the cache key
    // (otherwise out-of-range inputs would keep missing the cache).
    fl::shared_ptr<const Gamma8> ptr =
        fl::make_shared<Gamma8Impl>(gamma_clamped);
    sCachedKey = key;
    sCachedPtr = ptr;
    sCacheValid = true;
    return ptr;
}

} // namespace fl
