1 /* 2 * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER. 3 * 4 * This code is free software; you can redistribute it and/or modify it 5 * under the terms of the GNU General Public License version 2 only, as 6 * published by the Free Software Foundation. Oracle designates this 7 * particular file as subject to the "Classpath" exception as provided 8 * by Oracle in the LICENSE file that accompanied this code. 9 * 10 * This code is distributed in the hope that it will be useful, but WITHOUT 11 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 12 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 13 * version 2 for more details (a copy is included in the LICENSE file that 14 * accompanied this code). 15 * 16 * You should have received a copy of the GNU General Public License version 17 * 2 along with this work; if not, write to the Free Software Foundation, 18 * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. 19 * 20 * Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA 21 * or visit www.oracle.com if you need additional information or have any 22 * questions. 23 */ 24 25 // This file is available under and governed by the GNU General Public 26 // License version 2 only, as published by the Free Software Foundation. 27 // However, the following notice accompanied the original version of this 28 // file: 29 // 30 //--------------------------------------------------------------------------------- 31 // 32 // Little Color Management System 33 // Copyright (c) 1998-2010 Marti Maria Saguer 34 // 35 // Permission is hereby granted, free of charge, to any person obtaining 36 // a copy of this software and associated documentation files (the "Software"), 37 // to deal in the Software without restriction, including without limitation 38 // the rights to use, copy, modify, merge, publish, distribute, sublicense, 39 // and/or sell copies of the Software, and to permit persons to whom the Software 40 // is furnished to do so, subject to the following conditions: 41 // 42 // The above copyright notice and this permission notice shall be included in 43 // all copies or substantial portions of the Software. 44 // 45 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, 46 // EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO 47 // THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND 48 // NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE 49 // LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION 50 // OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION 51 // WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. 52 // 53 //--------------------------------------------------------------------------------- 54 // 55 56 #include "lcms2_internal.h" 57 58 59 60 // This file contains routines for resampling and LUT optimization, black point detection 61 // and black preservation. 62 63 // Black point detection ------------------------------------------------------------------------- 64 65 66 // PCS -> PCS round trip transform, always uses relative intent on the device -> pcs 67 static 68 cmsHTRANSFORM CreateRoundtripXForm(cmsHPROFILE hProfile, cmsUInt32Number nIntent) 69 { 70 cmsHPROFILE hLab = cmsCreateLab4Profile(NULL); 71 cmsHTRANSFORM xform; 72 cmsBool BPC[4] = { FALSE, FALSE, FALSE, FALSE }; 73 cmsFloat64Number States[4] = { 1.0, 1.0, 1.0, 1.0 }; 74 cmsHPROFILE hProfiles[4]; 75 cmsUInt32Number Intents[4]; 76 cmsContext ContextID = cmsGetProfileContextID(hProfile); 77 78 hProfiles[0] = hLab; hProfiles[1] = hProfile; hProfiles[2] = hProfile; hProfiles[3] = hLab; 79 Intents[0] = INTENT_RELATIVE_COLORIMETRIC; Intents[1] = nIntent; Intents[2] = INTENT_RELATIVE_COLORIMETRIC; Intents[3] = INTENT_RELATIVE_COLORIMETRIC; 80 81 xform = cmsCreateExtendedTransform(ContextID, 4, hProfiles, BPC, Intents, 82 States, NULL, 0, TYPE_Lab_DBL, TYPE_Lab_DBL, cmsFLAGS_NOCACHE|cmsFLAGS_NOOPTIMIZE); 83 84 cmsCloseProfile(hLab); 85 return xform; 86 } 87 88 // Use darker colorants to obtain black point. This works in the relative colorimetric intent and 89 // assumes more ink results in darker colors. No ink limit is assumed. 90 static 91 cmsBool BlackPointAsDarkerColorant(cmsHPROFILE hInput, 92 cmsUInt32Number Intent, 93 cmsCIEXYZ* BlackPoint, 94 cmsUInt32Number dwFlags) 95 { 96 cmsUInt16Number *Black; 97 cmsHTRANSFORM xform; 98 cmsColorSpaceSignature Space; 99 cmsUInt32Number nChannels; 100 cmsUInt32Number dwFormat; 101 cmsHPROFILE hLab; 102 cmsCIELab Lab; 103 cmsCIEXYZ BlackXYZ; 104 cmsContext ContextID = cmsGetProfileContextID(hInput); 105 106 // If the profile does not support input direction, assume Black point 0 107 if (!cmsIsIntentSupported(hInput, Intent, LCMS_USED_AS_INPUT)) { 108 109 BlackPoint -> X = BlackPoint ->Y = BlackPoint -> Z = 0.0; 110 return FALSE; 111 } 112 113 // Create a formatter which has n channels and floating point 114 dwFormat = cmsFormatterForColorspaceOfProfile(hInput, 2, FALSE); 115 116 // Try to get black by using black colorant 117 Space = cmsGetColorSpace(hInput); 118 119 // This function returns darker colorant in 16 bits for several spaces 120 if (!_cmsEndPointsBySpace(Space, NULL, &Black, &nChannels)) { 121 122 BlackPoint -> X = BlackPoint ->Y = BlackPoint -> Z = 0.0; 123 return FALSE; 124 } 125 126 if (nChannels != T_CHANNELS(dwFormat)) { 127 BlackPoint -> X = BlackPoint ->Y = BlackPoint -> Z = 0.0; 128 return FALSE; 129 } 130 131 // Lab will be used as the output space, but lab2 will avoid recursion 132 hLab = cmsCreateLab2ProfileTHR(ContextID, NULL); 133 if (hLab == NULL) { 134 BlackPoint -> X = BlackPoint ->Y = BlackPoint -> Z = 0.0; 135 return FALSE; 136 } 137 138 // Create the transform 139 xform = cmsCreateTransformTHR(ContextID, hInput, dwFormat, 140 hLab, TYPE_Lab_DBL, Intent, cmsFLAGS_NOOPTIMIZE|cmsFLAGS_NOCACHE); 141 cmsCloseProfile(hLab); 142 143 if (xform == NULL) { 144 // Something went wrong. Get rid of open resources and return zero as black 145 146 BlackPoint -> X = BlackPoint ->Y = BlackPoint -> Z = 0.0; 147 return FALSE; 148 } 149 150 // Convert black to Lab 151 cmsDoTransform(xform, Black, &Lab, 1); 152 153 // Force it to be neutral, clip to max. L* of 50 154 Lab.a = Lab.b = 0; 155 if (Lab.L > 50) Lab.L = 50; 156 157 // Free the resources 158 cmsDeleteTransform(xform); 159 160 // Convert from Lab (which is now clipped) to XYZ. 161 cmsLab2XYZ(NULL, &BlackXYZ, &Lab); 162 163 if (BlackPoint != NULL) 164 *BlackPoint = BlackXYZ; 165 166 return TRUE; 167 168 cmsUNUSED_PARAMETER(dwFlags); 169 } 170 171 // Get a black point of output CMYK profile, discounting any ink-limiting embedded 172 // in the profile. For doing that, we use perceptual intent in input direction: 173 // Lab (0, 0, 0) -> [Perceptual] Profile -> CMYK -> [Rel. colorimetric] Profile -> Lab 174 static 175 cmsBool BlackPointUsingPerceptualBlack(cmsCIEXYZ* BlackPoint, cmsHPROFILE hProfile) 176 177 { 178 cmsHTRANSFORM hRoundTrip; 179 cmsCIELab LabIn, LabOut; 180 cmsCIEXYZ BlackXYZ; 181 182 // Is the intent supported by the profile? 183 if (!cmsIsIntentSupported(hProfile, INTENT_PERCEPTUAL, LCMS_USED_AS_INPUT)) { 184 185 BlackPoint -> X = BlackPoint ->Y = BlackPoint -> Z = 0.0; 186 return TRUE; 187 } 188 189 hRoundTrip = CreateRoundtripXForm(hProfile, INTENT_PERCEPTUAL); 190 if (hRoundTrip == NULL) { 191 BlackPoint -> X = BlackPoint ->Y = BlackPoint -> Z = 0.0; 192 return FALSE; 193 } 194 195 LabIn.L = LabIn.a = LabIn.b = 0; 196 cmsDoTransform(hRoundTrip, &LabIn, &LabOut, 1); 197 198 // Clip Lab to reasonable limits 199 if (LabOut.L > 50) LabOut.L = 50; 200 LabOut.a = LabOut.b = 0; 201 202 cmsDeleteTransform(hRoundTrip); 203 204 // Convert it to XYZ 205 cmsLab2XYZ(NULL, &BlackXYZ, &LabOut); 206 207 if (BlackPoint != NULL) 208 *BlackPoint = BlackXYZ; 209 210 return TRUE; 211 } 212 213 // This function shouldn't exist at all -- there is such quantity of broken 214 // profiles on black point tag, that we must somehow fix chromaticity to 215 // avoid huge tint when doing Black point compensation. This function does 216 // just that. There is a special flag for using black point tag, but turned 217 // off by default because it is bogus on most profiles. The detection algorithm 218 // involves to turn BP to neutral and to use only L component. 219 cmsBool CMSEXPORT cmsDetectBlackPoint(cmsCIEXYZ* BlackPoint, cmsHPROFILE hProfile, cmsUInt32Number Intent, cmsUInt32Number dwFlags) 220 { 221 222 // Zero for black point 223 if (cmsGetDeviceClass(hProfile) == cmsSigLinkClass) { 224 225 BlackPoint -> X = BlackPoint ->Y = BlackPoint -> Z = 0.0; 226 return FALSE; 227 } 228 229 // v4 + perceptual & saturation intents does have its own black point, and it is 230 // well specified enough to use it. Black point tag is deprecated in V4. 231 232 if ((cmsGetEncodedICCversion(hProfile) >= 0x4000000) && 233 (Intent == INTENT_PERCEPTUAL || Intent == INTENT_SATURATION)) { 234 235 // Matrix shaper share MRC & perceptual intents 236 if (cmsIsMatrixShaper(hProfile)) 237 return BlackPointAsDarkerColorant(hProfile, INTENT_RELATIVE_COLORIMETRIC, BlackPoint, 0); 238 239 // Get Perceptual black out of v4 profiles. That is fixed for perceptual & saturation intents 240 BlackPoint -> X = cmsPERCEPTUAL_BLACK_X; 241 BlackPoint -> Y = cmsPERCEPTUAL_BLACK_Y; 242 BlackPoint -> Z = cmsPERCEPTUAL_BLACK_Z; 243 244 return TRUE; 245 } 246 247 248 #ifdef CMS_USE_PROFILE_BLACK_POINT_TAG 249 250 // v2, v4 rel/abs colorimetric 251 if (cmsIsTag(hProfile, cmsSigMediaBlackPointTag) && 252 Intent == INTENT_RELATIVE_COLORIMETRIC) { 253 254 cmsCIEXYZ *BlackPtr, BlackXYZ, UntrustedBlackPoint, TrustedBlackPoint, MediaWhite; 255 cmsCIELab Lab; 256 257 // If black point is specified, then use it, 258 259 BlackPtr = cmsReadTag(hProfile, cmsSigMediaBlackPointTag); 260 if (BlackPtr != NULL) { 261 262 BlackXYZ = *BlackPtr; 263 _cmsReadMediaWhitePoint(&MediaWhite, hProfile); 264 265 // Black point is absolute XYZ, so adapt to D50 to get PCS value 266 cmsAdaptToIlluminant(&UntrustedBlackPoint, &MediaWhite, cmsD50_XYZ(), &BlackXYZ); 267 268 // Force a=b=0 to get rid of any chroma 269 cmsXYZ2Lab(NULL, &Lab, &UntrustedBlackPoint); 270 Lab.a = Lab.b = 0; 271 if (Lab.L > 50) Lab.L = 50; // Clip to L* <= 50 272 cmsLab2XYZ(NULL, &TrustedBlackPoint, &Lab); 273 274 if (BlackPoint != NULL) 275 *BlackPoint = TrustedBlackPoint; 276 277 return TRUE; 278 } 279 } 280 #endif 281 282 // That is about v2 profiles. 283 284 // If output profile, discount ink-limiting and that's all 285 if (Intent == INTENT_RELATIVE_COLORIMETRIC && 286 (cmsGetDeviceClass(hProfile) == cmsSigOutputClass) && 287 (cmsGetColorSpace(hProfile) == cmsSigCmykData)) 288 return BlackPointUsingPerceptualBlack(BlackPoint, hProfile); 289 290 // Nope, compute BP using current intent. 291 return BlackPointAsDarkerColorant(hProfile, Intent, BlackPoint, dwFlags); 292 } 293 294 295 296 // --------------------------------------------------------------------------------------------------------- 297 298 // Least Squares Fit of a Quadratic Curve to Data 299 // http://www.personal.psu.edu/jhm/f90/lectures/lsq2.html 300 301 static 302 cmsFloat64Number RootOfLeastSquaresFitQuadraticCurve(int n, cmsFloat64Number x[], cmsFloat64Number y[]) 303 { 304 double sum_x = 0, sum_x2 = 0, sum_x3 = 0, sum_x4 = 0; 305 double sum_y = 0, sum_yx = 0, sum_yx2 = 0; 306 double disc; 307 int i; 308 cmsMAT3 m; 309 cmsVEC3 v, res; 310 311 if (n < 4) return 0; 312 313 for (i=0; i < n; i++) { 314 315 double xn = x[i]; 316 double yn = y[i]; 317 318 sum_x += xn; 319 sum_x2 += xn*xn; 320 sum_x3 += xn*xn*xn; 321 sum_x4 += xn*xn*xn*xn; 322 323 sum_y += yn; 324 sum_yx += yn*xn; 325 sum_yx2 += yn*xn*xn; 326 } 327 328 _cmsVEC3init(&m.v[0], n, sum_x, sum_x2); 329 _cmsVEC3init(&m.v[1], sum_x, sum_x2, sum_x3); 330 _cmsVEC3init(&m.v[2], sum_x2, sum_x3, sum_x4); 331 332 _cmsVEC3init(&v, sum_y, sum_yx, sum_yx2); 333 334 if (!_cmsMAT3solve(&res, &m, &v)) return 0; 335 336 // y = t x2 + u x + c 337 // x = ( - u + Sqrt( u^2 - 4 t c ) ) / ( 2 t ) 338 disc = res.n[1]*res.n[1] - 4.0 * res.n[0] * res.n[2]; 339 if (disc < 0) return -1; 340 341 return ( -1.0 * res.n[1] + sqrt( disc )) / (2.0 * res.n[0]); 342 } 343 344 static 345 cmsBool IsMonotonic(int n, const cmsFloat64Number Table[]) 346 { 347 int i; 348 cmsFloat64Number last; 349 350 last = Table[n-1]; 351 352 for (i = n-2; i >= 0; --i) { 353 354 if (Table[i] > last) 355 356 return FALSE; 357 else 358 last = Table[i]; 359 360 } 361 362 return TRUE; 363 } 364 365 // Calculates the black point of a destination profile. 366 // This algorithm comes from the Adobe paper disclosing its black point compensation method. 367 cmsBool CMSEXPORT cmsDetectDestinationBlackPoint(cmsCIEXYZ* BlackPoint, cmsHPROFILE hProfile, cmsUInt32Number Intent, cmsUInt32Number dwFlags) 368 { 369 cmsColorSpaceSignature ColorSpace; 370 cmsHTRANSFORM hRoundTrip = NULL; 371 cmsCIELab InitialLab, destLab, Lab; 372 373 cmsFloat64Number MinL, MaxL; 374 cmsBool NearlyStraightMidRange = FALSE; 375 cmsFloat64Number L; 376 cmsFloat64Number x[101], y[101]; 377 cmsFloat64Number lo, hi, NonMonoMin; 378 int n, l, i, NonMonoIndx; 379 380 381 // Make sure intent is adequate 382 if (Intent != INTENT_PERCEPTUAL && 383 Intent != INTENT_RELATIVE_COLORIMETRIC && 384 Intent != INTENT_SATURATION) { 385 BlackPoint -> X = BlackPoint ->Y = BlackPoint -> Z = 0.0; 386 return FALSE; 387 } 388 389 390 // v4 + perceptual & saturation intents does have its own black point, and it is 391 // well specified enough to use it. Black point tag is deprecated in V4. 392 if ((cmsGetEncodedICCversion(hProfile) >= 0x4000000) && 393 (Intent == INTENT_PERCEPTUAL || Intent == INTENT_SATURATION)) { 394 395 // Matrix shaper share MRC & perceptual intents 396 if (cmsIsMatrixShaper(hProfile)) 397 return BlackPointAsDarkerColorant(hProfile, INTENT_RELATIVE_COLORIMETRIC, BlackPoint, 0); 398 399 // Get Perceptual black out of v4 profiles. That is fixed for perceptual & saturation intents 400 BlackPoint -> X = cmsPERCEPTUAL_BLACK_X; 401 BlackPoint -> Y = cmsPERCEPTUAL_BLACK_Y; 402 BlackPoint -> Z = cmsPERCEPTUAL_BLACK_Z; 403 return TRUE; 404 } 405 406 407 // Check if the profile is lut based and gray, rgb or cmyk (7.2 in Adobe's document) 408 ColorSpace = cmsGetColorSpace(hProfile); 409 if (!cmsIsCLUT(hProfile, Intent, LCMS_USED_AS_OUTPUT ) || 410 (ColorSpace != cmsSigGrayData && 411 ColorSpace != cmsSigRgbData && 412 ColorSpace != cmsSigCmykData)) { 413 414 // In this case, handle as input case 415 return cmsDetectBlackPoint(BlackPoint, hProfile, Intent, dwFlags); 416 } 417 418 // It is one of the valid cases!, presto chargo hocus pocus, go for the Adobe magic 419 420 // Step 1 421 // ====== 422 423 // Set a first guess, that should work on good profiles. 424 if (Intent == INTENT_RELATIVE_COLORIMETRIC) { 425 426 cmsCIEXYZ IniXYZ; 427 428 // calculate initial Lab as source black point 429 if (!cmsDetectBlackPoint(&IniXYZ, hProfile, Intent, dwFlags)) { 430 return FALSE; 431 } 432 433 // convert the XYZ to lab 434 cmsXYZ2Lab(NULL, &InitialLab, &IniXYZ); 435 436 } else { 437 438 // set the initial Lab to zero, that should be the black point for perceptual and saturation 439 InitialLab.L = 0; 440 InitialLab.a = 0; 441 InitialLab.b = 0; 442 } 443 444 445 // Step 2 446 // ====== 447 448 // Create a roundtrip. Define a Transform BT for all x in L*a*b* 449 hRoundTrip = CreateRoundtripXForm(hProfile, Intent); 450 if (hRoundTrip == NULL) return FALSE; 451 452 // Calculate Min L* 453 Lab = InitialLab; 454 Lab.L = 0; 455 cmsDoTransform(hRoundTrip, &Lab, &destLab, 1); 456 MinL = destLab.L; 457 458 // Calculate Max L* 459 Lab = InitialLab; 460 Lab.L = 100; 461 cmsDoTransform(hRoundTrip, &Lab, &destLab, 1); 462 MaxL = destLab.L; 463 464 // Step 3 465 // ====== 466 467 // check if quadratic estimation needs to be done. 468 if (Intent == INTENT_RELATIVE_COLORIMETRIC) { 469 470 // Conceptually, this code tests how close the source l and converted L are to one another in the mid-range 471 // of the values. If the converted ramp of L values is close enough to a straight line y=x, then InitialLab 472 // is good enough to be the DestinationBlackPoint, 473 NearlyStraightMidRange = TRUE; 474 475 for (l=0; l <= 100; l++) { 476 477 Lab.L = l; 478 Lab.a = InitialLab.a; 479 Lab.b = InitialLab.b; 480 481 cmsDoTransform(hRoundTrip, &Lab, &destLab, 1); 482 483 L = destLab.L; 484 485 // Check the mid range in 20% after MinL 486 if (L > (MinL + 0.2 * (MaxL - MinL))) { 487 488 // Is close enough? 489 if (fabs(L - l) > 4.0) { 490 491 // Too far away, profile is buggy! 492 NearlyStraightMidRange = FALSE; 493 break; 494 } 495 } 496 } 497 } 498 else { 499 // Check is always performed for perceptual and saturation intents 500 NearlyStraightMidRange = FALSE; 501 } 502 503 504 // If no furter checking is needed, we are done 505 if (NearlyStraightMidRange) { 506 507 cmsLab2XYZ(NULL, BlackPoint, &InitialLab); 508 cmsDeleteTransform(hRoundTrip); 509 return TRUE; 510 } 511 512 // The round-trip curve normally looks like a nearly constant section at the black point, 513 // with a corner and a nearly straight line to the white point. 514 515 // STEP 4 516 // ======= 517 518 // find the black point using the least squares error quadratic curve fitting 519 520 if (Intent == INTENT_RELATIVE_COLORIMETRIC) { 521 lo = 0.1; 522 hi = 0.5; 523 } 524 else { 525 526 // Perceptual and saturation 527 lo = 0.03; 528 hi = 0.25; 529 } 530 531 // Capture points for the fitting. 532 n = 0; 533 for (l=0; l <= 100; l++) { 534 535 cmsFloat64Number ff; 536 537 Lab.L = (cmsFloat64Number) l; 538 Lab.a = InitialLab.a; 539 Lab.b = InitialLab.b; 540 541 cmsDoTransform(hRoundTrip, &Lab, &destLab, 1); 542 543 ff = (destLab.L - MinL)/(MaxL - MinL); 544 545 if (ff >= lo && ff < hi) { 546 547 x[n] = Lab.L; 548 y[n] = ff; 549 n++; 550 } 551 552 } 553 554 // This part is not on the Adobe paper, but I found is necessary for getting any result. 555 556 if (IsMonotonic(n, y)) { 557 558 // Monotonic means lower point is stil valid 559 cmsLab2XYZ(NULL, BlackPoint, &InitialLab); 560 cmsDeleteTransform(hRoundTrip); 561 return TRUE; 562 } 563 564 // No suitable points, regret and use safer algorithm 565 if (n == 0) { 566 cmsDeleteTransform(hRoundTrip); 567 return cmsDetectBlackPoint(BlackPoint, hProfile, Intent, dwFlags); 568 } 569 570 571 NonMonoMin = 100; 572 NonMonoIndx = 0; 573 for (i=0; i < n; i++) { 574 575 if (y[i] < NonMonoMin) { 576 NonMonoIndx = i; 577 NonMonoMin = y[i]; 578 } 579 } 580 581 Lab.L = x[NonMonoIndx]; 582 583 // fit and get the vertex of quadratic curve 584 Lab.L = RootOfLeastSquaresFitQuadraticCurve(n, x, y); 585 586 if (Lab.L < 0.0 || Lab.L > 50.0) { // clip to zero L* if the vertex is negative 587 Lab.L = 0; 588 } 589 590 Lab.a = InitialLab.a; 591 Lab.b = InitialLab.b; 592 593 cmsLab2XYZ(NULL, BlackPoint, &Lab); 594 595 cmsDeleteTransform(hRoundTrip); 596 return TRUE; 597 }