add colorcvt module

docs/modules/colorcvt.md

@@ -0,0 +1,10 @@
+---
+icon: package
+label: colorcvt
+---
+
+Provide color conversion functions.
+
+#### Source code
+
+:::code source="../../include/typings/colorcvt.pyi" :::

include/pocketpy/config.h

@@ -39,6 +39,12 @@
 // This is the maximum character length of a module path
 #define PK_MAX_MODULE_PATH_LEN      63
 
+// This is some math constants
+#define PK_M_PI                     3.1415926535897932384
+#define PK_M_E                      2.7182818284590452354
+#define PK_M_DEG2RAD                0.017453292519943295
+#define PK_M_RAD2DEG                57.29577951308232
+
 #ifdef _WIN32
     #define PK_PLATFORM_SEP '\\'
 #else

include/pocketpy/interpreter/modules.h

@@ -18,6 +18,7 @@ void pk__add_module_pickle();
 
 void pk__add_module_linalg();
 void pk__add_module_array2d();
+void pk__add_module_colorcvt();
 
 void pk__add_module_conio();
 void pk__add_module_lz4();

include/typings/colorcvt.pyi

@@ -0,0 +1,8 @@
+from linalg import vec3
+
+def linear_srgb_to_srgb(rgb: vec3) -> vec3: ...
+def srgb_to_linear_srgb(rgb: vec3) -> vec3: ...
+def srgb_to_hsv(rgb: vec3) -> vec3: ...
+def hsv_to_srgb(hsv: vec3) -> vec3: ...
+def oklch_to_linear_srgb(lch: vec3) -> vec3: ...
+def linear_srgb_to_oklch(rgb: vec3) -> vec3: ...

src/interpreter/vm.c

