| 1 | greg | 2.1 | /* Copyright (c) 1995 Regents of the University of California */ | 
| 2 |  |  |  | 
| 3 |  |  | #ifndef lint | 
| 4 |  |  | static char SCCSid[] = "$SunId$ LBL"; | 
| 5 |  |  | #endif | 
| 6 |  |  |  | 
| 7 |  |  | /* | 
| 8 |  |  | * Calibrate a scanned MacBeth Color Checker Chart | 
| 9 |  |  | * | 
| 10 |  |  | * Produce a .cal file suitable for use with pcomb. | 
| 11 |  |  | */ | 
| 12 |  |  |  | 
| 13 |  |  | #include <stdio.h> | 
| 14 |  |  | #ifdef MSDOS | 
| 15 |  |  | #include <fcntl.h> | 
| 16 |  |  | #endif | 
| 17 |  |  | #include "color.h" | 
| 18 |  |  | #include "resolu.h" | 
| 19 |  |  | #include "pmap.h" | 
| 20 |  |  |  | 
| 21 | greg | 2.4 | /* MacBeth colors */ | 
| 22 |  |  | #define DarkSkin        0 | 
| 23 |  |  | #define LightSkin       1 | 
| 24 |  |  | #define BlueSky         2 | 
| 25 |  |  | #define Foliage         3 | 
| 26 |  |  | #define BlueFlower      4 | 
| 27 |  |  | #define BluishGreen     5 | 
| 28 |  |  | #define Orange          6 | 
| 29 |  |  | #define PurplishBlue    7 | 
| 30 |  |  | #define ModerateRed     8 | 
| 31 |  |  | #define Purple          9 | 
| 32 |  |  | #define YellowGreen     10 | 
| 33 |  |  | #define OrangeYellow    11 | 
| 34 |  |  | #define Blue            12 | 
| 35 |  |  | #define Green           13 | 
| 36 |  |  | #define Red             14 | 
| 37 |  |  | #define Yellow          15 | 
| 38 |  |  | #define Magenta         16 | 
| 39 |  |  | #define Cyan            17 | 
| 40 |  |  | #define White           18 | 
| 41 |  |  | #define Neutral8        19 | 
| 42 |  |  | #define Neutral65       20 | 
| 43 |  |  | #define Neutral5        21 | 
| 44 |  |  | #define Neutral35       22 | 
| 45 |  |  | #define Black           23 | 
| 46 |  |  | /* computed from 5nm spectral measurements */ | 
| 47 |  |  | /* CIE 1931 2 degree obs, equal-energy white */ | 
| 48 | greg | 2.1 | float   mbxyY[24][3] = { | 
| 49 |  |  | {0.462, 0.3769, 0.0932961},     /* DarkSkin */ | 
| 50 |  |  | {0.4108, 0.3542, 0.410348},     /* LightSkin */ | 
| 51 |  |  | {0.2626, 0.267, 0.181554},      /* BlueSky */ | 
| 52 |  |  | {0.36, 0.4689, 0.108447},       /* Foliage */ | 
| 53 |  |  | {0.2977, 0.2602, 0.248407},     /* BlueFlower */ | 
| 54 |  |  | {0.2719, 0.3485, 0.401156},     /* BluishGreen */ | 
| 55 |  |  | {0.52, 0.4197, 0.357899},       /* Orange */ | 
| 56 |  |  | {0.229, 0.1866, 0.103911},      /* PurplishBlue */ | 
| 57 |  |  | {0.4909, 0.3262, 0.242615},     /* ModerateRed */ | 
| 58 |  |  | {0.3361, 0.2249, 0.0600102},    /* Purple */ | 
| 59 |  |  | {0.3855, 0.4874, 0.42963},      /* YellowGreen */ | 
| 60 |  |  | {0.4853, 0.4457, 0.476343},     /* OrangeYellow */ | 
| 61 |  |  | {0.2026, 0.1369, 0.0529249},    /* Blue */ | 
| 62 |  |  | {0.3007, 0.4822, 0.221226},     /* Green */ | 
| 63 |  |  | {0.5805, 0.3238, 0.162167},     /* Red */ | 
| 64 |  |  | {0.4617, 0.472, 0.64909},       /* Yellow */ | 
| 65 |  |  | {0.4178, 0.2625, 0.233662},     /* Magenta */ | 
| 66 |  |  | {0.2038, 0.2508, 0.167275},     /* Cyan */ | 
| 67 |  |  | {0.3358, 0.337, 0.916877},      /* White */ | 
| 68 |  |  | {0.3338, 0.3348, 0.604678},     /* Neutral.8 */ | 
| 69 |  |  | {0.3333, 0.3349, 0.364566},     /* Neutral.65 */ | 
| 70 |  |  | {0.3353, 0.3359, 0.200238},     /* Neutral.5 */ | 
| 71 |  |  | {0.3363, 0.336, 0.0878721},     /* Neutral.35 */ | 
| 72 |  |  | {0.3346, 0.3349, 0.0308383}     /* Black */ | 
| 73 |  |  | }; | 
| 74 |  |  |  | 
| 75 |  |  | COLOR   mbRGB[24];              /* MacBeth RGB values */ | 
| 76 |  |  |  | 
| 77 |  |  | #define NMBNEU          6       /* Number of MacBeth neutral colors */ | 
| 78 | greg | 2.4 | short   mbneu[NMBNEU] = {Black,Neutral35,Neutral5,Neutral65,Neutral8,White}; | 
| 79 | greg | 2.1 |  | 
| 80 | greg | 2.6 | #define NMBMOD          16      /* Number of MacBeth unsaturated colors */ | 
| 81 | greg | 2.4 | short   mbmod[NMBMOD] = { | 
| 82 |  |  | DarkSkin,LightSkin,BlueSky,Foliage,BlueFlower,BluishGreen, | 
| 83 | greg | 2.6 | PurplishBlue,ModerateRed,YellowGreen,OrangeYellow, | 
| 84 | greg | 2.4 | Black,Neutral35,Neutral5,Neutral65,Neutral8,White | 
| 85 |  |  | }; | 
| 86 |  |  |  | 
| 87 |  |  | #define NMBSAT          8       /* Number of MacBeth saturated colors */ | 
| 88 |  |  | short   mbsat[NMBSAT] = { | 
| 89 |  |  | Red,Green,Blue,Magenta,Yellow,Cyan, | 
| 90 | greg | 2.6 | Orange,Purple | 
| 91 | greg | 2.4 | }; | 
| 92 |  |  |  | 
| 93 | greg | 2.7 | #define  REQFLGS        (1L<<White|1L<<Neutral8|1L<<Neutral65| \ | 
| 94 |  |  | 1L<<Neutral5|1L<<Neutral35|1L<<Black) | 
| 95 |  |  |  | 
| 96 |  |  | #define  CENTCVG        0.3     /* measured coverage of square sample */ | 
| 97 |  |  | #define  FULLCVG        0.9     /* coverage of entire square */ | 
| 98 |  |  |  | 
| 99 | greg | 2.1 | int     xmax, ymax;             /* input image dimensions */ | 
| 100 |  |  | int     bounds[4][2];           /* image coordinates of chart corners */ | 
| 101 |  |  | double  imgxfm[3][3];           /* coordinate transformation matrix */ | 
| 102 |  |  |  | 
| 103 | greg | 2.7 | COLOR   inpRGB[24];             /* measured or scanned input colors */ | 
| 104 |  |  | long    inpflags = 0;           /* flags of which colors were input */ | 
| 105 | greg | 2.1 |  | 
| 106 | greg | 2.2 | COLOR   bramp[NMBNEU][2];       /* brightness ramp (per primary) */ | 
| 107 | greg | 2.1 | double  solmat[3][3];           /* color mapping matrix */ | 
| 108 |  |  |  | 
| 109 | greg | 2.4 | FILE    *debugfp = NULL;        /* debug output picture */ | 
| 110 | greg | 2.1 | char    *progname; | 
| 111 |  |  |  | 
| 112 |  |  | extern char     *malloc(); | 
| 113 |  |  |  | 
| 114 |  |  |  | 
| 115 |  |  | main(argc, argv) | 
| 116 |  |  | int     argc; | 
| 117 |  |  | char    **argv; | 
| 118 |  |  | { | 
| 119 | greg | 2.7 | int     inpispic = 1; | 
| 120 | greg | 2.1 | int     i; | 
| 121 |  |  |  | 
| 122 |  |  | progname = argv[0]; | 
| 123 | greg | 2.7 | for (i = 1; i < argc && argv[i][0] == '-'; i++) | 
| 124 |  |  | switch (argv[i][1]) { | 
| 125 |  |  | case 'd':                               /* debug output */ | 
| 126 |  |  | i++; | 
| 127 |  |  | if (badarg(argc-i, argv+i, "s")) | 
| 128 |  |  | goto userr; | 
| 129 |  |  | if ((debugfp = fopen(argv[i], "w")) == NULL) { | 
| 130 |  |  | perror(argv[i]); | 
| 131 |  |  | exit(1); | 
| 132 |  |  | } | 
| 133 | greg | 2.1 | #ifdef MSDOS | 
| 134 | greg | 2.7 | setmode(fileno(debugfp), O_BINARY); | 
| 135 | greg | 2.1 | #endif | 
| 136 | greg | 2.7 | newheader("RADIANCE", debugfp);         /* start */ | 
| 137 |  |  | printargs(argc, argv, debugfp);         /* header */ | 
| 138 |  |  | break; | 
| 139 |  |  | case 'p':                               /* picture position */ | 
| 140 |  |  | if (badarg(argc-i-1, argv+i+1, "iiiiiiii")) | 
| 141 |  |  | goto userr; | 
| 142 |  |  | bounds[0][0] = atoi(argv[++i]); | 
| 143 |  |  | bounds[0][1] = atoi(argv[++i]); | 
| 144 |  |  | bounds[1][0] = atoi(argv[++i]); | 
| 145 |  |  | bounds[1][1] = atoi(argv[++i]); | 
| 146 |  |  | bounds[2][0] = atoi(argv[++i]); | 
| 147 |  |  | bounds[2][1] = atoi(argv[++i]); | 
| 148 |  |  | bounds[3][0] = atoi(argv[++i]); | 
| 149 |  |  | bounds[3][1] = atoi(argv[++i]); | 
| 150 |  |  | inpispic = 2; | 
| 151 |  |  | break; | 
| 152 |  |  | case 'c':                               /* color input */ | 
| 153 |  |  | inpispic = 0; | 
| 154 |  |  | break; | 
| 155 |  |  | default: | 
| 156 |  |  | goto userr; | 
| 157 |  |  | } | 
| 158 |  |  | /* open files */ | 
| 159 |  |  | if (i < argc && freopen(argv[i], "r", stdin) == NULL) { | 
| 160 | greg | 2.1 | perror(argv[1]); | 
| 161 |  |  | exit(1); | 
| 162 |  |  | } | 
| 163 | greg | 2.7 | if (i+1 < argc && freopen(argv[i+1], "w", stdout) == NULL) { | 
| 164 | greg | 2.1 | perror(argv[2]); | 
| 165 |  |  | exit(1); | 
| 166 |  |  | } | 
| 167 | greg | 2.7 | if (inpispic) {                 /* load input picture header */ | 
| 168 | greg | 2.1 | #ifdef MSDOS | 
| 169 | greg | 2.7 | setmode(fileno(stdin), O_BINARY); | 
| 170 | greg | 2.1 | #endif | 
| 171 | greg | 2.7 | if (checkheader(stdin, COLRFMT, NULL) < 0 || | 
| 172 |  |  | fgetresolu(&xmax, &ymax, stdin) < 0) { | 
| 173 |  |  | fprintf(stderr, "%s: bad input picture\n", progname); | 
| 174 |  |  | exit(1); | 
| 175 |  |  | } | 
| 176 |  |  | } else {                        /* else set default xmax and ymax */ | 
| 177 |  |  | xmax = 512; | 
| 178 |  |  | ymax = 2*512/3; | 
| 179 | greg | 2.1 | } | 
| 180 | greg | 2.7 | if (inpispic != 2) {            /* use default boundaries */ | 
| 181 | greg | 2.1 | bounds[0][0] = bounds[2][0] = .029*xmax + .5; | 
| 182 |  |  | bounds[0][1] = bounds[1][1] = .956*ymax + .5; | 
| 183 |  |  | bounds[1][0] = bounds[3][0] = .971*xmax + .5; | 
| 184 |  |  | bounds[2][1] = bounds[3][1] = .056*ymax + .5; | 
| 185 |  |  | } | 
| 186 |  |  | init();                         /* initialize */ | 
| 187 | greg | 2.7 | if (inpispic)                   /* get picture colors */ | 
| 188 |  |  | getpicture(); | 
| 189 |  |  | else | 
| 190 |  |  | getcolors(); | 
| 191 | greg | 2.2 | compute();                      /* compute color mapping */ | 
| 192 | greg | 2.4 | /* print comment */ | 
| 193 | greg | 2.7 | printf("{ Color correction file computed by:\n\t"); | 
| 194 |  |  | printargs(argc, argv, stdout); | 
| 195 |  |  | printf("}\n"); | 
| 196 | greg | 2.1 | putmapping();                   /* put out color mapping */ | 
| 197 | greg | 2.7 | if (debugfp != NULL)            /* put out debug picture */ | 
| 198 |  |  | if (inpispic) | 
| 199 |  |  | picdebug(); | 
| 200 |  |  | else | 
| 201 |  |  | clrdebug(); | 
| 202 | greg | 2.1 | exit(0); | 
| 203 |  |  | userr: | 
| 204 | greg | 2.7 | fprintf(stderr, | 
| 205 |  |  | "Usage: %s [-d dbg.pic][-p xul yul xur yur xll yll xlr ylr] input.pic [output.cal]\n", | 
| 206 | greg | 2.1 | progname); | 
| 207 | greg | 2.7 | fprintf(stderr, "   or: %s [-d dbg.pic] -c [xyY.dat [output.cal]]\n", | 
| 208 |  |  | progname); | 
| 209 | greg | 2.1 | exit(1); | 
| 210 |  |  | } | 
| 211 |  |  |  | 
| 212 |  |  |  | 
| 213 |  |  | init()                          /* initialize */ | 
| 214 |  |  | { | 
| 215 |  |  | double  quad[4][2]; | 
| 216 | greg | 2.7 | register int    i; | 
| 217 | greg | 2.1 | /* make coordinate transformation */ | 
| 218 |  |  | quad[0][0] = bounds[0][0]; | 
| 219 |  |  | quad[0][1] = bounds[0][1]; | 
| 220 |  |  | quad[1][0] = bounds[1][0]; | 
| 221 |  |  | quad[1][1] = bounds[1][1]; | 
| 222 |  |  | quad[2][0] = bounds[3][0]; | 
| 223 |  |  | quad[2][1] = bounds[3][1]; | 
| 224 |  |  | quad[3][0] = bounds[2][0]; | 
| 225 |  |  | quad[3][1] = bounds[2][1]; | 
| 226 |  |  |  | 
| 227 |  |  | if (pmap_quad_rect(0., 0., 6., 4., quad, imgxfm) == PMAP_BAD) { | 
| 228 |  |  | fprintf(stderr, "%s: bad chart boundaries\n", progname); | 
| 229 |  |  | exit(1); | 
| 230 |  |  | } | 
| 231 | greg | 2.7 | /* map MacBeth colors to RGB space */ | 
| 232 |  |  | for (i = 0; i < 24; i++) | 
| 233 |  |  | xyY2RGB(mbRGB[i], mbxyY[i]); | 
| 234 | greg | 2.1 | } | 
| 235 |  |  |  | 
| 236 |  |  |  | 
| 237 |  |  | int | 
| 238 | greg | 2.7 | chartndx(x, y, cvg)                     /* find color number for position */ | 
| 239 | greg | 2.1 | int     x, y; | 
| 240 | greg | 2.7 | double  cvg; | 
| 241 | greg | 2.1 | { | 
| 242 |  |  | double  ipos[3], cpos[3]; | 
| 243 |  |  | int     ix, iy; | 
| 244 |  |  | double  fx, fy; | 
| 245 | greg | 2.7 | double  cmin, cmax; | 
| 246 | greg | 2.1 |  | 
| 247 |  |  | ipos[0] = x; | 
| 248 |  |  | ipos[1] = y; | 
| 249 |  |  | ipos[2] = 1; | 
| 250 |  |  | mx3d_transform(ipos, imgxfm, cpos); | 
| 251 |  |  | cpos[0] /= cpos[2]; | 
| 252 |  |  | cpos[1] /= cpos[2]; | 
| 253 |  |  | if (cpos[0] < 0. || cpos[0] >= 6. || cpos[1] < 0. || cpos[1] >= 4.) | 
| 254 |  |  | return(-1); | 
| 255 |  |  | ix = cpos[0]; | 
| 256 |  |  | iy = cpos[1]; | 
| 257 |  |  | fx = cpos[0] - ix; | 
| 258 |  |  | fy = cpos[1] - iy; | 
| 259 | greg | 2.7 | cmin = .5*(1.-cvg); | 
| 260 |  |  | cmax = 1. - cmin; | 
| 261 |  |  | if (fx < cmin || fx >= cmax || fy < cmin || fy >= cmax) | 
| 262 | greg | 2.1 | return(-1); | 
| 263 |  |  | return(iy*6 + ix); | 
| 264 |  |  | } | 
| 265 |  |  |  | 
| 266 |  |  |  | 
| 267 | greg | 2.7 | getpicture()                            /* load in picture colors */ | 
| 268 | greg | 2.1 | { | 
| 269 |  |  | COLR    *scanln; | 
| 270 |  |  | COLOR   pval; | 
| 271 |  |  | int     ccount[24]; | 
| 272 |  |  | double  d; | 
| 273 |  |  | int     y; | 
| 274 |  |  | register int    x, i; | 
| 275 |  |  |  | 
| 276 |  |  | scanln = (COLR *)malloc(xmax*sizeof(COLR)); | 
| 277 |  |  | if (scanln == NULL) { | 
| 278 |  |  | perror(progname); | 
| 279 |  |  | exit(1); | 
| 280 |  |  | } | 
| 281 |  |  | for (i = 0; i < 24; i++) { | 
| 282 | greg | 2.7 | setcolor(inpRGB[i], 0., 0., 0.); | 
| 283 | greg | 2.1 | ccount[i] = 0; | 
| 284 |  |  | } | 
| 285 |  |  | for (y = ymax-1; y >= 0; y--) { | 
| 286 |  |  | if (freadcolrs(scanln, xmax, stdin) < 0) { | 
| 287 |  |  | fprintf(stderr, "%s: error reading input picture\n", | 
| 288 |  |  | progname); | 
| 289 |  |  | exit(1); | 
| 290 |  |  | } | 
| 291 |  |  | for (x = 0; x < xmax; x++) { | 
| 292 | greg | 2.7 | i = chartndx(x, y, CENTCVG); | 
| 293 | greg | 2.1 | if (i >= 0) { | 
| 294 |  |  | colr_color(pval, scanln[x]); | 
| 295 | greg | 2.7 | addcolor(inpRGB[i], pval); | 
| 296 | greg | 2.1 | ccount[i]++; | 
| 297 |  |  | } | 
| 298 |  |  | } | 
| 299 |  |  | } | 
| 300 | greg | 2.7 | for (i = 0; i < 24; i++) {              /* compute averages */ | 
| 301 |  |  | if (ccount[i] == 0) | 
| 302 |  |  | continue; | 
| 303 |  |  | d = 1./ccount[i]; | 
| 304 |  |  | scalecolor(inpRGB[i], d); | 
| 305 |  |  | inpflags |= 1L<<i; | 
| 306 |  |  | } | 
| 307 |  |  | free((char *)scanln); | 
| 308 |  |  | } | 
| 309 |  |  |  | 
| 310 |  |  |  | 
| 311 |  |  | getcolors()                     /* get xyY colors from standard input */ | 
| 312 |  |  | { | 
| 313 |  |  | int     gotwhite = 0; | 
| 314 |  |  | COLOR   whiteclr; | 
| 315 |  |  | int     n; | 
| 316 |  |  | float   xyYin[3]; | 
| 317 |  |  |  | 
| 318 |  |  | while (fgetval(stdin, 'i', &n) == 1) {          /* read colors */ | 
| 319 |  |  | if (n < 0 | n > 24 || | 
| 320 |  |  | fgetval(stdin, 'f', &xyYin[0]) != 1 || | 
| 321 |  |  | fgetval(stdin, 'f', &xyYin[1]) != 1 || | 
| 322 |  |  | fgetval(stdin, 'f', &xyYin[2]) != 1 || | 
| 323 |  |  | xyYin[0] < 0. | xyYin[0] > 1. | | 
| 324 |  |  | xyYin[1] < 0. | xyYin[1] > 1.) { | 
| 325 |  |  | fprintf(stderr, "%s: bad color input data\n", | 
| 326 | greg | 2.1 | progname); | 
| 327 |  |  | exit(1); | 
| 328 |  |  | } | 
| 329 | greg | 2.7 | if (n == 0) {                           /* calibration white */ | 
| 330 |  |  | xyY2RGB(whiteclr, xyYin); | 
| 331 |  |  | gotwhite++; | 
| 332 |  |  | } else {                                /* standard color */ | 
| 333 |  |  | n--; | 
| 334 |  |  | xyY2RGB(inpRGB[n], xyYin); | 
| 335 |  |  | inpflags |= 1L<<n; | 
| 336 |  |  | } | 
| 337 | greg | 2.1 | } | 
| 338 | greg | 2.7 | /* normalize colors */ | 
| 339 |  |  | if (!gotwhite) { | 
| 340 |  |  | if (!(inpflags & 1L<<White)) { | 
| 341 |  |  | fprintf(stderr, "%s: missing input for White\n", | 
| 342 |  |  | progname); | 
| 343 |  |  | exit(1); | 
| 344 |  |  | } | 
| 345 |  |  | setcolor(whiteclr, | 
| 346 |  |  | colval(inpRGB[White],RED)/colval(mbRGB[White],RED), | 
| 347 |  |  | colval(inpRGB[White],GRN)/colval(mbRGB[White],GRN), | 
| 348 |  |  | colval(inpRGB[White],BLU)/colval(mbRGB[White],BLU)); | 
| 349 |  |  | } | 
| 350 |  |  | for (n = 0; n < 24; n++) | 
| 351 |  |  | if (inpflags & 1L<<n) | 
| 352 |  |  | setcolor(inpRGB[n], | 
| 353 |  |  | colval(inpRGB[n],RED)/colval(whiteclr,RED), | 
| 354 |  |  | colval(inpRGB[n],GRN)/colval(whiteclr,GRN), | 
| 355 |  |  | colval(inpRGB[n],BLU)/colval(whiteclr,BLU)); | 
| 356 | greg | 2.1 | } | 
| 357 |  |  |  | 
| 358 |  |  |  | 
| 359 | greg | 2.2 | bresp(y, x)             /* piecewise linear interpolation of primaries */ | 
| 360 |  |  | COLOR   y, x; | 
| 361 | greg | 2.1 | { | 
| 362 | greg | 2.3 | double  cv[3]; | 
| 363 | greg | 2.2 | register int    i, n; | 
| 364 | greg | 2.1 |  | 
| 365 | greg | 2.2 | for (i = 0; i < 3; i++) { | 
| 366 | greg | 2.5 | for (n = 0; n < NMBNEU-2; n++) | 
| 367 |  |  | if (colval(x,i) < colval(bramp[n+1][0],i)) | 
| 368 |  |  | break; | 
| 369 | greg | 2.3 | cv[i] = ((colval(bramp[n+1][0],i) - colval(x,i)) * | 
| 370 | greg | 2.2 | colval(bramp[n][1],i) + | 
| 371 |  |  | (colval(x,i) - colval(bramp[n][0],i)) * | 
| 372 |  |  | colval(bramp[n+1][1],i)) / | 
| 373 |  |  | (colval(bramp[n+1][0],i) - colval(bramp[n][0],i)); | 
| 374 | greg | 2.3 | if (cv[i] < 0.) cv[i] = 0.; | 
| 375 | greg | 2.2 | } | 
| 376 | greg | 2.3 | setcolor(y, cv[0], cv[1], cv[2]); | 
| 377 | greg | 2.1 | } | 
| 378 |  |  |  | 
| 379 |  |  |  | 
| 380 | greg | 2.2 | compute()                       /* compute color mapping */ | 
| 381 | greg | 2.1 | { | 
| 382 | greg | 2.2 | COLOR   clrin[NMBMOD], clrout[NMBMOD]; | 
| 383 | greg | 2.7 | register int    i, n; | 
| 384 |  |  | /* did we get what we need? */ | 
| 385 |  |  | if ((inpflags & REQFLGS) != REQFLGS) { | 
| 386 |  |  | fprintf(stderr, "%s: missing required input colors\n", | 
| 387 |  |  | progname); | 
| 388 |  |  | exit(1); | 
| 389 | greg | 2.2 | } | 
| 390 | greg | 2.1 | /* compute piecewise luminance curve */ | 
| 391 |  |  | for (i = 0; i < NMBNEU; i++) { | 
| 392 | greg | 2.7 | copycolor(bramp[i][0], inpRGB[mbneu[i]]); | 
| 393 | greg | 2.2 | copycolor(bramp[i][1], mbRGB[mbneu[i]]); | 
| 394 | greg | 2.1 | } | 
| 395 |  |  | /* compute color matrix */ | 
| 396 | greg | 2.7 | for (n = 0, i = 0; i < NMBMOD; i++) | 
| 397 |  |  | if (inpflags & 1L<<mbmod[i]) { | 
| 398 |  |  | bresp(clrin[n], inpRGB[mbmod[i]]); | 
| 399 |  |  | copycolor(clrout[n], mbRGB[mbmod[i]]); | 
| 400 |  |  | n++; | 
| 401 |  |  | } | 
| 402 |  |  | compsoln(clrin, clrout, n); | 
| 403 | greg | 2.2 | } | 
| 404 |  |  |  | 
| 405 |  |  |  | 
| 406 |  |  | putmapping()                    /* put out color mapping for pcomb -f */ | 
| 407 |  |  | { | 
| 408 |  |  | static char     cchar[3] = {'r', 'g', 'b'}; | 
| 409 |  |  | register int    i, j; | 
| 410 | greg | 2.1 | /* print brightness mapping */ | 
| 411 | greg | 2.2 | for (j = 0; j < 3; j++) { | 
| 412 |  |  | printf("%cxa(i) : select(i", cchar[j]); | 
| 413 |  |  | for (i = 0; i < NMBNEU; i++) | 
| 414 |  |  | printf(",%g", colval(bramp[i][0],j)); | 
| 415 |  |  | printf(");\n"); | 
| 416 |  |  | printf("%cya(i) : select(i", cchar[j]); | 
| 417 |  |  | for (i = 0; i < NMBNEU; i++) | 
| 418 |  |  | printf(",%g", colval(bramp[i][1],j)); | 
| 419 |  |  | printf(");\n"); | 
| 420 |  |  | printf("%c = %ci(1);\n", cchar[j], cchar[j]); | 
| 421 |  |  | printf("%cfi(n) = if(n-%g, %d, if(%cxa(n+1)-%c, n, %cfi(n+1)));\n", | 
| 422 |  |  | cchar[j], NMBNEU-1.5, NMBNEU-1, cchar[j], | 
| 423 |  |  | cchar[j], cchar[j]); | 
| 424 |  |  | printf("%cndx = %cfi(1);\n", cchar[j], cchar[j]); | 
| 425 |  |  | printf("%cn = ((%cxa(%cndx+1)-%c)*%cya(%cndx) + ", | 
| 426 |  |  | cchar[j], cchar[j], cchar[j], | 
| 427 |  |  | cchar[j], cchar[j], cchar[j]); | 
| 428 |  |  | printf("(%c-%cxa(%cndx))*%cya(%cndx+1)) /\n", | 
| 429 |  |  | cchar[j], cchar[j], cchar[j], | 
| 430 |  |  | cchar[j], cchar[j]); | 
| 431 |  |  | printf("\t\t(%cxa(%cndx+1) - %cxa(%cndx)) ;\n", | 
| 432 |  |  | cchar[j], cchar[j], cchar[j], cchar[j]); | 
| 433 |  |  | } | 
| 434 | greg | 2.1 | /* print color mapping */ | 
| 435 | greg | 2.2 | printf("ro = %g*rn + %g*gn + %g*bn ;\n", | 
| 436 | greg | 2.3 | solmat[0][0], solmat[0][1], solmat[0][2]); | 
| 437 | greg | 2.2 | printf("go = %g*rn + %g*gn + %g*bn ;\n", | 
| 438 | greg | 2.3 | solmat[1][0], solmat[1][1], solmat[1][2]); | 
| 439 | greg | 2.2 | printf("bo = %g*rn + %g*gn + %g*bn ;\n", | 
| 440 | greg | 2.3 | solmat[2][0], solmat[2][1], solmat[2][2]); | 
| 441 | greg | 2.1 | } | 
| 442 |  |  |  | 
| 443 |  |  |  | 
| 444 | greg | 2.4 | compsoln(cin, cout, n)          /* solve 3xN system using least-squares */ | 
| 445 | greg | 2.2 | COLOR   cin[], cout[]; | 
| 446 | greg | 2.1 | int     n; | 
| 447 |  |  | { | 
| 448 |  |  | extern double   mx3d_adjoint(), fabs(); | 
| 449 |  |  | double  mat[3][3], invmat[3][3]; | 
| 450 |  |  | double  det; | 
| 451 |  |  | double  colv[3], rowv[3]; | 
| 452 | greg | 2.4 | register int    i, j, k; | 
| 453 | greg | 2.1 |  | 
| 454 | greg | 2.4 | if (n < 3 | n > NMBMOD) { | 
| 455 | greg | 2.7 | fprintf(stderr, "%s: bad number of colors to match\n", progname); | 
| 456 | greg | 2.1 | exit(1); | 
| 457 |  |  | } | 
| 458 | greg | 2.4 | if (n == 3) | 
| 459 |  |  | for (i = 0; i < 3; i++) | 
| 460 |  |  | for (j = 0; j < 3; j++) | 
| 461 |  |  | mat[i][j] = colval(cin[j],i); | 
| 462 |  |  | else {                          /* compute A^t A */ | 
| 463 |  |  | for (i = 0; i < 3; i++) | 
| 464 |  |  | for (j = i; j < 3; j++) { | 
| 465 |  |  | mat[i][j] = 0.; | 
| 466 |  |  | for (k = 0; k < n; k++) | 
| 467 |  |  | mat[i][j] += colval(cin[k],i) * | 
| 468 |  |  | colval(cin[k],j); | 
| 469 |  |  | } | 
| 470 |  |  | for (i = 1; i < 3; i++)         /* using symmetry */ | 
| 471 |  |  | for (j = 0; j < i; j++) | 
| 472 |  |  | mat[i][j] = mat[j][i]; | 
| 473 |  |  | } | 
| 474 | greg | 2.1 | det = mx3d_adjoint(mat, invmat); | 
| 475 |  |  | if (fabs(det) < 1e-4) { | 
| 476 |  |  | fprintf(stderr, "%s: cannot compute color mapping\n", | 
| 477 |  |  | progname); | 
| 478 |  |  | solmat[0][0] = solmat[1][1] = solmat[2][2] = 1.; | 
| 479 |  |  | solmat[0][1] = solmat[0][2] = solmat[1][0] = | 
| 480 |  |  | solmat[1][2] = solmat[2][0] = solmat[2][1] = 0.; | 
| 481 |  |  | return; | 
| 482 |  |  | } | 
| 483 |  |  | for (i = 0; i < 3; i++) | 
| 484 |  |  | for (j = 0; j < 3; j++) | 
| 485 |  |  | invmat[i][j] /= det; | 
| 486 |  |  | for (i = 0; i < 3; i++) { | 
| 487 | greg | 2.4 | if (n == 3) | 
| 488 |  |  | for (j = 0; j < 3; j++) | 
| 489 |  |  | colv[j] = colval(cout[j],i); | 
| 490 |  |  | else | 
| 491 |  |  | for (j = 0; j < 3; j++) { | 
| 492 |  |  | colv[j] = 0.; | 
| 493 |  |  | for (k = 0; k < n; k++) | 
| 494 |  |  | colv[j] += colval(cout[k],i) * | 
| 495 |  |  | colval(cin[k],j); | 
| 496 |  |  | } | 
| 497 | greg | 2.3 | mx3d_transform(colv, invmat, rowv); | 
| 498 | greg | 2.1 | for (j = 0; j < 3; j++) | 
| 499 | greg | 2.3 | solmat[i][j] = rowv[j]; | 
| 500 | greg | 2.1 | } | 
| 501 |  |  | } | 
| 502 |  |  |  | 
| 503 | greg | 2.3 |  | 
| 504 | greg | 2.1 | cvtcolor(cout, cin)             /* convert color according to our mapping */ | 
| 505 |  |  | COLOR   cout, cin; | 
| 506 |  |  | { | 
| 507 |  |  | double  r, g, b; | 
| 508 |  |  |  | 
| 509 | greg | 2.2 | bresp(cout, cin); | 
| 510 | greg | 2.3 | r = colval(cout,0)*solmat[0][0] + colval(cout,1)*solmat[0][1] | 
| 511 |  |  | + colval(cout,2)*solmat[0][2]; | 
| 512 | greg | 2.2 | if (r < 0) r = 0; | 
| 513 | greg | 2.3 | g = colval(cout,0)*solmat[1][0] + colval(cout,1)*solmat[1][1] | 
| 514 |  |  | + colval(cout,2)*solmat[1][2]; | 
| 515 | greg | 2.2 | if (g < 0) g = 0; | 
| 516 | greg | 2.3 | b = colval(cout,0)*solmat[2][0] + colval(cout,1)*solmat[2][1] | 
| 517 | greg | 2.2 | + colval(cout,2)*solmat[2][2]; | 
| 518 |  |  | if (b < 0) b = 0; | 
| 519 |  |  | setcolor(cout, r, g, b); | 
| 520 | greg | 2.1 | } | 
| 521 |  |  |  | 
| 522 |  |  |  | 
| 523 | greg | 2.7 | xyY2RGB(rgbout, xyYin)          /* convert xyY to RGB */ | 
| 524 |  |  | COLOR   rgbout; | 
| 525 |  |  | register float  xyYin[3]; | 
| 526 | greg | 2.1 | { | 
| 527 | greg | 2.7 | COLOR   ctmp; | 
| 528 |  |  | double  d; | 
| 529 |  |  |  | 
| 530 |  |  | d = xyYin[2] / xyYin[1]; | 
| 531 |  |  | ctmp[0] = xyYin[0] * d; | 
| 532 |  |  | ctmp[1] = xyYin[2]; | 
| 533 |  |  | ctmp[2] = (1. - xyYin[0] - xyYin[1]) * d; | 
| 534 |  |  | cie_rgb(rgbout, ctmp); | 
| 535 |  |  | } | 
| 536 |  |  |  | 
| 537 |  |  |  | 
| 538 |  |  | picdebug()                      /* put out debugging picture */ | 
| 539 |  |  | { | 
| 540 | greg | 2.1 | COLOR   *scan; | 
| 541 |  |  | int     y; | 
| 542 |  |  | register int    x, i; | 
| 543 |  |  |  | 
| 544 |  |  | if (fseek(stdin, 0L, 0) == EOF) { | 
| 545 |  |  | fprintf(stderr, "%s: cannot seek on input picture\n", progname); | 
| 546 |  |  | exit(1); | 
| 547 |  |  | } | 
| 548 |  |  | getheader(stdin, NULL, NULL);           /* skip input header */ | 
| 549 |  |  | fgetresolu(&xmax, &ymax, stdin); | 
| 550 |  |  | /* allocate scanline */ | 
| 551 |  |  | scan = (COLOR *)malloc(xmax*sizeof(COLOR)); | 
| 552 |  |  | if (scan == NULL) { | 
| 553 |  |  | perror(progname); | 
| 554 |  |  | exit(1); | 
| 555 |  |  | } | 
| 556 |  |  | /* finish debug header */ | 
| 557 | greg | 2.5 | fputformat(COLRFMT, debugfp); | 
| 558 | greg | 2.1 | putc('\n', debugfp); | 
| 559 |  |  | fprtresolu(xmax, ymax, debugfp); | 
| 560 | greg | 2.7 | /* write debug picture */ | 
| 561 | greg | 2.1 | for (y = ymax-1; y >= 0; y--) { | 
| 562 |  |  | if (freadscan(scan, xmax, stdin) < 0) { | 
| 563 |  |  | fprintf(stderr, "%s: error rereading input picture\n", | 
| 564 |  |  | progname); | 
| 565 |  |  | exit(1); | 
| 566 |  |  | } | 
| 567 |  |  | for (x = 0; x < xmax; x++) { | 
| 568 | greg | 2.7 | i = chartndx(x, y, CENTCVG); | 
| 569 | greg | 2.1 | if (i < 0) | 
| 570 |  |  | cvtcolor(scan[x], scan[x]); | 
| 571 |  |  | else | 
| 572 |  |  | copycolor(scan[x], mbRGB[i]); | 
| 573 |  |  | } | 
| 574 |  |  | if (fwritescan(scan, xmax, debugfp) < 0) { | 
| 575 |  |  | fprintf(stderr, "%s: error writing debugging picture\n", | 
| 576 |  |  | progname); | 
| 577 |  |  | exit(1); | 
| 578 |  |  | } | 
| 579 |  |  | } | 
| 580 | greg | 2.7 | /* clean up */ | 
| 581 |  |  | fclose(debugfp); | 
| 582 |  |  | free((char *)scan); | 
| 583 |  |  | } | 
| 584 |  |  |  | 
| 585 |  |  |  | 
| 586 |  |  | clrdebug()                      /* put out debug picture from color input */ | 
| 587 |  |  | { | 
| 588 |  |  | static COLR     blkclr = BLKCOLR; | 
| 589 |  |  | COLR    mbclr[24], cvclr[24]; | 
| 590 |  |  | COLR    *scan; | 
| 591 |  |  | COLOR   ctmp; | 
| 592 |  |  | int     y; | 
| 593 |  |  | register int    i, x; | 
| 594 |  |  | /* convert colors */ | 
| 595 |  |  | for (i = 0; i < 24; i++) { | 
| 596 |  |  | setcolr(mbclr[i], colval(mbRGB[i],RED), | 
| 597 |  |  | colval(mbRGB[i],GRN), colval(mbRGB[i],BLU)); | 
| 598 |  |  | if (inpflags & 1L<<i) { | 
| 599 |  |  | cvtcolor(ctmp, inpRGB[i]); | 
| 600 |  |  | setcolr(cvclr[i], colval(ctmp,RED), | 
| 601 |  |  | colval(ctmp,GRN), colval(ctmp,BLU)); | 
| 602 |  |  | } | 
| 603 |  |  | } | 
| 604 |  |  | /* allocate scanline */ | 
| 605 |  |  | scan = (COLR *)malloc(xmax*sizeof(COLR)); | 
| 606 |  |  | if (scan == NULL) { | 
| 607 |  |  | perror(progname); | 
| 608 |  |  | exit(1); | 
| 609 |  |  | } | 
| 610 |  |  | /* finish debug header */ | 
| 611 |  |  | fputformat(COLRFMT, debugfp); | 
| 612 |  |  | putc('\n', debugfp); | 
| 613 |  |  | fprtresolu(xmax, ymax, debugfp); | 
| 614 |  |  | /* write debug picture */ | 
| 615 |  |  | for (y = ymax-1; y >= 0; y--) { | 
| 616 |  |  | for (x = 0; x < xmax; x++) | 
| 617 |  |  | if ((i = chartndx(x, y, CENTCVG)) >= 0) | 
| 618 |  |  | copycolr(scan[x], mbclr[i]); | 
| 619 |  |  | else if ((i = chartndx(x, y, FULLCVG)) >= 0 && | 
| 620 |  |  | inpflags & 1L<<i) | 
| 621 |  |  | copycolr(scan[x], cvclr[i]); | 
| 622 |  |  | else | 
| 623 |  |  | copycolr(scan[x], blkclr); | 
| 624 |  |  | if (fwritecolrs(scan, xmax, debugfp) < 0) { | 
| 625 |  |  | fprintf(stderr, "%s: error writing debugging picture\n", | 
| 626 |  |  | progname); | 
| 627 |  |  | exit(1); | 
| 628 |  |  | } | 
| 629 |  |  | } | 
| 630 |  |  | /* clean up */ | 
| 631 |  |  | fclose(debugfp); | 
| 632 | greg | 2.1 | free((char *)scan); | 
| 633 |  |  | } |