summaryrefslogtreecommitdiff
path: root/xicc/cam97test.c
diff options
context:
space:
mode:
Diffstat (limited to 'xicc/cam97test.c')
-rw-r--r--xicc/cam97test.c456
1 files changed, 456 insertions, 0 deletions
diff --git a/xicc/cam97test.c b/xicc/cam97test.c
new file mode 100644
index 0000000..ff71bc8
--- /dev/null
+++ b/xicc/cam97test.c
@@ -0,0 +1,456 @@
+
+/*
+ * International Color Consortium color transform expanded support
+ *
+ * Author: Graeme W. Gill
+ * Date: 1/5/01
+ * Version: 1.20
+ *
+ * Copyright 2000-2004 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.
+ *
+ */
+
+/*
+ * This is some test code to test the CIECAM97s3 functionality.
+ */
+
+
+#include <stdio.h>
+#include <math.h>
+#include "xcam.h"
+#include "cam97s3.h"
+
+#define DIAG /* Print internal value diagnostics for each conversion */
+#define SPOTTEST /* Test known spot colors */
+#undef INVTEST /* Lab cube to XYZ to Jab to XYZ */
+#undef INVTEST2 /* Jab cube to XYZ to Jab */
+#define TRES 33
+#define USE_HK 0 /* Use Helmholtz-Kohlraush */
+#define NOCAMCLIP 1 /* Don't clip before XYZ2Jab */
+
+#ifndef _isnan
+#define _isnan(x) ((x) != (x))
+#define _finite(x) ((x) == (x))
+#endif
+
+static void
+Lab2XYZ(double *out, double *in) {
+ double L = in[0], a = in[1], b = in[2];
+ double x,y,z,fx,fy,fz;
+
+ if (L > 8.0) {
+ fy = (L + 16.0)/116.0;
+ y = pow(fy,3.0);
+ } else {
+ y = L/903.2963058;
+ fy = 7.787036979 * y + 16.0/116.0;
+ }
+
+ fx = a/500.0 + fy;
+ if (fx > 24.0/116.0)
+ x = pow(fx,3.0);
+ else
+ x = (fx - 16.0/116.0)/7.787036979;
+
+ fz = fy - b/200.0;
+ if (fz > 24.0/116.0)
+ z = pow(fz,3.0);
+ else
+ z = (fz - 16.0/116.0)/7.787036979;
+
+ out[0] = x * 0.9642;
+ out[1] = y * 1.0000;
+ out[2] = z * 0.8249;
+}
+
+/* CIE XYZ to perceptual Lab */
+static void
+XYZ2Lab(double *out, double *in) {
+ double X = in[0], Y = in[1], Z = in[2];
+ double x,y,z,fx,fy,fz;
+
+ x = X/0.9642;
+ y = Y/1.0000;
+ z = Z/0.8249;
+
+ if (x > 0.008856451586)
+ fx = pow(x,1.0/3.0);
+ else
+ fx = 7.787036979 * x + 16.0/116.0;
+
+ if (y > 0.008856451586)
+ fy = pow(y,1.0/3.0);
+ else
+ fy = 7.787036979 * y + 16.0/116.0;
+
+ if (z > 0.008856451586)
+ fz = pow(z,1.0/3.0);
+ else
+ fz = 7.787036979 * z + 16.0/116.0;
+
+ out[0] = 116.0 * fy - 16.0;
+ out[1] = 500.0 * (fx - fy);
+ out[2] = 200.0 * (fy - fz);
+}
+
+
+/* Return maximum difference */
+double maxdiff(double in1[3], double in2[3]) {
+ int i;
+ double tt, rv = 0.0;
+
+ /* Deal with the possibility that we have nan's */
+ for (i = 0; i < 3; i++) {
+ tt = fabs(in1[i] - in2[i]);
+ if (_isnan(tt))
+ return tt;
+ if (tt > rv)
+ rv = tt;
+ }
+ return rv;
+}
+
+/* Return absolute difference */
+double absdiff(double in1[3], double in2[3]) {
+ double tt, rv = 0.0;
+ tt = in1[0] - in2[0];
+ rv += tt * tt;
+ tt = in1[1] - in2[1];
+ rv += tt * tt;
+ tt = in1[2] - in2[2];
+ rv += tt * tt;
+ return sqrt(rv);
+}
+
+/* Return maximum Lab difference of XYZ */
+double maxxyzdiff(double i1[3], double i2[3]) {
+ int i;
+ double tt, rv = 0.0;
+ double in1[3], in2[3];
+
+ XYZ2Lab(in1, i1);
+ XYZ2Lab(in2, i2);
+
+ /* Deal with the possibility that we have nan's */
+ for (i = 0; i < 3; i++) {
+ tt = fabs(in1[i] - in2[i]);
+ if (_isnan(tt))
+ return tt;
+ if (tt > rv)
+ rv = tt;
+ }
+ return rv;
+}
+
+int
+main(void) {
+ int ok = 1;
+ double white[6][3] = {
+ { 0.9505, 1.000, 1.0888 },
+ { 0.9505, 1.000, 1.0888 },
+ { 1.0985, 1.000, 0.3558 },
+ { 1.0985, 1.000, 0.3558 },
+ { 0.9505, 1.000, 1.0888 },
+ { 0.9642, 1.000, 0.8249 } /* D50 for inversion tests */
+ };
+ double La[6] = { 318.31, 31.83, 318.31, 31.83, 318.31, 150.0 };
+ double sample[5][3] = {
+ { 0.1901, 0.2000, 0.2178 },
+ { 0.5706, 0.4306, 0.3196 },
+ { 0.0353, 0.0656, 0.0214 },
+ { 0.1901, 0.2000, 0.2178 },
+ { 0.9505, 1.000, 1.0888 } /* Check white */
+ };
+
+#ifdef CIECAM97S3_HK
+#ifdef CIECAM97S3_SPLINE_E
+ double correct[5][4] = { /* Hue range is last two */
+ { 41.14, 0.06, 225.3, 251.9 },
+ { 69.04, 71.04, 19.4, 19.4 },
+ { 35.09, 87.15, 175.3, 175.3 },
+ { 55.64, 82.43, 252.5, 252.5 },
+ { 100.00000, 0.0, 0.0, 360.0 }
+ };
+#else
+ double correct[5][4] = { /* Hue range is last two */
+ { 41.14, 0.06, 225.3, 251.9 },
+ { 69.056048, 71.20763, 19.4, 19.4 },
+ { 35.365009, 88.65238, 175.3, 175.3 },
+ { 55.266184, 80.55128, 252.5, 252.5 },
+ { 100.00000, 0.0, 0.0, 360.0 }
+ };
+#endif
+#else
+ double correct[5][4] = { /* Hue range is last two */
+ { 41.13, 0.05, 225.3, 251.9 },
+ { 64.14, 71.22, 19.4, 19.4 },
+ { 19.18, 88.64, 175.3, 175.3 },
+ { 39.11, 80.55, 252.5, 252.5 },
+ { 100.0, 0.00, 0.0, 360.0 } /* Check white */
+ };
+#endif
+
+ double Jab[3];
+ double res[3];
+ cam97s3 *cam;
+ int c;
+
+ cam = new_cam97s3();
+
+#ifdef SPOTTEST
+ for (c = 0; c < 5; c++) {
+
+#ifdef DIAG
+ printf("Case %d:\n",c);
+#endif /* DIAG */
+ cam->set_view(
+ cam,
+ vc_average, /* Enumerated Viewing Condition */
+ white[c], /* Reference/Adapted White XYZ (Y range 0.0 .. 1.0) */
+ La[c], /* Adapting/Surround Luminance cd/m^2 */
+ 0.20, /* Relative Luminance of Background to reference white */
+ 0.0, /* Luminance of white in image - not used */
+ 0.00, /* Flare as a fraction of the reference white (Y range 0.0 .. 1.0) */
+ white[c], /* The Flare color coordinates (typically the Ambient color) */
+ USE_HK
+ );
+#ifdef DIAG
+ printf("\n");
+#endif /* DIAG */
+
+ {
+ double JCh[3], target[3];
+ cam->XYZ_to_cam(cam, Jab, sample[c]);
+
+ /* Convert to JCh for checking */
+ JCh[0] = Jab[0];
+
+ /* Compute hue angle */
+ JCh[2] = (180.0/3.14159265359) * atan2(Jab[2], Jab[1]);
+ JCh[2] = (JCh[2] < 0.0) ? JCh[2] + 360.0 : JCh[2];
+
+ /* Compute chroma value */
+ JCh[1] = sqrt(Jab[1] * Jab[1] + Jab[2] * Jab[2]);
+
+ target[0] = correct[c][0];
+ target[1] = correct[c][1];
+ if (JCh[2] >= correct[c][2] && JCh[2] <= correct[c][3]) {
+ target[2] = JCh[2];
+ } else {
+ if (fabs(JCh[2] - correct[c][2]) < fabs(JCh[2] - correct[c][3]))
+ target[2] = correct[c][2];
+ else
+ target[2] = correct[c][3];
+ }
+
+ cam->cam_to_XYZ(cam, res, Jab);
+
+ if (maxdiff(JCh, target) > 0.05) {
+ printf("Spottest: Excesive error in conversion to CAM %f\n",maxdiff(JCh, target));
+ ok = 0;
+ }
+ if (maxdiff(sample[c], res) > 1e-5) {
+ printf("Spottest: Excessive error in inversion %f\n",maxdiff(sample[c], res));
+ ok = 0;
+ }
+#ifdef DIAG
+ printf("Jab is %f %f %f, Jch is %f %f %f\n",
+ Jab[0], Jab[1], Jab[2],
+ JCh[0], JCh[1], JCh[2]);
+
+ printf("Error to expected value = %f\n",maxdiff(JCh, target));
+
+ printf("XYZ is %f %f %f\n",res[0], res[1], res[2]);
+ printf("Inversion error = %f\n",maxdiff(sample[c], res));
+#endif /* DIAG */
+ }
+ }
+#endif /* SPOTTEST */
+
+#ifdef INVTEST
+ for (c = 5; c < 6; c++) {
+ /* Get the range of Jab values */
+ double min[3] = { 1e6, 1e6, 1e6 };
+ double max[3] = { -1e6, -1e6, -1e6 };
+
+ cam->set_view(
+ cam,
+ vc_average, /* Enumerated Viewing Condition */
+ white[c], /* Reference/Adapted White XYZ (Y range 0.0 .. 1.0) = D50 */
+ 34.0, /* Adapting/Surround Luminance cd/m^2 */
+ 0.20, /* Relative Luminance of Background to reference white */
+ 0.0, /* Luminance of white in image - not used */
+ 0.01, /* Flare as a fraction of the reference white (Y range 0.0 .. 1.0) */
+ white[c], /* The Flare color coordinates (typically the Ambient color) */
+ USE_HK
+ );
+
+ {
+ int co0, co1, co2; /* (using co[3] triggers compiler bug) */
+ double merr = 0.0;
+ double xyz[3], Lab[3], Jab[3], checkxyz[3];
+
+
+#ifdef NEVER
+ /* Test case */
+ Lab[0] = 0.0;
+ Lab[1] = -180.0;
+ Lab[2] = -180.0;
+ Lab2XYZ(xyz, Lab);
+ cam->XYZ_to_cam(cam, Jab, xyz);
+ cam->cam_to_XYZ(cam, checkxyz, Jab);
+
+#else /* !NEVER */
+
+ /* itterate through Lab space */
+ for (co0 = 0; co0 < TRES; co0++) {
+ Lab[0] = co0/(TRES-1.0);
+ Lab[0] = Lab[0] * 100.0;
+ for (co1 = 0; co1 < TRES; co1++) {
+ Lab[1] = co1/(TRES-1.0);
+ Lab[1] = (Lab[1] - 0.5) * 256.0;
+ for (co2 = 0; co2 < TRES; co2++) {
+ double mxd;
+ int i;
+ Lab[2] = co2/(TRES-1.0);
+ Lab[2] = (Lab[2] - 0.5) * 256.0;
+
+ Lab2XYZ(xyz, Lab);
+ cam->XYZ_to_cam(cam, Jab, xyz);
+ for (i = 0; i < 3; i++) {
+ if (Jab[i] < min[i])
+ min[i] = Jab[i];
+ if (Jab[i] > max[i])
+ max[i] = Jab[i];
+ }
+ cam->cam_to_XYZ(cam, checkxyz, Jab);
+
+ /* Check the result */
+ mxd = maxxyzdiff(checkxyz, xyz);
+ if (_finite(merr) && (_isnan(mxd) || mxd > merr))
+ merr = mxd;
+#ifdef DIAG
+ if (_isnan(mxd) || mxd > 0.1) {
+ printf("\n");
+ printf("Lab = %f %f %f\n",Lab[0], Lab[1], Lab[2]);
+ printf("%f %f %f -> %f %f %f -> %f %f %f [%f]\n",
+ xyz[0],xyz[1],xyz[2],Jab[0],Jab[1],Jab[2],
+ checkxyz[0],checkxyz[1],checkxyz[2],
+ maxxyzdiff(checkxyz, xyz));
+ }
+#endif /* DIAG */
+ }
+ }
+ }
+ if (_isnan(merr) || merr > 0.1) {
+ printf("Excessive error in inversion check %f\n",merr);
+ ok = 0;
+ }
+#ifdef DIAG
+ printf("Inversion check complete, peak error = %f\n",merr);
+#endif /* !NEVER */
+
+#endif /* DIAG */
+ printf("\n");
+ printf("Range of Jab values was:\n");
+ printf("J: %f -> %f\n", min[0], max[0]);
+ printf("a: %f -> %f\n", min[1], max[1]);
+ printf("b: %f -> %f\n", min[2], max[2]);
+ }
+ }
+#endif /* INVTEST */
+
+#ifdef INVTEST2
+ for (c = 5; c < 6; c++) {
+ cam->set_view(
+ cam,
+ vc_average, /* Enumerated Viewing Condition */
+ white[c], /* Reference/Adapted White XYZ (Y range 0.0 .. 1.0) = D50 */
+ 34.0, /* Adapting/Surround Luminance cd/m^2 */
+ 0.20, /* Relative Luminance of Background to reference white */
+ 0.0, /* Luminance of white in image - not used */
+ 0.01, /* Flare as a fraction of the reference white (Y range 0.0 .. 1.0) */
+ white[c], /* The Flare color coordinates (typically the Ambient color) */
+ USE_HK
+ );
+
+ {
+ int co0, co1, co2; /* (using co[3] triggers compiler bug) */
+ double merr = 0.0;
+ double xyz[3], Jab[3], checkJab[3];
+
+#ifdef NEVER
+ /* Test case */
+#ifdef NEVER
+ Jab[0] = 0.0;
+ Jab[1] = -180.0;
+ Jab[2] = -180.0;
+ cam->cam_to_XYZ(cam, xyz, Jab);
+ cam->XYZ_to_cam(cam, checkJab, xyz);
+#else
+ xyz[0] = 555.365115;
+ xyz[1] = -518.091908;
+ xyz[2] = -1591.731701;
+ cam->XYZ_to_cam(cam, Jab, xyz);
+ cam->cam_to_XYZ(cam, xyz, Jab);
+#endif
+
+#else /* !NEVER */
+
+ /* itterate through Jab space */
+ for (co0 = 0; co0 < TRES; co0++) {
+ Jab[0] = co0/(TRES-1.0);
+ Jab[0] = 0.0 + Jab[0] * 100.0;
+ for (co1 = 0; co1 < TRES; co1++) {
+ Jab[1] = co1/(TRES-1.0);
+ Jab[1] = (Jab[1] - 0.5) * 256.0;
+ for (co2 = 0; co2 < TRES; co2++) {
+ double mxd;
+ Jab[2] = co2/(TRES-1.0);
+ Jab[2] = (Jab[2] - 0.5) * 256.0;
+
+ cam->cam_to_XYZ(cam, xyz, Jab);
+ cam->XYZ_to_cam(cam, checkJab, xyz);
+
+ /* Check the result */
+ mxd = maxdiff(checkJab, Jab);
+ if (_finite(merr) && (_isnan(mxd) || mxd > merr))
+ merr = mxd;
+#ifdef DIAG
+ if (_isnan(mxd) || mxd > 0.1) {
+ printf("\n");
+ printf("Jab = %f %f %f\n",Jab[0], Jab[1], Jab[2]);
+ printf("%f %f %f -> %f %f %f -> %f %f %f [%f]\n",
+ Jab[0],Jab[1],Jab[2], xyz[0],xyz[1],xyz[2],
+ checkJab[0],checkJab[1],checkJab[2],
+ maxdiff(checkJab, Jab));
+ }
+#endif /* DIAG */
+ }
+ }
+ }
+ if (_isnan(merr) || merr > 1.0) {
+ printf("Excessive error in inversion check 2 %f\n",merr);
+ ok = 0;
+ }
+#ifdef DIAG
+ printf("Inversion check 2 complete, peak error = %f\n",merr);
+#endif /* DIAG */
+#endif /* !NEVER */
+ }
+ }
+#endif /* INVTEST2 */
+
+ printf("\n");
+ if (ok == 0) {
+ printf("Cam testing FAILED\n");
+ } else {
+ printf("Cam testing OK\n");
+ }
+
+ return 0;
+}
+