@@ -201,6 +201,7 @@ void VM__ctor(VM* self) {
 
     pk__add_module_linalg();
     pk__add_module_array2d();
+    pk__add_module_colorcvt();
 
     // add modules
     pk__add_module_os();

src/modules/colorcvt.c

@@ -0,0 +1,230 @@
+#include "pocketpy/pocketpy.h"
+
+#include "pocketpy/common/utils.h"
+#include "pocketpy/objects/object.h"
+#include "pocketpy/common/sstream.h"
+#include "pocketpy/interpreter/vm.h"
+#include <math.h>
+
+// https://bottosson.github.io/posts/gamutclipping/#oklab-to-linear-srgb-conversion
+
+// clang-format off
+static c11_vec3 linear_srgb_to_oklab(c11_vec3 c)
+{
+	float l = 0.4122214708f * c.x + 0.5363325363f * c.y + 0.0514459929f * c.z;
+	float m = 0.2119034982f * c.x + 0.6806995451f * c.y + 0.1073969566f * c.z;
+	float s = 0.0883024619f * c.x + 0.2817188376f * c.y + 0.6299787005f * c.z;
+
+	float l_ = cbrtf(l);
+	float m_ = cbrtf(m);
+	float s_ = cbrtf(s);
+
+	return (c11_vec3){{
+		0.2104542553f * l_ + 0.7936177850f * m_ - 0.0040720468f * s_,
+		1.9779984951f * l_ - 2.4285922050f * m_ + 0.4505937099f * s_,
+		0.0259040371f * l_ + 0.7827717662f * m_ - 0.8086757660f * s_,
+	}};
+}
+
+static c11_vec3 oklab_to_linear_srgb(c11_vec3 c)
+{
+    float l_ = c.x + 0.3963377774f * c.y + 0.2158037573f * c.z;
+    float m_ = c.x - 0.1055613458f * c.y - 0.0638541728f * c.z;
+    float s_ = c.x - 0.0894841775f * c.y - 1.2914855480f * c.z;
+
+    float l = l_ * l_ * l_;
+    float m = m_ * m_ * m_;
+    float s = s_ * s_ * s_;
+
+    return (c11_vec3){{
+        +4.0767416621f * l - 3.3077115913f * m + 0.2309699292f * s,
+        -1.2684380046f * l + 2.6097574011f * m - 0.3413193965f * s,
+        -0.0041960863f * l - 0.7034186147f * m + 1.7076147010f * s,
+    }};
+}
+
+// clang-format on
+
+static float _gamma_correct_inv(float x) {
+    return (x <= 0.04045f) ? (x / 12.92f) : powf((x + 0.055f) / 1.055f, 2.4f);
+}
+
+static float _gamma_correct(float x) {
+    return (x <= 0.0031308f) ? (12.92f * x) : (1.055f * powf(x, 1.0f / 2.4f) - 0.055f);
+}
+
+static c11_vec3 srgb_to_linear_srgb(c11_vec3 c) {
+    c.x = _gamma_correct_inv(c.x);
+    c.y = _gamma_correct_inv(c.y);
+    c.z = _gamma_correct_inv(c.z);
+    return c;
+}
+
+static c11_vec3 linear_srgb_to_srgb(c11_vec3 c) {
+    c.x = _gamma_correct(c.x);
+    c.y = _gamma_correct(c.y);
+    c.z = _gamma_correct(c.z);
+    return c;
+}
+
+static c11_vec3 _oklab_to_oklch(c11_vec3 c) {
+    c11_vec3 res;
+    res.x = c.x;
+    res.y = sqrtf(c.y * c.y + c.z * c.z);
+    res.z = fmodf(atan2f(c.z, c.y), 2 * (float)PK_M_PI);
+    res.z = res.z * PK_M_RAD2DEG;
+    return res;
+}
+
+static c11_vec3 _oklch_to_oklab(c11_vec3 c) {
+    c11_vec3 res;
+    res.x = c.x;
+    res.y = c.y * cosf(c.z * PK_M_DEG2RAD);
+    res.z = c.y * sinf(c.z * PK_M_DEG2RAD);
+    return res;
+}
+
+static c11_vec3 linear_srgb_to_oklch(c11_vec3 c) {
+    return _oklab_to_oklch(linear_srgb_to_oklab(c));
+}
+
+static bool _is_valid_srgb(c11_vec3 c) {
+    return c.x >= 0.0f && c.x <= 1.0f && c.y >= 0.0f && c.y <= 1.0f && c.z >= 0.0f && c.z <= 1.0f;
+}
+
+static c11_vec3 oklch_to_linear_srgb(c11_vec3 c) {
+    c11_vec3 candidate = oklab_to_linear_srgb(_oklch_to_oklab(c));
+    if(_is_valid_srgb(candidate)) return candidate;
+
+    // try with chroma = 0
+    c11_vec3 clamped = {
+        {c.x, 0.0f, c.z}
+    };
+
+    // if not even chroma = 0 is displayable
+    // fall back to RGB clamping
+    candidate = oklab_to_linear_srgb(_oklch_to_oklab(clamped));
+    if(!_is_valid_srgb(candidate)) {
+        candidate.x = fmaxf(0.0f, fminf(1.0f, candidate.x));
+        candidate.y = fmaxf(0.0f, fminf(1.0f, candidate.y));
+        candidate.z = fmaxf(0.0f, fminf(1.0f, candidate.z));
+        return candidate;
+    }
+
+    // By this time we know chroma = 0 is displayable and our current chroma is not.
+    // Find the displayable chroma through the bisection method.
+    float start = 0.0f;
+    float end = c.y;
+    float range[2] = {0.0f, 0.4f};
+    float resolution = (range[1] - range[0]) / powf(2, 13);
+    float _last_good_c = clamped.y;
+
+    while(end - start > resolution) {
+        clamped.y = start + (end - start) * 0.5f;
+        candidate = oklab_to_linear_srgb(_oklch_to_oklab(clamped));
+        if(_is_valid_srgb(candidate)) {
+            _last_good_c = clamped.y;
+            start = clamped.y;
+        } else {
+            end = clamped.y;
+        }
+    }
+
+    candidate = oklab_to_linear_srgb(_oklch_to_oklab(clamped));
+    if(_is_valid_srgb(candidate)) return candidate;
+    clamped.y = _last_good_c;
+    return oklab_to_linear_srgb(_oklch_to_oklab(clamped));
+}
+
+// https://github.com/python/cpython/blob/3.13/Lib/colorsys.py
+static c11_vec3 srgb_to_hsv(c11_vec3 c) {
+    float r = c.x;
+    float g = c.y;
+    float b = c.z;
+
+    float maxc = fmaxf(r, fmaxf(g, b));
+    float minc = fminf(r, fminf(g, b));
+    float v = maxc;
+    if(minc == maxc) {
+        return (c11_vec3){
+            {0.0f, 0.0f, v}
+        };
+    }
+
+    float s = (maxc - minc) / maxc;
+    float rc = (maxc - r) / (maxc - minc);
+    float gc = (maxc - g) / (maxc - minc);
+    float bc = (maxc - b) / (maxc - minc);
+    float h;
+    if(r == maxc) {
+        h = bc - gc;
+    } else if(g == maxc) {
+        h = 2.0f + rc - bc;
+    } else {
+        h = 4.0f + gc - rc;
+    }
+    h = fmodf(h / 6.0f, 1.0f);
+    return (c11_vec3){
+        {h, s, v}
+    };
+}
+
+static c11_vec3 hsv_to_srgb(c11_vec3 c) {
+    float h = c.x;
+    float s = c.y;
+    float v = c.z;
+
+    if(s == 0.0f) {
+        return (c11_vec3){
+            {v, v, v}
+        };
+    }
+
+    int i = (int)(h * 6.0f);
+    float f = (h * 6.0f) - i;
+    float p = v * (1.0f - s);
+    float q = v * (1.0f - s * f);
+    float t = v * (1.0f - s * (1.0f - f));
+    i = i % 6;
+    switch(i) {
+        // clang-format off
+        case 0: return (c11_vec3){{v, t, p}};
+        case 1: return (c11_vec3){{q, v, p}};
+        case 2: return (c11_vec3){{p, v, t}};
+        case 3: return (c11_vec3){{p, q, v}};
+        case 4: return (c11_vec3){{t, p, v}};
+        case 5: return (c11_vec3){{v, p, q}};
+        // clang-format on
+        default: c11__unreachable();
+    }
+}
+
+#define DEF_VEC3_WRAPPER(F)                                                                        \
+    static bool colorcvt_##F(int argc, py_Ref argv);                                               \
+    static bool colorcvt_##F(int argc, py_Ref argv) {                                              \
+        PY_CHECK_ARGC(1);                                                                          \
+        PY_CHECK_ARG_TYPE(0, tp_vec3);                                                             \
+        c11_vec3 c = py_tovec3(argv);                                                              \
+        py_newvec3(py_retval(), F(c));                                                             \
+        return true;                                                                               \
+    }
+
+DEF_VEC3_WRAPPER(linear_srgb_to_srgb)
+DEF_VEC3_WRAPPER(srgb_to_linear_srgb)
+DEF_VEC3_WRAPPER(srgb_to_hsv)
+DEF_VEC3_WRAPPER(hsv_to_srgb)
+DEF_VEC3_WRAPPER(oklch_to_linear_srgb)
+DEF_VEC3_WRAPPER(linear_srgb_to_oklch)
+
+void pk__add_module_colorcvt() {
+    py_Ref mod = py_newmodule("colorcvt");
+
+    py_bindfunc(mod, "linear_srgb_to_srgb", colorcvt_linear_srgb_to_srgb);
+    py_bindfunc(mod, "srgb_to_linear_srgb", colorcvt_srgb_to_linear_srgb);
+    py_bindfunc(mod, "srgb_to_hsv", colorcvt_srgb_to_hsv);
+    py_bindfunc(mod, "hsv_to_srgb", colorcvt_hsv_to_srgb);
+    py_bindfunc(mod, "oklch_to_linear_srgb", colorcvt_oklch_to_linear_srgb);
+    py_bindfunc(mod, "linear_srgb_to_oklch", colorcvt_linear_srgb_to_oklch);
+}
+
+#undef DEF_VEC3_WRAPPER

src/modules/linalg.c

@@ -350,9 +350,8 @@ static bool vec2_angle_STATIC(int argc, py_Ref argv) {
     PY_CHECK_ARG_TYPE(0, tp_vec2);
     PY_CHECK_ARG_TYPE(1, tp_vec2);
     float val = atan2f(argv[1]._vec2.y, argv[1]._vec2.x) - atan2f(argv[0]._vec2.y, argv[0]._vec2.x);
-    const float PI = 3.1415926535897932384f;
-    if(val > PI) val -= 2 * PI;
-    if(val < -PI) val += 2 * PI;
+    if(val > PK_M_PI) val -= 2 * (float)PK_M_PI;
+    if(val < -PK_M_PI) val += 2 * (float)PK_M_PI;
     py_newfloat(py_retval(), val);
     return true;
 }

src/modules/math.c

@@ -118,7 +118,7 @@ static bool math_degrees(int argc, py_Ref argv) {
     PY_CHECK_ARGC(1);
     double x;
     if(!py_castfloat(py_arg(0), &x)) return false;
-    py_newfloat(py_retval(), x * 180 / 3.1415926535897932384);
+    py_newfloat(py_retval(), x * PK_M_RAD2DEG);
     return true;
 }
 
@@ -126,10 +126,12 @@ static bool math_radians(int argc, py_Ref argv) {
     PY_CHECK_ARGC(1);
     double x;
     if(!py_castfloat(py_arg(0), &x)) return false;
-    py_newfloat(py_retval(), x * 3.1415926535897932384 / 180);
+    py_newfloat(py_retval(), x * PK_M_DEG2RAD);
     return true;
 }
 
