--- ray/src/cv/bsdf2klems.c 2013/04/21 23:01:14 2.1 +++ ray/src/cv/bsdf2klems.c 2014/03/12 22:24:59 2.14 @@ -1,5 +1,5 @@ #ifndef lint -static const char RCSid[] = "$Id: bsdf2klems.c,v 2.1 2013/04/21 23:01:14 greg Exp $"; +static const char RCSid[] = "$Id: bsdf2klems.c,v 2.14 2014/03/12 22:24:59 greg Exp $"; #endif /* * Load measured BSDF interpolant and write out as XML file with Klems matrix. @@ -17,15 +17,58 @@ static const char RCSid[] = "$Id: bsdf2klems.c,v 2.1 2 #include "calcomp.h" #include "bsdfrep.h" #include "bsdf_m.h" + /* assumed maximum # Klems patches */ +#define MAXPATCHES 145 /* global argv[0] */ char *progname; /* selected basis function name */ static const char *kbasis = "LBNL/Klems Full"; /* number of BSDF samples per patch */ static int npsamps = 256; + /* limit on number of RBF lobes */ +static int lobe_lim = 15000; + /* progress bar length */ +static int do_prog = 79; + +/* Start new progress bar */ +#define prog_start(s) if (do_prog) fprintf(stderr, "%s: %s...\n", progname, s); else + +/* Draw progress bar of the appropriate length */ +static void +prog_show(double frac) +{ + char pbar[256]; + int nchars; + + if (do_prog <= 1) return; + if (do_prog > sizeof(pbar)-2) + do_prog = sizeof(pbar)-2; + if (frac < 0) frac = 0; + else if (frac > 1) frac = 1; + nchars = do_prog*frac + .5; + pbar[0] = '\r'; + memset(pbar+1, '*', nchars); + memset(pbar+1+nchars, '-', do_prog-nchars); + pbar[do_prog+1] = '\0'; + fputs(pbar, stderr); +} + +/* Finish progress bar */ +static void +prog_done(void) +{ + int n = do_prog; + + if (n <= 1) return; + fputc('\r', stderr); + while (n--) + fputc(' ', stderr); + fputc('\r', stderr); +} + /* Return angle basis corresponding to the given name */ -ANGLE_BASIS * +static ANGLE_BASIS * get_basis(const char *bn) { int n = nabases; @@ -36,17 +79,10 @@ get_basis(const char *bn) return NULL; } -/* Output XML prologue to stdout */ +/* Output XML header to stdout */ static void -xml_prologue(int ac, char *av[]) +xml_header(int ac, char *av[]) { - ANGLE_BASIS *abp = get_basis(kbasis); - int i; - - if (abp == NULL) { - fprintf(stderr, "%s: unknown angle basis '%s'\n", progname, kbasis); - exit(1); - } puts(""); puts(""); fputs(""); +} + +/* Output XML prologue to stdout */ +static void +xml_prologue(const SDData *sd) +{ + const char *matn = (sd && sd->matn[0]) ? sd->matn : + bsdf_name[0] ? bsdf_name : "Unknown"; + const char *makr = (sd && sd->makr[0]) ? sd->makr : + bsdf_manuf[0] ? bsdf_manuf : "Unknown"; + ANGLE_BASIS *abp = get_basis(kbasis); + int i; + + if (abp == NULL) { + fprintf(stderr, "%s: unknown angle basis '%s'\n", progname, kbasis); + exit(1); + } puts("System"); puts("BSDF"); puts(""); puts(""); puts("\t"); - puts("\t\tName"); - puts("\t\tManufacturer"); + printf("\t\t%s\n", matn); + printf("\t\t%s\n", makr); + if (sd && sd->dim[2] > .001) { + printf("\t\t%.3f\n", sd->dim[2]); + printf("\t\t%.3f\n", sd->dim[0]); + printf("\t\t%.3f\n", sd->dim[1]); + } puts("\t\tOther"); puts("\t"); + if (sd && sd->mgf != NULL) { + puts("\t"); + puts("\t\t"); + fputs(sd->mgf, stdout); + puts(""); + puts("\t"); + } puts("\t"); puts("\t\tColumns"); puts("\t\t"); @@ -73,7 +138,7 @@ xml_prologue(int ac, char *av[]) printf("\t\t\t%g\n", i ? .5*(abp->lat[i].tmin + abp->lat[i+1].tmin) : .0 ); - printf("\t\t\t%d", abp->lat[i].nphis); + printf("\t\t\t%d\n", abp->lat[i].nphis); puts("\t\t\t"); printf("\t\t\t\t%g\n", abp->lat[i].tmin); printf("\t\t\t\t%g\n", abp->lat[i+1].tmin); @@ -127,7 +192,7 @@ xml_epilogue(void) puts(""); } -/* Load and resample XML BSDF description */ +/* Load and resample XML BSDF description using Klems basis */ static void eval_bsdf(const char *fname) { @@ -142,6 +207,7 @@ eval_bsdf(const char *fname) SDclearBSDF(&bsd, fname); /* load BSDF file */ if ((ec = SDloadFile(&bsd, fname)) != SDEnone) goto err; + xml_prologue(&bsd); /* pass geometry */ /* front reflection */ if (bsd.rf != NULL || bsd.rLambFront.cieY > .002) { input_orient = 1; output_orient = 1; @@ -151,7 +217,7 @@ eval_bsdf(const char *fname) sum = 0; /* average over patches */ for (n = npsamps; n-- > 0; ) { fo_getvec(vout, j+(n+frandom())/npsamps, abp); - fi_getvec(vin, i+(n+frandom())/npsamps, abp); + fi_getvec(vin, i+urand(n), abp); ec = SDevalBSDF(&sv, vout, vin, &bsd); if (ec != SDEnone) goto err; @@ -172,7 +238,7 @@ eval_bsdf(const char *fname) sum = 0; /* average over patches */ for (n = npsamps; n-- > 0; ) { bo_getvec(vout, j+(n+frandom())/npsamps, abp); - bi_getvec(vin, i+(n+frandom())/npsamps, abp); + bi_getvec(vin, i+urand(n), abp); ec = SDevalBSDF(&sv, vout, vin, &bsd); if (ec != SDEnone) goto err; @@ -193,7 +259,7 @@ eval_bsdf(const char *fname) sum = 0; /* average over patches */ for (n = npsamps; n-- > 0; ) { bo_getvec(vout, j+(n+frandom())/npsamps, abp); - fi_getvec(vin, i+(n+frandom())/npsamps, abp); + fi_getvec(vin, i+urand(n), abp); ec = SDevalBSDF(&sv, vout, vin, &bsd); if (ec != SDEnone) goto err; @@ -206,15 +272,15 @@ eval_bsdf(const char *fname) data_epilogue(); } /* back transmission */ - if (bsd.tb != NULL) { + if ((bsd.tb != NULL) | (bsd.tf != NULL)) { input_orient = -1; output_orient = 1; data_prologue(); for (j = 0; j < abp->nangles; j++) { for (i = 0; i < abp->nangles; i++) { - sum = 0; /* average over patches */ + sum = 0; /* average over patches */ for (n = npsamps; n-- > 0; ) { fo_getvec(vout, j+(n+frandom())/npsamps, abp); - bi_getvec(vin, i+(n+frandom())/npsamps, abp); + bi_getvec(vin, i+urand(n), abp); ec = SDevalBSDF(&sv, vout, vin, &bsd); if (ec != SDEnone) goto err; @@ -233,15 +299,17 @@ err: exit(1); } -/* Interpolate and output a BSDF function */ +/* Interpolate and output a BSDF function using Klems basis */ static void eval_function(char *funame) { ANGLE_BASIS *abp = get_basis(kbasis); + int assignD = (fundefined(funame) < 6); double iovec[6]; double sum; int i, j, n; + initurand(npsamps); data_prologue(); /* begin output */ for (j = 0; j < abp->nangles; j++) { /* run through directions */ for (i = 0; i < abp->nangles; i++) { @@ -253,17 +321,25 @@ eval_function(char *funame) bo_getvec(iovec+3, j+(n+frandom())/npsamps, abp); if (input_orient > 0) - fi_getvec(iovec, j+(n+frandom())/npsamps, abp); + fi_getvec(iovec, i+urand(n), abp); else - bi_getvec(iovec, j+(n+frandom())/npsamps, abp); + bi_getvec(iovec, i+urand(n), abp); + if (assignD) { + varset("Dx", '=', -iovec[3]); + varset("Dy", '=', -iovec[4]); + varset("Dz", '=', -iovec[5]); + ++eclock; + } sum += funvalue(funame, 6, iovec); } printf("\t%.3e\n", sum/npsamps); } putchar('\n'); + prog_show((j+1.)/abp->nangles); } data_epilogue(); /* finish output */ + prog_done(); } /* Interpolate and output a radial basis function BSDF representation */ @@ -271,12 +347,50 @@ static void eval_rbf(void) { ANGLE_BASIS *abp = get_basis(kbasis); - double iovec[6]; + float bsdfarr[MAXPATCHES*MAXPATCHES]; + FVECT vin, vout; + RBFNODE *rbf; double sum; int i, j, n; + /* sanity check */ + if (abp->nangles > MAXPATCHES) { + fprintf(stderr, "%s: too many patches!\n", progname); + exit(1); + } + data_prologue(); /* begin output */ + for (i = 0; i < abp->nangles; i++) { + if (input_orient > 0) /* use incident patch center */ + fi_getvec(vin, i+.5*(i>0), abp); + else + bi_getvec(vin, i+.5*(i>0), abp); -fprintf(stder, "%s: RBF evaluation currently unimplemented\n", progname); -exit(1); + rbf = advect_rbf(vin, lobe_lim); /* compute radial basis func */ + + for (j = 0; j < abp->nangles; j++) { + sum = 0; /* sample over exiting patch */ + for (n = npsamps; n--; ) { + if (output_orient > 0) + fo_getvec(vout, j+(n+frandom())/npsamps, abp); + else + bo_getvec(vout, j+(n+frandom())/npsamps, abp); + + sum += eval_rbfrep(rbf, vout); + } + fo_getvec(vout, j+.5, abp); /* use centered secant */ + bsdfarr[j*abp->nangles + i] = sum / (npsamps*vout[2]); + } + if (rbf != NULL) + free(rbf); + prog_show((i+1.)/abp->nangles); + } + n = 0; /* write out our matrix */ + for (j = 0; j < abp->nangles; j++) { + for (i = 0; i < abp->nangles; i++) + printf("\t%.3e\n", bsdfarr[n++]); + putchar('\n'); + } + data_epilogue(); /* finish output */ + prog_done(); } /* Read in BSDF and interpolate as Klems matrix representation */ @@ -292,7 +406,7 @@ main(int argc, char *argv[]) esupport &= ~(E_INCHAN|E_OUTCHAN); scompile("PI:3.14159265358979323846", NULL, 0); biggerlib(); - for (i = 1; i < argc-1 && (argv[i][0] == '-') | (argv[i][0] == '+'); i++) + for (i = 1; i < argc && (argv[i][0] == '-') | (argv[i][0] == '+'); i++) switch (argv[i][1]) { /* get options */ case 'n': npsamps = atoi(argv[++i]); @@ -319,6 +433,12 @@ main(int argc, char *argv[]) case 'q': kbasis = "LBNL/Klems Quarter"; break; + case 'l': + lobe_lim = atoi(argv[++i]); + break; + case 'p': + do_prog = atoi(argv[i]+2); + break; default: goto userr; } @@ -327,29 +447,37 @@ main(int argc, char *argv[]) fprintf(stderr, "%s: need single function with 6 arguments: bsdf(ix,iy,iz,ox,oy,oz)\n", progname); + fprintf(stderr, "\tor 3 arguments using Dx,Dy,Dz: bsdf(ix,iy,iz)\n"); goto userr; } - xml_prologue(argc, argv); /* start XML output */ + ++eclock; + xml_header(argc, argv); /* start XML output */ + xml_prologue(NULL); if (dofwd) { input_orient = -1; output_orient = -1; - eval_function(argv[i]); /* outside reflectance */ + prog_start("Evaluating outside reflectance"); + eval_function(argv[i]); output_orient = 1; - eval_function(argv[i]); /* outside -> inside */ + prog_start("Evaluating outside->inside transmission"); + eval_function(argv[i]); } if (dobwd) { input_orient = 1; output_orient = 1; - eval_function(argv[i]); /* inside reflectance */ + prog_start("Evaluating inside reflectance"); + eval_function(argv[i]); output_orient = -1; - eval_function(argv[i]); /* inside -> outside */ + prog_start("Evaluating inside->outside transmission"); + eval_function(argv[i]); } xml_epilogue(); /* finish XML output & exit */ return(0); } - if (i == argc-1 && (cp = strstr(argv[i], ".xml")) != NULL && - strlen(cp) == 4) { /* XML input? */ - xml_prologue(argc, argv); /* start XML output */ + /* XML input? */ + if (i == argc-1 && (cp = argv[i]+strlen(argv[i])-4) > argv[i] && + !strcasecmp(cp, ".xml")) { + xml_header(argc, argv); /* start XML output */ eval_bsdf(argv[i]); /* load & resample BSDF */ xml_epilogue(); /* finish XML output & exit */ return(0); @@ -357,6 +485,7 @@ main(int argc, char *argv[]) if (i < argc) { /* open input files if given */ int nbsdf = 0; for ( ; i < argc; i++) { /* interpolate each component */ + char pbuf[256]; FILE *fpin = fopen(argv[i], "rb"); if (fpin == NULL) { fprintf(stderr, "%s: cannot open BSDF interpolant '%s'\n", @@ -366,8 +495,12 @@ main(int argc, char *argv[]) if (!load_bsdf_rep(fpin)) return(1); fclose(fpin); - if (!nbsdf++) /* start XML on first dist. */ - xml_prologue(argc, argv); + if (!nbsdf++) { /* start XML on first dist. */ + xml_header(argc, argv); + xml_prologue(NULL); + } + sprintf(pbuf, "Interpolating component '%s'", argv[i]); + prog_start(pbuf); eval_rbf(); } xml_epilogue(); /* finish XML output & exit */ @@ -376,17 +509,17 @@ main(int argc, char *argv[]) SET_FILE_BINARY(stdin); /* load from stdin */ if (!load_bsdf_rep(stdin)) return(1); - xml_prologue(argc, argv); /* start XML output */ + xml_header(argc, argv); /* start XML output */ + xml_prologue(NULL); + prog_start("Interpolating from standard input"); eval_rbf(); /* resample dist. */ xml_epilogue(); /* finish XML output & exit */ return(0); userr: fprintf(stderr, - "Usage: %s [-n spp][-h|-q][bsdf.sir ..] > bsdf.xml\n", - progname); + "Usage: %s [-n spp][-h|-q][-l maxlobes] [bsdf.sir ..] > bsdf.xml\n", progname); fprintf(stderr, - " or: %s [-n spp][-h|-q] bsdf_in.xml > bsdf_out.xml\n", - progname); + " or: %s [-n spp][-h|-q] bsdf_in.xml > bsdf_out.xml\n", progname); fprintf(stderr, " or: %s [-n spp][-h|-q][{+|-}for[ward]][{+|-}b[ackward]][-e expr][-f file] bsdf_func > bsdf.xml\n", progname);