1
2 /*
3 * Argyll Color Correction System
4 * Input device profile creator.
5 *
6 * Author: Graeme W. Gill
7 * Date: 11/10/00
8 *
9 * Copyright 2000 - 2011 Graeme W. Gill
10 * All rights reserved.
11 *
12 * This material is licenced under the GNU AFFERO GENERAL PUBLIC LICENSE Version 3 :-
13 * see the License.txt file for licencing details.
14 */
15
16 /*
17 * This program takes in the scattered test chart
18 * points, and interpolates them into a gridded
19 * forward ICC device profile.
20 * It also creates (at the moment) a limited backward
21 * profile based on the forward grid.
22 *
23 */
24
25 /*
26 * TTBD:
27 * Need to make this more of a library:
28 * ** By default limit matrix primaries to have +ve XYZ
29 * Add flag to override this.
30 * Fix error handling
31 * fix verbose output
32 * hand icc object back rather than writing file ?
33 */
34
35 #undef DEBUG
36
37 #define verbo stdout
38
39 #include <stdio.h>
40 #include "counters.h"
41 #include "numlib.h"
42 #include "icc.h"
43 #include "cgats.h"
44 #include "xicc.h"
45 #include "rspl.h"
46 #include "prof.h"
47
48 #define DOB2A /* Create B2A table as well (not implemented) */
49 #define NO_B2A_PCS_CURVES /* PCS curves seem to make B2A less accurate. Why ? */
50 #define USE_CAM_CLIP_OPT /* Clip out of gamut in CAM space rather than XYZ or L*a*b* */
51 #undef USE_EXTRA_FITTING /* Turn on data point error compensation */
52 #define USE_2ASS_SMOOTHING /* Turn on Gaussian smoothing */
53 #undef WARN_CLUT_CLIPPING /* Print warning if setting clut clips */
54 #define EXTRAP_MAXPNTS 10 /* Maximum number of extra extrapolated points per direction */
55 #define EXTRAP_WEIGHT 1.0 /* Extra extrapolated point weighting */
56
57 /*
58 Basic algorithm outline:
59
60 Scanner:
61
62 Figure out the input curves to give
63 the flattest grid.
64
65 Figure out the grid values.
66
67 Use them to generate the A2B table.
68
69 Do all the calculations in Lab space,
70 but represent the profile in XYZ space, so that
71 the white/black point normalisation doesn't cause
72 the clut values to be clipped.
73
74 This leads to a poorer accuracy as an XYZ profile,
75 but can then be compensated for using the ICX_MERGE_CLUT flag
76 together with a PCS override.
77
78 Note we're hard coded as RGB device space, so we're not coping
79 with grey scale or CMY.
80 */
81
82 #ifdef DEBUG
83 #undef DBG
84 #define DBG(xxx) printf xxx ;
85 #else
86 #undef DBG
87 #define DBG(xxx)
88 #endif
89
90 /* ---------------------------------------- */
91 #ifdef DOB2A
92
93 /* structure to support output icc B2A Lut initialisation calbacks */
94 /* Note that we don't cope with a LUT matrix - assume it's unity. */
95
96 typedef struct {
97 int verb;
98 int total, count, last; /* Progress count information */
99 int noPCScurves; /* Flag set if we don't want PCS curves */
100 icColorSpaceSignature pcsspace; /* The PCS colorspace */
101 icColorSpaceSignature devspace; /* The device colorspace */
102 icxLuLut *x; /* A2B icxLuLut we are inverting in std PCS */
103
104 double swxyz[3]; /* Source white point in XYZ */
105
106 int wantLab; /* 0 if want is XYZ PCS, 1 want is Lab PCS */
107 } in_b2a_callback;
108
109
110 /* --------------------------------------------------------- */
111
112 /* Extra non-linearity applied to BtoA XYZ PCS */
113 /* This distributes the LUT indexes more evenly in */
114 /* perceptual space, greatly improving the B2A accuracy of XYZ LUT */
xyzcurve(double * out,double * in)115 static void xyzcurve(double *out, double *in) {
116 int i;
117 double sc = 65535.0/32768.0;
118
119 /* Use an L* like curve, scaled to the maximum XYZ valu */
120 out[0] = in[0]/sc;
121 out[1] = in[1]/sc;
122 out[2] = in[2]/sc;
123 for (i = 0; i < 3; i++) {
124 if (out[i] > 0.08)
125 out[i] = pow((out[i] + 0.16)/1.16, 3.0);
126 else
127 out[i] = out[i]/9.032962896;
128 }
129 out[0] = out[0] * sc;
130 out[1] = out[1] * sc;
131 out[2] = out[2] * sc;
132 }
133
invxyzcurve(double * out,double * in)134 static void invxyzcurve(double *out, double *in) {
135 int i;
136 double sc = 65535.0/32768.0;
137
138 out[0] = in[0]/sc;
139 out[1] = in[1]/sc;
140 out[2] = in[2]/sc;
141 for (i = 0; i < 3; i++) {
142 if (out[i] > 0.008856451586)
143 out[i] = 1.16 * pow(out[i],1.0/3.0) - 0.16;
144 else
145 out[i] = 9.032962896 * out[i];
146 }
147 out[0] = out[0] * sc;
148 out[1] = out[1] * sc;
149 out[2] = out[2] * sc;
150 }
151
152 /* --------------------------------------------------------- */
153 /* NOTE :- the assumption that each stage of the BtoA is a mirror */
154 /* of the AtoB makes for inflexibility. */
155 /* Perhaps it would be better to remove this asumption from the */
156 /* in_b2a_clut processing ? */
157 /* To do this we then need inv_in_b2a_input(), and */
158 /* inv_in_b2a_output(), and we need to clearly distinguish between */
159 /* AtoB PCS' & DEV', and BtoA PCS' & DEV', since they are not */
160 /* necessarily the same... */
161
162
163 /* B2A Input table is the inverse of the AtoB output table */
164 /* Input PCS output PCS'' */
in_b2a_input(void * cntx,double out[3],double in[3])165 void in_b2a_input(void *cntx, double out[3], double in[3]) {
166 in_b2a_callback *p = (in_b2a_callback *)cntx;
167
168 DBG(("out_b2a_input got PCS %f %f %f\n",in[0],in[1],in[2]))
169
170 /* PCS to PCS' */
171 if (p->noPCScurves) {
172 out[0] = in[0];
173 out[1] = in[1];
174 out[2] = in[2];
175 } else {
176 if (p->x->inv_output(p->x, out, in) > 1)
177 error("%d, %s",p->x->pp->errc,p->x->pp->err);
178 }
179 /* PCS' to PCS'' */
180 if (p->pcsspace == icSigXYZData) /* Apply XYZ non-linearity curve */
181 invxyzcurve(out, out);
182
183 DBG(("in_b2a_input returning PCS'' %f %f %f\n",out[0],out[1],out[2]))
184 }
185
186 /* clut - multitable */
187 /* Input PCS' output Dev' */
188 /* We're applying any abstract profile after gamut mapping, */
189 /* on the assumption is is primarily being used to "correct" the */
190 /* output device. Ideally the gamut mapping should take the change */
191 /* the abstract profile has on the output device into account, but */
192 /* currently we're not doing this.. */
in_b2a_clut(void * cntx,double * out,double in[3])193 void in_b2a_clut(void *cntx, double *out, double in[3]) {
194 in_b2a_callback *p = (in_b2a_callback *)cntx;
195 double in1[3];
196
197 in1[0] = in[0]; /* in[] may be aliased with out[] */
198 in1[1] = in[1]; /* so take a copy. */
199 in1[2] = in[2];
200
201 DBG(("in_b2a_clut got PCS' %f %f %f\n",in[0],in[1],in[2]))
202
203 if (p->pcsspace == icSigXYZData) /* Undo effects of extra XYZ non-linearity curve */
204 xyzcurve(in1, in1);
205
206 if (p->noPCScurves) { /* We were given PCS or have converted to PCS */
207
208 /* PCS to PCS' */
209 if (p->x->inv_output(p->x, in1, in1) > 1)
210 error("%d, %s",p->x->pp->errc,p->x->pp->err);
211
212 DBG(("convert to PCS' got %f %f %f\n",in1[0],in1[1],in1[2]))
213 }
214
215 /* Invert AtoB clut (PCS' to Dev') Colorimetric */
216 /* to producte the colorimetric tables output. */
217 if (p->x->inv_clut(p->x, out, in1) > 1)
218 error("%d, %s",p->x->pp->errc,p->x->pp->err);
219
220 DBG(("convert PCS' to DEV' got %f %f %f %f\n",out[0],out[1],out[2],out[3]))
221 DBG(("in_b2a_clut returning DEV' %f %f %f\n",out[0],out[1],out[2]))
222
223 if (p->verb) { /* Output percent intervals */
224 int pc;
225 p->count++;
226 pc = (int)(p->count * 100.0/p->total + 0.5);
227 if (pc != p->last) {
228 printf("%c%2d%%",cr_char,pc); fflush(stdout);
229 p->last = pc;
230 }
231 }
232 }
233
234 /* Output table is the inverse of the AtoB input table */
235 /* Input Dev' output Dev */
in_b2a_output(void * cntx,double out[4],double in[4])236 void in_b2a_output(void *cntx, double out[4], double in[4]) {
237 in_b2a_callback *p = (in_b2a_callback *)cntx;
238
239 DBG(("in_b2a_output got DEV' %f %f %f\n",in[0],in[1],in[2]))
240
241 if (p->x->inv_input(p->x, out, in) > 1)
242 error("%d, %s",p->x->pp->errc,p->x->pp->err);
243
244 DBG(("in_b2a_output returning DEV %f %f %f\n",out[0],out[1],out[2]))
245 }
246
247 #endif /* DOB2A */
248 /* ---------------------------------------- */
249
250 /* Make an input device profile, where we create an A2B lut */
251 /* directly from the scattered input data. */
252 void
make_input_icc(prof_atype ptype,icmICCVersion iccver,int verb,int iquality,int oquality,int noisluts,int noipluts,int nooluts,int nocied,int verify,int autowpsc,int clipovwp,double wpscale,int dob2a,int extrap,int clipprims,char * in_name,char * file_name,cgats * icg,int emis,int spec,icxIllumeType illum,xspect * cust_illum,icxObserverType observ,double smooth,double avgdev,profxinf * xpi)253 make_input_icc(
254 prof_atype ptype, /* Profile algorithm type */
255 icmICCVersion iccver, /* ICC profile version to create */
256 int verb,
257 int iquality, /* A2B table quality, 0..3 */
258 int oquality, /* B2A table quality, 0..3 */
259 int noisluts, /* nz to supress creation of input (Device) shaper luts */
260 int noipluts, /* nz to supress creation of input (Device) position luts */
261 int nooluts, /* nz to supress creation of output (PCS) shaper luts */
262 int nocied, /* nz to supress inclusion of .ti3 data in profile */
263 int verify,
264 int autowpsc, /* nz for Auto scale the WP to prevent clipping above WP patch */
265 int clipovwp, /* nz for Clip cLUT values above WP */
266 double wpscale, /* >= 0.0 for media white point scale factor */
267 int dob2a, /* nz to create a B2A table as well */
268 int extrap, /* nz to create extra cLUT interpolation points */
269 int clipprims, /* Clip white, black and primaries */
270 char *in_name, /* input .ti3 file name */
271 char *file_name, /* output icc name */
272 cgats *icg, /* input cgats structure */
273 int emis, /* emissive reference data */
274 int spec, /* Use spectral data flag */
275 icxIllumeType illum, /* Spectral illuminant */
276 xspect *cust_illum, /* Possible custom illumination */
277 icxObserverType observ, /* Spectral observer */
278 double smooth, /* RSPL smoothing factor, -ve if raw */
279 double avgdev, /* reading Average Deviation as a proportion of the input range */
280 profxinf *xpi /* Optional Profile creation extra data */
281 ) {
282 icmFile *wr_fp;
283 icc *wr_icco;
284 int npat; /* Number of patches */
285 int npxpat = 0; /* Number of possible extrap extrapolation patches */
286 int nxpat = 0; /* Number of extrap extrapolation patches */
287 cow *tpat; /* Patch input values */
288 int i, rv = 0;
289 int isLab = 0; /* 0 if input is XYZ, 1 if input is Lab */
290 int wantLab = 0; /* 0 if want is XYZ, 1 want is Lab. */
291 /* Values will be wantLab after reading */
292 int isLut = 0; /* 0 if shaper+ matrix, 1 if lut type */
293 int isShTRC = 0; /* 0 if separate gamma/shaper TRC, 1 if shared */
294
295 if (ptype == prof_clutLab) { /* Lab lut */
296 wantLab = 1;
297 isLut = 1;
298 } else if (ptype == prof_clutXYZ) { /* XYZ lut */
299 wantLab = 0;
300 isLut = 1;
301 } else {
302 wantLab = 0; /* gamma/shaper + matrix profile must be XYZ */
303 isLut = 0;
304 extrap = 0;
305
306 if (ptype == prof_gam1mat
307 || ptype == prof_sha1mat
308 || ptype == prof_matonly) {
309 isShTRC = 1; /* Single curve */
310 }
311 }
312
313 /* Open up the file for writing */
314 if ((wr_fp = new_icmFileStd_name(file_name,"w")) == NULL)
315 error ("Write: Can't open file '%s'",file_name);
316
317 if ((wr_icco = new_icc()) == NULL)
318 error ("Write: Creation of ICC object failed");
319
320 /* Add all the tags required */
321
322 /* The header: */
323 {
324 icmHeader *wh = wr_icco->header;
325
326 /* Values that must be set before writing */
327 wh->deviceClass = icSigInputClass;
328 wh->colorSpace = icSigRgbData; /* It's an RGB profile */
329 if (wantLab)
330 wh->pcs = icSigLabData;
331 else
332 wh->pcs = icSigXYZData;
333
334 if (xpi->default_ri != icMaxEnumIntent)
335 wh->renderingIntent = xpi->default_ri;
336 else
337 wh->renderingIntent = icRelativeColorimetric;
338
339 /* Values that should be set before writing */
340 if (xpi != NULL && xpi->manufacturer != 0L)
341 wh->manufacturer = xpi->manufacturer;
342 else
343 wh->manufacturer = icmSigUnknownType;
344
345 if (xpi != NULL && xpi->model != 0L)
346 wh->model = xpi->model;
347 else
348 wh->model = icmSigUnknownType;
349
350 /* Values that may be set before writing */
351 if (xpi != NULL && xpi->creator != 0L)
352 wh->creator = xpi->creator;
353 #ifdef NT
354 wh->platform = icSigMicrosoft;
355 #endif
356 #ifdef __APPLE__
357 wh->platform = icSigMacintosh;
358 #endif
359 #if defined(UNIX) && !defined(__APPLE__)
360 wh->platform = icmSig_nix;
361 #endif
362
363 if (xpi != NULL && xpi->transparency)
364 wh->attributes.l |= icTransparency;
365 if (xpi != NULL && xpi->matte)
366 wh->attributes.l |= icMatte;
367 if (xpi != NULL && xpi->negative)
368 wh->attributes.l |= icNegative;
369 if (xpi != NULL && xpi->blackandwhite)
370 wh->attributes.l |= icBlackAndWhite;
371 }
372 /* Profile Description Tag: */
373 {
374 icmTextDescription *wo;
375 char *dst; /* description */
376
377 if (xpi != NULL && xpi->profDesc != NULL)
378 dst = xpi->profDesc;
379 else {
380 dst = "This is a Lut style RGB - XYZ Input Profile";
381 }
382
383 if ((wo = (icmTextDescription *)wr_icco->add_tag(
384 wr_icco, icSigProfileDescriptionTag, icSigTextDescriptionType)) == NULL)
385 error("add_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
386
387 wo->size = strlen(dst)+1; /* Allocated and used size of desc, inc null */
388 wo->allocate((icmBase *)wo);/* Allocate space */
389 strcpy(wo->desc, dst); /* Copy the string in */
390 }
391 /* Copyright Tag: */
392 {
393 icmText *wo;
394 char *crt;
395
396 if (xpi != NULL && xpi->copyright != NULL)
397 crt = xpi->copyright;
398 else
399 crt = "Copyright, the creator of this profile";
400
401 if ((wo = (icmText *)wr_icco->add_tag(
402 wr_icco, icSigCopyrightTag, icSigTextType)) == NULL)
403 error("add_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
404
405 wo->size = strlen(crt)+1; /* Allocated and used size of text, inc null */
406 wo->allocate((icmBase *)wo);/* Allocate space */
407 strcpy(wo->data, crt); /* Copy the text in */
408 }
409 /* Device Manufacturers Description Tag: */
410 if (xpi != NULL && xpi->deviceMfgDesc != NULL) {
411 icmTextDescription *wo;
412 char *dst = xpi->deviceMfgDesc;
413
414 if ((wo = (icmTextDescription *)wr_icco->add_tag(
415 wr_icco, icSigDeviceMfgDescTag, icSigTextDescriptionType)) == NULL)
416 error("add_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
417
418 wo->size = strlen(dst)+1; /* Allocated and used size of desc, inc null */
419 wo->allocate((icmBase *)wo);/* Allocate space */
420 strcpy(wo->desc, dst); /* Copy the string in */
421 }
422 /* Model Description Tag: */
423 if (xpi != NULL && xpi->modelDesc != NULL) {
424 icmTextDescription *wo;
425 char *dst = xpi->modelDesc;
426
427 if ((wo = (icmTextDescription *)wr_icco->add_tag(
428 wr_icco, icSigDeviceModelDescTag, icSigTextDescriptionType)) == NULL)
429 error("add_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
430
431 wo->size = strlen(dst)+1; /* Allocated and used size of desc, inc null */
432 wo->allocate((icmBase *)wo);/* Allocate space */
433 strcpy(wo->desc, dst); /* Copy the string in */
434 }
435 /* White Point Tag: */
436 {
437 icmXYZArray *wo;
438 /* Note that tag types icSigXYZType and icSigXYZArrayType are identical */
439 if ((wo = (icmXYZArray *)wr_icco->add_tag(
440 wr_icco, icSigMediaWhitePointTag, icSigXYZArrayType)) == NULL)
441 error("add_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
442
443 wo->size = 1;
444 wo->allocate((icmBase *)wo); /* Allocate space */
445 wo->data[0].X = 0.9642; /* Set a default value - D50 */
446 wo->data[0].Y = 1.0000;
447 wo->data[0].Z = 0.8249;
448 }
449 /* Black Point Tag: */
450 {
451 icmXYZArray *wo;
452 if ((wo = (icmXYZArray *)wr_icco->add_tag(
453 wr_icco, icSigMediaBlackPointTag, icSigXYZArrayType)) == NULL)
454 error("add_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
455
456 wo->size = 1;
457 wo->allocate((icmBase *)wo); /* Allocate space */
458 wo->data[0].X = 0.00; /* Set default perfect black */
459 wo->data[0].Y = 0.00;
460 wo->data[0].Z = 0.00;
461 }
462
463 if (isLut == 0) { /* shaper + matrix type */
464
465 /* Red, Green and Blue Colorant Tags: */
466 {
467 icmXYZArray *wor, *wog, *wob;
468 if ((wor = (icmXYZArray *)wr_icco->add_tag(
469 wr_icco, icSigRedColorantTag, icSigXYZArrayType)) == NULL)
470 error("add_tag failed: %d, %s",rv,wr_icco->err);
471 if ((wog = (icmXYZArray *)wr_icco->add_tag(
472 wr_icco, icSigGreenColorantTag, icSigXYZArrayType)) == NULL)
473 error("add_tag failed: %d, %s",rv,wr_icco->err);
474 if ((wob = (icmXYZArray *)wr_icco->add_tag(
475 wr_icco, icSigBlueColorantTag, icSigXYZArrayType)) == NULL)
476 error("add_tag failed: %d, %s",rv,wr_icco->err);
477
478 wor->size = wog->size = wob->size = 1;
479 wor->allocate((icmBase *)wor); /* Allocate space */
480 wog->allocate((icmBase *)wog);
481 wob->allocate((icmBase *)wob);
482
483 /* Setup some sane dummy values */
484 /* icxMatrix will override these later */
485 wor->data[0].X = 1.0; wor->data[0].Y = 0.0; wor->data[0].Z = 0.0;
486 wog->data[0].X = 0.0; wog->data[0].Y = 1.0; wog->data[0].Z = 0.0;
487 wob->data[0].X = 0.0; wob->data[0].Y = 0.0; wob->data[0].Z = 1.0;
488 }
489
490 /* Red, Green and Blue Tone Reproduction Curve Tags: */
491 {
492 icmCurve *wor, *wog, *wob;
493 if ((wor = (icmCurve *)wr_icco->add_tag(
494 wr_icco, icSigRedTRCTag, icSigCurveType)) == NULL)
495 error("add_tag failed: %d, %s",rv,wr_icco->err);
496
497 if (isShTRC) { /* Make all TRCs shared */
498 if ((wog = (icmCurve *)wr_icco->link_tag(
499 wr_icco, icSigGreenTRCTag, icSigRedTRCTag)) == NULL)
500 error("link_tag failed: %d, %s",rv,wr_icco->err);
501 if ((wob = (icmCurve *)wr_icco->link_tag(
502 wr_icco, icSigBlueTRCTag, icSigRedTRCTag)) == NULL)
503 error("link_tag failed: %d, %s",rv,wr_icco->err);
504
505 } else { /* Else individual */
506 if ((wog = (icmCurve *)wr_icco->add_tag(
507 wr_icco, icSigGreenTRCTag, icSigCurveType)) == NULL)
508 error("add_tag failed: %d, %s",rv,wr_icco->err);
509 if ((wob = (icmCurve *)wr_icco->add_tag(
510 wr_icco, icSigBlueTRCTag, icSigCurveType)) == NULL)
511 error("add_tag failed: %d, %s",rv,wr_icco->err);
512 }
513
514 if (ptype == prof_shamat || ptype == prof_sha1mat) { /* Shaper */
515 wor->flag = wog->flag = wob->flag = icmCurveSpec;
516 wor->size = wog->size = wob->size = 256; /* Number of entries */
517 } else { /* Gamma */
518 wor->flag = wog->flag = wob->flag = icmCurveGamma;
519 wor->size = wog->size = wob->size = 1; /* Must be 1 for gamma */
520 }
521 wor->allocate((icmBase *)wor); /* Allocate space */
522 wog->allocate((icmBase *)wog);
523 wob->allocate((icmBase *)wob);
524
525 /* icxMatrix will set curve values */
526 }
527
528 } else { /* Lut type profile */
529
530 /* 16 bit dev -> pcs lut: */
531 {
532 icmLut *wo;
533
534 /* Only A2B0, no intent */
535 if ((wo = (icmLut *)wr_icco->add_tag(
536 wr_icco, icSigAToB0Tag, icSigLut16Type)) == NULL)
537 error("add_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
538
539 wo->inputChan = 3;
540 wo->outputChan = 3;
541 if (iquality >= 3) {
542 wo->clutPoints = 45;
543 wo->inputEnt = 2048;
544 wo->outputEnt = 2048;
545 } else if (iquality == 2) {
546 wo->clutPoints = 33;
547 wo->inputEnt = 2048;
548 wo->outputEnt = 2048;
549 } else if (iquality == 1) {
550 wo->clutPoints = 17;
551 wo->inputEnt = 1024;
552 wo->outputEnt = 1024;
553 } else {
554 wo->clutPoints = 9;
555 wo->inputEnt = 512;
556 wo->outputEnt = 512;
557 }
558
559 wo->allocate((icmBase *)wo);/* Allocate space */
560
561 /* icxLuLut will set tables values */
562 }
563
564 #ifdef DOB2A
565 /* 16 bit pcs -> dev lut: */
566 if (dob2a) {
567 icmLut *wo;
568
569 /* Only B2A0, no intent */
570 if ((wo = (icmLut *)wr_icco->add_tag(
571 wr_icco, icSigBToA0Tag, icSigLut16Type)) == NULL)
572 error("add_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
573
574 wo->inputChan = 3;
575 wo->outputChan = 3;
576 if (oquality >= 3) {
577 wo->clutPoints = 45;
578 wo->inputEnt = 2048;
579 wo->outputEnt = 2048;
580 } else if (oquality == 2) {
581 wo->clutPoints = 33;
582 wo->inputEnt = 2048;
583 wo->outputEnt = 2048;
584 } else if (oquality == 1) {
585 wo->clutPoints = 17;
586 wo->inputEnt = 1024;
587 wo->outputEnt = 1024;
588 } else if (oquality >= 0) {
589 wo->clutPoints = 9;
590 wo->inputEnt = 512;
591 wo->outputEnt = 512;
592 } else { /* Special, Extremely low quality */
593 wo->clutPoints = 3;
594 wo->inputEnt = 64;
595 wo->outputEnt = 64;
596 }
597
598 wo->allocate((icmBase *)wo);/* Allocate space */
599
600 /* We set the tables below */
601 }
602 #endif /* DOB2A */
603
604 }
605
606 /* Sample data use to create profile: */
607 if (nocied == 0) {
608 icmText *wo;
609 char *crt;
610 FILE *fp;
611
612 if ((wo = (icmText *)wr_icco->add_tag(
613 wr_icco, icmMakeTag('t','a','r','g'), icSigTextType)) == NULL)
614 error("add_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
615
616 #if defined(O_BINARY) || defined(_O_BINARY)
617 if ((fp = fopen(in_name, "rb")) == NULL)
618 #else
619 if ((fp = fopen(in_name, "r")) == NULL)
620 #endif
621 error("Unable to open input file '%s' for reading",in_name);
622
623 if (fseek(fp, 0, SEEK_END))
624 error("Unable to seek to end of file '%s'",in_name);
625 wo->size = ftell(fp) + 1; /* Size needed + null */
626 wo->allocate((icmBase *)wo);/* Allocate space */
627
628 if (fseek(fp, 0, SEEK_SET))
629 error("Unable to seek to end of file '%s'",in_name);
630
631 if (fread(wo->data, 1, wo->size-1, fp) != wo->size-1)
632 error("Failed to read file '%s'",in_name);
633 wo->data[wo->size-1] = '\000';
634 fclose(fp);
635
636 /* Duplicate for compatibility */
637 if (wr_icco->link_tag(
638 wr_icco, icmMakeTag('D','e','v','D'), icmMakeTag('t','a','r','g')) == NULL)
639 error("link_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
640 if (wr_icco->link_tag(
641 wr_icco, icmMakeTag('C','I','E','D'), icmMakeTag('t','a','r','g')) == NULL)
642 error("link_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
643 }
644
645 if ((npat = icg->t[0].nsets) <= 0)
646 error ("No sets of data");
647
648 if (verb) {
649 fprintf(verbo,"No of test patches = %d\n",npat);
650 }
651
652 if (extrap) {
653 npxpat = 4 * EXTRAP_MAXPNTS; /* Allow for up to 20 extra patches */
654 }
655
656 /* Allocate arrays to hold test patch input and output values */
657 if ((tpat = (cow *)malloc(sizeof(cow) * (npat + npxpat))) == NULL)
658 error("Malloc failed - tpat[]");
659
660 /* Read in the CGATs fields */
661 {
662 int ti;
663 int Xi, Yi, Zi;
664 int ri, gi, bi;
665
666 /* Check that we handle the color space */
667 if ((ti = icg->find_kword(icg, 0, "COLOR_REP")) < 0)
668 error ("Input file doesn't contain keyword COLOR_REPS");
669 if (strcmp(icg->t[0].kdata[ti],"LAB_RGB") == 0) {
670 isLab = 1;
671 } else {
672 if (strcmp(icg->t[0].kdata[ti],"XYZ_RGB") == 0) {
673 isLab = 0;
674 } else {
675 error ("Input device input file has unhandled color representation");
676 }
677 }
678
679 if ((ri = icg->find_field(icg, 0, "RGB_R")) < 0)
680 error ("Input file doesn't contain field RGB_R");
681 if (icg->t[0].ftype[ri] != r_t)
682 error ("Field RGB_R is wrong type - expect float");
683 if ((gi = icg->find_field(icg, 0, "RGB_G")) < 0)
684 error ("Input file doesn't contain field RGB_G");
685 if (icg->t[0].ftype[gi] != r_t)
686 error ("Field RGB_G is wrong type - expect float");
687 if ((bi = icg->find_field(icg, 0, "RGB_B")) < 0)
688 error ("Input file doesn't contain field RGB_B");
689 if (icg->t[0].ftype[bi] != r_t)
690 error ("Field RGB_B is wrong type - expect float");
691
692 if (spec == 0) { /* Using instrument tristimulous value */
693
694 if (isLab) {
695 if ((Xi = icg->find_field(icg, 0, "LAB_L")) < 0)
696 error ("Input file doesn't contain field LAB_L");
697 if (icg->t[0].ftype[Xi] != r_t)
698 error ("Field LAB_L is wrong type - expect float");
699 if ((Yi = icg->find_field(icg, 0, "LAB_A")) < 0)
700 error ("Input file doesn't contain field LAB_A");
701 if (icg->t[0].ftype[Yi] != r_t)
702 error ("Field LAB_A is wrong type - expect float");
703 if ((Zi = icg->find_field(icg, 0, "LAB_B")) < 0)
704 error ("Input file doesn't contain field LAB_B");
705 if (icg->t[0].ftype[Zi] != r_t)
706 error ("Field LAB_B is wrong type - expect float");
707 } else {
708 if ((Xi = icg->find_field(icg, 0, "XYZ_X")) < 0)
709 error ("Input file doesn't contain field XYZ_X");
710 if (icg->t[0].ftype[Xi] != r_t)
711 error ("Field XYZ_X is wrong type - expect float");
712 if ((Yi = icg->find_field(icg, 0, "XYZ_Y")) < 0)
713 error ("Input file doesn't contain field XYZ_Y");
714 if (icg->t[0].ftype[Yi] != r_t)
715 error ("Field XYZ_Y is wrong type - expect float");
716 if ((Zi = icg->find_field(icg, 0, "XYZ_Z")) < 0)
717 error ("Input file doesn't contain field XYZ_Z");
718 if (icg->t[0].ftype[Zi] != r_t)
719 error ("Field XYZ_Z is wrong type - expect float");
720 }
721
722 for (i = 0; i < npat; i++) {
723 tpat[i].w = 1.0;
724 tpat[i].p[0] = *((double *)icg->t[0].fdata[i][ri]) / 100.0;
725 tpat[i].p[1] = *((double *)icg->t[0].fdata[i][gi]) / 100.0;
726 tpat[i].p[2] = *((double *)icg->t[0].fdata[i][bi]) / 100.0;
727 if (tpat[i].p[0] > 1.0
728 || tpat[i].p[1] > 1.0
729 || tpat[i].p[2] > 1.0) {
730 error("At %d device values %f %f %f field exceeds 100.0!",i,100.0 * tpat[i].p[0],100.0 * tpat[i].p[1],100.0 * tpat[i].p[2]);
731 }
732 tpat[i].v[0] = *((double *)icg->t[0].fdata[i][Xi]);
733 tpat[i].v[1] = *((double *)icg->t[0].fdata[i][Yi]);
734 tpat[i].v[2] = *((double *)icg->t[0].fdata[i][Zi]);
735 if (!isLab) {
736 tpat[i].v[0] /= 100.0; /* Normalise XYZ to range 0.0 - 1.0 */
737 tpat[i].v[1] /= 100.0;
738 tpat[i].v[2] /= 100.0;
739 }
740 if (!isLab && wantLab) { /* Convert test patch result XYZ to PCS (D50 Lab) */
741 icmXYZ2Lab(&icmD50, tpat[i].v, tpat[i].v);
742 } else if (isLab && !wantLab) {
743 icmLab2XYZ(&icmD50, tpat[i].v, tpat[i].v);
744 }
745 }
746
747 } else { /* Using spectral data */
748 int j, ii;
749 xspect sp;
750 char buf[100];
751 int spi[XSPECT_MAX_BANDS]; /* CGATS indexes for each wavelength */
752 xsp2cie *sp2cie; /* Spectral conversion object */
753
754 if ((ii = icg->find_kword(icg, 0, "SPECTRAL_BANDS")) < 0)
755 error ("Input file doesn't contain keyword SPECTRAL_BANDS");
756 sp.spec_n = atoi(icg->t[0].kdata[ii]);
757 if ((ii = icg->find_kword(icg, 0, "SPECTRAL_START_NM")) < 0)
758 error ("Input file doesn't contain keyword SPECTRAL_START_NM");
759 sp.spec_wl_short = atof(icg->t[0].kdata[ii]);
760 if ((ii = icg->find_kword(icg, 0, "SPECTRAL_END_NM")) < 0)
761 error ("Input file doesn't contain keyword SPECTRAL_END_NM");
762 sp.spec_wl_long = atof(icg->t[0].kdata[ii]);
763 sp.norm = 100.0;
764
765 /* Find the fields for spectral values */
766 for (j = 0; j < sp.spec_n; j++) {
767 int nm;
768
769 /* Compute nearest integer wavelength */
770 nm = (int)(sp.spec_wl_short + ((double)j/(sp.spec_n-1.0))
771 * (sp.spec_wl_long - sp.spec_wl_short) + 0.5);
772
773 sprintf(buf,"SPEC_%03d",nm);
774
775 if ((spi[j] = icg->find_field(icg, 0, buf)) < 0)
776 error("Input file doesn't contain field %s",buf);
777
778 if (icg->t[0].ftype[spi[j]] != r_t)
779 error("Field %s is wrong type - expect float",buf);
780 }
781
782 /* Create a spectral conversion object */
783 if (emis) {
784 illum = icxIT_none;
785 cust_illum = NULL;
786 }
787 if ((sp2cie = new_xsp2cie(illum, cust_illum, observ, NULL,
788 wantLab ? icSigLabData : icSigXYZData, icxClamp)) == NULL)
789 error("Creation of spectral conversion object failed");
790
791 for (i = 0; i < npat; i++) {
792 tpat[i].w = 1.0;
793 tpat[i].p[0] = *((double *)icg->t[0].fdata[i][ri]) / 100.0;
794 tpat[i].p[1] = *((double *)icg->t[0].fdata[i][gi]) / 100.0;
795 tpat[i].p[2] = *((double *)icg->t[0].fdata[i][bi]) / 100.0;
796 if (tpat[i].p[0] > 1.0
797 || tpat[i].p[1] > 1.0
798 || tpat[i].p[2] > 1.0) {
799 error("Patch %d device value field exceeds 100.0!",i+1);
800 }
801
802 /* Read the spectral values for this patch */
803 for (j = 0; j < sp.spec_n; j++) {
804 sp.spec[j] = *((double *)icg->t[0].fdata[i][spi[j]]);
805 }
806
807 /* Convert it to CIE space */
808 sp2cie->convert(sp2cie, tpat[i].v, &sp);
809 }
810
811 sp2cie->del(sp2cie); /* Done with this */
812
813 }
814 } /* End of reading in CGATs file */
815
816 if (isLut == 0) { /* Gamma/Shaper + matrix profile */
817 xicc *wr_xicc; /* extention object */
818 icxLuBase *xluo; /* Forward ixcLu */
819 int flags = 0;
820
821 /* Wrap with an expanded icc */
822 if ((wr_xicc = new_xicc(wr_icco)) == NULL)
823 error("Creation of xicc failed");
824
825 if (verb)
826 flags |= ICX_VERBOSE;
827
828 if (ptype == prof_matonly)
829 flags |= ICX_NO_IN_SHP_LUTS; /* Make it linear */
830
831 if (clipprims)
832 flags |= ICX_CLIP_WB | ICX_CLIP_PRIMS;
833
834 flags |= ICX_SET_BLACK; /* Compute & use black */
835 flags |= ICX_SET_WHITE; /* Compute & use white */
836
837 /* ICX_SET_WHITE_C isn't applicable to matrix profiles */
838 if (autowpsc)
839 flags |= ICX_SET_WHITE_US; /* Compute & use white without scaling to L */
840
841 flags |= ICX_WRITE_WBL; /* Matrix: write white/black/luminence */
842
843 /* Setup Device -> XYZ conversion (Fwd) object from scattered data. */
844 if ((xluo = wr_xicc->set_luobj(
845 // wr_xicc, icmFwd, icRelativeColorimetric,
846 wr_xicc, icmFwd, icmDefaultIntent,
847 icmLuOrdNorm,
848 flags, /* Flags */
849 npat, npat, tpat, NULL, 0.0, wpscale,
850 // NULL, /* bpo */
851 smooth, avgdev, 1.0,
852 NULL, NULL, NULL, iquality)) == NULL)
853 error("%d, %s",wr_xicc->errc, wr_xicc->err);
854
855 /* Free up xicc stuff */
856 xluo->del(xluo);
857 wr_xicc->del(wr_xicc);
858
859 } else { /* cLUT based profile */
860 int flags = 0;
861 icxMatrixModel *mm = NULL;
862
863 xicc *wr_xicc; /* extention object */
864 icxLuBase *AtoB; /* AtoB ixcLu */
865
866 if (extrap) {
867 cow *mpat;
868 int nmpat;
869 int j;
870
871 if (verb) printf("Creating extrapolation black and white points:\n");
872
873 if ((mpat = (cow *)malloc(sizeof(cow) * npat)) == NULL)
874 error("Malloc failed - mpat[]");
875
876 /* Weight points from full set to build matrix model */
877 /* to extrapolate the neutral axis */
878 for (nmpat = j = 0; j < npat; j++) {
879 double mnp, mxp;
880 int k;
881
882 icmCpy3(mpat[nmpat].p, tpat[j].p);
883 icmCpy3(mpat[nmpat].v, tpat[j].v);
884
885 /* Locate largest/smallest RGB value */
886 mxp = -1e6, mnp = 1e6;
887 for (k = 0; k < 3; k++) {
888 if (tpat[j].p[k] > mxp)
889 mxp = tpat[j].p[k];
890 if (tpat[j].p[k] < mnp)
891 mnp = tpat[j].p[k];
892 }
893 mxp -= mnp; /* Spread; 0 for R=G=B */
894
895 mpat[nmpat].w = pow(1.1 - mxp, 2.0);
896 //printf("~1 added value %d: %f %f %f -> %f %f %f wt %f\n",j, mpat[nmpat].p[0], mpat[nmpat].p[1], mpat[nmpat].p[2], mpat[nmpat].v[0], mpat[nmpat].v[1], mpat[nmpat].v[2],mpat[nmpat].w);
897 nmpat++;
898 }
899
900 /* Create gamma/matrix model to extrapolate with. */
901 /* (Use ofset & gain, gamma curve as 0th order with 1 harmonic, */
902 /* and smooth it.) */
903 if ((mm = new_MatrixModel(wr_icco, verb, nmpat, mpat, wantLab,
904 /* quality */ -1, /* isLinear */ ptype == prof_matonly,
905 /* isGamma */ 0, /* isShTRC */ 0,
906 /* shape0gam */ 1, /* clipbw */ 0, /* clipprims */ 0,
907 // /* smooth */ 1.0, /* scale */ 0.7)) == NULL) {
908 /* smooth */ 1.0, /* scale */ 1.0)) == NULL) {
909 error("Creating extrapolation matrix model failed - memory ?");
910 }
911
912 #ifdef NEVER /* Plot Lab of model */
913 {
914 #define XRES 100
915 double xx[XRES];
916 double y0[XRES];
917 double y1[XRES];
918 double y2[XRES];
919
920 /* Display the result fit */
921 for (i = 0; i < XRES; i++) {
922 double rgb[3], lab[3];
923 xx[i] = rgb[0] = rgb[1] = rgb[2] = i/(double)(XRES-1);
924 mm->lookup(mm, lab, rgb);
925 if (wantLab)
926 icmLab2XYZ(&icmD50,lab,lab);
927 y0[i] = lab[0];
928 y1[i] = lab[1];
929 y2[i] = lab[2];
930 }
931 do_plot(xx,y0,y1,y2,XRES);
932 }
933 #endif /* DEBUG_PLOT */
934 }
935
936 if (extrap) {
937 int ii, wix = 0, j;
938 int pcsy; /* Effective PCS L or Y chanel index */
939 double wpy = -1e60;
940 double dwhite[MXDI]; /* Device white */
941 double mxdw;
942 double avgdist; /* Average distance between points */
943
944 /* Figure out the device white point. */
945 /* Note that this is duplicating code in xicc/xmatrix.c */
946 /* and xfit.c */
947
948 if (wantLab)
949 pcsy = 0; /* L or Lab */
950 else
951 pcsy = 1; /* Y of XYZ */
952
953 for (i = 0; i < npat; i++) {
954 double labv[3], yv;
955
956 /* Create D50 Lab to allow some chromatic sensitivity */
957 /* in picking the white point */
958 if (wantLab)
959 icmCpy3(labv, tpat[i].v);
960 else
961 icmXYZ2Lab(&icmD50, labv, tpat[i].v);
962
963 /* Tilt things towards D50 neutral white patches */
964 yv = labv[0] - 0.3 * sqrt(labv[1] * labv[1] + labv[2] * labv[2]);
965 if (yv > wpy) {
966 for (j = 0; j < 3; j++)
967 dwhite[j] = tpat[i].p[j];
968 wpy = yv;
969 wix = i;
970 }
971 }
972
973 /* Fix extrapolation matrix to be perfect at white point */
974 mm->force(mm, tpat[wix].v, tpat[wix].p);
975
976 /* Scale the white point to make one dev value 1.0 */
977 mxdw = -1;
978 for (j = 0; j < 3; j++) {
979 if (dwhite[j] > mxdw)
980 mxdw = dwhite[j];
981 }
982 for (j = 0; j < 3; j++) {
983 dwhite[j] /= mxdw;
984 }
985
986 avgdist = pow(1.0/(double)npat, 1.0/3.0);
987 if (avgdist < 0.001)
988 avgdist = 0.001;
989 else if (avgdist > 0.3)
990 avgdist = 0.3;
991 //printf("~1 avgdist = %f\n",avgdist);
992
993 /* For points with white point device ratio, */
994 /* and points with R=G=B ratio, create extrapolation points. */
995 for (ii = 0; ii < 2; ii++) {
996
997 if (ii > 1)
998 dwhite[0] = dwhite[1] = dwhite[2] = 1.0;
999
1000 /* Create a series of black and white patch */
1001 for (i = 0; i < 2; i++) {
1002 int cix; /* Closest point index */
1003 int eix; /* End point index */
1004 double cde = 1e60; /* Closest point distance */
1005 double tt;
1006 double corr[3], cwt; /* Correction */
1007
1008 eix = npat + nxpat;
1009
1010 icmScale3(tpat[eix].p, dwhite, (double)i);
1011
1012 /* Locate closest point */
1013 for (j = 0; j < npat; j++) {
1014 double mnp, mxp;
1015 int k;
1016
1017 /* Locate largest/smallest RGB value */
1018 mxp = -1e6, mnp = 1e6;
1019 for (k = 0; k < 3; k++) {
1020 if (tpat[j].p[k] > mxp)
1021 mxp = tpat[j].p[k];
1022 if (tpat[j].p[k] < mnp)
1023 mnp = tpat[j].p[k];
1024 }
1025 mxp -= mnp; /* Spread; 0 for R=G=B */
1026
1027 tt = icmNorm33(tpat[eix].p, tpat[j].p);
1028 tt += mxp;
1029
1030 if (tt < cde) {
1031 cde = tt;
1032 cix = j;
1033 }
1034 }
1035
1036 //printf("~1 closest %d: de %f, %f %f %f -> %f %f %f\n",cix, cde, tpat[cix].p[0], tpat[cix].p[1], tpat[cix].p[2], tpat[cix].v[0], tpat[cix].v[1], tpat[cix].v[2]);
1037
1038 //{
1039 //double val[3];
1040 //mm->lookup(mm, val, tpat[cix].p);
1041 //printf("~1 closest gam/matrix -> %f %f %f\n",val[0],val[1],val[2]);
1042 //}
1043
1044 /* Lookup matrix value for our new point */
1045 mm->lookup(mm, tpat[eix].v, tpat[eix].p);
1046 //printf("~1 got value %d: %f %f %f -> %f %f %f\n",i, tpat[eix].p[0], tpat[eix].p[1], tpat[eix].p[2], tpat[eix].v[0], tpat[eix].v[1], tpat[eix].v[2]);
1047 /* Weight the extra point so that it doesn't overpower the */
1048 /* nearest real point to it too much. */
1049 tt = cde;
1050 if (tt > avgdist) /* Distance at which sythetic point has 100% weight */
1051 tt = avgdist;
1052 tpat[eix].w = 0.5 * EXTRAP_WEIGHT * tt/avgdist;
1053 //printf("~1 weight %f\n",tpat[eix].w);
1054 if (verb)
1055 printf("Added synthetic point @ %f %f %f, val %f %f %f, weight %f\n",tpat[eix].p[0], tpat[eix].p[1], tpat[eix].p[2], tpat[eix].v[0], tpat[eix].v[1], tpat[eix].v[2],tpat[eix].w);
1056 nxpat++;
1057
1058 /* If there is a lot of space, add a second intemediate point */
1059 //printf("~1 cde = %f, avgdist = %f\n",cde,avgdist);
1060 if (cde >= (0.5 * avgdist)) {
1061 int nxps; /* Number of extra points including end point */
1062 nxps = 1 + (int)(cde/(0.5 * avgdist));
1063 if (nxps > EXTRAP_MAXPNTS)
1064 nxps = EXTRAP_MAXPNTS;
1065
1066 //printf("~1 nxps = %d\n",nxps);
1067 for (j = 1; j < nxps; j++) {
1068 double bl, ipos;
1069
1070 bl = j/(nxps + 1.0);
1071
1072 ipos = (1.0 - bl) * tpat[eix].p[0]
1073 + bl * (tpat[cix].p[0] + tpat[cix].p[1] + tpat[cix].p[1])/3.0;
1074 icmScale3(tpat[npat + nxpat].p, dwhite, ipos);
1075
1076 /* Lookup matrix value for our new point */
1077 mm->lookup(mm, tpat[npat + nxpat].v, tpat[npat + nxpat].p);
1078
1079 /* Weight the extra point so that it doesn't overpower the */
1080 /* nearest real point to it too much. */
1081 cde = icmNorm33(tpat[cix].p, tpat[npat + nxpat].p);
1082
1083 if (cde > avgdist) /* Distance at which sythetic point has 100% weight */
1084 cde = avgdist;
1085 tpat[npat + nxpat].w = 0.5 * EXTRAP_WEIGHT * cde/avgdist;
1086 if (verb)
1087 printf("Added synthetic point @ %f %f %f, val %f %f %f, weight %f\n",tpat[npat + nxpat].p[0], tpat[npat + nxpat].p[1], tpat[npat + nxpat].p[2], tpat[npat + nxpat].v[0], tpat[npat + nxpat].v[1], tpat[npat + nxpat].v[2],tpat[npat + nxpat].w);
1088 nxpat++;
1089 }
1090 }
1091 }
1092 }
1093 }
1094
1095 /* Wrap with an expanded icc */
1096 if ((wr_xicc = new_xicc(wr_icco)) == NULL)
1097 error ("Creation of xicc failed");
1098
1099 flags |= ICX_CLIP_NEAREST; /* This will avoid clip caused rev setup */
1100
1101 if (noisluts)
1102 flags |= ICX_NO_IN_SHP_LUTS;
1103
1104 if (noipluts)
1105 flags |= ICX_NO_IN_POS_LUTS;
1106
1107 if (nooluts)
1108 flags |= ICX_NO_OUT_LUTS;
1109
1110 if (verb)
1111 flags |= ICX_VERBOSE;
1112
1113 if (clipprims)
1114 flags |= ICX_CLIP_WB;
1115
1116 flags |= ICX_SET_BLACK; /* Compute & use black */
1117 flags |= ICX_SET_WHITE; /* Compute & use white */
1118 if (clipovwp)
1119 flags |= ICX_SET_WHITE_C; /* Compute & use white and clip cLUT over D50 */
1120 else if (autowpsc)
1121 flags |= ICX_SET_WHITE_US; /* Compute & use white without scaling to L */
1122
1123 /* Setup RGB -> Lab conversion object from scattered data. */
1124 /* Note that we've layered it on a native XYZ icc profile. */
1125 /* (The skeleton model is not used - it doesn't seem to help) */
1126 if ((AtoB = wr_xicc->set_luobj(
1127 wr_xicc, icmFwd, icmDefaultIntent,
1128 icmLuOrdNorm,
1129 #ifdef USE_EXTRA_FITTING
1130 ICX_EXTRA_FIT |
1131 #endif
1132 #ifdef USE_2PASSSMTH
1133 ICX_2PASSSMTH |
1134 #endif
1135 flags, /* Flags */
1136 npat + nxpat, npat, tpat, NULL, 0.0, wpscale,
1137 // NULL, /* bpo */
1138 smooth, avgdev, 1.0,
1139 NULL, NULL, NULL, iquality)) == NULL)
1140 error ("%d, %s",wr_xicc->errc, wr_xicc->err);
1141
1142 if (mm != NULL)
1143 mm->del(mm);
1144
1145 /* Free up xicc stuff */
1146 AtoB->del(AtoB);
1147
1148 #ifdef DOB2A
1149 if (dob2a) {
1150 icmLut *wo;
1151
1152 in_b2a_callback cx;
1153
1154 if (verb)
1155 printf("Setting up B to A table lookup\n");
1156
1157 /* Get a suitable forward conversion object to invert. */
1158 /* By creating a separate one to the one created using scattered data, */
1159 /* we ge the chance to set ICX_CAM_CLIP. It is always set to Lab 'PCS' */
1160 {
1161 int flags = 0;
1162
1163 if (verb)
1164 flags |= ICX_VERBOSE;
1165
1166 flags |= ICX_CLIP_NEAREST; /* Not vector clip */
1167
1168 #ifdef USE_CAM_CLIP_OPT
1169 flags |= ICX_CAM_CLIP; /* Clip in CAM Jab space rather than Lab */
1170 #else
1171 warning("!!!! USE_CAM_CLIP_OPT in profout.c is off !!!!");
1172 #endif
1173 if ((AtoB = wr_xicc->get_luobj(wr_xicc, flags, icmFwd,
1174 icmDefaultIntent,
1175 wantLab ? icSigLabData : icSigXYZData,
1176 icmLuOrdNorm, NULL, NULL)) == NULL)
1177 error ("%d, %s",wr_xicc->errc, wr_xicc->err);
1178 }
1179
1180 /* setup context ready for B2A table setting */
1181 cx.verb = verb;
1182 cx.pcsspace = wantLab ? icSigLabData : icSigXYZData;
1183 cx.wantLab = wantLab; /* Copy PCS flag over */
1184 #ifdef NO_B2A_PCS_CURVES
1185 cx.noPCScurves = 1; /* Don't use PCS curves */
1186 #else
1187 cx.noPCScurves = 0;
1188 #endif
1189 cx.devspace = icSigRgbData;
1190 cx.x = (icxLuLut *)AtoB; /* A2B icxLuLut created from scattered data */
1191
1192 if ((wo = (icmLut *)wr_icco->read_tag(
1193 wr_icco, icSigBToA0Tag)) == NULL)
1194 error("read_tag failed: %d, %s",wr_icco->errc,wr_icco->err);
1195
1196 /* We now setup an exact inverse, colorimetric style */
1197 /* Use helper function to do the hard work. */
1198
1199 if (cx.verb) {
1200 unsigned int ui;
1201 int extra;
1202 cx.count = 0;
1203 cx.last = -1;
1204 for (cx.total = 1, ui = 0; ui < wo->inputChan; ui++, cx.total *= wo->clutPoints)
1205 ;
1206 /* Add in cell center points */
1207 for (extra = 1, ui = 0; ui < wo->inputChan; ui++, extra *= (wo->clutPoints-1))
1208 ;
1209 cx.total += extra;
1210 printf("Creating B to A tables\n");
1211 printf(" 0%%"); fflush(stdout);
1212 }
1213
1214 if (icmSetMultiLutTables(
1215 1,
1216 &wo,
1217 ICM_CLUT_SET_APXLS, /* Use least squared aprox. */
1218 &cx, /* Context */
1219 cx.pcsspace, /* Input color space */
1220 icSigRgbData, /* Output color space */
1221 in_b2a_input, /* Input transform PCS->PCS' */
1222 NULL, NULL, /* Use default Lab' range */
1223 in_b2a_clut, /* Lab' -> Device' transfer function */
1224 NULL, NULL, /* Use default Device' range */
1225 in_b2a_output, /* Output transfer function, Device'->Device */
1226 NULL, NULL) != 0) /* Use default APXLS range */
1227 error("Setting 16 bit PCS->Device Lut failed: %d, %s",wr_icco->errc,wr_icco->err);
1228 if (cx.verb) {
1229 printf("\n");
1230 }
1231 #ifdef WARN_CLUT_CLIPPING
1232 if (wr_icco->warnc)
1233 warning("Values clipped in setting LUT");
1234 #endif /* WARN_CLUT_CLIPPING */
1235
1236 if (verb)
1237 printf("Done B to A table\n");
1238 AtoB->del(AtoB);
1239 }
1240 #endif /* DOB2A */
1241 wr_xicc->del(wr_xicc);
1242
1243 }
1244
1245 /* Write the file (including all tags) out */
1246 if ((rv = wr_icco->write(wr_icco,wr_fp,0)) != 0)
1247 error ("Write file: %d, %s",rv,wr_icco->err);
1248
1249 /* Close the file */
1250 wr_icco->del(wr_icco);
1251 wr_fp->del(wr_fp);
1252
1253 /* Check the profile accuracy against the data points */
1254 if (verb || verify) {
1255 icmFile *rd_fp;
1256 icc *rd_icco;
1257 icmLuBase *luo;
1258 double merr = 0.0;
1259 double aerr = 0.0;
1260 double nsamps = 0.0;
1261
1262 /* Open up the file for reading */
1263 if ((rd_fp = new_icmFileStd_name(file_name,"r")) == NULL)
1264 error ("Write: Can't open file '%s'",file_name);
1265
1266 if ((rd_icco = new_icc()) == NULL)
1267 error ("Write: Creation of ICC object failed");
1268
1269 /* Read the header and tag list */
1270 if ((rv = rd_icco->read(rd_icco,rd_fp,0)) != 0)
1271 error ("Read: %d, %s",rv,rd_icco->err);
1272
1273 /* ~~ should use an xluobj with merge output ~~~ */
1274 /* Get the A2B table */
1275 if ((luo = rd_icco->get_luobj(rd_icco, icmFwd,
1276 icAbsoluteColorimetric, icSigLabData, icmLuOrdNorm)) == NULL) {
1277 error ("%d, %s",rd_icco->errc, rd_icco->err);
1278 }
1279
1280 for (i = 0; i < npat; i++) {
1281 double out[3], ref[3];
1282 double mxd;
1283
1284 if (luo->lookup(luo, out, tpat[i].p) > 1)
1285 error ("%d, %s",rd_icco->errc,rd_icco->err);
1286
1287 /* Our tpat data might be in XYZ, so generate an Lab ref value */
1288 if (!wantLab) { /* Convert test patch result XYZ to PCS (D50 Lab) */
1289 icmXYZ2Lab(&icmD50, ref, tpat[i].v);
1290
1291 } else {
1292 ref[0] = tpat[i].v[0];
1293 ref[1] = tpat[i].v[1];
1294 ref[2] = tpat[i].v[2];
1295 }
1296
1297 if (verb && verify) {
1298 printf("[%f] %f %f %f -> %f %f %f should be %f %f %f\n",
1299 icmLabDE(ref, out),
1300 tpat[i].p[0],tpat[i].p[1],tpat[i].p[2],
1301 out[0],out[1],out[2],
1302 ref[0],ref[1],ref[2]);
1303 }
1304
1305 /* Check the result */
1306 mxd = icmLabDE(ref, out);
1307 if (mxd > merr)
1308 merr = mxd;
1309
1310 aerr += mxd;
1311 nsamps++;
1312 }
1313 printf("Profile check complete, peak err = %f, avg err = %f\n",merr,aerr/nsamps);
1314
1315 /* Done with lookup object */
1316 luo->del(luo);
1317
1318 /* Close the file */
1319 rd_icco->del(rd_icco);
1320 rd_fp->del(rd_fp);
1321 }
1322
1323 free(tpat);
1324 }
1325
1326
1327