diff options
author | Niklas Haas <git@nand.wakku.to> | 2014-06-22 08:33:43 +0200 |
---|---|---|
committer | wm4 <wm4@nowhere> | 2014-06-22 19:02:06 +0200 |
commit | 17762a1919947db0e66e025bd2084de896eaa3fa (patch) | |
tree | 86aaf7c3c7bf37c9bf7a3c7f19eac9bacf7cd07e /video/csputils.c | |
parent | 204fed4d5b4aa20b5a6b5824f5d4e71ccbaf87fb (diff) | |
download | mpv-17762a1919947db0e66e025bd2084de896eaa3fa.tar.bz2 mpv-17762a1919947db0e66e025bd2084de896eaa3fa.tar.xz |
video: Generate an accurate CMS matrix instead of hard-coding
This also avoids an extra matrix multiplication when using :srgb, making
that path both more efficient and also eliminating more hard-coded
values.
In addition, the previously hard-coded XYZ to RGB matrix will be
dynamically generated.
Diffstat (limited to 'video/csputils.c')
-rw-r--r-- | video/csputils.c | 298 |
1 files changed, 216 insertions, 82 deletions
diff --git a/video/csputils.c b/video/csputils.c index 43dd17511a..8b36615de6 100644 --- a/video/csputils.c +++ b/video/csputils.c @@ -213,6 +213,195 @@ void mp_gen_gamma_map(uint8_t *map, int size, float gamma) } } +void mp_invert_matrix3x3(float m[3][3]) +{ + float m00 = m[0][0], m01 = m[0][1], m02 = m[0][2], + m10 = m[1][0], m11 = m[1][1], m12 = m[1][2], + m20 = m[2][0], m21 = m[2][1], m22 = m[2][2]; + + // calculate the adjoint + m[0][0] = (m11 * m22 - m21 * m12); + m[0][1] = -(m01 * m22 - m21 * m02); + m[0][2] = (m01 * m12 - m11 * m02); + m[1][0] = -(m10 * m22 - m20 * m12); + m[1][1] = (m00 * m22 - m20 * m02); + m[1][2] = -(m00 * m12 - m10 * m02); + m[2][0] = (m10 * m21 - m20 * m11); + m[2][1] = -(m00 * m21 - m20 * m01); + m[2][2] = (m00 * m11 - m10 * m01); + + // calculate the determinant (as inverse == 1/det * adjoint, + // adjoint * m == identity * det, so this calculates the det) + float det = m00 * m[0][0] + m10 * m[0][1] + m20 * m[0][2]; + det = 1.0f / det; + + for (int i = 0; i < 3; i++) { + for (int j = 0; j < 3; j++) + m[i][j] *= det; + } +} + +// A := A * B +void mp_mul_matrix3x3(float a[3][3], float b[3][3]) +{ + float a00 = a[0][0], a01 = a[0][1], a02 = a[0][2], + a10 = a[1][0], a11 = a[1][1], a12 = a[1][2], + a20 = a[2][0], a21 = a[2][1], a22 = a[2][2]; + + for (int i = 0; i < 3; i++) { + a[0][i] = a00 * b[0][i] + a01 * b[1][i] + a02 * b[2][i]; + a[1][i] = a10 * b[0][i] + a11 * b[1][i] + a12 * b[2][i]; + a[2][i] = a20 * b[0][i] + a21 * b[1][i] + a22 * b[2][i]; + } +} + +/** + * \brief return the primaries associated with a certain mp_csp_primaries val + * \param csp the colorspace for which to return the primaries + */ +struct mp_csp_primaries mp_get_csp_primaries(enum mp_csp_prim spc) +{ + /* + Values from: ITU-R Recommendations BT.601-7, BT.709-5, BT.2020-0 + + https://www.itu.int/dms_pubrec/itu-r/rec/bt/R-REC-BT.601-7-201103-I!!PDF-E.pdf + https://www.itu.int/dms_pubrec/itu-r/rec/bt/R-REC-BT.709-5-200204-I!!PDF-E.pdf + https://www.itu.int/dms_pubrec/itu-r/rec/bt/R-REC-BT.2020-0-201208-I!!PDF-E.pdf + */ + + static const struct mp_csp_col_xy d65 = {0.3127, 0.3290}; + + switch (spc) { + case MP_CSP_PRIM_BT_601_525: + return (struct mp_csp_primaries) { + .red = {0.630, 0.340}, + .green = {0.310, 0.595}, + .blue = {0.155, 0.070}, + .white = d65 + }; + case MP_CSP_PRIM_BT_601_625: + return (struct mp_csp_primaries) { + .red = {0.640, 0.330}, + .green = {0.290, 0.600}, + .blue = {0.150, 0.060}, + .white = d65 + }; + case MP_CSP_PRIM_BT_709: + return (struct mp_csp_primaries) { + .red = {0.640, 0.330}, + .green = {0.300, 0.600}, + .blue = {0.150, 0.060}, + .white = d65 + }; + case MP_CSP_PRIM_BT_2020: + return (struct mp_csp_primaries) { + .red = {0.708, 0.292}, + .green = {0.170, 0.797}, + .blue = {0.131, 0.046}, + .white = d65 + }; + default: + return (struct mp_csp_primaries) {{0}}; + } +} + +// Compute the RGB/XYZ matrix as described here: +// http://www.brucelindbloom.com/index.html?Eqn_RGB_XYZ_Matrix.html +void mp_get_rgb2xyz_matrix(struct mp_csp_primaries space, float m[3][3]) +{ + float S[3], X[4], Z[4]; + + // Convert from CIE xyY to XYZ. Note that Y=1 holds true for all primaries + X[0] = space.red.x / space.red.y; + X[1] = space.green.x / space.green.y; + X[2] = space.blue.x / space.blue.y; + X[3] = space.white.x / space.white.y; + + Z[0] = (1 - space.red.x - space.red.y) / space.red.y; + Z[1] = (1 - space.green.x - space.green.y) / space.green.y; + Z[2] = (1 - space.blue.x - space.blue.y) / space.blue.y; + Z[3] = (1 - space.white.x - space.white.y) / space.white.y; + + // S = XYZ^-1 * W + for (int i = 0; i < 3; i++) { + m[0][i] = X[i]; + m[1][i] = 1; + m[2][i] = Z[i]; + } + + mp_invert_matrix3x3(m); + + for (int i = 0; i < 3; i++) + S[i] = m[i][0] * X[3] + m[i][1] * 1 + m[i][2] * Z[3]; + + // M = [Sc * XYZc] + for (int i = 0; i < 3; i++) { + m[0][i] = S[i] * X[i]; + m[1][i] = S[i] * 1; + m[2][i] = S[i] * Z[i]; + } +} + +/** + * \brief get the coefficients of the source -> bt2020 cms matrix + * \param src primaries of the source gamut + * \param dest primaries of the destination gamut + * \param m array to store coefficients into + */ +void mp_get_cms_matrix(struct mp_csp_primaries src, struct mp_csp_primaries dest, float m[3][3]) +{ + // RGBd<-RGBs = RGBd<-XYZd * XYZd<-XYZs * XYZs<-RGBs + // Equations from: http://www.brucelindbloom.com/index.html?Math.html + float tmp[3][3] = {{0}}; + + // RGBd<-XYZd, inverted from XYZd<-RGBd + mp_get_rgb2xyz_matrix(dest, m); + mp_invert_matrix3x3(m); + + // Chromatic adaptation, only needed if the white point differs + if (fabs(src.white.x - dest.white.x) > 1e-6 || + fabs(src.white.y - dest.white.y) > 1e-6) { + // XYZd<-XYZs = Ma^-1 * (I*[Cd/Cs]) * Ma + // http://www.brucelindbloom.com/index.html?Eqn_ChromAdapt.html + float C[3][2]; + + // Ma = Bradford matrix, arguably most popular method in use today. + // This is derived experimentally and thus hard-coded. + float bradford[3][3] = { + { 0.8951, 0.2664, -0.1614 }, + { -0.7502, 1.7135, 0.0367 }, + { 0.0389, -0.0685, 1.0296 }, + }; + + for (int i = 0; i < 3; i++) { + // source cone + C[i][0] = bradford[i][0] * src.white.x / src.white.y + + bradford[i][1] * 1 + + bradford[i][2] * (1 - src.white.x - src.white.y) / src.white.y; + + // dest cone + C[i][1] = bradford[i][0] * dest.white.x / dest.white.y + + bradford[i][1] * 1 + + bradford[i][2] * (1 - dest.white.x - dest.white.y) / dest.white.y; + } + + // tmp := I * [Cd/Cs] * Ma + for (int i = 0; i < 3; i++) + tmp[i][i] = C[i][1] / C[i][0]; + + mp_mul_matrix3x3(tmp, bradford); + + // M := M * Ma^-1 * tmp + mp_invert_matrix3x3(bradford); + mp_mul_matrix3x3(m, bradford); + mp_mul_matrix3x3(m, tmp); + } + + // XYZs<-RGBs + mp_get_rgb2xyz_matrix(src, tmp); + mp_mul_matrix3x3(m, tmp); +} + /* Fill in the Y, U, V vectors of a yuv2rgb conversion matrix * based on the given luma weights of the R, G and B components (lr, lg, lb). * lr+lg+lb is assumed to equal 1. @@ -247,53 +436,6 @@ static void luma_coeffs(float m[3][4], float lr, float lg, float lb) } /** - * Get the coefficients of the source -> bt2020 gamut mapping matrix, - * given an identifier describing the source gamut primaries. - */ -void mp_get_cms_matrix(enum mp_csp_prim source, float m[3][3]) -{ - // Conversion matrices to BT.2020 primaries - // These were computed using: http://lpaste.net/101796 - switch (source) { - case MP_CSP_PRIM_BT_601_525: { - static const float from_525[3][3] = { - {0.5952542, 0.3493139, 0.0554319}, - {0.0812437, 0.8915033, 0.0272530}, - {0.0155123, 0.0819116, 0.9025760}, - }; - memcpy(m, from_525, sizeof(from_525)); - break; - } - case MP_CSP_PRIM_BT_601_625: { - static const float from_625[3][3] = { - {0.6550368, 0.3021610, 0.0428023}, - {0.0721406, 0.9166311, 0.0112283}, - {0.0171134, 0.0978535, 0.8850332}, - }; - memcpy(m, from_625, sizeof(from_625)); - break; - } - case MP_CSP_PRIM_BT_709: { - static const float from_709[3][3] = { - {0.6274039, 0.3292830, 0.0433131}, - {0.0690973, 0.9195404, 0.0113623}, - {0.0163914, 0.0880133, 0.8955953}, - }; - memcpy(m, from_709, sizeof(from_709)); - break; - } - case MP_CSP_PRIM_BT_2020: - default: { - // No conversion necessary. This matrix should not even be used - // in this case, but assign it to the identity just to be safe. - static const float ident[3][3] = {{1, 0, 0}, {0, 1, 0}, {0, 0, 1}}; - memcpy(m, ident, sizeof(ident)); - break; - } - } -} - -/** * \brief get the coefficients of the yuv -> rgb conversion matrix * \param params struct specifying the properties of the conversion like * brightness, ... @@ -329,12 +471,20 @@ void mp_get_yuv2rgb_coeffs(struct mp_csp_params *params, float m[3][4]) break; } case MP_CSP_XYZ: { - static const float xyz_to_rgb[3][4] = { - {3.2404542, -1.5371385, -0.4985314}, - {-0.9692660, 1.8760108, 0.0415560}, - {0.0556434, -0.2040259, 1.0572252}, - }; - memcpy(m, xyz_to_rgb, sizeof(xyz_to_rgb)); + // The vo should probably not be using a matrix generated by this + // function for XYZ sources, but if it does, let's just assume we + // want BT.709. + float xyz_to_rgb[3][3]; + mp_get_rgb2xyz_matrix(mp_get_csp_primaries(MP_CSP_PRIM_BT_709), xyz_to_rgb); + mp_invert_matrix3x3(xyz_to_rgb); + + for (int i = 0; i < 3; i++) { + for (int j = 0; j < 3; j++) + m[i][j] = xyz_to_rgb[i][j]; + + m[i][3] = 0; + } + levels_in = -1; break; } @@ -512,44 +662,28 @@ int mp_csp_equalizer_set(struct mp_csp_equalizer *eq, const char *property, void mp_invert_yuv2rgb(float out[3][4], float in[3][4]) { - float m00 = in[0][0], m01 = in[0][1], m02 = in[0][2], m03 = in[0][3], - m10 = in[1][0], m11 = in[1][1], m12 = in[1][2], m13 = in[1][3], - m20 = in[2][0], m21 = in[2][1], m22 = in[2][2], m23 = in[2][3]; + float tmp[3][3]; - // calculate the adjoint - out[0][0] = (m11 * m22 - m21 * m12); - out[0][1] = -(m01 * m22 - m21 * m02); - out[0][2] = (m01 * m12 - m11 * m02); - out[1][0] = -(m10 * m22 - m20 * m12); - out[1][1] = (m00 * m22 - m20 * m02); - out[1][2] = -(m00 * m12 - m10 * m02); - out[2][0] = (m10 * m21 - m20 * m11); - out[2][1] = -(m00 * m21 - m20 * m01); - out[2][2] = (m00 * m11 - m10 * m01); + for (int i = 0; i < 3; i++) { + for (int j = 0; j < 3; j++) + tmp[i][j] = in[i][j]; + } - // calculate the determinant (as inverse == 1/det * adjoint, - // adjoint * m == identity * det, so this calculates the det) - float det = m00 * out[0][0] + m10 * out[0][1] + m20 * out[0][2]; - det = 1.0f / det; + mp_invert_matrix3x3(tmp); - out[0][0] *= det; - out[0][1] *= det; - out[0][2] *= det; - out[1][0] *= det; - out[1][1] *= det; - out[1][2] *= det; - out[2][0] *= det; - out[2][1] *= det; - out[2][2] *= det; + for (int i = 0; i < 3; i++) { + for (int j = 0; j < 3; j++) + out[i][j] = tmp[i][j]; + } // fix the constant coefficient // rgb = M * yuv + C // M^-1 * rgb = yuv + M^-1 * C // yuv = M^-1 * rgb - M^-1 * C // ^^^^^^^^^^ - out[0][3] = -(out[0][0] * m03 + out[0][1] * m13 + out[0][2] * m23); - out[1][3] = -(out[1][0] * m03 + out[1][1] * m13 + out[1][2] * m23); - out[2][3] = -(out[2][0] * m03 + out[2][1] * m13 + out[2][2] * m23); + out[0][3] = -(out[0][0] * in[0][3] + out[0][1] * in[1][3] + out[0][2] * in[2][3]); + out[1][3] = -(out[1][0] * in[0][3] + out[1][1] * in[1][3] + out[1][2] * in[2][3]); + out[2][3] = -(out[2][0] * in[0][3] + out[2][1] * in[1][3] + out[2][2] * in[2][3]); } // Multiply the color in c with the given matrix. |