ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/pcond.c
Revision: 3.7
Committed: Wed Jan 8 21:40:38 1997 UTC (27 years, 3 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 3.6: +5 -5 lines
Log Message:
fixed bad calls to log() in fovhist()

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