ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/px/pmblur2.c
Revision: 2.5
Committed: Sat Oct 13 05:18:18 2012 UTC (11 years, 6 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad4R2P2, rad5R0, rad5R1, rad4R2, rad4R2P1
Changes since 2.4: +3 -4 lines
Log Message:
Reduced sampling artifacts substantially

File Contents

# User Rev Content
1 greg 2.1 #ifndef lint
2 greg 2.5 static const char RCSid[] = "$Id: pmblur2.c,v 2.4 2012/10/05 18:54:40 greg Exp $";
3 greg 2.1 #endif
4     /*
5     * pmblur2.c - program to computer better motion blur from ranimove frames.
6     *
7     * Created by Greg Ward on 10/3/12.
8     */
9    
10     #include "copyright.h"
11    
12     #include "standard.h"
13     #include "platform.h"
14     #include "color.h"
15     #include "view.h"
16    
17     #define ZGONE (2.*FHUGE) /* non-existent depth value */
18    
19     const char *hdrspec; /* HDR image specification string */
20     const char *zbfspec; /* depth buffer specification string */
21     const char *mvospec; /* motion vector specification string */
22    
23     RESOLU imres; /* image size and orientation */
24     double pixaspect = 1.; /* pixel aspect ratio */
25    
26     COLOR *imsum; /* blurred image in progress */
27     VIEW vwsum; /* averaged view for blurred image */
28     char imfmt[32] = PICFMT; /* image format */
29    
30     int nsamps; /* number of time samples */
31    
32     COLOR *imbuf; /* accumulation buffer for this time sample */
33     float *zbuf; /* associated depth buffer */
34    
35     COLR *imprev; /* previous loaded image frame */
36     float *zprev; /* loaded depth buffer */
37     double exprev; /* previous image exposure */
38     VIEW vwprev = STDVIEW; /* previous view */
39     int fnprev; /* frame number for previous view */
40    
41     char *progname; /* global argv[0] */
42    
43     typedef struct {
44     VIEW *vp; /* view pointer */
45     int gotview; /* got view parameters? */
46     double ev; /* exposure value */
47     } IMGHEAD; /* holder for individual header information */
48    
49    
50     /* Process line from picture header */
51     static int
52     headline(char *s, void *p)
53     {
54     IMGHEAD *ip = (IMGHEAD *)p;
55     char fmt[32];
56    
57     if (isview(s)) {
58     ip->gotview += (sscanview(ip->vp, s) > 0);
59     return(1);
60     }
61     if (isexpos(s)) {
62     ip->ev *= exposval(s);
63     return(1);
64     }
65     if (formatval(fmt, s)) {
66     if (globmatch(imfmt, fmt)) {
67     strcpy(imfmt, fmt);
68     return(1);
69     }
70     error(WARNING, "wrong picture format");
71     return(-1);
72     }
73     if (isaspect(s)) {
74     pixaspect = aspectval(s);
75     return(1);
76     }
77     return(0);
78     }
79    
80    
81     /* Load previous frame (allocates arrays on first call) */
82     static void
83     loadprev(int fno)
84     {
85     char *err;
86     IMGHEAD ih;
87     RESOLU rs;
88     char fname[256];
89     FILE *fp;
90     int y, hres;
91     int fd;
92    
93     if (fno == fnprev) /* already loaded? */
94     return;
95     sprintf(fname, hdrspec, fno); /* open picture file */
96     if ((fp = fopen(fname, "r")) == NULL) {
97     sprintf(errmsg, "cannot open picture file \"%s\"", fname);
98     error(SYSTEM, errmsg);
99     }
100     ih.vp = &vwprev;
101     ih.gotview = 0;
102     ih.ev = 1.;
103     if (getheader(fp, headline, &ih) < 0)
104     goto readerr;
105     if (!ih.gotview) {
106     sprintf(errmsg, "missing view in picture \"%s\"", fname);
107     error(USER, errmsg);
108     }
109     if ((err = setview(&vwprev)) != NULL)
110     error(USER, err);
111     exprev = ih.ev;
112     if (!fgetsresolu(&rs, fp))
113     goto readerr;
114 greg 2.3 if (rs.rt != PIXSTANDARD) {
115     sprintf(errmsg, "unsupported orientation in picture \"%s\"", fname);
116     error(USER, errmsg);
117     }
118 greg 2.1 if (!imres.xr) { /* allocate buffers */
119     imres = rs;
120     imsum = (COLOR *)ecalloc(imres.xr*imres.yr, sizeof(COLOR));
121     imbuf = (COLOR *)emalloc(sizeof(COLOR)*imres.xr*imres.yr);
122     zbuf = (float *)emalloc(sizeof(float)*imres.xr*imres.yr);
123     imprev = (COLR *)emalloc(sizeof(COLR)*rs.xr*rs.yr);
124     zprev = (float *)emalloc(sizeof(float)*rs.xr*rs.yr);
125     } else if ((imres.xr != rs.xr) | (imres.yr != rs.yr)) {
126     sprintf(errmsg, "resolution mismatch in picture \"%s\"", fname);
127     error(USER, errmsg);
128     }
129     hres = scanlen(&rs); /* load picture */
130     for (y = 0; y < numscans(&rs); y++)
131     if (freadcolrs(imprev + y*hres, hres, fp) < 0)
132     goto readerr;
133     fclose(fp);
134     sprintf(fname, zbfspec, fno); /* load depth buffer */
135     if ((fd = open(fname, O_RDONLY)) < 0) {
136     sprintf(errmsg, "cannot open depth buffer \"%s\"", fname);
137     error(SYSTEM, errmsg);
138     }
139     if (read(fd, zprev, sizeof(float)*rs.xr*rs.yr) != sizeof(float)*rs.xr*rs.yr)
140     goto readerr;
141     close(fd);
142     fnprev = fno;
143     return;
144     readerr:
145     sprintf(errmsg, "error loading file \"%s\"", fname);
146     error(USER, errmsg);
147     }
148    
149    
150     /* Interpolate between two image pixels */
151     /* XXX skipping expensive ray vector calculations for now */
152     static void
153 greg 2.3 interp_pixel(COLOR con, const VIEW *vwn, int xn, int yn, double zval,
154     double pos, int xp, int yp, int interpOK)
155 greg 2.1 {
156     const int hres = scanlen(&imres);
157     RREAL ipos;
158     FVECT wprev, wcoor, rdir;
159     int np, xd, yd, nd;
160     COLOR cpr;
161 greg 2.3 double sf;
162 greg 2.1 /* check if off image */
163     if ((xp < 0) | (xp >= hres))
164     return;
165     if ((yp < 0) | (yp >= numscans(&imres)))
166     return;
167     np = yp*hres + xp; /* get indexed destination */
168     xd = (1.-pos)*xp + pos*xn + .5;
169     yd = (1.-pos)*yp + pos*yn + .5;
170     nd = yd*hres + xd;
171 greg 2.3 /* check depth */
172     if (interpOK)
173     zval = (1.-pos)*zprev[np] + pos*zval;
174 greg 2.1 if (zval >= zbuf[nd])
175     return;
176 greg 2.3 zbuf[nd] = zval; /* assign new depth */
177     copycolor(imbuf[nd], con); /* assign new color */
178     if (!interpOK)
179     return;
180 greg 2.1 scalecolor(imbuf[nd], pos);
181     sf = (1.-pos)/exprev;
182     colr_color(cpr, imprev[np]);
183     scalecolor(cpr, sf);
184     addcolor(imbuf[nd], cpr);
185     }
186    
187    
188     /* Find neighbor with minimum depth from buffer */
189     static int
190     neigh_zmin(const float *zb, int n)
191     {
192     const int hres = scanlen(&imres);
193     const int vc = n / hres;
194     const int hc = n - vc*hres;
195     float zbest = ZGONE;
196     int nbest = 0;
197    
198     if (hc > 0 && zb[n-1] < zbest)
199     zbest = zb[nbest = n-1];
200     if (hc < hres-1 && zb[n+1] < zbest)
201     zbest = zb[nbest = n+1];
202     if (vc > 0 && zb[n-hres] < zbest)
203     zbest = zb[nbest = n-hres];
204     if (vc < numscans(&imres)-1 && zb[n+hres] < zbest)
205     return(n+hres);
206     return(nbest);
207     }
208    
209    
210 greg 2.5 /* Expand foreground pixels to mitigate aliasing/fill errors */
211 greg 2.1 static void
212     fill_missing(void)
213     {
214     int n, m;
215    
216     for (n = imres.xr*imres.yr; n--; )
217 greg 2.5 if (zbuf[n] > 1.25*zbuf[m = neigh_zmin(zbuf,n)])
218 greg 2.1 copycolor(imbuf[n], imbuf[m]);
219     }
220    
221    
222     /* Load a subframe */
223     static void
224     addframe(double fpos)
225     {
226     COLOR backg;
227     char *err;
228     VIEW fvw;
229     RESOLU rs;
230     IMGHEAD ihead;
231     char fname[256];
232     FILE *fpimg;
233     int fdzbf, fdmvo;
234     COLR *iscan;
235     float *zscan;
236     uint16 *mscan;
237     int hres, n, h;
238    
239     loadprev((int)fpos); /* make sure we have previous */
240     fpos -= floor(fpos);
241     if (fpos <= .02) { /* special case */
242     COLOR col;
243     for (n = imres.xr*imres.yr; n--; ) {
244     colr_color(col, imprev[n]);
245     addcolor(imsum[n], col);
246     }
247     return;
248     }
249     /* open input files */
250     sprintf(fname, hdrspec, fnprev+1);
251     if ((fpimg = fopen(fname, "r")) == NULL)
252     goto openerr;
253     fvw = vwprev;
254     ihead.vp = &fvw;
255     ihead.gotview = 0;
256     ihead.ev = 1.;
257     if (getheader(fpimg, headline, &ihead) < 0)
258     goto readerr;
259     if (!ihead.gotview) {
260     sprintf(errmsg, "missing view in picture \"%s\"", fname);
261     error(USER, errmsg);
262     }
263     if ((err = setview(&fvw)) != NULL)
264     error(USER, err);
265     if (!fgetsresolu(&rs, fpimg))
266     goto readerr;
267 greg 2.3 if (rs.rt != PIXSTANDARD) {
268     sprintf(errmsg, "unsupported orientation in picture \"%s\"", fname);
269     error(USER, errmsg);
270     }
271 greg 2.1 if ((rs.xr != imres.xr) | (rs.yr != imres.yr)) {
272     sprintf(errmsg, "resolution mismatch for picture \"%s\"", fname);
273     error(USER, errmsg);
274     }
275     vwsum.type = fvw.type; /* average in view */
276     for (n = 3; n--; ) {
277     vwsum.vp[n] += (1.-fpos)*vwprev.vp[n] + fpos*fvw.vp[n];
278     vwsum.vdir[n] += (1.-fpos)*vwprev.vdir[n] + fpos*fvw.vdir[n];
279     vwsum.vup[n] += (1.-fpos)*vwprev.vdir[n] + fpos*fvw.vdir[n];
280     }
281     vwsum.vdist += (1.-fpos)*vwprev.vdist + fpos*fvw.vdist;
282     vwsum.horiz += (1.-fpos)*vwprev.horiz + fpos*fvw.horiz;
283     vwsum.vert += (1.-fpos)*vwprev.vert + fpos*fvw.vert;
284     vwsum.voff += (1.-fpos)*vwprev.voff + fpos*fvw.voff;
285     vwsum.hoff += (1.-fpos)*vwprev.hoff + fpos*fvw.hoff;
286     vwsum.vfore += (1.-fpos)*vwprev.vfore + fpos*fvw.vfore;
287     vwsum.vaft += (1.-fpos)*vwprev.vaft + fpos*fvw.vaft;
288     sprintf(fname, zbfspec, fnprev+1);
289     if ((fdzbf = open(fname, O_RDONLY)) < 0)
290     goto openerr;
291     sprintf(fname, mvospec, fnprev+1);
292     if ((fdmvo = open(fname, O_RDONLY)) < 0)
293     goto openerr;
294     setcolor(backg, .0, .0, .0); /* clear image & depth buffers */
295     for (n = rs.xr*rs.yr; n--; ) {
296     if (zprev[n] >= .9*FHUGE &&
297     zprev[neigh_zmin(zprev,n)] >= .9*FHUGE)
298     colr_color(backg, imprev[n]);
299     copycolor(imbuf[n], backg);
300     zbuf[n] = ZGONE;
301     }
302     hres = scanlen(&rs); /* process a scanline at a time */
303     iscan = (COLR *)emalloc(sizeof(COLR)*hres);
304     zscan = (float *)emalloc(sizeof(float)*hres);
305     mscan = (uint16 *)emalloc(sizeof(uint16)*3*hres);
306     for (n = 0; n < numscans(&rs); n++) {
307     const double sf = 1./ihead.ev;
308     COLOR cval;
309     if (freadcolrs(iscan, hres, fpimg) < 0)
310     goto readerr;
311     if (read(fdzbf, zscan, sizeof(float)*hres) !=
312     sizeof(float)*hres)
313     goto readerr;
314     if (read(fdmvo, mscan, sizeof(uint16)*3*hres) !=
315     sizeof(uint16)*3*hres)
316     goto readerr;
317     for (h = hres; h--; ) {
318     if (!mscan[3*h])
319     continue; /* unknown motion */
320     colr_color(cval, iscan[h]);
321     scalecolor(cval, sf);
322     interp_pixel(cval, &fvw, h, n, zscan[h], fpos,
323     h + (int)mscan[3*h] - 32768,
324 greg 2.3 n - (int)mscan[3*h+1] + 32768,
325     mscan[3*h+2]);
326 greg 2.1 }
327     }
328     /* fill in missing pixels */
329     fill_missing();
330     /* add in results */
331     for (n = imres.xr*imres.yr; n--; )
332     addcolor(imsum[n], imbuf[n]);
333     /* clean up */
334     free(mscan);
335     free(zscan);
336     free(iscan);
337     close(fdmvo);
338     close(fdzbf);
339     fclose(fpimg);
340     return;
341     openerr:
342     sprintf(errmsg, "cannot open input \"%s\"", fname);
343     error(SYSTEM, errmsg);
344     readerr:
345     error(SYSTEM, "error reading input file");
346     }
347    
348    
349     /* Average and write our image out to the given file */
350     static void
351     write_average(FILE *fp)
352     {
353     const int hres = scanlen(&imres);
354     const double sf = 1./(double)nsamps;
355     const float csf = exprev*sf;
356     int n, h;
357     /* normalize view & image */
358     for (n = 3; n--; ) {
359     vwsum.vp[n] *= sf;
360     vwsum.vdir[n] *= sf;
361     vwsum.vup[n] *= sf;
362     }
363     vwsum.vdist *= sf;
364     vwsum.horiz *= sf;
365     vwsum.vert *= sf;
366     vwsum.voff *= sf;
367     vwsum.hoff *= sf;
368     vwsum.vfore *= sf;
369     vwsum.vaft *= sf;
370     for (n = imres.xr*imres.yr; n--; )
371     scalecolor(imsum[n], csf);
372     /* write it out */
373     fputs(VIEWSTR, fp);
374     fprintview(&vwsum, fp);
375     fputc('\n', fp);
376     if ((pixaspect < .98) | (pixaspect > 1.02))
377     fputaspect(pixaspect, fp);
378     fputexpos(exprev, fp);
379     if (strcmp(imfmt, PICFMT))
380     fputformat(imfmt, fp);
381     fputc('\n', fp); /* end of info. header */
382     fputsresolu(&imres, fp);
383     for (n = 0; n < numscans(&imres); n++)
384     if (fwritescan(imsum + n*hres, hres, fp) < 0)
385     error(SYSTEM, "write error on picture output");
386     }
387    
388     /* Load frame sequence and compute motion blur output */
389     int
390     main(int argc, char *argv[])
391     {
392     double fstart, fend, fstep, fcur;
393    
394     progname = argv[0];
395     SET_DEFAULT_BINARY();
396     if (argc != 5)
397     goto userr;
398 greg 2.2 /* get frame range & sampling */
399     switch (sscanf(argv[1], "%lf,%lf/%d", &fstart, &fend, &nsamps)) {
400 greg 2.1 case 1:
401 greg 2.4 fend = fstart;
402 greg 2.2 nsamps = 0;
403 greg 2.1 break;
404     case 2:
405 greg 2.2 nsamps = 0;
406     /* fall through */
407     case 3:
408 greg 2.1 if (fend < fstart)
409     goto userr;
410 greg 2.4 if (fend <= fstart+FTINY)
411     nsamps = 0;
412 greg 2.1 break;
413     default:
414     goto userr;
415     }
416     if (fstart < 1)
417     goto userr;
418 greg 2.2 if (nsamps <= 0)
419     nsamps = (fend - fstart)*17.;
420 greg 2.1 if (nsamps) {
421 greg 2.2 if (nsamps > 100) nsamps = 100;
422 greg 2.1 fstep = (fend - fstart)/nsamps;
423     } else {
424     fstep = 1.;
425     fstart -= .5;
426     nsamps = 1;
427     }
428 greg 2.2 hdrspec = argv[2];
429     zbfspec = argv[3];
430     mvospec = argv[4];
431 greg 2.1 /* load and filter each subframe */
432     for (fcur = fstart + .5*fstep; fcur <= fend; fcur += fstep)
433     addframe(fcur);
434     /* write results to stdout */
435     SET_FILE_BINARY(stdout);
436     newheader("RADIANCE", stdout);
437     printargs(argc, argv, stdout);
438     fputnow(stdout);
439     write_average(stdout);
440     return(fflush(stdout) == EOF);
441     userr:
442 greg 2.4 fprintf(stderr, "Usage: %s f0,f1[/n] HDRspec ZBUFspec MVOspec\n", progname);
443 greg 2.1 return(1);
444     }