ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/pcond.c
Revision: 3.8
Committed: Thu Jan 9 13:56:23 1997 UTC (27 years, 3 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 3.7: +33 -26 lines
Log Message:
changed to histogram truncation
made veiling reflection modify FOV sample image

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 * Condition Radiance picture for display/output
9 */
10
11 #include "pcond.h"
12
13
14 #define LDMAX 100 /* default max. display luminance */
15 #define LDMINF 0.01 /* default min. display lum. factor */
16
17 int what2do = 0; /* desired adjustments */
18
19 double ldmax = LDMAX; /* maximum output luminance */
20 double ldmin = 0.; /* minimum output luminance */
21 double Bldmin, Bldmax; /* Bl(ldmin) and Bl(ldmax) */
22
23 char *progname; /* global argv[0] */
24
25 char *infn; /* input file name */
26 FILE *infp; /* input stream */
27 FILE *mapfp = NULL; /* tone-mapping function stream */
28 VIEW ourview = STDVIEW; /* picture view */
29 int gotview = 0; /* picture has view */
30 double pixaspect = 1.0; /* pixel aspect ratio */
31 RESOLU inpres; /* input picture resolution */
32
33 COLOR *fovimg; /* foveal (1 degree) averaged image */
34 short fvxr, fvyr; /* foveal image resolution */
35 float bwhist[HISTRES]; /* luminance histogram */
36 double histot; /* total count of histogram */
37 double bwmin, bwmax; /* histogram limits */
38 double bwavg; /* mean brightness */
39
40 double scalef = 0.; /* linear scaling factor */
41
42
43 main(argc, argv)
44 int argc;
45 char *argv[];
46 {
47 static RGBPRIMS outprimS;
48 int i;
49 #define bool(flg) switch (argv[i][2]) { \
50 case '\0': what2do ^= flg; break; \
51 case 'y': case 'Y': case 't': case 'T': \
52 case '+': case '1': what2do |= flg; break; \
53 case 'n': case 'N': case 'f': case 'F': \
54 case '-': case '0': what2do &= ~(flg); break; \
55 default: goto userr; }
56
57 progname = argv[0];
58
59 for (i = 1; i < argc && argv[i][0] == '-'; i++)
60 switch (argv[i][1]) {
61 case 'h':
62 bool(DO_HUMAN);
63 break;
64 case 'a':
65 bool(DO_ACUITY);
66 break;
67 case 'v':
68 bool(DO_VEIL);
69 break;
70 case 's':
71 bool(DO_HSENS);
72 break;
73 case 'c':
74 bool(DO_COLOR);
75 break;
76 case 'w':
77 bool(DO_CWEIGHT);
78 break;
79 case 'l':
80 bool(DO_LINEAR);
81 break;
82 case 'p':
83 if (i+8 >= argc) goto userr;
84 outprimS[RED][CIEX] = atof(argv[++i]);
85 outprimS[RED][CIEY] = atof(argv[++i]);
86 outprimS[GRN][CIEX] = atof(argv[++i]);
87 outprimS[GRN][CIEY] = atof(argv[++i]);
88 outprimS[BLU][CIEX] = atof(argv[++i]);
89 outprimS[BLU][CIEY] = atof(argv[++i]);
90 outprimS[WHT][CIEX] = atof(argv[++i]);
91 outprimS[WHT][CIEY] = atof(argv[++i]);
92 outprims = outprimS;
93 break;
94 case 'e':
95 if (i+1 >= argc) goto userr;
96 scalef = atof(argv[++i]);
97 if (argv[i][0] == '+' | argv[i][0] == '-')
98 scalef = pow(2.0, scalef);
99 what2do |= DO_LINEAR;
100 break;
101 case 'f':
102 if (i+1 >= argc) goto userr;
103 mbcalfile = argv[++i];
104 break;
105 case 't':
106 if (i+1 >= argc) goto userr;
107 ldmax = atof(argv[++i]);
108 if (ldmax <= FTINY)
109 goto userr;
110 break;
111 case 'b':
112 if (i+1 >= argc) goto userr;
113 ldmin = atof(argv[++i]);
114 break;
115 case 'm':
116 if (i+1 >= argc) goto userr;
117 if ((mapfp = fopen(argv[++i], "w")) == NULL) {
118 fprintf(stderr,
119 "%s: cannot open for writing\n",
120 argv[i]);
121 exit(1);
122 }
123 break;
124 default:
125 goto userr;
126 }
127 if (mbcalfile != NULL & outprims != stdprims) {
128 fprintf(stderr, "%s: only one of -p or -f option supported\n",
129 progname);
130 exit(1);
131 }
132 if (outprims == stdprims & inprims != stdprims)
133 outprims = inprims;
134 if (ldmin <= FTINY)
135 ldmin = ldmax*LDMINF;
136 else if (ldmin >= ldmax) {
137 fprintf(stderr, "%s: Ldmin (%f) >= Ldmax (%f)!\n", progname,
138 ldmin, ldmax);
139 exit(1);
140 }
141 Bldmin = Bl(ldmin);
142 Bldmax = Bl(ldmax);
143 if (i >= argc || i+2 < argc)
144 goto userr;
145 if ((infp = fopen(infn=argv[i], "r")) == NULL)
146 syserror(infn);
147 if (i+2 == argc && freopen(argv[i+1], "w", stdout) == NULL)
148 syserror(argv[i+1]);
149 #ifdef MSDOS
150 setmode(fileno(infp), O_BINARY);
151 setmode(fileno(stdout), O_BINARY);
152 #endif
153 getahead(); /* load input header */
154 printargs(argc, argv, stdout); /* add to output header */
155 if (outprims != inprims)
156 fputprims(outprims, stdout);
157 mapimage(); /* map the picture */
158 if (mapfp != NULL) /* write out basic mapping */
159 putmapping(mapfp);
160 exit(0);
161 userr:
162 fprintf(stderr, "Usage: %s [-{h|a|v|s|c|l|w}[+-]][-e ev][-p xr yr xg yg xb yb xw yw|-f mbf.cal][-t Ldmax][-b Ldmin][-m mapfile] inpic [outpic]\n",
163 progname);
164 exit(1);
165 #undef bool
166 }
167
168
169 syserror(s) /* report system error and exit */
170 char *s;
171 {
172 fprintf(stderr, "%s: ", progname);
173 perror(s);
174 exit(2);
175 }
176
177
178 headline(s) /* process header line */
179 char *s;
180 {
181 static RGBPRIMS inprimS;
182 char fmt[32];
183
184 if (formatval(fmt, s)) { /* check if format string */
185 if (!strcmp(fmt,COLRFMT)) lumf = rgblum;
186 else if (!strcmp(fmt,CIEFMT)) lumf = cielum;
187 else lumf = NULL;
188 return; /* don't echo */
189 }
190 if (isprims(s)) { /* get input primaries */
191 primsval(inprimS, s);
192 inprims= inprimS;
193 return; /* don't echo */
194 }
195 if (isexpos(s)) { /* picture exposure */
196 inpexp *= exposval(s);
197 return; /* don't echo */
198 }
199 if (isaspect(s)) /* pixel aspect ratio */
200 pixaspect *= aspectval(s);
201 if (isview(s)) /* image view */
202 gotview += sscanview(&ourview, s);
203 fputs(s, stdout);
204 }
205
206
207 getahead() /* load picture header */
208 {
209 char *err;
210
211 getheader(infp, headline, NULL);
212 if (lumf == NULL || !fgetsresolu(&inpres, infp)) {
213 fprintf(stderr, "%s: %s: not a Radiance picture\n",
214 progname, infn);
215 exit(1);
216 }
217 if (lumf == rgblum)
218 comprgb2xyzmat(inrgb2xyz, inprims);
219 else if (mbcalfile != NULL) {
220 fprintf(stderr, "%s: macbethcal only works with RGB pictures\n",
221 progname);
222 exit(1);
223 }
224 if (!gotview || ourview.type == VT_PAR) {
225 copystruct(&ourview, &stdview);
226 ourview.type = VT_PER;
227 if (pixaspect*inpres.yr < inpres.xr) {
228 ourview.horiz = 40.0;
229 ourview.vert = 2.*180./PI *
230 atan(.3639702*pixaspect*inpres.yr/inpres.xr);
231 } else {
232 ourview.vert = 40.0;
233 ourview.horiz = 2.*180./PI *
234 atan(.3639702*inpres.xr/pixaspect/inpres.yr);
235 }
236 }
237 if ((err = setview(&ourview)) != NULL) {
238 fprintf(stderr, "%s: view error in picture \"%s\": %s\n",
239 progname, infn, err);
240 exit(1);
241 }
242 }
243
244
245 mapimage() /* map picture and send to stdout */
246 {
247 COLOR *scan;
248
249 #ifdef DEBUG
250 fprintf(stderr, "%s: generating histogram...", progname);
251 #endif
252 getfovimg(); /* get foveal sample image */
253 comphist(); /* generate adaptation histogram */
254 #ifdef DEBUG
255 fputs("done\n", stderr);
256 #endif
257 check2do(); /* modify what2do flags */
258 if (what2do&DO_VEIL) {
259 #ifdef DEBUG
260 fprintf(stderr, "%s: computing veiling...", progname);
261 #endif
262 compveil();
263 #ifdef DEBUG
264 fputs("done\n", stderr);
265 #endif
266 }
267 #ifdef DEBUG
268 fprintf(stderr, "%s: computing brightness mapping...", progname);
269 #endif
270 if (!(what2do&DO_LINEAR) && mkbrmap() < 0) { /* make tone map */
271 what2do |= DO_LINEAR; /* use linear scaling */
272 #ifdef DEBUG
273 fputs("failed!\n", stderr);
274 } else
275 fputs("done\n", stderr);
276 #else
277 }
278 #endif
279 if (what2do&DO_LINEAR) {
280 if (scalef <= FTINY) {
281 if (what2do&DO_HSENS)
282 scalef = htcontrs(Lb(0.5*(Bldmax+Bldmin))) /
283 htcontrs(Lb(bwavg));
284 else
285 scalef = Lb(0.5*(Bldmax+Bldmin)) / Lb(bwavg);
286 scalef *= WHTEFFICACY/(inpexp*ldmax);
287 }
288 #ifdef DEBUG
289 fprintf(stderr, "%s: linear scaling factor = %f\n",
290 progname, scalef);
291 #endif
292 fputexpos(inpexp*scalef, stdout); /* record exposure */
293 if (lumf == cielum) scalef /= WHTEFFICACY;
294 }
295 putchar('\n'); /* complete header */
296 fputsresolu(&inpres, stdout); /* resolution doesn't change */
297
298 for (scan = firstscan(); scan != NULL; scan = nextscan())
299 if (fwritescan(scan, scanlen(&inpres), stdout) < 0) {
300 fprintf(stderr, "%s: scanline write error\n",
301 progname);
302 exit(1);
303 }
304 }
305
306
307 double
308 centprob(x, y) /* center-weighting probability function */
309 int x, y;
310 {
311 double xr, yr, p;
312 /* paraboloid, 0 at 90 degrees from center */
313 xr = (x - .5*(fvxr-1))/90.; /* 180 degree fisheye has fv?r == 90 */
314 yr = (y - .5*(fvyr-1))/90.;
315 p = 1. - xr*xr - yr*yr;
316 return(p < 0. ? 0. : p);
317 }
318
319
320 getfovimg() /* load foveal sampled image */
321 {
322 extern FILE *popen();
323 char combuf[128];
324 FILE *fp;
325 int x, y;
326 /* compute image size */
327 fvxr = sqrt(ourview.hn2)/FOVDIA + 0.5;
328 if (fvxr < 2) fvxr = 2;
329 fvyr = sqrt(ourview.vn2)/FOVDIA + 0.5;
330 if (fvyr < 2) fvyr = 2;
331 if (!(inpres.or & YMAJOR)) { /* picture is rotated? */
332 y = fvyr;
333 fvyr = fvxr;
334 fvxr = y;
335 }
336 if ((fovimg = (COLOR *)malloc(fvxr*fvyr*sizeof(COLOR))) == NULL)
337 syserror("malloc");
338 sprintf(combuf, "pfilt -1 -b -x %d -y %d %s", fvxr, fvyr, infn);
339 if ((fp = popen(combuf, "r")) == NULL)
340 syserror("popen");
341 getheader(fp, NULL, NULL); /* skip header */
342 if (fgetresolu(&x, &y, fp) < 0 || x != fvxr | y != fvyr)
343 goto readerr;
344 for (y = 0; y < fvyr; y++)
345 if (freadscan(fovscan(y), fvxr, fp) < 0)
346 goto readerr;
347 pclose(fp);
348 return;
349 readerr:
350 fprintf(stderr, "%s: error reading from pfilt process in fovimage\n",
351 progname);
352 exit(1);
353 }
354
355
356 comphist() /* create foveal sampled image and histogram */
357 {
358 double l, b, lwmin, lwmax;
359 register int x, y;
360
361 lwmin = 1e10; /* find extrema */
362 lwmax = 0.;
363 for (y = 0; y < fvyr; y++)
364 for (x = 0; x < fvxr; x++) {
365 l = plum(fovscan(y)[x]);
366 if (l < lwmin) lwmin = l;
367 if (l > lwmax) lwmax = l;
368 }
369 lwmin -= FTINY;
370 lwmax += FTINY;
371 if (lwmin < LMIN) lwmin = LMIN;
372 if (lwmax > LMAX) lwmax = LMAX;
373 /* compute histogram */
374 bwmin = Bl(lwmin);
375 bwmax = Bl(lwmax);
376 bwavg = 0.;
377 for (y = 0; y < fvyr; y++)
378 for (x = 0; x < fvxr; x++) {
379 l = plum(fovscan(y)[x]);
380 if (l < lwmin) continue;
381 if (l > lwmax) continue;
382 b = Bl(l);
383 bwavg += b;
384 l = what2do&DO_CWEIGHT ? centprob(x,y) : 1.;
385 bwhist[bwhi(b)] += l;
386 histot += l;
387 }
388 bwavg /= histot;
389 }
390
391
392 check2do() /* check histogram to see what isn't worth doing */
393 {
394 double sum;
395 double b, l;
396 register int i;
397
398 /* check for within display range */
399 l = Lb(bwmax)/Lb(bwmin);
400 if (l <= ldmax/ldmin)
401 what2do |= DO_LINEAR;
402 /* determine if veiling significant */
403 if (l < 100.) /* heuristic */
404 what2do &= ~DO_VEIL;
405
406 if (!(what2do & (DO_ACUITY|DO_COLOR)))
407 return;
408 /* find 5th percentile */
409 sum = histot*0.05;
410 for (i = 0; i < HISTRES; i++)
411 if ((sum -= bwhist[i]) <= 0)
412 break;
413 b = (i+.5)*(bwmax-bwmin)/HISTRES + bwmin;
414 l = Lb(b);
415 /* determine if acuity adj. useful */
416 if (what2do&DO_ACUITY &&
417 hacuity(l) >= (inpres.xr/sqrt(ourview.hn2) +
418 inpres.yr/sqrt(ourview.vn2))/(2.*180./PI))
419 what2do &= ~DO_ACUITY;
420 /* color sensitivity loss? */
421 if (l >= TopMesopic)
422 what2do &= ~DO_COLOR;
423 }