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 // inter PCS conversions XYZ <-> CIE L* a* b*
30 /*
31
32
33 CIE 15:2004 CIELab is defined as:
34
35 L* = 116*f(Y/Yn) - 16 0 <= L* <= 100
36 a* = 500*[f(X/Xn) - f(Y/Yn)]
37 b* = 200*[f(Y/Yn) - f(Z/Zn)]
38
39 and
40
41 f(t) = t^(1/3) 1 >= t > (24/116)^3
42 (841/108)*t + (16/116) 0 <= t <= (24/116)^3
43
44
45 Reverse transform is:
46
47 X = Xn*[a* / 500 + (L* + 16) / 116] ^ 3 if (X/Xn) > (24/116)
48 = Xn*(a* / 500 + L* / 116) / 7.787 if (X/Xn) <= (24/116)
49
50
51
52 PCS in Lab2 is encoded as:
53
54 8 bit Lab PCS:
55
56 L* 0..100 into a 0..ff byte.
57 a* t + 128 range is -128.0 +127.0
58 b*
59
60 16 bit Lab PCS:
61
62 L* 0..100 into a 0..ff00 word.
63 a* t + 128 range is -128.0 +127.9961
64 b*
65
66
67
68 Interchange Space Component Actual Range Encoded Range
69 CIE XYZ X 0 -> 1.99997 0x0000 -> 0xffff
70 CIE XYZ Y 0 -> 1.99997 0x0000 -> 0xffff
71 CIE XYZ Z 0 -> 1.99997 0x0000 -> 0xffff
72
73 Version 2,3
74 -----------
75
76 CIELAB (16 bit) L* 0 -> 100.0 0x0000 -> 0xff00
77 CIELAB (16 bit) a* -128.0 -> +127.996 0x0000 -> 0x8000 -> 0xffff
78 CIELAB (16 bit) b* -128.0 -> +127.996 0x0000 -> 0x8000 -> 0xffff
79
80
81 Version 4
82 ---------
83
84 CIELAB (16 bit) L* 0 -> 100.0 0x0000 -> 0xffff
85 CIELAB (16 bit) a* -128.0 -> +127 0x0000 -> 0x8080 -> 0xffff
86 CIELAB (16 bit) b* -128.0 -> +127 0x0000 -> 0x8080 -> 0xffff
87
88 */
89
90 // Conversions
cmsXYZ2xyY(cmsCIExyY * Dest,const cmsCIEXYZ * Source)91 void CMSEXPORT cmsXYZ2xyY(cmsCIExyY* Dest, const cmsCIEXYZ* Source)
92 {
93 cmsFloat64Number ISum;
94
95 ISum = 1./(Source -> X + Source -> Y + Source -> Z);
96
97 Dest -> x = (Source -> X) * ISum;
98 Dest -> y = (Source -> Y) * ISum;
99 Dest -> Y = Source -> Y;
100 }
101
cmsxyY2XYZ(cmsCIEXYZ * Dest,const cmsCIExyY * Source)102 void CMSEXPORT cmsxyY2XYZ(cmsCIEXYZ* Dest, const cmsCIExyY* Source)
103 {
104 Dest -> X = (Source -> x / Source -> y) * Source -> Y;
105 Dest -> Y = Source -> Y;
106 Dest -> Z = ((1 - Source -> x - Source -> y) / Source -> y) * Source -> Y;
107 }
108
109 /*
110 The break point (24/116)^3 = (6/29)^3 is a very small amount of tristimulus
111 primary (0.008856). Generally, this only happens for
112 nearly ideal blacks and for some orange / amber colors in transmission mode.
113 For example, the Z value of the orange turn indicator lamp lens on an
114 automobile will often be below this value. But the Z does not
115 contribute to the perceived color directly.
116 */
117
118 static
f(cmsFloat64Number t)119 cmsFloat64Number f(cmsFloat64Number t)
120 {
121 const cmsFloat64Number Limit = (24.0/116.0) * (24.0/116.0) * (24.0/116.0);
122
123 if (t <= Limit)
124 return (841.0/108.0) * t + (16.0/116.0);
125 else
126 return pow(t, 1.0/3.0);
127 }
128
129 static
f_1(cmsFloat64Number t)130 cmsFloat64Number f_1(cmsFloat64Number t)
131 {
132 const cmsFloat64Number Limit = (24.0/116.0);
133
134 if (t <= Limit) {
135 return (108.0/841.0) * (t - (16.0/116.0));
136 }
137
138 return t * t * t;
139 }
140
141
142 // Standard XYZ to Lab. it can handle negative XZY numbers in some cases
cmsXYZ2Lab(const cmsCIEXYZ * WhitePoint,cmsCIELab * Lab,const cmsCIEXYZ * xyz)143 void CMSEXPORT cmsXYZ2Lab(const cmsCIEXYZ* WhitePoint, cmsCIELab* Lab, const cmsCIEXYZ* xyz)
144 {
145 cmsFloat64Number fx, fy, fz;
146
147 if (WhitePoint == NULL)
148 WhitePoint = cmsD50_XYZ();
149
150 fx = f(xyz->X / WhitePoint->X);
151 fy = f(xyz->Y / WhitePoint->Y);
152 fz = f(xyz->Z / WhitePoint->Z);
153
154 Lab->L = 116.0*fy - 16.0;
155 Lab->a = 500.0*(fx - fy);
156 Lab->b = 200.0*(fy - fz);
157 }
158
159
160 // Standard XYZ to Lab. It can return negative XYZ in some cases
cmsLab2XYZ(const cmsCIEXYZ * WhitePoint,cmsCIEXYZ * xyz,const cmsCIELab * Lab)161 void CMSEXPORT cmsLab2XYZ(const cmsCIEXYZ* WhitePoint, cmsCIEXYZ* xyz, const cmsCIELab* Lab)
162 {
163 cmsFloat64Number x, y, z;
164
165 if (WhitePoint == NULL)
166 WhitePoint = cmsD50_XYZ();
167
168 y = (Lab-> L + 16.0) / 116.0;
169 x = y + 0.002 * Lab -> a;
170 z = y - 0.005 * Lab -> b;
171
172 xyz -> X = f_1(x) * WhitePoint -> X;
173 xyz -> Y = f_1(y) * WhitePoint -> Y;
174 xyz -> Z = f_1(z) * WhitePoint -> Z;
175
176 }
177
178 static
L2float2(cmsUInt16Number v)179 cmsFloat64Number L2float2(cmsUInt16Number v)
180 {
181 return (cmsFloat64Number) v / 652.800;
182 }
183
184 // the a/b part
185 static
ab2float2(cmsUInt16Number v)186 cmsFloat64Number ab2float2(cmsUInt16Number v)
187 {
188 return ((cmsFloat64Number) v / 256.0) - 128.0;
189 }
190
191 static
L2Fix2(cmsFloat64Number L)192 cmsUInt16Number L2Fix2(cmsFloat64Number L)
193 {
194 return _cmsQuickSaturateWord(L * 652.8);
195 }
196
197 static
ab2Fix2(cmsFloat64Number ab)198 cmsUInt16Number ab2Fix2(cmsFloat64Number ab)
199 {
200 return _cmsQuickSaturateWord((ab + 128.0) * 256.0);
201 }
202
203
204 static
L2float4(cmsUInt16Number v)205 cmsFloat64Number L2float4(cmsUInt16Number v)
206 {
207 return (cmsFloat64Number) v / 655.35;
208 }
209
210 // the a/b part
211 static
ab2float4(cmsUInt16Number v)212 cmsFloat64Number ab2float4(cmsUInt16Number v)
213 {
214 return ((cmsFloat64Number) v / 257.0) - 128.0;
215 }
216
217
cmsLabEncoded2FloatV2(cmsCIELab * Lab,const cmsUInt16Number wLab[3])218 void CMSEXPORT cmsLabEncoded2FloatV2(cmsCIELab* Lab, const cmsUInt16Number wLab[3])
219 {
220 Lab->L = L2float2(wLab[0]);
221 Lab->a = ab2float2(wLab[1]);
222 Lab->b = ab2float2(wLab[2]);
223 }
224
225
cmsLabEncoded2Float(cmsCIELab * Lab,const cmsUInt16Number wLab[3])226 void CMSEXPORT cmsLabEncoded2Float(cmsCIELab* Lab, const cmsUInt16Number wLab[3])
227 {
228 Lab->L = L2float4(wLab[0]);
229 Lab->a = ab2float4(wLab[1]);
230 Lab->b = ab2float4(wLab[2]);
231 }
232
233 static
Clamp_L_doubleV2(cmsFloat64Number L)234 cmsFloat64Number Clamp_L_doubleV2(cmsFloat64Number L)
235 {
236 const cmsFloat64Number L_max = (cmsFloat64Number) (0xFFFF * 100.0) / 0xFF00;
237
238 if (L < 0) L = 0;
239 if (L > L_max) L = L_max;
240
241 return L;
242 }
243
244
245 static
Clamp_ab_doubleV2(cmsFloat64Number ab)246 cmsFloat64Number Clamp_ab_doubleV2(cmsFloat64Number ab)
247 {
248 if (ab < MIN_ENCODEABLE_ab2) ab = MIN_ENCODEABLE_ab2;
249 if (ab > MAX_ENCODEABLE_ab2) ab = MAX_ENCODEABLE_ab2;
250
251 return ab;
252 }
253
cmsFloat2LabEncodedV2(cmsUInt16Number wLab[3],const cmsCIELab * fLab)254 void CMSEXPORT cmsFloat2LabEncodedV2(cmsUInt16Number wLab[3], const cmsCIELab* fLab)
255 {
256 cmsCIELab Lab;
257
258 Lab.L = Clamp_L_doubleV2(fLab ->L);
259 Lab.a = Clamp_ab_doubleV2(fLab ->a);
260 Lab.b = Clamp_ab_doubleV2(fLab ->b);
261
262 wLab[0] = L2Fix2(Lab.L);
263 wLab[1] = ab2Fix2(Lab.a);
264 wLab[2] = ab2Fix2(Lab.b);
265 }
266
267
268 static
Clamp_L_doubleV4(cmsFloat64Number L)269 cmsFloat64Number Clamp_L_doubleV4(cmsFloat64Number L)
270 {
271 if (L < 0) L = 0;
272 if (L > 100.0) L = 100.0;
273
274 return L;
275 }
276
277 static
Clamp_ab_doubleV4(cmsFloat64Number ab)278 cmsFloat64Number Clamp_ab_doubleV4(cmsFloat64Number ab)
279 {
280 if (ab < MIN_ENCODEABLE_ab4) ab = MIN_ENCODEABLE_ab4;
281 if (ab > MAX_ENCODEABLE_ab4) ab = MAX_ENCODEABLE_ab4;
282
283 return ab;
284 }
285
286 static
L2Fix4(cmsFloat64Number L)287 cmsUInt16Number L2Fix4(cmsFloat64Number L)
288 {
289 return _cmsQuickSaturateWord(L * 655.35);
290 }
291
292 static
ab2Fix4(cmsFloat64Number ab)293 cmsUInt16Number ab2Fix4(cmsFloat64Number ab)
294 {
295 return _cmsQuickSaturateWord((ab + 128.0) * 257.0);
296 }
297
cmsFloat2LabEncoded(cmsUInt16Number wLab[3],const cmsCIELab * fLab)298 void CMSEXPORT cmsFloat2LabEncoded(cmsUInt16Number wLab[3], const cmsCIELab* fLab)
299 {
300 cmsCIELab Lab;
301
302 Lab.L = Clamp_L_doubleV4(fLab ->L);
303 Lab.a = Clamp_ab_doubleV4(fLab ->a);
304 Lab.b = Clamp_ab_doubleV4(fLab ->b);
305
306 wLab[0] = L2Fix4(Lab.L);
307 wLab[1] = ab2Fix4(Lab.a);
308 wLab[2] = ab2Fix4(Lab.b);
309 }
310
311 // Auxiliary: convert to Radians
312 static
RADIANS(cmsFloat64Number deg)313 cmsFloat64Number RADIANS(cmsFloat64Number deg)
314 {
315 return (deg * M_PI) / 180.;
316 }
317
318
319 // Auxiliary: atan2 but operating in degrees and returning 0 if a==b==0
320 static
atan2deg(cmsFloat64Number a,cmsFloat64Number b)321 cmsFloat64Number atan2deg(cmsFloat64Number a, cmsFloat64Number b)
322 {
323 cmsFloat64Number h;
324
325 if (a == 0 && b == 0)
326 h = 0;
327 else
328 h = atan2(a, b);
329
330 h *= (180. / M_PI);
331
332 while (h > 360.)
333 h -= 360.;
334
335 while ( h < 0)
336 h += 360.;
337
338 return h;
339 }
340
341
342 // Auxiliary: Square
343 static
Sqr(cmsFloat64Number v)344 cmsFloat64Number Sqr(cmsFloat64Number v)
345 {
346 return v * v;
347 }
348 // From cylindrical coordinates. No check is performed, then negative values are allowed
cmsLab2LCh(cmsCIELCh * LCh,const cmsCIELab * Lab)349 void CMSEXPORT cmsLab2LCh(cmsCIELCh* LCh, const cmsCIELab* Lab)
350 {
351 LCh -> L = Lab -> L;
352 LCh -> C = pow(Sqr(Lab ->a) + Sqr(Lab ->b), 0.5);
353 LCh -> h = atan2deg(Lab ->b, Lab ->a);
354 }
355
356
357 // To cylindrical coordinates. No check is performed, then negative values are allowed
cmsLCh2Lab(cmsCIELab * Lab,const cmsCIELCh * LCh)358 void CMSEXPORT cmsLCh2Lab(cmsCIELab* Lab, const cmsCIELCh* LCh)
359 {
360 cmsFloat64Number h = (LCh -> h * M_PI) / 180.0;
361
362 Lab -> L = LCh -> L;
363 Lab -> a = LCh -> C * cos(h);
364 Lab -> b = LCh -> C * sin(h);
365 }
366
367 // In XYZ All 3 components are encoded using 1.15 fixed point
368 static
XYZ2Fix(cmsFloat64Number d)369 cmsUInt16Number XYZ2Fix(cmsFloat64Number d)
370 {
371 return _cmsQuickSaturateWord(d * 32768.0);
372 }
373
cmsFloat2XYZEncoded(cmsUInt16Number XYZ[3],const cmsCIEXYZ * fXYZ)374 void CMSEXPORT cmsFloat2XYZEncoded(cmsUInt16Number XYZ[3], const cmsCIEXYZ* fXYZ)
375 {
376 cmsCIEXYZ xyz;
377
378 xyz.X = fXYZ -> X;
379 xyz.Y = fXYZ -> Y;
380 xyz.Z = fXYZ -> Z;
381
382 // Clamp to encodeable values.
383 if (xyz.Y <= 0) {
384
385 xyz.X = 0;
386 xyz.Y = 0;
387 xyz.Z = 0;
388 }
389
390 if (xyz.X > MAX_ENCODEABLE_XYZ)
391 xyz.X = MAX_ENCODEABLE_XYZ;
392
393 if (xyz.X < 0)
394 xyz.X = 0;
395
396 if (xyz.Y > MAX_ENCODEABLE_XYZ)
397 xyz.Y = MAX_ENCODEABLE_XYZ;
398
399 if (xyz.Y < 0)
400 xyz.Y = 0;
401
402 if (xyz.Z > MAX_ENCODEABLE_XYZ)
403 xyz.Z = MAX_ENCODEABLE_XYZ;
404
405 if (xyz.Z < 0)
406 xyz.Z = 0;
407
408
409 XYZ[0] = XYZ2Fix(xyz.X);
410 XYZ[1] = XYZ2Fix(xyz.Y);
411 XYZ[2] = XYZ2Fix(xyz.Z);
412 }
413
414
415 // To convert from Fixed 1.15 point to cmsFloat64Number
416 static
XYZ2float(cmsUInt16Number v)417 cmsFloat64Number XYZ2float(cmsUInt16Number v)
418 {
419 cmsS15Fixed16Number fix32;
420
421 // From 1.15 to 15.16
422 fix32 = v << 1;
423
424 // From fixed 15.16 to cmsFloat64Number
425 return _cms15Fixed16toDouble(fix32);
426 }
427
428
cmsXYZEncoded2Float(cmsCIEXYZ * fXYZ,const cmsUInt16Number XYZ[3])429 void CMSEXPORT cmsXYZEncoded2Float(cmsCIEXYZ* fXYZ, const cmsUInt16Number XYZ[3])
430 {
431 fXYZ -> X = XYZ2float(XYZ[0]);
432 fXYZ -> Y = XYZ2float(XYZ[1]);
433 fXYZ -> Z = XYZ2float(XYZ[2]);
434 }
435
436
437 // Returns dE on two Lab values
cmsDeltaE(const cmsCIELab * Lab1,const cmsCIELab * Lab2)438 cmsFloat64Number CMSEXPORT cmsDeltaE(const cmsCIELab* Lab1, const cmsCIELab* Lab2)
439 {
440 cmsFloat64Number dL, da, db;
441
442 dL = fabs(Lab1 -> L - Lab2 -> L);
443 da = fabs(Lab1 -> a - Lab2 -> a);
444 db = fabs(Lab1 -> b - Lab2 -> b);
445
446 return pow(Sqr(dL) + Sqr(da) + Sqr(db), 0.5);
447 }
448
449
450 // Return the CIE94 Delta E
cmsCIE94DeltaE(const cmsCIELab * Lab1,const cmsCIELab * Lab2)451 cmsFloat64Number CMSEXPORT cmsCIE94DeltaE(const cmsCIELab* Lab1, const cmsCIELab* Lab2)
452 {
453 cmsCIELCh LCh1, LCh2;
454 cmsFloat64Number dE, dL, dC, dh, dhsq;
455 cmsFloat64Number c12, sc, sh;
456
457 dL = fabs(Lab1 ->L - Lab2 ->L);
458
459 cmsLab2LCh(&LCh1, Lab1);
460 cmsLab2LCh(&LCh2, Lab2);
461
462 dC = fabs(LCh1.C - LCh2.C);
463 dE = cmsDeltaE(Lab1, Lab2);
464
465 dhsq = Sqr(dE) - Sqr(dL) - Sqr(dC);
466 if (dhsq < 0)
467 dh = 0;
468 else
469 dh = pow(dhsq, 0.5);
470
471 c12 = sqrt(LCh1.C * LCh2.C);
472
473 sc = 1.0 + (0.048 * c12);
474 sh = 1.0 + (0.014 * c12);
475
476 return sqrt(Sqr(dL) + Sqr(dC) / Sqr(sc) + Sqr(dh) / Sqr(sh));
477 }
478
479
480 // Auxiliary
481 static
ComputeLBFD(const cmsCIELab * Lab)482 cmsFloat64Number ComputeLBFD(const cmsCIELab* Lab)
483 {
484 cmsFloat64Number yt;
485
486 if (Lab->L > 7.996969)
487 yt = (Sqr((Lab->L+16)/116)*((Lab->L+16)/116))*100;
488 else
489 yt = 100 * (Lab->L / 903.3);
490
491 return (54.6 * (M_LOG10E * (log(yt + 1.5))) - 9.6);
492 }
493
494
495
496 // bfd - gets BFD(1:1) difference between Lab1, Lab2
cmsBFDdeltaE(const cmsCIELab * Lab1,const cmsCIELab * Lab2)497 cmsFloat64Number CMSEXPORT cmsBFDdeltaE(const cmsCIELab* Lab1, const cmsCIELab* Lab2)
498 {
499 cmsFloat64Number lbfd1,lbfd2,AveC,Aveh,dE,deltaL,
500 deltaC,deltah,dc,t,g,dh,rh,rc,rt,bfd;
501 cmsCIELCh LCh1, LCh2;
502
503
504 lbfd1 = ComputeLBFD(Lab1);
505 lbfd2 = ComputeLBFD(Lab2);
506 deltaL = lbfd2 - lbfd1;
507
508 cmsLab2LCh(&LCh1, Lab1);
509 cmsLab2LCh(&LCh2, Lab2);
510
511 deltaC = LCh2.C - LCh1.C;
512 AveC = (LCh1.C+LCh2.C)/2;
513 Aveh = (LCh1.h+LCh2.h)/2;
514
515 dE = cmsDeltaE(Lab1, Lab2);
516
517 if (Sqr(dE)>(Sqr(Lab2->L-Lab1->L)+Sqr(deltaC)))
518 deltah = sqrt(Sqr(dE)-Sqr(Lab2->L-Lab1->L)-Sqr(deltaC));
519 else
520 deltah =0;
521
522
523 dc = 0.035 * AveC / (1 + 0.00365 * AveC)+0.521;
524 g = sqrt(Sqr(Sqr(AveC))/(Sqr(Sqr(AveC))+14000));
525 t = 0.627+(0.055*cos((Aveh-254)/(180/M_PI))-
526 0.040*cos((2*Aveh-136)/(180/M_PI))+
527 0.070*cos((3*Aveh-31)/(180/M_PI))+
528 0.049*cos((4*Aveh+114)/(180/M_PI))-
529 0.015*cos((5*Aveh-103)/(180/M_PI)));
530
531 dh = dc*(g*t+1-g);
532 rh = -0.260*cos((Aveh-308)/(180/M_PI))-
533 0.379*cos((2*Aveh-160)/(180/M_PI))-
534 0.636*cos((3*Aveh+254)/(180/M_PI))+
535 0.226*cos((4*Aveh+140)/(180/M_PI))-
536 0.194*cos((5*Aveh+280)/(180/M_PI));
537
538 rc = sqrt((AveC*AveC*AveC*AveC*AveC*AveC)/((AveC*AveC*AveC*AveC*AveC*AveC)+70000000));
539 rt = rh*rc;
540
541 bfd = sqrt(Sqr(deltaL)+Sqr(deltaC/dc)+Sqr(deltah/dh)+(rt*(deltaC/dc)*(deltah/dh)));
542
543 return bfd;
544 }
545
546
547 // cmc - CMC(l:c) difference between Lab1, Lab2
cmsCMCdeltaE(const cmsCIELab * Lab1,const cmsCIELab * Lab2,cmsFloat64Number l,cmsFloat64Number c)548 cmsFloat64Number CMSEXPORT cmsCMCdeltaE(const cmsCIELab* Lab1, const cmsCIELab* Lab2, cmsFloat64Number l, cmsFloat64Number c)
549 {
550 cmsFloat64Number dE,dL,dC,dh,sl,sc,sh,t,f,cmc;
551 cmsCIELCh LCh1, LCh2;
552
553 if (Lab1 ->L == 0 && Lab2 ->L == 0) return 0;
554
555 cmsLab2LCh(&LCh1, Lab1);
556 cmsLab2LCh(&LCh2, Lab2);
557
558
559 dL = Lab2->L-Lab1->L;
560 dC = LCh2.C-LCh1.C;
561
562 dE = cmsDeltaE(Lab1, Lab2);
563
564 if (Sqr(dE)>(Sqr(dL)+Sqr(dC)))
565 dh = sqrt(Sqr(dE)-Sqr(dL)-Sqr(dC));
566 else
567 dh =0;
568
569 if ((LCh1.h > 164) && (LCh1.h < 345))
570 t = 0.56 + fabs(0.2 * cos(((LCh1.h + 168)/(180/M_PI))));
571 else
572 t = 0.36 + fabs(0.4 * cos(((LCh1.h + 35 )/(180/M_PI))));
573
574 sc = 0.0638 * LCh1.C / (1 + 0.0131 * LCh1.C) + 0.638;
575 sl = 0.040975 * Lab1->L /(1 + 0.01765 * Lab1->L);
576
577 if (Lab1->L<16)
578 sl = 0.511;
579
580 f = sqrt((LCh1.C * LCh1.C * LCh1.C * LCh1.C)/((LCh1.C * LCh1.C * LCh1.C * LCh1.C)+1900));
581 sh = sc*(t*f+1-f);
582 cmc = sqrt(Sqr(dL/(l*sl))+Sqr(dC/(c*sc))+Sqr(dh/sh));
583
584 return cmc;
585 }
586
587 // dE2000 The weightings KL, KC and KH can be modified to reflect the relative
588 // importance of lightness, chroma and hue in different industrial applications
cmsCIE2000DeltaE(const cmsCIELab * Lab1,const cmsCIELab * Lab2,cmsFloat64Number Kl,cmsFloat64Number Kc,cmsFloat64Number Kh)589 cmsFloat64Number CMSEXPORT cmsCIE2000DeltaE(const cmsCIELab* Lab1, const cmsCIELab* Lab2,
590 cmsFloat64Number Kl, cmsFloat64Number Kc, cmsFloat64Number Kh)
591 {
592 cmsFloat64Number L1 = Lab1->L;
593 cmsFloat64Number a1 = Lab1->a;
594 cmsFloat64Number b1 = Lab1->b;
595 cmsFloat64Number C = sqrt( Sqr(a1) + Sqr(b1) );
596
597 cmsFloat64Number Ls = Lab2 ->L;
598 cmsFloat64Number as = Lab2 ->a;
599 cmsFloat64Number bs = Lab2 ->b;
600 cmsFloat64Number Cs = sqrt( Sqr(as) + Sqr(bs) );
601
602 cmsFloat64Number G = 0.5 * ( 1 - sqrt(pow((C + Cs) / 2 , 7.0) / (pow((C + Cs) / 2, 7.0) + pow(25.0, 7.0) ) ));
603
604 cmsFloat64Number a_p = (1 + G ) * a1;
605 cmsFloat64Number b_p = b1;
606 cmsFloat64Number C_p = sqrt( Sqr(a_p) + Sqr(b_p));
607 cmsFloat64Number h_p = atan2deg(b_p, a_p);
608
609
610 cmsFloat64Number a_ps = (1 + G) * as;
611 cmsFloat64Number b_ps = bs;
612 cmsFloat64Number C_ps = sqrt(Sqr(a_ps) + Sqr(b_ps));
613 cmsFloat64Number h_ps = atan2deg(b_ps, a_ps);
614
615 cmsFloat64Number meanC_p =(C_p + C_ps) / 2;
616
617 cmsFloat64Number hps_plus_hp = h_ps + h_p;
618 cmsFloat64Number hps_minus_hp = h_ps - h_p;
619
620 cmsFloat64Number meanh_p = fabs(hps_minus_hp) <= 180.000001 ? (hps_plus_hp)/2 :
621 (hps_plus_hp) < 360 ? (hps_plus_hp + 360)/2 :
622 (hps_plus_hp - 360)/2;
623
624 cmsFloat64Number delta_h = (hps_minus_hp) <= -180.000001 ? (hps_minus_hp + 360) :
625 (hps_minus_hp) > 180 ? (hps_minus_hp - 360) :
626 (hps_minus_hp);
627 cmsFloat64Number delta_L = (Ls - L1);
628 cmsFloat64Number delta_C = (C_ps - C_p );
629
630
631 cmsFloat64Number delta_H =2 * sqrt(C_ps*C_p) * sin(RADIANS(delta_h) / 2);
632
633 cmsFloat64Number T = 1 - 0.17 * cos(RADIANS(meanh_p-30))
634 + 0.24 * cos(RADIANS(2*meanh_p))
635 + 0.32 * cos(RADIANS(3*meanh_p + 6))
636 - 0.2 * cos(RADIANS(4*meanh_p - 63));
637
638 cmsFloat64Number Sl = 1 + (0.015 * Sqr((Ls + L1) /2- 50) )/ sqrt(20 + Sqr( (Ls+L1)/2 - 50) );
639
640 cmsFloat64Number Sc = 1 + 0.045 * (C_p + C_ps)/2;
641 cmsFloat64Number Sh = 1 + 0.015 * ((C_ps + C_p)/2) * T;
642
643 cmsFloat64Number delta_ro = 30 * exp( -Sqr(((meanh_p - 275 ) / 25)));
644
645 cmsFloat64Number Rc = 2 * sqrt(( pow(meanC_p, 7.0) )/( pow(meanC_p, 7.0) + pow(25.0, 7.0)));
646
647 cmsFloat64Number Rt = -sin(2 * RADIANS(delta_ro)) * Rc;
648
649 cmsFloat64Number deltaE00 = sqrt( Sqr(delta_L /(Sl * Kl)) +
650 Sqr(delta_C/(Sc * Kc)) +
651 Sqr(delta_H/(Sh * Kh)) +
652 Rt*(delta_C/(Sc * Kc)) * (delta_H / (Sh * Kh)));
653
654 return deltaE00;
655 }
656
657 // This function returns a number of gridpoints to be used as LUT table. It assumes same number
658 // of gripdpoints in all dimensions. Flags may override the choice.
_cmsReasonableGridpointsByColorspace(cmsColorSpaceSignature Colorspace,cmsUInt32Number dwFlags)659 cmsUInt32Number _cmsReasonableGridpointsByColorspace(cmsColorSpaceSignature Colorspace, cmsUInt32Number dwFlags)
660 {
661 cmsUInt32Number nChannels;
662
663 // Already specified?
664 if (dwFlags & 0x00FF0000) {
665 // Yes, grab'em
666 return (dwFlags >> 16) & 0xFF;
667 }
668
669 nChannels = cmsChannelsOf(Colorspace);
670
671 // HighResPrecalc is maximum resolution
672 if (dwFlags & cmsFLAGS_HIGHRESPRECALC) {
673
674 if (nChannels > 4)
675 return 7; // 7 for Hifi
676
677 if (nChannels == 4) // 23 for CMYK
678 return 23;
679
680 return 49; // 49 for RGB and others
681 }
682
683
684 // LowResPrecal is lower resolution
685 if (dwFlags & cmsFLAGS_LOWRESPRECALC) {
686
687 if (nChannels > 4)
688 return 6; // 6 for more than 4 channels
689
690 if (nChannels == 1)
691 return 33; // For monochrome
692
693 return 17; // 17 for remaining
694 }
695
696 // Default values
697 if (nChannels > 4)
698 return 7; // 7 for Hifi
699
700 if (nChannels == 4)
701 return 17; // 17 for CMYK
702
703 return 33; // 33 for RGB
704 }
705
706
_cmsEndPointsBySpace(cmsColorSpaceSignature Space,cmsUInt16Number ** White,cmsUInt16Number ** Black,cmsUInt32Number * nOutputs)707 cmsBool _cmsEndPointsBySpace(cmsColorSpaceSignature Space,
708 cmsUInt16Number **White,
709 cmsUInt16Number **Black,
710 cmsUInt32Number *nOutputs)
711 {
712 // Only most common spaces
713
714 static cmsUInt16Number RGBblack[4] = { 0, 0, 0 };
715 static cmsUInt16Number RGBwhite[4] = { 0xffff, 0xffff, 0xffff };
716 static cmsUInt16Number CMYKblack[4] = { 0xffff, 0xffff, 0xffff, 0xffff }; // 400% of ink
717 static cmsUInt16Number CMYKwhite[4] = { 0, 0, 0, 0 };
718 static cmsUInt16Number LABblack[4] = { 0, 0x8080, 0x8080 }; // V4 Lab encoding
719 static cmsUInt16Number LABwhite[4] = { 0xFFFF, 0x8080, 0x8080 };
720 static cmsUInt16Number CMYblack[4] = { 0xffff, 0xffff, 0xffff };
721 static cmsUInt16Number CMYwhite[4] = { 0, 0, 0 };
722 static cmsUInt16Number Grayblack[4] = { 0 };
723 static cmsUInt16Number GrayWhite[4] = { 0xffff };
724
725 switch (Space) {
726
727 case cmsSigGrayData: if (White) *White = GrayWhite;
728 if (Black) *Black = Grayblack;
729 if (nOutputs) *nOutputs = 1;
730 return TRUE;
731
732 case cmsSigRgbData: if (White) *White = RGBwhite;
733 if (Black) *Black = RGBblack;
734 if (nOutputs) *nOutputs = 3;
735 return TRUE;
736
737 case cmsSigLabData: if (White) *White = LABwhite;
738 if (Black) *Black = LABblack;
739 if (nOutputs) *nOutputs = 3;
740 return TRUE;
741
742 case cmsSigCmykData: if (White) *White = CMYKwhite;
743 if (Black) *Black = CMYKblack;
744 if (nOutputs) *nOutputs = 4;
745 return TRUE;
746
747 case cmsSigCmyData: if (White) *White = CMYwhite;
748 if (Black) *Black = CMYblack;
749 if (nOutputs) *nOutputs = 3;
750 return TRUE;
751
752 default:;
753 }
754
755 return FALSE;
756 }
757
758
759
760 // Several utilities -------------------------------------------------------
761
762 // Translate from our colorspace to ICC representation
763
_cmsICCcolorSpace(int OurNotation)764 cmsColorSpaceSignature CMSEXPORT _cmsICCcolorSpace(int OurNotation)
765 {
766 switch (OurNotation) {
767
768 case 1:
769 case PT_GRAY: return cmsSigGrayData;
770
771 case 2:
772 case PT_RGB: return cmsSigRgbData;
773
774 case PT_CMY: return cmsSigCmyData;
775 case PT_CMYK: return cmsSigCmykData;
776 case PT_YCbCr:return cmsSigYCbCrData;
777 case PT_YUV: return cmsSigLuvData;
778 case PT_XYZ: return cmsSigXYZData;
779
780 case PT_LabV2:
781 case PT_Lab: return cmsSigLabData;
782
783 case PT_YUVK: return cmsSigLuvKData;
784 case PT_HSV: return cmsSigHsvData;
785 case PT_HLS: return cmsSigHlsData;
786 case PT_Yxy: return cmsSigYxyData;
787
788 case PT_MCH1: return cmsSigMCH1Data;
789 case PT_MCH2: return cmsSigMCH2Data;
790 case PT_MCH3: return cmsSigMCH3Data;
791 case PT_MCH4: return cmsSigMCH4Data;
792 case PT_MCH5: return cmsSigMCH5Data;
793 case PT_MCH6: return cmsSigMCH6Data;
794 case PT_MCH7: return cmsSigMCH7Data;
795 case PT_MCH8: return cmsSigMCH8Data;
796
797 case PT_MCH9: return cmsSigMCH9Data;
798 case PT_MCH10: return cmsSigMCHAData;
799 case PT_MCH11: return cmsSigMCHBData;
800 case PT_MCH12: return cmsSigMCHCData;
801 case PT_MCH13: return cmsSigMCHDData;
802 case PT_MCH14: return cmsSigMCHEData;
803 case PT_MCH15: return cmsSigMCHFData;
804
805 default: return (cmsColorSpaceSignature) 0;
806 }
807 }
808
809
_cmsLCMScolorSpace(cmsColorSpaceSignature ProfileSpace)810 int CMSEXPORT _cmsLCMScolorSpace(cmsColorSpaceSignature ProfileSpace)
811 {
812 switch (ProfileSpace) {
813
814 case cmsSigGrayData: return PT_GRAY;
815 case cmsSigRgbData: return PT_RGB;
816 case cmsSigCmyData: return PT_CMY;
817 case cmsSigCmykData: return PT_CMYK;
818 case cmsSigYCbCrData:return PT_YCbCr;
819 case cmsSigLuvData: return PT_YUV;
820 case cmsSigXYZData: return PT_XYZ;
821 case cmsSigLabData: return PT_Lab;
822 case cmsSigLuvKData: return PT_YUVK;
823 case cmsSigHsvData: return PT_HSV;
824 case cmsSigHlsData: return PT_HLS;
825 case cmsSigYxyData: return PT_Yxy;
826
827 case cmsSig1colorData:
828 case cmsSigMCH1Data: return PT_MCH1;
829
830 case cmsSig2colorData:
831 case cmsSigMCH2Data: return PT_MCH2;
832
833 case cmsSig3colorData:
834 case cmsSigMCH3Data: return PT_MCH3;
835
836 case cmsSig4colorData:
837 case cmsSigMCH4Data: return PT_MCH4;
838
839 case cmsSig5colorData:
840 case cmsSigMCH5Data: return PT_MCH5;
841
842 case cmsSig6colorData:
843 case cmsSigMCH6Data: return PT_MCH6;
844
845 case cmsSigMCH7Data:
846 case cmsSig7colorData:return PT_MCH7;
847
848 case cmsSigMCH8Data:
849 case cmsSig8colorData:return PT_MCH8;
850
851 case cmsSigMCH9Data:
852 case cmsSig9colorData:return PT_MCH9;
853
854 case cmsSigMCHAData:
855 case cmsSig10colorData:return PT_MCH10;
856
857 case cmsSigMCHBData:
858 case cmsSig11colorData:return PT_MCH11;
859
860 case cmsSigMCHCData:
861 case cmsSig12colorData:return PT_MCH12;
862
863 case cmsSigMCHDData:
864 case cmsSig13colorData:return PT_MCH13;
865
866 case cmsSigMCHEData:
867 case cmsSig14colorData:return PT_MCH14;
868
869 case cmsSigMCHFData:
870 case cmsSig15colorData:return PT_MCH15;
871
872 default: return (cmsColorSpaceSignature) 0;
873 }
874 }
875
876
cmsChannelsOf(cmsColorSpaceSignature ColorSpace)877 cmsUInt32Number CMSEXPORT cmsChannelsOf(cmsColorSpaceSignature ColorSpace)
878 {
879 switch (ColorSpace) {
880
881 case cmsSigMCH1Data:
882 case cmsSig1colorData:
883 case cmsSigGrayData: return 1;
884
885 case cmsSigMCH2Data:
886 case cmsSig2colorData: return 2;
887
888 case cmsSigXYZData:
889 case cmsSigLabData:
890 case cmsSigLuvData:
891 case cmsSigYCbCrData:
892 case cmsSigYxyData:
893 case cmsSigRgbData:
894 case cmsSigHsvData:
895 case cmsSigHlsData:
896 case cmsSigCmyData:
897 case cmsSigMCH3Data:
898 case cmsSig3colorData: return 3;
899
900 case cmsSigLuvKData:
901 case cmsSigCmykData:
902 case cmsSigMCH4Data:
903 case cmsSig4colorData: return 4;
904
905 case cmsSigMCH5Data:
906 case cmsSig5colorData: return 5;
907
908 case cmsSigMCH6Data:
909 case cmsSig6colorData: return 6;
910
911 case cmsSigMCH7Data:
912 case cmsSig7colorData: return 7;
913
914 case cmsSigMCH8Data:
915 case cmsSig8colorData: return 8;
916
917 case cmsSigMCH9Data:
918 case cmsSig9colorData: return 9;
919
920 case cmsSigMCHAData:
921 case cmsSig10colorData: return 10;
922
923 case cmsSigMCHBData:
924 case cmsSig11colorData: return 11;
925
926 case cmsSigMCHCData:
927 case cmsSig12colorData: return 12;
928
929 case cmsSigMCHDData:
930 case cmsSig13colorData: return 13;
931
932 case cmsSigMCHEData:
933 case cmsSig14colorData: return 14;
934
935 case cmsSigMCHFData:
936 case cmsSig15colorData: return 15;
937
938 default: return 3;
939 }
940 }
941