summaryrefslogtreecommitdiff
path: root/xicc/xicc.c
diff options
context:
space:
mode:
authorJörg Frings-Fürst <debian@jff-webhosting.net>2014-09-01 13:56:46 +0200
committerJörg Frings-Fürst <debian@jff-webhosting.net>2014-09-01 13:56:46 +0200
commit22f703cab05b7cd368f4de9e03991b7664dc5022 (patch)
tree6f4d50beaa42328e24b1c6b56b6ec059e4ef21a5 /xicc/xicc.c
Initial import of argyll version 1.5.1-8debian/1.5.1-8
Diffstat (limited to 'xicc/xicc.c')
-rw-r--r--xicc/xicc.c3607
1 files changed, 3607 insertions, 0 deletions
diff --git a/xicc/xicc.c b/xicc/xicc.c
new file mode 100644
index 0000000..9b6d867
--- /dev/null
+++ b/xicc/xicc.c
@@ -0,0 +1,3607 @@
+
+/*
+ * International Color Consortium color transform expanded support
+ *
+ * Author: Graeme W. Gill
+ * Date: 2/7/00
+ * Version: 1.00
+ *
+ * Copyright 2000, 2001 Graeme W. Gill
+ * All rights reserved.
+ * This material is licenced under the GNU AFFERO GENERAL PUBLIC LICENSE Version 3 :-
+ * see the License.txt file for licencing details.
+ *
+ * Based on the old iccXfm class.
+ */
+
+/*
+ * This module expands the basic icclib functionality,
+ * providing more functionality in exercising.
+ * The implementation for the three different types
+ * of profile representation, are in their own source files.
+ */
+
+/*
+ * TTBD:
+ * Some of the error handling is crude. Shouldn't use
+ * error(), should return status.
+ *
+ */
+
+#include <sys/types.h>
+#include <string.h>
+#include <ctype.h>
+#ifdef __sun
+#include <unistd.h>
+#endif
+#if defined(__IBMC__) && defined(_M_IX86)
+#include <float.h>
+#endif
+#include "numlib.h"
+#include "counters.h"
+#include "plot.h"
+#include "../h/sort.h"
+#include "xicc.h" /* definitions for this library */
+
+#define USE_CAM /* Use CIECAM02 for clipping and gamut mapping, else use Lab */
+
+#undef DEBUG /* Plot 1d Luts */
+
+#ifdef DEBUG
+#include "plot.h"
+#endif
+
+#define MAX_INVSOLN 4
+
+static void xicc_del(xicc *p);
+icxLuBase * xicc_get_luobj(xicc *p, int flags, icmLookupFunc func, icRenderingIntent intent,
+ icColorSpaceSignature pcsor, icmLookupOrder order,
+ icxViewCond *vc, icxInk *ink);
+static icxLuBase *xicc_set_luobj(xicc *p, icmLookupFunc func, icRenderingIntent intent,
+ icmLookupOrder order, int flags, int no, int nobw, cow *points,
+ icxMatrixModel *skm,
+ double dispLuminance, double wpscale, double smooth, double avgdev,
+ icxViewCond *vc, icxInk *ink, xcal *cal, int quality);
+static void icxLutSpaces(icxLuBase *p, icColorSpaceSignature *ins, int *inn,
+ icColorSpaceSignature *outs, int *outn,
+ icColorSpaceSignature *pcs);
+static void icxLuSpaces(icxLuBase *p, icColorSpaceSignature *ins, int *inn,
+ icColorSpaceSignature *outs, int *outn,
+ icmLuAlgType *alg, icRenderingIntent *intt,
+ icmLookupFunc *fnc, icColorSpaceSignature *pcs);
+static void icxLu_get_native_ranges (icxLuBase *p,
+ double *inmin, double *inmax, double *outmin, double *outmax);
+static void icxLu_get_ranges (icxLuBase *p,
+ double *inmin, double *inmax, double *outmin, double *outmax);
+static void icxLuEfv_wh_bk_points(icxLuBase *p, double *wht, double *blk, double *kblk);
+int xicc_get_viewcond(xicc *p, icxViewCond *vc);
+
+/* The different profile types are in their own source filesm */
+/* and are included to keep their functions private. (static) */
+#include "xmono.c"
+#include "xmatrix.c"
+#include "xlut.c" /* New xfit3 in & out optimising based profiles */
+//#include "xlut1.c" /* Old xfit1 device curve based profiles */
+
+#ifdef NT /* You'd think there might be some standards.... */
+# ifndef __BORLANDC__
+# define stricmp _stricmp
+# endif
+#else
+# define stricmp strcasecmp
+#endif
+/* Utilities */
+
+/* Return a string description of the given enumeration value */
+const char *icx2str(icmEnumType etype, int enumval) {
+
+ if (etype == icmColorSpaceSignature) {
+ if (((icColorSpaceSignature)enumval) == icxSigJabData)
+ return "Jab";
+ else if (((icColorSpaceSignature)enumval) == icxSigJChData)
+ return "JCh";
+ else if (((icColorSpaceSignature)enumval) == icxSigLChData)
+ return "LCh";
+ } else if (etype == icmRenderingIntent) {
+ if (((icRenderingIntent)enumval) == icxAppearance)
+ return "icxAppearance";
+ else if (((icRenderingIntent)enumval) == icxAbsAppearance)
+ return "icxAbsAppearance";
+ else if (((icRenderingIntent)enumval) == icxPerceptualAppearance)
+ return "icxPerceptualAppearance";
+ else if (((icRenderingIntent)enumval) == icxAbsPerceptualAppearance)
+ return "icxAbsPerceptualAppearance";
+ else if (((icRenderingIntent)enumval) == icxSaturationAppearance)
+ return "icxSaturationAppearance";
+ else if (((icRenderingIntent)enumval) == icxAbsSaturationAppearance)
+ return "icxAbsSaturationAppearance";
+ }
+ return icm2str(etype, enumval);
+}
+
+/* Common xicc stuff */
+
+/* Return information about the native lut in/out colorspaces. */
+/* Any pointer may be NULL if value is not to be returned */
+static void
+icxLutSpaces(
+ icxLuBase *p, /* This */
+ icColorSpaceSignature *ins, /* Return input color space */
+ int *inn, /* Return number of input components */
+ icColorSpaceSignature *outs, /* Return output color space */
+ int *outn, /* Return number of output components */
+ icColorSpaceSignature *pcs /* Return PCS color space */
+) {
+ p->plu->lutspaces(p->plu, ins, inn, outs, outn, pcs);
+}
+
+/* Return information about the overall lookup in/out colorspaces, */
+/* including allowance for any PCS override. */
+/* Any pointer may be NULL if value is not to be returned */
+static void
+icxLuSpaces(
+ icxLuBase *p, /* This */
+ icColorSpaceSignature *ins, /* Return input color space */
+ int *inn, /* Return number of input components */
+ icColorSpaceSignature *outs, /* Return output color space */
+ int *outn, /* Return number of output components */
+ icmLuAlgType *alg, /* Return type of lookup algorithm used */
+ icRenderingIntent *intt, /* Return the intent implemented */
+ icmLookupFunc *fnc, /* Return the profile function being implemented */
+ icColorSpaceSignature *pcs /* Return the effective PCS */
+) {
+ icmLookupFunc function;
+ icColorSpaceSignature npcs; /* Native PCS */
+
+ p->plu->spaces(p->plu, NULL, inn, NULL, outn, alg, NULL, &function, &npcs, NULL);
+
+ if (intt != NULL)
+ *intt = p->intent;
+
+ if (fnc != NULL)
+ *fnc = function;
+
+ if (ins != NULL)
+ *ins = p->ins;
+
+ if (outs != NULL)
+ *outs = p->outs;
+
+ if (pcs != NULL)
+ *pcs = p->pcs;
+}
+
+/* Return the native (internaly visible) colorspace value ranges */
+static void
+icxLu_get_native_ranges (
+icxLuBase *p,
+double *inmin, double *inmax, /* Return maximum range of inspace values */
+double *outmin, double *outmax /* Return maximum range of outspace values */
+) {
+ int i;
+ if (inmin != NULL) {
+ for (i = 0; i < p->inputChan; i++)
+ inmin[i] = p->ninmin[i];
+ }
+ if (inmax != NULL) {
+ for (i = 0; i < p->inputChan; i++)
+ inmax[i] = p->ninmax[i];
+ }
+ if (outmin != NULL) {
+ for (i = 0; i < p->outputChan; i++)
+ outmin[i] = p->noutmin[i];
+ }
+ if (outmax != NULL) {
+ for (i = 0; i < p->outputChan; i++)
+ outmax[i] = p->noutmax[i];
+ }
+}
+
+/* Return the effective (externaly visible) colorspace value ranges */
+static void
+icxLu_get_ranges (
+icxLuBase *p,
+double *inmin, double *inmax, /* Return maximum range of inspace values */
+double *outmin, double *outmax /* Return maximum range of outspace values */
+) {
+ int i;
+ if (inmin != NULL) {
+ for (i = 0; i < p->inputChan; i++)
+ inmin[i] = p->inmin[i];
+ }
+ if (inmax != NULL) {
+ for (i = 0; i < p->inputChan; i++)
+ inmax[i] = p->inmax[i];
+ }
+ if (outmin != NULL) {
+ for (i = 0; i < p->outputChan; i++)
+ outmin[i] = p->outmin[i];
+ }
+ if (outmax != NULL) {
+ for (i = 0; i < p->outputChan; i++)
+ outmax[i] = p->outmax[i];
+ }
+}
+
+/* - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+/* Routine to figure out a suitable black point for CMYK */
+
+/* Structure to hold optimisation information */
+typedef struct {
+ icmLuBase *p;
+ int kch; /* K channel, -1 if none */
+ double tlimit, klimit; /* Ink limit values */
+ int inn; /* Number of input channels */
+ icColorSpaceSignature outs; /* Output space */
+ double p1[3]; /* white pivot point in abs Lab */
+ double p2[3]; /* Point on vector towards black in abs Lab */
+ double toll; /* Tollerance of black direction */
+} bpfind;
+
+/* Optimise device values to minimise L, while remaining */
+/* within the ink limit, and staying in line between p1 (white) and p2 (black dir) */
+static double bpfindfunc(void *adata, double pv[]) {
+ bpfind *b = (bpfind *)adata;
+ double rv = 0.0;
+ double Lab[3];
+ double lr, ta, tb, terr; /* L ratio, target a, target b, target error */
+ double ovr = 0.0;
+ int e;
+
+ /* Compute amount outside total limit */
+ if (b->tlimit >= 0.0) {
+ double sum;
+ for (sum = 0.0, e = 0; e < b->inn; e++)
+ sum += pv[e];
+ if (sum > b->tlimit) {
+ ovr = sum - b->tlimit;
+#ifdef DEBUG
+ printf("~1 total ink ovr = %f\n",ovr);
+#endif
+ }
+ }
+
+ /* Compute amount outside black limit */
+ if (b->klimit >= 0.0 && b->kch >= 0) {
+ double kval = pv[b->kch] - b->klimit;
+ if (kval > ovr) {
+ ovr = kval;
+#ifdef DEBUG
+ printf("~1 black ink ovr = %f\n",ovr);
+#endif
+ }
+ }
+ /* Compute amount outside device value limits 0.0 - 1.0 */
+ {
+ double dval;
+ for (dval = -1.0, e = 0; e < b->inn; e++) {
+ if (pv[e] < 0.0) {
+ if (-pv[e] > dval)
+ dval = -pv[e];
+ } else if (pv[e] > 1.0) {
+ if ((pv[e] - 1.0) > dval)
+ dval = pv[e] - 1.0;
+ }
+ }
+ if (dval > ovr)
+ ovr = dval;
+ }
+
+ /* Compute the Lab value: */
+ b->p->lookup(b->p, Lab, pv);
+ if (b->outs == icSigXYZData)
+ icmXYZ2Lab(&icmD50, Lab, Lab);
+
+#ifdef DEBUG
+ printf("~1 p1 = %f %f %f, p2 = %f %f %f\n",b->p1[0],b->p1[1],b->p1[2],b->p2[0],b->p2[1],b->p2[2]);
+ printf("~1 device value %f %f %f %f, Lab = %f %f %f\n",pv[0],pv[1],pv[2],pv[3],Lab[0],Lab[1],Lab[2]);
+#endif
+
+ /* Primary aim is to minimise L value */
+ rv = Lab[0];
+
+ /* See how out of line from p1 to p2 we are */
+ lr = (Lab[0] - b->p1[0])/(b->p2[0] - b->p1[0]); /* Distance towards p2 from p1 */
+ ta = lr * (b->p2[1] - b->p1[1]) + b->p1[1]; /* Target a value */
+ tb = lr * (b->p2[2] - b->p1[2]) + b->p1[2]; /* Target b value */
+
+ terr = (ta - Lab[1]) * (ta - Lab[1])
+ + (tb - Lab[2]) * (tb - Lab[2]);
+
+ if (terr < b->toll) /* Tollerance error doesn't count until it's over tollerance */
+ terr = 0.0;
+
+#ifdef DEBUG
+ printf("~1 target error %f\n",terr);
+#endif
+ rv += XICC_BLACK_FIND_ABERR_WEIGHT * terr; /* Make ab match more important than min. L */
+
+#ifdef DEBUG
+ printf("~1 out of range error %f\n",ovr);
+#endif
+ rv += 200 * ovr;
+
+#ifdef DEBUG
+ printf("~1 black find tc ret %f\n",rv);
+#endif
+ return rv;
+}
+
+/* Try and compute a real black point in XYZ given an iccLu, */
+/* and also return the K only black or the normal black if the device doesn't have K */
+/* black[] will be unchanged if black cannot be computed. */
+/* Note that the black point will be in the space of the Lu */
+/* converted to XYZ, so will have the Lu's intent etc. */
+/* (Note that this is duplicated in xlut.c set_icxLuLut() !!!) */
+static void icxLu_comp_bk_point(
+icxLuBase *x,
+int gblk, /* If nz, compute black if possible. */
+double *white, /* XYZ Input, used for computing black */
+double *black, /* XYZ Input & Output. Set if gblk NZ and can be computed */
+double *kblack /* XYZ Output. Looked up if possible or set to black[] otherwise */
+) {
+ icmLuBase *p = x->plu;
+ icmLuBase *op = p; /* Original icmLu, in case we replace p */
+ icc *icco = p->icp;
+ icmHeader *h = icco->header;
+ icColorSpaceSignature ins, outs;
+ int inn, outn;
+ icmLuAlgType alg;
+ icRenderingIntent intt;
+ icmLookupFunc fnc;
+ icmLookupOrder ord;
+ int kch = -1;
+ double dblack[MAX_CHAN]; /* device black value */
+ int e;
+
+#ifdef DEBUG
+ printf("~1 icxLu_comp_bk_point() called, gblk %d, white = %s, black = %s\n",gblk,icmPdv(3, white),icmPdv(3,black));
+#endif
+ /* Default return incoming black as K only black */
+ kblack[0] = black[0];
+ kblack[1] = black[1];
+ kblack[2] = black[2];
+
+ /* Get the effective characteristics of the Lu */
+ p->spaces(p, &ins, &inn, &outs, &outn, &alg, &intt, &fnc, NULL, &ord);
+
+ if (fnc == icmBwd) { /* Hmm. We've got PCS to device, and we want device to PCS. */
+
+ /* Strictly speaking this is a dubious approach, since for a cLut profile */
+ /* the B2A table could make the effective white and black points */
+ /* anything it likes, and they don't have to match what the corresponding */
+ /* A2B table does. In our usage it's probably OK, since we tend */
+ /* to use colorimetric B2A */
+#ifdef DEBUG
+ printf("~1 getting icmFwd\n");
+#endif
+ if ((p = icco->get_luobj(icco, icmFwd, intt, ins, ord)) == NULL)
+ error("icxLu_comp_bk_point: assert: getting Fwd Lookup failed!");
+
+ p->spaces(p, &ins, &inn, &outs, &outn, &alg, &intt, &fnc, NULL, &ord);
+ }
+
+ if (outs != icSigXYZData && outs != icSigLabData) {
+ error("icxLu_comp_bk_point: assert: icc Lu output is not XYZ or Lab!, outs = 0x%x, ");
+ }
+
+#ifdef DEBUG
+ printf("~1 icxLu_comp_bk_point called for inn = %d, ins = %s\n", inn, icx2str(icmColorSpaceSignature,ins));
+#endif
+
+ switch (ins) {
+
+ case icSigXYZData:
+ case icSigLabData:
+ case icSigLuvData:
+ case icSigYxyData:
+#ifdef DEBUG
+ printf("~1 Assuming CIE colorspace black is 0.0\n");
+#endif
+ if (gblk) {
+ for (e = 0; e < inn; e++)
+ black[0] = 0.0;
+ }
+ kblack[0] = black[0];
+ kblack[1] = black[1];
+ kblack[2] = black[2];
+ return;
+
+ case icSigRgbData:
+#ifdef DEBUG
+ printf("~1 RGB:\n");
+#endif
+ for (e = 0; e < inn; e++)
+ dblack[e] = 0.0;
+ break;
+
+ case icSigGrayData: { /* Could be additive or subtractive */
+ double dval[1];
+ double minv[3], maxv[3];
+#ifdef DEBUG
+ printf("~1 Gray:\n");
+#endif
+ /* Check out 0 and 100% colorant */
+ dval[0] = 0.0;
+ p->lookup(p, minv, dval);
+ if (outs == icSigXYZData)
+ icmXYZ2Lab(&icmD50, minv, minv);
+ dval[0] = 1.0;
+ p->lookup(p, maxv, dval);
+ if (outs == icSigXYZData)
+ icmXYZ2Lab(&icmD50, maxv, maxv);
+
+ if (minv[0] < maxv[0])
+ dblack[0] = 0.0;
+ else
+ dblack[0] = 1.0;
+ }
+ break;
+
+ case icSigCmyData:
+ for (e = 0; e < inn; e++)
+ dblack[e] = 1.0;
+ break;
+
+ case icSigCmykData:
+#ifdef DEBUG
+ printf("~1 CMYK:\n");
+#endif
+ kch = 3;
+ dblack[0] = 0.0;
+ dblack[1] = 0.0;
+ dblack[2] = 0.0;
+ dblack[3] = 1.0;
+ if (alg == icmLutType) {
+ icxLuLut *pp = (icxLuLut *)x;
+
+ if (pp->ink.tlimit >= 0.0)
+ dblack[kch] = pp->ink.tlimit;
+ };
+ break;
+
+ /* Use a heursistic. */
+ /* This duplicates code in icxGetLimits() :-( */
+ /* Colorant guessing should go in icclib ? */
+ case icSig2colorData:
+ case icSig3colorData:
+ case icSig4colorData:
+ case icSig5colorData:
+ case icSig6colorData:
+ case icSig7colorData:
+ case icSig8colorData:
+ case icSig9colorData:
+ case icSig10colorData:
+ case icSig11colorData:
+ case icSig12colorData:
+ case icSig13colorData:
+ case icSig14colorData:
+ case icSig15colorData:
+ case icSigMch5Data:
+ case icSigMch6Data:
+ case icSigMch7Data:
+ case icSigMch8Data: {
+ double dval[MAX_CHAN];
+ double ncval[3];
+ double cvals[MAX_CHAN][3];
+ int nlighter, ndarker;
+
+ /* Decide if the colorspace is additive or subtractive */
+#ifdef DEBUG
+ printf("~1 N channel:\n");
+#endif
+
+ /* First the no colorant value */
+ for (e = 0; e < inn; e++)
+ dval[e] = 0.0;
+ p->lookup(p, ncval, dval);
+ if (outs == icSigXYZData)
+ icmXYZ2Lab(&icmD50, ncval, ncval);
+
+ /* Then all the colorants */
+ nlighter = ndarker = 0;
+ for (e = 0; e < inn; e++) {
+ dval[e] = 1.0;
+ p->lookup(p, cvals[e], dval);
+ if (outs == icSigXYZData)
+ icmXYZ2Lab(&icmD50, cvals[e], cvals[e]);
+ dval[e] = 0.0;
+ if (fabs(cvals[e][0] - ncval[0]) > 5.0) {
+ if (cvals[e][0] > ncval[0])
+ nlighter++;
+ else
+ ndarker++;
+ }
+ }
+ if (ndarker == 0 && nlighter > 0) { /* Assume additive */
+ for (e = 0; e < inn; e++)
+ dblack[e] = 0.0;
+#ifdef DEBUG
+ printf("~1 N channel is additive:\n");
+#endif
+
+ } else if (ndarker > 0 && nlighter == 0) { /* Assume subtractive. */
+ double pbk[3] = { 0.0,0.0,0.0 }; /* Perfect black */
+ double smd = 1e10; /* Smallest distance */
+
+#ifdef DEBUG
+ printf("~1 N channel is subtractive:\n");
+#endif
+ /* See if we can guess the black channel */
+ for (e = 0; e < inn; e++) {
+ double tt;
+ tt = icmNorm33sq(pbk, cvals[e]);
+ if (tt < smd) {
+ smd = tt;
+ kch = e;
+ }
+ }
+ /* See if the black seems sane */
+ if (cvals[kch][0] > 40.0
+ || fabs(cvals[kch][1]) > 10.0
+ || fabs(cvals[kch][2]) > 10.0) {
+ if (p != op)
+ p->del(p);
+#ifdef DEBUG
+ printf("~1 black doesn't look sanem so assume nothing\n");
+#endif
+ return; /* Assume nothing */
+ }
+
+ /* Chosen kch as black */
+ for (e = 0; e < inn; e++)
+ dblack[e] = 0.0;
+ dblack[kch] = 1.0;
+ if (alg == icmLutType) {
+ icxLuLut *pp = (icxLuLut *)x;
+
+ if (pp->ink.tlimit >= 0.0)
+ dblack[kch] = pp->ink.tlimit;
+ };
+#ifdef DEBUG
+ printf("~1 N channel K = chan %d\n",kch);
+#endif
+ } else {
+ if (p != op)
+ p->del(p);
+#ifdef DEBUG
+ printf("~1 can't figure if additive or subtractive, so assume nothing\n");
+#endif
+ return; /* Assume nothing */
+ }
+ }
+ break;
+
+ default:
+#ifdef DEBUG
+ printf("~1 unhandled colorspace, so assume nothing\n");
+#endif
+ if (p != op)
+ p->del(p);
+ return; /* Don't do anything */
+ }
+
+ /* Lookup the K only value */
+ if (kch >= 0) {
+ p->lookup(p, kblack, dblack);
+
+ /* We always return XYZ */
+ if (outs == icSigLabData)
+ icmLab2XYZ(&icmD50, kblack, kblack);
+ }
+
+ if (gblk == 0) { /* That's all we have to do */
+#ifdef DEBUG
+ printf("~1 gblk == 0, so only return kblack\n");
+#endif
+ if (p != op)
+ p->del(p);
+ return;
+ }
+
+ /* Lookup the device black or K only value as a default */
+ p->lookup(p, black, dblack); /* May be XYZ or Lab */
+#ifdef DEBUG
+ printf("~1 Got default lu black %f %f %f, kch = %d\n", black[0],black[1],black[2],kch);
+#endif
+
+ /* !!! Hmm. For CMY and RGB we are simply using the device */
+ /* combination values as the black point. In reality we might */
+ /* want to have the option of using a neutral black point, */
+ /* just like CMYK ?? */
+
+ if (kch >= 0) { /* The space is subtractive with a K channel. */
+ /* If XICC_NEUTRAL_CMYK_BLACK then locate the darkest */
+ /* CMYK within limits with the same chromaticity as the white point, */
+ /* otherwise locate the device value within the ink limits that is */
+ /* in the direction of the K channel */
+ bpfind bfs; /* Callback context */
+ double sr[MXDO]; /* search radius */
+ double tt[MXDO]; /* Temporary */
+ double rs0[MXDO], rs1[MXDO]; /* Random start candidates */
+ int trial;
+ double brv;
+
+ /* Setup callback function context */
+ bfs.p = p;
+ bfs.inn = inn;
+ bfs.outs = outs;
+
+ bfs.kch = kch;
+ bfs.tlimit = -1.0;
+ bfs.klimit = -1.0;
+ bfs.toll = XICC_BLACK_POINT_TOLL;
+
+ if (alg == icmLutType) {
+ icxLuLut *pp = (icxLuLut *)x;
+
+ pp->kch = kch;
+ bfs.tlimit = pp->ink.tlimit;
+ bfs.klimit = pp->ink.klimit;
+#ifdef DEBUG
+ printf("~1 tlimit = %f, klimit = %f\n",bfs.tlimit,bfs.klimit);
+#endif
+ };
+
+#ifdef XICC_NEUTRAL_CMYK_BLACK
+#ifdef DEBUG
+ printf("~1 Searching for neutral black\n");
+#endif
+ /* white has been given to us in XYZ */
+ icmXYZ2Lab(&icmD50, bfs.p1, white); /* pivot Lab */
+ icmCpy3(bfs.p2, white); /* temp white XYZ */
+ icmScale3(bfs.p2, bfs.p2, 0.02); /* Scale white XYZ towards 0,0,0 */
+ icmXYZ2Lab(&icmD50, bfs.p2, bfs.p2); /* Convert black direction to Lab */
+
+#else /* Use K directin black */
+#ifdef DEBUG
+ printf("~1 Searching for K direction black\n");
+#endif
+ icmXYZ2Lab(&icmD50, bfs.p1, white); /* Pivot */
+
+ /* Now figure abs Lab value of K only, as the direction */
+ /* to use for the rich black. */
+ for (e = 0; e < inn; e++)
+ dblack[e] = 0.0;
+ if (bfs.klimit < 0.0)
+ dblack[kch] = 1.0;
+ else
+ dblack[kch] = bfs.klimit; /* K value */
+
+ p->lookup(p, black, dblack);
+
+ if (outs == icSigXYZData) {
+ icmXYZ2Lab(&icmD50, bfs.p2, black); /* K direction */
+ } else {
+ icmAry2Ary(bfs.p2, black);
+ }
+#endif
+
+#ifdef DEBUG
+ printf("~1 Lab pivot %f %f %f, Lab K direction %f %f %f\n",bfs.p1[0],bfs.p1[1],bfs.p1[2],bfs.p2[0],bfs.p2[1],bfs.p2[2]);
+#endif
+ /* Start with the K only as the current best value */
+ brv = bpfindfunc((void *)&bfs, dblack);
+#ifdef DEBUG
+ printf("~1 initial brv for K only = %f\n",brv);
+#endif
+
+ /* Set the random start 0 location as 000K */
+ /* and the random start 1 location as CMY0 */
+ {
+ double tt;
+
+ for (e = 0; e < inn; e++)
+ dblack[e] = rs0[e] = 0.0;
+ if (bfs.klimit < 0.0)
+ dblack[kch] = rs0[kch] = 1.0;
+ else
+ dblack[kch] = rs0[kch] = bfs.klimit; /* K value */
+
+ if (bfs.tlimit < 0.0)
+ tt = 1.0;
+ else
+ tt = bfs.tlimit/(inn - 1.0);
+ for (e = 0; e < inn; e++)
+ rs1[e] = tt;
+ rs1[kch] = 0.0; /* K value */
+ }
+
+ /* Find the device black point using optimization */
+ /* Do several trials to avoid local minima. */
+ rand32(0x12345678); /* Make trial values deterministic */
+ for (trial = 0; trial < 200; trial++) {
+ double rv; /* Temporary */
+
+ /* Start first trial at 000K */
+ if (trial == 0) {
+ for (e = 0; e < inn; e++) {
+ tt[e] = rs0[e];
+ sr[e] = 0.1;
+ }
+
+ } else {
+ /* Base is random between 000K and CMY0: */
+ if (trial < 100) {
+ rv = d_rand(0.0, 1.0);
+ for (e = 0; e < inn; e++) {
+ tt[e] = rv * rs0[e] + (1.0 - rv) * rs1[e];
+ sr[e] = 0.1;
+ }
+ /* Base on current best */
+ } else {
+ for (e = 0; e < inn; e++) {
+ tt[e] = dblack[e];
+ sr[e] = 0.1;
+ }
+ }
+
+ /* Then add random start offset */
+ for (rv = 0.0, e = 0; e < inn; e++) {
+ tt[e] += d_rand(-0.5, 0.5);
+ if (tt[e] < 0.0)
+ tt[e] = 0.0;
+ else if (tt[e] > 1.0)
+ tt[e] = 1.0;
+ }
+ }
+
+ /* Clip black */
+ if (bfs.klimit >= 0.0 && tt[kch] > bfs.klimit)
+ tt[kch] = bfs.klimit;
+
+ /* Compute amount outside total limit */
+ if (bfs.tlimit >= 0.0) {
+ for (rv = 0.0, e = 0; e < inn; e++)
+ rv += tt[e];
+
+ if (rv > bfs.tlimit) {
+ rv /= (double)inn;
+ for (e = 0; e < inn; e++)
+ tt[e] -= rv;
+ }
+ }
+
+ if (powell(&rv, inn, tt, sr, 0.000001, 1000, bpfindfunc,
+ (void *)&bfs, NULL, NULL) == 0) {
+#ifdef DEBUG
+ printf("~1 trial %d, rv %f bp %f %f %f %f\n",trial,rv,tt[0],tt[1],tt[2],tt[3]);
+#endif
+ if (rv < brv) {
+#ifdef DEBUG
+ printf("~1 new best\n");
+#endif
+ brv = rv;
+ for (e = 0; e < inn; e++)
+ dblack[e] = tt[e];
+ }
+ }
+ }
+ if (brv > 1000.0)
+ error("icxLu_comp_bk_point: Black point powell failed");
+
+ for (e = 0; e < inn; e++) { /* Make sure device values are in range */
+ if (dblack[e] < 0.0)
+ dblack[e] = 0.0;
+ else if (dblack[e] > 1.0)
+ dblack[e] = 1.0;
+ }
+ /* Now have device black in dblack[] */
+#ifdef DEBUG
+ printf("~1 got device black %f %f %f %f\n",dblack[0], dblack[1], dblack[2], dblack[3]);
+#endif
+
+ p->lookup(p, black, dblack); /* Convert to PCS */
+ }
+
+ if (p != op)
+ p->del(p);
+
+ /* We always return XYZ */
+ if (outs == icSigLabData)
+ icmLab2XYZ(&icmD50, black, black);
+
+#ifdef DEBUG
+ printf("~1 returning %f %f %f\n", black[0], black[1], black[2]);
+#endif
+
+ return;
+}
+
+/* - - - - - - - - - - - - - - - - - - - - - - - */
+
+/* Return the media white and black points */
+/* in the xlu effective PCS colorspace. Pointers may be NULL. */
+/* (ie. these will be relative values for relative intent etc.) */
+static void icxLuEfv_wh_bk_points(
+icxLuBase *p,
+double *wht,
+double *blk,
+double *kblk /* K only black */
+) {
+ double white[3], black[3], kblack[3];
+
+ /* Get the Lu PCS converted to XYZ icc black and white points in XYZ */
+ if (p->plu->lu_wh_bk_points(p->plu, white, black)) {
+ /* Black point is assumed. We should determine one instead. */
+ /* Lookup K only black too */
+ icxLu_comp_bk_point(p, 1, white, black, kblack);
+
+ } else {
+ /* Lookup a possible K only black */
+ icxLu_comp_bk_point(p, 0, white, black, kblack);
+ }
+
+//printf("~1 white %f %f %f, black %f %f %f, kblack %f %f %f\n",white[0],white[1],white[2],black[0],black[1],black[2],kblack[0],kblack[1],kblack[2]);
+ /* Convert to possibl xicc override PCS */
+ switch (p->pcs) {
+ case icSigXYZData:
+ break; /* Don't have to do anyting */
+ case icSigLabData:
+ icmXYZ2Lab(&icmD50, white, white); /* Convert from XYZ to Lab */
+ icmXYZ2Lab(&icmD50, black, black);
+ icmXYZ2Lab(&icmD50, kblack, kblack);
+ break;
+ case icxSigJabData:
+ p->cam->XYZ_to_cam(p->cam, white, white); /* Convert from XYZ to Jab */
+ p->cam->XYZ_to_cam(p->cam, black, black);
+ p->cam->XYZ_to_cam(p->cam, kblack, kblack);
+ break;
+ default:
+ break;
+ }
+
+//printf("~1 icxLuEfv_wh_bk_points: pcsor %s White %f %f %f, Black %f %f %f\n", icx2str(icmColorSpaceSignature,p->pcs), white[0], white[1], white[2], black[0], black[1], black[2]);
+ if (wht != NULL) {
+ wht[0] = white[0];
+ wht[1] = white[1];
+ wht[2] = white[2];
+ }
+
+ if (blk != NULL) {
+ blk[0] = black[0];
+ blk[1] = black[1];
+ blk[2] = black[2];
+ }
+
+ if (kblk != NULL) {
+ kblk[0] = kblack[0];
+ kblk[1] = kblack[1];
+ kblk[2] = kblack[2];
+ }
+}
+
+/* Create an instance of an xicc object */
+xicc *new_xicc(
+icc *picc /* icc we are expanding */
+) {
+ xicc *p;
+ if ((p = (xicc *) calloc(1,sizeof(xicc))) == NULL)
+ return NULL;
+ p->pp = picc;
+ p->del = xicc_del;
+ p->get_luobj = xicc_get_luobj;
+ p->set_luobj = xicc_set_luobj;
+ p->get_viewcond = xicc_get_viewcond;
+
+ /* Create an xcal if there is the right tag in the profile */
+ p->cal = xiccReadCalTag(p->pp);
+ p->nodel_cal = 0; /* We created it, we will delete it */
+
+ return p;
+}
+
+/* Do away with the xicc (but not the icc!) */
+static void xicc_del(
+xicc *p
+) {
+ if (p->cal != NULL && p->nodel_cal == 0)
+ p->cal->del(p->cal);
+ free (p);
+}
+
+
+/* Return an expanded lookup object, initialised */
+/* from the icc. */
+/* Return NULL on error, check errc+err for reason. */
+/* Set the pcsor & intent to consistent and values if */
+/* Jab and/or icxAppearance has been requested. */
+/* Create the underlying icm lookup object that is used */
+/* to create and implement the icx one. The icm will be used */
+/* to translate from native to effective PCS, unless the */
+/* effective PCS is Jab, in which case the icm will be set to */
+/* have an effective PCS of XYZ. Since native<->effecive PCS conversion */
+/* is done at the to/from_abs() stage, none of this affects the individual */
+/* conversion steps, which will all talk the native PCS (unless merged). */
+icxLuBase *xicc_get_luobj(
+xicc *p, /* this */
+int flags, /* clip, merge flags */
+icmLookupFunc func, /* Functionality */
+icRenderingIntent intent, /* Intent */
+icColorSpaceSignature pcsor,/* PCS override (0 = def) */
+icmLookupOrder order, /* Search Order */
+icxViewCond *vc, /* Viewing Condition (may be NULL if pcsor is not CIECAM) */
+icxInk *ink /* inking details (NULL for default) */
+) {
+ icmLuBase *plu;
+ icxLuBase *xplu;
+ icmLuAlgType alg;
+ icRenderingIntent n_intent = intent; /* Native Intent to request */
+ icColorSpaceSignature n_pcs = icmSigDefaultData; /* Native PCS to request */
+
+//printf("~1 xicc_get_luobj got intent %s and pcsor %s\n",icx2str(icmRenderingIntent,intent),icx2str(icmColorSpaceSignature,pcsor));
+
+ /* Ensure that appropriate PCS is slected for an appearance intent */
+ if (intent == icxAppearance
+ || intent == icxAbsAppearance
+ || intent == icxPerceptualAppearance
+ || intent == icxAbsPerceptualAppearance
+ || intent == icxSaturationAppearance
+ || intent == icxAbsSaturationAppearance) {
+ pcsor = icxSigJabData;
+
+ /* Translate non-Jab intents to the equivalent appearance "intent" if pcsor == Jab. */
+ /* This is how we get these when the UI's don't list all the apperances intents, */
+ /* we select the analogous non-apperance intent with pcsor = Jab. */
+ /* Note that Abs/non-abs selects between Apperance and AbsAppearance. */
+ } else if (pcsor == icxSigJabData) {
+ if (intent == icRelativeColorimetric)
+ intent = icxAppearance;
+ else if (intent == icAbsoluteColorimetric)
+ intent = icxAbsAppearance;
+ else if (intent == icPerceptual)
+ intent = icxPerceptualAppearance;
+ else if (intent == icmAbsolutePerceptual)
+ intent = icxAbsPerceptualAppearance;
+ else if (intent == icSaturation)
+ intent = icxSaturationAppearance;
+ else if (intent == icmAbsoluteSaturation)
+ intent = icxAbsSaturationAppearance;
+ else
+ intent = icxAppearance;
+ }
+
+ /* Translate intent asked for into intent needed in icclib */
+ if (intent == icxAppearance
+ || intent == icxAbsAppearance)
+ n_intent = icAbsoluteColorimetric;
+ else if (intent == icxPerceptualAppearance
+ || intent == icxAbsPerceptualAppearance)
+ n_intent = icmAbsolutePerceptual;
+ else if (intent == icxSaturationAppearance
+ || intent == icxAbsSaturationAppearance)
+ n_intent = icmAbsoluteSaturation;
+
+ if (pcsor != icmSigDefaultData)
+ n_pcs = pcsor; /* There is an icclib override */
+
+ if (pcsor == icxSigJabData) /* xicc override */
+ n_pcs = icSigXYZData; /* Translate to XYZ */
+
+//printf("~1 xicc_get_luobj processed intent %s and pcsor %s\n",icx2str(icmRenderingIntent,intent),icx2str(icmColorSpaceSignature,pcsor));
+//printf("~1 xicc_get_luobj icclib intent %s and pcsor %s\n",icx2str(icmRenderingIntent,n_intent),icx2str(icmColorSpaceSignature,n_pcs));
+ /* Get icclib lookup object */
+ if ((plu = p->pp->get_luobj(p->pp, func, n_intent, n_pcs, order)) == NULL) {
+ p->errc = p->pp->errc; /* Copy error */
+ strcpy(p->err, p->pp->err);
+ return NULL;
+ }
+
+ /* Figure out what the algorithm is */
+ plu->spaces(plu, NULL, NULL, NULL, NULL, &alg, NULL, NULL, &n_pcs, NULL);
+
+ /* make sure its "Abs CAM" */
+ if (vc!= NULL
+ && (intent == icxAbsAppearance
+ || intent == icxAbsPerceptualAppearance
+ || intent == icxAbsSaturationAppearance)) { /* make sure its "Abs CAM" */
+ /* Set white point and flare color to D50 */
+ /* (Hmm. This doesn't match what happens within collink with absolute intent!!) */
+ vc->Wxyz[0] = icmD50.X/icmD50.Y;
+ vc->Wxyz[1] = icmD50.Y/icmD50.Y; // Normalise white reference to Y = 1 ?
+ vc->Wxyz[2] = icmD50.Z/icmD50.Y;
+
+ vc->Fxyz[0] = icmD50.X;
+ vc->Fxyz[1] = icmD50.Y;
+ vc->Fxyz[2] = icmD50.Z;
+ }
+
+ /* Call xiccLu wrapper creation */
+ switch (alg) {
+ case icmMonoFwdType:
+ xplu = new_icxLuMono(p, flags, plu, func, intent, pcsor, vc, 0);
+ break;
+ case icmMonoBwdType:
+ xplu = new_icxLuMono(p, flags, plu, func, intent, pcsor, vc, 1);
+ break;
+ case icmMatrixFwdType:
+ xplu = new_icxLuMatrix(p, flags, plu, func, intent, pcsor, vc, 0);
+ break;
+ case icmMatrixBwdType:
+ xplu = new_icxLuMatrix(p, flags, plu, func, intent, pcsor, vc, 1);
+ break;
+ case icmLutType:
+ xplu = new_icxLuLut(p, flags, plu, func, intent, pcsor, vc, ink);
+ break;
+ default:
+ xplu = NULL;
+ break;
+ }
+
+ return xplu;
+}
+
+
+/* Return an expanded lookup object, initialised */
+/* from the icc, and then overwritten by a conversion */
+/* created from the supplied scattered data points. */
+/* The Lut is assumed to be a device -> native PCS profile. */
+/* If the SET_WHITE and/or SET_BLACK flags are set, */
+/* discover the white/black point, set it in the icc, */
+/* and make the Lut relative to them. */
+/* Return NULL on error, check errc+err for reason */
+static icxLuBase *xicc_set_luobj(
+xicc *p, /* this */
+icmLookupFunc func, /* Functionality */
+icRenderingIntent intent, /* Intent */
+icmLookupOrder order, /* Search Order */
+int flags, /* white/black point, verbose flags etc. */
+int no, /* Number of points */
+int nobw, /* Number of points to look for white & black patches in */
+cow *points, /* Array of input points in target PCS space */
+icxMatrixModel *skm, /* Optional skeleton model (used for input profiles) */
+double dispLuminance, /* > 0.0 if display luminance value and is known */
+double wpscale, /* > 0.0 if input white point is to be scaled */
+double smooth, /* RSPL smoothing factor, -ve if raw */
+double avgdev, /* reading Average Deviation as a proportion of the input range */
+icxViewCond *vc, /* Viewing Condition (NULL if not using CAM) */
+icxInk *ink, /* inking details (NULL for default) */
+xcal *cal, /* Optional cal, will override any existing (not deleted with xicc)*/
+int quality /* Quality metric, 0..3 */
+) {
+ icmLuBase *plu;
+ icxLuBase *xplu = NULL;
+ icmLuAlgType alg;
+
+ if (cal != NULL) {
+ if (p->cal != NULL && p->nodel_cal == 0)
+ p->cal->del(p->cal);
+ p->cal = cal;
+ p->nodel_cal = 1; /* We were given it, so don't delete it */
+ }
+
+ if (func != icmFwd) {
+ p->errc = 1;
+ sprintf(p->err,"Can only create Device->PCS profiles from scattered data.");
+ xplu = NULL;
+ return xplu;
+ }
+
+ /* Get icclib lookup object */
+ if ((plu = p->pp->get_luobj(p->pp, func, intent, 0, order)) == NULL) {
+ p->errc = p->pp->errc; /* Copy error */
+ strcpy(p->err, p->pp->err);
+ return NULL;
+ }
+
+ /* Figure out what the algorithm is */
+ plu->spaces(plu, NULL, NULL, NULL, NULL, &alg, NULL, NULL, NULL, NULL);
+
+ /* Call xiccLu wrapper creation */
+ switch (alg) {
+ case icmMonoFwdType:
+ p->errc = 1;
+ sprintf(p->err,"Setting Monochrome Fwd profile from scattered data not supported.");
+ plu->del(plu);
+ xplu = NULL; /* Not supported yet */
+ break;
+
+ case icmMatrixFwdType:
+ if (smooth < 0.0)
+ smooth = -smooth;
+ xplu = set_icxLuMatrix(p, plu, flags, no, nobw, points, skm, dispLuminance, wpscale, quality, smooth);
+ break;
+
+ case icmLutType:
+ /* ~~~ Should add check that it is a fwd profile ~~~ */
+ xplu = set_icxLuLut(p, plu, func, intent, flags, no, nobw, points, skm, dispLuminance, wpscale, smooth, avgdev, vc, ink, quality);
+ break;
+
+ default:
+ break;
+ }
+
+ return xplu;
+}
+
+/* ------------------------------------------------------ */
+/* Viewing Condition Parameter stuff */
+
+#ifdef NEVER /* Not currently used */
+
+/* Guess viewing parameters from the technology signature */
+static void guess_from_techsig(
+icTechnologySignature tsig,
+double *Ybp
+) {
+double Yb = -1.0;
+
+ switch (tsig) {
+ /* These are all inputing either a representation of */
+ /* a natural scene captured on another medium, or are assuming */
+ /* that the medium is the original. A _good_ system would */
+ /* let the user indicate which is the case. */
+ case icSigReflectiveScanner:
+ case icSigFilmScanner:
+ Yb = 0.2;
+ break;
+
+ /* Direct scene to value devices. */
+ case icSigDigitalCamera:
+ case icSigVideoCamera:
+ Yb = 0.2;
+ break;
+
+ /* Emmisive displays. */
+ /* We could try tweaking the white point on the assumption */
+ /* that the viewer will be adapted to a combination of both */
+ /* the CRT white point, and the ambient light. */
+ case icSigVideoMonitor:
+ case icSigCRTDisplay:
+ case icSigPMDisplay:
+ case icSigAMDisplay:
+ Yb = 0.2;
+ break;
+
+ /* Photo CD has its own viewing definitions */
+ /* (It represents original scene colors) */
+ case icSigPhotoCD:
+ Yb = 0.2;
+ break;
+
+ /* Projection devices, either direct, or */
+ /* via another intermediate medium. */
+ case icSigProjectionTelevision:
+ Yb = 0.1; /* Assume darkened room, little background */
+ break;
+ case icSigFilmWriter:
+ Yb = 0.0; /* Assume a dark room - no background */
+ break;
+
+ /* Printed media devices. */
+ case icSigInkJetPrinter:
+ case icSigThermalWaxPrinter:
+ case icSigElectrophotographicPrinter:
+ case icSigElectrostaticPrinter:
+ case icSigDyeSublimationPrinter:
+ case icSigPhotographicPaperPrinter:
+ case icSigPhotoImageSetter:
+ case icSigGravure:
+ case icSigOffsetLithography:
+ case icSigSilkscreen:
+ case icSigFlexography:
+ Yb = 0.2;
+ break;
+
+ default:
+ Yb = 0.2;
+ }
+
+ if (Ybp != NULL)
+ *Ybp = Yb;
+}
+
+#endif /* NEVER */
+
+
+/* See if we can read or guess the viewing conditions */
+/* for an ICC profile. */
+/* Return value 0 if it is well defined */
+/* Return value 1 if it is a guess */
+/* Return value 2 if it is not possible/appropriate */
+int xicc_get_viewcond(
+xicc *p, /* Expanded profile we're working with */
+icxViewCond *vc /* Viewing parameters to return */
+) {
+ icc *pp = p->pp; /* Base ICC */
+
+ /* Numbers we're trying to find */
+ ViewingCondition Ev = vc_none;
+ double Wxyz[3] = {-1.0, -1.0, -1.0}; /* Adapting white color */
+ double La = -1.0; /* Adapting luminance */
+ double Ixyz[3] = {-1.0, -1.0, -1.0}; /* Illuminant color */
+ double Li = -1.0; /* Illuminant luminance */
+ double Lb = -1.0; /* Backgrount luminance */
+ double Yb = -1.0; /* Background relative luminance to Lv */
+ double Lve = -1.0; /* Emissive device image luminance */
+ double Lvr = -1.0; /* Reflective device image luminance */
+ double Lv = -1.0; /* device image luminance */
+ double Yf = -1.0; /* Flare relative luminance to Lv */
+ double Fxyz[3] = {-1.0, -1.0, -1.0}; /* Flare color */
+ icTechnologySignature tsig = icMaxEnumTechnology; /* Technology Signature */
+ icProfileClassSignature devc = icMaxEnumClass;
+ int trans = -1; /* Set to 0 if not transparency, 1 if it is */
+
+ /* Collect all the information we can find */
+
+ /* Emmisive devices image white luminance */
+ {
+ icmXYZArray *luminanceTag;
+
+ if ((luminanceTag = (icmXYZArray *)pp->read_tag(pp, icSigLuminanceTag)) != NULL
+ && luminanceTag->ttype == icSigXYZType && luminanceTag->size >= 1) {
+ Lve = luminanceTag->data[0].Y; /* Copy structure */
+ }
+ }
+
+ /* Flare: */
+ {
+ icmMeasurement *ro;
+
+ if ((ro = (icmMeasurement *)pp->read_tag(pp, icSigMeasurementTag)) != NULL
+ && ro->ttype == icSigMeasurementType) {
+
+ Yf = ro->flare;
+ /* ro->illuminant ie D50, D65, D93, A etc. */
+ }
+ }
+
+ /* Media White Point */
+ {
+ icmXYZArray *whitePointTag;
+
+ if ((whitePointTag = (icmXYZArray *)pp->read_tag(pp, icSigMediaWhitePointTag)) != NULL
+ && whitePointTag->ttype == icSigXYZType && whitePointTag->size >= 1) {
+ Wxyz[0] = whitePointTag->data[0].X;
+ Wxyz[1] = whitePointTag->data[0].Y;
+ Wxyz[2] = whitePointTag->data[0].Z;
+ }
+ }
+
+ /* ViewingConditions: */
+ {
+ icmViewingConditions *ro;
+
+ if ((ro = (icmViewingConditions *)pp->read_tag(pp, icSigViewingConditionsTag)) != NULL
+ && ro->ttype == icSigViewingConditionsType) {
+
+ /* ro->illuminant.X */
+ /* ro->illuminant.Z */
+
+ Li = ro->illuminant.Y;
+
+ /* Reflect illuminant off the media white */
+ Lvr = Li * Wxyz[1];
+
+ /* Illuminant color */
+ Ixyz[0] = ro->illuminant.X/ro->illuminant.Y;
+ Ixyz[1] = 1.0;
+ Ixyz[2] = ro->illuminant.Z/ro->illuminant.Y;
+
+ /* Assume ICC surround is CICAM97 background */
+ /* ro->surround.X */
+ /* ro->surround.Z */
+ La = ro->surround.Y;
+
+ /* ro->stdIlluminant ie D50, D65, D93, A etc. */
+ }
+ }
+
+ /* Stuff we might need */
+
+ /* Technology: */
+ {
+ icmSignature *ro;
+
+ /* Try and read the tag from the file */
+ if ((ro = (icmSignature *)pp->read_tag(pp, icSigTechnologyTag)) != NULL
+ && ro->ttype != icSigSignatureType) {
+
+ tsig = ro->sig;
+ }
+ }
+
+ devc = pp->header->deviceClass; /* Type of profile */
+ if (devc == icSigLinkClass
+ || devc == icSigAbstractClass
+ || devc == icSigColorSpaceClass
+ || devc == icSigNamedColorClass)
+ return 2;
+
+ /*
+ icSigInputClass
+ icSigDisplayClass
+ icSigOutputClass
+ */
+
+ if ((pp->header->flags & icTransparency) != 0)
+ trans = 1;
+ else
+ trans = 0;
+
+
+ /* figure Lv if we have the information */
+ if (Lve >= 0.0)
+ Lv = Lve; /* Emmisive image white luminance */
+ else
+ Lv = Lvr; /* Reflectance image white luminance */
+
+ /* Fudge the technology signature */
+ if (tsig == icMaxEnumTechnology) {
+ if (devc == icSigDisplayClass)
+ tsig = icSigCRTDisplay;
+ }
+
+#ifndef NEVER
+ printf("Enumeration = %d\n", Ev);
+ printf("Viewing Conditions:\n");
+ printf("White adaptation color %f %f %f\n",Wxyz[0], Wxyz[1], Wxyz[2]);
+ printf("Adapting Luminance La = %f\n",La);
+ printf("Illuminant color %f %f %f\n",Ixyz[0], Ixyz[1], Ixyz[2]);
+ printf("Illuminant Luminance Li = %f\n",Li);
+ printf("Background Luminance Lb = %f\n",Lb);
+ printf("Relative Background Yb = %f\n",Yb);
+ printf("Emissive Image White Lve = %f\n",Lve);
+ printf("Reflective Image White Lvr = %f\n",Lvr);
+ printf("Device Image White Lv = %f\n",Lv);
+ printf("Relative Flare Yf = %f\n",Yf);
+ printf("Flare color %f %f %f\n",Fxyz[0], Fxyz[1], Fxyz[2]);
+ printf("Technology = %s\n",tag2str(tsig));
+ printf("deviceClass = %s\n",tag2str(devc));
+ printf("Transparency = %d\n",trans);
+#endif
+
+ /* See if the viewing conditions are completely defined as ICC can do it */
+ if (Wxyz[0] >= 0.0 && Wxyz[1] >= 0.0 && Wxyz[2] >= 0.0
+ && La >= 0.0
+ && Yb >= 0.0
+ && Lv >= 0.0
+ && Yf >= 0.0
+ && Fxyz[0] >= 0.0 && Fxyz[1] >= 0.0 && Fxyz[2] >= 0.0) {
+
+ vc->Ev = vc_none;
+ vc->Wxyz[0] = Wxyz[0];
+ vc->Wxyz[1] = Wxyz[1];
+ vc->Wxyz[2] = Wxyz[2];
+ vc->La = La;
+ vc->Yb = Yb;
+ vc->Lv = Lv;
+ vc->Yf = Yf;
+ vc->Fxyz[0] = Fxyz[0];
+ vc->Fxyz[1] = Fxyz[1];
+ vc->Fxyz[2] = Fxyz[2];
+ return 0;
+ }
+
+ /* Hmm. We didn't get all the info an ICC can contain. */
+ /* We will try to guess some reasonable defaults */
+
+ /* Have we at least got an adaptation white point ? */
+ if (Wxyz[0] < 0.0 || Wxyz[1] < 0.0 || Wxyz[2] < 0.0)
+ return 2; /* No */
+
+ /* Have we got the technology ? */
+ if (tsig == icMaxEnumTechnology)
+ return 2; /* Hopeless */
+
+ /* Guess from the technology */
+ switch (tsig) {
+
+ /* This is inputing either a representation of */
+ /* a natural scene captured on another a print medium, or */
+ /* are is assuming that the medium is the original. */
+ /* We will assume that the print is the original. */
+ case icSigReflectiveScanner:
+ {
+ if (La < 0.0) /* No adapting luminance */
+ La = 34.0; /* Use a practical print evaluation number */
+ if (Yb < 0.0) /* No background relative luminance */
+ Yb = 0.2; /* Assume grey world */
+ if (Lv < 0.0) /* No device image luminance */
+ Ev = vc_average; /* Assume average viewing conditions */
+ if (Yf < 0.0) /* No flare figure */
+ Yf = 0.01; /* Assume 1% flare */
+ if (Fxyz[0] < 0.0 || Fxyz[1] < 0.0 || Fxyz[2] < 0.0) /* No flare color */
+ Fxyz[0] = Wxyz[0], Fxyz[1] = Wxyz[1], Fxyz[2] = Wxyz[2];
+ break;
+ }
+
+ /* This is inputing either a representation of */
+ /* a natural scene captured on another a photo medium, or */
+ /* are is assuming that the medium is the original. */
+ /* We will assume a compromise media original, natural scene */
+ case icSigFilmScanner:
+ {
+ if (La < 0.0) /* No adapting luminance */
+ La = 50.0; /* Use bright indoors, dull outdoors */
+ if (Yb < 0.0) /* No background relative luminance */
+ Yb = 0.2; /* Assume grey world */
+ if (Lv < 0.0) /* No device image luminance */
+ Ev = vc_average; /* Assume average viewing conditions */
+ if (Yf < 0.0) /* No flare figure */
+ Yf = 0.005; /* Assume 0.5% flare */
+ if (Fxyz[0] < 0.0 || Fxyz[1] < 0.0 || Fxyz[2] < 0.0) /* No flare color */
+ Fxyz[0] = Wxyz[0], Fxyz[1] = Wxyz[1], Fxyz[2] = Wxyz[2];
+ break;
+ }
+
+ /* Direct scene to value devices. */
+ case icSigDigitalCamera:
+ case icSigVideoCamera:
+ {
+ if (La < 0.0) /* No adapting luminance */
+ La = 110.0; /* Use very bright indoors, usual outdoors */
+ if (Yb < 0.0) /* No background relative luminance */
+ Yb = 0.2; /* Assume grey world */
+ if (Lv < 0.0) /* No device image luminance */
+ Ev = vc_average; /* Assume average viewing conditions */
+ if (Yf < 0.0) /* No flare figure */
+ Yf = 0.0; /* Assume 0% flare */
+ if (Fxyz[0] < 0.0 || Fxyz[1] < 0.0 || Fxyz[2] < 0.0) /* No flare color */
+ Fxyz[0] = Wxyz[0], Fxyz[1] = Wxyz[1], Fxyz[2] = Wxyz[2];
+ break;
+ }
+
+ /* Emmisive displays. */
+ /* Assume a video monitor is in a darker environment than a CRT */
+ case icSigVideoMonitor:
+ {
+ if (La < 0.0) /* No adapting luminance */
+ La = 4.0; /* Darkened work environment */
+ if (Yb < 0.0) /* No background relative luminance */
+ Yb = 0.2; /* Assume grey world */
+ if (Lv < 0.0) /* No device image luminance */
+ Ev = vc_dim; /* Assume dim viewing conditions */
+ if (Yf < 0.0) /* No flare figure */
+ Yf = 0.01; /* Assume 1% flare */
+ if (Fxyz[0] < 0.0 || Fxyz[1] < 0.0 || Fxyz[2] < 0.0) /* No flare color */
+ Fxyz[0] = Wxyz[0], Fxyz[1] = Wxyz[1], Fxyz[2] = Wxyz[2];
+ break;
+ }
+
+
+ /* Assume a typical work environment */
+ case icSigCRTDisplay:
+ case icSigPMDisplay:
+ case icSigAMDisplay:
+ {
+ if (La < 0.0) /* No adapting luminance */
+ La = 33.0; /* Typical work environment */
+ if (Yb < 0.0) /* No background relative luminance */
+ Yb = 0.2; /* Assume grey world */
+ if (Lv < 0.0) /* No device image luminance */
+ Ev = vc_average; /* Assume average viewing conditions */
+ if (Yf < 0.0) /* No flare figure */
+ Yf = 0.02; /* Assume 2% flare */
+ if (Fxyz[0] < 0.0 || Fxyz[1] < 0.0 || Fxyz[2] < 0.0) /* No flare color */
+ Fxyz[0] = Wxyz[0], Fxyz[1] = Wxyz[1], Fxyz[2] = Wxyz[2];
+ break;
+ }
+
+ /* Photo CD has its own viewing definitions */
+ /* (It represents original scene colors) */
+ case icSigPhotoCD:
+ {
+ if (La < 0.0) /* No adapting luminance */
+ La = 320.0; /* Bright outdoors */
+ if (Yb < 0.0) /* No background relative luminance */
+ Yb = 0.2; /* Assume grey world */
+ if (Lv < 0.0) /* No device image luminance */
+ Ev = vc_average; /* Assume average viewing conditions */
+ if (Yf < 0.0) /* No flare figure */
+ Yf = 0.00; /* Assume 0% flare */
+ if (Fxyz[0] < 0.0 || Fxyz[1] < 0.0 || Fxyz[2] < 0.0) /* No flare color */
+ Fxyz[0] = Wxyz[0], Fxyz[1] = Wxyz[1], Fxyz[2] = Wxyz[2];
+ break;
+ }
+
+ /* Projection devices, either direct, or */
+ /* via another intermediate medium. */
+ /* Assume darkened room, little background */
+ case icSigProjectionTelevision:
+ {
+ if (La < 0.0) /* No adapting luminance */
+ La = 7.0; /* Dark environment */
+ if (Yb < 0.0) /* No background relative luminance */
+ Yb = 0.1; /* Assume little background */
+ if (Lv < 0.0) /* No device image luminance */
+ Ev = vc_dim; /* Dim environment */
+ if (Yf < 0.0) /* No flare figure */
+ Yf = 0.01; /* Assume 1% flare */
+ if (Fxyz[0] < 0.0 || Fxyz[1] < 0.0 || Fxyz[2] < 0.0) /* No flare color */
+ Fxyz[0] = Wxyz[0], Fxyz[1] = Wxyz[1], Fxyz[2] = Wxyz[2];
+ break;
+ }
+ /* Assume very darkened room, no background */
+ case icSigFilmWriter:
+ {
+ if (La < 0.0) /* No adapting luminance */
+ La = 7.0; /* Dark environment */
+ if (Yb < 0.0) /* No background relative luminance */
+ Yb = 0.0; /* Assume no background */
+ if (Lv < 0.0) /* No device image luminance */
+ Ev = vc_dark; /* Dark environment */
+ if (Yf < 0.0) /* No flare figure */
+ Yf = 0.01; /* Assume 1% flare */
+ if (Fxyz[0] < 0.0 || Fxyz[1] < 0.0 || Fxyz[2] < 0.0) /* No flare color */
+ Fxyz[0] = Wxyz[0], Fxyz[1] = Wxyz[1], Fxyz[2] = Wxyz[2];
+ break;
+ }
+
+ /* Printed media devices. */
+ /* Assume a normal print viewing environment */
+ case icSigInkJetPrinter:
+ case icSigThermalWaxPrinter:
+ case icSigElectrophotographicPrinter:
+ case icSigElectrostaticPrinter:
+ case icSigDyeSublimationPrinter:
+ case icSigPhotographicPaperPrinter:
+ case icSigPhotoImageSetter:
+ case icSigGravure:
+ case icSigOffsetLithography:
+ case icSigSilkscreen:
+ case icSigFlexography:
+ {
+ if (La < 0.0) /* No adapting luminance */
+ La = 40.0; /* Use a practical print evaluation number */
+ if (Yb < 0.0) /* No background relative luminance */
+ Yb = 0.2; /* Assume grey world */
+ if (Lv < 0.0) /* No device image luminance */
+ Ev = vc_average; /* Assume average viewing conditions */
+ if (Yf < 0.0) /* No flare figure */
+ Yf = 0.01; /* Assume 1% flare */
+ if (Fxyz[0] < 0.0 || Fxyz[1] < 0.0 || Fxyz[2] < 0.0) /* No flare color */
+ Fxyz[0] = Wxyz[0], Fxyz[1] = Wxyz[1], Fxyz[2] = Wxyz[2];
+ break;
+ }
+
+ default:
+ {
+ return 2;
+ }
+ }
+
+ return 1;
+}
+
+/* Write our viewing conditions to the underlying ICC profile, */
+/* using a private tag. */
+void xicc_set_viewcond(
+xicc *p, /* Expanded profile we're working with */
+icxViewCond *vc /* Viewing parameters to return */
+) {
+ //icc *pp = p->pp; /* Base ICC */
+
+ // ~~1 Not implemented yet
+}
+
+
+
+/* Return an enumerated viewing condition */
+/* Return enumeration if OK, -999 if there is no such enumeration. */
+/* xicc may be NULL if just the description is wanted, */
+/* or an explicit white point is provided. */
+int xicc_enum_viewcond(
+xicc *p, /* Expanded profile to get white point (May be NULL if desc NZ) */
+icxViewCond *vc, /* Viewing parameters to return, May be NULL if desc is nz */
+int no, /* Enumeration to return, -1 for default, -2 for none */
+char *as, /* String alias to number, NULL if none */
+int desc, /* NZ - Just return a description of this enumeration in vc */
+double *wp /* Provide white point if xicc is NULL */
+) {
+
+ if (desc == 0) { /* We're setting the viewing condition */
+ icc *pp; /* Base ICC */
+ icmXYZArray *whitePointTag;
+
+ if (vc == NULL)
+ return -999;
+
+ if (p == NULL) {
+ if (wp == NULL)
+ return -999;
+ vc->Wxyz[0] = wp[0];
+ vc->Wxyz[1] = wp[1];
+ vc->Wxyz[2] = wp[2];
+ } else {
+
+ pp = p->pp;
+ if ((whitePointTag = (icmXYZArray *)pp->read_tag(pp, icSigMediaWhitePointTag)) != NULL
+ && whitePointTag->ttype == icSigXYZType && whitePointTag->size >= 1) {
+ vc->Wxyz[0] = whitePointTag->data[0].X;
+ vc->Wxyz[1] = whitePointTag->data[0].Y;
+ vc->Wxyz[2] = whitePointTag->data[0].Z;
+ } else {
+ if (wp == NULL) {
+ sprintf(p->err,"Enum VC: Failed to read Media White point");
+ p->errc = 2;
+ return -999;
+ }
+ vc->Wxyz[0] = wp[0];
+ vc->Wxyz[1] = wp[1];
+ vc->Wxyz[2] = wp[2];
+ }
+ }
+
+ /* Set a default flare color */
+ vc->Fxyz[0] = vc->Wxyz[0];
+ vc->Fxyz[1] = vc->Wxyz[1];
+ vc->Fxyz[2] = vc->Wxyz[2];
+ }
+
+ /*
+
+ Typical adapting field luminances and white luminance in reflective setup:
+
+ E = illuminance in Lux
+ Lv = White luminance assuming 100% reflectance
+ La = Adapting field luminance in cd/m^2, assuming 20% reflectance from surround
+
+ E La Lv Condition
+ 11 0.7 4 Twilight
+ 32 2 10 Subdued indoor lighting
+ 64 4 20 Less than typical office light; sometimes recommended for
+ display-only workplaces (sRGB)
+ 350 22 111 Typical Office (sRGB annex D)
+ 500 32 160 Practical print evaluationa (ISO-3664 P2)
+ 1000 64 318 Good Print evaluation (CIE 116-1995)
+ 1000 64 318 Television Studio lighting
+ 1000 64 318 Overcast Outdoors
+ 2000 127 637 Critical print evaluation (ISO-3664 P1)
+ 10000 637 3183 Typical outdoors, full daylight
+ 50000 3185 15915 Bright summers day
+
+ */
+
+ if (no == -1
+ || (as != NULL && stricmp(as,"d") == 0)) {
+
+ no = -1;
+ if (vc != NULL) {
+ vc->desc = " d - Default Viewing Condition";
+ vc->Ev = vc_average; /* Average viewing conditions */
+ vc->La = 50.0; /* Practical to Good lighting */
+ vc->Lv = 250.0; /* Average viewing conditions ratio */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.01; /* 1% flare */
+ }
+ }
+ else if (no == 0
+ || (as != NULL && stricmp(as,"pp") == 0)) {
+
+ no = 0;
+ if (vc != NULL) {
+ vc->desc = " pp - Practical Reflection Print (ISO-3664 P2)";
+ vc->Ev = vc_average; /* Average viewing conditions */
+ vc->La = 32.0; /* Use a practical print evaluation number */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.01; /* 1% flare */
+ }
+ }
+ else if (no == 1
+ || (as != NULL && stricmp(as,"pe") == 0)) {
+
+ no = 1;
+ if (vc != NULL) {
+ vc->desc = " pe - Print evaluation environment (CIE 116-1995)";
+ vc->Ev = vc_average; /* Average viewing conditions */
+ vc->La = 64.0; /* Good */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.01; /* 1% flare */
+ }
+ }
+ else if (no == 2
+ || (as != NULL && stricmp(as,"pc") == 0)) {
+
+ no = 2;
+ if (vc != NULL) {
+ vc->desc = " pc - Critical print evaluation environment (ISO-3664 P1)";
+ vc->Ev = vc_average; /* Average viewing conditions */
+ vc->La = 127.0; /* Critical */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.01; /* 1% flare */
+ }
+ }
+ else if (no == 3
+ || (as != NULL && stricmp(as,"mt") == 0)) {
+
+ no = 3;
+ if (vc != NULL) {
+ vc->desc = " mt - Monitor in typical work environment";
+ vc->Ev = vc_average; /* Average viewing conditions */
+ vc->La = 22.0; /* Typical work environment */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.02; /* 2% flare */
+ }
+ }
+ else if (no == 4
+ || (as != NULL && stricmp(as,"mb") == 0)) {
+
+ no = 4;
+ if (vc != NULL) {
+ vc->desc = " mb - Bright monitor in bright work environment";
+ vc->Ev = vc_average; /* Average viewing conditions */
+ vc->La = 42.0; /* Bright work environment */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.02; /* 2% flare */
+ }
+ }
+ else if (no == 5
+ || (as != NULL && stricmp(as,"md") == 0)) {
+
+ no = 5;
+ if (vc != NULL) {
+ vc->desc = " md - Monitor in darkened work environment";
+ vc->Ev = vc_dim; /* Dim viewing conditions */
+ vc->La = 4.0; /* Darkened work environment */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.01; /* 1% flare */
+ }
+ }
+ else if (no == 6
+ || (as != NULL && stricmp(as,"jm") == 0)) {
+
+ no = 6;
+ if (vc != NULL) {
+ vc->desc = " jm - Projector in dim environment";
+ vc->Ev = vc_dim; /* Dim viewing conditions */
+ vc->La = 10.0; /* Adaptation is from display */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.01; /* 1% flare */
+ }
+ }
+ else if (no == 7
+ || (as != NULL && stricmp(as,"jd") == 0)) {
+
+ no = 7;
+ if (vc != NULL) {
+ vc->desc = " jd - Projector in dark environment";
+ vc->Ev = vc_dark; /* Dark viewing conditions */
+ vc->La = 10.0; /* Adaptation is from display */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.01; /* 1% flare ? */
+ }
+ }
+ else if (no == 8
+ || (as != NULL && stricmp(as,"pcd") == 0)) {
+
+ no = 8;
+ if (vc != NULL) {
+ vc->desc = "pcd - Photo CD - original scene outdoors";
+ vc->Ev = vc_average; /* Average viewing conditions */
+ vc->La = 320.0; /* Typical outdoors, 1600 cd/m^2 */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.00; /* 0% flare */
+ }
+ }
+ else if (no == 9
+ || (as != NULL && stricmp(as,"ob") == 0)) {
+
+ no = 9;
+ if (vc != NULL) {
+ vc->desc = " ob - Original scene - Bright Outdoors";
+ vc->Ev = vc_average; /* Average viewing conditions */
+ vc->La = 2000.0; /* Bright Outdoors */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.00; /* 0% flare */
+ }
+ }
+ else if (no == 10
+ || (as != NULL && stricmp(as,"cx") == 0)) {
+
+ no = 10;
+ if (vc != NULL) {
+ vc->desc = " cx - Cut Sheet Transparencies on a viewing box";
+ vc->Ev = vc_cut_sheet; /* Cut sheet viewing conditions */
+ vc->La = 53.0; /* Dim, adapted to slide ? */
+ vc->Yb = 0.2; /* Grey world */
+ vc->Yf = 0.01; /* 1% flare ? */
+ }
+ }
+ else {
+ if (p != NULL) {
+ sprintf(p->err,"Enum VC: Unrecognised enumeration %d",no);
+ p->errc = 1;
+ }
+ return -999;
+ }
+
+ return no;
+}
+
+/* Debug: dump a Viewing Condition to standard out */
+void xicc_dump_viewcond(
+icxViewCond *vc
+) {
+ printf("Viewing Condition:\n");
+ if (vc->Ev == vc_dark)
+ printf(" Surround to Image: Dark\n");
+ else if (vc->Ev == vc_dim)
+ printf(" Surround to Image: Dim\n");
+ else if (vc->Ev == vc_average)
+ printf(" Surround to Image: Average\n");
+ else if (vc->Ev == vc_cut_sheet)
+ printf(" Transparency on Light box\n");
+ printf(" Adapted white = %f %f %f\n",vc->Wxyz[0], vc->Wxyz[1], vc->Wxyz[2]);
+ printf(" Adapted luminance = %f cd/m^2\n",vc->La);
+ printf(" Background to image ratio = %f\n",vc->Yb);
+ if (vc->Ev == vc_none)
+ printf(" Image luminance = %f cd/m^2\n",vc->Lv);
+ printf(" Flare to image ratio = %f\n",vc->Yf);
+ printf(" Flare color = %f %f %f\n",vc->Fxyz[0], vc->Fxyz[1], vc->Fxyz[2]);
+}
+
+
+/* Debug: dump an Inking setup to standard out */
+void xicc_dump_inking(icxInk *ik) {
+ printf("Inking settings:\n");
+ if (ik->tlimit < 0.0)
+ printf("No total limit\n");
+ else
+ printf("Total limit = %f%%\n",ik->tlimit * 100.0);
+
+ if (ik->klimit < 0.0)
+ printf("No black limit\n");
+ else
+ printf("Black limit = %f%%\n",ik->klimit * 100.0);
+
+ if (ik->KonlyLmin)
+ printf("K only black as locus Lmin\n");
+ else
+ printf("Normal black as locus Lmin\n");
+
+ if (ik->k_rule == icxKvalue) {
+ printf("Inking rule is a fixed K target\n");
+ } if (ik->k_rule == icxKlocus) {
+ printf("Inking rule is a fixed locus target\n");
+ } if (ik->k_rule == icxKluma5 || ik->k_rule == icxKluma5k) {
+ if (ik->k_rule == icxKluma5)
+ printf("Inking rule is a 5 parameter locus function of L\n");
+ else
+ printf("Inking rule is a 5 parameter K function of L\n");
+ printf("Ksmth = %f\n",ik->c.Ksmth);
+ printf("Kskew = %f\n",ik->c.Kskew);
+ printf("Kstle = %f\n",ik->c.Kstle);
+ printf("Kstpo = %f\n",ik->c.Kstpo);
+ printf("Kenpo = %f\n",ik->c.Kenpo);
+ printf("Kenle = %f\n",ik->c.Kenle);
+ printf("Kshap = %f\n",ik->c.Kshap);
+ } if (ik->k_rule == icxKl5l || ik->k_rule == icxKl5lk) {
+ if (ik->k_rule == icxKl5l)
+ printf("Inking rule is a 2x5 parameter locus function of L and K aux\n");
+ else
+ printf("Inking rule is a 2x5 parameter K function of L and K aux\n");
+ printf("Min Ksmth = %f\n",ik->c.Ksmth);
+ printf("Min Kskew = %f\n",ik->c.Kskew);
+ printf("Min Kstle = %f\n",ik->c.Kstle);
+ printf("Min Kstpo = %f\n",ik->c.Kstpo);
+ printf("Min Kenpo = %f\n",ik->c.Kenpo);
+ printf("Min Kenle = %f\n",ik->c.Kenle);
+ printf("Min Kshap = %f\n",ik->c.Kshap);
+ printf("Max Ksmth = %f\n",ik->x.Ksmth);
+ printf("Max Kskew = %f\n",ik->x.Kskew);
+ printf("Max Kstle = %f\n",ik->x.Kstle);
+ printf("Max Kstpo = %f\n",ik->x.Kstpo);
+ printf("Max Kenpo = %f\n",ik->x.Kenpo);
+ printf("Max Kenle = %f\n",ik->x.Kenle);
+ printf("Max Kshap = %f\n",ik->x.Kshap);
+ }
+}
+
+/* ------------------------------------------------------ */
+/* Gamut Mapping Intent stuff */
+
+/* Return an enumerated gamut mapping intent */
+/* Return enumeration if OK, icxIllegalGMIntent if there is no such enumeration. */
+int xicc_enum_gmapintent(
+icxGMappingIntent *gmi, /* Gamut Mapping parameters to return */
+int no, /* Enumeration selected, icxNoGMIntent for none */
+char *as /* Alias string selector, NULL for none */
+) {
+#ifdef USE_CAM
+ int colccas = 0x2; /* Use cas clipping for colorimetric style intents */
+ int perccas = 0x1; /* Use cas for perceptual style intents */
+#else
+ int colccas = 0x0; /* Use Lab for colorimetric style intents */
+ int perccas = 0x0; /* Use Lab for perceptual style intents */
+ fprintf(stderr,"!!!!!! Warning, USE_CAM is off in xicc.c !!!!!!\n");
+#endif
+
+ /* Assert default if no guidance given */
+ if (no == icxNoGMIntent && as == NULL)
+ no = icxDefaultGMIntent;
+
+ if (no == 0
+ || no == icxAbsoluteGMIntent
+ || (as != NULL && stricmp(as,"a") == 0)) {
+ /* Map Absolute appearance space Jab to Jab and clip out of gamut */
+ no = 0;
+ gmi->as = "a";
+ gmi->desc = " a - Absolute Colorimetric (in Jab) [ICC Absolute Colorimetric]";
+ gmi->icci = icAbsoluteColorimetric;
+ gmi->usecas = colccas; /* Use absolute appearance space */
+ gmi->usemap = 0; /* Don't use gamut mapping */
+ gmi->greymf = 0.0;
+ gmi->glumwcpf = 0.0;
+ gmi->glumwexf = 0.0;
+ gmi->glumbcpf = 0.0;
+ gmi->glumbexf = 0.0;
+ gmi->glumknf = 0.0;
+ gmi->gamcpf = 0.0;
+ gmi->gamexf = 0.0;
+ gmi->gamcknf = 0.0;
+ gmi->gamxknf = 0.0;
+ gmi->gampwf = 0.0;
+ gmi->gamswf = 0.0;
+ gmi->satenh = 0.0; /* No saturation enhancement */
+ }
+ else if (no == 1
+ || (as != NULL && stricmp(as,"aw") == 0)) {
+
+ /* I'm not sure how often this intent is useful. It's less likely than */
+ /* I though that a printer white point won't fit within the gamut */
+ /* of a display profile, since the display white always has Y = 1.0, */
+ /* and no paper has better than about 95% reflectance. */
+ /* Perhaps it may be more useful for targeting printer profiles ? */
+
+ /* Map Absolute Jab to Jab and scale source to avoid clipping the white point */
+ no = 1;
+ gmi->as = "aw";
+ gmi->desc = "aw - Absolute Colorimetric (in Jab) with scaling to fit white point";
+ gmi->icci = icAbsoluteColorimetric;
+ gmi->usecas = 0x100 | colccas; /* Absolute Appearance space with scaling */
+ /* to avoid clipping the source white point */
+ gmi->usemap = 0; /* Don't use gamut mapping */
+ gmi->greymf = 0.0;
+ gmi->glumwcpf = 0.0;
+ gmi->glumwexf = 0.0;
+ gmi->glumbcpf = 0.0;
+ gmi->glumbexf = 0.0;
+ gmi->glumknf = 0.0;
+ gmi->gamcpf = 0.0;
+ gmi->gamexf = 0.0;
+ gmi->gamcknf = 0.0;
+ gmi->gamxknf = 0.0;
+ gmi->gampwf = 0.0;
+ gmi->gamswf = 0.0;
+ gmi->satenh = 0.0; /* No saturation enhancement */
+ }
+ else if (no == 2
+ || (as != NULL && stricmp(as,"aa") == 0)) {
+
+ /* Map appearance space Jab to Jab and clip out of gamut */
+ no = 2;
+ gmi->as = "aa";
+ gmi->desc = "aa - Absolute Appearance";
+ gmi->icci = icRelativeColorimetric;
+ gmi->usecas = perccas; /* Appearance space */
+ gmi->usemap = 0; /* Don't use gamut mapping */
+ gmi->greymf = 0.0;
+ gmi->glumwcpf = 0.0;
+ gmi->glumwexf = 0.0;
+ gmi->glumbcpf = 0.0;
+ gmi->glumbexf = 0.0;
+ gmi->glumknf = 0.0;
+ gmi->gamcpf = 0.0;
+ gmi->gamexf = 0.0;
+ gmi->gamcknf = 0.0;
+ gmi->gamxknf = 0.0;
+ gmi->gampwf = 0.0;
+ gmi->gamswf = 0.0;
+ gmi->satenh = 0.0; /* No saturation enhancement */
+ }
+ else if (no == 3
+ || no == icxRelativeGMIntent
+ || (as != NULL && stricmp(as,"r") == 0)) {
+
+ /* Align neutral axes and linearly map white point, then */
+ /* map appearance space Jab to Jab and clip out of gamut */
+ no = 3;
+ gmi->as = "r";
+ gmi->desc = " r - White Point Matched Appearance [ICC Relative Colorimetric]";
+ gmi->icci = icRelativeColorimetric;
+ gmi->usecas = perccas; /* Appearance space */
+ gmi->usemap = 1; /* Use gamut mapping */
+ gmi->greymf = 1.0; /* Fully align grey axis */
+ gmi->glumwcpf = 1.0; /* Fully compress grey axis at white end */
+ gmi->glumwexf = 1.0; /* Fully expand grey axis at white end */
+ gmi->glumbcpf = 0.0; /* No compression at black end */
+ gmi->glumbexf = 0.0; /* No expansion at black end */
+ gmi->glumknf = 0.0;
+ gmi->gamcpf = 0.0;
+ gmi->gamexf = 0.0;
+ gmi->gamcknf = 0.0;
+ gmi->gamxknf = 0.0;
+ gmi->gampwf = 0.0;
+ gmi->gamswf = 0.0;
+ gmi->satenh = 0.0; /* No saturation enhancement */
+ }
+ else if (no == 4
+ || (as != NULL && stricmp(as,"la") == 0)) {
+
+ /* Align neutral axes and linearly map white and black points, then */
+ /* map appearance space Jab to Jab and clip out of gamut */
+ no = 4;
+ gmi->as = "la";
+ gmi->desc = "la - Luminance axis matched Appearance";
+ gmi->icci = icRelativeColorimetric;
+ gmi->usecas = perccas; /* Appearance space */
+ gmi->usemap = 1; /* Use gamut mapping */
+ gmi->greymf = 1.0; /* Fully align grey axis */
+ gmi->glumwcpf = 1.0; /* Fully compress grey axis at white end */
+ gmi->glumwexf = 1.0; /* Fully expand grey axis at white end */
+ gmi->glumbcpf = 1.0; /* Fully compress grey axis at black end */
+ gmi->glumbexf = 1.0; /* Fully expand grey axis at black end */
+ gmi->glumknf = 0.0; /* No knee on grey mapping */
+ gmi->gamcpf = 0.0; /* No gamut compression */
+ gmi->gamexf = 0.0; /* No gamut expansion */
+ gmi->gamcknf = 0.0; /* No knee in gamut compress */
+ gmi->gamxknf = 0.0; /* No knee in gamut expand */
+ gmi->gampwf = 0.0; /* No Perceptual surface weighting factor */
+ gmi->gamswf = 0.0; /* No Saturation surface weighting factor */
+ gmi->satenh = 0.0; /* No saturation enhancement */
+ }
+ else if (no == 5
+ || no == icxDefaultGMIntent
+ || no == icxPerceptualGMIntent
+ || (as != NULL && stricmp(as,"p") == 0)) {
+
+ /* Align neutral axes and perceptually map white and black points, */
+ /* perceptually compress out of gamut and map appearance space Jab to Jab. */
+ no = 5;
+ gmi->as = "p";
+ gmi->desc = " p - Perceptual (Preferred) (Default) [ICC Perceptual]";
+ gmi->icci = icPerceptual;
+ gmi->usecas = perccas; /* Appearance space */
+ gmi->usemap = 1; /* Use gamut mapping */
+ gmi->greymf = 1.0; /* Fully align grey axis */
+ gmi->glumwcpf = 1.0; /* Fully compress grey axis at white end */
+ gmi->glumwexf = 1.0; /* Fully expand grey axis at white end */
+ gmi->glumbcpf = 1.0; /* Fully compress grey axis at black end */
+ gmi->glumbexf = 1.0; /* Fully expand grey axis at black end */
+ gmi->glumknf = 1.0; /* Sigma knee in grey compress/expand */
+ gmi->gamcpf = 1.0; /* Full gamut compression */
+ gmi->gamexf = 0.0; /* No gamut expansion */
+ gmi->gamcknf = 0.8; /* High Sigma knee in gamut compress */
+ gmi->gamxknf = 0.0; /* No knee in gamut expand */
+ gmi->gampwf = 1.0; /* Full Perceptual surface weighting factor */
+ gmi->gamswf = 0.0; /* No Saturation surface weighting factor */
+ gmi->satenh = 0.0; /* No saturation enhancement */
+ }
+ else if (no == 6
+ || (as != NULL && stricmp(as,"pa") == 0)) {
+
+ /* Don't align neutral axes, but perceptually compress out of gamut */
+ /* and map appearance space Jab to Jab. */
+ no = 5;
+ gmi->as = "pa";
+ gmi->desc = "pa - Perceptual Apperance ";
+ gmi->icci = icPerceptual;
+ gmi->usecas = perccas; /* Appearance space */
+ gmi->usemap = 1; /* Use gamut mapping */
+ gmi->greymf = 0.0; /* Don't align grey axis */
+ gmi->glumwcpf = 1.0; /* Fully compress grey axis at white end */
+ gmi->glumwexf = 1.0; /* Fully expand grey axis at white end */
+ gmi->glumbcpf = 1.0; /* Fully compress grey axis at black end */
+ gmi->glumbexf = 1.0; /* Fully expand grey axis at black end */
+ gmi->glumknf = 1.0; /* Sigma knee in grey compress/expand */
+ gmi->gamcpf = 1.0; /* Full gamut compression */
+ gmi->gamexf = 0.0; /* No gamut expansion */
+ gmi->gamcknf = 0.8; /* High Sigma knee in gamut compress */
+ gmi->gamxknf = 0.0; /* No knee in gamut expand */
+ gmi->gampwf = 1.0; /* Full Perceptual surface weighting factor */
+ gmi->gamswf = 0.0; /* No Saturation surface weighting factor */
+ gmi->satenh = 0.0; /* No saturation enhancement */
+ }
+ else if (no == 7
+ || (as != NULL && stricmp(as,"ms") == 0)) {
+
+ /* Align neutral axes and perceptually map white and black points, */
+ /* perceptually compress and expand to match gamuts and map Jab to Jab. */
+ no = 6;
+ gmi->as = "ms";
+ gmi->desc = "ms - Saturation";
+ gmi->icci = icSaturation;
+ gmi->usecas = perccas; /* Appearance space */
+ gmi->usemap = 1; /* Use gamut mapping */
+ gmi->greymf = 1.0; /* Fully align grey axis */
+ gmi->glumwcpf = 1.0; /* Fully compress grey axis at white end */
+ gmi->glumwexf = 1.0; /* Fully expand grey axis at white end */
+ gmi->glumbcpf = 1.0; /* Fully compress grey axis at black end */
+ gmi->glumbexf = 1.0; /* Fully expand grey axis at black end */
+ gmi->glumknf = 1.0; /* Sigma knee in grey compress/expand */
+ gmi->gamcpf = 1.0; /* Full gamut compression */
+ gmi->gamexf = 1.0; /* Full gamut expansion */
+ gmi->gamcknf = 1.0; /* High Sigma knee in gamut compress/expand */
+ gmi->gamxknf = 0.4; /* Moderate Sigma knee in gamut compress/expand */
+ gmi->gampwf = 0.2; /* Slight perceptual surface weighting factor */
+ gmi->gamswf = 0.8; /* Most saturation surface weighting factor */
+ gmi->satenh = 0.0; /* No saturation enhancement */
+ }
+ else if (no == 8
+ || no == icxSaturationGMIntent
+ || (as != NULL && stricmp(as,"s") == 0)) {
+
+ /* Same as "ms" but enhance saturation */
+ no = 7;
+ gmi->as = "s";
+ gmi->desc = " s - Enhanced Saturation [ICC Saturation]";
+ gmi->icci = icSaturation;
+ gmi->usecas = perccas; /* Appearance space */
+ gmi->usemap = 1; /* Use gamut mapping */
+ gmi->greymf = 1.0; /* Fully align grey axis */
+ gmi->glumwcpf = 1.0; /* Fully compress grey axis at white end */
+ gmi->glumwexf = 1.0; /* Fully expand grey axis at white end */
+ gmi->glumbcpf = 1.0; /* Fully compress grey axis at black end */
+ gmi->glumbexf = 1.0; /* Fully expand grey axis at black end */
+ gmi->glumknf = 1.0; /* Sigma knee in grey compress/expand */
+ gmi->gamcpf = 1.0; /* Full gamut compression */
+ gmi->gamexf = 1.0; /* Full gamut expansion */
+ gmi->gamcknf = 1.0; /* High sigma knee in gamut compress */
+ gmi->gamxknf = 0.5; /* Moderate sigma knee in gamut expand */
+ gmi->gampwf = 0.0; /* No Perceptual surface weighting factor */
+ gmi->gamswf = 1.0; /* Full Saturation surface weighting factor */
+ gmi->satenh = 0.9; /* Medium saturation enhancement */
+ }
+ else if (no == 9
+ || (as != NULL && stricmp(as,"al") == 0)) {
+
+ /* Map absolute L*a*b* to L*a*b* and clip out of gamut */
+ no = 8;
+ gmi->as = "al";
+ gmi->desc = "al - Absolute Colorimetric (Lab)";
+ gmi->icci = icAbsoluteColorimetric;
+ gmi->usecas = 0x0; /* Don't use appearance space, use L*a*b* */
+ gmi->usemap = 0; /* Don't use gamut mapping */
+ gmi->greymf = 0.0;
+ gmi->glumwcpf = 0.0;
+ gmi->glumwexf = 0.0;
+ gmi->glumbcpf = 0.0;
+ gmi->glumbexf = 0.0;
+ gmi->glumknf = 0.0;
+ gmi->gamcpf = 0.0;
+ gmi->gamexf = 0.0;
+ gmi->gamcknf = 0.0;
+ gmi->gamxknf = 0.0;
+ gmi->gampwf = 0.0;
+ gmi->gamswf = 0.0;
+ gmi->satenh = 0.0; /* No saturation enhancement */
+ }
+ else if (no == 10
+ || (as != NULL && stricmp(as,"rl") == 0)) {
+
+ /* Align neutral axes and linearly map white point, then */
+ /* map L*a*b* to L*a*b* and clip out of gamut */
+ no = 3;
+ gmi->as = "rl";
+ gmi->desc = "rl - White Point Matched Appearance (Lab)";
+ gmi->icci = icRelativeColorimetric;
+ gmi->usecas = 0x0; /* Don't use appearance space, use L*a*b* */
+ gmi->usemap = 1; /* Use gamut mapping */
+ gmi->greymf = 1.0; /* And linearly map white point */
+ gmi->glumwcpf = 1.0;
+ gmi->glumwexf = 1.0;
+ gmi->glumbcpf = 0.0;
+ gmi->glumbexf = 0.0;
+ gmi->glumknf = 0.0;
+ gmi->gamcpf = 0.0;
+ gmi->gamexf = 0.0;
+ gmi->gamcknf = 0.0;
+ gmi->gamxknf = 0.0;
+ gmi->gampwf = 0.0;
+ gmi->gamswf = 0.0;
+ gmi->satenh = 0.0; /* No saturation enhancement */
+ }
+ else { /* icxIllegalGMIntent */
+ return icxIllegalGMIntent;
+ }
+
+ return no;
+}
+
+
+/* Debug: dump a Gamut Mapping specification */
+void xicc_dump_gmi(
+icxGMappingIntent *gmi /* Gamut Mapping parameters to return */
+) {
+ printf(" Gamut Mapping Specification:\n");
+ if (gmi->desc != NULL)
+ printf(" Description = '%s'\n",gmi->desc);
+ printf(" Closest ICC intent = '%s'\n",icm2str(icmRenderingIntent,gmi->icci));
+
+ if ((gmi->usecas & 0xff) == 0)
+ printf(" Not using Color Apperance Space\n");
+ else if ((gmi->usecas & 0xff) == 1)
+ printf(" Using Color Apperance Space\n");
+ else if ((gmi->usecas & 0xff) == 2)
+ printf(" Using Absolute Color Apperance Space\n");
+
+ if ((gmi->usecas & 0x100) != 0)
+ printf(" Scaling source to avoid white point clipping\n");
+
+ if (gmi->usemap == 0)
+ printf(" Not using Mapping\n");
+ else {
+ printf(" Using Mapping with parameters:\n");
+ printf(" Grey axis alignment factor %f\n", gmi->greymf);
+ printf(" Grey axis white compression factor %f\n", gmi->glumwcpf);
+ printf(" Grey axis white expansion factor %f\n", gmi->glumwexf);
+ printf(" Grey axis black compression factor %f\n", gmi->glumbcpf);
+ printf(" Grey axis black expansion factor %f\n", gmi->glumbexf);
+ printf(" Grey axis knee factor %f\n", gmi->glumknf);
+ printf(" Gamut compression factor %f\n", gmi->gamcpf);
+ printf(" Gamut expansion factor %f\n", gmi->gamexf);
+ printf(" Gamut compression knee factor %f\n", gmi->gamcknf);
+ printf(" Gamut expansion knee factor %f\n", gmi->gamxknf);
+ printf(" Gamut Perceptual mapping weighting factor %f\n", gmi->gampwf);
+ printf(" Gamut Saturation mapping weighting factor %f\n", gmi->gamswf);
+ printf(" Saturation enhancement factor %f\n", gmi->satenh);
+ }
+}
+
+/* ------------------------------------------------------ */
+/* Turn xicc xcal into limit calibration callback */
+
+/* Given an icc profile, try and create an xcal */
+/* Return NULL on error or no cal */
+xcal *xiccReadCalTag(icc *p) {
+ xcal *cal = NULL;
+ icTagSignature sig = icmMakeTag('t','a','r','g');
+ icmText *ro;
+ int oi, tab;
+
+//printf("~1 about to look for CAL in profile\n");
+ if ((ro = (icmText *)p->read_tag(p, sig)) != NULL) {
+ cgatsFile *cgf;
+ cgats *icg;
+
+ if (ro->ttype != icSigTextType)
+ return NULL;
+
+//printf("~1 found 'targ' tag\n");
+ if ((icg = new_cgats()) == NULL) {
+ return NULL;
+ }
+ if ((cgf = new_cgatsFileMem(ro->data, ro->size)) != NULL) {
+ icg->add_other(icg, "CTI3");
+ oi = icg->add_other(icg, "CAL");
+
+//printf("~1 created cgats object from 'targ' tag\n");
+ if (icg->read(icg, cgf) == 0) {
+
+ for (tab = 0; tab < icg->ntables; tab++) {
+ if (icg->t[tab].tt == tt_other && icg->t[tab].oi == oi) {
+ break;
+ }
+ }
+ if (tab < icg->ntables) {
+//printf("~1 found CAL table\n");
+
+ if ((cal = new_xcal()) == NULL) {
+ icg->del(icg);
+ cgf->del(cgf);
+ return NULL;
+ }
+ if (cal->read_cgats(cal, icg, tab, "'targ' tag") != 0) {
+#ifdef DEBUG
+ printf("read_cgats on cal tag failed\n");
+#endif
+ cal->del(cal);
+ cal = NULL;
+ }
+//else printf("~1 read CAL and creaded xcal object OK\n");
+ }
+ }
+ cgf->del(cgf);
+ }
+ icg->del(icg);
+ }
+ return cal;
+}
+
+/* A callback that uses an xcal, that can be used with icc get_tac */
+void xiccCalCallback(void *cntx, double *out, double *in) {
+ xcal *cal = (xcal *)cntx;
+
+ cal->interp(cal, out, in);
+}
+
+/* ---------------------------------------------- */
+
+/* Utility function - given an open icc profile, */
+/* guess which channel is the black. */
+/* Return -1 if there is no black channel or it can't be guessed */
+int icxGuessBlackChan(icc *p) {
+ int kch = -1;
+
+ switch (p->header->colorSpace) {
+ case icSigCmykData:
+ kch = 3;
+ break;
+
+ /* Use a heuristic to detect the black channel. */
+ /* This duplicates code in icxLu_comp_bk_point() :-( */
+ /* Colorant guessing should go in icclib ? */
+ case icSig2colorData:
+ case icSig3colorData:
+ case icSig4colorData:
+ case icSig5colorData:
+ case icSig6colorData:
+ case icSig7colorData:
+ case icSig8colorData:
+ case icSig9colorData:
+ case icSig10colorData:
+ case icSig11colorData:
+ case icSig12colorData:
+ case icSig13colorData:
+ case icSig14colorData:
+ case icSig15colorData:
+ case icSigMch5Data:
+ case icSigMch6Data:
+ case icSigMch7Data:
+ case icSigMch8Data: {
+ icmLuBase *lu;
+ double dval[MAX_CHAN];
+ double ncval[3];
+ double cvals[MAX_CHAN][3];
+ int inn, e, nlighter, ndarker;
+
+ /* Grab a lookup object */
+ if ((lu = p->get_luobj(p, icmFwd, icRelativeColorimetric, icSigLabData, icmLuOrdNorm)) == NULL)
+ error("icxGetLimits: assert: getting Fwd Lookup failed!");
+
+ lu->spaces(lu, NULL, &inn, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
+
+ /* Decide if the colorspace is aditive or subtractive */
+
+ /* First the no colorant value */
+ for (e = 0; e < inn; e++)
+ dval[e] = 0.0;
+ lu->lookup(lu, ncval, dval);
+
+ /* Then all the colorants */
+ nlighter = ndarker = 0;
+ for (e = 0; e < inn; e++) {
+ dval[e] = 1.0;
+ lu->lookup(lu, cvals[e], dval);
+ dval[e] = 0.0;
+ if (fabs(cvals[e][0] - ncval[0]) > 5.0) {
+ if (cvals[e][0] > ncval[0])
+ nlighter++;
+ else
+ ndarker++;
+ }
+ }
+
+ if (ndarker > 0 && nlighter == 0) { /* Assume subtractive. */
+ double pbk[3] = { 0.0,0.0,0.0 }; /* Perfect black */
+ double smd = 1e10; /* Smallest distance */
+
+ /* Guess the black channel */
+ for (e = 0; e < inn; e++) {
+ double tt;
+ tt = icmNorm33sq(pbk, cvals[e]);
+ if (tt < smd) {
+ smd = tt;
+ kch = e;
+ }
+ }
+ /* See if the black seems sane */
+ if (cvals[kch][0] > 40.0
+ || fabs(cvals[kch][1]) > 10.0
+ || fabs(cvals[kch][2]) > 10.0) {
+ kch = -1;
+ }
+ }
+ lu->del(lu);
+ }
+ break;
+
+ default:
+ break;
+ }
+
+ return kch;
+}
+
+/* Utility function - given an open icc profile, */
+/* estmate the total ink limit and black ink limit. */
+/* Note that this is rather rough, because ICC profiles */
+/* don't have a tag for this information, and ICC profiles */
+/* don't have any straightforward way of identifying particular */
+/* color channels for > 4 color. */
+/* If there are no limits, or they are not discoverable or */
+/* applicable, return values of -1.0 */
+
+void icxGetLimits(
+xicc *xp,
+double *tlimit,
+double *klimit
+) {
+ icc *p = xp->pp;
+ int nch;
+ double max[MAX_CHAN]; /* Max of each channel */
+ double total;
+
+ total = p->get_tac(p, max, xp->cal != NULL ? xiccCalCallback : NULL, (void *)xp->cal);
+
+ if (total < 0.0) { /* Not valid */
+ if (tlimit != NULL)
+ *tlimit = -1.0;
+ if (klimit != NULL)
+ *klimit = -1.0;
+ return;
+ }
+
+ nch = icmCSSig2nchan(p->header->colorSpace);
+
+ /* No effective limit */
+ if (tlimit != NULL) {
+ if (total >= (double)nch) {
+ *tlimit = -1.0;
+ } else {
+ *tlimit = total;
+ }
+ }
+
+ if (klimit != NULL) {
+ int kch;
+
+ kch = icxGuessBlackChan(p);
+
+ if (kch < 0 || max[kch] >= 1.0) {
+ *klimit = -1.0;
+ } else {
+ *klimit = max[kch];
+ }
+ }
+}
+
+/* Replace a non-set limit (ie. < 0.0) with the heuristic from */
+/* the given profile. */
+void icxDefaultLimits(
+xicc *xp,
+double *tlout,
+double tlin,
+double *klout,
+double klin
+) {
+ if (tlin < 0.0 || klin < 0.0) {
+ double tl, kl;
+
+ icxGetLimits(xp, &tl, &kl);
+
+ if (tlin < 0.0)
+ tlin = tl;
+
+ if (klin < 0.0)
+ klin = kl;
+ }
+
+ if (tlout != NULL)
+ *tlout = tlin;
+
+ if (klout != NULL)
+ *klout = klin;
+}
+
+/* Structure to hold optimisation information */
+typedef struct {
+ xcal *cal;
+ double ilimit;
+ double uilimit;
+} ulimctx;
+
+/* Callback to find equivalent underlying total limit */
+/* and try and maximize it while remaining within gamut */
+static double ulimitfunc(void *cntx, double pv[]) {
+ ulimctx *cx = (ulimctx *)cntx;
+ xcal *cal = cx->cal;
+ int devchan = cal->devchan;
+ int i;
+ double dv, odv;
+ double og = 0.0, rv = 0.0;
+
+ double usum = 0.0, sum = 0.0;
+
+ /* Comute calibrated sum of channels except last */
+ for (i = 0; i < (devchan-1); i++) {
+ double dv = pv[i]; /* Underlying (pre-calibration) device value */
+ usum += dv; /* Underlying sum */
+ if (dv < 0.0) {
+ og += -dv;
+ dv = 0.0;
+ } else if (dv > 1.0) {
+ og += dv - 1.0;
+ dv = 1.0;
+ } else
+ dv = cal->interp_ch(cal, i, dv); /* Calibrated device value */
+ sum += dv; /* Calibrated device sum */
+ }
+ /* Compute the omitted channel value */
+ dv = cx->ilimit - sum; /* Omitted calibrated device value */
+ if (dv < 0.0) {
+ og += -dv;
+ dv = 0.0;
+ } else if (dv > 1.0) {
+ og += dv - 1.0;
+ dv = 1.0;
+ } else
+ dv = cal->inv_interp_ch(cal, i, dv); /* Omitted underlying device value */
+ usum += dv; /* Underlying sum */
+ cx->uilimit = usum;
+
+ rv = 10000.0 * og - usum; /* Penalize out of gamut, maximize underlying sum */
+
+//printf("~1 returning %f from %f %f %f %f\n",rv,pv[0],pv[1],pv[2],dv);
+ return rv;
+}
+
+/* Given a calibrated total ink limit and an xcal, return the */
+/* equivalent underlying (pre-calibration) total ink limit. */
+/* This is the maximum equivalent, that makes sure that */
+/* the calibrated limit is met or exceeded. */
+double icxMaxUnderlyingLimit(xcal *cal, double ilimit) {
+ ulimctx cx;
+ int i;
+ double dv[MAX_CHAN];
+ double sr[MAX_CHAN];
+ double rv; /* Residual value */
+
+ if (cal->devchan <= 1) {
+ return cal->inv_interp_ch(cal, 0, ilimit);
+ }
+
+ cx.cal = cal;
+ cx.ilimit = ilimit;
+
+ for (i = 0; i < (cal->devchan-1); i++) {
+ sr[i] = 0.05;
+ dv[i] = 0.1;
+ }
+ if (powell(&rv, cal->devchan-1, dv, sr, 0.000001, 1000, ulimitfunc,
+ (void *)&cx, NULL, NULL) != 0) {
+ warning("icxUnderlyingLimit() failed for chan %d, ilimit %f\n",cal->devchan,ilimit);
+ return ilimit;
+ }
+ ulimitfunc((void *)&cx, dv);
+
+ return cx.uilimit;
+}
+
+/* ------------------------------------------------------ */
+/* Conversion and deltaE formular that include partial */
+/* derivatives, for use within fit parameter optimisations. */
+
+/* CIE XYZ to perceptual Lab with partial derivatives. */
+void icxdXYZ2Lab(icmXYZNumber *w, double *out, double dout[3][3], double *in) {
+ double wp[3], tin[3], dtin[3];
+ int i;
+
+ wp[0] = w->X, wp[1] = w->Y, wp[2] = w->Z;
+
+ for (i = 0; i < 3; i++) {
+ tin[i] = in[i]/wp[i];
+ dtin[i] = 1.0/wp[i];
+
+ if (tin[i] > 0.008856451586) {
+ dtin[i] *= pow(tin[i], -2.0/3.0) / 3.0;
+ tin[i] = pow(tin[i],1.0/3.0);
+ } else {
+ dtin[i] *= 7.787036979;
+ tin[i] = 7.787036979 * tin[i] + 16.0/116.0;
+ }
+ }
+
+ out[0] = 116.0 * tin[1] - 16.0;
+ dout[0][0] = 0.0;
+ dout[0][1] = 116.0 * dtin[1];
+ dout[0][2] = 0.0;
+
+ out[1] = 500.0 * (tin[0] - tin[1]);
+ dout[1][0] = 500.0 * dtin[0];
+ dout[1][1] = 500.0 * -dtin[1];
+ dout[1][2] = 0.0;
+
+ out[2] = 200.0 * (tin[1] - tin[2]);
+ dout[2][0] = 0.0;
+ dout[2][1] = 200.0 * dtin[1];
+ dout[2][2] = 200.0 * -dtin[2];
+}
+
+
+/* Return the normal Delta E squared, given two Lab values, */
+/* including partial derivatives. */
+double icxdLabDEsq(double dout[2][3], double *Lab0, double *Lab1) {
+ double rv = 0.0, tt;
+
+ tt = Lab0[0] - Lab1[0];
+ dout[0][0] = 2.0 * tt;
+ dout[1][0] = -2.0 * tt;
+ rv += tt * tt;
+ tt = Lab0[1] - Lab1[1];
+ dout[0][1] = 2.0 * tt;
+ dout[1][1] = -2.0 * tt;
+ rv += tt * tt;
+ tt = Lab0[2] - Lab1[2];
+ dout[0][2] = 2.0 * tt;
+ dout[1][2] = -2.0 * tt;
+ rv += tt * tt;
+ return rv;
+}
+
+/* Return the CIE94 Delta E color difference measure, squared */
+/* including partial derivatives. */
+double icxdCIE94sq(double dout[2][3], double Lab0[3], double Lab1[3]) {
+ double desq, _desq[2][3];
+ double dlsq;
+ double dcsq, _dcsq[2][2]; /* == [x][1,2] */
+ double c12, _c12[2][2]; /* == [x][1,2] */
+ double dhsq, _dhsq[2][2]; /* == [x][1,2] */
+ double rv;
+
+ {
+ double dl, da, db;
+
+ dl = Lab0[0] - Lab1[0];
+ dlsq = dl * dl; /* dl squared */
+ da = Lab0[1] - Lab1[1];
+ db = Lab0[2] - Lab1[2];
+
+
+ /* Compute normal Lab delta E squared */
+ desq = dlsq + da * da + db * db;
+ _desq[0][0] = 2.0 * dl;
+ _desq[1][0] = -2.0 * dl;
+ _desq[0][1] = 2.0 * da;
+ _desq[1][1] = -2.0 * da;
+ _desq[0][2] = 2.0 * db;
+ _desq[1][2] = -2.0 * db;
+ }
+
+ {
+ double c1, c2, dc, tt;
+
+ /* Compute chromanance for the two colors */
+ c1 = sqrt(Lab0[1] * Lab0[1] + Lab0[2] * Lab0[2]);
+ c2 = sqrt(Lab1[1] * Lab1[1] + Lab1[2] * Lab1[2]);
+ c12 = sqrt(c1 * c2); /* Symetric chromanance */
+
+ tt = 0.5 * (pow(c2, 0.5) + 1e-12)/(pow(c1, 1.5) + 1e-12);
+ _c12[0][0] = Lab0[1] * tt;
+ _c12[0][1] = Lab0[2] * tt;
+ tt = 0.5 * (pow(c1, 0.5) + 1e-12)/(pow(c2, 1.5) + 1e-12);
+ _c12[1][0] = Lab1[1] * tt;
+ _c12[1][1] = Lab1[2] * tt;
+
+ /* delta chromanance squared */
+ dc = c2 - c1;
+ dcsq = dc * dc;
+ if (c1 < 1e-12 || c2 < 1e-12) {
+ c1 += 1e-12;
+ c2 += 1e-12;
+ }
+ _dcsq[0][0] = -2.0 * Lab0[1] * (c2 - c1)/c1;
+ _dcsq[0][1] = -2.0 * Lab0[2] * (c2 - c1)/c1;
+ _dcsq[1][0] = 2.0 * Lab1[1] * (c2 - c1)/c2;
+ _dcsq[1][1] = 2.0 * Lab1[2] * (c2 - c1)/c2;
+ }
+
+ /* Compute delta hue squared */
+ dhsq = desq - dlsq - dcsq;
+ if (dhsq >= 0.0) {
+ _dhsq[0][0] = _desq[0][1] - _dcsq[0][0];
+ _dhsq[0][1] = _desq[0][2] - _dcsq[0][1];
+ _dhsq[1][0] = _desq[1][1] - _dcsq[1][0];
+ _dhsq[1][1] = _desq[1][2] - _dcsq[1][1];
+ } else {
+ dhsq = 0.0;
+ _dhsq[0][0] = 0.0;
+ _dhsq[0][1] = 0.0;
+ _dhsq[1][0] = 0.0;
+ _dhsq[1][1] = 0.0;
+ }
+
+ {
+ double sc, scsq, scf;
+ double sh, shsq, shf;
+
+ /* Weighting factors for delta chromanance & delta hue */
+ sc = 1.0 + 0.048 * c12;
+ scsq = sc * sc;
+
+ sh = 1.0 + 0.014 * c12;
+ shsq = sh * sh;
+
+ rv = dlsq + dcsq/scsq + dhsq/shsq;
+
+ scf = 0.048 * -2.0 * dcsq/(scsq * sc);
+ shf = 0.014 * -2.0 * dhsq/(shsq * sh);
+ dout[0][0] = _desq[0][0];
+ dout[0][1] = _dcsq[0][0]/scsq + _c12[0][0] * scf
+ + _dhsq[0][0]/shsq + _c12[0][0] * shf;
+ dout[0][2] = _dcsq[0][1]/scsq + _c12[0][1] * scf
+ + _dhsq[0][1]/shsq + _c12[0][1] * shf;
+ dout[1][0] = _desq[1][0];
+ dout[1][1] = _dcsq[1][0]/scsq + _c12[1][0] * scf
+ + _dhsq[1][0]/shsq + _c12[1][0] * shf;
+ dout[1][2] = _dcsq[1][1]/scsq + _c12[1][1] * scf
+ + _dhsq[1][1]/shsq + _c12[1][1] * shf;
+ return rv;
+ }
+}
+
+// ~~99 not sure if these are correct:
+
+/* Return the normal Delta E given two Lab values, */
+/* including partial derivatives. */
+double icxdLabDE(double dout[2][3], double *Lab0, double *Lab1) {
+ double rv = 0.0, tt;
+
+ tt = Lab0[0] - Lab1[0];
+ dout[0][0] = 1.0 * tt;
+ dout[1][0] = -1.0 * tt;
+ rv += tt * tt;
+ tt = Lab0[1] - Lab1[1];
+ dout[0][1] = 1.0 * tt;
+ dout[1][1] = -1.0 * tt;
+ rv += tt * tt;
+ tt = Lab0[2] - Lab1[2];
+ dout[0][2] = 1.0 * tt;
+ dout[1][2] = -1.0 * tt;
+ rv += tt * tt;
+ return sqrt(rv);
+}
+
+/* Return the CIE94 Delta E color difference measure */
+/* including partial derivatives. */
+double icxdCIE94(double dout[2][3], double Lab0[3], double Lab1[3]) {
+ double desq, _desq[2][3];
+ double dlsq;
+ double dcsq, _dcsq[2][2]; /* == [x][1,2] */
+ double c12, _c12[2][2]; /* == [x][1,2] */
+ double dhsq, _dhsq[2][2]; /* == [x][1,2] */
+ double rv;
+
+ {
+ double dl, da, db;
+
+ dl = Lab0[0] - Lab1[0];
+ dlsq = dl * dl; /* dl squared */
+ da = Lab0[1] - Lab1[1];
+ db = Lab0[2] - Lab1[2];
+
+
+ /* Compute normal Lab delta E squared */
+ desq = dlsq + da * da + db * db;
+ _desq[0][0] = 1.0 * dl;
+ _desq[1][0] = -1.0 * dl;
+ _desq[0][1] = 1.0 * da;
+ _desq[1][1] = -1.0 * da;
+ _desq[0][2] = 1.0 * db;
+ _desq[1][2] = -1.0 * db;
+ }
+
+ {
+ double c1, c2, dc, tt;
+
+ /* Compute chromanance for the two colors */
+ c1 = sqrt(Lab0[1] * Lab0[1] + Lab0[2] * Lab0[2]);
+ c2 = sqrt(Lab1[1] * Lab1[1] + Lab1[2] * Lab1[2]);
+ c12 = sqrt(c1 * c2); /* Symetric chromanance */
+
+ tt = 0.5 * (pow(c2, 0.5) + 1e-12)/(pow(c1, 1.5) + 1e-12);
+ _c12[0][0] = Lab0[1] * tt;
+ _c12[0][1] = Lab0[2] * tt;
+ tt = 0.5 * (pow(c1, 0.5) + 1e-12)/(pow(c2, 1.5) + 1e-12);
+ _c12[1][0] = Lab1[1] * tt;
+ _c12[1][1] = Lab1[2] * tt;
+
+ /* delta chromanance squared */
+ dc = c2 - c1;
+ dcsq = dc * dc;
+ if (c1 < 1e-12 || c2 < 1e-12) {
+ c1 += 1e-12;
+ c2 += 1e-12;
+ }
+ _dcsq[0][0] = -1.0 * Lab0[1] * (c2 - c1)/c1;
+ _dcsq[0][1] = -1.0 * Lab0[2] * (c2 - c1)/c1;
+ _dcsq[1][0] = 1.0 * Lab1[1] * (c2 - c1)/c2;
+ _dcsq[1][1] = 1.0 * Lab1[2] * (c2 - c1)/c2;
+ }
+
+ /* Compute delta hue squared */
+ dhsq = desq - dlsq - dcsq;
+ if (dhsq >= 0.0) {
+ _dhsq[0][0] = _desq[0][1] - _dcsq[0][0];
+ _dhsq[0][1] = _desq[0][2] - _dcsq[0][1];
+ _dhsq[1][0] = _desq[1][1] - _dcsq[1][0];
+ _dhsq[1][1] = _desq[1][2] - _dcsq[1][1];
+ } else {
+ dhsq = 0.0;
+ _dhsq[0][0] = 0.0;
+ _dhsq[0][1] = 0.0;
+ _dhsq[1][0] = 0.0;
+ _dhsq[1][1] = 0.0;
+ }
+
+ {
+ double sc, scsq, scf;
+ double sh, shsq, shf;
+
+ /* Weighting factors for delta chromanance & delta hue */
+ sc = 1.0 + 0.048 * c12;
+ scsq = sc * sc;
+
+ sh = 1.0 + 0.014 * c12;
+ shsq = sh * sh;
+
+ rv = dlsq + dcsq/scsq + dhsq/shsq;
+
+ scf = 0.048 * -1.0 * dcsq/(scsq * sc);
+ shf = 0.014 * -1.0 * dhsq/(shsq * sh);
+ dout[0][0] = _desq[0][0];
+ dout[0][1] = _dcsq[0][0]/scsq + _c12[0][0] * scf
+ + _dhsq[0][0]/shsq + _c12[0][0] * shf;
+ dout[0][2] = _dcsq[0][1]/scsq + _c12[0][1] * scf
+ + _dhsq[0][1]/shsq + _c12[0][1] * shf;
+ dout[1][0] = _desq[1][0];
+ dout[1][1] = _dcsq[1][0]/scsq + _c12[1][0] * scf
+ + _dhsq[1][0]/shsq + _c12[1][0] * shf;
+ dout[1][2] = _dcsq[1][1]/scsq + _c12[1][1] * scf
+ + _dhsq[1][1]/shsq + _c12[1][1] * shf;
+ return sqrt(rv);
+ }
+}
+
+/* ------------------------------------------------------ */
+/* A power-like function, based on Graphics Gems adjustment curve. */
+/* Avoids "toe" problem of pure power. */
+/* Adjusted so that "power" 2 and 0.5 agree with real power at 0.5 */
+
+double icx_powlike(double vv, double pp) {
+ double tt, g;
+
+ if (pp >= 1.0) {
+ g = 2.0 * (pp - 1.0);
+ vv = vv/(g - g * vv + 1.0);
+ } else {
+ g = 2.0 - 2.0/pp;
+ vv = (vv - g * vv)/(1.0 - g * vv);
+ }
+
+ return vv;
+}
+
+/* Compute the necessary aproximate power, to transform */
+/* the given value from src to dst. They are assumed to be */
+/* in the range 0.0 .. 1.0 */
+double icx_powlike_needed(double src, double dst) {
+ double pp, g;
+
+ if (dst <= src) {
+ g = -((src - dst)/(dst * src - dst));
+ pp = (0.5 * g) + 1.0;
+ } else {
+ g = -((src - dst)/((dst - 1.0) * src));
+ pp = 1.0/(1.0 - 0.5 * g);
+ }
+
+ return pp;
+}
+
+/* ------------------------------------------------------ */
+/* Parameterized transfer/dot gain function. */
+/* Used for device modelling. Including partial */
+/* derivative for input and parameters. */
+
+/* NOTE that clamping the input values seems to cause */
+/* conjgrad() problems. */
+
+/* Transfer function */
+double icxTransFunc(
+double *v, /* Pointer to first parameter */
+int luord, /* Number of parameters */
+double vv /* Source of value */
+) {
+ double g;
+ int ord;
+
+ /* Process all the shaper orders from low to high. */
+ /* [These shapers were inspired by a Graphics Gem idea */
+ /* (Gems IV, VI.3, "Fast Alternatives to Perlin's Bias and */
+ /* Gain Functions, pp 401). */
+ /* They have the nice properties that they are smooth, and */
+ /* are monotonic. The control parameter has been */
+ /* altered to have a range from -oo to +oo rather than 0.0 to 1.0 */
+ /* so that the search space is less non-linear. */
+ for (ord = 0; ord < luord; ord++) {
+ int nsec; /* Number of sections */
+ double sec; /* Section */
+
+ g = v[ord]; /* Parameter */
+
+ nsec = ord + 1; /* Increase sections for each order */
+
+ vv *= (double)nsec;
+
+ sec = floor(vv);
+ if (((int)sec) & 1)
+ g = -g; /* Alternate action in each section */
+ vv -= sec;
+ if (g >= 0.0) {
+ vv = vv/(g - g * vv + 1.0);
+ } else {
+ vv = (vv - g * vv)/(1.0 - g * vv);
+ }
+ vv += sec;
+ vv /= (double)nsec;
+ }
+
+ return vv;
+}
+
+/* Inverse transfer function */
+double icxInvTransFunc(
+double *v, /* Pointer to first parameter */
+int luord, /* Number of parameters */
+double vv /* Source of value */
+) {
+ double g;
+ int ord;
+
+ /* Process the shaper orders in reverse from high to low. */
+ /* [These shapers were inspired by a Graphics Gem idea */
+ /* (Gems IV, VI.3, "Fast Alternatives to Perlin's Bias and */
+ /* Gain Functions, pp 401). */
+ /* They have the nice properties that they are smooth, and */
+ /* are monotonic. The control parameter has been */
+ /* altered to have a range from -oo to +oo rather than 0.0 to 1.0 */
+ /* so that the search space is less non-linear. */
+ for (ord = luord-1; ord >= 0; ord--) {
+ int nsec; /* Number of sections */
+ double sec; /* Section */
+
+ g = -v[ord]; /* Inverse parameter */
+
+ nsec = ord + 1; /* Increase sections for each order */
+
+ vv *= (double)nsec;
+
+ sec = floor(vv);
+ if (((int)sec) & 1)
+ g = -g; /* Alternate action in each section */
+ vv -= sec;
+ if (g >= 0.0) {
+ vv = vv/(g - g * vv + 1.0);
+ } else {
+ vv = (vv - g * vv)/(1.0 - g * vv);
+ }
+ vv += sec;
+ vv /= (double)nsec;
+ }
+
+ return vv;
+}
+
+/* Transfer function with offset and scale */
+double icxSTransFunc(
+double *v, /* Pointer to first parameter */
+int luord, /* Number of parameters */
+double vv, /* Source of value */
+double min, /* Scale values */
+double max
+) {
+ max -= min;
+
+ vv = (vv - min)/max;
+ vv = icxTransFunc(v, luord, vv);
+ vv = (vv * max) + min;
+ return vv;
+}
+
+/* Inverse Transfer function with offset and scale */
+double icxInvSTransFunc(
+double *v, /* Pointer to first parameter */
+int luord, /* Number of parameters */
+double vv, /* Source of value */
+double min, /* Scale values */
+double max
+) {
+ max -= min;
+
+ vv = (vv - min)/max;
+ vv = icxInvTransFunc(v, luord, vv);
+ vv = (vv * max) + min;
+ return vv;
+}
+
+/* Transfer function with partial derivative */
+/* with respect to the parameters. */
+double icxdpTransFunc(
+double *v, /* Pointer to first parameter */
+double *dv, /* Return derivative wrt each parameter */
+int luord, /* Number of parameters */
+double vv /* Source of value */
+) {
+ double g;
+ int i, ord;
+
+ /* Process all the shaper orders from high to low. */
+ for (ord = 0; ord < luord; ord++) {
+ double dsv; /* del for del in g */
+ double ddv; /* del for del in vv */
+ int nsec; /* Number of sections */
+ double sec; /* Section */
+
+ g = v[ord]; /* Parameter */
+
+ nsec = ord + 1; /* Increase sections for each order */
+
+ vv *= (double)nsec;
+
+ sec = floor(vv);
+ if (((int)sec) & 1) {
+ g = -g; /* Alternate action in each section */
+ }
+ vv -= sec;
+ if (g >= 0.0) {
+ double tt = g - g * vv + 1.0;
+ dsv = (vv * vv - vv)/(tt * tt);
+ ddv = (g + 1.0)/(tt * tt);
+ vv = vv/tt;
+ } else {
+ double tt = 1.0 - g * vv;
+ dsv = (vv * vv - vv)/(tt * tt);
+ ddv = (1.0 - g)/(tt * tt);
+ vv = (vv - g * vv)/tt;
+ }
+
+ vv += sec;
+ vv /= (double)nsec;
+ dsv /= (double)nsec;
+ if (((int)sec) & 1)
+ dsv = -dsv;
+
+ dv[ord] = dsv;
+ for (i = ord - 1; i >= 0; i--)
+ dv[i] *= ddv;
+ }
+
+ return vv;
+}
+
+/* Transfer function with offset and scale, and */
+/* partial derivative with respect to the parameters. */
+double icxdpSTransFunc(
+double *v, /* Pointer to first parameter */
+double *dv, /* Return derivative wrt each parameter */
+int luord, /* Number of parameters */
+double vv, /* Source of value */
+double min, /* Scale values */
+double max
+) {
+ int i;
+ max -= min;
+
+ vv = (vv - min)/max;
+ vv = icxdpTransFunc(v, dv, luord, vv);
+ vv = (vv * max) + min;
+ for (i = 0; i < luord; i++)
+ dv[i] *= max;
+ return vv;
+}
+
+
+/* Transfer function with partial derivative */
+/* with respect to the input value. */
+double icxdiTransFunc(
+double *v, /* Pointer to first parameter */
+double *pdin, /* Return derivative wrt source value */
+int luord, /* Number of parameters */
+double vv /* Source of value */
+) {
+ double g, din;
+ int ord;
+
+#ifdef NEVER
+ if (vv < 0.0 || vv > 1.0) {
+ if (vv < 0.0)
+ vv = 0.0;
+ else
+ vv = 1.0;
+
+ *pdin = 0.0;
+ return vv;
+ }
+#endif
+ din = 1.0;
+
+ /* Process all the shaper orders from high to low. */
+ for (ord = 0; ord < luord; ord++) {
+ double ddv; /* del for del in vv */
+ int nsec; /* Number of sections */
+ double sec; /* Section */
+
+ g = v[ord]; /* Parameter */
+
+ nsec = ord + 1; /* Increase sections for each order */
+
+ vv *= (double)nsec;
+
+ sec = floor(vv);
+ if (((int)sec) & 1) {
+ g = -g; /* Alternate action in each section */
+ }
+ vv -= sec;
+ if (g >= 0.0) {
+ double tt = g - g * vv + 1.0;
+ ddv = (g + 1.0)/(tt * tt);
+ vv = vv/tt;
+ } else {
+ double tt = 1.0 - g * vv;
+ ddv = (1.0 - g)/(tt * tt);
+ vv = (vv - g * vv)/tt;
+ }
+
+ vv += sec;
+ vv /= (double)nsec;
+ din *= ddv;
+ }
+
+ *pdin = din;
+ return vv;
+}
+
+/* Transfer function with offset and scale, and */
+/* partial derivative with respect to the input value. */
+double icxdiSTransFunc(
+double *v, /* Pointer to first parameter */
+double *pdv, /* Return derivative wrt source value */
+int luord, /* Number of parameters */
+double vv, /* Source of value */
+double min, /* Scale values */
+double max
+) {
+ max -= min;
+
+ vv = (vv - min)/max;
+ vv = icxdiTransFunc(v, pdv, luord, vv);
+ vv = (vv * max) + min;
+ return vv;
+}
+
+
+/* Transfer function with partial derivative */
+/* with respect to the parameters and the input value. */
+double icxdpdiTransFunc(
+double *v, /* Pointer to first parameter */
+double *dv, /* Return derivative wrt each parameter */
+double *pdin, /* Return derivative wrt source value */
+int luord, /* Number of parameters */
+double vv /* Source of value */
+) {
+ double g, din;
+ int i, ord;
+
+#ifdef NEVER
+ if (vv < 0.0 || vv > 1.0) {
+ if (vv < 0.0)
+ vv = 0.0;
+ else
+ vv = 1.0;
+
+ for (ord = 0; ord < luord; ord++)
+ dv[ord] = 0.0;
+ *pdin = 0.0;
+ return vv;
+ }
+#endif
+ din = 1.0;
+
+ /* Process all the shaper orders from high to low. */
+ for (ord = 0; ord < luord; ord++) {
+ double dsv; /* del for del in g */
+ double ddv; /* del for del in vv */
+ int nsec; /* Number of sections */
+ double sec; /* Section */
+
+ g = v[ord]; /* Parameter */
+
+ nsec = ord + 1; /* Increase sections for each order */
+
+ vv *= (double)nsec;
+
+ sec = floor(vv);
+ if (((int)sec) & 1) {
+ g = -g; /* Alternate action in each section */
+ }
+ vv -= sec;
+ if (g >= 0.0) {
+ double tt = g - g * vv + 1.0;
+ dsv = (vv * vv - vv)/(tt * tt);
+ ddv = (g + 1.0)/(tt * tt);
+ vv = vv/tt;
+ } else {
+ double tt = 1.0 - g * vv;
+ dsv = (vv * vv - vv)/(tt * tt);
+ ddv = (1.0 - g)/(tt * tt);
+ vv = (vv - g * vv)/tt;
+ }
+
+ vv += sec;
+ vv /= (double)nsec;
+ dsv /= (double)nsec;
+ if (((int)sec) & 1)
+ dsv = -dsv;
+
+ dv[ord] = dsv;
+ for (i = ord - 1; i >= 0; i--)
+ dv[i] *= ddv;
+ din *= ddv;
+ }
+
+ *pdin = din;
+ return vv;
+}
+
+/* Transfer function with offset and scale, and */
+/* partial derivative with respect to the */
+/* parameters and the input value. */
+double icxdpdiSTransFunc(
+double *v, /* Pointer to first parameter */
+double *dv, /* Return derivative wrt each parameter */
+double *pdin, /* Return derivative wrt source value */
+int luord, /* Number of parameters */
+double vv, /* Source of value */
+double min, /* Scale values */
+double max
+) {
+ int i;
+ max -= min;
+
+ vv = (vv - min)/max;
+ vv = icxdpdiTransFunc(v, dv, pdin, luord, vv);
+ vv = (vv * max) + min;
+ for (i = 0; i < luord; i++)
+ dv[i] *= max;
+ return vv;
+}
+
+/* ------------------------------------------------------ */
+/* Multi-plane interpolation, used for device modelling. */
+/* Including partial derivative for input and parameters. */
+/* A simple flat plane is used for each output. */
+
+/* Multi-plane interpolation - uses base + di slope values. */
+/* Parameters are assumed to be fdi groups of di + 1 parameters. */
+void icxPlaneInterp(
+double *v, /* Pointer to first parameter [fdi * (di + 1)] */
+int fdi, /* Number of output channels */
+int di, /* Number of input channels */
+double *out, /* Resulting fdi values */
+double *in /* Input di values */
+) {
+ int e, f;
+
+ for (f = 0; f < fdi; f++) {
+ for (out[f] = 0.0, e = 0; e < di; e++, v++) {
+ out[f] += in[e] * *v;
+ }
+ out[f] += *v;
+ }
+}
+
+
+/* Multii-plane interpolation with partial derivative */
+/* with respect to the input and parameters. */
+void icxdpdiPlaneInterp(
+double *v, /* Pointer to first parameter value [fdi * (di + 1)] */
+double *dv, /* Return [1 + di] deriv. wrt each parameter v */
+double *din, /* Return [fdi * di] deriv. wrt each input value */
+int fdi, /* Number of output channels */
+int di, /* Number of input channels */
+double *out, /* Resulting fdi values */
+double *in /* Input di values */
+) {
+ int e, ee, f, g;
+ int dip2 = (di + 1); /* Output dim increment through parameters */
+
+ /* Compute the output values */
+ for (f = 0; f < fdi; f++) {
+ for (out[f] = 0.0, e = 0; e < di; e++)
+ out[f] += in[e] * v[f * dip2 + e];
+ out[f] += v[f * dip2 + e];
+ }
+
+ /* Since interpolation is verys simple, derivative are also simple */
+
+ /* Copy del for parameter to return array */
+ for (e = 0; e < di; e++)
+ dv[e] = in[e];
+ dv[e] = 1.0;
+
+ /* Compute del of out[] from in[] */
+ for (f = 0; f < fdi; f++) {
+ for (e = 0; e < di; e++) {
+ din[f * di + e] = v[f * dip2 + e];
+ }
+ }
+}
+
+/* ------------------------------------------------------ */
+/* Matrix cube interpolation, used for device modelling. */
+/* Including partial derivative for input and parameters. */
+
+/* Matrix cube interpolation - interpolate between 2^di output corner values. */
+/* Parameters are assumed to be fdi groups of 2^di parameters. */
+void icxCubeInterp(
+double *v, /* Pointer to first parameter */
+int fdi, /* Number of output channels */
+int di, /* Number of input channels */
+double *out, /* Resulting fdi values */
+double *in /* Input di values */
+) {
+ int e, f, g;
+ double gw[1 << MXDI]; /* weight for each matrix grid cube corner */
+
+ /* Compute corner weights needed for interpolation */
+ gw[0] = 1.0;
+ for (e = 0, g = 1; e < di; e++, g *= 2) {
+ int i;
+ for (i = 0; i < g; i++) {
+ gw[g+i] = gw[i] * in[e];
+ gw[i] *= (1.0 - in[e]);
+ }
+ }
+
+ /* Now compute the output values */
+ for (f = 0; f < fdi; f++) {
+ out[f] = 0.0; /* For each output value */
+ for (e = 0; e < (1 << di); e++) { /* For all corners of cube */
+ out[f] += gw[e] * *v;
+ v++;
+ }
+ }
+}
+
+
+/* Matrix cube interpolation. with partial derivative */
+/* with respect to the input and parameters. */
+void icxdpdiCubeInterp(
+double *v, /* Pointer to first parameter value [fdi * 2^di] */
+double *dv, /* Return [2^di] deriv. wrt each parameter v */
+double *din, /* Return [fdi * di] deriv. wrt each input value */
+int fdi, /* Number of output channels */
+int di, /* Number of input channels */
+double *out, /* Resulting fdi values */
+double *in /* Input di values */
+) {
+ int e, ee, f, g;
+ int dip2 = (1 << di);
+ double gw[1 << MXDI]; /* weight for each matrix grid cube corner */
+
+ /* Compute corner weights needed for interpolation */
+ gw[0] = 1.0;
+ for (e = 0, g = 1; e < di; e++, g *= 2) {
+ int i;
+ for (i = 0; i < g; i++) {
+ gw[g+i] = gw[i] * in[e];
+ gw[i] *= (1.0 - in[e]);
+ }
+ }
+
+ /* Now compute the output values */
+ for (f = 0; f < fdi; f++) {
+ out[f] = 0.0; /* For each output value */
+ for (ee = 0; ee < dip2; ee++) { /* For all corners of cube */
+ out[f] += gw[ee] * v[f * dip2 + ee];
+ }
+ }
+
+ /* Copy del for parameter to return array */
+ for (ee = 0; ee < dip2; ee++) { /* For all other corners of cube */
+ dv[ee] = gw[ee]; /* del from parameter */
+ }
+
+ /* Compute del from in[] value we want */
+ for (e = 0; e < di; e++) { /* For input we want del wrt */
+
+ for (f = 0; f < fdi; f++)
+ din[f * di + e] = 0.0; /* Zero del ready of accumulation */
+
+ for (ee = 0; ee < dip2; ee++) { /* For all corners of cube weights, */
+ int e2; /* accumulate del from in[] we want. */
+ double vv = 1.0;
+
+ /* Compute in[] weighted cube corners for all except del of in[] we want */
+ for (e2 = 0; e2 < di; e2++) { /* const from non del inputs */
+ if (e2 == e)
+ continue;
+ if (ee & (1 << e2))
+ vv *= in[e2];
+ else
+ vv *= (1.0 - in[e2]);
+ }
+
+ /* Accumulate contribution of in[] we want for corner to out[] we want */
+ if (ee & (1 << e)) {
+ for (f = 0; f < fdi; f++)
+ din[f * di + e] += v[f * dip2 + ee] * vv;
+ } else {
+ for (f = 0; f < fdi; f++)
+ din[f * di + e] -= v[f * dip2 + ee] * vv;
+ }
+ }
+ }
+}
+
+/* ------------------------------------------------------ */
+/* Matrix cube simplex interpolation, used for device modelling. */
+
+/* Matrix cube simplex interpolation - interpolate between 2^di output corner values. */
+/* Parameters are assumed to be fdi groups of 2^di parameters. */
+void icxCubeSxInterp(
+double *v, /* Pointer to first parameter */
+int fdi, /* Number of output channels */
+int di, /* Number of input channels */
+double *out, /* Resulting fdi values */
+double *in /* Input di values */
+) {
+ int si[MAX_CHAN]; /* in[] Sort index, [0] = smalest */
+
+//{
+// double tout[MXDO];
+//
+// icxCubeInterp(v, fdi, di, tout, in);
+//printf("\n~1 Cube interp result = %f\n",tout[0]);
+//}
+
+//printf("~1 icxCubeSxInterp: %f %f %f\n", in[0], in[1], in[2]);
+ /* Do insertion sort on coordinates, smallest to largest. */
+ {
+ int ff, vf;
+ unsigned int e;
+ double v;
+ for (e = 0; e < di; e++)
+ si[e] = e; /* Initial unsorted indexes */
+
+ for (e = 1; e < di; e++) {
+ ff = e;
+ v = in[si[ff]];
+ vf = ff;
+ while (ff > 0 && in[si[ff-1]] > v) {
+ si[ff] = si[ff-1];
+ ff--;
+ }
+ si[ff] = vf;
+ }
+ }
+//printf("~1 sort order %d %d %d\n", si[0], si[1], si[2]);
+//printf(" from %f %f %f\n", in[si[0]], in[si[1]], in[si[2]]);
+
+ /* Now compute the weightings, simplex vertices and output values */
+ {
+ unsigned int e, f;
+ double w; /* Current vertex weight */
+
+ w = 1.0 - in[si[di-1]]; /* Vertex at base of cell */
+ for (f = 0; f < fdi; f++) {
+ out[f] = w * v[f * (1 << di)];
+//printf("~1 out[%d] = %f = %f * %f\n",f,out[f],w,v[f * (1 << di)]);
+ }
+
+ for (e = di-1; e > 0; e--) { /* Middle verticies */
+ w = in[si[e]] - in[si[e-1]];
+ v += (1 << si[e]); /* Move to top of cell in next largest dimension */
+ for (f = 0; f < fdi; f++) {
+ out[f] += w * v[f * (1 << di)];
+//printf("~1 out[%d] = %f += %f * %f\n",f,out[f],w,v[f * (1 << di)]);
+ }
+ }
+
+ w = in[si[0]];
+ v += (1 << si[0]); /* Far corner from base of cell */
+ for (f = 0; f < fdi; f++) {
+ out[f] += w * v[f * (1 << di)];
+//printf("~1 out[%d] = %f += %f * %f\n",f,out[f],w,v[f * (1 << di)]);
+ }
+ }
+}
+
+/* ------------------------------------------------------ */
+/* Matrix multiplication, used for device modelling. */
+/* Including partial derivative for input and parameters. */
+
+
+/* 3x3 matrix multiplication, with the matrix in a 1D array */
+/* with respect to the input and parameters. */
+void icxMulBy3x3Parm(
+ double out[3], /* Return input multiplied by matrix */
+ double mat[9], /* Matrix organised in [slow][fast] order */
+ double in[3] /* Input values */
+) {
+ double *v = mat, ov[3];
+ int e, f;
+
+ /* Compute the output values */
+ for (f = 0; f < 3; f++) {
+ ov[f] = 0.0; /* For each output value */
+ for (e = 0; e < 3; e++) {
+ ov[f] += *v++ * in[e];
+ }
+ }
+ out[0] = ov[0];
+ out[1] = ov[1];
+ out[2] = ov[2];
+}
+
+
+/* 3x3 matrix multiplication, with partial derivatives */
+/* with respect to the input and parameters. */
+void icxdpdiMulBy3x3Parm(
+ double out[3], /* Return input multiplied by matrix */
+ double dv[3][9], /* Return deriv for each [output] with respect to [param] */
+ double din[3][3], /* Return deriv for each [output] with respect to [input] */
+ double mat[9], /* Matrix organised in [slow][fast] order */
+ double in[3] /* Input values */
+) {
+ double *v, ov[3];
+ int e, f;
+
+ /* Compute the output values */
+ v = mat;
+ for (f = 0; f < 3; f++) {
+ ov[f] = 0.0; /* For each output value */
+ for (e = 0; e < 3; e++) {
+ ov[f] += *v++ * in[e];
+ }
+ }
+
+ /* Compute deriv. with respect to the matrix parameter % 3 */
+ /* This is pretty simple for a matrix ... */
+ for (f = 0; f < 3; f++) {
+ for (e = 0; e < 9; e++) {
+ if (e/3 == f)
+ dv[f][e] = in[e % 3];
+ else
+ dv[f][e] = 0.0;
+ }
+ }
+
+ /* Compute deriv. with respect to the input values */
+ /* This is pretty simple for a matrix ... */
+ v = mat;
+ for (f = 0; f < 3; f++)
+ for (e = 0; e < 3; e++)
+ din[f][e] = *v++;
+
+ out[0] = ov[0];
+ out[1] = ov[1];
+ out[2] = ov[2];
+}
+
+#undef stricmp
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+