ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/dielectric.c
Revision: 2.18
Committed: Tue Mar 30 16:13:01 2004 UTC (20 years, 1 month ago) by schorsch
Content type: text/plain
Branch: MAIN
Changes since 2.17: +30 -19 lines
Log Message:
Continued ANSIfication. There are only bits and pieces left now.

File Contents

# User Rev Content
1 greg 1.1 #ifndef lint
2 schorsch 2.18 static const char RCSid[] = "$Id: dielectric.c,v 2.17 2003/07/27 22:12:03 schorsch Exp $";
3 greg 1.1 #endif
4     /*
5     * dielectric.c - shading function for transparent materials.
6 greg 2.15 */
7    
8 greg 2.16 #include "copyright.h"
9 greg 1.1
10     #include "ray.h"
11     #include "otypes.h"
12 schorsch 2.18 #include "rtotypes.h"
13 greg 1.1
14     #ifdef DISPERSE
15     #include "source.h"
16 schorsch 2.18 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 greg 1.1 #endif
19    
20 schorsch 2.18 static double mylog(double x);
21    
22    
23 greg 1.1 /*
24     * Explicit calculations for Fresnel's equation are performed,
25     * but only one square root computation is necessary.
26     * The index of refraction is given as a Hartmann equation
27     * with lambda0 equal to zero. If the slope of Hartmann's
28     * equation is non-zero, the material disperses light upon
29     * refraction. This condition is examined on rays traced to
30     * light sources. If a ray is exiting a dielectric material, we
31     * check the sources to see if any would cause bright color to be
32     * directed to the viewer due to dispersion. This gives colorful
33     * sparkle to crystals, etc. (Only if DISPERSE is defined!)
34     *
35     * Arguments for MAT_DIELECTRIC are:
36     * red grn blu rndx Hartmann
37     *
38     * Arguments for MAT_INTERFACE are:
39     * red1 grn1 blu1 rndx1 red2 grn2 blu2 rndx2
40     *
41     * The primaries are material transmission per unit length.
42     * MAT_INTERFACE uses dielectric1 for inside and dielectric2 for
43     * outside.
44     */
45    
46    
47     #define MLAMBDA 500 /* mean lambda */
48     #define MAXLAMBDA 779 /* maximum lambda */
49     #define MINLAMBDA 380 /* minimum lambda */
50    
51     #define MINCOS 0.997 /* minimum dot product for dispersion */
52    
53 greg 2.9
54 greg 2.10 static double
55 schorsch 2.18 mylog( /* special log for extinction coefficients */
56     double x
57     )
58 greg 2.10 {
59     if (x < 1e-40)
60     return(-100.);
61     if (x >= 1.)
62     return(0.);
63     return(log(x));
64     }
65    
66    
67 schorsch 2.18 extern int
68     m_dielectric( /* color a ray which hit a dielectric interface */
69     OBJREC *m,
70     register RAY *r
71     )
72 greg 1.1 {
73     double cos1, cos2, nratio;
74 greg 2.9 COLOR ctrans;
75 greg 2.11 COLOR talb;
76 gwlarson 2.14 int hastexture;
77 greg 1.5 double refl, trans;
78 greg 1.1 FVECT dnorm;
79     double d1, d2;
80     RAY p;
81     register int i;
82    
83     if (m->oargs.nfargs != (m->otype==MAT_DIELECTRIC ? 5 : 8))
84     objerror(m, USER, "bad arguments");
85    
86     raytexture(r, m->omod); /* get modifiers */
87    
88 schorsch 2.17 if ( (hastexture = DOT(r->pert,r->pert) > FTINY*FTINY) )
89 gwlarson 2.14 cos1 = raynormal(dnorm, r); /* perturb normal */
90     else {
91     VCOPY(dnorm, r->ron);
92     cos1 = r->rod;
93     }
94 greg 1.1 /* index of refraction */
95     if (m->otype == MAT_DIELECTRIC)
96     nratio = m->oargs.farg[3] + m->oargs.farg[4]/MLAMBDA;
97     else
98     nratio = m->oargs.farg[3] / m->oargs.farg[7];
99    
100     if (cos1 < 0.0) { /* inside */
101 gwlarson 2.14 hastexture = -hastexture;
102 greg 1.1 cos1 = -cos1;
103     dnorm[0] = -dnorm[0];
104     dnorm[1] = -dnorm[1];
105     dnorm[2] = -dnorm[2];
106 greg 2.10 setcolor(r->cext, -mylog(m->oargs.farg[0]*colval(r->pcol,RED)),
107     -mylog(m->oargs.farg[1]*colval(r->pcol,GRN)),
108     -mylog(m->oargs.farg[2]*colval(r->pcol,BLU)));
109 greg 2.11 setcolor(r->albedo, 0., 0., 0.);
110 greg 2.9 r->gecc = 0.;
111     if (m->otype == MAT_INTERFACE) {
112     setcolor(ctrans,
113 greg 2.10 -mylog(m->oargs.farg[4]*colval(r->pcol,RED)),
114     -mylog(m->oargs.farg[5]*colval(r->pcol,GRN)),
115     -mylog(m->oargs.farg[6]*colval(r->pcol,BLU)));
116 greg 2.11 setcolor(talb, 0., 0., 0.);
117 greg 2.9 } else {
118     copycolor(ctrans, cextinction);
119 greg 2.11 copycolor(talb, salbedo);
120 greg 2.9 }
121 greg 1.1 } else { /* outside */
122     nratio = 1.0 / nratio;
123 greg 2.9
124 greg 2.10 setcolor(ctrans, -mylog(m->oargs.farg[0]*colval(r->pcol,RED)),
125     -mylog(m->oargs.farg[1]*colval(r->pcol,GRN)),
126     -mylog(m->oargs.farg[2]*colval(r->pcol,BLU)));
127 greg 2.11 setcolor(talb, 0., 0., 0.);
128 greg 2.9 if (m->otype == MAT_INTERFACE) {
129     setcolor(r->cext,
130 greg 2.10 -mylog(m->oargs.farg[4]*colval(r->pcol,RED)),
131     -mylog(m->oargs.farg[5]*colval(r->pcol,GRN)),
132     -mylog(m->oargs.farg[6]*colval(r->pcol,BLU)));
133 greg 2.11 setcolor(r->albedo, 0., 0., 0.);
134 greg 2.9 r->gecc = 0.;
135     }
136 greg 1.1 }
137    
138     d2 = 1.0 - nratio*nratio*(1.0 - cos1*cos1); /* compute cos theta2 */
139    
140     if (d2 < FTINY) /* total reflection */
141    
142     refl = 1.0;
143    
144     else { /* refraction occurs */
145     /* compute Fresnel's equations */
146     cos2 = sqrt(d2);
147     d1 = cos1;
148     d2 = nratio*cos2;
149     d1 = (d1 - d2) / (d1 + d2);
150     refl = d1 * d1;
151    
152     d1 = 1.0 / cos1;
153     d2 = nratio / cos2;
154     d1 = (d1 - d2) / (d1 + d2);
155     refl += d1 * d1;
156    
157 greg 2.9 refl *= 0.5;
158 greg 1.1 trans = 1.0 - refl;
159 greg 2.15
160     trans *= nratio*nratio; /* solid angle ratio */
161 greg 1.1
162 gwlarson 2.13 if (rayorigin(&p, r, REFRACTED, trans) == 0) {
163 greg 1.1
164     /* compute refracted ray */
165     d1 = nratio*cos1 - cos2;
166     for (i = 0; i < 3; i++)
167     p.rdir[i] = nratio*r->rdir[i] + d1*dnorm[i];
168 gwlarson 2.14 /* accidental reflection? */
169     if (hastexture &&
170     DOT(p.rdir,r->ron)*hastexture >= -FTINY) {
171     d1 *= (double)hastexture;
172     for (i = 0; i < 3; i++) /* ignore texture */
173     p.rdir[i] = nratio*r->rdir[i] +
174     d1*r->ron[i];
175     normalize(p.rdir); /* not exact */
176     }
177 greg 1.1 #ifdef DISPERSE
178     if (m->otype != MAT_DIELECTRIC
179     || r->rod > 0.0
180     || r->crtype & SHADOW
181 greg 2.3 || !directvis
182 greg 1.1 || m->oargs.farg[4] == 0.0
183 greg 2.12 || !disperse(m, r, p.rdir,
184     trans, ctrans, talb))
185 greg 1.1 #endif
186     {
187 greg 2.9 copycolor(p.cext, ctrans);
188 greg 2.11 copycolor(p.albedo, talb);
189 greg 1.1 rayvalue(&p);
190     scalecolor(p.rcol, trans);
191     addcolor(r->rcol, p.rcol);
192 greg 2.4 if (nratio >= 1.0-FTINY && nratio <= 1.0+FTINY)
193     r->rt = r->rot + p.rt;
194 greg 1.1 }
195     }
196     }
197    
198     if (!(r->crtype & SHADOW) &&
199 gwlarson 2.13 rayorigin(&p, r, REFLECTED, refl) == 0) {
200 greg 1.1
201     /* compute reflected ray */
202     for (i = 0; i < 3; i++)
203     p.rdir[i] = r->rdir[i] + 2.0*cos1*dnorm[i];
204 gwlarson 2.14 /* accidental penetration? */
205     if (hastexture && DOT(p.rdir,r->ron)*hastexture <= FTINY)
206     for (i = 0; i < 3; i++) /* ignore texture */
207     p.rdir[i] = r->rdir[i] + 2.0*r->rod*r->ron[i];
208 greg 1.1
209     rayvalue(&p); /* reflected ray value */
210    
211     scalecolor(p.rcol, refl); /* color contribution */
212     addcolor(r->rcol, p.rcol);
213     }
214 greg 2.9 /* rayvalue() computes absorption */
215 greg 2.7 return(1);
216 greg 1.1 }
217    
218    
219     #ifdef DISPERSE
220    
221 schorsch 2.18 static int
222     disperse( /* check light sources for dispersion */
223     OBJREC *m,
224     RAY *r,
225     FVECT vt,
226     double tr,
227     COLOR cet,
228     COLOR abt
229     )
230 greg 1.1 {
231     RAY sray, *entray;
232     FVECT v1, v2, n1, n2;
233     FVECT dv, v2Xdv;
234     double v2Xdvv2Xdv;
235 greg 1.7 int success = 0;
236     SRCINDEX si;
237 greg 1.1 FVECT vtmp1, vtmp2;
238     double dtmp1, dtmp2;
239     int l1, l2;
240     COLOR ctmp;
241     int i;
242    
243     /*
244     * This routine computes dispersion to the first order using
245     * the following assumptions:
246     *
247     * 1) The dependency of the index of refraction on wavelength
248     * is approximated by Hartmann's equation with lambda0
249     * equal to zero.
250     * 2) The entry and exit locations are constant with respect
251     * to dispersion.
252     *
253     * The second assumption permits us to model dispersion without
254     * having to sample refracted directions. We assume that the
255     * geometry inside the material is constant, and concern ourselves
256     * only with the relationship between the entering and exiting ray.
257     * We compute the first derivatives of the entering and exiting
258     * refraction with respect to the index of refraction. This
259     * is then used in a first order Taylor series to determine the
260     * index of refraction necessary to send the exiting ray to each
261     * light source.
262     * If an exiting ray hits a light source within the refraction
263     * boundaries, we sum all the frequencies over the disc of the
264     * light source to determine the resulting color. A smaller light
265     * source will therefore exhibit a sharper spectrum.
266     */
267    
268     if (!(r->crtype & REFRACTED)) { /* ray started in material */
269     VCOPY(v1, r->rdir);
270     n1[0] = -r->rdir[0]; n1[1] = -r->rdir[1]; n1[2] = -r->rdir[2];
271     } else {
272     /* find entry point */
273     for (entray = r; entray->rtype != REFRACTED;
274     entray = entray->parent)
275     ;
276     entray = entray->parent;
277     if (entray->crtype & REFRACTED) /* too difficult */
278     return(0);
279     VCOPY(v1, entray->rdir);
280     VCOPY(n1, entray->ron);
281     }
282     VCOPY(v2, vt); /* exiting ray */
283     VCOPY(n2, r->ron);
284    
285     /* first order dispersion approx. */
286     dtmp1 = DOT(n1, v1);
287     dtmp2 = DOT(n2, v2);
288     for (i = 0; i < 3; i++)
289     dv[i] = v1[i] + v2[i] - n1[i]/dtmp1 - n2[i]/dtmp2;
290    
291     if (DOT(dv, dv) <= FTINY) /* null effect */
292     return(0);
293     /* compute plane normal */
294     fcross(v2Xdv, v2, dv);
295     v2Xdvv2Xdv = DOT(v2Xdv, v2Xdv);
296    
297     /* check sources */
298 greg 1.7 initsrcindex(&si);
299     while (srcray(&sray, r, &si)) {
300 greg 1.1
301 greg 1.7 if (DOT(sray.rdir, v2) < MINCOS)
302 greg 1.1 continue; /* bad source */
303     /* adjust source ray */
304    
305     dtmp1 = DOT(v2Xdv, sray.rdir) / v2Xdvv2Xdv;
306     sray.rdir[0] -= dtmp1 * v2Xdv[0];
307     sray.rdir[1] -= dtmp1 * v2Xdv[1];
308     sray.rdir[2] -= dtmp1 * v2Xdv[2];
309    
310     l1 = lambda(m, v2, dv, sray.rdir); /* mean lambda */
311    
312     if (l1 > MAXLAMBDA || l1 < MINLAMBDA) /* not visible */
313     continue;
314     /* trace source ray */
315 greg 2.12 copycolor(sray.cext, cet);
316     copycolor(sray.albedo, abt);
317 greg 1.1 normalize(sray.rdir);
318     rayvalue(&sray);
319 greg 1.2 if (bright(sray.rcol) <= FTINY) /* missed it */
320 greg 1.1 continue;
321    
322     /*
323     * Compute spectral sum over diameter of source.
324     * First find directions for rays going to opposite
325     * sides of source, then compute wavelengths for each.
326     */
327    
328     fcross(vtmp1, v2Xdv, sray.rdir);
329 greg 1.7 dtmp1 = sqrt(si.dom / v2Xdvv2Xdv / PI);
330 greg 1.1
331     /* compute first ray */
332     for (i = 0; i < 3; i++)
333     vtmp2[i] = sray.rdir[i] + dtmp1*vtmp1[i];
334    
335     l1 = lambda(m, v2, dv, vtmp2); /* first lambda */
336     if (l1 < 0)
337     continue;
338     /* compute second ray */
339     for (i = 0; i < 3; i++)
340     vtmp2[i] = sray.rdir[i] - dtmp1*vtmp1[i];
341    
342     l2 = lambda(m, v2, dv, vtmp2); /* second lambda */
343     if (l2 < 0)
344     continue;
345     /* compute color from spectrum */
346     if (l1 < l2)
347     spec_rgb(ctmp, l1, l2);
348     else
349     spec_rgb(ctmp, l2, l1);
350     multcolor(ctmp, sray.rcol);
351     scalecolor(ctmp, tr);
352     addcolor(r->rcol, ctmp);
353     success++;
354     }
355     return(success);
356     }
357    
358    
359     static int
360 schorsch 2.18 lambda( /* compute lambda for material */
361     register OBJREC *m,
362     FVECT v2,
363     FVECT dv,
364     FVECT lr
365     )
366 greg 1.1 {
367     FVECT lrXdv, v2Xlr;
368     double dtmp, denom;
369     int i;
370    
371     fcross(lrXdv, lr, dv);
372     for (i = 0; i < 3; i++)
373     if (lrXdv[i] > FTINY || lrXdv[i] < -FTINY)
374     break;
375     if (i >= 3)
376     return(-1);
377    
378     fcross(v2Xlr, v2, lr);
379    
380     dtmp = m->oargs.farg[4] / MLAMBDA;
381     denom = dtmp + v2Xlr[i]/lrXdv[i] * (m->oargs.farg[3] + dtmp);
382    
383     if (denom < FTINY)
384     return(-1);
385    
386     return(m->oargs.farg[4] / denom);
387     }
388    
389     #endif /* DISPERSE */