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-2023 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 // Allocates an empty multi profile element
  60 cmsStage* CMSEXPORT _cmsStageAllocPlaceholder(cmsContext ContextID,
  61                                 cmsStageSignature Type,
  62                                 cmsUInt32Number InputChannels,
  63                                 cmsUInt32Number OutputChannels,
  64                                 _cmsStageEvalFn     EvalPtr,
  65                                 _cmsStageDupElemFn  DupElemPtr,
  66                                 _cmsStageFreeElemFn FreePtr,
  67                                 void*             Data)
  68 {
  69     cmsStage* ph = (cmsStage*) _cmsMallocZero(ContextID, sizeof(cmsStage));
  70 
  71     if (ph == NULL) return NULL;
  72 
  73 
  74     ph ->ContextID = ContextID;
  75 
  76     ph ->Type       = Type;
  77     ph ->Implements = Type;   // By default, no clue on what is implementing
  78 
  79     ph ->InputChannels  = InputChannels;
  80     ph ->OutputChannels = OutputChannels;
  81     ph ->EvalPtr        = EvalPtr;
  82     ph ->DupElemPtr     = DupElemPtr;
  83     ph ->FreePtr        = FreePtr;
  84     ph ->Data           = Data;
  85 
  86     return ph;
  87 }
  88 
  89 
  90 static
  91 void EvaluateIdentity(const cmsFloat32Number In[],
  92                             cmsFloat32Number Out[],
  93                       const cmsStage *mpe)
  94 {
  95     memmove(Out, In, mpe ->InputChannels * sizeof(cmsFloat32Number));
  96 }
  97 
  98 
  99 cmsStage* CMSEXPORT cmsStageAllocIdentity(cmsContext ContextID, cmsUInt32Number nChannels)
 100 {
 101     return _cmsStageAllocPlaceholder(ContextID,
 102                                    cmsSigIdentityElemType,
 103                                    nChannels, nChannels,
 104                                    EvaluateIdentity,
 105                                    NULL,
 106                                    NULL,
 107                                    NULL);
 108  }
 109 
 110 // Conversion functions. From floating point to 16 bits
 111 static
 112 void FromFloatTo16(const cmsFloat32Number In[], cmsUInt16Number Out[], cmsUInt32Number n)
 113 {
 114     cmsUInt32Number i;
 115 
 116     for (i=0; i < n; i++) {
 117         Out[i] = _cmsQuickSaturateWord(In[i] * 65535.0);
 118     }
 119 }
 120 
 121 // From 16 bits to floating point
 122 static
 123 void From16ToFloat(const cmsUInt16Number In[], cmsFloat32Number Out[], cmsUInt32Number n)
 124 {
 125     cmsUInt32Number i;
 126 
 127     for (i=0; i < n; i++) {
 128         Out[i] = (cmsFloat32Number) In[i] / 65535.0F;
 129     }
 130 }
 131 
 132 
 133 // This function is quite useful to analyze the structure of a LUT and retrieve the MPE elements
 134 // that conform the LUT. It should be called with the LUT, the number of expected elements and
 135 // then a list of expected types followed with a list of cmsFloat64Number pointers to MPE elements. If
 136 // the function founds a match with current pipeline, it fills the pointers and returns TRUE
 137 // if not, returns FALSE without touching anything. Setting pointers to NULL does bypass
 138 // the storage process.
 139 cmsBool  CMSEXPORT cmsPipelineCheckAndRetreiveStages(const cmsPipeline* Lut, cmsUInt32Number n, ...)
 140 {
 141     va_list args;
 142     cmsUInt32Number i;
 143     cmsStage* mpe;
 144     cmsStageSignature Type;
 145     void** ElemPtr;
 146 
 147     // Make sure same number of elements
 148     if (cmsPipelineStageCount(Lut) != n) return FALSE;
 149 
 150     va_start(args, n);
 151 
 152     // Iterate across asked types
 153     mpe = Lut ->Elements;
 154     for (i=0; i < n; i++) {
 155 
 156         // Get asked type. cmsStageSignature is promoted to int by compiler
 157         Type  = (cmsStageSignature)va_arg(args, int);
 158         if (mpe ->Type != Type) {
 159 
 160             va_end(args);       // Mismatch. We are done.
 161             return FALSE;
 162         }
 163         mpe = mpe ->Next;
 164     }
 165 
 166     // Found a combination, fill pointers if not NULL
 167     mpe = Lut ->Elements;
 168     for (i=0; i < n; i++) {
 169 
 170         ElemPtr = va_arg(args, void**);
 171         if (ElemPtr != NULL)
 172             *ElemPtr = mpe;
 173 
 174         mpe = mpe ->Next;
 175     }
 176 
 177     va_end(args);
 178     return TRUE;
 179 }
 180 
 181 // Below there are implementations for several types of elements. Each type may be implemented by a
 182 // evaluation function, a duplication function, a function to free resources and a constructor.
 183 
 184 // *************************************************************************************************
 185 // Type cmsSigCurveSetElemType (curves)
 186 // *************************************************************************************************
 187 
 188 cmsToneCurve** _cmsStageGetPtrToCurveSet(const cmsStage* mpe)
 189 {
 190     _cmsStageToneCurvesData* Data = (_cmsStageToneCurvesData*) mpe ->Data;
 191 
 192     return Data ->TheCurves;
 193 }
 194 
 195 static
 196 void EvaluateCurves(const cmsFloat32Number In[],
 197                     cmsFloat32Number Out[],
 198                     const cmsStage *mpe)
 199 {
 200     _cmsStageToneCurvesData* Data;
 201     cmsUInt32Number i;
 202 
 203     _cmsAssert(mpe != NULL);
 204 
 205     Data = (_cmsStageToneCurvesData*) mpe ->Data;
 206     if (Data == NULL) return;
 207 
 208     if (Data ->TheCurves == NULL) return;
 209 
 210     for (i=0; i < Data ->nCurves; i++) {
 211         Out[i] = cmsEvalToneCurveFloat(Data ->TheCurves[i], In[i]);
 212     }
 213 }
 214 
 215 static
 216 void CurveSetElemTypeFree(cmsStage* mpe)
 217 {
 218     _cmsStageToneCurvesData* Data;
 219     cmsUInt32Number i;
 220 
 221     _cmsAssert(mpe != NULL);
 222 
 223     Data = (_cmsStageToneCurvesData*) mpe ->Data;
 224     if (Data == NULL) return;
 225 
 226     if (Data ->TheCurves != NULL) {
 227         for (i=0; i < Data ->nCurves; i++) {
 228             if (Data ->TheCurves[i] != NULL)
 229                 cmsFreeToneCurve(Data ->TheCurves[i]);
 230         }
 231     }
 232     _cmsFree(mpe ->ContextID, Data ->TheCurves);
 233     _cmsFree(mpe ->ContextID, Data);
 234 }
 235 
 236 
 237 static
 238 void* CurveSetDup(cmsStage* mpe)
 239 {
 240     _cmsStageToneCurvesData* Data = (_cmsStageToneCurvesData*) mpe ->Data;
 241     _cmsStageToneCurvesData* NewElem;
 242     cmsUInt32Number i;
 243 
 244     NewElem = (_cmsStageToneCurvesData*) _cmsMallocZero(mpe ->ContextID, sizeof(_cmsStageToneCurvesData));
 245     if (NewElem == NULL) return NULL;
 246 
 247     NewElem ->nCurves   = Data ->nCurves;
 248     NewElem ->TheCurves = (cmsToneCurve**) _cmsCalloc(mpe ->ContextID, NewElem ->nCurves, sizeof(cmsToneCurve*));
 249 
 250     if (NewElem ->TheCurves == NULL) goto Error;
 251 
 252     for (i=0; i < NewElem ->nCurves; i++) {
 253 
 254         // Duplicate each curve. It may fail.
 255         NewElem ->TheCurves[i] = cmsDupToneCurve(Data ->TheCurves[i]);
 256         if (NewElem ->TheCurves[i] == NULL) goto Error;
 257 
 258 
 259     }
 260     return (void*) NewElem;
 261 
 262 Error:
 263 
 264     if (NewElem ->TheCurves != NULL) {
 265         for (i=0; i < NewElem ->nCurves; i++) {
 266             if (NewElem ->TheCurves[i])
 267                 cmsFreeToneCurve(NewElem ->TheCurves[i]);
 268         }
 269     }
 270     _cmsFree(mpe ->ContextID, NewElem ->TheCurves);
 271     _cmsFree(mpe ->ContextID, NewElem);
 272     return NULL;
 273 }
 274 
 275 
 276 // Curves == NULL forces identity curves
 277 cmsStage* CMSEXPORT cmsStageAllocToneCurves(cmsContext ContextID, cmsUInt32Number nChannels, cmsToneCurve* const Curves[])
 278 {
 279     cmsUInt32Number i;
 280     _cmsStageToneCurvesData* NewElem;
 281     cmsStage* NewMPE;
 282 
 283 
 284     NewMPE = _cmsStageAllocPlaceholder(ContextID, cmsSigCurveSetElemType, nChannels, nChannels,
 285                                      EvaluateCurves, CurveSetDup, CurveSetElemTypeFree, NULL );
 286     if (NewMPE == NULL) return NULL;
 287 
 288     NewElem = (_cmsStageToneCurvesData*) _cmsMallocZero(ContextID, sizeof(_cmsStageToneCurvesData));
 289     if (NewElem == NULL) {
 290         cmsStageFree(NewMPE);
 291         return NULL;
 292     }
 293 
 294     NewMPE ->Data  = (void*) NewElem;
 295 
 296     NewElem ->nCurves   = nChannels;
 297     NewElem ->TheCurves = (cmsToneCurve**) _cmsCalloc(ContextID, nChannels, sizeof(cmsToneCurve*));
 298     if (NewElem ->TheCurves == NULL) {
 299         cmsStageFree(NewMPE);
 300         return NULL;
 301     }
 302 
 303     for (i=0; i < nChannels; i++) {
 304 
 305         if (Curves == NULL) {
 306             NewElem ->TheCurves[i] = cmsBuildGamma(ContextID, 1.0);
 307         }
 308         else {
 309             NewElem ->TheCurves[i] = cmsDupToneCurve(Curves[i]);
 310         }
 311 
 312         if (NewElem ->TheCurves[i] == NULL) {
 313             cmsStageFree(NewMPE);
 314             return NULL;
 315         }
 316 
 317     }
 318 
 319    return NewMPE;
 320 }
 321 
 322 
 323 // Create a bunch of identity curves
 324 cmsStage* CMSEXPORT _cmsStageAllocIdentityCurves(cmsContext ContextID, cmsUInt32Number nChannels)
 325 {
 326     cmsStage* mpe = cmsStageAllocToneCurves(ContextID, nChannels, NULL);
 327 
 328     if (mpe == NULL) return NULL;
 329     mpe ->Implements = cmsSigIdentityElemType;
 330     return mpe;
 331 }
 332 
 333 
 334 // *************************************************************************************************
 335 // Type cmsSigMatrixElemType (Matrices)
 336 // *************************************************************************************************
 337 
 338 
 339 // Special care should be taken here because precision loss. A temporary cmsFloat64Number buffer is being used
 340 static
 341 void EvaluateMatrix(const cmsFloat32Number In[],
 342                     cmsFloat32Number Out[],
 343                     const cmsStage *mpe)
 344 {
 345     cmsUInt32Number i, j;
 346     _cmsStageMatrixData* Data = (_cmsStageMatrixData*) mpe ->Data;
 347     cmsFloat64Number Tmp;
 348 
 349     // Input is already in 0..1.0 notation
 350     for (i=0; i < mpe ->OutputChannels; i++) {
 351 
 352         Tmp = 0;
 353         for (j=0; j < mpe->InputChannels; j++) {
 354             Tmp += In[j] * Data->Double[i*mpe->InputChannels + j];
 355         }
 356 
 357         if (Data ->Offset != NULL)
 358             Tmp += Data->Offset[i];
 359 
 360         Out[i] = (cmsFloat32Number) Tmp;
 361     }
 362 
 363 
 364     // Output in 0..1.0 domain
 365 }
 366 
 367 
 368 // Duplicate a yet-existing matrix element
 369 static
 370 void* MatrixElemDup(cmsStage* mpe)
 371 {
 372     _cmsStageMatrixData* Data = (_cmsStageMatrixData*) mpe ->Data;
 373     _cmsStageMatrixData* NewElem;
 374     cmsUInt32Number sz;
 375 
 376     NewElem = (_cmsStageMatrixData*) _cmsMallocZero(mpe ->ContextID, sizeof(_cmsStageMatrixData));
 377     if (NewElem == NULL) return NULL;
 378 
 379     sz = mpe ->InputChannels * mpe ->OutputChannels;
 380 
 381     NewElem ->Double = (cmsFloat64Number*) _cmsDupMem(mpe ->ContextID, Data ->Double, sz * sizeof(cmsFloat64Number)) ;
 382 
 383     if (Data ->Offset)
 384         NewElem ->Offset = (cmsFloat64Number*) _cmsDupMem(mpe ->ContextID,
 385                                                 Data ->Offset, mpe -> OutputChannels * sizeof(cmsFloat64Number)) ;
 386 
 387     return (void*) NewElem;
 388 }
 389 
 390 
 391 static
 392 void MatrixElemTypeFree(cmsStage* mpe)
 393 {
 394     _cmsStageMatrixData* Data = (_cmsStageMatrixData*) mpe ->Data;
 395     if (Data == NULL)
 396         return;
 397     if (Data ->Double)
 398         _cmsFree(mpe ->ContextID, Data ->Double);
 399 
 400     if (Data ->Offset)
 401         _cmsFree(mpe ->ContextID, Data ->Offset);
 402 
 403     _cmsFree(mpe ->ContextID, mpe ->Data);
 404 }
 405 
 406 
 407 
 408 cmsStage*  CMSEXPORT cmsStageAllocMatrix(cmsContext ContextID, cmsUInt32Number Rows, cmsUInt32Number Cols,
 409                                      const cmsFloat64Number* Matrix, const cmsFloat64Number* Offset)
 410 {
 411     cmsUInt32Number i, n;
 412     _cmsStageMatrixData* NewElem;
 413     cmsStage* NewMPE;
 414 
 415     n = Rows * Cols;
 416 
 417     // Check for overflow
 418     if (n == 0) return NULL;
 419     if (n >= UINT_MAX / Cols) return NULL;
 420     if (n >= UINT_MAX / Rows) return NULL;
 421     if (n < Rows || n < Cols) return NULL;
 422 
 423     NewMPE = _cmsStageAllocPlaceholder(ContextID, cmsSigMatrixElemType, Cols, Rows,
 424                                      EvaluateMatrix, MatrixElemDup, MatrixElemTypeFree, NULL );
 425     if (NewMPE == NULL) return NULL;
 426 
 427 
 428     NewElem = (_cmsStageMatrixData*) _cmsMallocZero(ContextID, sizeof(_cmsStageMatrixData));
 429     if (NewElem == NULL) goto Error;
 430     NewMPE->Data = (void*)NewElem;
 431 
 432     NewElem ->Double = (cmsFloat64Number*) _cmsCalloc(ContextID, n, sizeof(cmsFloat64Number));
 433     if (NewElem->Double == NULL) goto Error;
 434 
 435     for (i=0; i < n; i++) {
 436         NewElem ->Double[i] = Matrix[i];
 437     }
 438 
 439     if (Offset != NULL) {
 440 
 441         NewElem ->Offset = (cmsFloat64Number*) _cmsCalloc(ContextID, Rows, sizeof(cmsFloat64Number));
 442         if (NewElem->Offset == NULL) goto Error;
 443 
 444         for (i=0; i < Rows; i++) {
 445                 NewElem ->Offset[i] = Offset[i];
 446         }
 447     }
 448 
 449     return NewMPE;
 450 
 451 Error:
 452     cmsStageFree(NewMPE);
 453     return NULL;
 454 }
 455 
 456 
 457 // *************************************************************************************************
 458 // Type cmsSigCLutElemType
 459 // *************************************************************************************************
 460 
 461 
 462 // Evaluate in true floating point
 463 static
 464 void EvaluateCLUTfloat(const cmsFloat32Number In[], cmsFloat32Number Out[], const cmsStage *mpe)
 465 {
 466     _cmsStageCLutData* Data = (_cmsStageCLutData*) mpe ->Data;
 467 
 468     Data -> Params ->Interpolation.LerpFloat(In, Out, Data->Params);
 469 }
 470 
 471 
 472 // Convert to 16 bits, evaluate, and back to floating point
 473 static
 474 void EvaluateCLUTfloatIn16(const cmsFloat32Number In[], cmsFloat32Number Out[], const cmsStage *mpe)
 475 {
 476     _cmsStageCLutData* Data = (_cmsStageCLutData*) mpe ->Data;
 477     cmsUInt16Number In16[MAX_STAGE_CHANNELS], Out16[MAX_STAGE_CHANNELS];
 478 
 479     _cmsAssert(mpe ->InputChannels  <= MAX_STAGE_CHANNELS);
 480     _cmsAssert(mpe ->OutputChannels <= MAX_STAGE_CHANNELS);
 481 
 482     FromFloatTo16(In, In16, mpe ->InputChannels);
 483     Data -> Params ->Interpolation.Lerp16(In16, Out16, Data->Params);
 484     From16ToFloat(Out16, Out,  mpe ->OutputChannels);
 485 }
 486 
 487 
 488 // Given an hypercube of b dimensions, with Dims[] number of nodes by dimension, calculate the total amount of nodes
 489 static
 490 cmsUInt32Number CubeSize(const cmsUInt32Number Dims[], cmsUInt32Number b)
 491 {
 492     cmsUInt32Number rv, dim;
 493 
 494     _cmsAssert(Dims != NULL);
 495 
 496     for (rv = 1; b > 0; b--) {
 497 
 498         dim = Dims[b-1];
 499         if (dim <= 1) return 0;  // Error
 500 
 501         rv *= dim;
 502 
 503         // Check for overflow
 504         if (rv > UINT_MAX / dim) return 0;
 505     }
 506 
 507     // Again, prevent overflow
 508     if (rv > UINT_MAX / 15) return 0;
 509 
 510     return rv;
 511 }
 512 
 513 static
 514 void* CLUTElemDup(cmsStage* mpe)
 515 {
 516     _cmsStageCLutData* Data = (_cmsStageCLutData*) mpe ->Data;
 517     _cmsStageCLutData* NewElem;
 518 
 519 
 520     NewElem = (_cmsStageCLutData*) _cmsMallocZero(mpe ->ContextID, sizeof(_cmsStageCLutData));
 521     if (NewElem == NULL) return NULL;
 522 
 523     NewElem ->nEntries       = Data ->nEntries;
 524     NewElem ->HasFloatValues = Data ->HasFloatValues;
 525 
 526     if (Data ->Tab.T) {
 527 
 528         if (Data ->HasFloatValues) {
 529             NewElem ->Tab.TFloat = (cmsFloat32Number*) _cmsDupMem(mpe ->ContextID, Data ->Tab.TFloat, Data ->nEntries * sizeof (cmsFloat32Number));
 530             if (NewElem ->Tab.TFloat == NULL)
 531                 goto Error;
 532         } else {
 533             NewElem ->Tab.T = (cmsUInt16Number*) _cmsDupMem(mpe ->ContextID, Data ->Tab.T, Data ->nEntries * sizeof (cmsUInt16Number));
 534             if (NewElem ->Tab.T == NULL)
 535                 goto Error;
 536         }
 537     }
 538 
 539     NewElem ->Params   = _cmsComputeInterpParamsEx(mpe ->ContextID,
 540                                                    Data ->Params ->nSamples,
 541                                                    Data ->Params ->nInputs,
 542                                                    Data ->Params ->nOutputs,
 543                                                    NewElem ->Tab.T,
 544                                                    Data ->Params ->dwFlags);
 545     if (NewElem->Params != NULL)
 546         return (void*) NewElem;
 547  Error:
 548     if (NewElem->Tab.T)
 549         // This works for both types
 550         _cmsFree(mpe ->ContextID, NewElem -> Tab.T);
 551     _cmsFree(mpe ->ContextID, NewElem);
 552     return NULL;
 553 }
 554 
 555 
 556 static
 557 void CLutElemTypeFree(cmsStage* mpe)
 558 {
 559 
 560     _cmsStageCLutData* Data = (_cmsStageCLutData*) mpe ->Data;
 561 
 562     // Already empty
 563     if (Data == NULL) return;
 564 
 565     // This works for both types
 566     if (Data -> Tab.T)
 567         _cmsFree(mpe ->ContextID, Data -> Tab.T);
 568 
 569     _cmsFreeInterpParams(Data ->Params);
 570     _cmsFree(mpe ->ContextID, mpe ->Data);
 571 }
 572 
 573 
 574 // Allocates a 16-bit multidimensional CLUT. This is evaluated at 16-bit precision. Table may have different
 575 // granularity on each dimension.
 576 cmsStage* CMSEXPORT cmsStageAllocCLut16bitGranular(cmsContext ContextID,
 577                                          const cmsUInt32Number clutPoints[],
 578                                          cmsUInt32Number inputChan,
 579                                          cmsUInt32Number outputChan,
 580                                          const cmsUInt16Number* Table)
 581 {
 582     cmsUInt32Number i, n;
 583     _cmsStageCLutData* NewElem;
 584     cmsStage* NewMPE;
 585 
 586     _cmsAssert(clutPoints != NULL);
 587 
 588     if (inputChan > MAX_INPUT_DIMENSIONS) {
 589         cmsSignalError(ContextID, cmsERROR_RANGE, "Too many input channels (%d channels, max=%d)", inputChan, MAX_INPUT_DIMENSIONS);
 590         return NULL;
 591     }
 592 
 593     NewMPE = _cmsStageAllocPlaceholder(ContextID, cmsSigCLutElemType, inputChan, outputChan,
 594                                      EvaluateCLUTfloatIn16, CLUTElemDup, CLutElemTypeFree, NULL );
 595 
 596     if (NewMPE == NULL) return NULL;
 597 
 598     NewElem = (_cmsStageCLutData*) _cmsMallocZero(ContextID, sizeof(_cmsStageCLutData));
 599     if (NewElem == NULL) {
 600         cmsStageFree(NewMPE);
 601         return NULL;
 602     }
 603 
 604     NewMPE ->Data  = (void*) NewElem;
 605 
 606     NewElem -> nEntries = n = outputChan * CubeSize(clutPoints, inputChan);
 607     NewElem -> HasFloatValues = FALSE;
 608 
 609     if (n == 0) {
 610         cmsStageFree(NewMPE);
 611         return NULL;
 612     }
 613 
 614 
 615     NewElem ->Tab.T  = (cmsUInt16Number*) _cmsCalloc(ContextID, n, sizeof(cmsUInt16Number));
 616     if (NewElem ->Tab.T == NULL) {
 617         cmsStageFree(NewMPE);
 618         return NULL;
 619     }
 620 
 621     if (Table != NULL) {
 622         for (i=0; i < n; i++) {
 623             NewElem ->Tab.T[i] = Table[i];
 624         }
 625     }
 626 
 627     NewElem ->Params = _cmsComputeInterpParamsEx(ContextID, clutPoints, inputChan, outputChan, NewElem ->Tab.T, CMS_LERP_FLAGS_16BITS);
 628     if (NewElem ->Params == NULL) {
 629         cmsStageFree(NewMPE);
 630         return NULL;
 631     }
 632 
 633     return NewMPE;
 634 }
 635 
 636 cmsStage* CMSEXPORT cmsStageAllocCLut16bit(cmsContext ContextID,
 637                                     cmsUInt32Number nGridPoints,
 638                                     cmsUInt32Number inputChan,
 639                                     cmsUInt32Number outputChan,
 640                                     const cmsUInt16Number* Table)
 641 {
 642     cmsUInt32Number Dimensions[MAX_INPUT_DIMENSIONS];
 643     int i;
 644 
 645    // Our resulting LUT would be same gridpoints on all dimensions
 646     for (i=0; i < MAX_INPUT_DIMENSIONS; i++)
 647         Dimensions[i] = nGridPoints;
 648 
 649     return cmsStageAllocCLut16bitGranular(ContextID, Dimensions, inputChan, outputChan, Table);
 650 }
 651 
 652 
 653 cmsStage* CMSEXPORT cmsStageAllocCLutFloat(cmsContext ContextID,
 654                                        cmsUInt32Number nGridPoints,
 655                                        cmsUInt32Number inputChan,
 656                                        cmsUInt32Number outputChan,
 657                                        const cmsFloat32Number* Table)
 658 {
 659    cmsUInt32Number Dimensions[MAX_INPUT_DIMENSIONS];
 660    int i;
 661 
 662     // Our resulting LUT would be same gridpoints on all dimensions
 663     for (i=0; i < MAX_INPUT_DIMENSIONS; i++)
 664         Dimensions[i] = nGridPoints;
 665 
 666     return cmsStageAllocCLutFloatGranular(ContextID, Dimensions, inputChan, outputChan, Table);
 667 }
 668 
 669 
 670 
 671 cmsStage* CMSEXPORT cmsStageAllocCLutFloatGranular(cmsContext ContextID, const cmsUInt32Number clutPoints[], cmsUInt32Number inputChan, cmsUInt32Number outputChan, const cmsFloat32Number* Table)
 672 {
 673     cmsUInt32Number i, n;
 674     _cmsStageCLutData* NewElem;
 675     cmsStage* NewMPE;
 676 
 677     _cmsAssert(clutPoints != NULL);
 678 
 679     if (inputChan > MAX_INPUT_DIMENSIONS) {
 680         cmsSignalError(ContextID, cmsERROR_RANGE, "Too many input channels (%d channels, max=%d)", inputChan, MAX_INPUT_DIMENSIONS);
 681         return NULL;
 682     }
 683 
 684     NewMPE = _cmsStageAllocPlaceholder(ContextID, cmsSigCLutElemType, inputChan, outputChan,
 685                                              EvaluateCLUTfloat, CLUTElemDup, CLutElemTypeFree, NULL);
 686     if (NewMPE == NULL) return NULL;
 687 
 688 
 689     NewElem = (_cmsStageCLutData*) _cmsMallocZero(ContextID, sizeof(_cmsStageCLutData));
 690     if (NewElem == NULL) {
 691         cmsStageFree(NewMPE);
 692         return NULL;
 693     }
 694 
 695     NewMPE ->Data  = (void*) NewElem;
 696 
 697     // There is a potential integer overflow on conputing n and nEntries.
 698     NewElem -> nEntries = n = outputChan * CubeSize(clutPoints, inputChan);
 699     NewElem -> HasFloatValues = TRUE;
 700 
 701     if (n == 0) {
 702         cmsStageFree(NewMPE);
 703         return NULL;
 704     }
 705 
 706     NewElem ->Tab.TFloat  = (cmsFloat32Number*) _cmsCalloc(ContextID, n, sizeof(cmsFloat32Number));
 707     if (NewElem ->Tab.TFloat == NULL) {
 708         cmsStageFree(NewMPE);
 709         return NULL;
 710     }
 711 
 712     if (Table != NULL) {
 713         for (i=0; i < n; i++) {
 714             NewElem ->Tab.TFloat[i] = Table[i];
 715         }
 716     }
 717 
 718     NewElem ->Params = _cmsComputeInterpParamsEx(ContextID, clutPoints,  inputChan, outputChan, NewElem ->Tab.TFloat, CMS_LERP_FLAGS_FLOAT);
 719     if (NewElem ->Params == NULL) {
 720         cmsStageFree(NewMPE);
 721         return NULL;
 722     }
 723 
 724     return NewMPE;
 725 }
 726 
 727 
 728 static
 729 int IdentitySampler(CMSREGISTER const cmsUInt16Number In[], CMSREGISTER cmsUInt16Number Out[], CMSREGISTER void * Cargo)
 730 {
 731     int nChan = *(int*) Cargo;
 732     int i;
 733 
 734     for (i=0; i < nChan; i++)
 735         Out[i] = In[i];
 736 
 737     return 1;
 738 }
 739 
 740 // Creates an MPE that just copies input to output
 741 cmsStage* CMSEXPORT _cmsStageAllocIdentityCLut(cmsContext ContextID, cmsUInt32Number nChan)
 742 {
 743     cmsUInt32Number Dimensions[MAX_INPUT_DIMENSIONS];
 744     cmsStage* mpe ;
 745     int i;
 746 
 747     for (i=0; i < MAX_INPUT_DIMENSIONS; i++)
 748         Dimensions[i] = 2;
 749 
 750     mpe = cmsStageAllocCLut16bitGranular(ContextID, Dimensions, nChan, nChan, NULL);
 751     if (mpe == NULL) return NULL;
 752 
 753     if (!cmsStageSampleCLut16bit(mpe, IdentitySampler, &nChan, 0)) {
 754         cmsStageFree(mpe);
 755         return NULL;
 756     }
 757 
 758     mpe ->Implements = cmsSigIdentityElemType;
 759     return mpe;
 760 }
 761 
 762 
 763 
 764 // Quantize a value 0 <= i < MaxSamples to 0..0xffff
 765 cmsUInt16Number CMSEXPORT _cmsQuantizeVal(cmsFloat64Number i, cmsUInt32Number MaxSamples)
 766 {
 767     cmsFloat64Number x;
 768 
 769     x = ((cmsFloat64Number) i * 65535.) / (cmsFloat64Number) (MaxSamples - 1);
 770     return _cmsQuickSaturateWord(x);
 771 }
 772 
 773 
 774 // This routine does a sweep on whole input space, and calls its callback
 775 // function on knots. returns TRUE if all ok, FALSE otherwise.
 776 cmsBool CMSEXPORT cmsStageSampleCLut16bit(cmsStage* mpe, cmsSAMPLER16 Sampler, void * Cargo, cmsUInt32Number dwFlags)
 777 {
 778     int i, t, index, rest;
 779     cmsUInt32Number nTotalPoints;
 780     cmsUInt32Number nInputs, nOutputs;
 781     cmsUInt32Number* nSamples;
 782     cmsUInt16Number In[MAX_INPUT_DIMENSIONS+1], Out[MAX_STAGE_CHANNELS];
 783     _cmsStageCLutData* clut;
 784 
 785     if (mpe == NULL) return FALSE;
 786 
 787     clut = (_cmsStageCLutData*) mpe->Data;
 788 
 789     if (clut == NULL) return FALSE;
 790 
 791     nSamples = clut->Params ->nSamples;
 792     nInputs  = clut->Params ->nInputs;
 793     nOutputs = clut->Params ->nOutputs;
 794 
 795     if (nInputs <= 0) return FALSE;
 796     if (nOutputs <= 0) return FALSE;
 797     if (nInputs > MAX_INPUT_DIMENSIONS) return FALSE;
 798     if (nOutputs >= MAX_STAGE_CHANNELS) return FALSE;
 799 
 800     memset(In, 0, sizeof(In));
 801     memset(Out, 0, sizeof(Out));
 802 
 803     nTotalPoints = CubeSize(nSamples, nInputs);
 804     if (nTotalPoints == 0) return FALSE;
 805 
 806     index = 0;
 807     for (i = 0; i < (int) nTotalPoints; i++) {
 808 
 809         rest = i;
 810         for (t = (int)nInputs - 1; t >= 0; --t) {
 811 
 812             cmsUInt32Number  Colorant = rest % nSamples[t];
 813 
 814             rest /= nSamples[t];
 815 
 816             In[t] = _cmsQuantizeVal(Colorant, nSamples[t]);
 817         }
 818 
 819         if (clut ->Tab.T != NULL) {
 820             for (t = 0; t < (int)nOutputs; t++)
 821                 Out[t] = clut->Tab.T[index + t];
 822         }
 823 
 824         if (!Sampler(In, Out, Cargo))
 825             return FALSE;
 826 
 827         if (!(dwFlags & SAMPLER_INSPECT)) {
 828 
 829             if (clut ->Tab.T != NULL) {
 830                 for (t=0; t < (int) nOutputs; t++)
 831                     clut->Tab.T[index + t] = Out[t];
 832             }
 833         }
 834 
 835         index += nOutputs;
 836     }
 837 
 838     return TRUE;
 839 }
 840 
 841 // Same as anterior, but for floating point
 842 cmsBool CMSEXPORT cmsStageSampleCLutFloat(cmsStage* mpe, cmsSAMPLERFLOAT Sampler, void * Cargo, cmsUInt32Number dwFlags)
 843 {
 844     int i, t, index, rest;
 845     cmsUInt32Number nTotalPoints;
 846     cmsUInt32Number nInputs, nOutputs;
 847     cmsUInt32Number* nSamples;
 848     cmsFloat32Number In[MAX_INPUT_DIMENSIONS+1], Out[MAX_STAGE_CHANNELS];
 849     _cmsStageCLutData* clut;
 850 
 851     if (mpe == NULL) return FALSE;
 852 
 853     clut = (_cmsStageCLutData*)mpe->Data;
 854 
 855     if (clut == NULL) return FALSE;
 856 
 857     nSamples = clut->Params ->nSamples;
 858     nInputs  = clut->Params ->nInputs;
 859     nOutputs = clut->Params ->nOutputs;
 860 
 861     if (nInputs <= 0) return FALSE;
 862     if (nOutputs <= 0) return FALSE;
 863     if (nInputs  > MAX_INPUT_DIMENSIONS) return FALSE;
 864     if (nOutputs >= MAX_STAGE_CHANNELS) return FALSE;
 865 
 866     nTotalPoints = CubeSize(nSamples, nInputs);
 867     if (nTotalPoints == 0) return FALSE;
 868 
 869     index = 0;
 870     for (i = 0; i < (int)nTotalPoints; i++) {
 871 
 872         rest = i;
 873         for (t = (int) nInputs-1; t >=0; --t) {
 874 
 875             cmsUInt32Number  Colorant = rest % nSamples[t];
 876 
 877             rest /= nSamples[t];
 878 
 879             In[t] =  (cmsFloat32Number) (_cmsQuantizeVal(Colorant, nSamples[t]) / 65535.0);
 880         }
 881 
 882         if (clut ->Tab.TFloat != NULL) {
 883             for (t=0; t < (int) nOutputs; t++)
 884                 Out[t] = clut->Tab.TFloat[index + t];
 885         }
 886 
 887         if (!Sampler(In, Out, Cargo))
 888             return FALSE;
 889 
 890         if (!(dwFlags & SAMPLER_INSPECT)) {
 891 
 892             if (clut ->Tab.TFloat != NULL) {
 893                 for (t=0; t < (int) nOutputs; t++)
 894                     clut->Tab.TFloat[index + t] = Out[t];
 895             }
 896         }
 897 
 898         index += nOutputs;
 899     }
 900 
 901     return TRUE;
 902 }
 903 
 904 
 905 
 906 // This routine does a sweep on whole input space, and calls its callback
 907 // function on knots. returns TRUE if all ok, FALSE otherwise.
 908 cmsBool CMSEXPORT cmsSliceSpace16(cmsUInt32Number nInputs, const cmsUInt32Number clutPoints[],
 909                                          cmsSAMPLER16 Sampler, void * Cargo)
 910 {
 911     int i, t, rest;
 912     cmsUInt32Number nTotalPoints;
 913     cmsUInt16Number In[cmsMAXCHANNELS];
 914 
 915     if (nInputs >= cmsMAXCHANNELS) return FALSE;
 916 
 917     nTotalPoints = CubeSize(clutPoints, nInputs);
 918     if (nTotalPoints == 0) return FALSE;
 919 
 920     for (i = 0; i < (int) nTotalPoints; i++) {
 921 
 922         rest = i;
 923         for (t = (int) nInputs-1; t >=0; --t) {
 924 
 925             cmsUInt32Number  Colorant = rest % clutPoints[t];
 926 
 927             rest /= clutPoints[t];
 928             In[t] = _cmsQuantizeVal(Colorant, clutPoints[t]);
 929 
 930         }
 931 
 932         if (!Sampler(In, NULL, Cargo))
 933             return FALSE;
 934     }
 935 
 936     return TRUE;
 937 }
 938 
 939 cmsInt32Number CMSEXPORT cmsSliceSpaceFloat(cmsUInt32Number nInputs, const cmsUInt32Number clutPoints[],
 940                                             cmsSAMPLERFLOAT Sampler, void * Cargo)
 941 {
 942     int i, t, rest;
 943     cmsUInt32Number nTotalPoints;
 944     cmsFloat32Number In[cmsMAXCHANNELS];
 945 
 946     if (nInputs >= cmsMAXCHANNELS) return FALSE;
 947 
 948     nTotalPoints = CubeSize(clutPoints, nInputs);
 949     if (nTotalPoints == 0) return FALSE;
 950 
 951     for (i = 0; i < (int) nTotalPoints; i++) {
 952 
 953         rest = i;
 954         for (t = (int) nInputs-1; t >=0; --t) {
 955 
 956             cmsUInt32Number  Colorant = rest % clutPoints[t];
 957 
 958             rest /= clutPoints[t];
 959             In[t] =  (cmsFloat32Number) (_cmsQuantizeVal(Colorant, clutPoints[t]) / 65535.0);
 960 
 961         }
 962 
 963         if (!Sampler(In, NULL, Cargo))
 964             return FALSE;
 965     }
 966 
 967     return TRUE;
 968 }
 969 
 970 // ********************************************************************************
 971 // Type cmsSigLab2XYZElemType
 972 // ********************************************************************************
 973 
 974 
 975 static
 976 void EvaluateLab2XYZ(const cmsFloat32Number In[],
 977                      cmsFloat32Number Out[],
 978                      const cmsStage *mpe)
 979 {
 980     cmsCIELab Lab;
 981     cmsCIEXYZ XYZ;
 982     const cmsFloat64Number XYZadj = MAX_ENCODEABLE_XYZ;
 983 
 984     // V4 rules
 985     Lab.L = In[0] * 100.0;
 986     Lab.a = In[1] * 255.0 - 128.0;
 987     Lab.b = In[2] * 255.0 - 128.0;
 988 
 989     cmsLab2XYZ(NULL, &XYZ, &Lab);
 990 
 991     // From XYZ, range 0..19997 to 0..1.0, note that 1.99997 comes from 0xffff
 992     // encoded as 1.15 fixed point, so 1 + (32767.0 / 32768.0)
 993 
 994     Out[0] = (cmsFloat32Number) ((cmsFloat64Number) XYZ.X / XYZadj);
 995     Out[1] = (cmsFloat32Number) ((cmsFloat64Number) XYZ.Y / XYZadj);
 996     Out[2] = (cmsFloat32Number) ((cmsFloat64Number) XYZ.Z / XYZadj);
 997     return;
 998 
 999     cmsUNUSED_PARAMETER(mpe);
1000 }
1001 
1002 
1003 // No dup or free routines needed, as the structure has no pointers in it.
1004 cmsStage* CMSEXPORT _cmsStageAllocLab2XYZ(cmsContext ContextID)
1005 {
1006     return _cmsStageAllocPlaceholder(ContextID, cmsSigLab2XYZElemType, 3, 3, EvaluateLab2XYZ, NULL, NULL, NULL);
1007 }
1008 
1009 // ********************************************************************************
1010 
1011 // v2 L=100 is supposed to be placed on 0xFF00. There is no reasonable
1012 // number of gridpoints that would make exact match. However, a prelinearization
1013 // of 258 entries, would map 0xFF00 exactly on entry 257, and this is good to avoid scum dot.
1014 // Almost all what we need but unfortunately, the rest of entries should be scaled by
1015 // (255*257/256) and this is not exact.
1016 
1017 cmsStage* _cmsStageAllocLabV2ToV4curves(cmsContext ContextID)
1018 {
1019     cmsStage* mpe;
1020     cmsToneCurve* LabTable[3];
1021     int i, j;
1022 
1023     LabTable[0] = cmsBuildTabulatedToneCurve16(ContextID, 258, NULL);
1024     LabTable[1] = cmsBuildTabulatedToneCurve16(ContextID, 258, NULL);
1025     LabTable[2] = cmsBuildTabulatedToneCurve16(ContextID, 258, NULL);
1026 
1027     for (j=0; j < 3; j++) {
1028 
1029         if (LabTable[j] == NULL) {
1030             cmsFreeToneCurveTriple(LabTable);
1031             return NULL;
1032         }
1033 
1034         // We need to map * (0xffff / 0xff00), that's same as (257 / 256)
1035         // So we can use 258-entry tables to do the trick (i / 257) * (255 * 257) * (257 / 256);
1036         for (i=0; i < 257; i++)  {
1037 
1038             LabTable[j]->Table16[i] = (cmsUInt16Number) ((i * 0xffff + 0x80) >> 8);
1039         }
1040 
1041         LabTable[j] ->Table16[257] = 0xffff;
1042     }
1043 
1044     mpe = cmsStageAllocToneCurves(ContextID, 3, LabTable);
1045     cmsFreeToneCurveTriple(LabTable);
1046 
1047     if (mpe == NULL) return NULL;
1048     mpe ->Implements = cmsSigLabV2toV4;
1049     return mpe;
1050 }
1051 
1052 // ********************************************************************************
1053 
1054 // Matrix-based conversion, which is more accurate, but slower and cannot properly be saved in devicelink profiles
1055 cmsStage* CMSEXPORT _cmsStageAllocLabV2ToV4(cmsContext ContextID)
1056 {
1057     static const cmsFloat64Number V2ToV4[] = { 65535.0/65280.0, 0, 0,
1058                                      0, 65535.0/65280.0, 0,
1059                                      0, 0, 65535.0/65280.0
1060                                      };
1061 
1062     cmsStage *mpe = cmsStageAllocMatrix(ContextID, 3, 3, V2ToV4, NULL);
1063 
1064     if (mpe == NULL) return mpe;
1065     mpe ->Implements = cmsSigLabV2toV4;
1066     return mpe;
1067 }
1068 
1069 
1070 // Reverse direction
1071 cmsStage* CMSEXPORT _cmsStageAllocLabV4ToV2(cmsContext ContextID)
1072 {
1073     static const cmsFloat64Number V4ToV2[] = { 65280.0/65535.0, 0, 0,
1074                                      0, 65280.0/65535.0, 0,
1075                                      0, 0, 65280.0/65535.0
1076                                      };
1077 
1078      cmsStage *mpe = cmsStageAllocMatrix(ContextID, 3, 3, V4ToV2, NULL);
1079 
1080     if (mpe == NULL) return mpe;
1081     mpe ->Implements = cmsSigLabV4toV2;
1082     return mpe;
1083 }
1084 
1085 
1086 // To Lab to float. Note that the MPE gives numbers in normal Lab range
1087 // and we need 0..1.0 range for the formatters
1088 // L* : 0...100 => 0...1.0  (L* / 100)
1089 // ab* : -128..+127 to 0..1  ((ab* + 128) / 255)
1090 
1091 cmsStage* _cmsStageNormalizeFromLabFloat(cmsContext ContextID)
1092 {
1093     static const cmsFloat64Number a1[] = {
1094         1.0/100.0, 0, 0,
1095         0, 1.0/255.0, 0,
1096         0, 0, 1.0/255.0
1097     };
1098 
1099     static const cmsFloat64Number o1[] = {
1100         0,
1101         128.0/255.0,
1102         128.0/255.0
1103     };
1104 
1105     cmsStage *mpe = cmsStageAllocMatrix(ContextID, 3, 3, a1, o1);
1106 
1107     if (mpe == NULL) return mpe;
1108     mpe ->Implements = cmsSigLab2FloatPCS;
1109     return mpe;
1110 }
1111 
1112 // Fom XYZ to floating point PCS
1113 cmsStage* _cmsStageNormalizeFromXyzFloat(cmsContext ContextID)
1114 {
1115 #define n (32768.0/65535.0)
1116     static const cmsFloat64Number a1[] = {
1117         n, 0, 0,
1118         0, n, 0,
1119         0, 0, n
1120     };
1121 #undef n
1122 
1123     cmsStage *mpe =  cmsStageAllocMatrix(ContextID, 3, 3, a1, NULL);
1124 
1125     if (mpe == NULL) return mpe;
1126     mpe ->Implements = cmsSigXYZ2FloatPCS;
1127     return mpe;
1128 }
1129 
1130 cmsStage* _cmsStageNormalizeToLabFloat(cmsContext ContextID)
1131 {
1132     static const cmsFloat64Number a1[] = {
1133         100.0, 0, 0,
1134         0, 255.0, 0,
1135         0, 0, 255.0
1136     };
1137 
1138     static const cmsFloat64Number o1[] = {
1139         0,
1140         -128.0,
1141         -128.0
1142     };
1143 
1144     cmsStage *mpe =  cmsStageAllocMatrix(ContextID, 3, 3, a1, o1);
1145     if (mpe == NULL) return mpe;
1146     mpe ->Implements = cmsSigFloatPCS2Lab;
1147     return mpe;
1148 }
1149 
1150 cmsStage* _cmsStageNormalizeToXyzFloat(cmsContext ContextID)
1151 {
1152 #define n (65535.0/32768.0)
1153 
1154     static const cmsFloat64Number a1[] = {
1155         n, 0, 0,
1156         0, n, 0,
1157         0, 0, n
1158     };
1159 #undef n
1160 
1161     cmsStage *mpe = cmsStageAllocMatrix(ContextID, 3, 3, a1, NULL);
1162     if (mpe == NULL) return mpe;
1163     mpe ->Implements = cmsSigFloatPCS2XYZ;
1164     return mpe;
1165 }
1166 
1167 // Clips values smaller than zero
1168 static
1169 void Clipper(const cmsFloat32Number In[], cmsFloat32Number Out[], const cmsStage *mpe)
1170 {
1171        cmsUInt32Number i;
1172        for (i = 0; i < mpe->InputChannels; i++) {
1173 
1174               cmsFloat32Number n = In[i];
1175               Out[i] = n < 0 ? 0 : n;
1176        }
1177 }
1178 
1179 cmsStage*  _cmsStageClipNegatives(cmsContext ContextID, cmsUInt32Number nChannels)
1180 {
1181        return _cmsStageAllocPlaceholder(ContextID, cmsSigClipNegativesElemType,
1182               nChannels, nChannels, Clipper, NULL, NULL, NULL);
1183 }
1184 
1185 // ********************************************************************************
1186 // Type cmsSigXYZ2LabElemType
1187 // ********************************************************************************
1188 
1189 static
1190 void EvaluateXYZ2Lab(const cmsFloat32Number In[], cmsFloat32Number Out[], const cmsStage *mpe)
1191 {
1192     cmsCIELab Lab;
1193     cmsCIEXYZ XYZ;
1194     const cmsFloat64Number XYZadj = MAX_ENCODEABLE_XYZ;
1195 
1196     // From 0..1.0 to XYZ
1197 
1198     XYZ.X = In[0] * XYZadj;
1199     XYZ.Y = In[1] * XYZadj;
1200     XYZ.Z = In[2] * XYZadj;
1201 
1202     cmsXYZ2Lab(NULL, &Lab, &XYZ);
1203 
1204     // From V4 Lab to 0..1.0
1205 
1206     Out[0] = (cmsFloat32Number) (Lab.L / 100.0);
1207     Out[1] = (cmsFloat32Number) ((Lab.a + 128.0) / 255.0);
1208     Out[2] = (cmsFloat32Number) ((Lab.b + 128.0) / 255.0);
1209     return;
1210 
1211     cmsUNUSED_PARAMETER(mpe);
1212 }
1213 
1214 cmsStage* CMSEXPORT _cmsStageAllocXYZ2Lab(cmsContext ContextID)
1215 {
1216     return _cmsStageAllocPlaceholder(ContextID, cmsSigXYZ2LabElemType, 3, 3, EvaluateXYZ2Lab, NULL, NULL, NULL);
1217 
1218 }
1219 
1220 // ********************************************************************************
1221 
1222 // For v4, S-Shaped curves are placed in a/b axis to increase resolution near gray
1223 
1224 cmsStage* _cmsStageAllocLabPrelin(cmsContext ContextID)
1225 {
1226     cmsToneCurve* LabTable[3];
1227     cmsFloat64Number Params[1] =  {2.4} ;
1228 
1229     LabTable[0] = cmsBuildGamma(ContextID, 1.0);
1230     LabTable[1] = cmsBuildParametricToneCurve(ContextID, 108, Params);
1231     LabTable[2] = cmsBuildParametricToneCurve(ContextID, 108, Params);
1232 
1233     return cmsStageAllocToneCurves(ContextID, 3, LabTable);
1234 }
1235 
1236 
1237 // Free a single MPE
1238 void CMSEXPORT cmsStageFree(cmsStage* mpe)
1239 {
1240     if (mpe ->FreePtr)
1241         mpe ->FreePtr(mpe);
1242 
1243     _cmsFree(mpe ->ContextID, mpe);
1244 }
1245 
1246 
1247 cmsUInt32Number  CMSEXPORT cmsStageInputChannels(const cmsStage* mpe)
1248 {
1249     return mpe ->InputChannels;
1250 }
1251 
1252 cmsUInt32Number  CMSEXPORT cmsStageOutputChannels(const cmsStage* mpe)
1253 {
1254     return mpe ->OutputChannels;
1255 }
1256 
1257 cmsStageSignature CMSEXPORT cmsStageType(const cmsStage* mpe)
1258 {
1259     return mpe -> Type;
1260 }
1261 
1262 void* CMSEXPORT cmsStageData(const cmsStage* mpe)
1263 {
1264     return mpe -> Data;
1265 }
1266 
1267 cmsContext CMSEXPORT cmsGetStageContextID(const cmsStage* mpe)
1268 {
1269     return mpe -> ContextID;
1270 }
1271 
1272 cmsStage*  CMSEXPORT cmsStageNext(const cmsStage* mpe)
1273 {
1274     return mpe -> Next;
1275 }
1276 
1277 
1278 // Duplicates an MPE
1279 cmsStage* CMSEXPORT cmsStageDup(cmsStage* mpe)
1280 {
1281     cmsStage* NewMPE;
1282 
1283     if (mpe == NULL) return NULL;
1284     NewMPE = _cmsStageAllocPlaceholder(mpe ->ContextID,
1285                                      mpe ->Type,
1286                                      mpe ->InputChannels,
1287                                      mpe ->OutputChannels,
1288                                      mpe ->EvalPtr,
1289                                      mpe ->DupElemPtr,
1290                                      mpe ->FreePtr,
1291                                      NULL);
1292     if (NewMPE == NULL) return NULL;
1293 
1294     NewMPE ->Implements = mpe ->Implements;
1295 
1296     if (mpe ->DupElemPtr) {
1297 
1298         NewMPE ->Data = mpe ->DupElemPtr(mpe);
1299 
1300         if (NewMPE->Data == NULL) {
1301 
1302             cmsStageFree(NewMPE);
1303             return NULL;
1304         }
1305 
1306     } else {
1307 
1308         NewMPE ->Data       = NULL;
1309     }
1310 
1311     return NewMPE;
1312 }
1313 
1314 
1315 // ***********************************************************************************************************
1316 
1317 // This function sets up the channel count
1318 static
1319 cmsBool BlessLUT(cmsPipeline* lut)
1320 {
1321     // We can set the input/output channels only if we have elements.
1322     if (lut ->Elements != NULL) {
1323 
1324         cmsStage* prev;
1325         cmsStage* next;
1326         cmsStage* First;
1327         cmsStage* Last;
1328 
1329         First  = cmsPipelineGetPtrToFirstStage(lut);
1330         Last   = cmsPipelineGetPtrToLastStage(lut);
1331 
1332         if (First == NULL || Last == NULL) return FALSE;
1333 
1334         lut->InputChannels = First->InputChannels;
1335         lut->OutputChannels = Last->OutputChannels;
1336 
1337         // Check chain consistency
1338         prev = First;
1339         next = prev->Next;
1340 
1341         while (next != NULL)
1342         {
1343             if (next->InputChannels != prev->OutputChannels)
1344                 return FALSE;
1345 
1346             next = next->Next;
1347             prev = prev->Next;
1348     }
1349 }
1350 
1351     return TRUE;
1352 }
1353 
1354 
1355 // Default to evaluate the LUT on 16 bit-basis. Precision is retained.
1356 static
1357 void _LUTeval16(CMSREGISTER const cmsUInt16Number In[], CMSREGISTER cmsUInt16Number Out[],  CMSREGISTER const void* D)
1358 {
1359     cmsPipeline* lut = (cmsPipeline*) D;
1360     cmsStage *mpe;
1361     cmsFloat32Number Storage[2][MAX_STAGE_CHANNELS];
1362     int Phase = 0, NextPhase;
1363 
1364     From16ToFloat(In, &Storage[Phase][0], lut ->InputChannels);
1365 
1366     for (mpe = lut ->Elements;
1367          mpe != NULL;
1368          mpe = mpe ->Next) {
1369 
1370              NextPhase = Phase ^ 1;
1371              mpe ->EvalPtr(&Storage[Phase][0], &Storage[NextPhase][0], mpe);
1372              Phase = NextPhase;
1373     }
1374 
1375 
1376     FromFloatTo16(&Storage[Phase][0], Out, lut ->OutputChannels);
1377 }
1378 
1379 
1380 
1381 // Does evaluate the LUT on cmsFloat32Number-basis.
1382 static
1383 void _LUTevalFloat(const cmsFloat32Number In[], cmsFloat32Number Out[], const void* D)
1384 {
1385     cmsPipeline* lut = (cmsPipeline*) D;
1386     cmsStage *mpe;
1387     cmsFloat32Number Storage[2][MAX_STAGE_CHANNELS];
1388     int Phase = 0, NextPhase;
1389 
1390     memmove(&Storage[Phase][0], In, lut ->InputChannels  * sizeof(cmsFloat32Number));
1391 
1392     for (mpe = lut ->Elements;
1393          mpe != NULL;
1394          mpe = mpe ->Next) {
1395 
1396               NextPhase = Phase ^ 1;
1397               mpe ->EvalPtr(&Storage[Phase][0], &Storage[NextPhase][0], mpe);
1398               Phase = NextPhase;
1399     }
1400 
1401     memmove(Out, &Storage[Phase][0], lut ->OutputChannels * sizeof(cmsFloat32Number));
1402 }
1403 
1404 
1405 // LUT Creation & Destruction
1406 cmsPipeline* CMSEXPORT cmsPipelineAlloc(cmsContext ContextID, cmsUInt32Number InputChannels, cmsUInt32Number OutputChannels)
1407 {
1408        cmsPipeline* NewLUT;
1409 
1410        // A value of zero in channels is allowed as placeholder
1411        if (InputChannels >= cmsMAXCHANNELS ||
1412            OutputChannels >= cmsMAXCHANNELS) return NULL;
1413 
1414        NewLUT = (cmsPipeline*) _cmsMallocZero(ContextID, sizeof(cmsPipeline));
1415        if (NewLUT == NULL) return NULL;
1416 
1417        NewLUT -> InputChannels  = InputChannels;
1418        NewLUT -> OutputChannels = OutputChannels;
1419 
1420        NewLUT ->Eval16Fn    = _LUTeval16;
1421        NewLUT ->EvalFloatFn = _LUTevalFloat;
1422        NewLUT ->DupDataFn   = NULL;
1423        NewLUT ->FreeDataFn  = NULL;
1424        NewLUT ->Data        = NewLUT;
1425        NewLUT ->ContextID   = ContextID;
1426 
1427        if (!BlessLUT(NewLUT))
1428        {
1429            _cmsFree(ContextID, NewLUT);
1430            return NULL;
1431        }
1432 
1433        return NewLUT;
1434 }
1435 
1436 cmsContext CMSEXPORT cmsGetPipelineContextID(const cmsPipeline* lut)
1437 {
1438     _cmsAssert(lut != NULL);
1439     return lut ->ContextID;
1440 }
1441 
1442 cmsUInt32Number CMSEXPORT cmsPipelineInputChannels(const cmsPipeline* lut)
1443 {
1444     _cmsAssert(lut != NULL);
1445     return lut ->InputChannels;
1446 }
1447 
1448 cmsUInt32Number CMSEXPORT cmsPipelineOutputChannels(const cmsPipeline* lut)
1449 {
1450     _cmsAssert(lut != NULL);
1451     return lut ->OutputChannels;
1452 }
1453 
1454 // Free a profile elements LUT
1455 void CMSEXPORT cmsPipelineFree(cmsPipeline* lut)
1456 {
1457     cmsStage *mpe, *Next;
1458 
1459     if (lut == NULL) return;
1460 
1461     for (mpe = lut ->Elements;
1462         mpe != NULL;
1463         mpe = Next) {
1464 
1465             Next = mpe ->Next;
1466             cmsStageFree(mpe);
1467     }
1468 
1469     if (lut ->FreeDataFn) lut ->FreeDataFn(lut ->ContextID, lut ->Data);
1470 
1471     _cmsFree(lut ->ContextID, lut);
1472 }
1473 
1474 
1475 // Default to evaluate the LUT on 16 bit-basis.
1476 void CMSEXPORT cmsPipelineEval16(const cmsUInt16Number In[], cmsUInt16Number Out[],  const cmsPipeline* lut)
1477 {
1478     _cmsAssert(lut != NULL);
1479     lut ->Eval16Fn(In, Out, lut->Data);
1480 }
1481 
1482 
1483 // Does evaluate the LUT on cmsFloat32Number-basis.
1484 void CMSEXPORT cmsPipelineEvalFloat(const cmsFloat32Number In[], cmsFloat32Number Out[], const cmsPipeline* lut)
1485 {
1486     _cmsAssert(lut != NULL);
1487     lut ->EvalFloatFn(In, Out, lut);
1488 }
1489 
1490 
1491 
1492 // Duplicates a LUT
1493 cmsPipeline* CMSEXPORT cmsPipelineDup(const cmsPipeline* lut)
1494 {
1495     cmsPipeline* NewLUT;
1496     cmsStage *NewMPE, *Anterior = NULL, *mpe;
1497     cmsBool  First = TRUE;
1498 
1499     if (lut == NULL) return NULL;
1500 
1501     NewLUT = cmsPipelineAlloc(lut ->ContextID, lut ->InputChannels, lut ->OutputChannels);
1502     if (NewLUT == NULL) return NULL;
1503 
1504     for (mpe = lut ->Elements;
1505          mpe != NULL;
1506          mpe = mpe ->Next) {
1507 
1508              NewMPE = cmsStageDup(mpe);
1509 
1510              if (NewMPE == NULL) {
1511                  cmsPipelineFree(NewLUT);
1512                  return NULL;
1513              }
1514 
1515              if (First) {
1516                  NewLUT ->Elements = NewMPE;
1517                  First = FALSE;
1518              }
1519              else {
1520                 if (Anterior != NULL)
1521                     Anterior ->Next = NewMPE;
1522              }
1523 
1524             Anterior = NewMPE;
1525     }
1526 
1527     NewLUT ->Eval16Fn    = lut ->Eval16Fn;
1528     NewLUT ->EvalFloatFn = lut ->EvalFloatFn;
1529     NewLUT ->DupDataFn   = lut ->DupDataFn;
1530     NewLUT ->FreeDataFn  = lut ->FreeDataFn;
1531 
1532     if (NewLUT ->DupDataFn != NULL)
1533         NewLUT ->Data = NewLUT ->DupDataFn(lut ->ContextID, lut->Data);
1534 
1535 
1536     NewLUT ->SaveAs8Bits    = lut ->SaveAs8Bits;
1537 
1538     if (!BlessLUT(NewLUT))
1539     {
1540         _cmsFree(lut->ContextID, NewLUT);
1541         return NULL;
1542     }
1543 
1544     return NewLUT;
1545 }
1546 
1547 
1548 int CMSEXPORT cmsPipelineInsertStage(cmsPipeline* lut, cmsStageLoc loc, cmsStage* mpe)
1549 {
1550     cmsStage* Anterior = NULL, *pt;
1551 
1552     if (lut == NULL || mpe == NULL)
1553         return FALSE;
1554 
1555     switch (loc) {
1556 
1557         case cmsAT_BEGIN:
1558             mpe ->Next = lut ->Elements;
1559             lut ->Elements = mpe;
1560             break;
1561 
1562         case cmsAT_END:
1563 
1564             if (lut ->Elements == NULL)
1565                 lut ->Elements = mpe;
1566             else {
1567 
1568                 for (pt = lut ->Elements;
1569                      pt != NULL;
1570                      pt = pt -> Next) Anterior = pt;
1571 
1572                 Anterior ->Next = mpe;
1573                 mpe ->Next = NULL;
1574             }
1575             break;
1576         default:;
1577             return FALSE;
1578     }
1579 
1580     return BlessLUT(lut);
1581 }
1582 
1583 // Unlink an element and return the pointer to it
1584 void CMSEXPORT cmsPipelineUnlinkStage(cmsPipeline* lut, cmsStageLoc loc, cmsStage** mpe)
1585 {
1586     cmsStage *Anterior, *pt, *Last;
1587     cmsStage *Unlinked = NULL;
1588 
1589 
1590     // If empty LUT, there is nothing to remove
1591     if (lut ->Elements == NULL) {
1592         if (mpe) *mpe = NULL;
1593         return;
1594     }
1595 
1596     // On depending on the strategy...
1597     switch (loc) {
1598 
1599         case cmsAT_BEGIN:
1600             {
1601                 cmsStage* elem = lut ->Elements;
1602 
1603                 lut ->Elements = elem -> Next;
1604                 elem ->Next = NULL;
1605                 Unlinked = elem;
1606 
1607             }
1608             break;
1609 
1610         case cmsAT_END:
1611             Anterior = Last = NULL;
1612             for (pt = lut ->Elements;
1613                 pt != NULL;
1614                 pt = pt -> Next) {
1615                     Anterior = Last;
1616                     Last = pt;
1617             }
1618 
1619             Unlinked = Last;  // Next already points to NULL
1620 
1621             // Truncate the chain
1622             if (Anterior)
1623                 Anterior ->Next = NULL;
1624             else
1625                 lut ->Elements = NULL;
1626             break;
1627         default:;
1628     }
1629 
1630     if (mpe)
1631         *mpe = Unlinked;
1632     else
1633         cmsStageFree(Unlinked);
1634 
1635     // May fail, but we ignore it
1636     BlessLUT(lut);
1637 }
1638 
1639 
1640 // Concatenate two LUT into a new single one
1641 cmsBool  CMSEXPORT cmsPipelineCat(cmsPipeline* l1, const cmsPipeline* l2)
1642 {
1643     cmsStage* mpe;
1644 
1645     // If both LUTS does not have elements, we need to inherit
1646     // the number of channels
1647     if (l1 ->Elements == NULL && l2 ->Elements == NULL) {
1648         l1 ->InputChannels  = l2 ->InputChannels;
1649         l1 ->OutputChannels = l2 ->OutputChannels;
1650     }
1651 
1652     // Cat second
1653     for (mpe = l2 ->Elements;
1654          mpe != NULL;
1655          mpe = mpe ->Next) {
1656 
1657             // We have to dup each element
1658             if (!cmsPipelineInsertStage(l1, cmsAT_END, cmsStageDup(mpe)))
1659                 return FALSE;
1660     }
1661 
1662     return BlessLUT(l1);
1663 }
1664 
1665 
1666 cmsBool CMSEXPORT cmsPipelineSetSaveAs8bitsFlag(cmsPipeline* lut, cmsBool On)
1667 {
1668     cmsBool Anterior = lut ->SaveAs8Bits;
1669 
1670     lut ->SaveAs8Bits = On;
1671     return Anterior;
1672 }
1673 
1674 
1675 cmsStage* CMSEXPORT cmsPipelineGetPtrToFirstStage(const cmsPipeline* lut)
1676 {
1677     return lut ->Elements;
1678 }
1679 
1680 cmsStage* CMSEXPORT cmsPipelineGetPtrToLastStage(const cmsPipeline* lut)
1681 {
1682     cmsStage *mpe, *Anterior = NULL;
1683 
1684     for (mpe = lut ->Elements; mpe != NULL; mpe = mpe ->Next)
1685         Anterior = mpe;
1686 
1687     return Anterior;
1688 }
1689 
1690 cmsUInt32Number CMSEXPORT cmsPipelineStageCount(const cmsPipeline* lut)
1691 {
1692     cmsStage *mpe;
1693     cmsUInt32Number n;
1694 
1695     for (n=0, mpe = lut ->Elements; mpe != NULL; mpe = mpe ->Next)
1696             n++;
1697 
1698     return n;
1699 }
1700 
1701 // This function may be used to set the optional evaluator and a block of private data. If private data is being used, an optional
1702 // duplicator and free functions should also be specified in order to duplicate the LUT construct. Use NULL to inhibit such functionality.
1703 void CMSEXPORT _cmsPipelineSetOptimizationParameters(cmsPipeline* Lut,
1704                                         _cmsPipelineEval16Fn Eval16,
1705                                         void* PrivateData,
1706                                         _cmsFreeUserDataFn FreePrivateDataFn,
1707                                         _cmsDupUserDataFn  DupPrivateDataFn)
1708 {
1709 
1710     Lut ->Eval16Fn = Eval16;
1711     Lut ->DupDataFn = DupPrivateDataFn;
1712     Lut ->FreeDataFn = FreePrivateDataFn;
1713     Lut ->Data = PrivateData;
1714 }
1715 
1716 
1717 // ----------------------------------------------------------- Reverse interpolation
1718 // Here's how it goes. The derivative Df(x) of the function f is the linear
1719 // transformation that best approximates f near the point x. It can be represented
1720 // by a matrix A whose entries are the partial derivatives of the components of f
1721 // with respect to all the coordinates. This is know as the Jacobian
1722 //
1723 // The best linear approximation to f is given by the matrix equation:
1724 //
1725 // y-y0 = A (x-x0)
1726 //
1727 // So, if x0 is a good "guess" for the zero of f, then solving for the zero of this
1728 // linear approximation will give a "better guess" for the zero of f. Thus let y=0,
1729 // and since y0=f(x0) one can solve the above equation for x. This leads to the
1730 // Newton's method formula:
1731 //
1732 // xn+1 = xn - A-1 f(xn)
1733 //
1734 // where xn+1 denotes the (n+1)-st guess, obtained from the n-th guess xn in the
1735 // fashion described above. Iterating this will give better and better approximations
1736 // if you have a "good enough" initial guess.
1737 
1738 
1739 #define JACOBIAN_EPSILON            0.001f
1740 #define INVERSION_MAX_ITERATIONS    30
1741 
1742 // Increment with reflexion on boundary
1743 static
1744 void IncDelta(cmsFloat32Number *Val)
1745 {
1746     if (*Val < (1.0 - JACOBIAN_EPSILON))
1747 
1748         *Val += JACOBIAN_EPSILON;
1749 
1750     else
1751         *Val -= JACOBIAN_EPSILON;
1752 
1753 }
1754 
1755 
1756 
1757 // Euclidean distance between two vectors of n elements each one
1758 static
1759 cmsFloat32Number EuclideanDistance(cmsFloat32Number a[], cmsFloat32Number b[], int n)
1760 {
1761     cmsFloat32Number sum = 0;
1762     int i;
1763 
1764     for (i=0; i < n; i++) {
1765         cmsFloat32Number dif = b[i] - a[i];
1766         sum +=  dif * dif;
1767     }
1768 
1769     return sqrtf(sum);
1770 }
1771 
1772 
1773 // Evaluate a LUT in reverse direction. It only searches on 3->3 LUT. Uses Newton method
1774 //
1775 // x1 <- x - [J(x)]^-1 * f(x)
1776 //
1777 // lut: The LUT on where to do the search
1778 // Target: LabK, 3 values of Lab plus destination K which is fixed
1779 // Result: The obtained CMYK
1780 // Hint:   Location where begin the search
1781 
1782 cmsBool CMSEXPORT cmsPipelineEvalReverseFloat(cmsFloat32Number Target[],
1783                                               cmsFloat32Number Result[],
1784                                               cmsFloat32Number Hint[],
1785                                               const cmsPipeline* lut)
1786 {
1787     cmsUInt32Number  i, j;
1788     cmsFloat64Number  error, LastError = 1E20;
1789     cmsFloat32Number  fx[4], x[4], xd[4], fxd[4];
1790     cmsVEC3 tmp, tmp2;
1791     cmsMAT3 Jacobian;
1792 
1793     // Only 3->3 and 4->3 are supported
1794     if (lut ->InputChannels != 3 && lut ->InputChannels != 4) return FALSE;
1795     if (lut ->OutputChannels != 3) return FALSE;
1796 
1797     // Take the hint as starting point if specified
1798     if (Hint == NULL) {
1799 
1800         // Begin at any point, we choose 1/3 of CMY axis
1801         x[0] = x[1] = x[2] = 0.3f;
1802     }
1803     else {
1804 
1805         // Only copy 3 channels from hint...
1806         for (j=0; j < 3; j++)
1807             x[j] = Hint[j];
1808     }
1809 
1810     // If Lut is 4-dimensions, then grab target[3], which is fixed
1811     if (lut ->InputChannels == 4) {
1812         x[3] = Target[3];
1813     }
1814     else x[3] = 0; // To keep lint happy
1815 
1816 
1817     // Iterate
1818     for (i = 0; i < INVERSION_MAX_ITERATIONS; i++) {
1819 
1820         // Get beginning fx
1821         cmsPipelineEvalFloat(x, fx, lut);
1822 
1823         // Compute error
1824         error = EuclideanDistance(fx, Target, 3);
1825 
1826         // If not convergent, return last safe value
1827         if (error >= LastError)
1828             break;
1829 
1830         // Keep latest values
1831         LastError     = error;
1832         for (j=0; j < lut ->InputChannels; j++)
1833                 Result[j] = x[j];
1834 
1835         // Found an exact match?
1836         if (error <= 0)
1837             break;
1838 
1839         // Obtain slope (the Jacobian)
1840         for (j = 0; j < 3; j++) {
1841 
1842             xd[0] = x[0];
1843             xd[1] = x[1];
1844             xd[2] = x[2];
1845             xd[3] = x[3];  // Keep fixed channel
1846 
1847             IncDelta(&xd[j]);
1848 
1849             cmsPipelineEvalFloat(xd, fxd, lut);
1850 
1851             Jacobian.v[0].n[j] = ((fxd[0] - fx[0]) / JACOBIAN_EPSILON);
1852             Jacobian.v[1].n[j] = ((fxd[1] - fx[1]) / JACOBIAN_EPSILON);
1853             Jacobian.v[2].n[j] = ((fxd[2] - fx[2]) / JACOBIAN_EPSILON);
1854         }
1855 
1856         // Solve system
1857         tmp2.n[0] = fx[0] - Target[0];
1858         tmp2.n[1] = fx[1] - Target[1];
1859         tmp2.n[2] = fx[2] - Target[2];
1860 
1861         if (!_cmsMAT3solve(&tmp, &Jacobian, &tmp2))
1862             return FALSE;
1863 
1864         // Move our guess
1865         x[0] -= (cmsFloat32Number) tmp.n[0];
1866         x[1] -= (cmsFloat32Number) tmp.n[1];
1867         x[2] -= (cmsFloat32Number) tmp.n[2];
1868 
1869         // Some clipping....
1870         for (j=0; j < 3; j++) {
1871             if (x[j] < 0) x[j] = 0;
1872             else
1873                 if (x[j] > 1.0) x[j] = 1.0;
1874         }
1875     }
1876 
1877     return TRUE;
1878 }
1879 
1880