Calculating RGB↔XYZ matrix
Posted by Michał ‘mina86’ Nazarewicz on 3rd of February 2019 | (cite)
I’ve recently found myself in need of an RGB↔XYZ transformation matrix expressed to the maximum possible precision. Sources on the Internet typically limit the precision to a handful decimal places so I’ve performed do the calculations myself.
What we’re looking for is a 3-by-3 matrix
Derivation
RGB primary colours are typically provided as chromaticity usually expressed as x and y coordinates. For sRGB they are defined in IEC 61966 standard (and also Rec. 709 document which is 170 francs cheaper). Converting them to XYZ space is simple:
Red | ||
---|---|---|
Green | ||
Blue | ||
White (D65) |
That’s where reference white point comes into play. Its coordinates in linear RGB space,
For each colour
All quantities on the right-hand side are known, therefore
Final formula
Given chromaticity of primary colours of an RGB space (
which can also be written as
Matrix converting coordinates in the opposite direction is the inverse of
Implementation
Having the theoretical part covered, it’s time to put the equations into practice, which brings up a question of language to use. Rust programmers may simply reach for the rgb_derivation crate which provides all the necessary algorithms and supports rational arithmetic using big integers. However, for the sake of brevity, rather than presenting the Rust implementation, we’re going to use Python instead. The code begins with an overview of the algorithm:
Chromaticity = collections.namedtuple('Chromaticity', 'x y') def calculate_rgb_matrix(primaries, white): M_prime = (tuple(c.x / c.y for c in primaries), tuple(1 for _ in primaries), tuple((1 - c.x - c.y) / c.y for c in primaries)) W = (white.x / white.y, 1, (1 - white.x - white.y) / white.y) Y = mul_matrix_by_column(inverse_3x3_matrix(M_prime), W) return mul_matrix_by_diag(M_prime, Y)
The function first constructs
All operations on matrices are delegated to separate functions. Since the matrices are small there is no need for any sophisticated algorithms. Instead, the most straightforward multiplication algorithms are used:
def mul_matrix_by_column(matrix, column): return tuple( sum(row[i] * column[i] for i in range(len(row))) for row in matrix) def mul_matrix_by_diag(matrix, column): return tuple( tuple(row[c] * column[c] for c in range(len(column))) for row in matrix)
Only the function inverting a 3-by-3 matrix is somewhat more complex:
def inverse_3x3_matrix(matrix): def cofactor(row, col): minor = [matrix[r][c] for r in (0, 1, 2) if r != row for c in (0, 1, 2) if c != col] a, b, c, d = minor det_minor = a * d - b * c return det_minor if (row ^ col) & 1 == 0 else -det_minor comatrix = tuple( tuple(cofactor(row, col) for col in (0, 1, 2)) for row in (0, 1, 2)) det = sum(matrix[0][col] * comatrix[0][col] for col in (0, 1, 2)) return tuple( tuple(comatrix[col][row] / det for col in (0, 1, 2)) for row in (0, 1, 2))
It first constructs a matrix of cofactors of the input (i.e. comatrix
). Because function’s argument is always a 3-by-3 matrix, each minor’s determinant can be trivially calculated using
The above code works for any RGB system. To get result for sRGB its sRGB’s primaries and white point chromaticities need to be passed:
def calculate_srgb_matrix(): primaries = (Chromaticity(fractions.Fraction(64, 100), fractions.Fraction(33, 100)), Chromaticity(fractions.Fraction(30, 100), fractions.Fraction(60, 100)), Chromaticity(fractions.Fraction(15, 100), fractions.Fraction( 6, 100))) white = Chromaticity(fractions.Fraction(312713, 1000000), fractions.Fraction(329016, 1000000)) return calculate_rgb_matrix(primaries, white)
Full implementation with other bells and whistles can be found inside of the ansi_colours
repository.
The matrix
Finally, there’s the matrix itself. Rust programmers are in luck again. The srgb crate provides the matrix along with other values and functions needed to work with sRGB colour space. If the values are needed in other programming languages, they can be copied from the following listing:
⎡ 0.4124108464885388 0.3575845678529519 0.18045380393360833 ⎤ M ≈ ⎢ 0.21264934272065283 0.7151691357059038 0.07218152157344333 ⎥ ⎣ 0.019331758429150258 0.11919485595098397 0.9503900340503373 ⎦ ⎡ 33786752 / 81924984 29295110 / 81924984 14783675 / 81924984 ⎤ M = ⎢ 8710647 / 40962492 29295110 / 40962492 2956735 / 40962492 ⎥ ⎣ 4751262 / 245774952 29295110 / 245774952 233582065 / 245774952 ⎦ ⎡ 3.240812398895283 -1.5373084456298136 -0.4985865229069666 ⎤ M⁻¹ ≈ ⎢ -0.9692430170086407 1.8759663029085742 0.04155503085668564 ⎥ ⎣ 0.055638398436112804 -0.20400746093241362 1.0571295702861434 ⎦ ⎡ 4277208 / 1319795 -2028932 / 1319795 -658032 / 1319795 ⎤ M⁻¹ = ⎢ -70985202 / 73237775 137391598 / 73237775 3043398 / 73237775 ⎥ ⎣ 164508 / 2956735 -603196 / 2956735 3125652 / 2956735 ⎦
It’s important to remember this is not all that is needed to perform conversion between sRGB and XYZ colour spaces. The matrix assumes a linear RGB coordinates normalised to [0, 1] range which requires dealing with gamma correction. I’ve written a follow-up article which describes how to deal with that.
Updated in March 2021 with more precise value for the D65 standard illuminant. This affected the resulting numbers slightly.