ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/cv/bsdf2klems.c
(Generate patch)

Comparing ray/src/cv/bsdf2klems.c (file contents):
Revision 2.20 by greg, Wed Feb 3 01:57:06 2016 UTC vs.
Revision 2.36 by greg, Fri Feb 23 03:47:57 2024 UTC

# Line 8 | Line 8 | static const char RCSid[] = "$Id$";
8   */
9  
10   #define _USE_MATH_DEFINES
11 #include <stdio.h>
11   #include <stdlib.h>
13 #include <string.h>
12   #include <math.h>
13 + #include <ctype.h>
14   #include "random.h"
15   #include "platform.h"
16   #include "paths.h"
# Line 31 | Line 30 | static const char      klems_half[] = "LBNL/Klems Half";
30   static const char       klems_quarter[] = "LBNL/Klems Quarter";
31   static const char       *kbasis = klems_full;
32                                  /* number of BSDF samples per patch */
33 < static int              npsamps = 256;
33 > static int              npsamps = 1024;
34                                  /* limit on number of RBF lobes */
35   static int              lobe_lim = 15000;
36                                  /* progress bar length */
# Line 41 | Line 40 | static int             do_prog = 79;
40   static char             *wrapBSDF[MAXCARG] = {"wrapBSDF", "-W", "-UU"};
41   static int              wbsdfac = 3;
42  
43 < /* Add argument to wrapBSDF, allocating space if isstatic */
43 > /* Add argument to wrapBSDF, allocating space if !isstatic */
44   static void
45   add_wbsdf(const char *arg, int isstatic)
46   {
# Line 197 | Line 196 | eval_bsdf(const char *fname)
196                      for (n = npsamps; n-- > 0; ) {
197                          fo_getvec(vout, j+(n+frandom())/npsamps, abp);
198                          fi_getvec(vin, i+urand(n), abp);
199 <                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
199 >                        ec = SDevalBSDF(&sdv, vin, vout, &bsd);
200                          if (ec != SDEnone)
201                                  goto err;
202                          sum += sdv.cieY;
203                          if (rbf_colorimetry == RBCtristimulus) {
204 <                                c_ccvt(&sdv.spec, C_CSXY);
205 <                                xsum += sdv.cieY*sdv.spec.cx;
207 <                                ysum += sdv.cieY*sdv.spec.cy;
204 >                                xsum += sdv.cieY * sdv.spec.cx;
205 >                                ysum += sdv.cieY * sdv.spec.cy;
206                          }
207                      }
208                      fprintf(cfp[CIE_Y], "\t%.3e\n", sum/npsamps);
209                      if (rbf_colorimetry == RBCtristimulus) {
210 <                        fprintf(cfp[CIE_X], "\t%3e\n", xsum*sum/(npsamps*ysum));
211 <                        fprintf(cfp[CIE_Z], "\t%3e\n",
210 >                        fprintf(cfp[CIE_X], "\t%.3e\n", xsum*sum/(npsamps*ysum));
211 >                        fprintf(cfp[CIE_Z], "\t%.3e\n",
212                                  (sum - xsum - ysum)*sum/(npsamps*ysum));
213                      }
214                  }
# Line 247 | Line 245 | eval_bsdf(const char *fname)
245                      for (n = npsamps; n-- > 0; ) {
246                          bo_getvec(vout, j+(n+frandom())/npsamps, abp);
247                          bi_getvec(vin, i+urand(n), abp);
248 <                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
248 >                        ec = SDevalBSDF(&sdv, vin, vout, &bsd);
249                          if (ec != SDEnone)
250                                  goto err;
251                          sum += sdv.cieY;
252                          if (rbf_colorimetry == RBCtristimulus) {
253 <                                c_ccvt(&sdv.spec, C_CSXY);
254 <                                xsum += sdv.cieY*sdv.spec.cx;
257 <                                ysum += sdv.cieY*sdv.spec.cy;
253 >                                xsum += sdv.cieY * sdv.spec.cx;
254 >                                ysum += sdv.cieY * sdv.spec.cy;
255                          }
256                      }
257                      fprintf(cfp[CIE_Y], "\t%.3e\n", sum/npsamps);
258                      if (rbf_colorimetry == RBCtristimulus) {
259 <                        fprintf(cfp[CIE_X], "\t%3e\n", xsum*sum/(npsamps*ysum));
260 <                        fprintf(cfp[CIE_Z], "\t%3e\n",
259 >                        fprintf(cfp[CIE_X], "\t%.3e\n", xsum*sum/(npsamps*ysum));
260 >                        fprintf(cfp[CIE_Z], "\t%.3e\n",
261                                  (sum - xsum - ysum)*sum/(npsamps*ysum));
262                      }
263                  }
# Line 280 | Line 277 | eval_bsdf(const char *fname)
277              }
278          }
279                                                  /* front transmission */
280 <        if (bsd.tf != NULL || bsd.tLamb.cieY > .002) {
280 >        if (bsd.tf != NULL || bsd.tLambFront.cieY > .002) {
281              input_orient = 1; output_orient = -1;
282              cfp[CIE_Y] = open_component_file(CIE_Y);
283              if (bsd.tf != NULL && bsd.tf->comp[0].cspec[2].flags) {
# Line 296 | Line 293 | eval_bsdf(const char *fname)
293                      for (n = npsamps; n-- > 0; ) {
294                          bo_getvec(vout, j+(n+frandom())/npsamps, abp);
295                          fi_getvec(vin, i+urand(n), abp);
296 <                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
296 >                        ec = SDevalBSDF(&sdv, vin, vout, &bsd);
297                          if (ec != SDEnone)
298                                  goto err;
299                          sum += sdv.cieY;
300                          if (rbf_colorimetry == RBCtristimulus) {
301 <                                c_ccvt(&sdv.spec, C_CSXY);
302 <                                xsum += sdv.cieY*sdv.spec.cx;
306 <                                ysum += sdv.cieY*sdv.spec.cy;
301 >                                xsum += sdv.cieY * sdv.spec.cx;
302 >                                ysum += sdv.cieY * sdv.spec.cy;
303                          }
304                      }
305                      fprintf(cfp[CIE_Y], "\t%.3e\n", sum/npsamps);
306                      if (rbf_colorimetry == RBCtristimulus) {
307 <                        fprintf(cfp[CIE_X], "\t%3e\n", xsum*sum/(npsamps*ysum));
308 <                        fprintf(cfp[CIE_Z], "\t%3e\n",
307 >                        fprintf(cfp[CIE_X], "\t%.3e\n", xsum*sum/(npsamps*ysum));
308 >                        fprintf(cfp[CIE_Z], "\t%.3e\n",
309                                  (sum - xsum - ysum)*sum/(npsamps*ysum));
310                      }
311                  }
# Line 332 | Line 328 | eval_bsdf(const char *fname)
328          if ((bsd.tb != NULL) | (bsd.tf != NULL)) {
329              input_orient = -1; output_orient = 1;
330              cfp[CIE_Y] = open_component_file(CIE_Y);
331 <            if (bsd.tb != NULL && bsd.tb->comp[0].cspec[2].flags) {
332 <                rbf_colorimetry = RBCtristimulus;
331 >            if (bsd.tb != NULL)
332 >                rbf_colorimetry = bsd.tb->comp[0].cspec[2].flags
333 >                                        ? RBCtristimulus : RBCphotopic ;
334 >            if (rbf_colorimetry == RBCtristimulus) {
335                  cfp[CIE_X] = open_component_file(CIE_X);
336                  cfp[CIE_Z] = open_component_file(CIE_Z);
337 <            } else
340 <                rbf_colorimetry = RBCphotopic;
337 >            }
338              for (j = 0; j < abp->nangles; j++) {
339                  for (i = 0; i < abp->nangles; i++) {
340                      sum = 0;            /* average over patches */
# Line 345 | Line 342 | eval_bsdf(const char *fname)
342                      for (n = npsamps; n-- > 0; ) {
343                          fo_getvec(vout, j+(n+frandom())/npsamps, abp);
344                          bi_getvec(vin, i+urand(n), abp);
345 <                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
345 >                        ec = SDevalBSDF(&sdv, vin, vout, &bsd);
346                          if (ec != SDEnone)
347                                  goto err;
348                          sum += sdv.cieY;
349                          if (rbf_colorimetry == RBCtristimulus) {
350 <                                c_ccvt(&sdv.spec, C_CSXY);
351 <                                xsum += sdv.cieY*sdv.spec.cx;
355 <                                ysum += sdv.cieY*sdv.spec.cy;
350 >                                xsum += sdv.cieY * sdv.spec.cx;
351 >                                ysum += sdv.cieY * sdv.spec.cy;
352                          }
353                      }
354                      fprintf(cfp[CIE_Y], "\t%.3e\n", sum/npsamps);
355                      if (rbf_colorimetry == RBCtristimulus) {
356 <                        fprintf(cfp[CIE_X], "\t%3e\n", xsum*sum/(npsamps*ysum));
357 <                        fprintf(cfp[CIE_Z], "\t%3e\n",
356 >                        fprintf(cfp[CIE_X], "\t%.3e\n", xsum*sum/(npsamps*ysum));
357 >                        fprintf(cfp[CIE_Z], "\t%.3e\n",
358                                  (sum - xsum - ysum)*sum/(npsamps*ysum));
359                      }
360                  }
# Line 434 | Line 430 | eval_function(char *funame)
430   static void
431   eval_rbf(void)
432   {
433 <        ANGLE_BASIS     *abp = get_basis(kbasis);
434 <        float           (*XZarr)[2] = NULL;
435 <        float           bsdfarr[MAXPATCHES*MAXPATCHES];
436 <        FILE            *cfp[3];
437 <        FVECT           vin, vout;
438 <        double          sum, xsum, ysum;
439 <        int             i, j, n;
440 <                                                /* sanity check */
441 <        if (abp->nangles > MAXPATCHES) {
442 <                fprintf(stderr, "%s: too many patches!\n", progname);
443 <                exit(1);
444 <        }
445 <        if (rbf_colorimetry == RBCtristimulus)
446 <                XZarr = (float (*)[2])malloc(sizeof(float)*2*abp->nangles*abp->nangles);
447 <        for (i = 0; i < abp->nangles; i++) {
448 <            RBFNODE     *rbf;
449 <            if (input_orient > 0)               /* use incident patch center */
450 <                fi_getvec(vin, i+.5*(i>0), abp);
451 <            else
452 <                bi_getvec(vin, i+.5*(i>0), abp);
433 >    ANGLE_BASIS *abp = get_basis(kbasis);
434 >    float       (*XZarr)[2] = NULL;
435 >    float       bsdfarr[MAXPATCHES*MAXPATCHES];
436 >    FILE        *cfp[3];
437 >    FVECT       vin, vout;
438 >    double      sum, xsum, ysum, normf;
439 >    int         i, j, ni, no, nisamps, nosamps;
440 >                                        /* sanity check */
441 >    if (abp->nangles > MAXPATCHES) {
442 >        fprintf(stderr, "%s: too many patches!\n", progname);
443 >        exit(1);
444 >    }
445 >    memset(bsdfarr, 0, sizeof(bsdfarr));
446 >    if (rbf_colorimetry == RBCtristimulus)
447 >        XZarr = (float (*)[2])calloc(abp->nangles*abp->nangles, 2*sizeof(float));
448 >    nosamps = (int)(pow((double)npsamps, 0.67) + .5);
449 >    nisamps = (npsamps + (nosamps>>1)) / nosamps;
450 >    normf = 1./(double)(nisamps*nosamps);
451 >    for (i = 0; i < abp->nangles; i++) {
452 >        for (ni = nisamps; ni--; ) {            /* sample over incident patch */
453 >            RBFNODE     *rbf;
454 >            if (input_orient > 0)               /* vary incident patch loc. */
455 >                fi_getvec(vin, i+urand(ni), abp);
456 >            else
457 >                bi_getvec(vin, i+urand(ni), abp);
458  
459 <            rbf = advect_rbf(vin, lobe_lim);    /* compute radial basis func */
459 >            rbf = advect_rbf(vin, lobe_lim);    /* compute radial basis func */
460  
461 <            for (j = 0; j < abp->nangles; j++) {
462 <                sum = 0;                        /* sample over exiting patch */
461 >            for (j = 0; j < abp->nangles; j++) {
462 >                sum = 0;                        /* sample over exiting patch */
463                  xsum = ysum = 0;
464 <                for (n = npsamps; n--; ) {
464 >                for (no = nosamps; no--; ) {
465                      SDValue     sdv;
466                      if (output_orient > 0)
467 <                        fo_getvec(vout, j+(n+frandom())/npsamps, abp);
467 >                        fo_getvec(vout, j+(no+frandom())/nosamps, abp);
468                      else
469 <                        bo_getvec(vout, j+(n+frandom())/npsamps, abp);
469 >                        bo_getvec(vout, j+(no+frandom())/nosamps, abp);
470  
471                      eval_rbfcol(&sdv, rbf, vout);
472                      sum += sdv.cieY;
473 <                    if (XZarr != NULL) {
474 <                        c_ccvt(&sdv.spec, C_CSXY);
475 <                        xsum += sdv.cieY*sdv.spec.cx;
476 <                        ysum += sdv.cieY*sdv.spec.cy;
476 <                    }
473 >                    if (rbf_colorimetry == RBCtristimulus) {
474 >                        xsum += sdv.cieY * sdv.spec.cx;
475 >                        ysum += sdv.cieY * sdv.spec.cy;
476 >                    }
477                  }
478 <                n = j*abp->nangles + i;
479 <                bsdfarr[n] = sum / (double)npsamps;
480 <                if (XZarr != NULL) {
481 <                    XZarr[n][0] = xsum*sum/(npsamps*ysum);
482 <                    XZarr[n][1] = (sum - xsum - ysum)*sum/(npsamps*ysum);
478 >                no = j*abp->nangles + i;
479 >                bsdfarr[no] += sum * normf;
480 >                if (rbf_colorimetry == RBCtristimulus) {
481 >                    XZarr[no][0] += xsum*sum*normf/ysum;
482 >                    XZarr[no][1] += (sum - xsum - ysum)*sum*normf/ysum;
483                  }
484              }
485 <            if (rbf != NULL)
485 >            if (rbf != NULL)
486                  free(rbf);
487            prog_show((i+1.)/abp->nangles);
487          }
488 <                                                /* write out our matrix */
489 <        cfp[CIE_Y] = open_component_file(CIE_Y);
490 <        n = 0;
491 <        for (j = 0; j < abp->nangles; j++) {
492 <            for (i = 0; i < abp->nangles; i++, n++)
493 <                fprintf(cfp[CIE_Y], "\t%.3e\n", bsdfarr[n]);
494 <            fputc('\n', cfp[CIE_Y]);
495 <        }
496 <        prog_done();
497 <        if (fclose(cfp[CIE_Y])) {
498 <                fprintf(stderr, "%s: error writing Y output\n", progname);
499 <                exit(1);
500 <        }
501 <        if (XZarr == NULL)                      /* no color? */
502 <                return;
503 <        cfp[CIE_X] = open_component_file(CIE_X);
504 <        cfp[CIE_Z] = open_component_file(CIE_Z);
505 <        n = 0;
506 <        for (j = 0; j < abp->nangles; j++) {
507 <            for (i = 0; i < abp->nangles; i++, n++) {
508 <                fprintf(cfp[CIE_X], "\t%.3e\n", XZarr[n][0]);
509 <                fprintf(cfp[CIE_Z], "\t%.3e\n", XZarr[n][1]);
510 <            }
511 <            fputc('\n', cfp[CIE_X]);
512 <            fputc('\n', cfp[CIE_Z]);
513 <        }
514 <        free(XZarr);
515 <        if (fclose(cfp[CIE_X]) || fclose(cfp[CIE_Z])) {
516 <                fprintf(stderr, "%s: error writing X/Z output\n", progname);
517 <                exit(1);
518 <        }
488 >        prog_show((i+1.)/abp->nangles);
489 >    }
490 >                                        /* write out our matrix */
491 >    cfp[CIE_Y] = open_component_file(CIE_Y);
492 >    no = 0;
493 >    for (j = 0; j < abp->nangles; j++) {
494 >        for (i = 0; i < abp->nangles; i++, no++)
495 >        fprintf(cfp[CIE_Y], "\t%.3e\n", bsdfarr[no]);
496 >        fputc('\n', cfp[CIE_Y]);
497 >    }
498 >    prog_done();
499 >    if (fclose(cfp[CIE_Y])) {
500 >        fprintf(stderr, "%s: error writing Y output\n", progname);
501 >        exit(1);
502 >    }
503 >    if (XZarr == NULL)                  /* no color? */
504 >        return;
505 >    cfp[CIE_X] = open_component_file(CIE_X);
506 >    cfp[CIE_Z] = open_component_file(CIE_Z);
507 >    no = 0;
508 >    for (j = 0; j < abp->nangles; j++) {
509 >        for (i = 0; i < abp->nangles; i++, no++) {
510 >        fprintf(cfp[CIE_X], "\t%.3e\n", XZarr[no][0]);
511 >        fprintf(cfp[CIE_Z], "\t%.3e\n", XZarr[no][1]);
512 >        }
513 >        fputc('\n', cfp[CIE_X]);
514 >        fputc('\n', cfp[CIE_Z]);
515 >    }
516 >    free(XZarr);
517 >    if (fclose(cfp[CIE_X]) || fclose(cfp[CIE_Z])) {
518 >        fprintf(stderr, "%s: error writing X/Z output\n", progname);
519 >        exit(1);
520 >    }
521   }
522  
523 < #ifdef _WIN32
523 > #if defined(_WIN32) || defined(_WIN64)
524   /* Execute wrapBSDF command (may never return) */
525   static int
526   wrap_up(void)
527   {
528 <        char    cmd[8192];
528 >        char    cmd[32700];
529  
530          if (bsdf_manuf[0]) {
531                  add_wbsdf("-f", 1);
# Line 575 | Line 576 | wrap_up(void)
576   }
577   #endif
578  
579 + #define HEAD_BUFLEN     10240
580 + static char     head_buf[HEAD_BUFLEN];
581 + static int      cur_headlen = 0;
582 +
583 + /* Record header line as comment associated with this SIR input */
584 + static int
585 + record2header(char *s)
586 + {
587 +        int     len = strlen(s);
588 +
589 +        if (cur_headlen+len >= HEAD_BUFLEN-6)
590 +                return(0);
591 +                                        /* includes EOL */
592 +        strcpy(head_buf+cur_headlen, s);
593 +        cur_headlen += len;
594 +
595 + #if defined(_WIN32) || defined(_WIN64)
596 +        if (head_buf[cur_headlen-1] == '\n')
597 +                head_buf[cur_headlen-1] = '\t';
598 + #endif
599 +        return(1);
600 + }
601 +
602 + /* Finish off header for this file */
603 + static void
604 + done_header(void)
605 + {
606 +        while (cur_headlen > 0 && isspace(head_buf[cur_headlen-1]))
607 +                --cur_headlen;
608 +        head_buf[cur_headlen] = '\0';
609 +        if (!cur_headlen)
610 +                return;
611 +        add_wbsdf("-C", 1);
612 +        add_wbsdf(head_buf, 0);
613 +        head_buf[cur_headlen=0] = '\0';
614 + }
615 +
616   /* Read in BSDF and interpolate as Klems matrix representation */
617   int
618   main(int argc, char *argv[])
619   {
620          int     dofwd = 0, dobwd = 1;
621 <        char    buf[2048];
621 >        char    buf[1024];
622          char    *cp;
623          int     i, na;
624  
# Line 601 | Line 639 | main(int argc, char *argv[])
639                          single_plane_incident = 0;
640                          break;
641                  case 'f':
642 <                        if (!argv[i][2]) {
642 >                        if ((argv[i][0] == '-') & !argv[i][2]) {
643                                  if (strchr(argv[++i], '=') != NULL) {
644                                          add_wbsdf("-f", 1);
645                                          add_wbsdf(argv[i], 1);
646                                  } else {
647 <                                        fcompile(argv[i]);
647 >                                        char    *fpath = getpath(argv[i],
648 >                                                            getrlibpath(), 0);
649 >                                        if (fpath == NULL) {
650 >                                                fprintf(stderr,
651 >                                                "%s: cannot find file '%s'\n",
652 >                                                        argv[0], argv[i]);
653 >                                                return(1);
654 >                                        }
655 >                                        fcompile(fpath);
656                                          single_plane_incident = 0;
657                                  }
658                          } else
# Line 654 | Line 700 | main(int argc, char *argv[])
700                          fprintf(stderr, "\tor 3 arguments using Dx,Dy,Dz: bsdf(ix,iy,iz)\n");
701                          goto userr;
702                  }
703 +                doptimize(1);                   /* optimize definitions */
704                  ++eclock;
705                  if (dofwd) {
706                          input_orient = -1;
# Line 684 | Line 731 | main(int argc, char *argv[])
731          if (i < argc) {                         /* open input files if given */
732                  int     nbsdf = 0;
733                  for ( ; i < argc; i++) {        /* interpolate each component */
687                        char    pbuf[256];
734                          FILE    *fpin = fopen(argv[i], "rb");
735                          if (fpin == NULL) {
736                                  fprintf(stderr, "%s: cannot open BSDF interpolant '%s'\n",
737                                                  progname, argv[i]);
738                                  return(1);
739                          }
740 +                        sprintf(buf, "%s:\n", argv[i]);
741 +                        record2header(buf);
742 +                        sir_headshare = &record2header;
743                          if (!load_bsdf_rep(fpin))
744                                  return(1);
745                          fclose(fpin);
746 <                        sprintf(pbuf, "Interpolating component '%s'", argv[i]);
747 <                        prog_start(pbuf);
746 >                        done_header();
747 >                        sprintf(buf, "Interpolating component '%s'", argv[i]);
748 >                        prog_start(buf);
749                          eval_rbf();
750                  }
751                  return(wrap_up());
752          }
753          SET_FILE_BINARY(stdin);                 /* load from stdin */
754 +        record2header("<stdin>:\n");
755 +        sir_headshare = &record2header;
756          if (!load_bsdf_rep(stdin))
757                  return(1);
758 +        done_header();
759          prog_start("Interpolating from standard input");
760          eval_rbf();                             /* resample dist. */
761          return(wrap_up());

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines