ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/macbethcal.c
(Generate patch)

Comparing ray/src/px/macbethcal.c (file contents):
Revision 2.2 by greg, Wed Oct 11 11:57:41 1995 UTC vs.
Revision 2.6 by greg, Tue Oct 17 11:05:17 1995 UTC

# Line 18 | Line 18 | static char SCCSid[] = "$SunId$ LBL";
18   #include "resolu.h"
19   #include "pmap.h"
20  
21 <                                /* MacBeth colors (CIE 1931, absolute white) */
22 <                                /* computed from spectral measurements */
21 >                                /* 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   float   mbxyY[24][3] = {
49                  {0.462, 0.3769, 0.0932961},     /* DarkSkin */
50                  {0.4108, 0.3542, 0.410348},     /* LightSkin */
# Line 50 | Line 75 | float  mbxyY[24][3] = {
75   COLOR   mbRGB[24];              /* MacBeth RGB values */
76  
77   #define NMBNEU          6       /* Number of MacBeth neutral colors */
78 < short   mbneu[NMBNEU] = {23,22,21,20,19,18};
54 < #define NMBMOD          3       /* Number of MacBeth moderate colors */
55 < short   mbmod[NMBMOD] = {1,2,21};
56 < #define NMBSAT          6       /* Number of MacBeth saturated colors */
57 < short   mbsat[NMBSAT] = {14,12,13,15,16,17};
78 > short   mbneu[NMBNEU] = {Black,Neutral35,Neutral5,Neutral65,Neutral8,White};
79  
80 + #define NMBMOD          16      /* Number of MacBeth unsaturated colors */
81 + short   mbmod[NMBMOD] = {
82 +                DarkSkin,LightSkin,BlueSky,Foliage,BlueFlower,BluishGreen,
83 +                PurplishBlue,ModerateRed,YellowGreen,OrangeYellow,
84 +                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 +                Orange,Purple
91 +        };
92 +
93   int     xmax, ymax;             /* input image dimensions */
94   int     bounds[4][2];           /* image coordinates of chart corners */
95   double  imgxfm[3][3];           /* coordinate transformation matrix */
# Line 65 | Line 99 | COLOR  picRGB[24];             /* picture colors */
99   COLOR   bramp[NMBNEU][2];       /* brightness ramp (per primary) */
100   double  solmat[3][3];           /* color mapping matrix */
101  
102 < FILE    *debugfp = NULL;
102 > FILE    *debugfp = NULL;        /* debug output picture */
103   char    *progname;
104  
105   extern char     *malloc();
# Line 126 | Line 160 | char   **argv;
160          init();                         /* initialize */
161          getcolors();                    /* get picture colors */
162          compute();                      /* compute color mapping */
163 +                                        /* print comment */
164 +        printf("{ Color correction file computed by %s }\n", progname);
165 +        printf("{ from scanned MacBetch color chart %s }\n", argv[1]);
166          putmapping();                   /* put out color mapping */
167          putdebug();                     /* put out debug picture */
168          exit(0);
# Line 231 | Line 268 | getcolors()                            /* load in picture colors */
268   bresp(y, x)             /* piecewise linear interpolation of primaries */
269   COLOR   y, x;
270   {
271 +        double  cv[3];
272          register int    i, n;
273  
274          for (i = 0; i < 3; i++) {
275 <                n = NMBNEU;
276 <                while (n > 0 && colval(x,i) < colval(bramp[--n][0],i))
277 <                        ;
278 <                colval(y,i) = ((colval(bramp[n+1][0],i) - colval(x,i)) *
275 >                for (n = 0; n < NMBNEU-2; n++)
276 >                        if (colval(x,i) < colval(bramp[n+1][0],i))
277 >                                break;
278 >                cv[i] = ((colval(bramp[n+1][0],i) - colval(x,i)) *
279                                                  colval(bramp[n][1],i) +
280                                  (colval(x,i) - colval(bramp[n][0],i)) *
281                                                  colval(bramp[n+1][1],i)) /
282                          (colval(bramp[n+1][0],i) - colval(bramp[n][0],i));
283 +                if (cv[i] < 0.) cv[i] = 0.;
284          }
285 +        setcolor(y, cv[0], cv[1], cv[2]);
286   }
287  
288  
# Line 304 | Line 344 | putmapping()                   /* put out color mapping for pcomb -f *
344          }
345                                          /* print color mapping */
346          printf("ro = %g*rn + %g*gn + %g*bn ;\n",
347 <                        solmat[0][0], solmat[1][0], solmat[2][0]);
347 >                        solmat[0][0], solmat[0][1], solmat[0][2]);
348          printf("go = %g*rn + %g*gn + %g*bn ;\n",
349 <                        solmat[0][1], solmat[1][1], solmat[2][1]);
349 >                        solmat[1][0], solmat[1][1], solmat[1][2]);
350          printf("bo = %g*rn + %g*gn + %g*bn ;\n",
351 <                        solmat[0][2], solmat[1][2], solmat[2][2]);
351 >                        solmat[2][0], solmat[2][1], solmat[2][2]);
352   }
353  
354  
355 < compsoln(cin, cout, n)          /* solve 3x3 system */
355 > compsoln(cin, cout, n)          /* solve 3xN system using least-squares */
356   COLOR   cin[], cout[];
357   int     n;
358   {
# Line 320 | Line 360 | int    n;
360          double  mat[3][3], invmat[3][3];
361          double  det;
362          double  colv[3], rowv[3];
363 <        register int    i, j;
363 >        register int    i, j, k;
364  
365 <        if (n != 3) {
365 >        if (n < 3 | n > NMBMOD) {
366                  fprintf(stderr, "%s: inconsistent code!\n", progname);
367                  exit(1);
368          }
369 <        for (i = 0; i < 3; i++)
370 <                for (j = 0; j < 3; j++)
371 <                        mat[i][j] = colval(cin[j],i);
369 >        if (n == 3)
370 >                for (i = 0; i < 3; i++)
371 >                        for (j = 0; j < 3; j++)
372 >                                mat[i][j] = colval(cin[j],i);
373 >        else {                          /* compute A^t A */
374 >                for (i = 0; i < 3; i++)
375 >                        for (j = i; j < 3; j++) {
376 >                                mat[i][j] = 0.;
377 >                                for (k = 0; k < n; k++)
378 >                                        mat[i][j] += colval(cin[k],i) *
379 >                                                        colval(cin[k],j);
380 >                        }
381 >                for (i = 1; i < 3; i++)         /* using symmetry */
382 >                        for (j = 0; j < i; j++)
383 >                                mat[i][j] = mat[j][i];
384 >        }
385          det = mx3d_adjoint(mat, invmat);
386          if (fabs(det) < 1e-4) {
387                  fprintf(stderr, "%s: cannot compute color mapping\n",
# Line 342 | Line 395 | int    n;
395                  for (j = 0; j < 3; j++)
396                          invmat[i][j] /= det;
397          for (i = 0; i < 3; i++) {
398 +                if (n == 3)
399 +                        for (j = 0; j < 3; j++)
400 +                                colv[j] = colval(cout[j],i);
401 +                else
402 +                        for (j = 0; j < 3; j++) {
403 +                                colv[j] = 0.;
404 +                                for (k = 0; k < n; k++)
405 +                                        colv[j] += colval(cout[k],i) *
406 +                                                        colval(cin[k],j);
407 +                        }
408 +                mx3d_transform(colv, invmat, rowv);
409                  for (j = 0; j < 3; j++)
410 <                        rowv[j] = colval(cout[j],i);
347 <                mx3d_transform(rowv, invmat, colv);
348 <                for (j = 0; j < 3; j++)
349 <                        solmat[j][i] = colv[j];
410 >                        solmat[i][j] = rowv[j];
411          }
412   }
413  
# Line 357 | Line 418 | COLOR  cout, cin;
418          double  r, g, b;
419  
420          bresp(cout, cin);
421 <        r = colval(cout,0)*solmat[0][0] + colval(cout,1)*solmat[1][0]
422 <                        + colval(cout,2)*solmat[2][0];
421 >        r = colval(cout,0)*solmat[0][0] + colval(cout,1)*solmat[0][1]
422 >                        + colval(cout,2)*solmat[0][2];
423          if (r < 0) r = 0;
424 <        g = colval(cout,0)*solmat[0][1] + colval(cout,1)*solmat[1][1]
425 <                        + colval(cout,2)*solmat[2][1];
424 >        g = colval(cout,0)*solmat[1][0] + colval(cout,1)*solmat[1][1]
425 >                        + colval(cout,2)*solmat[1][2];
426          if (g < 0) g = 0;
427 <        b = colval(cout,0)*solmat[0][2] + colval(cout,1)*solmat[1][2]
427 >        b = colval(cout,0)*solmat[2][0] + colval(cout,1)*solmat[2][1]
428                          + colval(cout,2)*solmat[2][2];
429          if (b < 0) b = 0;
430          setcolor(cout, r, g, b);
# Line 391 | Line 452 | putdebug()                     /* put out debugging picture */
452                  exit(1);
453          }
454                                                  /* finish debug header */
455 +        fputformat(COLRFMT, debugfp);
456          putc('\n', debugfp);
457          fprtresolu(xmax, ymax, debugfp);
458          for (y = ymax-1; y >= 0; y--) {

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines