[css-color-4] Linear transformation matrices are very slightly inaccurate

[nex3]

The sample code for color conversions in Color Level 4 lists a number of matrices to use to convert between color values. I've been pre-computing the linear transformations from one colorspace to another for Sass using lossless numeric logic to avoid floating-point rounding errors, and in doing so I noticed some very very slight inaccuracies in the listed matrices. These are all at the level of floating-point rounding errors, so unlikely to be particularly significant in the real world, but I thought it might be desirable to have the example code be as accurate as possible up to the limits of the floating-point representation.

Here are the most precise 64-bit floating-point representations of the various transformations I've looked at so far, generated by this script:

// srgb to xyz-d65
var M = [
         0.41239079926595950,  0.35758433938387796,  0.18048078840183430,
         0.21263900587151036,  0.71516867876775590,  0.07219231536073371,
         0.01933081871559185,  0.11919477979462598,  0.95053215224966060
];

// xyz-d65 to srgb
var M = [
         3.24096994190452130, -1.53738317757009350, -0.49861076029300330,
        -0.96924363628087980,  1.87596750150772060,  0.04155505740717561,
         0.05563007969699360, -0.20397695888897657,  1.05697151424287860
];

// display-p3 to xyz-d65
var M = [
         0.48657094864821626,  0.26566769316909294,  0.19821728523436250,
         0.22897456406974884,  0.69173852183650620,  0.07928691409374500,
         0.00000000000000000,  0.04511338185890257,  1.04394436890097570
];

// xyz-d65 to display-p3
var M = [
         2.49349691194142450, -0.93138361791912360, -0.40271078445071684,
        -0.82948896956157490,  1.76266406031834680,  0.02362468584194359,
         0.03584583024378433, -0.07617238926804170,  0.95688452400768730
];

// a98-rgb to xyz-d65
var M = [
         0.57666904291013080,  0.18555823790654627,  0.18822864623499472,
         0.29734497525053616,  0.62736356625546600,  0.07529145849399789,
         0.02703136138641237,  0.07068885253582714,  0.99133753683763890
];

// xyz-d65 to a98-rgb
var M = [
         2.04158790381074600, -0.56500697427885960, -0.34473135077832950,
        -0.96924363628087980,  1.87596750150772060,  0.04155505740717561,
         0.01344428063203102, -0.11836239223101823,  1.01517499439120540
];

// rec2020 to xyz-d65
var M = [
         0.63695804830129130,  0.14461690358620838,  0.16888097516417205,
         0.26270021201126703,  0.67799807151887100,  0.05930171646986194,
         0.00000000000000000,  0.02807269304908750,  1.06098505771079090
];

// xyz-d65 to rec2020
var M = [
         1.71665118797126760, -0.35567078377639240, -0.25336628137365980,
        -0.66668435183248900,  1.61648123663493900,  0.01576854581391113,
         0.01763985744531091, -0.04277061325780865,  0.94210312123547400
];

// prophoto-rgb to xyz-d50
var M = [
         0.79776664490064230,  0.13518129740053308,  0.03134773412839220,
         0.28807482881940130,  0.71183523424187300,  0.00008993693872564,
         0.00000000000000000,  0.00000000000000000,  0.82510460251046020
];

// xyz-d50 to prophoto-rgb
var M = [
         1.34578688164715830, -0.25557208737979464, -0.05110186497554526,
        -0.54463070512490190,  1.50824774284514680,  0.02052744743642139,
         0.00000000000000000,  0.00000000000000000,  1.21196754563894520
];

// xyz-d65 to xyz-d50
var M = [
         1.04792979254499660,  0.02294687060160952, -0.05019226628920519,
         0.02962780877005567,  0.99043442675388000, -0.01707379906341879,
        -0.00924304064620452,  0.01505519149029816,  0.75187428142813700
];

// xyz-d50 to xyz-d65
var M = [
         0.95547342148807520, -0.02309845494876452,  0.06325924320057065,
        -0.02836970933386358,  1.00999539808130410,  0.02104144119191730,
         0.01231401486448199, -0.02050764929889898,  1.33036592624212400
];

[svgeesus]

Thanks for reporting this. I agree the values should be as accurate as possible. When a lot of conversions are being done, small error residuals can really build up.

Recently these matrices were updated as a result of this pull request to use rational numbers where possible. You can see these in the editors draft. We didn't update the ProPhoto ones, because the rational forms exceed the precision of JavaScript math. And that pull request did not update the Bradford CAT matrices. So it would be good to get the most accurate versions of these.

I wonder if you would mind checking your results above against the new rational forms?

[nex3]

At a quick spot check, it looks like my results match up exactly with the values from #7320 where available. I've also confirmed that the XYZ-to-XYZ whitepoint-changing matrices have rational forms that are too large to represent as rationals.

[facelessuser]

If these are updated, Oklab should also be updated as its matrices are affected when the XYZ to sRGB matrix changes. Just an FYI as it is easy to forget.

[nex3]

@facelessuser Can you elaborate? I thought the OKLab matrices just came verbatim from https://bottosson.github.io/posts/oklab/#converting-from-xyz-to-oklab.

[facelessuser]

@nex3 Sure. The matrices found on that page are calculated from the RGB transforms from XYZ assuming some white point. That whitepoint did not initially match what the CSS spec was using and caused conversions to not be as good. So, we recalculated, based on the author's guidance.

Now, if we are changing the RGB to XYZ transforms again, conversions will be off again as the Oklab matrices will transform the color from XYZ to sRGB Linear with a matrix other than what was done to get the color from sRGB Linear. They must be in sync.

Ensuring rational numbers for sRGB and then calculating Oklab matrices, this is what I get:

sRGB

===== sRGB =====
--- rgb -> xyz ---
[[0.4123907992659595, 0.35758433938387796, 0.1804807884018343],
 [0.21263900587151036, 0.7151686787677559, 0.07219231536073371],
 [0.01933081871559185, 0.11919477979462599, 0.9505321522496606]]
--- xyz -> rgb ---
[[3.2409699419045213, -1.5373831775700935, -0.4986107602930033],
 [-0.9692436362808798, 1.8759675015077206, 0.04155505740717561],
 [0.05563007969699361, -0.20397695888897657, 1.0569715142428786]]

Oklab

===== sRGB Linear -> lms =====
[[0.412221469470763, 0.5363325372617348, 0.05144599326750221],
 [0.21190349581782522, 0.6806995506452342, 0.1073969535369405],
 [0.08830245919005643, 0.2817188391361215, 0.6299787016738221]]
===== lms -> sRGB Linear =====
[[4.076741636075958, -3.307711539258062, 0.23096990318210417],
 [-1.268437973285032, 2.609757349287689, -0.3413193760026571],
 [-0.0041960761386755155, -0.7034186179359364, 1.707614694074612]]
===== XYZ D65 Linear -> lms =====
[[0.819022437996703, 0.3619062600528904, -0.12887378152098788],
 [0.03298365393238847, 0.9292868615863433, 0.03614466635064236],
 [0.0481771893596242, 0.2642395317527308, 0.6335478284694308]]
===== lms -> XYZ D65 =====
[[1.2268798758459243, -0.5578149944602171, 0.2813910456659646],
 [-0.04057574521480084, 1.1122868032803173, -0.07171105806551635],
 [-0.07637293667466007, -0.42149333240224324, 1.5869240198367818]]
===== lms ** 1/3 -> Oklab =====
[[0.2104542553, 0.793617785, -0.0040720468],
 [1.9779984951, -2.428592205, 0.4505937099],
 [0.0259040371, 0.7827717662, -0.808675766]]
===== Oklab -> lms ** 1/3 =====
[[0.9999999984505198, 0.39633779217376786, 0.2158037580607588],
 [1.0000000088817609, -0.10556134232365635, -0.06385417477170591],
 [1.0000000546724108, -0.08948418209496575, -1.2914855378640917]]

[nex3]

I'm still a little confused because the sRGB linear -> LMS matrix and its inverse don't seem to appear in the sample code, and I believe the other LMS-related matrices don't depend on anything related to sRGB.

[facelessuser]

You are correct that he doesn't show exactly how to calculate them, he only mentions:

Depending on use case you might want to use the sRGB matrices your application uses directly instead of the ones provided here.

I'd have to try and dig up the thread here in the CSS issues where Björn Ottosson specified how to calculate the matrices, but the basic algorithm is here https://github.com/facelessuser/coloraide/blob/main/tools/calc_oklab_matrices.py.

[facelessuser]

Found the Oklab author's original response. Just so it doesn't look like I'm making stuff up 🙃. #6642 (comment)

[nex3]

I'm afraid I'm still a bit lost 😅. I see that you're deriving M1 (a.k.a XYZ_TO_LMS) from M0 in your calc_oklabl_matrices.py, but it still seems like that derivation is independent of the sRGB matrices. Which specific matrix in the color conversion sample code needs to be updated based on the sRGB changes?

[facelessuser]

Forgive me, the SRGBL_TO_LMS conversion (for Okhsl and Okhsv) needs it, not the XYZ to LMS conversion 🤦🏻.

With that said, this whole step, if using rational numbers, can be improved:

XYZ_TO_LMS = alg.divide(M0, alg.outer(alg.dot(M0, xyzt.white_d65), alg.ones(3)))

Now, I am currently a bit skeptical as to how much a difference this will help in the end as the differences are so far out (in decimal places), and most color spaces when being developed don't even bother going through this effort either, but I don't think is bad either. I'm just not sure in real-world cases whether the ever slightly better matrices will be noticed 🙂 .

[facelessuser]

The only reason I would suggest updating the Oklab matrices is they are indeed calculated for CSS just like the RGB matrices. If it is important enough to do it for RGB matrices, then it is probably important enough to do it for Oklab matrices as well 🤷🏻.

[svgeesus]

I see that you're deriving M1 (a.k.a XYZ_TO_LMS) from M0 in your calc_oklabl_matrices.py, but it still seems like that derivation is independent of the sRGB matrices.

All XYZ values assume adaptation to a particular white point.

The D65 white point assumed in the @bottosson post is the ASTM value which is not identical to the value used throughout CSS Color 4. The difference is only in the 5th decimal place, but this was already enough to give significant deviation from the OKLab neutral axis for neutral greys (in sRGB, display-p3, and any other space with a D65 whitepoint) until we fixed it.

Which we did by taking the published sRGB to LMS matrix and factoring out the (revised) sRGB to XYZ-D65 part, giving a new M1.

[nex3]

Related to this: the floating-point errors in the sample code are substantial enough that they're actually visible in the test cases. For example, the xyz-d50-001 test case lists the conversion for #8000 as color(xyz-d50 0.08312 0.154746 0.020961). Running the transformation with lossless rationals yields color(xyz-d50 0.08313899747975508 0.15474758917173437 0.02095602320031179), which would round to color(xyz-d50 0.08314 0.154748 0.020956). The rounding errors will accumulate even more dramatically for transformations that require more computations—I think the lms-001 test is incorrect to the second decimal digit in one channel.

Obviously it's unrealistic to expect implementations to do these transformations losslessly, but it seems like the test cases maybe should reflect the true mathematical values with the expectation that individual implementations will be less accurate than that given the realities of floating-point arithmetic.

[facelessuser]

I'm not sure if that test case was implying that the XYZ D50 value is meant to be the exact translation of #008000 or not. I also don't know if that test was generated with earlier matrices or not (the RGB matrices have changed at least once or twice).

One thing that can be noted is that the XYZ value is not off due to rounding errors as the non-rational matrices are no where near that lossy. They actually have pretty decent round tripping.

>>> xyz = Color('#008000').convert('xyz-d50')
>>> xyz
color(xyz-d50 0.08314 0.15475 0.02096 / 1)
>>> xyz.convert('srgb').convert('xyz-d50')
color(xyz-d50 0.08314 0.15475 0.02096 / 1)

We can also see the test holds true as well since you have to round things off to get the hex values.

>>> Color('color(xyz-d50 0.08312 0.154746 0.020961)').convert('srgb').to_string(hex=True)
#008000

---

EDIT: This isn't arguing against using rational numbers, just some objective observations. I was mainly curious to see how much of a difference using the matrices calculated using rational numbers made with other conversions.

On a side note. I did some experiments using rational numbers for the RGB matrices in my own library.

Using approximate values for the RGB matrices, I get some pretty clean conversions to lab:

>>> Color('white').convert('lab-d65')[:]
[100.0, 0.0, 0.0, 1.0]
>>> Color('white').convert('lab')[:]
[100.0, 0.0, 0.0, 1.0]

But when I used rational matrices for the RGB matrices, I get:

>>> Color('white').convert('lab-d65')[:]
[100.0, 0.0, 0.0, 1.0]
>>> Color('white').convert('lab')[:]
[100.0, 1.1102230246251565e-13, 0.0, 1.0]

Now, does that mean the rational-matrices aren't as good? No, they definitely ensure less error is introduced during the matrix calculation, they just don't prevent all the other errors introduced during the rest of the calculations.

Since most of these color spaces are designed using approximate values anyways, it's probably just dumb luck that the approximate matrices give me such nice values. Would I get even better values for the rational-matrix examples if I ensured the chromatic adaptation matrices also used rational numbers and ensured that all other calculations used rational numbers? Maybe, but it's not realistic to ensure the whole conversion path conforms to using rational numbers.

In general, I really didn't see significant improvements in overall calculations to warrant using rational numbers in my personal library. I also didn't see any significantly worse calculations to argue against using them either. Since the difference between the rational matrices and the approximate matrices is so slight, I guess this is to be expected.

If someone was building an application, I'd personally argue it probably isn't going to matter which is used.

[svgeesus]

For example, the xyz-d50-001 test case

Converting color(xyz-d50 0.08312 0.154746 0.020961) I get an OOG for sRGB: rgb(-0.07681% 50.19777% 0.007727%)

Converting #008000 I get color(xyz-d50 0.083139 0.154748 0.020956) which, back-converted, is rgb(-0.00088% 50.19617% -0.00018%) - closer, still shows round-off error though. I suspect the chromatic adaptation is the culprit here.

Running the transformation with lossless rationals yields color(xyz-d50 0.08313899747975508 0.15474758917173437 0.02095602320031179), which would round to color(xyz-d50 0.08314 0.154748 0.020956).

That is very close to what I got with the current color.js, which is not using the rational form.

[svgeesus]

The Oklab matrices are now 64bit accurate and were copied over from this color.js PR

Bradford matrices were updated as a result of

• [css-color-4] Linear Bradford matrices are inexact #9607

[svgeesus]

All matrices are now using either the rational form, if it exists, or the 64 bit accurate form.