Back to index

radiance  4R0+20100331
normal.c
Go to the documentation of this file.
00001 #ifndef lint
00002 static const char RCSid[] = "$Id: normal.c,v 2.51 2005/12/08 17:56:38 greg Exp $";
00003 #endif
00004 /*
00005  *  normal.c - shading function for normal materials.
00006  *
00007  *     8/19/85
00008  *     12/19/85 - added stuff for metals.
00009  *     6/26/87 - improved specular model.
00010  *     9/28/87 - added model for translucent materials.
00011  *     Later changes described in delta comments.
00012  */
00013 
00014 #include "copyright.h"
00015 
00016 #include  "ray.h"
00017 #include  "ambient.h"
00018 #include  "source.h"
00019 #include  "otypes.h"
00020 #include  "rtotypes.h"
00021 #include  "random.h"
00022 
00023 #ifndef  MAXITER
00024 #define  MAXITER     10            /* maximum # specular ray attempts */
00025 #endif
00026                                    /* estimate of Fresnel function */
00027 #define  FRESNE(ci)  (exp(-5.85*(ci)) - 0.00287989916)
00028 #define  FRESTHRESH  0.017999      /* minimum specularity for approx. */
00029 
00030 
00031 /*
00032  *     This routine implements the isotropic Gaussian
00033  *  model described by Ward in Siggraph `92 article.
00034  *     We orient the surface towards the incoming ray, so a single
00035  *  surface can be used to represent an infinitely thin object.
00036  *
00037  *  Arguments for MAT_PLASTIC and MAT_METAL are:
00038  *     red    grn    blu    specular-frac.       facet-slope
00039  *
00040  *  Arguments for MAT_TRANS are:
00041  *     red    grn    blu    rspec  rough  trans  tspec
00042  */
00043 
00044                             /* specularity flags */
00045 #define  SP_REFL     01            /* has reflected specular component */
00046 #define  SP_TRAN     02            /* has transmitted specular */
00047 #define  SP_PURE     04            /* purely specular (zero roughness) */
00048 #define  SP_FLAT     010           /* flat reflecting surface */
00049 #define  SP_RBLT     020           /* reflection below sample threshold */
00050 #define  SP_TBLT     040           /* transmission below threshold */
00051 
00052 typedef struct {
00053        OBJREC  *mp;         /* material pointer */
00054        RAY  *rp;            /* ray pointer */
00055        short  specfl;              /* specularity flags, defined above */
00056        COLOR  mcolor;              /* color of this material */
00057        COLOR  scolor;              /* color of specular component */
00058        FVECT  vrefl;        /* vector in direction of reflected ray */
00059        FVECT  prdir;        /* vector in transmitted direction */
00060        double  alpha2;             /* roughness squared */
00061        double  rdiff, rspec;       /* reflected specular, diffuse */
00062        double  trans;              /* transmissivity */
00063        double  tdiff, tspec;       /* transmitted specular, diffuse */
00064        FVECT  pnorm;        /* perturbed surface normal */
00065        double  pdot;        /* perturbed dot product */
00066 }  NORMDAT;          /* normal material data */
00067 
00068 static srcdirf_t dirnorm;
00069 static void gaussamp(RAY  *r, NORMDAT  *np);
00070 
00071 
00072 static void
00073 dirnorm(             /* compute source contribution */
00074        COLOR  cval,                /* returned coefficient */
00075        void  *nnp,          /* material data */
00076        FVECT  ldir,                /* light source direction */
00077        double  omega               /* light source size */
00078 )
00079 {
00080        register NORMDAT *np = nnp;
00081        double  ldot;
00082        double  lrdiff, ltdiff;
00083        double  dtmp, d2;
00084        FVECT  vtmp;
00085        COLOR  ctmp;
00086 
00087        setcolor(cval, 0.0, 0.0, 0.0);
00088 
00089        ldot = DOT(np->pnorm, ldir);
00090 
00091        if (ldot < 0.0 ? np->trans <= FTINY : np->trans >= 1.0-FTINY)
00092               return;              /* wrong side */
00093 
00094                             /* Fresnel estimate */
00095        lrdiff = np->rdiff;
00096        ltdiff = np->tdiff;
00097        if (np->specfl & SP_PURE && np->rspec >= FRESTHRESH &&
00098                      (lrdiff > FTINY) | (ltdiff > FTINY)) {
00099               dtmp = 1. - FRESNE(fabs(ldot));
00100               lrdiff *= dtmp;
00101               ltdiff *= dtmp;
00102        }
00103 
00104        if (ldot > FTINY && lrdiff > FTINY) {
00105               /*
00106                *  Compute and add diffuse reflected component to returned
00107                *  color.  The diffuse reflected component will always be
00108                *  modified by the color of the material.
00109                */
00110               copycolor(ctmp, np->mcolor);
00111               dtmp = ldot * omega * lrdiff * (1.0/PI);
00112               scalecolor(ctmp, dtmp);
00113               addcolor(cval, ctmp);
00114        }
00115        if (ldot > FTINY && (np->specfl&(SP_REFL|SP_PURE)) == SP_REFL) {
00116               /*
00117                *  Compute specular reflection coefficient using
00118                *  gaussian distribution model.
00119                */
00120                                           /* roughness */
00121               dtmp = np->alpha2;
00122                                           /* + source if flat */
00123               if (np->specfl & SP_FLAT)
00124                      dtmp += omega * (0.25/PI);
00125                                           /* half vector */
00126               vtmp[0] = ldir[0] - np->rp->rdir[0];
00127               vtmp[1] = ldir[1] - np->rp->rdir[1];
00128               vtmp[2] = ldir[2] - np->rp->rdir[2];
00129               d2 = DOT(vtmp, np->pnorm);
00130               d2 *= d2;
00131               d2 = (DOT(vtmp,vtmp) - d2) / d2;
00132                                           /* gaussian */
00133               dtmp = exp(-d2/dtmp)/(4.*PI * np->pdot * dtmp);
00134                                           /* worth using? */
00135               if (dtmp > FTINY) {
00136                      copycolor(ctmp, np->scolor);
00137                      dtmp *= omega;
00138                      scalecolor(ctmp, dtmp);
00139                      addcolor(cval, ctmp);
00140               }
00141        }
00142        if (ldot < -FTINY && ltdiff > FTINY) {
00143               /*
00144                *  Compute diffuse transmission.
00145                */
00146               copycolor(ctmp, np->mcolor);
00147               dtmp = -ldot * omega * ltdiff * (1.0/PI);
00148               scalecolor(ctmp, dtmp);
00149               addcolor(cval, ctmp);
00150        }
00151        if (ldot < -FTINY && (np->specfl&(SP_TRAN|SP_PURE)) == SP_TRAN) {
00152               /*
00153                *  Compute specular transmission.  Specular transmission
00154                *  is always modified by material color.
00155                */
00156                                           /* roughness + source */
00157               dtmp = np->alpha2 + omega*(1.0/PI);
00158                                           /* gaussian */
00159               dtmp = exp((2.*DOT(np->prdir,ldir)-2.)/dtmp) /
00160                                    (PI*np->pdot*dtmp);
00161                                           /* worth using? */
00162               if (dtmp > FTINY) {
00163                      copycolor(ctmp, np->mcolor);
00164                      dtmp *= np->tspec * omega;
00165                      scalecolor(ctmp, dtmp);
00166                      addcolor(cval, ctmp);
00167               }
00168        }
00169 }
00170 
00171 
00172 extern int
00173 m_normal(                   /* color a ray that hit something normal */
00174        register OBJREC  *m,
00175        register RAY  *r
00176 )
00177 {
00178        NORMDAT  nd;
00179        double  fest;
00180        double  transtest, transdist;
00181        double mirtest, mirdist;
00182        int    hastexture;
00183        double d;
00184        COLOR  ctmp;
00185        register int  i;
00186                                           /* easy shadow test */
00187        if (r->crtype & SHADOW && m->otype != MAT_TRANS)
00188               return(1);
00189 
00190        if (m->oargs.nfargs != (m->otype == MAT_TRANS ? 7 : 5))
00191               objerror(m, USER, "bad number of arguments");
00192                                           /* check for back side */
00193        if (r->rod < 0.0) {
00194               if (!backvis && m->otype != MAT_TRANS) {
00195                      raytrans(r);
00196                      return(1);
00197               }
00198               raytexture(r, m->omod);
00199               flipsurface(r);                    /* reorient if backvis */
00200        } else
00201               raytexture(r, m->omod);
00202        nd.mp = m;
00203        nd.rp = r;
00204                                           /* get material color */
00205        setcolor(nd.mcolor, m->oargs.farg[0],
00206                         m->oargs.farg[1],
00207                         m->oargs.farg[2]);
00208                                           /* get roughness */
00209        nd.specfl = 0;
00210        nd.alpha2 = m->oargs.farg[4];
00211        if ((nd.alpha2 *= nd.alpha2) <= FTINY)
00212               nd.specfl |= SP_PURE;
00213 
00214        if ( (hastexture = (DOT(r->pert,r->pert) > FTINY*FTINY)) ) {
00215               nd.pdot = raynormal(nd.pnorm, r);  /* perturb normal */
00216        } else {
00217               VCOPY(nd.pnorm, r->ron);
00218               nd.pdot = r->rod;
00219        }
00220        if (r->ro != NULL && isflat(r->ro->otype))
00221               nd.specfl |= SP_FLAT;
00222        if (nd.pdot < .001)
00223               nd.pdot = .001;                    /* non-zero for dirnorm() */
00224        multcolor(nd.mcolor, r->pcol);            /* modify material color */
00225        mirtest = transtest = 0;
00226        mirdist = transdist = r->rot;
00227        nd.rspec = m->oargs.farg[3];
00228                                           /* compute Fresnel approx. */
00229        if (nd.specfl & SP_PURE && nd.rspec >= FRESTHRESH) {
00230               fest = FRESNE(r->rod);
00231               nd.rspec += fest*(1. - nd.rspec);
00232        } else
00233               fest = 0.;
00234                                           /* compute transmission */
00235        if (m->otype == MAT_TRANS) {
00236               nd.trans = m->oargs.farg[5]*(1.0 - nd.rspec);
00237               nd.tspec = nd.trans * m->oargs.farg[6];
00238               nd.tdiff = nd.trans - nd.tspec;
00239               if (nd.tspec > FTINY) {
00240                      nd.specfl |= SP_TRAN;
00241                                                  /* check threshold */
00242                      if (!(nd.specfl & SP_PURE) &&
00243                                    specthresh >= nd.tspec-FTINY)
00244                             nd.specfl |= SP_TBLT;
00245                      if (!hastexture || r->crtype & SHADOW) {
00246                             VCOPY(nd.prdir, r->rdir);
00247                             transtest = 2;
00248                      } else {
00249                             for (i = 0; i < 3; i++)            /* perturb */
00250                                    nd.prdir[i] = r->rdir[i] - r->pert[i];
00251                             if (DOT(nd.prdir, r->ron) < -FTINY)
00252                                    normalize(nd.prdir); /* OK */
00253                             else
00254                                    VCOPY(nd.prdir, r->rdir);
00255                      }
00256               }
00257        } else
00258               nd.tdiff = nd.tspec = nd.trans = 0.0;
00259                                           /* transmitted ray */
00260        if ((nd.specfl&(SP_TRAN|SP_PURE|SP_TBLT)) == (SP_TRAN|SP_PURE)) {
00261               RAY  lr;
00262               copycolor(lr.rcoef, nd.mcolor);    /* modified by color */
00263               scalecolor(lr.rcoef, nd.tspec);
00264               if (rayorigin(&lr, TRANS, r, lr.rcoef) == 0) {
00265                      VCOPY(lr.rdir, nd.prdir);
00266                      rayvalue(&lr);
00267                      multcolor(lr.rcol, lr.rcoef);
00268                      addcolor(r->rcol, lr.rcol);
00269                      transtest *= bright(lr.rcol);
00270                      transdist = r->rot + lr.rt;
00271               }
00272        } else
00273               transtest = 0;
00274 
00275        if (r->crtype & SHADOW) {          /* the rest is shadow */
00276               r->rt = transdist;
00277               return(1);
00278        }
00279                                           /* get specular reflection */
00280        if (nd.rspec > FTINY) {
00281               nd.specfl |= SP_REFL;
00282                                           /* compute specular color */
00283               if (m->otype != MAT_METAL) {
00284                      setcolor(nd.scolor, nd.rspec, nd.rspec, nd.rspec);
00285               } else if (fest > FTINY) {
00286                      d = nd.rspec*(1. - fest);
00287                      for (i = 0; i < 3; i++)
00288                             nd.scolor[i] = fest + nd.mcolor[i]*d;
00289               } else {
00290                      copycolor(nd.scolor, nd.mcolor);
00291                      scalecolor(nd.scolor, nd.rspec);
00292               }
00293                                           /* check threshold */
00294               if (!(nd.specfl & SP_PURE) && specthresh >= nd.rspec-FTINY)
00295                      nd.specfl |= SP_RBLT;
00296                                           /* compute reflected ray */
00297               for (i = 0; i < 3; i++)
00298                      nd.vrefl[i] = r->rdir[i] + 2.*nd.pdot*nd.pnorm[i];
00299                                           /* penetration? */
00300               if (hastexture && DOT(nd.vrefl, r->ron) <= FTINY)
00301                      for (i = 0; i < 3; i++)            /* safety measure */
00302                             nd.vrefl[i] = r->rdir[i] + 2.*r->rod*r->ron[i];
00303        }
00304                                           /* reflected ray */
00305        if ((nd.specfl&(SP_REFL|SP_PURE|SP_RBLT)) == (SP_REFL|SP_PURE)) {
00306               RAY  lr;
00307               if (rayorigin(&lr, REFLECTED, r, nd.scolor) == 0) {
00308                      VCOPY(lr.rdir, nd.vrefl);
00309                      rayvalue(&lr);
00310                      multcolor(lr.rcol, lr.rcoef);
00311                      addcolor(r->rcol, lr.rcol);
00312                      if (!hastexture && nd.specfl & SP_FLAT) {
00313                             mirtest = 2.*bright(lr.rcol);
00314                             mirdist = r->rot + lr.rt;
00315                      }
00316               }
00317        }
00318                                           /* diffuse reflection */
00319        nd.rdiff = 1.0 - nd.trans - nd.rspec;
00320 
00321        if (nd.specfl & SP_PURE && nd.rdiff <= FTINY && nd.tdiff <= FTINY)
00322               return(1);                  /* 100% pure specular */
00323 
00324        if (!(nd.specfl & SP_PURE))
00325               gaussamp(r, &nd);           /* checks *BLT flags */
00326 
00327        if (nd.rdiff > FTINY) {            /* ambient from this side */
00328               copycolor(ctmp, nd.mcolor); /* modified by material color */
00329               if (nd.specfl & SP_RBLT)
00330                      scalecolor(ctmp, 1.0-nd.trans);
00331               else
00332                      scalecolor(ctmp, nd.rdiff);
00333               multambient(ctmp, r, hastexture ? nd.pnorm : r->ron);
00334               addcolor(r->rcol, ctmp);    /* add to returned color */
00335        }
00336        if (nd.tdiff > FTINY) {            /* ambient from other side */
00337               copycolor(ctmp, nd.mcolor); /* modified by color */
00338               if (nd.specfl & SP_TBLT)
00339                      scalecolor(ctmp, nd.trans);
00340               else
00341                      scalecolor(ctmp, nd.tdiff);
00342               flipsurface(r);
00343               if (hastexture) {
00344                      FVECT  bnorm;
00345                      bnorm[0] = -nd.pnorm[0];
00346                      bnorm[1] = -nd.pnorm[1];
00347                      bnorm[2] = -nd.pnorm[2];
00348                      multambient(ctmp, r, bnorm);
00349               } else
00350                      multambient(ctmp, r, r->ron);
00351               addcolor(r->rcol, ctmp);
00352               flipsurface(r);
00353        }
00354                                    /* add direct component */
00355        direct(r, dirnorm, &nd);
00356                                    /* check distance */
00357        d = bright(r->rcol);
00358        if (transtest > d)
00359               r->rt = transdist;
00360        else if (mirtest > d)
00361               r->rt = mirdist;
00362 
00363        return(1);
00364 }
00365 
00366 
00367 static void
00368 gaussamp(                   /* sample gaussian specular */
00369        RAY  *r,
00370        register NORMDAT  *np
00371 )
00372 {
00373        RAY  sr;
00374        FVECT  u, v, h;
00375        double  rv[2];
00376        double  d, sinp, cosp;
00377        int  niter;
00378        register int  i;
00379                                    /* quick test */
00380        if ((np->specfl & (SP_REFL|SP_RBLT)) != SP_REFL &&
00381                      (np->specfl & (SP_TRAN|SP_TBLT)) != SP_TRAN)
00382               return;
00383                                    /* set up sample coordinates */
00384        v[0] = v[1] = v[2] = 0.0;
00385        for (i = 0; i < 3; i++)
00386               if (np->pnorm[i] < 0.6 && np->pnorm[i] > -0.6)
00387                      break;
00388        v[i] = 1.0;
00389        fcross(u, v, np->pnorm);
00390        normalize(u);
00391        fcross(v, np->pnorm, u);
00392                                    /* compute reflection */
00393        if ((np->specfl & (SP_REFL|SP_RBLT)) == SP_REFL &&
00394                      rayorigin(&sr, SPECULAR, r, np->scolor) == 0) {
00395               dimlist[ndims++] = (int)np->mp;
00396               for (niter = 0; niter < MAXITER; niter++) {
00397                      if (niter)
00398                             d = frandom();
00399                      else
00400                             d = urand(ilhash(dimlist,ndims)+samplendx);
00401                      multisamp(rv, 2, d);
00402                      d = 2.0*PI * rv[0];
00403                      cosp = tcos(d);
00404                      sinp = tsin(d);
00405                      rv[1] = 1.0 - specjitter*rv[1];
00406                      if (rv[1] <= FTINY)
00407                             d = 1.0;
00408                      else
00409                             d = sqrt( np->alpha2 * -log(rv[1]) );
00410                      for (i = 0; i < 3; i++)
00411                             h[i] = np->pnorm[i] + d*(cosp*u[i] + sinp*v[i]);
00412                      d = -2.0 * DOT(h, r->rdir) / (1.0 + d*d);
00413                      for (i = 0; i < 3; i++)
00414                             sr.rdir[i] = r->rdir[i] + d*h[i];
00415                      if (DOT(sr.rdir, r->ron) > FTINY) {
00416                             rayvalue(&sr);
00417                             multcolor(sr.rcol, sr.rcoef);
00418                             addcolor(r->rcol, sr.rcol);
00419                             break;
00420                      }
00421               }
00422               ndims--;
00423        }
00424                                    /* compute transmission */
00425        copycolor(sr.rcoef, np->mcolor);   /* modified by color */
00426        scalecolor(sr.rcoef, np->tspec);
00427        if ((np->specfl & (SP_TRAN|SP_TBLT)) == SP_TRAN &&
00428                      rayorigin(&sr, SPECULAR, r, sr.rcoef) == 0) {
00429               dimlist[ndims++] = (int)np->mp;
00430               for (niter = 0; niter < MAXITER; niter++) {
00431                      if (niter)
00432                             d = frandom();
00433                      else
00434                             d = urand(ilhash(dimlist,ndims)+1823+samplendx);
00435                      multisamp(rv, 2, d);
00436                      d = 2.0*PI * rv[0];
00437                      cosp = tcos(d);
00438                      sinp = tsin(d);
00439                      rv[1] = 1.0 - specjitter*rv[1];
00440                      if (rv[1] <= FTINY)
00441                             d = 1.0;
00442                      else
00443                             d = sqrt( np->alpha2 * -log(rv[1]) );
00444                      for (i = 0; i < 3; i++)
00445                             sr.rdir[i] = np->prdir[i] + d*(cosp*u[i] + sinp*v[i]);
00446                      if (DOT(sr.rdir, r->ron) < -FTINY) {
00447                             normalize(sr.rdir);  /* OK, normalize */
00448                             rayvalue(&sr);
00449                             multcolor(sr.rcol, sr.rcoef);
00450                             addcolor(r->rcol, sr.rcol);
00451                             break;
00452                      }
00453               }
00454               ndims--;
00455        }
00456 }