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.26 by greg, Tue Mar 19 22:03:24 2019 UTC vs.
Revision 2.37 by greg, Tue Jun 3 21:31:51 2025 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 23 | Line 22 | static const char RCSid[] = "$Id$";
22   enum {CIE_X, CIE_Y, CIE_Z};
23                                  /* assumed maximum # Klems patches */
24   #define MAXPATCHES      145
26 –                                /* global argv[0] */
27 – char                    *progname;
25                                  /* selected basis function name */
26   static const char       klems_full[] = "LBNL/Klems Full";
27   static const char       klems_half[] = "LBNL/Klems Half";
# Line 197 | Line 194 | eval_bsdf(const char *fname)
194                      for (n = npsamps; n-- > 0; ) {
195                          fo_getvec(vout, j+(n+frandom())/npsamps, abp);
196                          fi_getvec(vin, i+urand(n), abp);
197 <                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
197 >                        ec = SDevalBSDF(&sdv, vin, vout, &bsd);
198                          if (ec != SDEnone)
199                                  goto err;
200                          sum += sdv.cieY;
# Line 246 | Line 243 | eval_bsdf(const char *fname)
243                      for (n = npsamps; n-- > 0; ) {
244                          bo_getvec(vout, j+(n+frandom())/npsamps, abp);
245                          bi_getvec(vin, i+urand(n), abp);
246 <                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
246 >                        ec = SDevalBSDF(&sdv, vin, vout, &bsd);
247                          if (ec != SDEnone)
248                                  goto err;
249                          sum += sdv.cieY;
# Line 278 | Line 275 | eval_bsdf(const char *fname)
275              }
276          }
277                                                  /* front transmission */
278 <        if (bsd.tf != NULL || bsd.tLamb.cieY > .002) {
278 >        if (bsd.tf != NULL || bsd.tLambFront.cieY > .002) {
279              input_orient = 1; output_orient = -1;
280              cfp[CIE_Y] = open_component_file(CIE_Y);
281              if (bsd.tf != NULL && bsd.tf->comp[0].cspec[2].flags) {
# Line 294 | Line 291 | eval_bsdf(const char *fname)
291                      for (n = npsamps; n-- > 0; ) {
292                          bo_getvec(vout, j+(n+frandom())/npsamps, abp);
293                          fi_getvec(vin, i+urand(n), abp);
294 <                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
294 >                        ec = SDevalBSDF(&sdv, vin, vout, &bsd);
295                          if (ec != SDEnone)
296                                  goto err;
297                          sum += sdv.cieY;
# Line 343 | Line 340 | eval_bsdf(const char *fname)
340                      for (n = npsamps; n-- > 0; ) {
341                          fo_getvec(vout, j+(n+frandom())/npsamps, abp);
342                          bi_getvec(vin, i+urand(n), abp);
343 <                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
343 >                        ec = SDevalBSDF(&sdv, vin, vout, &bsd);
344                          if (ec != SDEnone)
345                                  goto err;
346                          sum += sdv.cieY;
# Line 431 | Line 428 | eval_function(char *funame)
428   static void
429   eval_rbf(void)
430   {
431 <        ANGLE_BASIS     *abp = get_basis(kbasis);
432 <        float           (*XZarr)[2] = NULL;
433 <        float           bsdfarr[MAXPATCHES*MAXPATCHES];
434 <        FILE            *cfp[3];
435 <        FVECT           vin, vout;
436 <        double          sum, xsum, ysum;
437 <        int             i, j, n;
438 <                                                /* sanity check */
439 <        if (abp->nangles > MAXPATCHES) {
440 <                fprintf(stderr, "%s: too many patches!\n", progname);
441 <                exit(1);
442 <        }
443 <        if (rbf_colorimetry == RBCtristimulus)
444 <                XZarr = (float (*)[2])malloc(sizeof(float)*2*abp->nangles*abp->nangles);
445 <        for (i = 0; i < abp->nangles; i++) {
446 <            RBFNODE     *rbf;
447 <            if (input_orient > 0)               /* use incident patch center */
448 <                fi_getvec(vin, i+.5*(i>0), abp);
449 <            else
450 <                bi_getvec(vin, i+.5*(i>0), abp);
431 >    ANGLE_BASIS *abp = get_basis(kbasis);
432 >    float       (*XZarr)[2] = NULL;
433 >    float       bsdfarr[MAXPATCHES*MAXPATCHES];
434 >    FILE        *cfp[3];
435 >    FVECT       vin, vout;
436 >    double      sum, xsum, ysum, normf;
437 >    int         i, j, ni, no, nisamps, nosamps;
438 >                                        /* sanity check */
439 >    if (abp->nangles > MAXPATCHES) {
440 >        fprintf(stderr, "%s: too many patches!\n", progname);
441 >        exit(1);
442 >    }
443 >    memset(bsdfarr, 0, sizeof(bsdfarr));
444 >    if (rbf_colorimetry == RBCtristimulus)
445 >        XZarr = (float (*)[2])calloc(abp->nangles*abp->nangles, 2*sizeof(float));
446 >    nosamps = (int)(pow((double)npsamps, 0.67) + .5);
447 >    nisamps = (npsamps + (nosamps>>1)) / nosamps;
448 >    normf = 1./(double)(nisamps*nosamps);
449 >    for (i = 0; i < abp->nangles; i++) {
450 >        for (ni = nisamps; ni--; ) {            /* sample over incident patch */
451 >            RBFNODE     *rbf;
452 >            if (input_orient > 0)               /* vary incident patch loc. */
453 >                fi_getvec(vin, i+urand(ni), abp);
454 >            else
455 >                bi_getvec(vin, i+urand(ni), abp);
456  
457 <            rbf = advect_rbf(vin, lobe_lim);    /* compute radial basis func */
457 >            rbf = advect_rbf(vin, lobe_lim);    /* compute radial basis func */
458  
459 <            for (j = 0; j < abp->nangles; j++) {
460 <                sum = 0;                        /* sample over exiting patch */
459 >            for (j = 0; j < abp->nangles; j++) {
460 >                sum = 0;                        /* sample over exiting patch */
461                  xsum = ysum = 0;
462 <                for (n = npsamps; n--; ) {
462 >                for (no = nosamps; no--; ) {
463                      SDValue     sdv;
464                      if (output_orient > 0)
465 <                        fo_getvec(vout, j+(n+frandom())/npsamps, abp);
465 >                        fo_getvec(vout, j+(no+frandom())/nosamps, abp);
466                      else
467 <                        bo_getvec(vout, j+(n+frandom())/npsamps, abp);
467 >                        bo_getvec(vout, j+(no+frandom())/nosamps, abp);
468  
469                      eval_rbfcol(&sdv, rbf, vout);
470                      sum += sdv.cieY;
471                      if (rbf_colorimetry == RBCtristimulus) {
472                          xsum += sdv.cieY * sdv.spec.cx;
473                          ysum += sdv.cieY * sdv.spec.cy;
474 <                    }
474 >                    }
475                  }
476 <                n = j*abp->nangles + i;
477 <                bsdfarr[n] = sum / npsamps;
476 >                no = j*abp->nangles + i;
477 >                bsdfarr[no] += sum * normf;
478                  if (rbf_colorimetry == RBCtristimulus) {
479 <                    XZarr[n][0] = xsum*sum/(npsamps*ysum);
480 <                    XZarr[n][1] = (sum - xsum - ysum)*sum/(npsamps*ysum);
479 >                    XZarr[no][0] += xsum*sum*normf/ysum;
480 >                    XZarr[no][1] += (sum - xsum - ysum)*sum*normf/ysum;
481                  }
482              }
483 <            if (rbf != NULL)
483 >            if (rbf != NULL)
484                  free(rbf);
483 –            prog_show((i+1.)/abp->nangles);
485          }
486 <                                                /* write out our matrix */
487 <        cfp[CIE_Y] = open_component_file(CIE_Y);
488 <        n = 0;
489 <        for (j = 0; j < abp->nangles; j++) {
490 <            for (i = 0; i < abp->nangles; i++, n++)
491 <                fprintf(cfp[CIE_Y], "\t%.3e\n", bsdfarr[n]);
492 <            fputc('\n', cfp[CIE_Y]);
493 <        }
494 <        prog_done();
495 <        if (fclose(cfp[CIE_Y])) {
496 <                fprintf(stderr, "%s: error writing Y output\n", progname);
497 <                exit(1);
498 <        }
499 <        if (XZarr == NULL)                      /* no color? */
500 <                return;
501 <        cfp[CIE_X] = open_component_file(CIE_X);
502 <        cfp[CIE_Z] = open_component_file(CIE_Z);
503 <        n = 0;
504 <        for (j = 0; j < abp->nangles; j++) {
505 <            for (i = 0; i < abp->nangles; i++, n++) {
506 <                fprintf(cfp[CIE_X], "\t%.3e\n", XZarr[n][0]);
507 <                fprintf(cfp[CIE_Z], "\t%.3e\n", XZarr[n][1]);
508 <            }
509 <            fputc('\n', cfp[CIE_X]);
510 <            fputc('\n', cfp[CIE_Z]);
511 <        }
512 <        free(XZarr);
513 <        if (fclose(cfp[CIE_X]) || fclose(cfp[CIE_Z])) {
514 <                fprintf(stderr, "%s: error writing X/Z output\n", progname);
515 <                exit(1);
516 <        }
486 >        prog_show((i+1.)/abp->nangles);
487 >    }
488 >                                        /* write out our matrix */
489 >    cfp[CIE_Y] = open_component_file(CIE_Y);
490 >    no = 0;
491 >    for (j = 0; j < abp->nangles; j++) {
492 >        for (i = 0; i < abp->nangles; i++, no++)
493 >        fprintf(cfp[CIE_Y], "\t%.3e\n", bsdfarr[no]);
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 >    no = 0;
506 >    for (j = 0; j < abp->nangles; j++) {
507 >        for (i = 0; i < abp->nangles; i++, no++) {
508 >        fprintf(cfp[CIE_X], "\t%.3e\n", XZarr[no][0]);
509 >        fprintf(cfp[CIE_Z], "\t%.3e\n", XZarr[no][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 >    }
519   }
520  
521   #if defined(_WIN32) || defined(_WIN64)
# Line 520 | Line 523 | eval_rbf(void)
523   static int
524   wrap_up(void)
525   {
526 <        char    cmd[8192];
526 >        char    cmd[32700];
527  
528          if (bsdf_manuf[0]) {
529                  add_wbsdf("-f", 1);
# Line 571 | Line 574 | wrap_up(void)
574   }
575   #endif
576  
577 + #define HEAD_BUFLEN     10240
578 + static char     head_buf[HEAD_BUFLEN];
579 + static int      cur_headlen = 0;
580 +
581 + /* Record header line as comment associated with this SIR input */
582 + static int
583 + record2header(char *s)
584 + {
585 +        int     len = strlen(s);
586 +
587 +        if (cur_headlen+len >= HEAD_BUFLEN-6)
588 +                return(0);
589 +                                        /* includes EOL */
590 +        strcpy(head_buf+cur_headlen, s);
591 +        cur_headlen += len;
592 +
593 + #if defined(_WIN32) || defined(_WIN64)
594 +        if (head_buf[cur_headlen-1] == '\n')
595 +                head_buf[cur_headlen-1] = '\t';
596 + #endif
597 +        return(1);
598 + }
599 +
600 + /* Finish off header for this file */
601 + static void
602 + done_header(void)
603 + {
604 +        while (cur_headlen > 0 && isspace(head_buf[cur_headlen-1]))
605 +                --cur_headlen;
606 +        head_buf[cur_headlen] = '\0';
607 +        if (!cur_headlen)
608 +                return;
609 +        add_wbsdf("-C", 1);
610 +        add_wbsdf(head_buf, 0);
611 +        head_buf[cur_headlen=0] = '\0';
612 + }
613 +
614   /* Read in BSDF and interpolate as Klems matrix representation */
615   int
616   main(int argc, char *argv[])
617   {
618          int     dofwd = 0, dobwd = 1;
619 <        char    buf[2048];
619 >        char    buf[1024];
620          char    *cp;
621          int     i, na;
622 <
623 <        progname = argv[0];
622 >                                                /* set global progname */
623 >        fixargv0(argv[0]);
624          esupport |= E_VARIABLE|E_FUNCTION|E_RCONST;
625          esupport &= ~(E_INCHAN|E_OUTCHAN);
626          scompile("PI:3.14159265358979323846", NULL, 0);
# Line 597 | Line 637 | main(int argc, char *argv[])
637                          single_plane_incident = 0;
638                          break;
639                  case 'f':
640 <                        if (!argv[i][2]) {
640 >                        if ((argv[i][0] == '-') & !argv[i][2]) {
641                                  if (strchr(argv[++i], '=') != NULL) {
642                                          add_wbsdf("-f", 1);
643                                          add_wbsdf(argv[i], 1);
# Line 658 | Line 698 | main(int argc, char *argv[])
698                          fprintf(stderr, "\tor 3 arguments using Dx,Dy,Dz: bsdf(ix,iy,iz)\n");
699                          goto userr;
700                  }
701 +                doptimize(1);                   /* optimize definitions */
702                  ++eclock;
703                  if (dofwd) {
704                          input_orient = -1;
# Line 688 | Line 729 | main(int argc, char *argv[])
729          if (i < argc) {                         /* open input files if given */
730                  int     nbsdf = 0;
731                  for ( ; i < argc; i++) {        /* interpolate each component */
691 –                        char    pbuf[256];
732                          FILE    *fpin = fopen(argv[i], "rb");
733                          if (fpin == NULL) {
734                                  fprintf(stderr, "%s: cannot open BSDF interpolant '%s'\n",
735                                                  progname, argv[i]);
736                                  return(1);
737                          }
738 +                        sprintf(buf, "%s:\n", argv[i]);
739 +                        record2header(buf);
740 +                        sir_headshare = &record2header;
741                          if (!load_bsdf_rep(fpin))
742                                  return(1);
743                          fclose(fpin);
744 <                        sprintf(pbuf, "Interpolating component '%s'", argv[i]);
745 <                        prog_start(pbuf);
744 >                        done_header();
745 >                        sprintf(buf, "Interpolating component '%s'", argv[i]);
746 >                        prog_start(buf);
747                          eval_rbf();
748                  }
749                  return(wrap_up());
750          }
751          SET_FILE_BINARY(stdin);                 /* load from stdin */
752 +        record2header("<stdin>:\n");
753 +        sir_headshare = &record2header;
754          if (!load_bsdf_rep(stdin))
755                  return(1);
756 +        done_header();
757          prog_start("Interpolating from standard input");
758          eval_rbf();                             /* resample dist. */
759          return(wrap_up());

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines