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

Comparing ray/src/gen/mkillum2.c (file contents):
Revision 2.29 by greg, Wed Mar 4 00:12:25 2009 UTC vs.
Revision 2.37 by greg, Tue Aug 16 18:09:53 2011 UTC

# Line 11 | Line 11 | static const char      RCSid[] = "$Id$";
11   #include  "face.h"
12   #include  "cone.h"
13   #include  "source.h"
14 + #include  "paths.h"
15  
16 + #ifndef NBSDFSAMPS
17 + #define NBSDFSAMPS      256             /* BSDF resampling count */
18 + #endif
19  
20   COLORV *        distarr = NULL;         /* distribution array */
21   int             distsiz = 0;
22   COLORV *        direct_discount = NULL; /* amount to take off direct */
23  
24 +
25   void
26   newdist(                        /* allocate & clear distribution array */
27          int siz
# Line 24 | Line 29 | newdist(                       /* allocate & clear distribution array */
29   {
30          if (siz <= 0) {
31                  if (distsiz > 0)
32 <                        free((void *)distarr);
32 >                        free(distarr);
33                  distarr = NULL;
34                  distsiz = 0;
35                  return;
36          }
37          if (distsiz < siz) {
38                  if (distsiz > 0)
39 <                        free((void *)distarr);
39 >                        free(distarr);
40                  distarr = (COLORV *)malloc(sizeof(COLOR)*siz);
41                  if (distarr == NULL)
42                          error(SYSTEM, "out of memory in newdist");
# Line 57 | Line 62 | done_discount()                        /* clear off direct contrib. record
62   {
63          if (direct_discount == NULL)
64                  return;
65 <        free((void *)direct_discount);
65 >        free(direct_discount);
66          direct_discount = NULL;
67   }
68  
# Line 103 | Line 108 | raysamp(                       /* queue a ray sample */
108          VCOPY(myRay.rorg, org);
109          VCOPY(myRay.rdir, dir);
110          myRay.rmax = .0;
111 <        rayorigin(&myRay, PRIMARY, NULL, NULL);
111 >        rayorigin(&myRay, PRIMARY|SPECULAR, NULL, NULL);
112          myRay.rno = ndx;
113                                          /* queue ray, check result */
114          process_ray(&myRay, ray_pqueue(&myRay));
# Line 118 | Line 123 | srcsamps(                      /* sample sources from this surface positi
123          MAT4 ixfm
124   )
125   {
126 <        int  nalt, nazi;
126 >        int  nalt=1, nazi=1;
127          SRCINDEX  si;
128          RAY  sr;
129          FVECT   v;
130          double  d;
131          int  i, j;
132                                                  /* get sampling density */
133 <        if (il->sampdens <= 0) {
129 <                nalt = nazi = 1;
130 <        } else {
133 >        if (il->sd == NULL && il->sampdens > 0) {
134                  i = PI * il->sampdens;
135                  nalt = sqrt(i/PI) + .5;
136                  nazi = PI*nalt + .5;
# Line 151 | Line 154 | srcsamps(                      /* sample sources from this surface positi
154                          d = -1.0001*il->thick - 5.*FTINY;
155                  else
156                          d = 5.*FTINY;
157 <                for (i = 3; i--; )
158 <                        sr.rorg[i] += d*nrm[i];
157 >                VSUM(sr.rorg, sr.rorg, nrm, d);
158 >                samplendx++;                    /* increment sample counter */
159                  if (!srcray(&sr, NULL, &si))
160                          break;                  /* end of sources */
161                                                  /* index direction */
# Line 243 | Line 246 | flatdir(               /* compute uniform hemispherical direction *
246          dv[2] = sqrt(1. - alt);
247   }
248  
249 +
250   int
251   flatindex(              /* compute index for hemispherical direction */
252          FVECT   dv,
# Line 264 | Line 268 | flatindex(             /* compute index for hemispherical directi
268  
269  
270   int
271 + printgeom(              /* print out detailed geometry for BSDF */
272 +        struct BSDF_data *sd,
273 +        char *xfrot,
274 +        FVECT ctr,
275 +        double s1,
276 +        double s2
277 + )
278 + {
279 +        static char     mgftemp[] = TEMPLATE;
280 +        char            cmdbuf[64];
281 +        FILE            *fp;
282 +        double          sca;
283 +
284 +        if (sd == NULL || sd->mgf == NULL)
285 +                return(0);
286 +        if (sd->dim[0] <= FTINY || sd->dim[1] <= FTINY)
287 +                return(0);
288 +        if ((s1 > s2) ^ (sd->dim[0] > sd->dim[1])) {
289 +                sca = s1; s1 = s2; s2 = sca;
290 +        }
291 +        s1 /= sd->dim[0];
292 +        s2 /= sd->dim[1];
293 +        sca = s1 > s2 ? s1 : s2;
294 +        strcpy(mgftemp, TEMPLATE);
295 +        if ((fp = fopen(mktemp(mgftemp), "w")) == NULL)
296 +                error(SYSTEM, "cannot create temporary file for MGF");
297 +                                        /* prepend our transform */
298 +        fprintf(fp, "xf%s -s %.5f -t %.5g %.5g %.5g\n",
299 +                        xfrot, sca, ctr[0], ctr[1], ctr[2]);
300 +                                        /* output given MGF description */
301 +        fputs(sd->mgf, fp);
302 +        fputs("\nxf\n", fp);
303 +        if (fclose(fp) == EOF)
304 +                error(SYSTEM, "error writing MGF temporary file");
305 +                                        /* execute mgf2rad to convert MGF */
306 +        strcpy(cmdbuf, "mgf2rad ");
307 +        strcpy(cmdbuf+8, mgftemp);
308 +        fflush(stdout);
309 +        system(cmdbuf);
310 +        unlink(mgftemp);                /* clean up */
311 +        return(1);
312 + }
313 +
314 +
315 + int
316   my_default(     /* default illum action */
317          OBJREC  *ob,
318          struct illum_args  *il,
# Line 293 | Line 342 | my_face(               /* make an illum face */
342          FVECT  u, v;
343          double  ur[2], vr[2];
344          MAT4  xfm;
345 +        char  xfrot[64];
346          int  nallow;
347          FACE  *fa;
348          int  i, j;
# Line 304 | Line 354 | my_face(               /* make an illum face */
354          }
355                                  /* set up sampling */
356          if (il->sd != NULL) {
357 <                if (!getBSDF_xfm(xfm, fa->norm, il->udir)) {
357 >                if (!getBSDF_xfm(xfm, fa->norm, il->udir, xfrot)) {
358                          objerror(ob, WARNING, "illegal up direction");
359                          freeface(ob);
360                          return(my_default(ob, il, nm));
# Line 346 | Line 396 | my_face(               /* make an illum face */
396                  if (r2 < vr[0]) vr[0] = r2;
397                  if (r2 > vr[1]) vr[1] = r2;
398          }
399 +                                /* output detailed geometry? */
400 +        if (!(il->flags & IL_LIGHT) && il->sd != NULL && il->sd->mgf != NULL &&
401 +                        il->thick <= FTINY) {
402 +                for (j = 3; j--; )
403 +                        org[j] = .5*(ur[0]+ur[1])*u[j] +
404 +                                        .5*(vr[0]+vr[1])*v[j] +
405 +                                        fa->offset*fa->norm[j];
406 +                printgeom(il->sd, xfrot, org, ur[1]-ur[0], vr[1]-vr[0]);
407 +        }
408          dim[0] = random();
409                                  /* sample polygon */
410          nallow = 5*n*il->nsamps;
# Line 390 | Line 449 | my_face(               /* make an illum face */
449                      raysamp(dim[1], org, dir);
450                  }
451                                  /* add in direct component? */
452 <        if (!directvis && (il->flags & IL_LIGHT || il->sd != NULL)) {
452 >        if (il->flags & IL_LIGHT || il->sd != NULL) {
453                  MAT4    ixfm;
454                  if (il->sd == NULL) {
455                          for (i = 3; i--; ) {
# Line 544 | Line 603 | my_ring(               /* make an illum ring */
603          co = getcone(ob, 0);
604                                  /* set up sampling */
605          if (il->sd != NULL) {
606 <                if (!getBSDF_xfm(xfm, co->ad, il->udir)) {
606 >                if (!getBSDF_xfm(xfm, co->ad, il->udir, NULL)) {
607                          objerror(ob, WARNING, "illegal up direction");
608                          freecone(ob);
609                          return(my_default(ob, il, nm));
# Line 598 | Line 657 | my_ring(               /* make an illum ring */
657                      raysamp(dim[1], org, dir);
658                  }
659                                  /* add in direct component? */
660 <        if (!directvis && (il->flags & IL_LIGHT || il->sd != NULL)) {
660 >        if (il->flags & IL_LIGHT || il->sd != NULL) {
661                  MAT4    ixfm;
662                  if (il->sd == NULL) {
663                          for (i = 3; i--; ) {
# Line 652 | Line 711 | my_ring(               /* make an illum ring */
711                                  /* clean up */
712          freecone(ob);
713          return(1);
714 + }
715 +
716 +
717 + void
718 + redistribute(           /* pass distarr ray sums through BSDF */
719 +        struct BSDF_data *b,
720 +        int nalt,
721 +        int nazi,
722 +        FVECT u,
723 +        FVECT v,
724 +        FVECT w,
725 +        MAT4 xm
726 + )
727 + {
728 +        int     nout = 0;
729 +        MAT4    mymat, inmat;
730 +        COLORV  *idist;
731 +        COLORV  *cp;
732 +        FVECT   dv;
733 +        double  wt;
734 +        int     i, j, k, c, o;
735 +        COLOR   col, cinc;
736 +                                        /* copy incoming distribution */
737 +        if (b->ninc > distsiz)
738 +                error(INTERNAL, "error 1 in redistribute");
739 +        idist = (COLORV *)malloc(sizeof(COLOR)*b->ninc);
740 +        if (idist == NULL)
741 +                error(SYSTEM, "out of memory in redistribute");
742 +        memcpy(idist, distarr, sizeof(COLOR)*b->ninc);
743 +                                        /* compose direction transform */
744 +        for (i = 3; i--; ) {
745 +                mymat[i][0] = u[i];
746 +                mymat[i][1] = v[i];
747 +                mymat[i][2] = w[i];
748 +                mymat[i][3] = 0.;
749 +        }
750 +        mymat[3][0] = mymat[3][1] = mymat[3][2] = 0.;
751 +        mymat[3][3] = 1.;
752 +        if (xm != NULL)
753 +                multmat4(mymat, xm, mymat);
754 +        for (i = 3; i--; ) {            /* make sure it's normalized */
755 +                wt = 1./sqrt(   mymat[0][i]*mymat[0][i] +
756 +                                mymat[1][i]*mymat[1][i] +
757 +                                mymat[2][i]*mymat[2][i] );
758 +                for (j = 3; j--; )
759 +                        mymat[j][i] *= wt;
760 +        }
761 +        if (!invmat4(inmat, mymat))     /* need inverse as well */
762 +                error(INTERNAL, "cannot invert BSDF transform");
763 +        newdist(nalt*nazi);             /* resample distribution */
764 +        for (i = b->ninc; i--; ) {
765 +                int     direct_out = -1;
766 +                COLOR   cdir;
767 +                getBSDF_incvec(dv, b, i);       /* compute incident irrad. */
768 +                multv3(dv, dv, mymat);
769 +                if (dv[2] < 0.0) {
770 +                        dv[0] = -dv[0]; dv[1] = -dv[1]; dv[2] = -dv[2];
771 +                        direct_out += (direct_discount != NULL);
772 +                }
773 +                wt = getBSDF_incohm(b, i);
774 +                wt *= dv[2];                    /* solid_angle*cosine(theta) */
775 +                cp = &idist[3*i];
776 +                copycolor(cinc, cp);
777 +                scalecolor(cinc, wt);
778 +                if (!direct_out) {              /* discount direct contr. */
779 +                        cp = &direct_discount[3*i];
780 +                        copycolor(cdir, cp);
781 +                        scalecolor(cdir, -wt);
782 +                        if (b->nout != b->ninc)
783 +                                direct_out = flatindex(dv, nalt, nazi);
784 +                        else
785 +                                direct_out = i; /* assumes dist. mirroring */
786 +                }
787 +                for (k = nalt; k--; )           /* loop over distribution */
788 +                  for (j = nazi; j--; ) {
789 +                    int rstart = random();
790 +                    for (c = NBSDFSAMPS; c--; ) {
791 +                        double  sp[2];
792 +                        multisamp(sp, 2, urand(rstart+c));
793 +                        flatdir(dv, (k + sp[0])/nalt,
794 +                                        (j + .5 - sp[1])/nazi);
795 +                        multv3(dv, dv, inmat);
796 +                                                /* evaluate BSDF @ outgoing */
797 +                        o = getBSDF_outndx(b, dv);
798 +                        if (o < 0) {
799 +                                nout++;
800 +                                continue;
801 +                        }
802 +                        wt = BSDF_value(b, i, o) * (1./NBSDFSAMPS);
803 +                        copycolor(col, cinc);
804 +                        if (b->nout != b->ninc)
805 +                                o = k*nazi + j;
806 +                        if (o == direct_out)
807 +                                addcolor(col, cdir);    /* minus direct */
808 +                        scalecolor(col, wt);
809 +                        cp = &distarr[3*(k*nazi + j)];
810 +                        addcolor(cp, col);      /* sum into distribution */
811 +                    }
812 +                  }
813 +        }
814 +        free(idist);                    /* free temp space */
815 +        if (nout) {
816 +                sprintf(errmsg, "missing %.1f%% of BSDF directions",
817 +                                100.*nout/(b->ninc*nalt*nazi*NBSDFSAMPS));
818 +                error(WARNING, errmsg);
819 +        }
820   }

Diff Legend

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