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

Comparing ray/src/common/spec_rgb.c (file contents):
Revision 2.4 by greg, Thu Jul 7 15:21:25 1994 UTC vs.
Revision 2.13 by greg, Mon Jun 30 19:04:29 2003 UTC

# Line 1 | Line 1
1 /* Copyright (c) 1990 Regents of the University of California */
2
1   #ifndef lint
2 < static char SCCSid[] = "$SunId$ LBL";
2 > static const char       RCSid[] = "$Id$";
3   #endif
6
4   /*
5   * Convert colors and spectral ranges.
6 + * Added von Kries white-balance calculations 10/01 (GW).
7 + *
8 + * Externals declared in color.h
9   */
10  
11 + #include "copyright.h"
12 +
13 + #include <stdio.h>
14 + #include <string.h>
15   #include "color.h"
16  
17 + #define CEPS    1e-4                    /* color epsilon */
18 +
19 + #define CEQ(v1,v2)      ((v1) <= (v2)+CEPS && (v2) <= (v1)+CEPS)
20 +
21 + #define XYEQ(c1,c2)     (CEQ((c1)[CIEX],(c2)[CIEX]) && CEQ((c1)[CIEY],(c2)[CIEY]))
22 +
23 +
24 + RGBPRIMS  stdprims = STDPRIMS;  /* standard primary chromaticities */
25 +
26 + COLOR  cblack = BLKCOLOR;               /* global black color */
27 + COLOR  cwhite = WHTCOLOR;               /* global white color */
28 +
29 + float  xyneu[2] = {1./3., 1./3.};       /* neutral xy chromaticities */
30 +
31   /*
32   *      The following table contains the CIE tristimulus integrals
33   *  for X, Y, and Z.  The table is cumulative, so that
# Line 39 | Line 57 | static BYTE  chroma[3][NINC] = {
57          }
58   };
59  
60 < static float  xyz2rgbmat[3][3] = {      /* XYZ to RGB */
60 > COLORMAT  xyz2rgbmat = {                /* XYZ to RGB (no white balance) */
61          {(CIE_y_g - CIE_y_b - CIE_x_b*CIE_y_g + CIE_y_b*CIE_x_g)/CIE_C_rD,
62           (CIE_x_b - CIE_x_g - CIE_x_b*CIE_y_g + CIE_x_g*CIE_y_b)/CIE_C_rD,
63           (CIE_x_g*CIE_y_b - CIE_x_b*CIE_y_g)/CIE_C_rD},
# Line 51 | Line 69 | static float  xyz2rgbmat[3][3] = {     /* XYZ to RGB */
69           (CIE_x_r*CIE_y_g - CIE_x_g*CIE_y_r)/CIE_C_bD}
70   };
71  
72 < static float  rgb2xyzmat[3][3] = {      /* RGB to XYZ */
72 > COLORMAT  rgb2xyzmat = {                /* RGB to XYZ (no white balance) */
73          {CIE_x_r*CIE_C_rD/CIE_D,CIE_x_g*CIE_C_gD/CIE_D,CIE_x_b*CIE_C_bD/CIE_D},
74          {CIE_y_r*CIE_C_rD/CIE_D,CIE_y_g*CIE_C_gD/CIE_D,CIE_y_b*CIE_C_bD/CIE_D},
75          {(1.-CIE_x_r-CIE_y_r)*CIE_C_rD/CIE_D,
# Line 59 | Line 77 | static float  rgb2xyzmat[3][3] = {     /* RGB to XYZ */
77           (1.-CIE_x_b-CIE_y_b)*CIE_C_bD/CIE_D}
78   };
79  
80 + COLORMAT  vkmat = {             /* Sharp primary matrix */
81 +        { 1.2694, -0.0988, -0.1706},
82 +        {-0.8364,  1.8006,  0.0357},
83 +        { 0.0297, -0.0315,  1.0018}
84 + };
85  
86 + COLORMAT  ivkmat = {            /* inverse Sharp primary matrix */
87 +        { 0.8156,  0.0472,  0.1372},
88 +        { 0.3791,  0.5769,  0.0440},
89 +        {-0.0123,  0.0167,  0.9955}
90 + };
91  
92 +
93 + void
94   spec_rgb(col, s, e)             /* compute RGB color from spectral range */
95   COLOR  col;
96   int  s, e;
# Line 72 | Line 102 | int  s, e;
102   }
103  
104  
105 + void
106   spec_cie(col, s, e)             /* compute a color from a spectral range */
107   COLOR  col;             /* returned color */
108   int  s, e;              /* starting and ending wavelengths */
# Line 84 | Line 115 | int  s, e;             /* starting and ending wavelengths */
115  
116          e -= STARTWL;
117          if (e <= s) {
118 <                col[RED] = col[GRN] = col[BLU] = 0.0;
118 >                col[CIEX] = col[CIEY] = col[CIEZ] = 0.0;
119                  return;
120          }
121          if (e >= INCWL*(NINC - 1))
# Line 100 | Line 131 | int  s, e;             /* starting and ending wavelengths */
131          for (i = 0; i < 3; i++)
132                  col[i] -= chroma[i][d]*(INCWL - r) - chroma[i][d + 1]*r;
133  
134 <        col[RED] = (col[RED] + 0.5) / (256*INCWL);
135 <        col[GRN] = (col[GRN] + 0.5) / (256*INCWL);
136 <        col[BLU] = (col[BLU] + 0.5) / (256*INCWL);
134 >        col[CIEX] = (col[CIEX] + 0.5) * (1./(256*INCWL));
135 >        col[CIEY] = (col[CIEY] + 0.5) * (1./(256*INCWL));
136 >        col[CIEZ] = (col[CIEZ] + 0.5) * (1./(256*INCWL));
137   }
138  
139  
140 < cie_rgb(rgbcolor, ciecolor)             /* convert CIE to RGB */
141 < register COLOR  rgbcolor, ciecolor;
140 > void
141 > cie_rgb(rgb, xyz)               /* convert CIE color to standard RGB */
142 > COLOR   rgb;
143 > COLOR  xyz;
144   {
145 <        register int  i;
145 >        colortrans(rgb, xyz2rgbmat, xyz);
146 >        clipgamut(rgb, xyz[CIEY], CGAMUT_LOWER, cblack, cwhite);
147 > }
148  
149 <        for (i = 0; i < 3; i++) {
150 <                rgbcolor[i] =   xyz2rgbmat[i][0]*ciecolor[0] +
151 <                                xyz2rgbmat[i][1]*ciecolor[1] +
152 <                                xyz2rgbmat[i][2]*ciecolor[2] ;
153 <                if (rgbcolor[i] < 0.0)
154 <                        rgbcolor[i] = 0.0;
149 >
150 > int
151 > clipgamut(col, brt, gamut, lower, upper)        /* clip to gamut cube */
152 > COLOR  col;
153 > double  brt;
154 > int  gamut;
155 > COLOR  lower, upper;
156 > {
157 >        int  rflags = 0;
158 >        double  brtmin, brtmax, v, vv;
159 >        COLOR  cgry;
160 >        register int  i;
161 >                                        /* check for no check */
162 >        if (gamut == 0) return(0);
163 >                                        /* check brightness limits */
164 >        brtmin = 1./3.*(lower[0]+lower[1]+lower[2]);
165 >        if (gamut & CGAMUT_LOWER && brt < brtmin) {
166 >                copycolor(col, lower);
167 >                return(CGAMUT_LOWER);
168          }
169 +        brtmax = 1./3.*(upper[0]+upper[1]+upper[2]);
170 +        if (gamut & CGAMUT_UPPER && brt > brtmax) {
171 +                copycolor(col, upper);
172 +                return(CGAMUT_UPPER);
173 +        }
174 +                                        /* compute equivalent grey */
175 +        v = (brt - brtmin)/(brtmax - brtmin);
176 +        for (i = 0; i < 3; i++)
177 +                cgry[i] = v*upper[i] + (1.-v)*lower[i];
178 +        vv = 1.;                        /* check each limit */
179 +        for (i = 0; i < 3; i++)
180 +                if (gamut & CGAMUT_LOWER && col[i] < lower[i]) {
181 +                        v = (lower[i]+CEPS - cgry[i])/(col[i] - cgry[i]);
182 +                        if (v < vv) vv = v;
183 +                        rflags |= CGAMUT_LOWER;
184 +                } else if (gamut & CGAMUT_UPPER && col[i] > upper[i]) {
185 +                        v = (upper[i]-CEPS - cgry[i])/(col[i] - cgry[i]);
186 +                        if (v < vv) vv = v;
187 +                        rflags |= CGAMUT_UPPER;
188 +                }
189 +        if (rflags)                     /* desaturate to cube face */
190 +                for (i = 0; i < 3; i++)
191 +                        col[i] = vv*col[i] + (1.-vv)*cgry[i];
192 +        return(rflags);
193   }
194  
195  
196 < rgb_cie(ciecolor, rgbcolor)             /* convert RGB to CIE */
197 < register COLOR  ciecolor, rgbcolor;
196 > void
197 > colortrans(c2, mat, c1)         /* convert c1 by mat and put into c2 */
198 > register COLOR  c2;
199 > register COLORMAT  mat;
200 > register COLOR  c1;
201   {
202 <        register int  i;
202 >        COLOR   cout;
203  
204 +        cout[0] = mat[0][0]*c1[0] + mat[0][1]*c1[1] + mat[0][2]*c1[2];
205 +        cout[1] = mat[1][0]*c1[0] + mat[1][1]*c1[1] + mat[1][2]*c1[2];
206 +        cout[2] = mat[2][0]*c1[0] + mat[2][1]*c1[1] + mat[2][2]*c1[2];
207 +
208 +        copycolor(c2, cout);
209 + }
210 +
211 +
212 + void
213 + multcolormat(m3, m2, m1)        /* multiply m1 by m2 and put into m3 */
214 + COLORMAT  m3;                   /* m3 can be either m1 or m2 w/o harm */
215 + COLORMAT  m2, m1;
216 + {
217 +        COLORMAT  mt;
218 +        register int  i, j;
219 +
220          for (i = 0; i < 3; i++)
221 <                ciecolor[i] =   rgb2xyzmat[i][0]*rgbcolor[0] +
222 <                                rgb2xyzmat[i][1]*rgbcolor[1] +
223 <                                rgb2xyzmat[i][2]*rgbcolor[2] ;
221 >                for (j = 0; j < 3; j++)
222 >                        mt[i][j] =      m1[i][0]*m2[0][j] +
223 >                                        m1[i][1]*m2[1][j] +
224 >                                        m1[i][2]*m2[2][j] ;
225 >        cpcolormat(m3, mt);
226 > }
227 >
228 >
229 > void
230 > compxyz2rgbmat(mat, pr)         /* compute conversion from CIE to RGB space */
231 > COLORMAT  mat;
232 > register RGBPRIMS  pr;
233 > {
234 >        double  C_rD, C_gD, C_bD;
235 >
236 >        if (pr == stdprims) {   /* can use xyz2rgbmat */
237 >                cpcolormat(mat, xyz2rgbmat);
238 >                return;
239 >        }
240 >        C_rD = (1./pr[WHT][CIEY]) *
241 >                        ( pr[WHT][CIEX]*(pr[GRN][CIEY] - pr[BLU][CIEY]) -
242 >                          pr[WHT][CIEY]*(pr[GRN][CIEX] - pr[BLU][CIEX]) +
243 >                  pr[GRN][CIEX]*pr[BLU][CIEY] - pr[BLU][CIEX]*pr[GRN][CIEY] ) ;
244 >        C_gD = (1./pr[WHT][CIEY]) *
245 >                        ( pr[WHT][CIEX]*(pr[BLU][CIEY] - pr[RED][CIEY]) -
246 >                          pr[WHT][CIEY]*(pr[BLU][CIEX] - pr[RED][CIEX]) -
247 >                  pr[RED][CIEX]*pr[BLU][CIEY] + pr[BLU][CIEX]*pr[RED][CIEY] ) ;
248 >        C_bD = (1./pr[WHT][CIEY]) *
249 >                        ( pr[WHT][CIEX]*(pr[RED][CIEY] - pr[GRN][CIEY]) -
250 >                          pr[WHT][CIEY]*(pr[RED][CIEX] - pr[GRN][CIEX]) +
251 >                  pr[RED][CIEX]*pr[GRN][CIEY] - pr[GRN][CIEX]*pr[RED][CIEY] ) ;
252 >
253 >        mat[0][0] = (pr[GRN][CIEY] - pr[BLU][CIEY] -
254 >                        pr[BLU][CIEX]*pr[GRN][CIEY] +
255 >                        pr[BLU][CIEY]*pr[GRN][CIEX])/C_rD ;
256 >        mat[0][1] = (pr[BLU][CIEX] - pr[GRN][CIEX] -
257 >                        pr[BLU][CIEX]*pr[GRN][CIEY] +
258 >                        pr[GRN][CIEX]*pr[BLU][CIEY])/C_rD ;
259 >        mat[0][2] = (pr[GRN][CIEX]*pr[BLU][CIEY] -
260 >                        pr[BLU][CIEX]*pr[GRN][CIEY])/C_rD ;
261 >        mat[1][0] = (pr[BLU][CIEY] - pr[RED][CIEY] -
262 >                        pr[BLU][CIEY]*pr[RED][CIEX] +
263 >                        pr[RED][CIEY]*pr[BLU][CIEX])/C_gD ;
264 >        mat[1][1] = (pr[RED][CIEX] - pr[BLU][CIEX] -
265 >                        pr[RED][CIEX]*pr[BLU][CIEY] +
266 >                        pr[BLU][CIEX]*pr[RED][CIEY])/C_gD ;
267 >        mat[1][2] = (pr[BLU][CIEX]*pr[RED][CIEY] -
268 >                        pr[RED][CIEX]*pr[BLU][CIEY])/C_gD ;
269 >        mat[2][0] = (pr[RED][CIEY] - pr[GRN][CIEY] -
270 >                        pr[RED][CIEY]*pr[GRN][CIEX] +
271 >                        pr[GRN][CIEY]*pr[RED][CIEX])/C_bD ;
272 >        mat[2][1] = (pr[GRN][CIEX] - pr[RED][CIEX] -
273 >                        pr[GRN][CIEX]*pr[RED][CIEY] +
274 >                        pr[RED][CIEX]*pr[GRN][CIEY])/C_bD ;
275 >        mat[2][2] = (pr[RED][CIEX]*pr[GRN][CIEY] -
276 >                        pr[GRN][CIEX]*pr[RED][CIEY])/C_bD ;
277 > }
278 >
279 >
280 > void
281 > comprgb2xyzmat(mat, pr)         /* compute conversion from RGB to CIE space */
282 > COLORMAT  mat;
283 > register RGBPRIMS  pr;
284 > {
285 >        double  C_rD, C_gD, C_bD, D;
286 >
287 >        if (pr == stdprims) {   /* can use rgb2xyzmat */
288 >                cpcolormat(mat, rgb2xyzmat);
289 >                return;
290 >        }
291 >        C_rD = (1./pr[WHT][CIEY]) *
292 >                        ( pr[WHT][CIEX]*(pr[GRN][CIEY] - pr[BLU][CIEY]) -
293 >                          pr[WHT][CIEY]*(pr[GRN][CIEX] - pr[BLU][CIEX]) +
294 >                  pr[GRN][CIEX]*pr[BLU][CIEY] - pr[BLU][CIEX]*pr[GRN][CIEY] ) ;
295 >        C_gD = (1./pr[WHT][CIEY]) *
296 >                        ( pr[WHT][CIEX]*(pr[BLU][CIEY] - pr[RED][CIEY]) -
297 >                          pr[WHT][CIEY]*(pr[BLU][CIEX] - pr[RED][CIEX]) -
298 >                  pr[RED][CIEX]*pr[BLU][CIEY] + pr[BLU][CIEX]*pr[RED][CIEY] ) ;
299 >        C_bD = (1./pr[WHT][CIEY]) *
300 >                        ( pr[WHT][CIEX]*(pr[RED][CIEY] - pr[GRN][CIEY]) -
301 >                          pr[WHT][CIEY]*(pr[RED][CIEX] - pr[GRN][CIEX]) +
302 >                  pr[RED][CIEX]*pr[GRN][CIEY] - pr[GRN][CIEX]*pr[RED][CIEY] ) ;
303 >        D = pr[RED][CIEX]*(pr[GRN][CIEY] - pr[BLU][CIEY]) +
304 >                        pr[GRN][CIEX]*(pr[BLU][CIEY] - pr[RED][CIEY]) +
305 >                        pr[BLU][CIEX]*(pr[RED][CIEY] - pr[GRN][CIEY]) ;
306 >        mat[0][0] = pr[RED][CIEX]*C_rD/D;
307 >        mat[0][1] = pr[GRN][CIEX]*C_gD/D;
308 >        mat[0][2] = pr[BLU][CIEX]*C_bD/D;
309 >        mat[1][0] = pr[RED][CIEY]*C_rD/D;
310 >        mat[1][1] = pr[GRN][CIEY]*C_gD/D;
311 >        mat[1][2] = pr[BLU][CIEY]*C_bD/D;
312 >        mat[2][0] = (1.-pr[RED][CIEX]-pr[RED][CIEY])*C_rD/D;
313 >        mat[2][1] = (1.-pr[GRN][CIEX]-pr[GRN][CIEY])*C_gD/D;
314 >        mat[2][2] = (1.-pr[BLU][CIEX]-pr[BLU][CIEY])*C_bD/D;
315 > }
316 >
317 >
318 > void
319 > comprgb2rgbmat(mat, pr1, pr2)   /* compute conversion from RGB1 to RGB2 */
320 > COLORMAT  mat;
321 > RGBPRIMS  pr1, pr2;
322 > {
323 >        COLORMAT  pr1toxyz, xyztopr2;
324 >
325 >        if (pr1 == pr2) {
326 >                mat[0][0] = mat[1][1] = mat[2][2] = 1.0;
327 >                mat[0][1] = mat[0][2] = mat[1][0] =
328 >                mat[1][2] = mat[2][0] = mat[2][1] = 0.0;
329 >                return;
330 >        }
331 >        comprgb2xyzmat(pr1toxyz, pr1);
332 >        compxyz2rgbmat(xyztopr2, pr2);
333 >                                /* combine transforms */
334 >        multcolormat(mat, pr1toxyz, xyztopr2);
335 > }
336 >
337 >
338 > void
339 > compxyzWBmat(mat, wht1, wht2)   /* CIE von Kries transform from wht1 to wht2 */
340 > COLORMAT  mat;
341 > float  wht1[2], wht2[2];
342 > {
343 >        COLOR   cw1, cw2;
344 >        if (XYEQ(wht1,wht2)) {
345 >                mat[0][0] = mat[1][1] = mat[2][2] = 1.0;
346 >                mat[0][1] = mat[0][2] = mat[1][0] =
347 >                mat[1][2] = mat[2][0] = mat[2][1] = 0.0;
348 >                return;
349 >        }
350 >        cw1[RED] = wht1[CIEX]/wht1[CIEY];
351 >        cw1[GRN] = 1.;
352 >        cw1[BLU] = (1. - wht1[CIEX] - wht1[CIEY])/wht1[CIEY];
353 >        colortrans(cw1, vkmat, cw1);
354 >        cw2[RED] = wht2[CIEX]/wht2[CIEY];
355 >        cw2[GRN] = 1.;
356 >        cw2[BLU] = (1. - wht2[CIEX] - wht2[CIEY])/wht2[CIEY];
357 >        colortrans(cw2, vkmat, cw2);
358 >        mat[0][0] = cw2[RED]/cw1[RED];
359 >        mat[1][1] = cw2[GRN]/cw1[GRN];
360 >        mat[2][2] = cw2[BLU]/cw1[BLU];
361 >        mat[0][1] = mat[0][2] = mat[1][0] =
362 >        mat[1][2] = mat[2][0] = mat[2][1] = 0.0;
363 >        multcolormat(mat, vkmat, mat);
364 >        multcolormat(mat, mat, ivkmat);
365 > }
366 >
367 >
368 > void
369 > compxyz2rgbWBmat(mat, pr)       /* von Kries conversion from CIE to RGB space */
370 > COLORMAT  mat;
371 > RGBPRIMS  pr;
372 > {
373 >        COLORMAT        wbmat;
374 >
375 >        compxyz2rgbmat(mat, pr);
376 >        if (XYEQ(pr[WHT],xyneu))
377 >                return;
378 >        compxyzWBmat(wbmat, xyneu, pr[WHT]);
379 >        multcolormat(mat, wbmat, mat);
380 > }
381 >
382 > void
383 > comprgb2xyzWBmat(mat, pr)       /* von Kries conversion from RGB to CIE space */
384 > COLORMAT  mat;
385 > RGBPRIMS  pr;
386 > {
387 >        COLORMAT        wbmat;
388 >        
389 >        comprgb2xyzmat(mat, pr);
390 >        if (XYEQ(pr[WHT],xyneu))
391 >                return;
392 >        compxyzWBmat(wbmat, pr[WHT], xyneu);
393 >        multcolormat(mat, mat, wbmat);
394 > }
395 >
396 > void
397 > comprgb2rgbWBmat(mat, pr1, pr2) /* von Kries conversion from RGB1 to RGB2 */
398 > COLORMAT  mat;
399 > RGBPRIMS  pr1, pr2;
400 > {
401 >        COLORMAT  pr1toxyz, xyztopr2, wbmat;
402 >
403 >        if (pr1 == pr2) {
404 >                mat[0][0] = mat[1][1] = mat[2][2] = 1.0;
405 >                mat[0][1] = mat[0][2] = mat[1][0] =
406 >                mat[1][2] = mat[2][0] = mat[2][1] = 0.0;
407 >                return;
408 >        }
409 >        comprgb2xyzmat(pr1toxyz, pr1);
410 >        compxyzWBmat(wbmat, pr1[WHT], pr2[WHT]);
411 >        compxyz2rgbmat(xyztopr2, pr2);
412 >                                /* combine transforms */
413 >        multcolormat(mat, pr1toxyz, wbmat);
414 >        multcolormat(mat, mat, xyztopr2);
415   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines