libmonet/include/monet/color.h

#pragma once

// Colour spaces and tone solving. No dependencies beyond the standard library.
//
// Colours are handled in OKLab/OKLCH (Björn Ottosson, 2020) and a "tone" is the
// CIE L* of the result, like Material 3's HCT tone. This is a light weight
// stand in for HCT: it keeps the idea (a hue, a chroma and a tone that really
// is perceptual lightness, so the same tone gives the same contrast whatever
// the hue) without the CAM16 machinery. It is not bit for bit what Material's
// own libraries produce.

#include <algorithm>
#include <cmath>
#include <cstdint>
#include <string>

namespace monet {

struct Rgb {
    std::uint8_t r = 0;
    std::uint8_t g = 0;
    std::uint8_t b = 0;
    friend bool operator==(const Rgb &, const Rgb &) = default;
};

struct Oklch {
    double L = 0;  // 0..1
    double C = 0;  // chroma, 0 for greys, about 0.37 at most in sRGB
    double h = 0;  // hue in degrees, [0, 360)
};

namespace detail {

struct Linear {
    double r, g, b;
};

constexpr double kPi = 3.14159265358979323846;

inline double toLinear(double c) {
    return c <= 0.04045 ? c / 12.92 : std::pow((c + 0.055) / 1.055, 2.4);
}
inline double fromLinear(double c) {
    c = std::clamp(c, 0.0, 1.0);
    return c <= 0.0031308 ? c * 12.92 : 1.055 * std::pow(c, 1.0 / 2.4) - 0.055;
}

inline Linear linearOf(Rgb c) {
    return {toLinear(c.r / 255.0), toLinear(c.g / 255.0), toLinear(c.b / 255.0)};
}

inline Rgb rgbOf(Linear c) {
    auto q = [](double v) {
        return static_cast<std::uint8_t>(std::lround(std::clamp(fromLinear(v), 0.0, 1.0) * 255.0));
    };
    return {q(c.r), q(c.g), q(c.b)};
}

inline double luminanceOf(Linear c) { return 0.2126729 * c.r + 0.7151522 * c.g + 0.0721750 * c.b; }

inline Linear linearFromOklab(double L, double a, double b) {
    const double l_ = L + 0.3963377774 * a + 0.2158037573 * b;
    const double m_ = L - 0.1055613458 * a - 0.0638541728 * b;
    const double s_ = L - 0.0894841775 * a - 1.2914855480 * b;
    const double l = l_ * l_ * l_, m = m_ * m_ * m_, s = s_ * s_ * s_;
    return {+4.0767416621 * l - 3.3077115913 * m + 0.2309699292 * s,
            -1.2684380046 * l + 2.6097574011 * m - 0.3413193965 * s,
            -0.0041960863 * l - 0.7034186147 * m + 1.7076147010 * s};
}

inline Linear linearFromOklch(Oklch c) {
    const double rad = c.h * kPi / 180.0;
    return linearFromOklab(c.L, c.C * std::cos(rad), c.C * std::sin(rad));
}

inline bool inGamut(Linear c) {
    constexpr double eps = 1e-5;
    return c.r >= -eps && c.r <= 1 + eps && c.g >= -eps && c.g <= 1 + eps && c.b >= -eps &&
           c.b <= 1 + eps;
}

// Lowers the chroma until the colour fits in sRGB, keeping lightness and hue.
inline Linear gamutMapped(double L, double C, double h) {
    Linear full = linearFromOklch({L, C, h});
    if (inGamut(full)) {
        return full;
    }
    double lo = 0, hi = C;
    for (int i = 0; i < 28; i++) {
        const double mid = (lo + hi) / 2;
        if (inGamut(linearFromOklch({L, mid, h}))) {
            lo = mid;
        } else {
            hi = mid;
        }
    }
    return linearFromOklch({L, lo, h});
}

} // namespace detail

inline Oklch toOklch(Rgb c) {
    const auto lin = detail::linearOf(c);
    const double l = std::cbrt(0.4122214708 * lin.r + 0.5363325363 * lin.g + 0.0514459929 * lin.b);
    const double m = std::cbrt(0.2119034982 * lin.r + 0.6806995451 * lin.g + 0.1073969566 * lin.b);
    const double s = std::cbrt(0.0883024619 * lin.r + 0.2817188376 * lin.g + 0.6299787005 * lin.b);
    const double L = 0.2104542553 * l + 0.7936177850 * m - 0.0040720468 * s;
    const double a = 1.9779984951 * l - 2.4285922050 * m + 0.4505937099 * s;
    const double b = 0.0259040371 * l + 0.7827717662 * m - 0.8086757660 * s;
    double h = std::atan2(b, a) * 180.0 / detail::kPi;
    if (h < 0) {
        h += 360.0;
    }
    return {L, std::hypot(a, b), h};
}

inline Rgb toRgb(Oklch c) { return detail::rgbOf(detail::gamutMapped(c.L, c.C, c.h)); }

// WCAG relative luminance, 0..1
inline double luminance(Rgb c) { return detail::luminanceOf(detail::linearOf(c)); }

// CIE L*, 0..100
inline double toneOf(Rgb c) {
    const double y = luminance(c);
    return y > 216.0 / 24389.0 ? 116.0 * std::cbrt(y) - 16.0 : y * 24389.0 / 27.0;
}

inline double luminanceOfTone(double tone) {
    tone = std::clamp(tone, 0.0, 100.0);
    return tone > 8.0 ? std::pow((tone + 16.0) / 116.0, 3.0) : tone * 27.0 / 24389.0;
}

// WCAG contrast ratio, 1..21
inline double contrast(Rgb a, Rgb b) {
    const double la = luminance(a), lb = luminance(b);
    return (std::max(la, lb) + 0.05) / (std::min(la, lb) + 0.05);
}

// The colour with the given hue, at most the given chroma, and CIE L* equal to
// `tone`. If the chroma is not available at that tone it is reduced, which is
// what makes the extremes (tone near 0 or 100) grey.
inline Rgb toneColor(double hue, double chroma, double tone) {
    if (tone <= 0.0) {
        return {0, 0, 0};
    }
    if (tone >= 100.0) {
        return {255, 255, 255};
    }
    const double target = luminanceOfTone(tone);
    double lo = 0.0, hi = 1.0;
    for (int i = 0; i < 42; i++) {
        const double mid = (lo + hi) / 2;
        if (detail::luminanceOf(detail::gamutMapped(mid, chroma, hue)) < target) {
            lo = mid;
        } else {
            hi = mid;
        }
    }
    return detail::rgbOf(detail::gamutMapped((lo + hi) / 2, chroma, hue));
}

inline std::string hex(Rgb c) {
    static const char *digits = "0123456789ABCDEF";
    std::string out = "#";
    for (std::uint8_t v : {c.r, c.g, c.b}) {
        out += digits[v >> 4];
        out += digits[v & 15];
    }
    return out;
}

} // namespace monet