1 /* gstyle-color-convert.c
2  *
3  * Copyright 2016 sebastien lafargue <slafargue@gnome.org>
4  *
5  * This program is free software: you can redistribute it and/or modify
6  * it under the terms of the GNU General Public License as published by
7  * the Free Software Foundation, either version 3 of the License, or
8  * (at your option) any later version.
9  *
10  * This program is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
13  * GNU General Public License for more details.
14  *
15  * You should have received a copy of the GNU General Public License
16  * along with this program.  If not, see <http://www.gnu.org/licenses/>.
17  *
18  * SPDX-License-Identifier: GPL-3.0-or-later
19  */
20 
21 #include <stdio.h>
22 #include <stdlib.h>
23 #include <math.h>
24 
25 #include "gstyle-color-convert.h"
26 
27 #define _2PI          6.2831853071795864769252867665590057683943387987502
28 #define PI_d_6        0.523598775598298873077107230546583814032861566562516
29 #define PI_d_30       0.104719755119659774615421446109316762806572313312503
30 #define _63PI_d_180   1.099557428756427633461925184147826009469009289781285
31 #define _180_d_PI     57.29577951308232087679815481410517033240547246656442
32 #define _25_pow_7     6103515625
33 #define _30PI_d_180   0.523598775598298873077107230546583814032861566562516
34 #define _16_d_116     0.137931034
35 
36 #define one_third     0.333333333333333333333333333333333333333333333333333
37 #define two_third     0.666666666666666666666666666666666666666666666666666
38 
39 #define D65_xref      0.95047
40 #define D65_yref      1.0
41 #define D65_zref      1.08883
42 
43 /* DeltaE algorithm described at:
44  * http://www.ece.rochester.edu/~gsharma/ciede2000/ciede2000noteCRNA.pdf
45  *
46  */
47 
48 /* pow_1_24 and pow_24 are adapted version from babl, published under LGPL */
49 /* babl - dynamically extendable universal pixel conversion library.
50  * Copyright 2012, Red Hat, Inc.
51  */
52 
53 /* Chebychev polynomial terms for x^(5/12) expanded around x=1.5
54  * Non-zero terms calculated via
55  * NIntegrate[(2/Pi)*ChebyshevT[N,u]/Sqrt [1-u^2]*((u+3)/2)^(5/12),{u,-1,1}, PrecisionGoal->20, WorkingPrecision -> 100]
56  * Zeroth term is similar except it uses 1/pi rather than 2/pi.
57  */
58 static const double Cn[] = {
59   1.1758200232996901923,
60   0.16665763094889061230,
61   -0.0083154894939042125035,
62   0.00075187976780420279038,
63   -0.000083240178519391795367,
64   0.000010229209410070008679,
65   -1.3401001466409860246e-6,
66   1.8333422241635376682e-7,
67   -2.5878596761348859722e-8
68 };
69 
70 /* Returns x^(5/12) for x in [1,2) */
71 static inline gdouble
pow512norm(gdouble x)72 pow512norm (gdouble x)
73 {
74   gdouble Tn[9];
75   gdouble u;
76 
77   u = 2.0*x - 3.0;
78   Tn[0] = 1.0;
79   Tn[1] = u;
80   Tn[2] = 2*u*Tn[2-1] - Tn[2-2];
81   Tn[3] = 2*u*Tn[3-1] - Tn[3-2];
82   Tn[4] = 2*u*Tn[4-1] - Tn[4-2];
83   Tn[5] = 2*u*Tn[5-1] - Tn[5-2];
84   Tn[6] = 2*u*Tn[6-1] - Tn[6-2];
85 
86   return Cn[0]*Tn[0] + Cn[1]*Tn[1] + Cn[2]*Tn[2] + Cn[3]*Tn[3] + Cn[4]*Tn[4] + Cn[5]*Tn[5] + Cn[6]*Tn[6];
87 }
88 
89 /* Precalculated (2^N) ^ (5 / 12) */
90 static const gdouble pow2_512[12] =
91 {
92   1.0,
93   1.3348398541700343678,
94   1.7817974362806785482,
95   2.3784142300054420538,
96   3.1748021039363991669,
97   4.2378523774371812394,
98   5.6568542494923805819,
99   7.5509945014535482244,
100   1.0079368399158985525e1,
101   1.3454342644059433809e1,
102   1.7959392772949968275e1,
103   2.3972913230026907883e1
104 };
105 
106 /* Returns x^(1/2.4) == x^(5/12) */
107 static inline gdouble
pow_1_24(gdouble x)108 pow_1_24 (gdouble x)
109 {
110   gdouble s;
111   gint iexp;
112   div_t qr = {0};
113 
114   s = frexp (x, &iexp);
115   s *= 2.0;
116   iexp -= 1;
117 
118   qr = div (iexp, 12);
119   if (qr.rem < 0)
120     {
121       qr.quot -= 1;
122       qr.rem += 12;
123     }
124 
125   return ldexp (pow512norm (s) * pow2_512[qr.rem], 5 * qr.quot);
126 }
127 
128 /* Chebychev polynomial terms for x^(7/5) expanded around x=1.5
129  * Non-zero terms calculated via
130  * NIntegrate[(2/Pi)*ChebyshevT[N,u]/Sqrt [1-u^2]*((u+3)/2)^(7/5),{u,-1,1}, PrecisionGoal->20, WorkingPrecision -> 100]
131  * Zeroth term is similar except it uses 1/pi rather than 2/pi.
132  */
133 static const gdouble iCn[] =
134 {
135   1.7917488588043277509,
136   0.82045614371976854984,
137   0.027694100686325412819,
138   -0.00094244335181762134018,
139   0.000064355540911469709545,
140   -5.7224404636060757485e-6,
141   5.8767669437311184313e-7,
142   -6.6139920053589721168e-8,
143   7.9323242696227458163e-9
144 };
145 
146 /* Returns x^(7/5) for x in [1,2) */
147 static inline gdouble
pow75norm(gdouble x)148 pow75norm (gdouble x)
149 {
150   gdouble Tn[9];
151   gdouble u;
152 
153   u = 2.0*x - 3.0;
154   Tn[0] = 1.0;
155   Tn[1] = u;
156   Tn[2] = 2*u*Tn[2-1] - Tn[2-2];
157   Tn[3] = 2*u*Tn[3-1] - Tn[3-2];
158   Tn[4] = 2*u*Tn[4-1] - Tn[4-2];
159   Tn[5] = 2*u*Tn[5-1] - Tn[5-2];
160 
161   return iCn[0]*Tn[0] + iCn[1]*Tn[1] + iCn[2]*Tn[2] + iCn[3]*Tn[3] + iCn[4]*Tn[4] + iCn[5]*Tn[5];
162 }
163 
164 /* Precalculated (2^N) ^ (7 / 5) */
165 static const gdouble pow2_75[5] =
166 {
167   1.0,
168   2.6390158215457883983,
169   6.9644045063689921093,
170   1.8379173679952558018e+1,
171   4.8502930128332728543e+1
172 };
173 
174 /* Returns x^2.4 == x * x ^1.4 == x * x^(7/5) */
175 static inline gdouble
pow_24(gdouble x)176 pow_24 (gdouble x)
177 {
178   gdouble s;
179   gint iexp;
180   div_t qr = {0};
181 
182   s = frexp (x, &iexp);
183   s *= 2.0;
184   iexp -= 1;
185 
186   qr = div (iexp, 5);
187   if (qr.rem < 0)
188     {
189       qr.quot -= 1;
190       qr.rem += 5;
191     }
192 
193   return x * ldexp (pow75norm (s) * pow2_75[qr.rem], 7 * qr.quot);
194 }
195 
196 static inline gboolean
fix_rgb_bounds(GdkRGBA * rgba)197 fix_rgb_bounds (GdkRGBA *rgba)
198 {
199   gdouble r, g, b;
200   gboolean res = TRUE;
201 
202   r = rgba->red;
203   g = rgba->green;
204   b = rgba->blue;
205 
206   if (r < 0.0)
207     {
208       rgba->red = 0.0;
209       res = FALSE;
210     }
211   else if (r > 1.0)
212     {
213       rgba->red = 1.0;
214       res = FALSE;
215     }
216 
217   if (g < 0.0)
218     {
219       rgba->green = 0.0;
220       res = FALSE;
221     }
222   else if (g > 1.0)
223     {
224       rgba->green = 1.0;
225       res = FALSE;
226     }
227 
228     if (b < 0.0)
229     {
230       rgba->blue = 0.0;
231       res = FALSE;
232     }
233   else if (b > 1.0)
234     {
235       rgba->blue = 1.0;
236       res = FALSE;
237     }
238 
239   return res;
240 }
241 
242 static inline void
gstyle_color_convert_rgb_to_srgb(GdkRGBA * rgba,gdouble * red,gdouble * green,gdouble * blue)243 gstyle_color_convert_rgb_to_srgb (GdkRGBA *rgba,
244                                   gdouble *red,
245                                   gdouble *green,
246                                   gdouble *blue)
247 {
248   /* rgba and srgb values range [0, 1] */
249 
250   *red = (rgba->red > 0.04045) ? pow_24 (((rgba->red + 0.055) / 1.055)) : rgba->red / 12.92;
251   *green = (rgba->green > 0.04045) ? pow_24 (((rgba->green + 0.055) / 1.055)) : rgba->green / 12.92;
252   *blue = (rgba->blue > 0.04045) ? pow_24 (((rgba->blue + 0.055) / 1.055)) : rgba->blue / 12.92;
253 }
254 
255 static inline void
gstyle_color_convert_srgb_to_rgb(gdouble red,gdouble green,gdouble blue,GdkRGBA * rgba)256 gstyle_color_convert_srgb_to_rgb (gdouble  red,
257                                   gdouble  green,
258                                   gdouble  blue,
259                                   GdkRGBA *rgba)
260 {
261   /* rgba and srgb values range [0, 1] */
262 
263   rgba->red = (red > 0.0031308) ? (pow_1_24 (red) * 1.055) - 0.055 : red * 12.92;
264   rgba->green = (green > 0.0031308) ? (pow_1_24 (green) * 1.055) - 0.055 : green * 12.92;
265   rgba->blue = (blue > 0.0031308) ? (pow_1_24 (blue) * 1.055) - 0.055 : blue * 12.92;
266 
267   fix_rgb_bounds (rgba);
268 }
269 
270 static inline void
gstyle_color_convert_srgb_to_xyz(gdouble red,gdouble green,gdouble blue,GstyleXYZ * xyz)271 gstyle_color_convert_srgb_to_xyz (gdouble    red,
272                                   gdouble    green,
273                                   gdouble    blue,
274                                   GstyleXYZ *xyz)
275 {
276   /* srgb range [0, 1] x [0, 0.9505] y [0, 1] z [0, 1.08883] Observer= 2°, Illuminant= D65 */
277 
278   xyz->x = (red * 0.4124564 + green * 0.3575761 + blue * 0.1804375);
279   xyz->y = (red * 0.2126729 + green * 0.7151522 + blue * 0.0721750);
280   xyz->z = (red * 0.0193339 + green * 0.1191920 + blue * 0.9503041);
281 
282   //fix_xyz_bounds (xyz);
283 }
284 
285 static inline void
gstyle_color_convert_xyz_to_srgb(GstyleXYZ * xyz,gdouble * red,gdouble * green,gdouble * blue)286 gstyle_color_convert_xyz_to_srgb (GstyleXYZ *xyz,
287                                   gdouble   *red,
288                                   gdouble   *green,
289                                   gdouble   *blue)
290 {
291   /* srgb range [0, 1] x [0, 0.9505] y [0, 1] z [0, 1.08883] Observer= 2°, Illuminant= D65 */
292 
293   *red   = xyz->x *  3.2404542 + xyz->y * -1.5371385 + xyz->z * -0.4985314;
294   *green = xyz->x * -0.9692660 + xyz->y *  1.8760108 + xyz->z *  0.0415560;
295   *blue  = xyz->x *  0.0556434 + xyz->y * -0.2040259 + xyz->z *  1.0572252;
296 }
297 
298 inline void
gstyle_color_convert_cielab_to_xyz(GstyleCielab * lab,GstyleXYZ * xyz)299 gstyle_color_convert_cielab_to_xyz (GstyleCielab *lab,
300                                     GstyleXYZ    *xyz)
301 {
302   gdouble tmp_x, tmp_y, tmp_z;
303   gdouble pow3_x, pow3_y, pow3_z;
304 
305   tmp_y = (lab->l + 16.0 ) / 116.0;
306   tmp_x = lab->a / 500.0 + tmp_y;
307   tmp_z = tmp_y - lab->b / 200.0;
308 
309   /* far faster than pow (x, 3) */
310   pow3_x = tmp_x * tmp_x * tmp_x;
311   pow3_y = tmp_y * tmp_y * tmp_y;
312   pow3_z = tmp_z * tmp_z * tmp_z;
313 
314   tmp_x = (pow3_x > 0.008856) ? pow3_x : (tmp_x - _16_d_116) / 7.787;
315   tmp_y = (pow3_y > 0.008856) ? pow3_y : (tmp_y - _16_d_116) / 7.787;
316   tmp_z = (pow3_z > 0.008856) ? pow3_z : (tmp_z - _16_d_116) / 7.787;
317 
318   xyz->x = tmp_x * D65_xref;
319   xyz->y = tmp_y * D65_yref;
320   xyz->z = tmp_z * D65_zref;
321 }
322 
323 /* fastpow (x, 0.333333333) */
324 inline void
gstyle_color_convert_xyz_to_cielab(GstyleXYZ * xyz,GstyleCielab * lab)325 gstyle_color_convert_xyz_to_cielab (GstyleXYZ    *xyz,
326                                     GstyleCielab *lab)
327 {
328   /* Observer= 2°, Illuminant= D65 */
329   gdouble x, y, z;
330 
331   x = xyz->x / D65_xref;
332   y = xyz->y / D65_yref;
333   z = xyz->z / D65_zref;
334 
335   x = (x > 0.008856) ? cbrt (x) : (x * 7.787) + _16_d_116;
336   y = (y > 0.008856) ? cbrt (y) : (y * 7.787) + _16_d_116;
337   z = (z > 0.008856) ? cbrt (z) : (z * 7.787) + _16_d_116;
338 
339   lab->l = y * 116.0 - 16.0;
340   lab->a = (x - y) * 500.0;
341   lab->b = (y - z) * 200.0;
342 }
343 
344 /**
345  * gstyle_color_convert_rgb_to_hsl:
346  * @rgba: a #GdkRGBA struct.
347  * @hue: (out): The hue component of a hsl color in range  [0.0-360.0[
348  * @saturation: (out): The saturation component of a hsl color in range [0.0-100.0]
349  * @lightness: (out): The lightness component of a hsl color in range [0.0-100.0]
350  *
351  * Convert rgb components to HSL ones.
352  * The alpha component is not used because it doesn't change in the conversion.
353  *
354  */
355 void
gstyle_color_convert_rgb_to_hsl(GdkRGBA * rgba,gdouble * hue,gdouble * saturation,gdouble * lightness)356 gstyle_color_convert_rgb_to_hsl (GdkRGBA *rgba,
357                                  gdouble *hue,
358                                  gdouble *saturation,
359                                  gdouble *lightness)
360 {
361   gdouble tmp_hue = 0.0;
362   gdouble tmp_saturation = 0.0;
363   gdouble tmp_lightness;
364 
365   gdouble red = rgba->red;
366   gdouble green = rgba->green;
367   gdouble blue = rgba->blue;
368 
369   gdouble min;
370   gdouble max;
371   gdouble d;
372   gdouble max_min;
373 
374   if (red > green)
375     {
376       max = (red > blue) ? red : blue;
377       min = (green < blue) ? green : blue;
378     }
379   else
380     {
381       max = (green > blue) ? green : blue;
382       min = (red < blue) ? red : blue;
383     }
384 
385   max_min = max + min;
386   tmp_lightness = max_min / 2.0;
387   if (max != min)
388     {
389       d = max - min;
390       tmp_saturation = (tmp_lightness > 0.5) ? d / (2.0 - max - min) : d / max_min;
391       if (max == red)
392         tmp_hue = (green - blue) / d + (green < blue ? 6.0 : 0.0);
393       else if (max == green)
394         tmp_hue = (blue - red) / d + 2.0;
395       else
396         tmp_hue = (red - green) / d + 4.0;
397     }
398 
399   if (hue != NULL)
400     *hue = tmp_hue * 60.0;
401 
402   if (saturation != NULL)
403     *saturation = tmp_saturation * 100.0;
404 
405   if (lightness != NULL)
406     *lightness = tmp_lightness * 100.0;
407 }
408 
409 static inline gdouble
hue2rgb(gdouble m1,gdouble m2,gdouble hue)410 hue2rgb (gdouble m1,
411          gdouble m2,
412          gdouble hue)
413 {
414   while (hue < 0.0)
415     hue += 360.0;
416 
417   while (hue > 360.0)
418     hue -= 360.0;
419 
420   if (hue < 60.0)
421     return m1 + (m2 - m1) * hue / 60.0;
422 
423   if (hue < 180.0)
424     return m2;
425 
426   if (hue < 240.0)
427     return m1 + (m2 - m1) * (240.0 - hue) / 60.0;
428 
429   return m1;
430 }
431 
432 /**
433  * gstyle_color_convert_hsl_to_rgb:
434  * @hue: The hue component of a hsl color in range  [0.0-360.0[
435  * @saturation: The saturation component of a hsl color in range [0.0-100.0]
436  * @lightness: The lightness component of a hsl color in range [0.0-100.0]
437  * @rgba: a #GdkRGBA.
438  *
439  * Convert RGB components to HSL ones.
440  * The alpha component is not used because it doesn't change in the conversion.
441  *
442  */
443 void
gstyle_color_convert_hsl_to_rgb(gdouble hue,gdouble saturation,gdouble lightness,GdkRGBA * rgba)444 gstyle_color_convert_hsl_to_rgb (gdouble   hue,
445                                  gdouble   saturation,
446                                  gdouble   lightness,
447                                  GdkRGBA  *rgba)
448 {
449   gdouble m1;
450   gdouble m2;
451 
452   if (saturation == 0.0)
453     rgba->red = rgba->green = rgba->blue = lightness;
454   else
455     {
456       m2 = (lightness > 0.5) ? lightness + saturation - (lightness * saturation) : lightness * (1.0 + saturation);
457       m1 = 2.0 * lightness - m2;
458 
459       rgba->red = hue2rgb (m1, m2, hue + 120.0);
460       rgba->green = hue2rgb (m1, m2, hue);
461       rgba->blue = hue2rgb (m1, m2, hue - 120.0);
462     }
463 }
464 
465 /**
466  * gstyle_color_convert_hsv_to_rgb:
467  * @hue: The hue component of a hsv color in range  [0.0-1.0[
468  * @saturation: The saturation component of a hsv color in range [0.0-1.0]
469  * @value: The value component of a hsv color in range [0.0-1.0]
470  * @rgba: a #GdkRGBA.
471  *
472  * Convert HSV components to RGB ones.
473  * The alpha component is not used because it doesn't change in the conversion.
474  *
475  */
476 void
gstyle_color_convert_hsv_to_rgb(gdouble hue,gdouble saturation,gdouble value,GdkRGBA * rgba)477 gstyle_color_convert_hsv_to_rgb (gdouble   hue,
478                                  gdouble   saturation,
479                                  gdouble   value,
480                                  GdkRGBA  *rgba)
481 {
482   gdouble f;
483   gdouble p;
484   gdouble q;
485   gdouble t;
486   gint i;
487 
488   if (saturation == 0.0)
489     rgba->red = rgba->green = rgba->blue = value;
490   else
491     {
492       hue *= 6.0;
493       if (hue == 6.0)
494         hue = 0.0;
495 
496       i = (int)hue;
497       f = hue - i;
498       p = value * (1.0 - saturation);
499       q = value * (1.0 - saturation * f);
500       t = value * (1.0 - saturation * (1.0 - f));
501 
502       switch (i)
503         {
504         case 0:
505           rgba->red = value;
506           rgba->green = t;
507           rgba->blue = p;
508           break;
509 
510         case 1:
511           rgba->red = q;
512           rgba->green = value;
513           rgba->blue = p;
514           break;
515 
516         case 2:
517           rgba->red = p;
518           rgba->green = value;
519           rgba->blue = t;
520           break;
521 
522         case 3:
523           rgba->red = p;
524           rgba->green = q;
525           rgba->blue = value;
526           break;
527 
528         case 4:
529           rgba->red = t;
530           rgba->green = p;
531           rgba->blue = value;
532           break;
533 
534         case 5:
535           rgba->red = value;
536           rgba->green = p;
537           rgba->blue = q;
538           break;
539 
540         default:
541           g_assert_not_reached ();
542         }
543     }
544 }
545 
546 /**
547  * gstyle_color_convert_rgb_to_xyz:
548  * @rgba: An #GdkRGBA.
549  * @xyz: a #GstyleXYZ.
550  *
551  * Convert RGB components to XYZ ones.
552  * The alpha component is not used because it doesn't change in the conversion.
553  *
554  */
555 void
gstyle_color_convert_rgb_to_xyz(GdkRGBA * rgba,GstyleXYZ * xyz)556 gstyle_color_convert_rgb_to_xyz (GdkRGBA   *rgba,
557                                  GstyleXYZ *xyz)
558 {
559   gdouble srgb_red, srgb_green, srgb_blue;
560 
561   gstyle_color_convert_rgb_to_srgb (rgba, &srgb_red, &srgb_green, &srgb_blue);
562   gstyle_color_convert_srgb_to_xyz (srgb_red, srgb_green, srgb_blue, xyz);
563 }
564 
565 /**
566  * gstyle_color_convert_rgb_to_hsv:
567  * @rgba: An #GdkRGBA.
568  * @hue: (out): The hue component of a hsv color in range  [0.0-1.0[
569  * @saturation: (out): The saturation component of a hsv color in range [0.0-1.0]
570  * @value: (out): The value component of a hsv color in range [0.0-1.0]
571  *
572  * Convert RGB components to HSV ones.
573  * The alpha component is not used because it doesn't change in the conversion.
574  *
575  */
576 void
gstyle_color_convert_rgb_to_hsv(GdkRGBA * rgba,gdouble * hue,gdouble * saturation,gdouble * value)577 gstyle_color_convert_rgb_to_hsv (GdkRGBA *rgba,
578                                  gdouble *hue,
579                                  gdouble *saturation,
580                                  gdouble *value)
581 {
582   gdouble vmin, vmax, delta;
583   gdouble d_red, d_green, d_blue;
584   gdouble delta_d_2;
585 
586   if (rgba->red > rgba->green)
587     {
588       vmax = (rgba->red > rgba->blue) ? rgba->red : rgba->blue;
589       vmin = (rgba->green < rgba->blue) ? rgba->green : rgba->blue;
590     }
591   else
592     {
593       vmax = (rgba->green > rgba->blue) ? rgba->green : rgba->blue;
594       vmin = (rgba->red < rgba->blue) ? rgba->red : rgba->blue;
595     }
596 
597   delta = vmax - vmin;
598   delta_d_2 = delta / 2.0;
599 
600   *value = vmax;
601 
602   if (delta < 1e-20 )
603     *hue = *saturation = 0.0;
604   else
605     {
606       *saturation = delta / vmax;
607 
608       d_red   = ((vmax - rgba->red)   / 6.0 + delta_d_2) / delta;
609       d_green = ((vmax - rgba->green) / 6.0 + delta_d_2) / delta;
610       d_blue  = ((vmax - rgba->blue)  / 6.0 + delta_d_2) / delta;
611 
612       if (vmax == rgba->red)
613          *hue = d_blue - d_green;
614       else if (vmax == rgba->green)
615          *hue = one_third + d_red - d_blue;
616       else if (vmax == rgba->blue)
617         *hue = two_third + d_green - d_red;
618 
619       if (*hue < 0.0)
620         *hue += 1.0;
621       else if (*hue > 1.0 )
622         *hue -= 1.0;
623     }
624 }
625 
626 /**
627  * gstyle_color_convert_rgb_to_cielab:
628  * @rgba: An #GdkRGBA.
629  * @lab: (out): a #GstyleCieLab struct.
630  *
631  * Convert RGB components to CIELAB ones.
632  * The alpha component is not used because it doesn't change in the conversion.
633  *
634  */
635 void
gstyle_color_convert_rgb_to_cielab(GdkRGBA * rgba,GstyleCielab * lab)636 gstyle_color_convert_rgb_to_cielab (GdkRGBA      *rgba,
637                                     GstyleCielab *lab)
638 {
639   gdouble srgb_red, srgb_green, srgb_blue;
640   GstyleXYZ xyz;
641 
642   gstyle_color_convert_rgb_to_srgb (rgba, &srgb_red, &srgb_green, &srgb_blue);
643   gstyle_color_convert_srgb_to_xyz (srgb_red, srgb_green, srgb_blue, &xyz);
644   gstyle_color_convert_xyz_to_cielab (&xyz, lab);
645 }
646 
647 /**
648  * gstyle_color_convert_cielab_to_rgb:
649  * @lab: a #GstyleCieLab struct.
650  * @rgba: (out): An #GdkRGBA.
651  *
652  * Convert CIELAB components to RGB ones.
653  * The alpha component is not used because it doesn't change in the conversion.
654  *
655  */
656 void
gstyle_color_convert_cielab_to_rgb(GstyleCielab * lab,GdkRGBA * rgba)657 gstyle_color_convert_cielab_to_rgb (GstyleCielab *lab,
658                                     GdkRGBA      *rgba)
659 {
660   gdouble srgb_red, srgb_green, srgb_blue;
661   GstyleXYZ xyz;
662 
663   gstyle_color_convert_cielab_to_xyz (lab, &xyz);
664   gstyle_color_convert_xyz_to_srgb (&xyz, &srgb_red, &srgb_green, &srgb_blue);
665   gstyle_color_convert_srgb_to_rgb (srgb_red, srgb_green, srgb_blue, rgba);
666 }
667 
668 inline void
gstyle_color_convert_xyz_to_rgb(GstyleXYZ * xyz,GdkRGBA * rgba)669 gstyle_color_convert_xyz_to_rgb (GstyleXYZ *xyz,
670                                  GdkRGBA   *rgba)
671 {
672   gdouble srgb_red, srgb_green, srgb_blue;
673 
674   gstyle_color_convert_xyz_to_srgb (xyz, &srgb_red, &srgb_green, &srgb_blue);
675   gstyle_color_convert_srgb_to_rgb (srgb_red, srgb_green, srgb_blue, rgba);
676 }
677 
678 /**
679  * gstyle_color_convert_hsv_to_xyz:
680  * @hue:  The hue component of a hsv color in range  [0.0-1.0[
681  * @saturation: The saturation component of a hsv color in range [0.0-1.0]
682  * @value: The value component of a hsv color in range [0.0-1.0]
683  * @xyz: (out): An #GstyleXYZ.
684  *
685  * Convert HSV components to XYZ ones.
686  *
687  */
688 void
gstyle_color_convert_hsv_to_xyz(gdouble hue,gdouble saturation,gdouble value,GstyleXYZ * xyz)689 gstyle_color_convert_hsv_to_xyz (gdouble    hue,
690                                  gdouble    saturation,
691                                  gdouble    value,
692                                  GstyleXYZ *xyz)
693 {
694   gdouble srgb_red, srgb_green, srgb_blue;
695   GdkRGBA rgba;
696 
697   gstyle_color_convert_hsv_to_rgb (hue, saturation, value, &rgba);
698   gstyle_color_convert_rgb_to_srgb (&rgba, &srgb_red, &srgb_green, &srgb_blue);
699   gstyle_color_convert_srgb_to_xyz (srgb_red, srgb_green, srgb_blue, xyz);
700 }
701 
702 /**
703  * gstyle_color_convert_xyz_to_hsv:
704  * @xyz: An #GstyleXYZ.
705  * @hue: (out): The hue component of a hsv color in range  [0.0-1.0[
706  * @saturation: (out): The saturation component of a hsv color in range [0.0-1.0]
707  * @value: (out): The value component of a hsv color in range [0.0-1.0]
708  *
709  * Convert XYZ components to HSV ones.
710  *
711  */
712 void
gstyle_color_convert_xyz_to_hsv(GstyleXYZ * xyz,gdouble * hue,gdouble * saturation,gdouble * value)713 gstyle_color_convert_xyz_to_hsv (GstyleXYZ *xyz,
714                                  gdouble   *hue,
715                                  gdouble   *saturation,
716                                  gdouble   *value)
717 {
718   GdkRGBA rgba;
719   gdouble srgb_red, srgb_green, srgb_blue;
720 
721   gstyle_color_convert_xyz_to_srgb (xyz, &srgb_red, &srgb_green, &srgb_blue);
722   gstyle_color_convert_srgb_to_rgb (srgb_red, srgb_green, srgb_blue, &rgba);
723   gstyle_color_convert_rgb_to_hsv (&rgba, hue, saturation, value);
724 }
725 
726 /**
727  * gstyle_color_delta_e:
728  * @lab1: a #GstyleCielab.
729  * @lab2: a #GstyleCielab.
730  *
731  * Compute the color difference between lab1 and lab2,
732  * based on the deltaE CIEDE2000 formula.
733  *
734  * Returns: the deltaE value.
735  *
736  */
737 gdouble
gstyle_color_delta_e(GstyleCielab * lab1,GstyleCielab * lab2)738 gstyle_color_delta_e (GstyleCielab *lab1,
739                       GstyleCielab *lab2)
740 {
741   gdouble ap1, Cp1, hp1, Cab1;
742   gdouble ap2, Cp2, hp2, Cab2;
743   gdouble Cab;
744 
745   gdouble Cp1Cp2;
746   gdouble Cab_pow_7, G;
747 
748   gdouble DLp, DCp, DHp;
749   gdouble dhp, Lp, Cp, hp;
750 
751   gdouble T;
752   gdouble Dtheta_rad;
753   gdouble Rc, RT;
754 
755   gdouble _50Lp_pow2;
756   gdouble SL, SC, SH;
757 
758   gdouble lab1_bb = lab1->b * lab1->b;
759   gdouble lab2_bb = lab2->b * lab2->b;
760 
761   Cab1 = sqrt (lab1->a * lab1->a + lab1_bb);
762   Cab2 = sqrt (lab2->a * lab2->a + lab2_bb);
763   Cab = (Cab1 + Cab2) / 2.0;
764   Cab_pow_7 = pow (Cab, 7);
765 
766   G = 0.5 * (1.0 - sqrt (Cab_pow_7 / (Cab_pow_7 + _25_pow_7)));
767 
768   ap1 = (1.0 + G) * lab1->a;
769   ap2 = (1.0 + G) * lab2->a;
770 
771   Cp1 = sqrt (ap1 * ap1 + lab1_bb);
772   Cp2 = sqrt (ap2 * ap2 + lab2_bb);
773   Cp1Cp2 = (Cp1 * Cp2);
774 
775   if (ap1 == 0 && lab1->b == 0)
776     hp1 = 0.0;
777   else
778     {
779       hp1 = atan2 (lab1->b, ap1);
780       if (hp1 < 0)
781         hp1 += _2PI;
782     }
783 
784   if (ap2 == 0 && lab2->b == 0)
785     hp2 = 0.0;
786   else
787     {
788       hp2 = atan2 (lab2->b, ap2);
789       if (hp2 < 0)
790         hp2 += _2PI;
791     }
792 
793   DLp = (lab2->l - lab1->l);
794   DCp = (Cp2 - Cp1);
795 
796   if (Cp1Cp2 == 0.0)
797     {
798       dhp = 0.0;
799       DHp = 0.0;
800 
801       hp = hp1 + hp2;
802     }
803   else
804     {
805       dhp = (hp2 - hp1);
806       if (dhp > G_PI)
807         dhp -= _2PI;
808 
809       if (dhp < -G_PI)
810         dhp += _2PI;
811 
812       DHp = 2.0 * sqrt (Cp1Cp2) * sin (dhp / 2.0);
813 
814       hp = (hp1 + hp2) / 2.0;
815       if (fabs (hp1 - hp2) > G_PI)
816         hp -= G_PI;
817 
818       if (hp < 0)
819         hp += _2PI;
820     }
821 
822   Lp = (lab1->l + lab2->l) / 2.0;
823   Cp = (Cp1 + Cp2) / 2.0;
824 
825   T = 1.0 - 0.17 * cos (hp - PI_d_6) +
826       0.24 * cos (2.0 * hp) +
827       0.32 * cos (3.0 * hp + PI_d_30) -
828       0.20 * cos (4.0 * hp - _63PI_d_180);
829 
830   Dtheta_rad = _30PI_d_180 * exp (-pow (((_180_d_PI * hp - 275.0) / 25.0), 2.0));
831 
832   Rc = 2.0 * sqrt (pow (Cp, 7.0) / (pow (Cp, 7.0) + _25_pow_7));
833 
834   _50Lp_pow2 = (Lp - 50.0) * (Lp - 50.0);
835   SL = 1.0 + (0.015 * _50Lp_pow2 / sqrt (20.0 + _50Lp_pow2));
836   SC = 1.0 + 0.045 * Cp;
837   SH = 1.0 + 0.015 * Cp * T;
838 
839   RT = -sin (2.0 * Dtheta_rad) * Rc;
840 
841   return sqrt( pow ((DLp / SL), 2.0) +
842          pow((DCp / SC), 2.0) +
843          pow((DHp / SH), 2.0) +
844          RT * (DCp / SC) * (DHp / SH));
845 }
846