+TWO_ARG_FUNC(fmod, fmod)
+
 static bool math_modf(int argc, py_Ref argv) {
     PY_CHECK_ARGC(1);
     double i;
@@ -157,8 +159,8 @@ static bool math_factorial(int argc, py_Ref argv) {
 void pk__add_module_math() {
     py_Ref mod = py_newmodule("math");
 
-    py_newfloat(py_emplacedict(mod, py_name("pi")), 3.1415926535897932384);
-    py_newfloat(py_emplacedict(mod, py_name("e")), 2.7182818284590452354);
+    py_newfloat(py_emplacedict(mod, py_name("pi")), PK_M_PI);
+    py_newfloat(py_emplacedict(mod, py_name("e")), PK_M_E);
     py_newfloat(py_emplacedict(mod, py_name("inf")), INFINITY);
     py_newfloat(py_emplacedict(mod, py_name("nan")), NAN);
 
@@ -196,6 +198,7 @@ void pk__add_module_math() {
     py_bindfunc(mod, "degrees", math_degrees);
     py_bindfunc(mod, "radians", math_radians);
 
+    py_bindfunc(mod, "fmod", math_fmod);
     py_bindfunc(mod, "modf", math_modf);
     py_bindfunc(mod, "factorial", math_factorial);
 }

tests/70_math.py

@@ -39,6 +39,12 @@ assert math.gcd(10, 7) == 1
 assert math.gcd(10, 10) == 10
 assert math.gcd(-10, 10) == 10
 
+# test fmod
+assert math.fmod(-2.0, 3.0) == -2.0
+assert math.fmod(2.0, 3.0) == 2.0
+assert math.fmod(4.0, 3.0) == 1.0
+assert math.fmod(-4.0, 3.0) == -1.0
+
 # test modf
 x, y = math.modf(1.5)
 assert isclose(x, 0.5)

tests/90_colorcvt.py

@@ -0,0 +1,47 @@
+import colorcvt
+from linalg import vec3
+
+def oklch(expr: str) -> vec3:
+    # oklch(82.33% 0.37 153)
+    expr = expr[6:-1]
+    l, c, h = expr.split()
+    l = float(l[:-1]) / 100
+    return vec3(l, float(c), float(h))
+
+def srgb32(expr: str) -> vec3:
+    # rgb(0, 239, 115)
+    expr = expr[4:-1]
+    r, g, b = expr.split(", ")
+    r, g, b = int(r), int(g), int(b)
+    c = vec3(r, g, b) / 255
+    return colorcvt.srgb_to_linear_srgb(c)
+
+def assert_equal(title: str, a: vec3, b: vec3) -> None:
+    epsilon = 1e-3
+    try:
+        assert abs(a.x - b.x) < epsilon
+        assert abs(a.y - b.y) < epsilon
+        assert abs(a.z - b.z) < epsilon
+    except AssertionError:
+        raise AssertionError(f"{title}\nExpected: {b}, got: {a}")
+
+def test(oklch_expr: str, srgb32_expr: str) -> None:
+    oklch_color = oklch(oklch_expr)
+    srgb32_color = srgb32(srgb32_expr)
+    assert_equal('oklch_to_linear_srgb', colorcvt.oklch_to_linear_srgb(oklch_color), srgb32_color)
+    # in range check
+    oklch_color = colorcvt.linear_srgb_to_oklch(srgb32_color)
+    assert_equal('oklch_to_linear_srgb+', colorcvt.oklch_to_linear_srgb(oklch_color), srgb32_color)
+    assert_equal('linear_srgb_to_oklch+', colorcvt.linear_srgb_to_oklch(srgb32_color), oklch_color)
+
+test("oklch(71.32% 0.1381 153)", "rgb(83, 187, 120)")
+test("oklch(45.15% 0.037 153)", "rgb(70, 92, 76)")
+test("oklch(22.5% 0.0518 153)", "rgb(4, 34, 16)")
+
+test("oklch(100% 0.37 153)", "rgb(255, 255, 255)")
+test("oklch(0% 0.0395 283.24)", "rgb(0, 0, 0)")
+
+# hard samples
+# test("oklch(95% 0.2911 264.18)", "rgb(224, 239, 255)")
+# test("oklch(28.09% 0.2245 153)", "rgb(0, 54, 12)")
+# test("oklch(82.33% 0.37 153)", "rgb(0, 239, 115)")