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

Comparing ray/src/rt/dielectric.c (file contents):
Revision 1.3 by greg, Tue Mar 27 11:39:58 1990 UTC vs.
Revision 2.23 by greg, Wed Aug 7 05:10:09 2013 UTC

# Line 1 | Line 1
1 – /* Copyright (c) 1986 Regents of the University of California */
2 –
1   #ifndef lint
2 < static char SCCSid[] = "$SunId$ LBL";
2 > static const char       RCSid[] = "$Id$";
3   #endif
6 –
4   /*
5   *  dielectric.c - shading function for transparent materials.
9 – *
10 – *     9/6/85
6   */
7  
8 < #include  "ray.h"
8 > #include "copyright.h"
9  
10 + #include  "ray.h"
11   #include  "otypes.h"
12 + #include  "rtotypes.h"
13  
14   #ifdef  DISPERSE
15   #include  "source.h"
16 + static int disperse(OBJREC *m,RAY *r,FVECT vt,double tr,COLOR cet,COLOR abt);
17 + static int lambda(OBJREC  *m, FVECT  v2, FVECT  dv, FVECT  lr);
18   #endif
19  
20 + static double mylog(double  x);
21 +
22 +
23   /*
24   *     Explicit calculations for Fresnel's equation are performed,
25   *  but only one square root computation is necessary.
# Line 49 | Line 51 | static char SCCSid[] = "$SunId$ LBL";
51   #define  MINCOS         0.997           /* minimum dot product for dispersion */
52  
53  
54 < m_dielectric(m, r)      /* color a ray which hit something transparent */
55 < OBJREC  *m;
56 < register RAY  *r;
54 > static double
55 > mylog(          /* special log for extinction coefficients */
56 >        double  x
57 > )
58   {
59 <        double  sqrt(), pow();
59 >        if (x < 1e-40)
60 >                return(-100.);
61 >        if (x >= 1.)
62 >                return(0.);
63 >        return(log(x));
64 > }
65 >
66 >
67 > extern int
68 > m_dielectric(   /* color a ray which hit a dielectric interface */
69 >        OBJREC  *m,
70 >        register RAY  *r
71 > )
72 > {
73          double  cos1, cos2, nratio;
74 <        COLOR  mcolor;
75 <        double  mabsorp;
74 >        COLOR  ctrans;
75 >        COLOR  talb;
76 >        int  hastexture;
77 >        double  transdist=0, transtest=0;
78 >        double  mirdist=0, mirtest=0;
79 >        int     flatsurface;
80          double  refl, trans;
81          FVECT  dnorm;
82          double  d1, d2;
# Line 68 | Line 88 | register RAY  *r;
88  
89          raytexture(r, m->omod);                 /* get modifiers */
90  
91 <        cos1 = raynormal(dnorm, r);             /* cosine of theta1 */
91 >        if ( (hastexture = DOT(r->pert,r->pert) > FTINY*FTINY) )
92 >                cos1 = raynormal(dnorm, r);     /* perturb normal */
93 >        else {
94 >                VCOPY(dnorm, r->ron);
95 >                cos1 = r->rod;
96 >        }
97 >        flatsurface = !hastexture && r->ro != NULL && isflat(r->ro->otype);
98 >
99                                                  /* index of refraction */
100          if (m->otype == MAT_DIELECTRIC)
101                  nratio = m->oargs.farg[3] + m->oargs.farg[4]/MLAMBDA;
# Line 76 | Line 103 | register RAY  *r;
103                  nratio = m->oargs.farg[3] / m->oargs.farg[7];
104          
105          if (cos1 < 0.0) {                       /* inside */
106 +                hastexture = -hastexture;
107                  cos1 = -cos1;
108                  dnorm[0] = -dnorm[0];
109                  dnorm[1] = -dnorm[1];
110                  dnorm[2] = -dnorm[2];
111 <                setcolor(mcolor, pow(m->oargs.farg[0], r->rot),
112 <                                 pow(m->oargs.farg[1], r->rot),
113 <                                 pow(m->oargs.farg[2], r->rot));
111 >                setcolor(r->cext, -mylog(m->oargs.farg[0]*colval(r->pcol,RED)),
112 >                                 -mylog(m->oargs.farg[1]*colval(r->pcol,GRN)),
113 >                                 -mylog(m->oargs.farg[2]*colval(r->pcol,BLU)));
114 >                setcolor(r->albedo, 0., 0., 0.);
115 >                r->gecc = 0.;
116 >                if (m->otype == MAT_INTERFACE) {
117 >                        setcolor(ctrans,
118 >                                -mylog(m->oargs.farg[4]*colval(r->pcol,RED)),
119 >                                -mylog(m->oargs.farg[5]*colval(r->pcol,GRN)),
120 >                                -mylog(m->oargs.farg[6]*colval(r->pcol,BLU)));
121 >                        setcolor(talb, 0., 0., 0.);
122 >                } else {
123 >                        copycolor(ctrans, cextinction);
124 >                        copycolor(talb, salbedo);
125 >                }
126          } else {                                /* outside */
127                  nratio = 1.0 / nratio;
128 <                if (m->otype == MAT_INTERFACE)
129 <                        setcolor(mcolor, pow(m->oargs.farg[4], r->rot),
130 <                                         pow(m->oargs.farg[5], r->rot),
131 <                                         pow(m->oargs.farg[6], r->rot));
132 <                else
133 <                        setcolor(mcolor, 1.0, 1.0, 1.0);
128 >
129 >                setcolor(ctrans, -mylog(m->oargs.farg[0]*colval(r->pcol,RED)),
130 >                                 -mylog(m->oargs.farg[1]*colval(r->pcol,GRN)),
131 >                                 -mylog(m->oargs.farg[2]*colval(r->pcol,BLU)));
132 >                setcolor(talb, 0., 0., 0.);
133 >                if (m->otype == MAT_INTERFACE) {
134 >                        setcolor(r->cext,
135 >                                -mylog(m->oargs.farg[4]*colval(r->pcol,RED)),
136 >                                -mylog(m->oargs.farg[5]*colval(r->pcol,GRN)),
137 >                                -mylog(m->oargs.farg[6]*colval(r->pcol,BLU)));
138 >                        setcolor(r->albedo, 0., 0., 0.);
139 >                        r->gecc = 0.;
140 >                }
141          }
95 –        mabsorp = bright(mcolor);
142  
143          d2 = 1.0 - nratio*nratio*(1.0 - cos1*cos1);     /* compute cos theta2 */
144  
# Line 113 | Line 159 | register RAY  *r;
159                  d1 = (d1 - d2) / (d1 + d2);
160                  refl += d1 * d1;
161  
162 <                refl /= 2.0;
162 >                refl *= 0.5;
163                  trans = 1.0 - refl;
164  
165 <                if (rayorigin(&p, r, REFRACTED, mabsorp*trans) == 0) {
165 >                trans *= nratio*nratio;         /* solid angle ratio */
166  
167 +                setcolor(p.rcoef, trans, trans, trans);
168 +
169 +                if (rayorigin(&p, REFRACTED, r, p.rcoef) == 0) {
170 +
171                                                  /* compute refracted ray */
172                          d1 = nratio*cos1 - cos2;
173                          for (i = 0; i < 3; i++)
174                                  p.rdir[i] = nratio*r->rdir[i] + d1*dnorm[i];
175 <
175 >                                                /* accidental reflection? */
176 >                        if (hastexture &&
177 >                                DOT(p.rdir,r->ron)*hastexture >= -FTINY) {
178 >                                d1 *= (double)hastexture;
179 >                                for (i = 0; i < 3; i++) /* ignore texture */
180 >                                        p.rdir[i] = nratio*r->rdir[i] +
181 >                                                        d1*r->ron[i];
182 >                                normalize(p.rdir);      /* not exact */
183 >                        } else
184 >                                checknorm(p.rdir);
185   #ifdef  DISPERSE
186                          if (m->otype != MAT_DIELECTRIC
187                                          || r->rod > 0.0
188                                          || r->crtype & SHADOW
189 +                                        || !directvis
190                                          || m->oargs.farg[4] == 0.0
191 <                                        || !disperse(m, r, p.rdir, trans))
191 >                                        || !disperse(m, r, p.rdir,
192 >                                                        trans, ctrans, talb))
193   #endif
194                          {
195 +                                copycolor(p.cext, ctrans);
196 +                                copycolor(p.albedo, talb);
197                                  rayvalue(&p);
198 <                                multcolor(mcolor, r->pcol);     /* modify */
136 <                                scalecolor(p.rcol, trans);
198 >                                multcolor(p.rcol, p.rcoef);
199                                  addcolor(r->rcol, p.rcol);
200 <                                r->rt = r->rot + p.rt;
200 >                                                /* virtual distance */
201 >                                if (flatsurface ||
202 >                                        (1.-FTINY <= nratio &&
203 >                                                nratio <= 1.+FTINY)) {
204 >                                        transtest = 2*bright(p.rcol);
205 >                                        transdist = r->rot + p.rt;
206 >                                }
207                          }
208                  }
209          }
210 <                
210 >        setcolor(p.rcoef, refl, refl, refl);
211 >
212          if (!(r->crtype & SHADOW) &&
213 <                        rayorigin(&p, r, REFLECTED, mabsorp*refl) == 0) {
213 >                        rayorigin(&p, REFLECTED, r, p.rcoef) == 0) {
214  
215                                          /* compute reflected ray */
216 <                for (i = 0; i < 3; i++)
217 <                        p.rdir[i] = r->rdir[i] + 2.0*cos1*dnorm[i];
218 <
216 >                VSUM(p.rdir, r->rdir, dnorm, 2.*cos1);
217 >                                        /* accidental penetration? */
218 >                if (hastexture && DOT(p.rdir,r->ron)*hastexture <= FTINY)
219 >                        VSUM(p.rdir, r->rdir, r->ron, 2.*r->rod);
220 >                checknorm(p.rdir);
221                  rayvalue(&p);                   /* reflected ray value */
222  
223 <                scalecolor(p.rcol, refl);       /* color contribution */
223 >                multcolor(p.rcol, p.rcoef);     /* color contribution */
224                  addcolor(r->rcol, p.rcol);
225 <                if (refl > trans)
226 <                        r->rt = r->rot + p.rt;
225 >                                                /* virtual distance */
226 >                if (flatsurface) {
227 >                        mirtest = 2*bright(p.rcol);
228 >                        mirdist = r->rot + p.rt;
229 >                }
230          }
231 <
232 <        multcolor(r->rcol, mcolor);             /* multiply by transmittance */
231 >                                /* check distance to return */
232 >        d1 = bright(r->rcol);
233 >        if (transtest > d1)
234 >                r->rt = transdist;
235 >        else if (mirtest > d1)
236 >                r->rt = mirdist;
237 >                                /* rayvalue() computes absorption */
238 >        return(1);
239   }
240  
241  
242   #ifdef  DISPERSE
243  
244 < static
245 < disperse(m, r, vt, tr)          /* check light sources for dispersion */
246 < OBJREC  *m;
247 < RAY  *r;
248 < FVECT  vt;
249 < double  tr;
244 > static int
245 > disperse(  /* check light sources for dispersion */
246 >        OBJREC  *m,
247 >        RAY  *r,
248 >        FVECT  vt,
249 >        double  tr,
250 >        COLOR  cet,
251 >        COLOR  abt
252 > )
253   {
254 <        double  sqrt();
255 <        RAY  sray, *entray;
254 >        RAY  sray;
255 >        const RAY  *entray;
256          FVECT  v1, v2, n1, n2;
257          FVECT  dv, v2Xdv;
258          double  v2Xdvv2Xdv;
259 <        int  sn, success = 0;
260 <        double  omega;
259 >        int  success = 0;
260 >        SRCINDEX  si;
261          FVECT  vtmp1, vtmp2;
262          double  dtmp1, dtmp2;
263          int  l1, l2;
# Line 224 | Line 307 | double  tr;
307          VCOPY(n2, r->ron);
308  
309                                          /* first order dispersion approx. */
310 <        dtmp1 = DOT(n1, v1);
311 <        dtmp2 = DOT(n2, v2);
310 >        dtmp1 = 1./DOT(n1, v1);
311 >        dtmp2 = 1./DOT(n2, v2);
312          for (i = 0; i < 3; i++)
313 <                dv[i] = v1[i] + v2[i] - n1[i]/dtmp1 - n2[i]/dtmp2;
313 >                dv[i] = v1[i] + v2[i] - n1[i]*dtmp1 - n2[i]*dtmp2;
314                  
315          if (DOT(dv, dv) <= FTINY)       /* null effect */
316                  return(0);
# Line 236 | Line 319 | double  tr;
319          v2Xdvv2Xdv = DOT(v2Xdv, v2Xdv);
320  
321                                          /* check sources */
322 <        for (sn = 0; sn < nsources; sn++) {
322 >        initsrcindex(&si);
323 >        while (srcray(&sray, r, &si)) {
324  
325 <                if ((omega = srcray(&sray, r, sn)) == 0.0 ||
242 <                                DOT(sray.rdir, v2) < MINCOS)
325 >                if (DOT(sray.rdir, v2) < MINCOS)
326                          continue;                       /* bad source */
244 –                
327                                                  /* adjust source ray */
328  
329                  dtmp1 = DOT(v2Xdv, sray.rdir) / v2Xdvv2Xdv;
# Line 254 | Line 336 | double  tr;
336                  if (l1 > MAXLAMBDA || l1 < MINLAMBDA)   /* not visible */
337                          continue;
338                                                  /* trace source ray */
339 +                copycolor(sray.cext, cet);
340 +                copycolor(sray.albedo, abt);
341                  normalize(sray.rdir);
342                  rayvalue(&sray);
343                  if (bright(sray.rcol) <= FTINY) /* missed it */
# Line 266 | Line 350 | double  tr;
350                   */
351                  
352                  fcross(vtmp1, v2Xdv, sray.rdir);
353 <                dtmp1 = sqrt(omega  / v2Xdvv2Xdv / PI);
353 >                dtmp1 = sqrt(si.dom  / v2Xdvv2Xdv / PI);
354  
355                                                          /* compute first ray */
356 <                for (i = 0; i < 3; i++)
273 <                        vtmp2[i] = sray.rdir[i] + dtmp1*vtmp1[i];
356 >                VSUM(vtmp2, sray.rdir, vtmp1, dtmp1);
357  
358                  l1 = lambda(m, v2, dv, vtmp2);          /* first lambda */
359                  if (l1 < 0)
360                          continue;
361                                                          /* compute second ray */
362 <                for (i = 0; i < 3; i++)
280 <                        vtmp2[i] = sray.rdir[i] - dtmp1*vtmp1[i];
362 >                VSUM(vtmp2, sray.rdir, vtmp1, -dtmp1);
363  
364                  l2 = lambda(m, v2, dv, vtmp2);          /* second lambda */
365                  if (l2 < 0)
# Line 297 | Line 379 | double  tr;
379  
380  
381   static int
382 < lambda(m, v2, dv, lr)                   /* compute lambda for material */
383 < register OBJREC  *m;
384 < FVECT  v2, dv, lr;
382 > lambda(                 /* compute lambda for material */
383 >        register OBJREC  *m,
384 >        FVECT  v2,
385 >        FVECT  dv,
386 >        FVECT  lr
387 > )
388   {
389          FVECT  lrXdv, v2Xlr;
390          double  dtmp, denom;
# Line 307 | Line 392 | FVECT  v2, dv, lr;
392  
393          fcross(lrXdv, lr, dv);
394          for (i = 0; i < 3; i++)
395 <                if (lrXdv[i] > FTINY || lrXdv[i] < -FTINY)
395 >                if ((lrXdv[i] > FTINY) | (lrXdv[i] < -FTINY))
396                          break;
397          if (i >= 3)
398                  return(-1);

Diff Legend

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