CIEDE2000 implementation in C++

Function version: v1.0.0
Site statistics
Number of visits617
Number of files viewed493 + 309

This page presents a reference implementation of the CIEDE2000 color difference formula in C++. If you wish to obtain an exact match with third-party implementations up to 10 decimal places, you may need to make some changes to the source code, including commenting and uncommenting a few lines, which can be applied automatically via the link below.

Diagram of the full-form CIEDE2000 formula with L*a*b* components and adjustments

The ΔE2000 function in C++

Let’s consider the more common and academic (Sharma, 2005) of the two formulations.

// This function written in C++ is not affiliated with the CIE (International Commission on Illumination),
// and is released into the public domain. It is provided "as is" without any warranty, express or implied.

#include <cmath>

// Expressly defining pi ensures that the code works on different platforms.
#ifndef M_PI
#define M_PI 3.14159265358979323846264338327950288419716939937511
#endif

// The classic CIE Ξ”E2000 implementation, which operates on two L*a*b* colors, and returns their difference.
// "l" ranges from 0 to 100, while "a" and "b" are unbounded and commonly clamped to the range of -128 to 127.
template<typename T>
static T ciede_2000(const T l_1, const T a_1, const T b_1, const T l_2, const T a_2, const T b_2) {
	// Working in C++ with the CIEDE2000 color-difference formula.
	// k_l, k_c, k_h are parametric factors to be adjusted according to
	// different viewing parameters such as textures, backgrounds...
	const T k_l = T(1.0);
	const T k_c = T(1.0);
	const T k_h = T(1.0);
	T n = (std::sqrt(a_1 * a_1 + b_1 * b_1) + std::sqrt(a_2 * a_2 + b_2 * b_2)) * T(0.5);
	n = n * n * n * n * n * n * n;
	// A factor involving chroma raised to the power of 7 designed to make
	// the influence of chroma on the total color difference more accurate.
	n = T(1.0) + T(0.5) * (T(1.0) - std::sqrt(n / (n + T(6103515625.0))));
	// Application of the chroma correction factor.
	const T c_1 = std::sqrt(a_1 * a_1 * n * n + b_1 * b_1);
	const T c_2 = std::sqrt(a_2 * a_2 * n * n + b_2 * b_2);
	// atan2 is preferred over atan because it accurately computes the angle of
	// a point (x, y) in all quadrants, handling the signs of both coordinates.
	T h_1 = std::atan2(b_1, a_1 * n);
	T h_2 = std::atan2(b_2, a_2 * n);
	h_1 += (h_1 < T(0.0)) * T(2.0) * T(M_PI);
	h_2 += (h_2 < T(0.0)) * T(2.0) * T(M_PI);
	n = std::fabs(h_2 - h_1);
	// Cross-implementation consistent rounding.
	if (T(M_PI) - T(1E-14) < n && n < T(M_PI) + T(1E-14))
		n = T(M_PI);
	// When the hue angles lie in different quadrants, the straightforward
	// average can produce a mean that incorrectly suggests a hue angle in
	// the wrong quadrant, the next lines handle this issue.
	T h_m = (h_1 + h_2) * T(0.5);
	T h_d = (h_2 - h_1) * T(0.5);
	h_d += (T(M_PI) < n) * T(M_PI);
	// πŸ“œ Sharma’s formulation doesn’t use the next line, but the one after it,
	// and these two variants differ by Β±0.0003 on the final color differences.
	h_m += (T(M_PI) < n) * T(M_PI);
	// h_m += (T(M_PI) < n) * ((h_m < T(M_PI)) - (T(M_PI) <= h_m)) * T(M_PI);
	const T p = T(36.0) * h_m - T(55.0) * T(M_PI);
	n = (c_1 + c_2) * T(0.5);
	n = n * n * n * n * n * n * n;
	// The hue rotation correction term is designed to account for the
	// non-linear behavior of hue differences in the blue region.
	const T r_t = T(-2.0) * std::sqrt(n / (n + T(6103515625.0)))
			* std::sin(T(M_PI) / T(3.0) * std::exp(p * p / (T(-25.0) * T(M_PI) * T(M_PI))));
	n = (l_1 + l_2) * T(0.5);
	n = (n - T(50.0)) * (n - T(50.0));
	// Lightness.
	const T l = (l_2 - l_1) / (k_l * (T(1.0) + T(3.0) / T(200.0) * n / std::sqrt(T(20.0) + n)));
	// These coefficients adjust the impact of different harmonic
	// components on the hue difference calculation.
	const T t = T(1.0) 	+ T(6.0) / T(25.0) * std::sin(T(2.0) * h_m + T(M_PI) / T(2.0))
				+ T(8.0) / T(25.0) * std::sin(T(3.0) * h_m + T(8.0) * T(M_PI) / T(15.0))
				- T(17.0) / T(100.0) * std::sin(h_m + T(M_PI) / T(3.0))
				- T(1.0) / T(5.0) * std::sin(T(4.0) * h_m + T(3.0) * T(M_PI) / T(20.0));
	n = c_1 + c_2;
	// Hue.
	const T h = T(2.0) * std::sqrt(c_1 * c_2) * std::sin(h_d) / (k_h * (T(1.0) + T(3.0) / T(400.0) * n * t));
	// Chroma.
	const T c = (c_2 - c_1) / (k_c * (T(1.0) + T(9.0) / T(400.0) * n));
	// Returning the square root ensures that dE00 accurately reflects the
	// geometric distance in color space, which can range from 0 to around 185.
	return std::sqrt(l * l + h * h + c * c + c * h * r_t);
}

// GitHub Project : https://github.com/michel-leonard/ciede2000-color-matching
//   Online Tests : https://michel-leonard.github.io/ciede2000-color-matching

// L1 = 96.5   a1 = 47.8   b1 = 4.6
// L2 = 96.8   a2 = 53.2   b2 = -4.1
// CIE Ξ”E00 = 4.6680978034 (Bruce Lindbloom, Netflix’s VMAF, ...)
// CIE Ξ”E00 = 4.6680847226 (Gaurav Sharma, OpenJDK, ...)
// Deviation between implementations β‰ˆ 1.3e-5

// See the source code comments for easy switching between these two widely used Ξ”E*00 implementation variants.

Source code accuracy and reliability

The difference between Sharma and Lindbloom formulations never exceeds Β±0.0003 on the final Ξ”E2000, which corresponds to the usual difference measured between two 32-bit implementations and is imperceptible to the human eye. Our 64-bit implementations, all consistent with each other, guarantee at least 10 correct decimal places, so choosing one formulation over another mainly depends on the desired interoperability. The formulation that appears by default on this page is the most commonly used (its micro-advantage lies in its community anchoring and its greater lightness than its analog when vectorized).

✎ If you find a comment in the source code that does not correspond to another language, please inform the author of the site, who will study your suggestion and incorporate it into the source code.

How do you convert RGB colors to L*a*b*?

Go to the AWK, C, Dart, Java, JavaScript, Kotlin, Lua, PHP, Python, Ruby or Rust page where such a converter (using D65 illuminant) is already implemented in addition to the color comparison function.

CIELAB value ranges and interpretation of the Ξ”E2000

In the CIELAB color space, the L* component represents lightness and typically ranges from 0 (black) to 100 (white). The a* and b* components represent color axes: a* goes from green to red, while b* goes from blue to yellow. In practice, a* and b* values usually fall between -128 and +127, although they can slightly exceed these limits depending on the color conversion.

Example of two colors presenting a just-noticeable difference (JND) according to CIEDE2000
Color 1Color 2Value of Ξ”E2000
1
2
3
Examples of CIEDE2000 values calculated between two distinct colors
Color 1Color 2Value of Ξ”E2000
5
10
15

k_l, k_c and k_h Parameters

The parameters k_l, k_c, and k_h are weighting factors applied to the lightness (Ξ”L*), chroma (Ξ”C*), and hue (Ξ”H*) terms in the CIEDE2000 formula. Their default value is 1, which corresponds to the standard viewing conditions recommended by the International Commission on Illumination. In practice, these coefficients are adjusted to reflect specific conditions: for example, k_l = 2 is sometimes used to give more weight to lightness differences (common in printing), while k_c or k_h may be reduced to increase tolerance to saturation or hue variations depending on quality control requirements. Depending on the context, these coefficients typically range between 0.5 and 2.

Ξ”E2000 (CIEDE2000) measures the perceived difference between two colors: 0 means the colors are identical, and higher values (up to around 185 in extreme cases) indicate a larger difference. For example, a Ξ”E2000 value around 5 means the colors are close, while a value around 15 means they are clearly different.

Example of use in C++

// Compute the Delta E (CIEDE2000) color difference between two L*a*b* colors in C++
// L1 = 75.5        a1 = 22.5        b1 = -2.5
// L2 = 76.5        a2 = 16.5        b2 = 2.25

const auto delta_e_32_bits = ciede_2000<float>(L1, a1, b1, L2, a2, b2);
const auto delta_e_64_bits = ciede_2000<double>(L1, a1, b1, L2, a2, b2);

std::printf("DeltaE 2000 (float):  %.8g\n", delta_e_32_bits);
std::printf("DeltaE 2000 (double): %.8g\n", delta_e_64_bits);

// .................................................. This shows a ΔE2000 of 4.8786078559
// As explained in the comments, compliance with Gaurav Sharma would display 4.8785929856

Test results

The driver, written in the C99 language and featuring 250 precise static tests, has proved that this C++ function is interoperable with the CIEDE2000 function available in other programming languages.

CIEDE2000 Verification Summary :
  First Verified Line : 45,65.26,-37.27,78.78,-80,-108,102.75551666659558236
             Duration : 13.66 s
            Successes : 10000000
               Errors : 0
      Average Delta E : 62.9626
    Average Deviation : 0
    Maximum Deviation : 0

Files to download

Feel free to use these files provided by Michel, even for commercial purposes.

Site statistics : downloads
FileSizeNumber of clicks
ciede-2000.cpp4 KB82
ciede-2000-driver.cpp6 KB76
ciede-2000-identity.cpp9 KB73
ciede-2000-random.cpp7 KB75
identity.yml4 KB64
test-cpp.yml3 KB61
vs-dvisvgm.yml5 KB62
reference-dataset.txt4 KB309
Click on cpp.zip to receive all these files in an archive.

Community

If you’d like to leave your opinion on this C++ source code or CIEDE2000 in general, the guestbook already contains 1 messages in English, and 9 messages in total, so let us know what you think.