ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/common/color.c
Revision: 2.42
Committed: Tue May 6 20:45:31 2025 UTC (2 months, 2 weeks ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad6R0, HEAD
Changes since 2.41: +2 -2 lines
Log Message:
fix: Fixed spectral resampling in previously untested case

File Contents

# User Rev Content
1 greg 1.1 #ifndef lint
2 greg 2.42 static const char RCSid[] = "$Id: color.c,v 2.41 2025/03/21 20:04:15 greg Exp $";
3 greg 1.1 #endif
4     /*
5     * color.c - routines for color calculations.
6     *
7 greg 2.9 * Externals declared in color.h
8     */
9    
10 greg 2.10 #include "copyright.h"
11 greg 1.1
12     #include <stdio.h>
13 greg 2.9 #include <stdlib.h>
14 greg 2.7 #include <math.h>
15 greg 1.1 #include "color.h"
16 greg 2.14
17     #ifdef getc_unlocked /* avoid horrendous overhead of flockfile */
18 greg 2.15 #undef getc
19     #undef putc
20 greg 2.19 #undef ferror
21 greg 2.14 #define getc getc_unlocked
22     #define putc putc_unlocked
23 greg 2.18 #define ferror ferror_unlocked
24 greg 2.14 #endif
25 greg 1.1
26 greg 2.34 #define MINELEN 17 /* minimum scanline length for encoding */
27 greg 2.6 #define MAXELEN 0x7fff /* maximum scanline length for encoding */
28 greg 1.14 #define MINRUN 4 /* minimum run length */
29    
30 greg 2.40 const float cxponent[256] = { /* exponent look-up */
31     .0f, 2.2958874e-41f, 4.59177481e-41f, 9.18354962e-41f,
32     1.83670992e-40f, 3.67341985e-40f, 7.34683969e-40f, 1.46936794e-39f,
33     2.93873588e-39f, 5.87747175e-39f, 1.17549435e-38f, 2.3509887e-38f,
34     4.7019774e-38f, 9.40395481e-38f, 1.88079096e-37f, 3.76158192e-37f,
35     7.52316385e-37f, 1.50463277e-36f, 3.00926554e-36f, 6.01853108e-36f,
36     1.20370622e-35f, 2.40741243e-35f, 4.81482486e-35f, 9.62964972e-35f,
37     1.92592994e-34f, 3.85185989e-34f, 7.70371978e-34f, 1.54074396e-33f,
38     3.08148791e-33f, 6.16297582e-33f, 1.23259516e-32f, 2.46519033e-32f,
39     4.93038066e-32f, 9.86076132e-32f, 1.97215226e-31f, 3.94430453e-31f,
40     7.88860905e-31f, 1.57772181e-30f, 3.15544362e-30f, 6.31088724e-30f,
41     1.26217745e-29f, 2.5243549e-29f, 5.04870979e-29f, 1.00974196e-28f,
42     2.01948392e-28f, 4.03896783e-28f, 8.07793567e-28f, 1.61558713e-27f,
43     3.23117427e-27f, 6.46234854e-27f, 1.29246971e-26f, 2.58493941e-26f,
44     5.16987883e-26f, 1.03397577e-25f, 2.06795153e-25f, 4.13590306e-25f,
45     8.27180613e-25f, 1.65436123e-24f, 3.30872245e-24f, 6.6174449e-24f,
46     1.32348898e-23f, 2.64697796e-23f, 5.29395592e-23f, 1.05879118e-22f,
47     2.11758237e-22f, 4.23516474e-22f, 8.47032947e-22f, 1.69406589e-21f,
48     3.38813179e-21f, 6.77626358e-21f, 1.35525272e-20f, 2.71050543e-20f,
49     5.42101086e-20f, 1.08420217e-19f, 2.16840434e-19f, 4.33680869e-19f,
50     8.67361738e-19f, 1.73472348e-18f, 3.46944695e-18f, 6.9388939e-18f,
51     1.38777878e-17f, 2.77555756e-17f, 5.55111512e-17f, 1.11022302e-16f,
52     2.22044605e-16f, 4.4408921e-16f, 8.8817842e-16f, 1.77635684e-15f,
53     3.55271368e-15f, 7.10542736e-15f, 1.42108547e-14f, 2.84217094e-14f,
54     5.68434189e-14f, 1.13686838e-13f, 2.27373675e-13f, 4.54747351e-13f,
55     9.09494702e-13f, 1.8189894e-12f, 3.63797881e-12f, 7.27595761e-12f,
56     1.45519152e-11f, 2.91038305e-11f, 5.82076609e-11f, 1.16415322e-10f,
57     2.32830644e-10f, 4.65661287e-10f, 9.31322575e-10f, 1.86264515e-09f,
58     3.7252903e-09f, 7.4505806e-09f, 1.49011612e-08f, 2.98023224e-08f,
59     5.96046448e-08f, 1.1920929e-07f, 2.38418579e-07f, 4.76837158e-07f,
60     9.53674316e-07f, 1.90734863e-06f, 3.81469727e-06f, 7.62939453e-06f,
61     1.52587891e-05f, 3.05175781e-05f, 6.10351562e-05f, 0.000122070312f,
62     0.000244140625f, 0.00048828125f, 0.0009765625f, 0.001953125f,
63     0.00390625f, 0.0078125f, 0.015625f, 0.03125f, 0.0625f, 0.125f,
64     0.25f, 0.5f, 1.f, 2.f, 4.f, 8.f, 16.f, 32.f, 64.f, 128.f, 256.f,
65     512.f, 1024.f, 2048.f, 4096.f, 8192.f, 16384.f, 32768.f, 65536.f,
66     131072.f, 262144.f, 524288.f, 1048576.f, 2097152.f, 4194304.f,
67     8388608.f, 16777216.f, 33554432.f, 67108864.f, 134217728.f,
68     268435456.f, 536870912.f, 1.07374182e+09f, 2.14748365e+09f,
69     4.2949673e+09f, 8.58993459e+09f, 1.71798692e+10f, 3.43597384e+10f,
70     6.87194767e+10f, 1.37438953e+11f, 2.74877907e+11f, 5.49755814e+11f,
71     1.09951163e+12f, 2.19902326e+12f, 4.39804651e+12f, 8.79609302e+12f,
72     1.7592186e+13f, 3.51843721e+13f, 7.03687442e+13f, 1.40737488e+14f,
73     2.81474977e+14f, 5.62949953e+14f, 1.12589991e+15f, 2.25179981e+15f,
74     4.50359963e+15f, 9.00719925e+15f, 1.80143985e+16f, 3.6028797e+16f,
75     7.2057594e+16f, 1.44115188e+17f, 2.88230376e+17f, 5.76460752e+17f,
76     1.1529215e+18f, 2.30584301e+18f, 4.61168602e+18f, 9.22337204e+18f,
77     1.84467441e+19f, 3.68934881e+19f, 7.37869763e+19f, 1.47573953e+20f,
78     2.95147905e+20f, 5.9029581e+20f, 1.18059162e+21f, 2.36118324e+21f,
79     4.72236648e+21f, 9.44473297e+21f, 1.88894659e+22f, 3.77789319e+22f,
80     7.55578637e+22f, 1.51115727e+23f, 3.02231455e+23f, 6.0446291e+23f,
81     1.20892582e+24f, 2.41785164e+24f, 4.83570328e+24f, 9.67140656e+24f,
82     1.93428131e+25f, 3.86856262e+25f, 7.73712525e+25f, 1.54742505e+26f,
83     3.0948501e+26f, 6.1897002e+26f, 1.23794004e+27f, 2.47588008e+27f,
84     4.95176016e+27f, 9.90352031e+27f, 1.98070406e+28f, 3.96140813e+28f,
85     7.92281625e+28f, 1.58456325e+29f, 3.1691265e+29f, 6.338253e+29f,
86     1.2676506e+30f, 2.5353012e+30f, 5.0706024e+30f, 1.01412048e+31f,
87     2.02824096e+31f, 4.05648192e+31f, 8.11296384e+31f, 1.62259277e+32f,
88     3.24518554e+32f, 6.49037107e+32f, 1.29807421e+33f, 2.59614843e+33f,
89     5.19229686e+33f, 1.03845937e+34f, 2.07691874e+34f, 4.15383749e+34f,
90     8.30767497e+34f, 1.66153499e+35f, 3.32306999e+35f, 6.64613998e+35f
91     };
92 greg 2.4
93 greg 2.27 int CNDX[4] = {0,1,2,3}; /* RGBE indices for SCOLOR, SCOLR */
94     float WLPART[4] = {780,588,480,380}; /* RGB wavelength limits+partitions (nm) */
95    
96    
97     int
98     setspectrsamp( /* assign spectral sampling, 1 if good, -1 if bad */
99     int cn[4], /* input cn[3]=nsamps */
100     float wlpt[4] /* input wlpt[0],wlpt[3]=extrema */
101     )
102     {
103     static const float PKWL[3] = {607, 553, 469};
104     int i, j;
105    
106     if (cn[3] < 3)
107     return(-1); /* reject this */
108    
109     if (wlpt[0] < wlpt[3]) {
110     float tf = wlpt[0];
111     wlpt[0] = wlpt[3]; wlpt[3] = tf;
112     }
113     if (wlpt[0] - wlpt[3] < 50.f)
114     return(-1); /* also reject */
115    
116     if (cn[3] > MAXCSAMP)
117     cn[3] = MAXCSAMP;
118    
119     if ((wlpt[3] >= PKWL[2]) | (wlpt[0] <= PKWL[0])) {
120     wlpt[1] = wlpt[0] + 0.333333f*(wlpt[3]-wlpt[0]);
121     wlpt[2] = wlpt[0] + 0.666667f*(wlpt[3]-wlpt[0]);
122     cn[0] = 0; cn[1] = cn[3]/3; cn[2] = cn[3]*2/3;
123     return(0); /* unhappy but non-fatal return value */
124     }
125     wlpt[1] = 588.f; /* tuned for standard green channel */
126     wlpt[2] = 480.f;
127     if (cn[3] == 3) { /* nothing to tune? */
128     cn[0] = 0; cn[1] = 1; cn[2] = 2;
129     } else { /* else find nearest color indices */
130     double curwl[3];
131     memset(curwl, 0, sizeof(curwl));
132     for (i = cn[3]; i--; ) {
133     const float cwl = (i+.5f)/cn[3]*(wlpt[3]-wlpt[0]) + wlpt[0];
134     for (j = 3; j--; )
135     if (fabs(cwl - PKWL[j]) < fabs(curwl[j] - PKWL[j])) {
136     curwl[j] = cwl;
137     cn[j] = i;
138     }
139     }
140     }
141     return(1); /* happy return value */
142     }
143    
144    
145     void
146     setscolor( /* assign spectral color from RGB */
147     SCOLOR scol,
148     double r,
149     double g,
150     double b
151     )
152     {
153 greg 2.39 double step, cwl;
154 greg 2.27 int i;
155    
156 greg 2.39 if (NCSAMP == 3) {
157     setcolor(scol, r, g, b);
158     return;
159     }
160     step = (WLPART[3] - WLPART[0])/(double)NCSAMP;
161     cwl = WLPART[0] + .5*step;
162 greg 2.27 for (i = 0; i < NCSAMP; i++) {
163     if (cwl >= WLPART[1])
164     scol[i] = r;
165     else if (cwl >= WLPART[2])
166     scol[i] = g;
167     else
168     scol[i] = b;
169     cwl += step;
170     }
171     }
172    
173    
174     void
175 greg 2.38 setscolr( /* assign common-exponent spectral from RGB */
176     SCOLR sclr,
177     double r,
178     double g,
179     double b
180     )
181     {
182     SCOLOR scol;
183    
184     setscolor(scol, r, g, b);
185     scolor2scolr(sclr, scol, NCSAMP);
186     }
187    
188    
189     void
190 greg 2.27 scolor2color( /* assign RGB color from spectrum */
191     COLOR col,
192 greg 2.35 const SCOLOR scol, /* uses average over bands */
193 greg 2.27 int ncs,
194 greg 2.30 const float wlpt[4]
195 greg 2.27 )
196     {
197     const double step = (wlpt[3] - wlpt[0])/(double)ncs;
198     double cwl = wlpt[0] + .5*step;
199     int i, j=0, n=0;
200    
201     setcolor(col, 0, 0, 0);
202     for (i = 0; i < ncs; i++) {
203     if (cwl < wlpt[j+1]) {
204     if (n > 1) col[j] /= (COLORV)n;
205     j++;
206     n = 0;
207     }
208     col[j] += scol[i];
209     n++;
210     cwl += step;
211     }
212     if (n > 1) col[j] /= (COLORV)n;
213     }
214    
215    
216     void
217 greg 2.36 scolr2colr( /* assign RGBE from common-exponent spectrum */
218     COLR clr,
219     const SCOLR sclr,
220     int ncs,
221     const float wlpt[4]
222     )
223 greg 2.37 #if 0 /* fancier method seems to be slower(!) */
224 greg 2.36 {
225     const double step = (wlpt[3] - wlpt[0])/(double)ncs;
226     double cwl;
227     int csum[3], cnt[3], eshft;
228     int i, j;
229    
230     csum[0] = csum[1] = csum[2] = 0;
231     cnt[0] = cnt[1] = cnt[2] = 0;
232     cwl = wlpt[j=0] + .5*step;
233     for (i = 0; i < ncs; i++) {
234     csum[j] += sclr[i];
235     ++cnt[j];
236     j += ((cwl += step) < wlpt[j+1]);
237     }
238     eshft = 7; /* compute exponent shift */
239     for (j = 3; (eshft > 0) & (j-- > 0); ) {
240     i = 0;
241     while (csum[j] < 128*cnt[j] >> i)
242     if (++i >= eshft)
243     break;
244     if (eshft > i)
245     eshft = i;
246     }
247     if (sclr[ncs] <= eshft) {
248     clr[RED] = clr[GRN] = clr[BLU] = 0;
249     clr[EXP] = 0;
250     return;
251     }
252     for (j = 3; j--; )
253     clr[j] = (csum[j]<<eshft)/cnt[j];
254    
255     clr[EXP] = sclr[ncs] - eshft;
256     }
257 greg 2.37 #else
258     {
259     SCOLOR scol;
260     COLOR col;
261    
262     scolr2scolor(scol, sclr, ncs);
263     scolor2color(col, scol, ncs, wlpt);
264     setcolr(clr, col[RED], col[GRN], col[BLU]);
265     }
266     #endif
267 greg 2.36
268    
269     void
270 greg 2.27 scolor2colr( /* assign RGBE from spectral color */
271     COLR clr,
272 greg 2.35 const SCOLOR scol, /* uses average over bands */
273 greg 2.27 int ncs,
274 greg 2.30 const float wlpt[4]
275 greg 2.27 )
276     {
277     COLOR col;
278    
279     scolor2color(col, scol, ncs, wlpt);
280     setcolr(clr, col[RED], col[GRN], col[BLU]);
281     }
282    
283    
284     void
285 greg 2.33 scolr2color( /* assign RGB from common exponent */
286     COLOR col,
287 greg 2.35 const SCOLR sclr,
288 greg 2.33 int ncs,
289     const float wlpt[4]
290     )
291     {
292     SCOLOR scol;
293    
294     scolr2scolor(scol, sclr, ncs);
295     scolor2color(col, scol, ncs, wlpt);
296     }
297    
298    
299     void
300 greg 2.27 scolor2scolr( /* float spectrum to common exponent */
301     SCOLR sclr,
302 greg 2.35 const SCOLOR scol,
303 greg 2.27 int ncs
304     )
305     {
306     int i = ncs;
307     COLORV p = scol[--i];
308    
309     while (i)
310     if (scol[--i] > p)
311     p = scol[i];
312     if (p <= 1e-32) {
313     memset(sclr, 0, ncs+1);
314     return;
315     }
316     p = frexp(p, &i) * 256.0 / p;
317     sclr[ncs] = i + COLXS;
318     for (i = ncs; i--; )
319     sclr[i] = (scol[i] > 0) * (int)(scol[i]*p);
320     }
321    
322    
323     void
324     scolr2scolor( /* common exponent to float spectrum */
325     SCOLOR scol,
326 greg 2.35 const SCOLR sclr,
327 greg 2.27 int ncs
328     )
329     {
330 greg 2.41 const COLORV f = cxponent[sclr[ncs]];
331     int i = ncs;
332 greg 2.27
333 greg 2.41 while (i--)
334     scol[i] = (sclr[i] + (COLORV)0.5)*f;
335 greg 2.27 }
336    
337    
338     double
339     scolor_mean( /* compute average for spectral color */
340 greg 2.35 const SCOLOR scol
341 greg 2.27 )
342     {
343     int i = NCSAMP;
344     double sum = 0;
345    
346     while (i--)
347     sum += scol[i];
348    
349     return sum/(double)NCSAMP;
350     }
351    
352    
353     double
354     sintens( /* find maximum value from spectrum */
355 greg 2.35 const SCOLOR scol
356 greg 2.27 )
357     {
358     int i = NCSAMP;
359     COLORV peak = scol[--i];
360    
361     while (i)
362     if (scol[--i] > peak)
363     peak = scol[i];
364    
365     return peak;
366     }
367    
368    
369     void
370     convertscolor( /* spectrum conversion, zero-fill ends */
371     SCOLOR dst, /* destination spectrum */
372     int dnc, /* destination # of spectral samples/intervals */
373     double dwl0, /* starting destination wavelength (longer) */
374     double dwl1, /* ending destination wavelength (shorter) */
375     const COLORV src[], /* source spectrum array */
376     int snc,
377     double swl0, /* long/short wavelengths may be reversed */
378     double swl1
379     )
380     {
381     const int sdir = 1 - 2*(swl0 < swl1);
382     const double sstp = (swl1 - swl0)/(double)snc;
383     const double dstp = (dwl1 - dwl0)/(double)dnc;
384     const double rdstp = 1./dstp;
385     int si, ssi, di;
386     double wl;
387    
388     if ((dnc < 3) | (dwl0 <= dwl1) | (dst == src))
389     return; /* invalid destination */
390    
391     if (dnc == snc && (dwl0-swl0)*(dwl0-swl0) + (dwl1-swl1)*(dwl1-swl1) <= .5) {
392     memcpy(dst, src, sizeof(COLORV)*dnc);
393     return; /* same spectral sampling */
394     }
395     memset(dst, 0, sizeof(COLORV)*dnc);
396     /* set starting positions */
397     if ((sdir>0 ? swl0 : swl1) <= dwl0) {
398     if (sdir > 0) {
399     wl = swl0;
400     ssi = 0;
401     } else {
402     wl = swl1;
403     ssi = snc-1;
404     }
405     si = 0;
406     di = (wl - dwl0)*rdstp;
407     } else {
408     wl = dwl0;
409 greg 2.29 if (sdir > 0) {
410     ssi = si = (wl - swl0)/sstp;
411     } else {
412 greg 2.42 si = (wl - swl1)/-sstp;
413 greg 2.29 ssi = snc-1 - si;
414     }
415 greg 2.27 di = 0;
416     }
417     swl0 += (sdir < 0)*sstp;
418     /* step through intervals */
419     while ((si < snc) & (di < dnc)) {
420     double intvl;
421     if (swl0 + (ssi+sdir)*sstp < dwl0 + (di+1)*dstp) {
422     intvl = dwl0 + (di+1)*dstp - wl;
423     dst[di++] += src[ssi]*intvl*rdstp;
424     } else {
425     intvl = swl0 + (ssi+sdir)*sstp - wl;
426     dst[di] += src[ssi]*intvl*rdstp;
427     ssi += sdir;
428     si++;
429     }
430     wl += intvl;
431     }
432     }
433    
434    
435 greg 2.22 void *
436 greg 2.17 tempbuffer( /* get a temporary buffer */
437 greg 2.27 size_t len
438 greg 2.17 )
439 greg 1.14 {
440 greg 2.27 static void *tempbuf = NULL;
441     static size_t tempbuflen = 0;
442 greg 1.14
443 greg 2.22 if (!len) { /* call to free */
444     if (tempbuflen) {
445 greg 2.20 free(tempbuf);
446 greg 2.22 tempbuf = NULL;
447     tempbuflen = 0;
448     }
449     return(NULL);
450 greg 1.14 }
451 greg 2.22 if (len <= tempbuflen) /* big enough already? */
452     return(tempbuf);
453     /* else free & reallocate */
454     if (tempbuflen)
455     free(tempbuf);
456     tempbuf = malloc(len);
457     tempbuflen = len*(tempbuf != NULL);
458 greg 1.14 return(tempbuf);
459     }
460    
461    
462 greg 2.9 int
463 greg 2.17 fwritecolrs( /* write out a colr scanline */
464     COLR *scanline,
465     int len,
466     FILE *fp
467     )
468 greg 1.1 {
469 greg 2.17 int i, j, beg, cnt = 1;
470 greg 1.14 int c2;
471 greg 1.1
472 schorsch 2.13 if ((len < MINELEN) | (len > MAXELEN)) /* OOBs, write out flat */
473 greg 1.14 return(fwrite((char *)scanline,sizeof(COLR),len,fp) - len);
474 greg 2.2 /* put magic header */
475     putc(2, fp);
476 greg 1.14 putc(2, fp);
477     putc(len>>8, fp);
478 greg 2.21 putc(len&0xff, fp);
479 greg 1.14 /* put components seperately */
480     for (i = 0; i < 4; i++) {
481     for (j = 0; j < len; j += cnt) { /* find next run */
482     for (beg = j; beg < len; beg += cnt) {
483 greg 2.20 for (cnt = 1; (cnt < 127) & (beg+cnt < len) &&
484 greg 1.14 scanline[beg+cnt][i] == scanline[beg][i]; cnt++)
485     ;
486     if (cnt >= MINRUN)
487     break; /* long enough */
488 greg 1.1 }
489 greg 2.20 if ((beg-j > 1) & (beg-j < MINRUN)) {
490 greg 1.15 c2 = j+1;
491     while (scanline[c2++][i] == scanline[j][i])
492     if (c2 == beg) { /* short run */
493     putc(128+beg-j, fp);
494     putc(scanline[j][i], fp);
495     j = beg;
496     break;
497     }
498     }
499     while (j < beg) { /* write out non-run */
500 greg 1.14 if ((c2 = beg-j) > 128) c2 = 128;
501     putc(c2, fp);
502     while (c2--)
503     putc(scanline[j++][i], fp);
504     }
505     if (cnt >= MINRUN) { /* write out run */
506     putc(128+cnt, fp);
507     putc(scanline[beg][i], fp);
508     } else
509     cnt = 0;
510     }
511 greg 1.1 }
512     return(ferror(fp) ? -1 : 0);
513     }
514    
515 greg 2.23 /*
516     * An old-format scanline is either a stream of valid RGBE or XYZE real
517     * pixels or at least one real pixel followed by some number of
518     * invalid real pixels of the form (1,1,1,n), where n is a count.
519     * These can themselves be repeated to create a multibyte repeat
520     * count, with the least significant byte first (little-endian order.)
521     * Repeat counts are limited by the size of an int; if a repetition
522     * leads to an overrun, the rest of the the repetition will be
523     * silently ignored.
524     */
525 greg 2.9 static int
526 greg 2.19 oldreadcolrs( /* read in an old-style colr scanline */
527 greg 2.17 COLR *scanline,
528     int len,
529     FILE *fp
530     )
531 greg 2.9 {
532 greg 2.25 int rshift = 0;
533 greg 2.17 int i;
534 greg 2.9
535     while (len > 0) {
536     scanline[0][RED] = getc(fp);
537     scanline[0][GRN] = getc(fp);
538     scanline[0][BLU] = getc(fp);
539 greg 2.18 scanline[0][EXP] = i = getc(fp);
540     if (i == EOF)
541 greg 2.9 return(-1);
542 greg 2.20 if (scanline[0][GRN] == 1 &&
543     (scanline[0][RED] == 1) &
544 greg 2.25 (scanline[0][BLU] == 1)) {
545 greg 2.20 i = scanline[0][EXP] << rshift;
546     while (i--) {
547 greg 2.9 copycolr(scanline[0], scanline[-1]);
548 greg 2.20 if (--len <= 0)
549     return(0);
550 greg 2.9 scanline++;
551     }
552     rshift += 8;
553     } else {
554     scanline++;
555     len--;
556     rshift = 0;
557     }
558     }
559     return(0);
560     }
561    
562 greg 2.23 /*
563     * There are two scanline formats: old and new. The old format
564     * compresses runs of RGBE or XYZE four-byte real pixels; the new
565     * format breaks the pixels into R, G, B, and E lines (or XYZE lines)
566     * which are individually run-length encoded.
567     *
568     * An old-format scanline always begins with a valid real pixel; at
569     * least one of the RGB (or XYZ) values will have its high-order bit
570     * set. A new-format scanline begins with four bytes which are not a
571     * valid real pixel: (2, 2, lenhigh, lenlow) where lenhigh is always
572     * less than 128 and hence never has a high-order bit set.
573     *
574     * A new-format scanline is broken into its RGBE or XYZE components.
575     * Each is output and run-length encoded separately so that a scanline
576     * is broken into four records. In turn, each record is organized
577     * into chunks of up to 128 characters, which begin with a count byte.
578     * If the count byte is greater than 128, the following data byte is
579     * repeated (count-128) times. If not, the count byte is followed by
580     * that many data bytes.
581     */
582 greg 2.9 int
583 greg 2.17 freadcolrs( /* read in an encoded colr scanline */
584     COLR *scanline,
585     int len,
586     FILE *fp
587     )
588 greg 1.1 {
589 greg 2.17 int i, j;
590 greg 2.6 int code, val;
591 greg 1.14 /* determine scanline type */
592 greg 2.26 if (len <= 0)
593     return(0);
594 greg 1.14 if ((i = getc(fp)) == EOF)
595     return(-1);
596 greg 2.26 scanline[0][RED] = i;
597 greg 1.14 scanline[0][GRN] = getc(fp);
598     scanline[0][BLU] = getc(fp);
599     if ((i = getc(fp)) == EOF)
600     return(-1);
601 greg 2.26 if ((scanline[0][RED] != 2) | (scanline[0][GRN] != 2) |
602     (scanline[0][BLU] & 0x80)) {
603 greg 1.14 scanline[0][EXP] = i;
604     return(oldreadcolrs(scanline+1, len-1, fp));
605     }
606     if ((scanline[0][BLU]<<8 | i) != len)
607     return(-1); /* length mismatch! */
608     /* read each component */
609     for (i = 0; i < 4; i++)
610     for (j = 0; j < len; ) {
611     if ((code = getc(fp)) == EOF)
612     return(-1);
613     if (code > 128) { /* run */
614 greg 2.6 code &= 127;
615 greg 2.9 if ((val = getc(fp)) == EOF)
616     return -1;
617 greg 2.16 if (j + code > len)
618     return -1; /* overrun */
619 greg 2.6 while (code--)
620     scanline[j++][i] = val;
621 greg 2.16 } else { /* non-run */
622     if (j + code > len)
623     return -1; /* overrun */
624 greg 2.9 while (code--) {
625     if ((val = getc(fp)) == EOF)
626     return -1;
627     scanline[j++][i] = val;
628     }
629 greg 2.16 }
630 greg 1.14 }
631 greg 1.1 return(0);
632     }
633    
634    
635 greg 2.30 /* read an nc-component common-exponent color scanline */
636     int
637 greg 2.35 freadscolrs(COLRV *scanline, int nc, int len, FILE *fp)
638 greg 2.30 {
639 greg 2.32 if (nc < 3)
640     return(-1);
641     if (nc == 3)
642     return(freadcolrs((COLR *)scanline, len, fp));
643    
644 greg 2.30 if (fread(scanline, nc+1, len, fp) != len)
645     return(-1);
646     return(0);
647     }
648    
649    
650 greg 2.36 /* read nc-component common-exponent color scan and convert to COLR's */
651     int
652     fread2colrs(COLR *scanline, int len, FILE *fp, int nc, const float wlpt[4])
653     {
654     COLRV *sclrscan;
655     int n;
656    
657     if (nc < 3)
658     return(-1);
659     if (nc == 3)
660     return(freadcolrs(scanline, len, fp));
661    
662     sclrscan = (COLRV *)tempbuffer(sizeof(COLRV)*(nc+1)*len);
663     if (sclrscan == NULL || freadscolrs(sclrscan, nc, len, fp) < 0)
664     return(-1);
665     for (n = len; n--; ) {
666     scolr2colr(*scanline++, sclrscan, nc, wlpt);
667     sclrscan += nc+1;
668     }
669     return(0);
670     }
671    
672    
673 greg 2.31 /* write an common-exponent spectral color scanline */
674 greg 2.30 int
675 greg 2.35 fwritescolrs(const COLRV *sscanline, int nc, int len, FILE *fp)
676 greg 2.30 {
677 greg 2.32 if (nc < 3)
678     return(-1);
679     if (nc == 3)
680     return(fwritecolrs((COLR *)sscanline, len, fp));
681    
682 greg 2.31 if (fwrite(sscanline, nc+1, len, fp) != len)
683 greg 2.30 return(-1);
684     return(0);
685     }
686    
687    
688 greg 2.9 int
689 greg 2.32 fwritescan( /* write out an RGB or XYZ scanline */
690 greg 2.17 COLOR *scanline,
691     int len,
692     FILE *fp
693     )
694 greg 1.1 {
695 greg 1.14 COLR *clrscan;
696     int n;
697 greg 2.17 COLR *sp;
698 greg 1.14 /* get scanline buffer */
699     if ((sp = (COLR *)tempbuffer(len*sizeof(COLR))) == NULL)
700     return(-1);
701     clrscan = sp;
702     /* convert scanline */
703     n = len;
704     while (n-- > 0) {
705     setcolr(sp[0], scanline[0][RED],
706 greg 1.1 scanline[0][GRN],
707     scanline[0][BLU]);
708     scanline++;
709 greg 1.14 sp++;
710 greg 1.1 }
711 greg 1.14 return(fwritecolrs(clrscan, len, fp));
712 greg 1.1 }
713    
714    
715 greg 2.9 int
716 greg 2.32 freadscan( /* read in an RGB or XYZ scanline */
717 greg 2.17 COLOR *scanline,
718     int len,
719     FILE *fp
720     )
721 greg 1.1 {
722 greg 2.17 COLR *clrscan;
723 greg 1.14
724     if ((clrscan = (COLR *)tempbuffer(len*sizeof(COLR))) == NULL)
725     return(-1);
726     if (freadcolrs(clrscan, len, fp) < 0)
727     return(-1);
728     /* convert scanline */
729 greg 2.40 while (len-- > 0) {
730     colr_color(scanline[0], clrscan[0]);
731 greg 1.14 scanline++; clrscan++;
732 greg 1.1 }
733     return(0);
734     }
735    
736    
737 greg 2.30 /* read an nc-component color scanline */
738     int
739     freadsscan(COLORV *sscanline, int nc, int len, FILE *fp)
740     {
741 greg 2.35 COLRV *tscn = (COLRV *)tempbuffer((nc+1)*len);
742 greg 2.41 int i, j;
743 greg 2.30
744     if (tscn == NULL || freadscolrs(tscn, nc, len, fp) < 0)
745     return(-1);
746 greg 2.41 for (i = len; i-- > 0; ) { /* inline scolr2scolor() */
747     const COLORV f = cxponent[tscn[nc]];
748     for (j = nc; j--; )
749     sscanline[j] = (tscn[j] + (COLORV)0.5)*f;
750 greg 2.30 sscanline += nc;
751     tscn += nc+1;
752     }
753     return(0);
754     }
755    
756    
757 greg 2.36 /* read an nc-component color scanline and return as RGB */
758     int
759     fread2scan(COLOR *scanline, int len, FILE *fp, int nc, const float wlpt[4])
760     {
761     COLRV *tscn;
762     int i;
763    
764     if (nc < 3)
765     return(-1);
766     if (nc == 3)
767     return(freadscan(scanline, len, fp));
768    
769     tscn = (COLRV *)tempbuffer((nc+1)*len);
770     if (tscn == NULL || freadscolrs(tscn, nc, len, fp) < 0)
771     return(-1);
772     for (i = len; i-- > 0; ) {
773     scolr2color(*scanline++, tscn, nc, wlpt);
774     tscn += nc+1;
775     }
776     return(0);
777     }
778    
779    
780 greg 2.32 /* write an nc-component spectral color scanline */
781 greg 2.30 int
782 greg 2.35 fwritesscan(const COLORV *sscanline, int nc, int len, FILE *fp)
783 greg 2.30 {
784 greg 2.35 COLRV *tscn = (COLRV *)tempbuffer((nc+1)*len);
785 greg 2.30 int i;
786    
787     if (tscn == NULL)
788     return(-1);
789     for (i = 0; i < len; i++) {
790 greg 2.31 scolor2scolr(tscn+i*(nc+1), sscanline, nc);
791     sscanline += nc;
792 greg 2.30 }
793 greg 2.31 return(fwritescolrs(tscn, nc, len, fp));
794 greg 2.30 }
795    
796    
797 greg 2.9 void
798 greg 2.17 setcolr( /* assign a short color value */
799     COLR clr,
800     double r,
801     double g,
802     double b
803     )
804 greg 1.1 {
805     double d;
806     int e;
807    
808     d = r > g ? r : g;
809     if (b > d) d = b;
810    
811 greg 1.4 if (d <= 1e-32) {
812 greg 1.1 clr[RED] = clr[GRN] = clr[BLU] = 0;
813     clr[EXP] = 0;
814     return;
815     }
816    
817 greg 2.21 d = frexp(d, &e) * 256.0 / d;
818 greg 1.1
819 greg 2.27 clr[RED] = (r > 0) * (int)(r*d);
820     clr[GRN] = (g > 0) * (int)(g*d);
821     clr[BLU] = (b > 0) * (int)(b*d);
822 greg 1.1 clr[EXP] = e + COLXS;
823     }
824    
825    
826 greg 2.9 int
827 greg 2.17 bigdiff( /* c1 delta c2 > md? */
828 greg 2.35 const COLOR c1,
829     const COLOR c2,
830 greg 2.17 double md
831     )
832 greg 1.7 {
833 greg 2.17 int i;
834 greg 1.7
835     for (i = 0; i < 3; i++)
836 greg 2.27 if ((colval(c1,i)-colval(c2,i) > md*colval(c2,i)) |
837     (colval(c2,i)-colval(c1,i) > md*colval(c1,i)))
838     return(1);
839     return(0);
840     }
841    
842    
843     int
844     sbigsdiff( /* sc1 delta sc2 > md? */
845 greg 2.35 const SCOLOR c1,
846     const SCOLOR c2,
847 greg 2.27 double md
848     )
849     {
850     int i = NCSAMP;
851    
852     while (i--)
853     if ((c1[i]-c2[i] > md*c2[i]) | (c2[i]-c1[i] > md*c1[i]))
854 greg 1.7 return(1);
855     return(0);
856     }