summaryrefslogtreecommitdiff
path: root/xicc/xlut.c
diff options
context:
space:
mode:
Diffstat (limited to 'xicc/xlut.c')
-rw-r--r--xicc/xlut.c669
1 files changed, 512 insertions, 157 deletions
diff --git a/xicc/xlut.c b/xicc/xlut.c
index 827dde8..a18956c 100644
--- a/xicc/xlut.c
+++ b/xicc/xlut.c
@@ -54,7 +54,7 @@
/*
- ~~~~~!!!!! This has all been fixed ?
+ !!!!! This has all been fixed ? !!!!!!
NOTE :- an alternative to the way display profile absolute is handled here
would be to always chromatically adapt the illuminant to D50, and encode
@@ -160,6 +160,9 @@
#undef DISABLE_KCURVE_FILTER /* [Undef] don't filter the Kcurve */
#undef REPORT_LOCUS_SEGMENTS /* [Undef[ Examine how many segments there are in aux inversion */
+#undef FASTREVSETUP_NON_CAM /* [Undef] Use fast setup on innerm non-CAM lookup, if we're */
+ /* going to use CAM clip for nn lookup */
+
#define XYZ_EXTRA_SMOOTH 20.0 /* Extra smoothing factor for XYZ profiles */
/* !!! Note this is mainly due to smoothing being */
/* scaled by data range in rspl code !!! */
@@ -171,7 +174,10 @@
/* Should this be smaller ? */
#undef USECAMCLIPSPLINE /* [Und] use spline blend between PCS and Jab */
-#define CCJSCALE 2.0 /* [2.0] Amount to emphasize J delta E in in computing clip */
+#define USELCHWEIGHT /* [def] Use LCh weighting for nn clip if possible */
+#define JCCWEIGHT 2.0 /* [2.0] Amount to emphasize J delta E in in computing clip */
+#define CCCWEIGHT 1.0 /* [1.0] Amount to emphasize C delta E in in computing clip */
+#define HCCWEIGHT 2.2 /* [2.2] Amount to emphasize H delta E in in computing clip */
/*
* TTBD:
@@ -429,7 +435,6 @@ icColorSpaceSignature is, /* Input space, XYZ or Lab */
double *out, double *in) {
icxLuLut *p = (icxLuLut *)pp;
-
/* Convert to the ICC PCS */
if (is == icSigLabData && p->natpcs == icSigXYZData) {
DBS(("fwd_relpcs_outpcs: Lab in = %s\n", icmPdv(p->inputChan, in)));
@@ -444,7 +449,7 @@ double *out, double *in) {
icmCpy3(out, in);
}
- /* Convert to absolute */
+ /* Convert to final PCS (i.e. absolute intent is set that way) */
((icmLuLut *)p->plu)->out_abs((icmLuLut *)p->plu, out, out);
DBS(("fwd_relpcs_outpcs: abs PCS = %s\n", icmPdv(p->inputChan, out)));
@@ -461,51 +466,98 @@ double *out, double *in) {
/* - - - - - - - - - - - - - - - - - - - - - - - - - - */
/* Components of inverse lookup, in order */
+static double sigfunc(double in, double pp) {
+ if (in < 0.5)
+ return pow(2 * in, pp) * 0.5;
+ else
+ return 1.0 - (pow(2 * (1.0 - in), pp) * 0.5);
+}
+
/* Utility function - compute the clip vector direction. */
/* return NULL if vector clip isn't used. */
double *icxClipVector(
icxClip *p, /* Clipping setup information */
double *in, /* Target point */
-double *cdirv /* Space for returned clip vector */
+double *cdirv, /* Returned clip vector */
+int safe /* Flag - return safe vector */
) {
int f;
if (p->nearclip != 0)
return NULL; /* Doing nearest clipping, not vector */
- /* Default is simple vector clip */
- for (f = 0; f < p->fdi; f++)
- cdirv[f] = p->ocent[f] - in[f]; /* Clip towards output gamut center */
+ if (p->cm == NULL) {
+ /* Default is simple vector clip */
+ for (f = 0; f < p->fdi; f++)
+ cdirv[f] = p->ocent[f] - in[f]; /* Clip towards output gamut center */
+
+ /* Else use CuspMap for more targeted vector */
+ } else {
+ double Cpow = 2.0; /* [2.0] 4.0 for less L change */
+ double Lsig = 2.5; /* [2.5] 1.6 for less L change */
+ double Lpow = 0.5; /* [0.5] 0.8 for less L change */
+ double Cratio = 0.9; /* [0.9] see cusptest.c */
+ double ratio;
+ double ss[3];
+ double inC;
+ double targ[3], cuspLCh[3];
- if (p->ocentl != 0.0) { /* Graduated vector clip */
- double cvl, nll;
+ inC = sqrt(in[1] * in[1] + in[2] * in[2]);
- /* Compute (negative) clip vector length */
- for (cvl = 0.0, f = 0; f < p->fdi; f++) {
- cvl += cdirv[f] * cdirv[f];
+ p->cm->getCusp(p->cm, cuspLCh, in);
+
+
+//icmLCh2Lab(targ, cuspLCh);
+//printf("\n~1 in %f %f %f -> cusp %f %f %f\n", in[0], in[1], in[2], targ[0], targ[1], targ[2]);
+//printf("~1 in %f %f %f -> cuspLCh %f %f %f\n", in[0], in[1], in[2], cuspLCh[0], cuspLCh[1], cuspLCh[2]);
+
+ /* Constrain clip vector to always be inwards pointing */
+ if (cuspLCh[1] > 0.90 * inC) {
+ cuspLCh[1] = 0.90 * inC;
+//printf("~1 Constrain cusp C: cuspLCh %f %f %f\n", cuspLCh[0], cuspLCh[1], cuspLCh[2]);
}
- cvl = sqrt(cvl);
- if (cvl > 1e-8) {
- /* Dot product of clip vector and clip center line */
- for (nll = 0.0, f = 0; f < p->fdi; f++)
- nll -= cdirv[f] * p->ocentv[f]; /* (Fix -ve clip vector) */
- nll /= (p->ocentl * p->ocentl); /* Normalised location along line */
-
- /* Limit to line */
- if (nll < 0.0)
- nll = 0.0;
- else if (nll > 1.0)
- nll = 1.0;
-
- if (p->LabLike) {
- /* Aim more towards center for saturated targets */
- double sat = sqrt(in[1] * in[1] + in[2] * in[2]);
- nll += sat/150.0 * (0.5 - nll);
- }
- /* Compute target clip direction */
- for (f = 0; f < p->fdi; f++)
- cdirv[f] = p->ocent[f] + nll * p->ocentv[f] - in[f];
+ ss[0] = in[0];
+ ss[1] = in[1];
+ ss[2] = in[2];
+
+ if (ss[0] < p->cm->Lmin[0])
+ ss[0] = p->cm->Lmin[0];
+ if (ss[0] > p->cm->Lmax[0])
+ ss[0] = p->cm->Lmax[0];
+
+ if (safe) {
+ targ[0] = cuspLCh[0]; /* L target is cusp L */
+ targ[1] = 0.0; /* Right on the neutral axis */
+
+ } else {
+ if (ss[0] >= cuspLCh[0]) {
+ ratio = (p->cm->Lmax[0] - ss[0])/(p->cm->Lmax[0] - cuspLCh[0]);
+ targ[0] = p->cm->Lmax[0] - sigfunc(pow(ratio, Lpow), Lsig)
+ * (p->cm->Lmax[0] - cuspLCh[0]);
+ targ[1] = pow(ratio, Cpow) * Cratio * cuspLCh[1];
+ } else {
+ ratio = (ss[0] - p->cm->Lmin[0])/(cuspLCh[0] - p->cm->Lmin[0]);
+ targ[0] = p->cm->Lmin[0] + sigfunc(pow(ratio, Lpow), Lsig)
+ * (cuspLCh[0] - p->cm->Lmin[0]);
+ targ[1] = pow(ratio, Cpow) * Cratio * cuspLCh[1];
+ }
}
+ targ[2] = cuspLCh[2]; /* h of in */
+
+//printf("~1 targLCH %f %f %f\n", targ[0], targ[1], targ[2]);
+ icmLCh2Lab(targ, targ);
+//printf("~1 targ %f %f %f\n", targ[0], targ[1], targ[2]);
+
+ /* line target point up with white and black point ? */
+ ratio = (ss[0] - p->cm->Lmin[0])/(p->cm->Lmax[0] - p->cm->Lmin[0]);
+ targ[1] += (1.0 - ratio) * p->cm->Lmin[1] + ratio * p->cm->Lmax[1];
+ targ[2] += (1.0 - ratio) * p->cm->Lmin[2] + ratio * p->cm->Lmax[2];
+
+//printf("~1 in %f %f %f -> targ %f %f %f\n", in[0], in[1], in[2], targ[0], targ[1], targ[2]);
+
+ /* Compute target clip direction */
+ for (f = 0; f < p->fdi; f++)
+ cdirv[f] = (targ[f] - in[f]);
}
return cdirv;
@@ -855,7 +907,7 @@ static double icxKcurveNF(double L, icxInkCurve *c) {
#define DBK(xxx)
#endif
-/* Same as above, but implement transition filters accross inflection points. */
+/* Same as above, but implement transition filters across inflection points. */
/* (The convolultion filter previously used could be */
/* re-instituted if something was done about compressing */
/* the filter at the boundaries so that the levels are met.) */
@@ -999,12 +1051,14 @@ double *clipd, /* If not NULL, return DE to gamut on clipi, 0 for not clip */
double *in /* Function input values to invert (== clut output' values) */
) {
co pp[MAX_INVSOLN]; /* Room for all the solutions found */
+ co upp; /* pp[0] value sent to rev_interp() for replay. */
int nsoln; /* Number of solutions found */
- double *cdir, cdirv[MXDO]; /* Clip vector direction and length */
+ double *cdir, cdirv[MXDO]; /* Clip vector direction and length/LCh weighting */
int e,f,i;
int fdi = p->clutTable->fdi;
int flags = 0; /* reverse interp flags */
int xflags = 0; /* extra clip/noclip flags */
+ int uflags = 0; /* flags value sent to rev_interp() for replay */
double tin[MXDO]; /* PCS value to be inverted */
double cdist = 0.0; /* clip DE */
int crv = 0; /* Return value - set to 1 if clipped */
@@ -1026,10 +1080,10 @@ double *in /* Function input values to invert (== clut output' values) */
/* Setup for reverse lookup */
for (f = 0; f < fdi; f++)
- pp[0].v[f] = in[f]; /* Target value */
+ upp.v[f] = pp[0].v[f] = in[f]; /* Target value */
/* Compute clip vector, if any */
- cdir = icxClipVector(&p->clip, in, cdirv);
+ cdir = icxClipVector(&p->clip, in, cdirv, 0);
if (p->clutTable->di > fdi) { /* ie. CMYK->Lab, there will be ambiguity */
double min[MXDI], max[MXDI]; /* Auxiliary locus range */
@@ -1084,7 +1138,7 @@ double *in /* Function input values to invert (== clut output' values) */
/* even though it won't have any effect. */
for (e = 0; e < p->clutTable->di; e++) {
if (p->auxm[e] != 0) {
- pp[0].p[e] = 0.5;
+ upp.p[e] = pp[0].p[e] = 0.5;
}
}
@@ -1123,7 +1177,7 @@ double *in /* Function input values to invert (== clut output' values) */
iv = min[e];
else if (iv > max[e])
iv = max[e];
- pp[0].p[e] = iv;
+ upp.p[e] = pp[0].p[e] = iv;
}
}
DBR(("inv_clut_aux: aux %f from auxt[] %f\n",pp[0].p[3],auxt[0]))
@@ -1137,7 +1191,7 @@ double *in /* Function input values to invert (== clut output' values) */
iv = min[e];
else if (iv > max[e])
iv = max[e];
- pp[0].p[e] = iv;
+ upp.p[e] = pp[0].p[e] = iv;
}
}
DBR(("inv_clut_aux: aux %f from out[0] K target %f min %f max %f\n",pp[0].p[3],out[3],min[3],max[3]))
@@ -1152,7 +1206,7 @@ double *in /* Function input values to invert (== clut output' values) */
iv = min[e];
else if (iv > max[e])
iv = max[e];
- pp[0].p[e] = iv;
+ upp.p[e] = pp[0].p[e] = iv;
}
}
DBR(("inv_clut_aux: aux %f from out[0] locus %f min %f max %f\n",pp[0].p[3],out[3],min[3],max[3]))
@@ -1195,7 +1249,7 @@ double *in /* Function input values to invert (== clut output' values) */
for (e = 0; e < p->clutTable->di; e++) {
if (p->auxm[e] != 0) {
- pp[0].p[e] = min[e] + rv * (max[e] - min[e]);
+ upp.p[e] = pp[0].p[e] = min[e] + rv * (max[e] - min[e]);
}
}
DBR(("inv_clut_aux: aux %f from locus %f min %f max %f\n",pp[0].p[3],rv,min[3],max[3]))
@@ -1209,7 +1263,7 @@ double *in /* Function input values to invert (== clut output' values) */
iv = min[e];
else if (iv > max[e])
iv = max[e];
- pp[0].p[e] = iv;
+ upp.p[e] = pp[0].p[e] = iv;
}
}
DBR(("inv_clut_aux: aux %f from out[0] K target %f min %f max %f\n",pp[0].p[3],rv,min[3],max[3]))
@@ -1241,7 +1295,7 @@ double *in /* Function input values to invert (== clut output' values) */
ii = 1.0;
ii = (1.0 - ii) * rv + ii * rv2;/* Blend between locus rule curves */
/* Out ink from output locus */
- pp[0].p[e] = min[e] + ii * (max[e] - min[e]);
+ upp.p[e] = pp[0].p[e] = min[e] + ii * (max[e] - min[e]);
} else {
double iv;
iv = out[e]; /* Input K level */
@@ -1249,7 +1303,7 @@ double *in /* Function input values to invert (== clut output' values) */
iv = rv;
else if (iv > rv2)
iv = rv2;
- pp[0].p[e] = iv;
+ upp.p[e] = pp[0].p[e] = iv;
}
}
}
@@ -1270,7 +1324,7 @@ double *in /* Function input values to invert (== clut output' values) */
tv = max[e];
tc.p[0] = tv;
p->inputTable[e]->interp(p->inputTable[e], &tc);
- pp[0].p[e] = tc.v[0];
+ upp.p[e] = pp[0].p[e] = tc.v[0];
}
}
@@ -1291,16 +1345,18 @@ double *in /* Function input values to invert (== clut output' values) */
tin[f] = pp[0].v[f];
}
+ uflags = RSPL_MAXAUX | flags | xflags; /* Combine all the flags */
+
/* Find reverse solution with target auxiliaries */
/* We choose the closest aux at or above the target */
/* to try and avoid glitches near black due to */
/* possible forked black locuses. */
nsoln = p->clutTable->rev_interp(
p->clutTable, /* rspl object */
- RSPL_MAXAUX | flags | xflags, /* Combine all the flags */
+ uflags,
MAX_INVSOLN, /* Maxumum solutions to return */
p->auxm, /* Auxiliary input chanel mask */
- cdir, /* Clip vector direction and length */
+ cdir, /* Clip vector direction/LCh weighting */
pp); /* Input target and output solutions */
/* returned solutions in pp[0..retval-1].p[] */
@@ -1312,13 +1368,15 @@ double *in /* Function input values to invert (== clut output' values) */
tin[f] = pp[0].v[f];
}
+ uflags = flags; /* No extra flags */
+
/* Color spaces don't need auxiliaries to choose from solution locus */
nsoln = p->clutTable->rev_interp(
p->clutTable, /* rspl object */
- flags, /* No extra flags */
+ uflags,
MAX_INVSOLN, /* Maxumum solutions to return */
NULL, /* No auxiliary input targets */
- cdir, /* Clip vector direction and length */
+ cdir, /* Clip vector direction/LCh weighting */
pp); /* Input target and output solutions */
/* returned solutions in pp[0..retval-1].p[] */
}
@@ -1343,8 +1401,9 @@ double *in /* Function input values to invert (== clut output' values) */
DBR(("inv_clut_aux got %d rev_interp solutions, clipflag = %d\n",nsoln,crv))
/* If we clipped and we should clip in CAM Jab space, compute reverse */
- /* clip solution using our additional CAM space. */
- /* (Note that we don't support vector clip in CAM space at the moment) */
+ /* clip solution using our additional CAM space rspl. */
+ /* Note that we don't support vector clip in CAM space at the moment */
+ /* - icxLuLut_init_clut_camclip() doesn't setup CAM rspl for vector clip. */
if (crv != 0 && p->camclip && p->nearclip) {
co cpp; /* Alternate CAM space solution */
double bf; /* Blend factor */
@@ -1375,7 +1434,6 @@ double *in /* Function input values to invert (== clut output' values) */
for (f = 0; f < fdi; f++) /* Transfer CAM targ */
cpp.v[f] = tin[f];
- cpp.v[0] *= CCJSCALE;
/* Make sure that the auxiliar value is initialized, */
/* even though it shouldn't have any effect, since should clipp. */
@@ -1391,7 +1449,7 @@ double *in /* Function input values to invert (== clut output' values) */
flags | xflags | RSPL_WILLCLIP, /* Combine all the flags + clip ?? */
1, /* Maximum solutions to return */
p->auxm, /* Auxiliary input chanel mask */
- cdir, /* Clip vector direction and length */
+ cdir, /* Clip vector direction/ LCh weighting */
&cpp); /* Input target and output solutions */
} else {
@@ -1400,7 +1458,7 @@ double *in /* Function input values to invert (== clut output' values) */
flags | RSPL_WILLCLIP, /* Because we know it will clip ?? */
1, /* Maximum solutions to return */
NULL, /* No auxiliary input targets */
- cdir, /* Clip vector direction and length */
+ cdir, /* Clip vector direction/LCh weighting */
&cpp); /* Input target and output solutions */
}
@@ -1411,7 +1469,6 @@ double *in /* Function input values to invert (== clut output' values) */
}
/* Compute the CAM clip distances */
- cpp.v[0] /= CCJSCALE;
for (cdist = 0.0, f = 0; f < fdi; f++) {
double tt;
tt = cpp.v[f] - tin[f];
@@ -1453,12 +1510,59 @@ double *in /* Function input values to invert (== clut output' values) */
/* Not CAM clip case */
} else {
+ /* If vector clip, replay with simple vector */
+ if (nsoln == 0 && p->nearclip == 0) {
+//printf("~1 !!!!!!!!!!!!!! Vector clip failed - trying safe replay !!!!!!!!!!!!!!!1\n");
+
+ // Reset pp[0] target values and auiliaries
+ for (e = 0; e < p->clutTable->di; e++)
+ pp[0].p[e] = upp.p[e];
+ for (f = 0; f < fdi; f++)
+ pp[0].v[f] = upp.v[f];
+
+ /* Compute safe clip vector */
+ cdir = icxClipVector(&p->clip, in, cdirv, 1);
+
+ /* Replay the lookup using safer vector */
+ nsoln = p->clutTable->rev_interp(
+ p->clutTable, /* rspl object */
+ uflags,
+ MAX_INVSOLN, /* Maxumum solutions to return */
+ NULL, /* No auxiliary input targets */
+ cdir, /* Clip vector direction and length */
+ pp); /* Input target and output solutions */
+ /* returned solutions in pp[0..retval-1].p[] */
+ nsoln &= RSPL_NOSOLNS; /* Get number of solutions */
+
+ if (nsoln == 0) { /* Hmm */
+//printf("~1 !!!!!!!!!!!!!! Vector clip failed again - using nn replay !!!!!!!!!!!!!!!1\n");
+ // Reset pp[0] target values and auiliaries
+ for (e = 0; e < p->clutTable->di; e++)
+ pp[0].p[e] = upp.p[e];
+ for (f = 0; f < fdi; f++)
+ pp[0].v[f] = upp.v[f];
+
+ /* Replay the lookup as a nearclip, to guarantee a solution */
+ nsoln = p->clutTable->rev_interp(
+ p->clutTable, /* rspl object */
+ RSPL_NEARCLIP | RSPL_NONNSETUP,
+ MAX_INVSOLN, /* Maxumum solutions to return */
+ NULL, /* No auxiliary input targets */
+ NULL, /* No LCh weighting */
+ pp); /* Input target and output solutions */
+ /* returned solutions in pp[0..retval-1].p[] */
+ nsoln &= RSPL_NOSOLNS; /* Get number of solutions */
+ }
+ }
+
if (nsoln == 1) { /* Exactly one solution */
i = 0;
} else if (nsoln == 0) { /* Zero solutions. This is unexpected. */
double in_v[MXDO];
p->output(p, in_v, pp[0].v); /* Get ICC inverse input values */
p->out_abs(p, in_v, in_v);
+ if (p->nearclip == 0)
+ a1logd(g_log,0,"Clip dst %f %f %f\n",pp[0].v[0]+cdir[0], pp[0].v[1]+cdir[1], pp[0].v[2]+cdir[2]);
error("Unexpected failure to find reverse solution for input to output table for value %f %f %f (ICC input %f %f %f)",pp[0].v[0],pp[0].v[1],pp[0].v[2], in_v[0], in_v[1], in_v[2]);
return 2;
} else { /* Multiple solutions */
@@ -1625,21 +1729,22 @@ int icxLuLut_inv_input(icxLuLut *p, double *out, double *in) {
int i,j;
int nsoln; /* Number of solutions found */
co pp[MAX_INVSOLN]; /* Room for all the solutions found */
- double cdir;
+// double cdir;
DBR(("inv_input got DEV' %f %f %f %f\n",in[0],in[1],in[2],in[3]))
for (i = 0; i < p->inputChan; i++) {
pp[0].p[0] = p->inputClipc[i];
pp[0].v[0] = in[i];
- cdir = p->inputClipc[i] - in[i]; /* Clip towards output range */
+// cdir = p->inputClipc[i] - in[i]; /* Clip towards output range */
nsoln = p->inputTable[i]->rev_interp (
p->inputTable[i], /* this */
RSPL_NEARCLIP, /* Clip to nearest (faster than vector) */
MAX_INVSOLN, /* Maximum number of solutions allowed for */
NULL, /* No auxiliary input targets */
- &cdir, /* Clip vector direction and length */
+ NULL, /* No LCH weight because this is 1D */
+// &cdir, /* Clip vector direction and length */
pp); /* Input and output values */
if (nsoln & RSPL_DIDCLIP)
@@ -1837,6 +1942,7 @@ icxLuBase *pp
/* - - - - - - - - - - - - - - - - - - - - - - - - - - */
static gamut *icxLuLutGamut(icxLuBase *plu, double detail);
+static icxCuspMap *icxLuLutCuspMap(icxLuBase *plu, int res);
/* Do the basic icxLuLut creation and initialisation */
static icxLuLut *
@@ -1860,6 +1966,7 @@ alloc_icxLuLut(
p->get_ranges = icxLu_get_ranges;
p->efv_wh_bk_points = icxLuEfv_wh_bk_points;
p->get_gamut = icxLuLutGamut;
+ p->get_cuspmap = icxLuLutCuspMap;
p->fwd_relpcs_outpcs = icxLuLut_fwd_relpcs_outpcs;
p->bwd_outpcs_relpcs = icxLuLut_bwd_outpcs_relpcs;
p->nearclip = 0; /* Set flag defaults */
@@ -2077,14 +2184,8 @@ icxLuLut *p /* Object being initialised */
p->clip.nearclip = 1;
} else { /* Vector clip */
- /* !!!! NOTE NOTE NOTE !!!! */
- /* We should re-write this to avoid calling p->clutTable->rev_interp(), */
- /* since this sets up all the rev acceleration tables for two calls, */
- /* if this lut is being setup from scattered data, and never used */
- /* for rev lookup */
icColorSpaceSignature clutos = p->natos;
- fprintf(stderr, "!!!!! setup_clip_icxLuLut with vector clip - possibly unnecessary rev setup !!!!\n");
p->clip.nearclip = 0;
p->clip.LabLike = 0;
p->clip.fdi = p->clutTable->fdi;
@@ -2092,95 +2193,25 @@ icxLuLut *p /* Object being initialised */
switch(clutos) {
case icxSigJabData:
case icSigLabData: {
- co pp; /* Room for all the solutions found */
- int nsoln; /* Number of solutions found */
- double cdir[MXDO]; /* Clip vector direction and length */
p->clip.LabLike = 1;
- /* Find high clip target */
- for (i = 0; i < p->inputChan; i++)
- pp.p[i] = 0.0; /* Set aux values */
- pp.v[0] = 105.0; pp.v[1] = pp.v[2] = 0.0; /* PCS Target value */
- cdir[0] = cdir[1] = cdir[2] = 0.0; /* Clip Target */
-
- p->inv_output(p, pp.v, pp.v); /* Compensate for output curve */
- p->inv_output(p, cdir, cdir);
-
- cdir[0] -= pp.v[0]; /* Clip vector */
- cdir[1] -= pp.v[1];
- cdir[2] -= pp.v[2];
-
- /* PCS -> Device with clipping */
- nsoln = p->clutTable->rev_interp(
- p->clutTable, /* rspl object */
- 0, /* No hint flags - might be in gamut, might vector clip */
- 1, /* Maxumum solutions to return */
- p->auxm, /* Auxiliary input targets */
- cdir, /* Clip vector direction and length */
- &pp); /* Input target and output solutions */
- /* returned solutions in pp[0..retval-1].p[] */
- nsoln &= RSPL_NOSOLNS; /* Get number of solutions */
-
- if (nsoln != 1) {
- p->pp->errc = 2;
- sprintf(p->pp->err,"Failed to find high clip target for Lab space");
- return p->pp->errc;
- }
-
- p->clip.ocent[0] = pp.v[0] - 0.001; /* Got main target */
- p->clip.ocent[1] = pp.v[1];
- p->clip.ocent[2] = pp.v[2];
-
- /* Find low clip target */
- pp.v[0] = -5.0; pp.v[1] = pp.v[2] = 0.0; /* PCS Target value */
- cdir[0] = 100.0; cdir[1] = cdir[2] = 0.0; /* Clip Target */
-
- p->inv_output(p, pp.v, pp.v); /* Compensate for output curve */
- p->inv_output(p, cdir, cdir);
- cdir[0] -= pp.v[0]; /* Clip vector */
- cdir[1] -= pp.v[1];
- cdir[2] -= pp.v[2];
-
- /* PCS -> Device with clipping */
- nsoln = p->clutTable->rev_interp(
- p->clutTable, /* rspl object */
- RSPL_WILLCLIP, /* Since there was no locus, we expect to have to clip */
- 1, /* Maxumum solutions to return */
- NULL, /* No auxiliary input targets */
- cdir, /* Clip vector direction and length */
- &pp); /* Input target and output solutions */
- /* returned solutions in pp[0..retval-1].p[] */
- nsoln &= RSPL_NOSOLNS; /* Get number of solutions */
- if (nsoln != 1) {
- p->pp->errc = 2;
- sprintf(p->pp->err,"Failed to find low clip target for Lab space");
- return p->pp->errc;
- }
-
- p->clip.ocentv[0] = pp.v[0] + 0.001 - p->clip.ocent[0]; /* Raw line vector */
- p->clip.ocentv[1] = pp.v[1] - p->clip.ocent[1];
- p->clip.ocentv[2] = pp.v[2] - p->clip.ocent[2];
-
- /* Compute vectors length */
- for (p->clip.ocentl = 0.0, i = 0; i < 3; i++)
- p->clip.ocentl += p->clip.ocentv[i] * p->clip.ocentv[i];
- p->clip.ocentl = sqrt(p->clip.ocentl);
- if (p->clip.ocentl <= 1e-8)
- p->clip.ocentl = 0.0;
-
+ /* Create a CuspMap to point vectors towards */
+ /* (Don't make it too fine, or there will be dips) */
+ p->clip.cm = p->get_cuspmap((icxLuBase *)p, 30);
break;
}
case icSigXYZData:
// ~~~~~~1 need to add this.
+ warning("xlut.c: setup_clip_icxLuLut() icSigXYZData case not implemented!");
+ /* Fall through */
default:
/* Do a crude approximation, that may not work. */
p->clutTable->get_out_range(p->clutTable, tmin, tmax);
- for (i = 0; i < p->clutTable->fdi; i++) {
+ for (i = 0; i < p->clutTable->fdi; i++)
p->clip.ocent[i] = (tmin[i] + tmax[i])/2.0;
- }
- p->clip.ocentl = 0.0;
+
break;
}
}
@@ -2306,13 +2337,13 @@ fprintf(stderr,"~1 Internal optimised 4D separations not yet implemented!\n");
/* Init the CAM model if it will be used */
if (pcsor == icxSigJabData || p->camclip) {
- if (vc != NULL) /* One has been provided */
+ if (vc != NULL) /* One has been provided */
p->vc = *vc; /* Copy the structure */
else
xicc_enum_viewcond(xicp, &p->vc, -1, NULL, 0, NULL); /* Use a default */
p->cam = new_icxcam(cam_default);
p->cam->set_view(p->cam, p->vc.Ev, p->vc.Wxyz, p->vc.La, p->vc.Yb, p->vc.Lv,
- p->vc.Yf, p->vc.Yg, p->vc.Gxyz, XICC_USE_HK);
+ p->vc.Yf, p->vc.Yg, p->vc.Gxyz, XICC_USE_HK, p->vc.hkscale);
} else {
p->cam = NULL;
}
@@ -2330,8 +2361,8 @@ fprintf(stderr,"~1 Internal optimised 4D separations not yet implemented!\n");
p->pcs = pcsor;
if (xicp->pp->header->deviceClass == icSigAbstractClass) {
- p->ins = pcsor;
- p->outs = pcsor;
+ p->ins = pcsor;
+ p->outs = pcsor;
} else if (xicp->pp->header->deviceClass != icSigLinkClass) {
if (func == icmBwd || func == icmGamut || func == icmPreview)
@@ -2412,13 +2443,21 @@ fprintf(stderr,"~1 Internal optimised 4D separations not yet implemented!\n");
/* ------------------------------- */
{
int gres[MXDI];
+ int xflags = 0;
for (i = 0; i < p->inputChan; i++)
gres[i] = p->lut->clutPoints;
+#ifdef FASTREVSETUP_NON_CAM
+ /* Don't fill in nnrev array if we aren't going to use it */
+ if (p->camclip && p->nearclip)
+ xflags = RSPL_FASTREVSETUP;
+#endif
+
/* Create rspl based multi-dim table */
if ((p->clutTable = new_rspl((p->fastsetup ? RSPL_FASTREVSETUP : RSPL_NOFLAGS)
- | (flags & ICX_VERBOSE ? RSPL_VERBOSE : RSPL_NOFLAGS),
+ | (flags & ICX_VERBOSE ? RSPL_VERBOSE : RSPL_NOFLAGS)
+ | xflags,
p->inputChan, p->outputChan)) == NULL) {
p->pp->errc = 2;
sprintf(p->pp->err,"Creation of clut table rspl failed");
@@ -2438,6 +2477,17 @@ fprintf(stderr,"~1 Internal optimised 4D separations not yet implemented!\n");
}
+#ifdef USELCHWEIGHT
+ /* If we are not doing camclip, but our output is an Lab like space, */
+ /* then apply lchw weighting anyway. */
+ if (!p->camclip && (p->outs == icSigLabData || p->outs == icxSigJabData)) {
+ double lchw[MXRO] = { JCCWEIGHT, CCCWEIGHT, HCCWEIGHT };
+
+ /* Set the Nearest Neighbor clipping Weighting */
+ p->clutTable->rev_set_lchw(p->clutTable, lchw);
+ }
+#endif /* USELCHWEIGHT */
+
/* clut clipping is setup separately */
}
@@ -2500,22 +2550,24 @@ icxLuLut_clut_camclip_func(
luluto->output(luluto, out, out);
luluto->out_abs(luluto, out, out);
p->cam->XYZ_to_cam(p->cam, out, out);
- out[0] *= CCJSCALE;
}
/* Initialise the additional CAM space clut rspl, used to compute */
/* reverse lookup CAM clipping results when the camclip flag is set. */
-/* We scale J by CCJSCALE, to give a more L* preserving clip direction */
-/* return error code. */
+/* We weight the CAM nn clipping, to give a more L* and H* preserving clip direction. */
+/* Return error code. */
+/* (We are assuming nearest clipping - we aren't setting up properly for */
+/* vector clipping) */
static int
icxLuLut_init_clut_camclip(
icxLuLut *p) {
int e, gres[MXDI];
+ double lchw[MXRO] = { JCCWEIGHT, CCCWEIGHT, HCCWEIGHT };
/* Setup so clut contains transform to CAM Jab */
/* (camclip is only used in fwd or invfwd direction lookup) */
double cmin[3], cmax[3];
- cmin[0] = 0.0; cmax[0] = CCJSCALE * 100.0; /* Nominal Jab output ranges */
+ cmin[0] = 0.0; cmax[0] = 100.0; /* Nominal Jab output ranges */
cmin[1] = -128.0; cmax[1] = 128.0;
cmin[2] = -128.0; cmax[2] = 128.0;
@@ -2536,6 +2588,11 @@ icxLuLut *p) {
return p->pp->errc;
}
+#ifdef USELCHWEIGHT
+ /* Set the Nearest Neighbor clipping Weighting */
+ p->cclutTable->rev_set_lchw(p->cclutTable, lchw);
+#endif /* USELCHWEIGHT */
+
for (e = 0; e < p->inputChan; e++)
gres[e] = p->lut->clutPoints;
@@ -3899,7 +3956,7 @@ int quality /* Quality metric, 0..3 */
xicc_enum_viewcond(xicp, &p->vc, -1, NULL, 0, NULL); /* Use a default */
p->cam = new_icxcam(cam_default);
p->cam->set_view(p->cam, p->vc.Ev, p->vc.Wxyz, p->vc.La, p->vc.Yb, p->vc.Lv,
- p->vc.Yf, p->vc.Yg, p->vc.Gxyz, XICC_USE_HK);
+ p->vc.Yf, p->vc.Yg, p->vc.Gxyz, XICC_USE_HK, p->vc.hkscale);
if (flags & ICX_VERBOSE)
printf("Done A to B table creation\n");
@@ -3912,7 +3969,7 @@ int quality /* Quality metric, 0..3 */
/* ========================================================== */
-/* Context for creating gamut boundary points fro, xicc */
+/* Context for creating gamut boundary points from, xicc */
typedef struct {
gamut *g; /* Gamut being created */
icxLuLut *x; /* xLut we are working from */
@@ -4096,9 +4153,9 @@ double detail /* gamut detail level, 0.0 = def */
}
/* Scan only device surface */
- for (m1 = 0; m1 < inn; m1++) { /* Choose first coord to scan */
+ for (m1 = 0; m1 < (inn-1); m1++) { /* Choose first coord to scan */
if (co[m1] != 0)
- continue; /* Not at lower corner */
+ continue; /* Not at lower corner */
for (m2 = m1 + 1; m2 < inn; m2++) { /* Choose second coord to scan */
int x, y;
@@ -4216,7 +4273,7 @@ double detail /* gamut detail level, 0.0 = def */
default:
break;
}
- if ((cx.flu = p->get_luobj(p, ICX_CLIP_NEAREST , icmFwd, intent, pcs, icmLuOrdNorm,
+ if ((cx.flu = p->get_luobj(p, ICX_CLIP_NEAREST, icmFwd, intent, pcs, icmLuOrdNorm,
&plu->vc, NULL)) == NULL) {
return NULL; /* oops */
}
@@ -4314,6 +4371,304 @@ double detail /* gamut detail level, 0.0 = def */
return gam;
}
+/* ========================================================== */
+/* Cusp Map finding code. */
+/* ========================================================== */
+
+/* Context for creating Cusp Map points from, xicc */
+typedef struct {
+ icxCuspMap *cm; /* Cusp Map being created */
+ icxLuLut *x; /* xLut we are working from */
+ icxLuBase *flu; /* Forward xlookup */
+ double in[MAX_CHAN]; /* Device input values */
+} lutcuspmapctx;
+
+/* Function to pass to rspl to create cusp map from */
+/* forward xLut transform grid points. */
+static void
+lutfwdcuspmap_func(
+ void *pp, /* lutcuspmapctx structure */
+ double *out, /* output' value at clut grid point (ie. PCS' value) */
+ double *in /* input' value at clut grid point (ie. device' value) */
+) {
+ lutcuspmapctx *p = (lutcuspmapctx *)pp;
+ double pcso[3]; /* PCS output value */
+
+ /* Figure if we are over the ink limit. */
+ if ( (p->x->ink.tlimit >= 0.0 || p->x->ink.klimit >= 0.0)
+ && icxLimitD(p->x, in) > 0.0) {
+ int i;
+ double sf;
+
+ /* We are, so use the bracket search to discover a scale */
+ /* for the clut input' value that will put us on the ink limit. */
+
+ for (i = 0; i < p->x->inputChan; i++)
+ p->in[i] = in[i];
+
+ if (zbrent(&sf, 0.0, 1.0, 1e-4, icxLimitFind, pp) != 0) {
+ return; /* Give up */
+ }
+
+ /* Compute ink limit value */
+ for (i = 0; i < p->x->inputChan; i++)
+ p->in[i] = sf * in[i];
+
+ /* Compute the clut output for this clut input */
+ p->x->clut(p->x, pcso, p->in);
+ p->x->output(p->x, pcso, pcso);
+ p->x->out_abs(p->x, pcso, pcso);
+ } else { /* No ink limiting */
+ /* Convert the clut PCS' values to PCS output values */
+ p->x->output(p->x, pcso, out);
+ p->x->out_abs(p->x, pcso, pcso);
+ }
+
+ /* Expand the usp map surface with this point */
+ p->cm->expand(p->cm, pcso);
+
+ /* Leave out[] unchanged */
+}
+
+/* Expand cusp map with given point */
+/* (We use a segmentned maxima approach) */
+static void cuspmap_expand(icxCuspMap *p, double lab[3]) {
+ double h, C;
+ int ix;
+
+ /* Hue angle 0.0 .. 1.0 */
+ h = (0.5/3.14159265359) * atan2(lab[2], lab[1]);
+ h = (h < 0.0) ? h + 1.0 : h;
+
+ /* Slot index 0..res-1 */
+ ix = (int)floor(h * p->res + 0.5);
+ if (ix >= p->res)
+ ix -= p->res;
+
+ /* Chromanance */
+ C = sqrt(lab[1] * lab[1] + lab[2] * lab[2]);
+
+ /* New point at this angle with largest chromanance */
+ if (C > p->C[ix]) {
+ p->C[ix] = C;
+ p->L[ix] = lab[0];
+ }
+
+ /* Tracl min * max L */
+ if (lab[0] > p->Lmax[0])
+ icmCpy3(p->Lmax, lab);
+ if (lab[0] < p->Lmin[0])
+ icmCpy3(p->Lmin, lab);
+}
+
+/* Interpolate over any gaps in map */
+static void cuspmap_complete(icxCuspMap *p) {
+ int i, j, k;
+
+// printf("cuspmap list before fixups:\n");
+// for (i = 0; i < p->res; i++) {
+// printf(" %d: C = %f, L = %f\n",i, p->C[i], p->L[i]);
+// }
+
+ /* First check if there are any entries at all */
+ for (i = 0; i < p->res; i++) {
+ if (p->C[i] >= 0.0)
+ break;
+ }
+
+ if (i >= p->res) /* Nothing there - give up */
+ return;
+
+ /* See if there are any gaps */
+ j = -1;
+ for (i = 0; i < p->res; i++) {
+//printf("~1 checking %d\n",i);
+ if (p->C[i] >= 0.0) {
+ j = i; /* Last valid slot */
+//printf("~1 last valid %d\n",j);
+
+ /* Got a gap */
+ } else {
+ int ii = i;
+//printf("~1 found gap at %d\n",i);
+ if (j < 0) { /* No previous */
+//printf("~1 no previous\n");
+ for (j = p->res-1; j >= 0; j--) {
+ if (p->C[j] >= 0.0)
+ break;
+ }
+ }
+//printf("~1 got previous %d\n",j);
+ /* Find next, even if it's behind us */
+ for (k = i+1; k != i; k = (k+1) % p->res) {
+ if (p->C[k] >= 0.0)
+ break;
+ }
+//printf("~1 got next %d\n",k);
+
+ /* Interpolate between them */
+ for (; i != k; i = (i+1) % p->res) {
+ double prop;
+ int bb, tt;
+ bb = k > i ? k - i : k + p->res - i;
+ tt = k > j ? k - j : k + p->res - j;
+ prop = (double)bb/(double)tt;
+ p->C[i] = prop * p->C[j] + (1.0 - prop) * p->C[k];
+ p->L[i] = prop * p->L[j] + (1.0 - prop) * p->L[k];
+//printf("~1 interpolating %d with weight %f from %d and %d\n",i,prop,j,k);
+ }
+ if (k > ii) {
+ i = k; /* Continue from next valid */
+//printf("~1 continuing from %d\n",i);
+ } else {
+//printf("~1 we've looped back\n");
+ break; /* We've looped back */
+ }
+ }
+ }
+
+// printf("cuspmap list after fixups:\n");
+// for (i = 0; i < p->res; i++) {
+// printf(" %d: C = %f, L = %f\n",i, p->C[i], p->L[i]);
+// }
+}
+
+/* Return the corresponding cusp location, given the source point */
+static void cuspmap_getCusp(icxCuspMap *p,double cuspLCh[3], double srcLab[3]) {
+ double h, C;
+ int ix, ix0, ix1;
+
+ /* Hue angle 0.0 .. 1.0 */
+ h = (0.5/3.14159265359) * atan2(srcLab[2], srcLab[1]);
+ h = (h < 0.0) ? h + 1.0 : h;
+
+//printf("~1 getCusp in %f %f %f, h = %f\n", srcLab[0], srcLab[1], srcLab[2],h);
+
+ /* Slot index 0..res-1 */
+ ix = (int)floor(h * p->res + 0.5);
+ if (ix >= p->res)
+ ix -= p->res;
+
+ /* Indexes each side */
+ ix0 = ix > 0 ? ix-1 : p->res-1;
+ ix1 = ix < (p->res-1) ? ix+1 : 0;
+
+ cuspLCh[0] = p->L[ix];
+ cuspLCh[1] = p->C[ix];
+ if (cuspLCh[1] > p->C[ix0]) /* Be conservative with C value */
+ cuspLCh[1] = p->C[ix0];
+ if (cuspLCh[1] > p->C[ix1])
+ cuspLCh[1] = p->C[ix1];
+ cuspLCh[2] = 360.0 * h;
+
+//printf("~1 index %d, L %f C %f h %f\n", ix, cuspLCh[0], cuspLCh[1], cuspLCh[2]);
+}
+
+/* We're done with CuspMap */
+static void cuspmap_del(icxCuspMap *p) {
+ if (p != NULL) {
+ if (p->L != NULL)
+ free(p->L);
+ if (p->C != NULL)
+ free(p->C);
+ free(p);
+ }
+}
+
+/* Given an xicc lookup object, return an icxCuspMap object. */
+/* Note that the PCS must be Lab or Jab. */
+/* An icxLuLut type must be icmFwd, and the ink limit (if supplied) */
+/* will be applied. */
+/* Return NULL on error, check errc+err for reason */
+static icxCuspMap *icxLuLutCuspMap(
+icxLuBase *plu, /* this */
+int res /* Hue resolution */
+) {
+ xicc *p = plu->pp; /* parent xicc */
+ icxLuLut *luluto = (icxLuLut *)plu; /* Lookup xLut type object */
+ icColorSpaceSignature ins, pcs, outs;
+ icmLookupFunc func;
+ icRenderingIntent intent;
+ int inn, outn;
+ lutcuspmapctx cx;
+ icxCuspMap *cm;
+ int i;
+
+ /* get some details */
+ plu->spaces(plu, &ins, &inn, &outs, &outn, NULL, &intent, &func, &pcs);
+
+ if (func != icmFwd) {
+ p->errc = 1;
+ sprintf(p->err,"Creating CuspMap for anything other than Device -> PCS is not supported.");
+ return NULL;
+ }
+
+ if (pcs != icSigLabData && pcs != icxSigJabData) {
+ p->errc = 1;
+ sprintf(p->err,"Creating CuspMap PCS of other than Lab or Jab is not supported.");
+ return NULL;
+ }
+
+ cx.cm = cm = (icxCuspMap *) calloc(1, sizeof(icxCuspMap));
+ cx.x = luluto;
+
+ if (cx.cm == NULL) {
+ p->errc = 2;
+ sprintf(p->err,"Malloc of icxCuspMap failed");
+ return NULL;
+ }
+
+ if ((cm->L = (double *)malloc(sizeof(double) * res)) == NULL) {
+ free(cm);
+ p->errc = 2;
+ sprintf(p->err,"Malloc of icxCuspMap failed");
+ return NULL;
+ }
+
+ if ((cm->C = (double *)malloc(sizeof(double) * res)) == NULL) {
+ free(cm->L);
+ free(cm);
+ p->errc = 2;
+ sprintf(p->err,"Malloc of icxCuspMap failed");
+ return NULL;
+ }
+
+ cm->res = res;
+ cm->expand = cuspmap_expand;
+ cm->getCusp = cuspmap_getCusp;
+ cm->del = cuspmap_del;
+
+ for (i = 0; i < res; i++) {
+ cm->L[i] = -1.0;
+ cm->C[i] = -1.0;
+ }
+ cm->Lmax[0] = -1.0;
+ cm->Lmin[0] = 101.0;
+
+ /* Scan through grid, expanding the CuspMap. */
+ luluto->clutTable->scan_rspl(
+ luluto->clutTable, /* this */
+ RSPL_NOFLAGS, /* Combination of flags */
+ (void *)&cx, /* Opaque function context */
+ lutfwdcuspmap_func /* Function to set from */
+ );
+
+ cm->Lmax[0] -= 0.1; /* Make sure tips intersect */
+ cm->Lmin[0] += 0.1;
+
+ for (i = 0; i < res; i++) {
+ if (cm->L[i] > cm->Lmax[0])
+ cm->L[i] = cm->Lmax[0];
+ if (cm->L[i] < cm->Lmin[0])
+ cm->L[i] = cm->Lmin[0];
+ }
+
+ /* Fill in any gaps */
+ cuspmap_complete(cm);
+
+ return cm;
+}
+
/* ----------------------------------------------- */
#ifdef DEBUG
#undef DEBUG /* Limit extent to this file */