--- ray/src/px/pinterp.c 1989/12/12 11:21:25 1.2 +++ ray/src/px/pinterp.c 1990/01/04 17:47:37 1.14 @@ -17,9 +17,18 @@ static char SCCSid[] = "$SunId$ LBL"; #define pscan(y) (ourpict+(y)*ourview.hresolu) #define zscan(y) (ourzbuf+(y)*ourview.hresolu) +#define F_FORE 1 /* fill foreground */ +#define F_BACK 2 /* fill background */ + +#define PACKSIZ 42 /* calculation packet size */ + +#define RTCOM "rtrace -h -ovl -fff -x %d %s" + +#define ABS(x) ((x)>0?(x):-(x)) + VIEW ourview = STDVIEW(512); /* desired view */ -double zeps = 0.001; /* allowed z epsilon */ +double zeps = .02; /* allowed z epsilon */ COLR *ourpict; /* output picture */ float *ourzbuf; /* corresponding z-buffer */ @@ -29,13 +38,29 @@ char *progname; VIEW theirview = STDVIEW(512); /* input view */ int gotview; /* got input view? */ +int fill = F_FORE|F_BACK; /* selected fill algorithm */ +extern int backfill(), rcalfill(); /* fill functions */ +int (*deffill)() = backfill; /* selected fill function */ +COLR backcolr = BLKCOLR; /* background color */ +double backz = 0.0; /* background z value */ +double theirs2ours[4][4]; /* transformation matrix */ +int normdist = 1; /* normalized distance? */ + +int childpid = -1; /* id of fill process */ +FILE *psend, *precv; /* pipes to/from fill calculation */ +int queue[PACKSIZ][2]; /* pending pixels */ +int queuesiz; /* number of pixels pending */ + + main(argc, argv) /* interpolate pictures */ int argc; char *argv[]; { #define check(olen,narg) if (argv[i][olen] || narg >= argc-i) goto badopt + extern double atof(); int gotvfile = 0; + char *zfile = NULL; char *err; int i; @@ -47,6 +72,61 @@ char *argv[]; check(2,1); zeps = atof(argv[++i]); break; + case 'n': /* dist. normalized? */ + check(2,0); + normdist = !normdist; + break; + case 'f': /* fill type */ + switch (argv[i][2]) { + case '0': /* none */ + check(3,0); + fill = 0; + break; + case 'f': /* foreground */ + check(3,0); + fill = F_FORE; + break; + case 'b': /* background */ + check(3,0); + fill = F_BACK; + break; + case 'a': /* all */ + check(3,0); + fill = F_FORE|F_BACK; + break; + case 'c': /* color */ + check(3,3); + deffill = backfill; + setcolr(backcolr, atof(argv[i+1]), + atof(argv[i+2]), atof(argv[i+3])); + i += 3; + break; + case 'z': /* z value */ + check(3,1); + deffill = backfill; + backz = atof(argv[++i]); + break; + case 'r': /* rtrace */ + check(3,1); + deffill = rcalfill; + calstart(RTCOM, argv[++i]); + break; + default: + goto badopt; + } + break; + case 'z': /* z file */ + check(2,1); + zfile = argv[++i]; + break; + case 'x': /* x resolution */ + check(2,1); + ourview.hresolu = atoi(argv[++i]); + break; + case 'y': /* y resolution */ + check(2,1); + ourview.vresolu = atoi(argv[++i]); + break; case 'v': /* view */ switch (argv[i][2]) { case 't': /* type */ @@ -97,23 +177,20 @@ char *argv[]; break; default: badopt: - fprintf(stderr, "%s: unknown option '%s'\n", + fprintf(stderr, "%s: command line error at '%s'\n", progname, argv[i]); - exit(1); + goto userr; } /* check arguments */ - if (argc-i < 2 || (argc-i)%2) { - fprintf(stderr, "Usage: %s [view args] pfile zfile ..\n", - progname); - exit(1); - } + if ((argc-i)%2) + goto userr; /* set view */ if (err = setview(&ourview)) { fprintf(stderr, "%s: %s\n", progname, err); exit(1); } /* allocate frame */ - ourpict = (COLR *)calloc(ourview.hresolu*ourview.vresolu,sizeof(COLR)); + ourpict = (COLR *)malloc(ourview.hresolu*ourview.vresolu*sizeof(COLR)); ourzbuf = (float *)calloc(ourview.hresolu*ourview.vresolu,sizeof(float)); if (ourpict == NULL || ourzbuf == NULL) { perror(progname); @@ -123,7 +200,12 @@ char *argv[]; for ( ; i < argc; i += 2) addpicture(argv[i], argv[i+1]); /* fill in spaces */ - fillpicture(); + if (fill&F_BACK) + backpicture(); + else + fillpicture(); + /* close calculation */ + caldone(); /* add to header */ printargs(argc, argv, stdout); if (gotvfile) { @@ -132,10 +214,18 @@ char *argv[]; printf("\n"); } printf("\n"); - /* write output */ + /* write picture */ writepicture(); + /* write z file */ + if (zfile != NULL) + writedistance(zfile); exit(0); +userr: + fprintf(stderr, + "Usage: %s [view opts][-t zthresh][-z zout][-fT][-n] pfile zspec ..\n", + progname); + exit(1); #undef check } @@ -157,99 +247,202 @@ char *s; } -addpicture(pfile, zfile) /* add picture to output */ -char *pfile, *zfile; +addpicture(pfile, zspec) /* add picture to output */ +char *pfile, *zspec; { + extern double atof(); FILE *pfp, *zfp; - COLR *scanin; - float *zin; char *err; - int xres, yres; + COLR *scanin; + float *zin, *zlast; + int *plast; int y; - /* open input files */ + /* open picture file */ if ((pfp = fopen(pfile, "r")) == NULL) { perror(pfile); exit(1); } - if ((zfp = fopen(zfile, "r")) == NULL) { - perror(zfile); - exit(1); - } /* get header and view */ printf("%s:\n", pfile); gotview = 0; getheader(pfp, headline); - if (!gotview || fgetresolu(&xres, &yres, pfp) != (YMAJOR|YDECR)) { + if (!gotview || fgetresolu(&theirview.hresolu, &theirview.vresolu, pfp) + != (YMAJOR|YDECR)) { fprintf(stderr, "%s: picture view error\n", pfile); exit(1); } - theirview.hresolu = xres; - theirview.vresolu = yres; if (err = setview(&theirview)) { fprintf(stderr, "%s: %s\n", pfile, err); exit(1); } + /* compute transformation */ + pixform(theirs2ours, &theirview, &ourview); /* allocate scanlines */ - scanin = (COLR *)malloc(xres*sizeof(COLR)); - zin = (float *)malloc(xres*sizeof(float)); - if (scanin == NULL || zin == NULL) { + scanin = (COLR *)malloc(theirview.hresolu*sizeof(COLR)); + zin = (float *)malloc(theirview.hresolu*sizeof(float)); + plast = (int *)calloc(theirview.hresolu, sizeof(int)); + zlast = (float *)calloc(theirview.hresolu, sizeof(float)); + if (scanin == NULL || zin == NULL || plast == NULL || zlast == NULL) { perror(progname); exit(1); } + /* get z specification or file */ + if ((zfp = fopen(zspec, "r")) == NULL) { + double zvalue; + register int x; + if (!isfloat(zspec) || (zvalue = atof(zspec)) <= 0.0) { + perror(zspec); + exit(1); + } + for (x = 0; x < theirview.hresolu; x++) + zin[x] = zvalue; + } /* load image */ - for (y = yres-1; y >= 0; y--) { - if (freadcolrs(scanin, xres, pfp) < 0) { + for (y = theirview.vresolu-1; y >= 0; y--) { + if (freadcolrs(scanin, theirview.hresolu, pfp) < 0) { fprintf(stderr, "%s: read error\n", pfile); exit(1); } - if (fread(zin, sizeof(float), xres, zfp) < xres) { - fprintf(stderr, "%s: read error\n", zfile); + if (zfp != NULL + && fread(zin,sizeof(float),theirview.hresolu,zfp) + < theirview.hresolu) { + fprintf(stderr, "%s: read error\n", zspec); exit(1); } - addscanline(y, scanin, zin); + addscanline(y, scanin, zin, plast, zlast); } /* clean up */ free((char *)scanin); free((char *)zin); + free((char *)plast); + free((char *)zlast); fclose(pfp); - fclose(zfp); + if (zfp != NULL) + fclose(zfp); } -addscanline(y, pline, zline) /* add scanline to output */ +pixform(xfmat, vw1, vw2) /* compute view1 to view2 matrix */ +register double xfmat[4][4]; +register VIEW *vw1, *vw2; +{ + double m4t[4][4]; + + setident4(xfmat); + xfmat[0][0] = vw1->vhinc[0]; + xfmat[0][1] = vw1->vhinc[1]; + xfmat[0][2] = vw1->vhinc[2]; + xfmat[1][0] = vw1->vvinc[0]; + xfmat[1][1] = vw1->vvinc[1]; + xfmat[1][2] = vw1->vvinc[2]; + xfmat[2][0] = vw1->vdir[0]; + xfmat[2][1] = vw1->vdir[1]; + xfmat[2][2] = vw1->vdir[2]; + xfmat[3][0] = vw1->vp[0]; + xfmat[3][1] = vw1->vp[1]; + xfmat[3][2] = vw1->vp[2]; + setident4(m4t); + m4t[0][0] = vw2->vhinc[0]/vw2->vhn2; + m4t[1][0] = vw2->vhinc[1]/vw2->vhn2; + m4t[2][0] = vw2->vhinc[2]/vw2->vhn2; + m4t[3][0] = -DOT(vw2->vp,vw2->vhinc)/vw2->vhn2; + m4t[0][1] = vw2->vvinc[0]/vw2->vvn2; + m4t[1][1] = vw2->vvinc[1]/vw2->vvn2; + m4t[2][1] = vw2->vvinc[2]/vw2->vvn2; + m4t[3][1] = -DOT(vw2->vp,vw2->vvinc)/vw2->vvn2; + m4t[0][2] = vw2->vdir[0]; + m4t[1][2] = vw2->vdir[1]; + m4t[2][2] = vw2->vdir[2]; + m4t[3][2] = -DOT(vw2->vp,vw2->vdir); + multmat4(xfmat, xfmat, m4t); +} + + +addscanline(y, pline, zline, lasty, lastyz) /* add scanline to output */ int y; COLR *pline; float *zline; +int *lasty; /* input/output */ +float *lastyz; /* input/output */ { - FVECT p, dir; - double xnew, ynew, znew; + extern double sqrt(); + double pos[3]; + int lastx = 0; + double lastxz = 0; + double zt; + int xpos, ypos; register int x; - register int xpos, ypos; - for (x = 0; x < theirview.hresolu; x++) { - rayview(p, dir, &theirview, x+.5, y+.5); - p[0] += zline[x]*dir[0]; - p[1] += zline[x]*dir[1]; - p[2] += zline[x]*dir[2]; - pixelview(&xnew, &ynew, &znew, &ourview, p); - if (znew <= 0.0 || xnew < 0 || xnew > ourview.hresolu - || ynew < 0 || ynew > ourview.vresolu) + for (x = theirview.hresolu-1; x >= 0; x--) { + pos[0] = x - .5*(theirview.hresolu-1); + pos[1] = y - .5*(theirview.vresolu-1); + pos[2] = zline[x]; + if (theirview.type == VT_PER) { + if (normdist) /* adjust for eye-ray distance */ + pos[2] /= sqrt( 1. + + pos[0]*pos[0]*theirview.vhn2 + + pos[1]*pos[1]*theirview.vvn2 ); + pos[0] *= pos[2]; + pos[1] *= pos[2]; + } + multp3(pos, pos, theirs2ours); + if (pos[2] <= 0) continue; - /* check current value at position */ - xpos = xnew; - ypos = ynew; - if (zscan(ypos)[xpos] <= 0.0 - || zscan(ypos)[xpos] - znew - > zeps*zscan(ypos)[xpos]) { - zscan(ypos)[xpos] = znew; - copycolr(pscan(ypos)[xpos], pline[x]); + if (ourview.type == VT_PER) { + pos[0] /= pos[2]; + pos[1] /= pos[2]; } + pos[0] += .5*ourview.hresolu; + pos[1] += .5*ourview.vresolu; + if (pos[0] < 0 || (xpos = pos[0]) >= ourview.hresolu + || pos[1] < 0 || (ypos = pos[1]) >= ourview.vresolu) + continue; + /* add pixel to our image */ + zt = 2.*zeps*zline[x]; + addpixel(xpos, ypos, + (fill&F_FORE && ABS(zline[x]-lastxz) <= zt) + ? lastx - xpos : 1, + (fill&F_FORE && ABS(zline[x]-lastyz[x]) <= zt) + ? lasty[x] - ypos : 1, + pline[x], pos[2]); + lastx = xpos; + lasty[x] = ypos; + lastxz = lastyz[x] = zline[x]; } } -fillpicture() /* fill in empty spaces */ +addpixel(xstart, ystart, width, height, pix, z) /* fill in area for pixel */ +int xstart, ystart; +int width, height; +COLR pix; +double z; { + register int x, y; + /* make width and height positive */ + if (width < 0) { + width = -width; + xstart = xstart-width+1; + } else if (width == 0) + width = 1; + if (height < 0) { + height = -height; + ystart = ystart-height+1; + } else if (height == 0) + height = 1; + /* fill pixel(s) within rectangle */ + for (y = ystart; y < ystart+height; y++) + for (x = xstart; x < xstart+width; x++) + if (zscan(y)[x] <= 0 + || zscan(y)[x]-z > zeps*zscan(y)[x]) { + zscan(y)[x] = z; + copycolr(pscan(y)[x], pix); + } +} + + +backpicture() /* background fill algorithm */ +{ int *yback, xback; int y; COLR pfill; @@ -261,50 +454,88 @@ fillpicture() /* fill in empty spaces */ return; } for (x = 0; x < ourview.hresolu; x++) - yback[x] = -1; + yback[x] = -2; + /* + * Xback and yback are the pixel locations of suitable + * background values in each direction. + * A value of -2 means unassigned, and -1 means + * that there is no suitable background in this direction. + */ /* fill image */ - for (y = 0; y < ourview.vresolu; y++) + for (y = 0; y < ourview.vresolu; y++) { + xback = -2; for (x = 0; x < ourview.hresolu; x++) - if (zscan(y)[x] <= 0.0) { /* found hole */ - xback = x-1; - do { /* find boundary */ - if (yback[x] < 0) { - for (i = y+1; - i < ourview.vresolu; - i++) - if (zscan(i)[x] > 0.0) - break; - if (i < ourview.vresolu + if (zscan(y)[x] <= 0) { /* empty pixel */ + /* + * First, find background from above or below. + * (farthest assigned pixel) + */ + if (yback[x] == -2) { + for (i = y+1; i < ourview.vresolu; i++) + if (zscan(i)[x] > 0) + break; + if (i < ourview.vresolu && (y <= 0 || zscan(y-1)[x] < zscan(i)[x])) - yback[x] = i; - else - yback[x] = y-1; - } - } while (++x < ourview.hresolu - && zscan(y)[x] <= 0.0); - i = xback; /* pick background */ - if (x < ourview.hresolu - && (i < 0 || zscan(y)[i] < zscan(y)[x])) - xback = x; - /* fill hole */ - if (xback < 0) { - while (++i < x) - if (yback[i] >= 0) - copycolr(pscan(y)[i],pscan(yback[i])[i]); + yback[x] = i; + else + yback[x] = y-1; + } + /* + * Next, find background from left or right. + */ + if (xback == -2) { + for (i = x+1; i < ourview.hresolu; i++) + if (zscan(y)[i] > 0) + break; + if (i < ourview.hresolu + && (x <= 0 || zscan(y)[x-1] < zscan(y)[i])) + xback = i; + else + xback = x-1; + } + /* + * Check to see if we have no background for + * this pixel. If not, use background color. + */ + if (xback < 0 && yback[x] < 0) { + (*deffill)(x,y); + continue; + } + /* + * Compare, and use the background that is + * farther, unless one of them is next to us. + */ + if ( yback[x] < 0 + || (xback >= 0 && ABS(x-xback) <= 1) + || ( ABS(y-yback[x]) > 1 + && zscan(yback[x])[x] + < zscan(y)[xback] ) ) { + copycolr(pscan(y)[x],pscan(y)[xback]); + zscan(y)[x] = zscan(y)[xback]; } else { - while (++i < x) - if (yback[i] < 0 || - zscan(yback[i])[i] < zscan(y)[xback]) - copycolr(pscan(y)[i],pscan(y)[xback]); - else - copycolr(pscan(y)[i],pscan(yback[i])[i]); + copycolr(pscan(y)[x],pscan(yback[x])[x]); + zscan(y)[x] = zscan(yback[x])[x]; } - } else - yback[x] = -1; /* clear boundary */ + } else { /* full pixel */ + yback[x] = -2; + xback = -2; + } + } free((char *)yback); } +fillpicture() /* paint in empty pixels with default */ +{ + register int x, y; + + for (y = 0; y < ourview.vresolu; y++) + for (x = 0; x < ourview.hresolu; x++) + if (zscan(y)[x] <= 0) + (*deffill)(x,y); +} + + writepicture() /* write out picture */ { int y; @@ -315,4 +546,174 @@ writepicture() /* write out picture */ perror(progname); exit(1); } +} + + +writedistance(fname) /* write out z file */ +char *fname; +{ + extern double sqrt(); + int donorm = normdist && ourview.type == VT_PER; + FILE *fp; + int y; + float *zout; + + if ((fp = fopen(fname, "w")) == NULL) { + perror(fname); + exit(1); + } + if (donorm + && (zout = (float *)malloc(ourview.hresolu*sizeof(float))) == NULL) { + perror(progname); + exit(1); + } + for (y = ourview.vresolu-1; y >= 0; y--) { + if (donorm) { + double vx, yzn2; + register int x; + yzn2 = y - .5*(ourview.vresolu-1); + yzn2 = 1. + yzn2*yzn2*ourview.vvn2; + for (x = 0; x < ourview.hresolu; x++) { + vx = x - .5*(ourview.hresolu-1); + zout[x] = zscan(y)[x] + * sqrt(vx*vx*ourview.vhn2 + yzn2); + } + } else + zout = zscan(y); + if (fwrite(zout, sizeof(float), ourview.hresolu, fp) + < ourview.hresolu) { + perror(fname); + exit(1); + } + } + if (donorm) + free((char *)zout); + fclose(fp); +} + + +isfloat(s) /* see if string is floating number */ +register char *s; +{ + for ( ; *s; s++) + if ((*s < '0' || *s > '9') && *s != '.' && *s != '-' + && *s != 'e' && *s != 'E' && *s != '+') + return(0); + return(1); +} + + +backfill(x, y) /* fill pixel with background */ +int x, y; +{ + register BYTE *dest = pscan(y)[x]; + + copycolr(dest, backcolr); + zscan(y)[x] = backz; +} + + +calstart(prog, args) /* start fill calculation */ +char *prog, *args; +{ + char combuf[512]; + int p0[2], p1[2]; + + if (childpid != -1) { + fprintf(stderr, "%s: too many calculations\n", progname); + exit(1); + } + sprintf(combuf, prog, PACKSIZ, args); + if (pipe(p0) < 0 || pipe(p1) < 0) + syserror(); + if ((childpid = vfork()) == 0) { /* fork calculation */ + close(p0[1]); + close(p1[0]); + if (p0[0] != 0) { + dup2(p0[0], 0); + close(p0[0]); + } + if (p1[1] != 1) { + dup2(p1[1], 1); + close(p1[1]); + } + execl("/bin/sh", "sh", "-c", combuf, 0); + perror("/bin/sh"); + _exit(127); + } + if (childpid == -1) + syserror(); + close(p0[0]); + close(p1[1]); + if ((psend = fdopen(p0[1], "w")) == NULL) + syserror(); + if ((precv = fdopen(p1[0], "r")) == NULL) + syserror(); + queuesiz = 0; +} + + +caldone() /* done with calculation */ +{ + int pid; + + if (childpid == -1) + return; + if (fclose(psend) == EOF) + syserror(); + clearqueue(); + fclose(precv); + while ((pid = wait(0)) != -1 && pid != childpid) + ; + childpid = -1; +} + + +rcalfill(x, y) /* fill with ray-calculated pixel */ +int x, y; +{ + FVECT orig, dir; + float outbuf[6]; + + if (queuesiz >= PACKSIZ) { /* flush queue */ + if (fflush(psend) == EOF) + syserror(); + clearqueue(); + } + /* send new ray */ + rayview(orig, dir, &ourview, x+.5, y+.5); + outbuf[0] = orig[0]; outbuf[1] = orig[1]; outbuf[2] = orig[2]; + outbuf[3] = dir[0]; outbuf[4] = dir[1]; outbuf[5] = dir[2]; + if (fwrite(outbuf, sizeof(float), 6, psend) < 6) + syserror(); + /* remember it */ + queue[queuesiz][0] = x; + queue[queuesiz][1] = y; + queuesiz++; +} + + +clearqueue() /* get results from queue */ +{ + float inbuf[4]; + register int i; + + for (i = 0; i < queuesiz; i++) { + if (fread(inbuf, sizeof(float), 4, precv) < 4) { + fprintf(stderr, "%s: read error in clearqueue\n", + progname); + exit(1); + } + setcolr(pscan(queue[i][1])[queue[i][0]], + inbuf[0], inbuf[1], inbuf[2]); + zscan(queue[i][1])[queue[i][0]] = inbuf[3]; + } + queuesiz = 0; +} + + +syserror() /* report error and exit */ +{ + perror(progname); + exit(1); }