ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/common/spec_rgb.c
Revision: 2.10
Committed: Wed Oct 7 17:11:40 1998 UTC (25 years, 6 months ago) by gwlarson
Content type: text/plain
Branch: MAIN
Changes since 2.9: +4 -2 lines
Log Message:
added conservative measures to make sure clipgamut() kept in bounds

File Contents

# Content
1 /* Copyright (c) 1997 Regents of the University of California */
2
3 #ifndef lint
4 static char SCCSid[] = "$SunId$ LBL";
5 #endif
6
7 /*
8 * Convert colors and spectral ranges.
9 */
10
11 #include "color.h"
12
13 #define CEPS 1e-7 /* color epsilon */
14
15
16 RGBPRIMS stdprims = STDPRIMS; /* standard primary chromaticities */
17
18 COLOR cblack = BLKCOLOR; /* global black color */
19 COLOR cwhite = WHTCOLOR; /* global white color */
20
21 /*
22 * The following table contains the CIE tristimulus integrals
23 * for X, Y, and Z. The table is cumulative, so that
24 * each color coordinate integrates to 1.
25 */
26
27 #define STARTWL 380 /* starting wavelength (nanometers) */
28 #define INCWL 10 /* wavelength increment */
29 #define NINC 40 /* # of values */
30
31 static BYTE chroma[3][NINC] = {
32 { /* X */
33 0, 0, 0, 2, 6, 13, 22, 30, 36, 41,
34 42, 43, 43, 44, 46, 52, 60, 71, 87, 106,
35 128, 153, 178, 200, 219, 233, 243, 249, 252, 254,
36 255, 255, 255, 255, 255, 255, 255, 255, 255, 255
37 }, { /* Y */
38 0, 0, 0, 0, 0, 1, 2, 4, 7, 11,
39 17, 24, 34, 48, 64, 84, 105, 127, 148, 169,
40 188, 205, 220, 232, 240, 246, 250, 253, 254, 255,
41 255, 255, 255, 255, 255, 255, 255, 255, 255, 255
42 }, { /* Z */
43 0, 0, 2, 10, 32, 66, 118, 153, 191, 220,
44 237, 246, 251, 253, 254, 255, 255, 255, 255, 255,
45 255, 255, 255, 255, 255, 255, 255, 255, 255, 255,
46 255, 255, 255, 255, 255, 255, 255, 255, 255, 255
47 }
48 };
49
50 COLORMAT xyz2rgbmat = { /* XYZ to RGB */
51 {(CIE_y_g - CIE_y_b - CIE_x_b*CIE_y_g + CIE_y_b*CIE_x_g)/CIE_C_rD,
52 (CIE_x_b - CIE_x_g - CIE_x_b*CIE_y_g + CIE_x_g*CIE_y_b)/CIE_C_rD,
53 (CIE_x_g*CIE_y_b - CIE_x_b*CIE_y_g)/CIE_C_rD},
54 {(CIE_y_b - CIE_y_r - CIE_y_b*CIE_x_r + CIE_y_r*CIE_x_b)/CIE_C_gD,
55 (CIE_x_r - CIE_x_b - CIE_x_r*CIE_y_b + CIE_x_b*CIE_y_r)/CIE_C_gD,
56 (CIE_x_b*CIE_y_r - CIE_x_r*CIE_y_b)/CIE_C_gD},
57 {(CIE_y_r - CIE_y_g - CIE_y_r*CIE_x_g + CIE_y_g*CIE_x_r)/CIE_C_bD,
58 (CIE_x_g - CIE_x_r - CIE_x_g*CIE_y_r + CIE_x_r*CIE_y_g)/CIE_C_bD,
59 (CIE_x_r*CIE_y_g - CIE_x_g*CIE_y_r)/CIE_C_bD}
60 };
61
62 COLORMAT rgb2xyzmat = { /* RGB to XYZ */
63 {CIE_x_r*CIE_C_rD/CIE_D,CIE_x_g*CIE_C_gD/CIE_D,CIE_x_b*CIE_C_bD/CIE_D},
64 {CIE_y_r*CIE_C_rD/CIE_D,CIE_y_g*CIE_C_gD/CIE_D,CIE_y_b*CIE_C_bD/CIE_D},
65 {(1.-CIE_x_r-CIE_y_r)*CIE_C_rD/CIE_D,
66 (1.-CIE_x_g-CIE_y_g)*CIE_C_gD/CIE_D,
67 (1.-CIE_x_b-CIE_y_b)*CIE_C_bD/CIE_D}
68 };
69
70
71
72 spec_rgb(col, s, e) /* compute RGB color from spectral range */
73 COLOR col;
74 int s, e;
75 {
76 COLOR ciecolor;
77
78 spec_cie(ciecolor, s, e);
79 cie_rgb(col, ciecolor);
80 }
81
82
83 spec_cie(col, s, e) /* compute a color from a spectral range */
84 COLOR col; /* returned color */
85 int s, e; /* starting and ending wavelengths */
86 {
87 register int i, d, r;
88
89 s -= STARTWL;
90 if (s < 0)
91 s = 0;
92
93 e -= STARTWL;
94 if (e <= s) {
95 col[CIEX] = col[CIEY] = col[CIEZ] = 0.0;
96 return;
97 }
98 if (e >= INCWL*(NINC - 1))
99 e = INCWL*(NINC - 1) - 1;
100
101 d = e / INCWL; /* interpolate values */
102 r = e % INCWL;
103 for (i = 0; i < 3; i++)
104 col[i] = chroma[i][d]*(INCWL - r) + chroma[i][d + 1]*r;
105
106 d = s / INCWL;
107 r = s % INCWL;
108 for (i = 0; i < 3; i++)
109 col[i] -= chroma[i][d]*(INCWL - r) - chroma[i][d + 1]*r;
110
111 col[CIEX] = (col[CIEX] + 0.5) * (1./(256*INCWL));
112 col[CIEY] = (col[CIEY] + 0.5) * (1./(256*INCWL));
113 col[CIEZ] = (col[CIEZ] + 0.5) * (1./(256*INCWL));
114 }
115
116
117 cie_rgb(rgb, xyz) /* convert CIE color to standard RGB */
118 COLOR rgb, xyz;
119 {
120 colortrans(rgb, xyz2rgbmat, xyz);
121 clipgamut(rgb, xyz[CIEY], CGAMUT_LOWER, cblack, cwhite);
122 }
123
124
125 int
126 clipgamut(col, brt, gamut, lower, upper) /* clip to gamut cube */
127 COLOR col;
128 double brt;
129 int gamut;
130 COLOR lower, upper;
131 {
132 int rflags = 0;
133 double brtmin, brtmax, v, vv;
134 COLOR cgry;
135 register int i;
136 /* check for no check */
137 if (gamut == 0) return(0);
138 /* check brightness limits */
139 brtmin = 1./3.*(lower[0]+lower[1]+lower[2]);
140 if (gamut & CGAMUT_LOWER && brt < brtmin) {
141 copycolor(col, lower);
142 return(CGAMUT_LOWER);
143 }
144 brtmax = 1./3.*(upper[0]+upper[1]+upper[2]);
145 if (gamut & CGAMUT_UPPER && brt > brtmax) {
146 copycolor(col, upper);
147 return(CGAMUT_UPPER);
148 }
149 /* compute equivalent grey */
150 v = (brt - brtmin)/(brtmax - brtmin);
151 for (i = 0; i < 3; i++)
152 cgry[i] = v*upper[i] + (1.-v)*lower[i];
153 vv = 1.; /* check each limit */
154 for (i = 0; i < 3; i++)
155 if (gamut & CGAMUT_LOWER && col[i] < lower[i]) {
156 v = (lower[i]+CEPS - cgry[i])/(col[i] - cgry[i]);
157 if (v < vv) vv = v;
158 rflags |= CGAMUT_LOWER;
159 } else if (gamut & CGAMUT_UPPER && col[i] > upper[i]) {
160 v = (upper[i]-CEPS - cgry[i])/(col[i] - cgry[i]);
161 if (v < vv) vv = v;
162 rflags |= CGAMUT_UPPER;
163 }
164 if (rflags) /* desaturate to cube face */
165 for (i = 0; i < 3; i++)
166 col[i] = vv*col[i] + (1.-vv)*cgry[i];
167 return(rflags);
168 }
169
170
171 colortrans(c2, mat, c1) /* convert c1 by mat and put into c2 */
172 register COLORMAT mat;
173 register COLOR c1, c2;
174 {
175 COLOR cout;
176
177 cout[0] = mat[0][0]*c1[0] + mat[0][1]*c1[1] + mat[0][2]*c1[2];
178 cout[1] = mat[1][0]*c1[0] + mat[1][1]*c1[1] + mat[1][2]*c1[2];
179 cout[2] = mat[2][0]*c1[0] + mat[2][1]*c1[1] + mat[2][2]*c1[2];
180
181 copycolor(c2, cout);
182 }
183
184
185 multcolormat(m3, m2, m1) /* multiply m1 by m2 and put into m3 */
186 COLORMAT m1, m2, m3; /* m3 can be either m1 or m2 w/o harm */
187 {
188 COLORMAT mt;
189 register int i, j;
190
191 for (i = 0; i < 3; i++)
192 for (j = 0; j < 3; j++)
193 mt[i][j] = m1[i][0]*m2[0][j] +
194 m1[i][1]*m2[1][j] +
195 m1[i][2]*m2[2][j] ;
196 cpcolormat(m3, mt);
197 }
198
199
200 compxyz2rgbmat(mat, pr) /* compute conversion from CIE to RGB space */
201 COLORMAT mat;
202 register RGBPRIMS pr;
203 {
204 double C_rD, C_gD, C_bD;
205
206 if (pr == stdprims) { /* can use xyz2rgbmat */
207 cpcolormat(mat, xyz2rgbmat);
208 return;
209 }
210 C_rD = (1./pr[WHT][CIEY]) *
211 ( pr[WHT][CIEX]*(pr[GRN][CIEY] - pr[BLU][CIEY]) -
212 pr[WHT][CIEY]*(pr[GRN][CIEX] - pr[BLU][CIEX]) +
213 pr[GRN][CIEX]*pr[BLU][CIEY] - pr[BLU][CIEX]*pr[GRN][CIEY] ) ;
214 C_gD = (1./pr[WHT][CIEY]) *
215 ( pr[WHT][CIEX]*(pr[BLU][CIEY] - pr[RED][CIEY]) -
216 pr[WHT][CIEY]*(pr[BLU][CIEX] - pr[RED][CIEX]) -
217 pr[RED][CIEX]*pr[BLU][CIEY] + pr[BLU][CIEX]*pr[RED][CIEY] ) ;
218 C_bD = (1./pr[WHT][CIEY]) *
219 ( pr[WHT][CIEX]*(pr[RED][CIEY] - pr[GRN][CIEY]) -
220 pr[WHT][CIEY]*(pr[RED][CIEX] - pr[GRN][CIEX]) +
221 pr[RED][CIEX]*pr[GRN][CIEY] - pr[GRN][CIEX]*pr[RED][CIEY] ) ;
222
223 mat[0][0] = (pr[GRN][CIEY] - pr[BLU][CIEY] -
224 pr[BLU][CIEX]*pr[GRN][CIEY] +
225 pr[BLU][CIEY]*pr[GRN][CIEX])/C_rD ;
226 mat[0][1] = (pr[BLU][CIEX] - pr[GRN][CIEX] -
227 pr[BLU][CIEX]*pr[GRN][CIEY] +
228 pr[GRN][CIEX]*pr[BLU][CIEY])/C_rD ;
229 mat[0][2] = (pr[GRN][CIEX]*pr[BLU][CIEY] -
230 pr[BLU][CIEX]*pr[GRN][CIEY])/C_rD ;
231 mat[1][0] = (pr[BLU][CIEY] - pr[RED][CIEY] -
232 pr[BLU][CIEY]*pr[RED][CIEX] +
233 pr[RED][CIEY]*pr[BLU][CIEX])/C_gD ;
234 mat[1][1] = (pr[RED][CIEX] - pr[BLU][CIEX] -
235 pr[RED][CIEX]*pr[BLU][CIEY] +
236 pr[BLU][CIEX]*pr[RED][CIEY])/C_gD ;
237 mat[1][2] = (pr[BLU][CIEX]*pr[RED][CIEY] -
238 pr[RED][CIEX]*pr[BLU][CIEY])/C_gD ;
239 mat[2][0] = (pr[RED][CIEY] - pr[GRN][CIEY] -
240 pr[RED][CIEY]*pr[GRN][CIEX] +
241 pr[GRN][CIEY]*pr[RED][CIEX])/C_bD ;
242 mat[2][1] = (pr[GRN][CIEX] - pr[RED][CIEX] -
243 pr[GRN][CIEX]*pr[RED][CIEY] +
244 pr[RED][CIEX]*pr[GRN][CIEY])/C_bD ;
245 mat[2][2] = (pr[RED][CIEX]*pr[GRN][CIEY] -
246 pr[GRN][CIEX]*pr[RED][CIEY])/C_bD ;
247 }
248
249
250 comprgb2xyzmat(mat, pr) /* compute conversion from RGB to CIE space */
251 COLORMAT mat;
252 register RGBPRIMS pr;
253 {
254 double C_rD, C_gD, C_bD, D;
255
256 if (pr == stdprims) { /* can use rgb2xyzmat */
257 cpcolormat(mat, rgb2xyzmat);
258 return;
259 }
260 C_rD = (1./pr[WHT][CIEY]) *
261 ( pr[WHT][CIEX]*(pr[GRN][CIEY] - pr[BLU][CIEY]) -
262 pr[WHT][CIEY]*(pr[GRN][CIEX] - pr[BLU][CIEX]) +
263 pr[GRN][CIEX]*pr[BLU][CIEY] - pr[BLU][CIEX]*pr[GRN][CIEY] ) ;
264 C_gD = (1./pr[WHT][CIEY]) *
265 ( pr[WHT][CIEX]*(pr[BLU][CIEY] - pr[RED][CIEY]) -
266 pr[WHT][CIEY]*(pr[BLU][CIEX] - pr[RED][CIEX]) -
267 pr[RED][CIEX]*pr[BLU][CIEY] + pr[BLU][CIEX]*pr[RED][CIEY] ) ;
268 C_bD = (1./pr[WHT][CIEY]) *
269 ( pr[WHT][CIEX]*(pr[RED][CIEY] - pr[GRN][CIEY]) -
270 pr[WHT][CIEY]*(pr[RED][CIEX] - pr[GRN][CIEX]) +
271 pr[RED][CIEX]*pr[GRN][CIEY] - pr[GRN][CIEX]*pr[RED][CIEY] ) ;
272 D = pr[RED][CIEX]*(pr[GRN][CIEY] - pr[BLU][CIEY]) +
273 pr[GRN][CIEX]*(pr[BLU][CIEY] - pr[RED][CIEY]) +
274 pr[BLU][CIEX]*(pr[RED][CIEY] - pr[GRN][CIEY]) ;
275 mat[0][0] = pr[RED][CIEX]*C_rD/D;
276 mat[0][1] = pr[GRN][CIEX]*C_gD/D;
277 mat[0][2] = pr[BLU][CIEX]*C_bD/D;
278 mat[1][0] = pr[RED][CIEY]*C_rD/D;
279 mat[1][1] = pr[GRN][CIEY]*C_gD/D;
280 mat[1][2] = pr[BLU][CIEY]*C_bD/D;
281 mat[2][0] = (1.-pr[RED][CIEX]-pr[RED][CIEY])*C_rD/D;
282 mat[2][1] = (1.-pr[GRN][CIEX]-pr[GRN][CIEY])*C_gD/D;
283 mat[2][2] = (1.-pr[BLU][CIEX]-pr[BLU][CIEY])*C_bD/D;
284 }
285
286
287 comprgb2rgbmat(mat, pr1, pr2) /* compute conversion from RGB1 to RGB2 */
288 COLORMAT mat;
289 RGBPRIMS pr1, pr2;
290 {
291 COLORMAT pr1toxyz, xyztopr2;
292
293 if (pr1 == pr2) {
294 mat[0][0] = mat[1][1] = mat[2][2] = 1.0;
295 mat[0][1] = mat[0][2] = mat[1][0] =
296 mat[1][2] = mat[2][0] = mat[2][1] = 0.0;
297 return;
298 }
299 comprgb2xyzmat(pr1toxyz, pr1);
300 compxyz2rgbmat(xyztopr2, pr2);
301 /* combine transforms */
302 multcolormat(mat, pr1toxyz, xyztopr2);
303 }