1//---------------------------------------------------------------------------------
2//
3// Little Color Management System
4// Copyright (c) 1998-2017 Marti Maria Saguer
5//
6// Permission is hereby granted, free of charge, to any person obtaining
7// a copy of this software and associated documentation files (the "Software"),
8// to deal in the Software without restriction, including without limitation
9// the rights to use, copy, modify, merge, publish, distribute, sublicense,
10// and/or sell copies of the Software, and to permit persons to whom the Software
11// is furnished to do so, subject to the following conditions:
12//
13// The above copyright notice and this permission notice shall be included in
14// all copies or substantial portions of the Software.
15//
16// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
17// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO
18// THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
19// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE
20// LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION
21// OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
22// WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
23//
24//---------------------------------------------------------------------------------
25//
26
27#include "lcms2_internal.h"
28
29// CIECAM 02 appearance model. Many thanks to Jordi Vilar for the debugging.
30
31// ---------- Implementation --------------------------------------------
32
33typedef struct {
34
35 cmsFloat64Number XYZ[3];
36 cmsFloat64Number RGB[3];
37 cmsFloat64Number RGBc[3];
38 cmsFloat64Number RGBp[3];
39 cmsFloat64Number RGBpa[3];
40 cmsFloat64Number a, b, h, e, H, A, J, Q, s, t, C, M;
41 cmsFloat64Number abC[2];
42 cmsFloat64Number abs[2];
43 cmsFloat64Number abM[2];
44
45} CAM02COLOR;
46
47typedef struct {
48
49 CAM02COLOR adoptedWhite;
50 cmsFloat64Number LA, Yb;
51 cmsFloat64Number F, c, Nc;
52 cmsUInt32Number surround;
53 cmsFloat64Number n, Nbb, Ncb, z, FL, D;
54} cmsCIECAM02;
55
56
57static
58cmsFloat64Number compute_n(cmsCIECAM02* pMod)
59{
60 return (pMod -> Yb / pMod -> adoptedWhite.XYZ[1]);
61}
62
63static
64cmsFloat64Number compute_z(cmsCIECAM02* pMod)
65{
66 return (1.48 + pow(pMod -> n, 0.5));
67}
68
69static
70cmsFloat64Number computeNbb(cmsCIECAM02* pMod)
71{
72 return (0.725 * pow((1.0 / pMod -> n), 0.2));
73}
74
75static
76cmsFloat64Number computeFL(cmsCIECAM02* pMod)
77{
78 cmsFloat64Number k, FL;
79
80 k = 1.0 / ((5.0 * pMod->LA) + 1.0);
81 FL = 0.2 * pow(k, 4.0) * (5.0 * pMod->LA) + 0.1 *
82 (pow((1.0 - pow(k, 4.0)), 2.0)) *
83 (pow((5.0 * pMod->LA), (1.0 / 3.0)));
84
85 return FL;
86}
87
88static
89cmsFloat64Number computeD(cmsCIECAM02* pMod)
90{
91 cmsFloat64Number D;
92
93 D = pMod->F - (1.0/3.6)*(exp(((-pMod ->LA-42) / 92.0)));
94
95 return D;
96}
97
98
99static
100CAM02COLOR XYZtoCAT02(CAM02COLOR clr)
101{
102 clr.RGB[0] = (clr.XYZ[0] * 0.7328) + (clr.XYZ[1] * 0.4296) + (clr.XYZ[2] * -0.1624);
103 clr.RGB[1] = (clr.XYZ[0] * -0.7036) + (clr.XYZ[1] * 1.6975) + (clr.XYZ[2] * 0.0061);
104 clr.RGB[2] = (clr.XYZ[0] * 0.0030) + (clr.XYZ[1] * 0.0136) + (clr.XYZ[2] * 0.9834);
105
106 return clr;
107}
108
109static
110CAM02COLOR ChromaticAdaptation(CAM02COLOR clr, cmsCIECAM02* pMod)
111{
112 cmsUInt32Number i;
113
114 for (i = 0; i < 3; i++) {
115 clr.RGBc[i] = ((pMod -> adoptedWhite.XYZ[1] *
116 (pMod->D / pMod -> adoptedWhite.RGB[i])) +
117 (1.0 - pMod->D)) * clr.RGB[i];
118 }
119
120 return clr;
121}
122
123
124static
125CAM02COLOR CAT02toHPE(CAM02COLOR clr)
126{
127 cmsFloat64Number M[9];
128
129 M[0] =(( 0.38971 * 1.096124) + (0.68898 * 0.454369) + (-0.07868 * -0.009628));
130 M[1] =(( 0.38971 * -0.278869) + (0.68898 * 0.473533) + (-0.07868 * -0.005698));
131 M[2] =(( 0.38971 * 0.182745) + (0.68898 * 0.072098) + (-0.07868 * 1.015326));
132 M[3] =((-0.22981 * 1.096124) + (1.18340 * 0.454369) + ( 0.04641 * -0.009628));
133 M[4] =((-0.22981 * -0.278869) + (1.18340 * 0.473533) + ( 0.04641 * -0.005698));
134 M[5] =((-0.22981 * 0.182745) + (1.18340 * 0.072098) + ( 0.04641 * 1.015326));
135 M[6] =(-0.009628);
136 M[7] =(-0.005698);
137 M[8] =( 1.015326);
138
139 clr.RGBp[0] = (clr.RGBc[0] * M[0]) + (clr.RGBc[1] * M[1]) + (clr.RGBc[2] * M[2]);
140 clr.RGBp[1] = (clr.RGBc[0] * M[3]) + (clr.RGBc[1] * M[4]) + (clr.RGBc[2] * M[5]);
141 clr.RGBp[2] = (clr.RGBc[0] * M[6]) + (clr.RGBc[1] * M[7]) + (clr.RGBc[2] * M[8]);
142
143 return clr;
144}
145
146static
147CAM02COLOR NonlinearCompression(CAM02COLOR clr, cmsCIECAM02* pMod)
148{
149 cmsUInt32Number i;
150 cmsFloat64Number temp;
151
152 for (i = 0; i < 3; i++) {
153 if (clr.RGBp[i] < 0) {
154
155 temp = pow((-1.0 * pMod->FL * clr.RGBp[i] / 100.0), 0.42);
156 clr.RGBpa[i] = (-1.0 * 400.0 * temp) / (temp + 27.13) + 0.1;
157 }
158 else {
159 temp = pow((pMod->FL * clr.RGBp[i] / 100.0), 0.42);
160 clr.RGBpa[i] = (400.0 * temp) / (temp + 27.13) + 0.1;
161 }
162 }
163
164 clr.A = (((2.0 * clr.RGBpa[0]) + clr.RGBpa[1] +
165 (clr.RGBpa[2] / 20.0)) - 0.305) * pMod->Nbb;
166
167 return clr;
168}
169
170static
171CAM02COLOR ComputeCorrelates(CAM02COLOR clr, cmsCIECAM02* pMod)
172{
173 cmsFloat64Number a, b, temp, e, t, r2d, d2r;
174
175 a = clr.RGBpa[0] - (12.0 * clr.RGBpa[1] / 11.0) + (clr.RGBpa[2] / 11.0);
176 b = (clr.RGBpa[0] + clr.RGBpa[1] - (2.0 * clr.RGBpa[2])) / 9.0;
177
178 r2d = (180.0 / 3.141592654);
179 if (a == 0) {
180 if (b == 0) clr.h = 0;
181 else if (b > 0) clr.h = 90;
182 else clr.h = 270;
183 }
184 else if (a > 0) {
185 temp = b / a;
186 if (b > 0) clr.h = (r2d * atan(temp));
187 else if (b == 0) clr.h = 0;
188 else clr.h = (r2d * atan(temp)) + 360;
189 }
190 else {
191 temp = b / a;
192 clr.h = (r2d * atan(temp)) + 180;
193 }
194
195 d2r = (3.141592654 / 180.0);
196 e = ((12500.0 / 13.0) * pMod->Nc * pMod->Ncb) *
197 (cos((clr.h * d2r + 2.0)) + 3.8);
198
199 if (clr.h < 20.14) {
200 temp = ((clr.h + 122.47)/1.2) + ((20.14 - clr.h)/0.8);
201 clr.H = 300 + (100*((clr.h + 122.47)/1.2)) / temp;
202 }
203 else if (clr.h < 90.0) {
204 temp = ((clr.h - 20.14)/0.8) + ((90.00 - clr.h)/0.7);
205 clr.H = (100*((clr.h - 20.14)/0.8)) / temp;
206 }
207 else if (clr.h < 164.25) {
208 temp = ((clr.h - 90.00)/0.7) + ((164.25 - clr.h)/1.0);
209 clr.H = 100 + ((100*((clr.h - 90.00)/0.7)) / temp);
210 }
211 else if (clr.h < 237.53) {
212 temp = ((clr.h - 164.25)/1.0) + ((237.53 - clr.h)/1.2);
213 clr.H = 200 + ((100*((clr.h - 164.25)/1.0)) / temp);
214 }
215 else {
216 temp = ((clr.h - 237.53)/1.2) + ((360 - clr.h + 20.14)/0.8);
217 clr.H = 300 + ((100*((clr.h - 237.53)/1.2)) / temp);
218 }
219
220 clr.J = 100.0 * pow((clr.A / pMod->adoptedWhite.A),
221 (pMod->c * pMod->z));
222
223 clr.Q = (4.0 / pMod->c) * pow((clr.J / 100.0), 0.5) *
224 (pMod->adoptedWhite.A + 4.0) * pow(pMod->FL, 0.25);
225
226 t = (e * pow(((a * a) + (b * b)), 0.5)) /
227 (clr.RGBpa[0] + clr.RGBpa[1] +
228 ((21.0 / 20.0) * clr.RGBpa[2]));
229
230 clr.C = pow(t, 0.9) * pow((clr.J / 100.0), 0.5) *
231 pow((1.64 - pow(0.29, pMod->n)), 0.73);
232
233 clr.M = clr.C * pow(pMod->FL, 0.25);
234 clr.s = 100.0 * pow((clr.M / clr.Q), 0.5);
235
236 return clr;
237}
238
239
240static
241CAM02COLOR InverseCorrelates(CAM02COLOR clr, cmsCIECAM02* pMod)
242{
243
244 cmsFloat64Number t, e, p1, p2, p3, p4, p5, hr, d2r;
245 d2r = 3.141592654 / 180.0;
246
247 t = pow( (clr.C / (pow((clr.J / 100.0), 0.5) *
248 (pow((1.64 - pow(0.29, pMod->n)), 0.73)))),
249 (1.0 / 0.9) );
250 e = ((12500.0 / 13.0) * pMod->Nc * pMod->Ncb) *
251 (cos((clr.h * d2r + 2.0)) + 3.8);
252
253 clr.A = pMod->adoptedWhite.A * pow(
254 (clr.J / 100.0),
255 (1.0 / (pMod->c * pMod->z)));
256
257 p1 = e / t;
258 p2 = (clr.A / pMod->Nbb) + 0.305;
259 p3 = 21.0 / 20.0;
260
261 hr = clr.h * d2r;
262
263 if (fabs(sin(hr)) >= fabs(cos(hr))) {
264 p4 = p1 / sin(hr);
265 clr.b = (p2 * (2.0 + p3) * (460.0 / 1403.0)) /
266 (p4 + (2.0 + p3) * (220.0 / 1403.0) *
267 (cos(hr) / sin(hr)) - (27.0 / 1403.0) +
268 p3 * (6300.0 / 1403.0));
269 clr.a = clr.b * (cos(hr) / sin(hr));
270 }
271 else {
272 p5 = p1 / cos(hr);
273 clr.a = (p2 * (2.0 + p3) * (460.0 / 1403.0)) /
274 (p5 + (2.0 + p3) * (220.0 / 1403.0) -
275 ((27.0 / 1403.0) - p3 * (6300.0 / 1403.0)) *
276 (sin(hr) / cos(hr)));
277 clr.b = clr.a * (sin(hr) / cos(hr));
278 }
279
280 clr.RGBpa[0] = ((460.0 / 1403.0) * p2) +
281 ((451.0 / 1403.0) * clr.a) +
282 ((288.0 / 1403.0) * clr.b);
283 clr.RGBpa[1] = ((460.0 / 1403.0) * p2) -
284 ((891.0 / 1403.0) * clr.a) -
285 ((261.0 / 1403.0) * clr.b);
286 clr.RGBpa[2] = ((460.0 / 1403.0) * p2) -
287 ((220.0 / 1403.0) * clr.a) -
288 ((6300.0 / 1403.0) * clr.b);
289
290 return clr;
291}
292
293static
294CAM02COLOR InverseNonlinearity(CAM02COLOR clr, cmsCIECAM02* pMod)
295{
296 cmsUInt32Number i;
297 cmsFloat64Number c1;
298
299 for (i = 0; i < 3; i++) {
300 if ((clr.RGBpa[i] - 0.1) < 0) c1 = -1;
301 else c1 = 1;
302 clr.RGBp[i] = c1 * (100.0 / pMod->FL) *
303 pow(((27.13 * fabs(clr.RGBpa[i] - 0.1)) /
304 (400.0 - fabs(clr.RGBpa[i] - 0.1))),
305 (1.0 / 0.42));
306 }
307
308 return clr;
309}
310
311static
312CAM02COLOR HPEtoCAT02(CAM02COLOR clr)
313{
314 cmsFloat64Number M[9];
315
316 M[0] = (( 0.7328 * 1.910197) + (0.4296 * 0.370950));
317 M[1] = (( 0.7328 * -1.112124) + (0.4296 * 0.629054));
318 M[2] = (( 0.7328 * 0.201908) + (0.4296 * 0.000008) - 0.1624);
319 M[3] = ((-0.7036 * 1.910197) + (1.6975 * 0.370950));
320 M[4] = ((-0.7036 * -1.112124) + (1.6975 * 0.629054));
321 M[5] = ((-0.7036 * 0.201908) + (1.6975 * 0.000008) + 0.0061);
322 M[6] = (( 0.0030 * 1.910197) + (0.0136 * 0.370950));
323 M[7] = (( 0.0030 * -1.112124) + (0.0136 * 0.629054));
324 M[8] = (( 0.0030 * 0.201908) + (0.0136 * 0.000008) + 0.9834);;
325
326 clr.RGBc[0] = (clr.RGBp[0] * M[0]) + (clr.RGBp[1] * M[1]) + (clr.RGBp[2] * M[2]);
327 clr.RGBc[1] = (clr.RGBp[0] * M[3]) + (clr.RGBp[1] * M[4]) + (clr.RGBp[2] * M[5]);
328 clr.RGBc[2] = (clr.RGBp[0] * M[6]) + (clr.RGBp[1] * M[7]) + (clr.RGBp[2] * M[8]);
329 return clr;
330}
331
332
333static
334CAM02COLOR InverseChromaticAdaptation(CAM02COLOR clr, cmsCIECAM02* pMod)
335{
336 cmsUInt32Number i;
337 for (i = 0; i < 3; i++) {
338 clr.RGB[i] = clr.RGBc[i] /
339 ((pMod->adoptedWhite.XYZ[1] * pMod->D / pMod->adoptedWhite.RGB[i]) + 1.0 - pMod->D);
340 }
341 return clr;
342}
343
344
345static
346CAM02COLOR CAT02toXYZ(CAM02COLOR clr)
347{
348 clr.XYZ[0] = (clr.RGB[0] * 1.096124) + (clr.RGB[1] * -0.278869) + (clr.RGB[2] * 0.182745);
349 clr.XYZ[1] = (clr.RGB[0] * 0.454369) + (clr.RGB[1] * 0.473533) + (clr.RGB[2] * 0.072098);
350 clr.XYZ[2] = (clr.RGB[0] * -0.009628) + (clr.RGB[1] * -0.005698) + (clr.RGB[2] * 1.015326);
351
352 return clr;
353}
354
355
356cmsHANDLE CMSEXPORT cmsCIECAM02Init(cmsContext ContextID, const cmsViewingConditions* pVC)
357{
358 cmsCIECAM02* lpMod;
359
360 _cmsAssert(pVC != NULL);
361
362 if((lpMod = (cmsCIECAM02*) _cmsMallocZero(ContextID, sizeof(cmsCIECAM02))) == NULL) {
363 return NULL;
364 }
365
366 lpMod ->adoptedWhite.XYZ[0] = pVC ->whitePoint.X;
367 lpMod ->adoptedWhite.XYZ[1] = pVC ->whitePoint.Y;
368 lpMod ->adoptedWhite.XYZ[2] = pVC ->whitePoint.Z;
369
370 lpMod -> LA = pVC ->La;
371 lpMod -> Yb = pVC ->Yb;
372 lpMod -> D = pVC ->D_value;
373 lpMod -> surround = pVC ->surround;
374
375 switch (lpMod -> surround) {
376
377
378 case CUTSHEET_SURROUND:
379 lpMod->F = 0.8;
380 lpMod->c = 0.41;
381 lpMod->Nc = 0.8;
382 break;
383
384 case DARK_SURROUND:
385 lpMod -> F = 0.8;
386 lpMod -> c = 0.525;
387 lpMod -> Nc = 0.8;
388 break;
389
390 case DIM_SURROUND:
391 lpMod -> F = 0.9;
392 lpMod -> c = 0.59;
393 lpMod -> Nc = 0.95;
394 break;
395
396 default:
397 // Average surround
398 lpMod -> F = 1.0;
399 lpMod -> c = 0.69;
400 lpMod -> Nc = 1.0;
401 }
402
403 lpMod -> n = compute_n(lpMod);
404 lpMod -> z = compute_z(lpMod);
405 lpMod -> Nbb = computeNbb(lpMod);
406 lpMod -> FL = computeFL(lpMod);
407
408 if (lpMod -> D == D_CALCULATE) {
409 lpMod -> D = computeD(lpMod);
410 }
411
412 lpMod -> Ncb = lpMod -> Nbb;
413
414 lpMod -> adoptedWhite = XYZtoCAT02(lpMod -> adoptedWhite);
415 lpMod -> adoptedWhite = ChromaticAdaptation(lpMod -> adoptedWhite, lpMod);
416 lpMod -> adoptedWhite = CAT02toHPE(lpMod -> adoptedWhite);
417 lpMod -> adoptedWhite = NonlinearCompression(lpMod -> adoptedWhite, lpMod);
418
419 return (cmsHANDLE) lpMod;
420
421}
422
423void CMSEXPORT cmsCIECAM02Done(cmsContext ContextID, cmsHANDLE hModel)
424{
425 cmsCIECAM02* lpMod = (cmsCIECAM02*) hModel;
426
427 if (lpMod) _cmsFree(ContextID, lpMod);
428}
429
430
431void CMSEXPORT cmsCIECAM02Forward(cmsContext ContextID, cmsHANDLE hModel, const cmsCIEXYZ* pIn, cmsJCh* pOut)
432{
433 CAM02COLOR clr;
434 cmsCIECAM02* lpMod = (cmsCIECAM02*) hModel;
435 cmsUNUSED_PARAMETER(ContextID);
436
437 _cmsAssert(lpMod != NULL);
438 _cmsAssert(pIn != NULL);
439 _cmsAssert(pOut != NULL);
440
441 memset(&clr, 0, sizeof(clr));
442
443 clr.XYZ[0] = pIn ->X;
444 clr.XYZ[1] = pIn ->Y;
445 clr.XYZ[2] = pIn ->Z;
446
447 clr = XYZtoCAT02(clr);
448 clr = ChromaticAdaptation(clr, lpMod);
449 clr = CAT02toHPE(clr);
450 clr = NonlinearCompression(clr, lpMod);
451 clr = ComputeCorrelates(clr, lpMod);
452
453 pOut ->J = clr.J;
454 pOut ->C = clr.C;
455 pOut ->h = clr.h;
456}
457
458void CMSEXPORT cmsCIECAM02Reverse(cmsContext ContextID, cmsHANDLE hModel, const cmsJCh* pIn, cmsCIEXYZ* pOut)
459{
460 CAM02COLOR clr;
461 cmsCIECAM02* lpMod = (cmsCIECAM02*) hModel;
462 cmsUNUSED_PARAMETER(ContextID);
463
464 _cmsAssert(lpMod != NULL);
465 _cmsAssert(pIn != NULL);
466 _cmsAssert(pOut != NULL);
467
468 memset(&clr, 0, sizeof(clr));
469
470 clr.J = pIn -> J;
471 clr.C = pIn -> C;
472 clr.h = pIn -> h;
473
474 clr = InverseCorrelates(clr, lpMod);
475 clr = InverseNonlinearity(clr, lpMod);
476 clr = HPEtoCAT02(clr);
477 clr = InverseChromaticAdaptation(clr, lpMod);
478 clr = CAT02toXYZ(clr);
479
480 pOut ->X = clr.XYZ[0];
481 pOut ->Y = clr.XYZ[1];
482 pOut ->Z = clr.XYZ[2];
483}
484