ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/normal.c
Revision: 2.48
Committed: Mon Sep 20 17:32:04 2004 UTC (19 years, 7 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad3R6, rad3R6P1
Changes since 2.47: +10 -9 lines
Log Message:
Corrected Gaussian reflectance model normalization (cosine factor)

File Contents

# Content
1 #ifndef lint
2 static const char RCSid[] = "$Id: normal.c,v 2.47 2004/03/30 16:13:01 schorsch Exp $";
3 #endif
4 /*
5 * normal.c - shading function for normal materials.
6 *
7 * 8/19/85
8 * 12/19/85 - added stuff for metals.
9 * 6/26/87 - improved specular model.
10 * 9/28/87 - added model for translucent materials.
11 * Later changes described in delta comments.
12 */
13
14 #include "copyright.h"
15
16 #include "ray.h"
17 #include "ambient.h"
18 #include "source.h"
19 #include "otypes.h"
20 #include "rtotypes.h"
21 #include "random.h"
22
23 #ifndef MAXITER
24 #define MAXITER 10 /* maximum # specular ray attempts */
25 #endif
26 /* estimate of Fresnel function */
27 #define FRESNE(ci) (exp(-5.85*(ci)) - 0.00287989916)
28
29
30 /*
31 * This routine implements the isotropic Gaussian
32 * model described by Ward in Siggraph `92 article.
33 * We orient the surface towards the incoming ray, so a single
34 * surface can be used to represent an infinitely thin object.
35 *
36 * Arguments for MAT_PLASTIC and MAT_METAL are:
37 * red grn blu specular-frac. facet-slope
38 *
39 * Arguments for MAT_TRANS are:
40 * red grn blu rspec rough trans tspec
41 */
42
43 /* specularity flags */
44 #define SP_REFL 01 /* has reflected specular component */
45 #define SP_TRAN 02 /* has transmitted specular */
46 #define SP_PURE 04 /* purely specular (zero roughness) */
47 #define SP_FLAT 010 /* flat reflecting surface */
48 #define SP_RBLT 020 /* reflection below sample threshold */
49 #define SP_TBLT 040 /* transmission below threshold */
50
51 typedef struct {
52 OBJREC *mp; /* material pointer */
53 RAY *rp; /* ray pointer */
54 short specfl; /* specularity flags, defined above */
55 COLOR mcolor; /* color of this material */
56 COLOR scolor; /* color of specular component */
57 FVECT vrefl; /* vector in direction of reflected ray */
58 FVECT prdir; /* vector in transmitted direction */
59 double alpha2; /* roughness squared */
60 double rdiff, rspec; /* reflected specular, diffuse */
61 double trans; /* transmissivity */
62 double tdiff, tspec; /* transmitted specular, diffuse */
63 FVECT pnorm; /* perturbed surface normal */
64 double pdot; /* perturbed dot product */
65 } NORMDAT; /* normal material data */
66
67 static srcdirf_t dirnorm;
68 static void gaussamp(RAY *r, NORMDAT *np);
69
70
71 static void
72 dirnorm( /* compute source contribution */
73 COLOR cval, /* returned coefficient */
74 void *nnp, /* material data */
75 FVECT ldir, /* light source direction */
76 double omega /* light source size */
77 )
78 {
79 register NORMDAT *np = nnp;
80 double ldot;
81 double ldiff;
82 double dtmp, d2;
83 FVECT vtmp;
84 COLOR ctmp;
85
86 setcolor(cval, 0.0, 0.0, 0.0);
87
88 ldot = DOT(np->pnorm, ldir);
89
90 if (ldot < 0.0 ? np->trans <= FTINY : np->trans >= 1.0-FTINY)
91 return; /* wrong side */
92
93 /* Fresnel estimate */
94 ldiff = np->rdiff;
95 if (np->specfl & SP_PURE && (np->rspec > FTINY) & (ldiff > FTINY))
96 ldiff *= 1. - FRESNE(fabs(ldot));
97
98 if (ldot > FTINY && ldiff > FTINY) {
99 /*
100 * Compute and add diffuse reflected component to returned
101 * color. The diffuse reflected component will always be
102 * modified by the color of the material.
103 */
104 copycolor(ctmp, np->mcolor);
105 dtmp = ldot * omega * ldiff * (1.0/PI);
106 scalecolor(ctmp, dtmp);
107 addcolor(cval, ctmp);
108 }
109 if (ldot > FTINY && (np->specfl&(SP_REFL|SP_PURE)) == SP_REFL) {
110 /*
111 * Compute specular reflection coefficient using
112 * gaussian distribution model.
113 */
114 /* roughness */
115 dtmp = np->alpha2;
116 /* + source if flat */
117 if (np->specfl & SP_FLAT)
118 dtmp += omega * (0.25/PI);
119 /* half vector */
120 vtmp[0] = ldir[0] - np->rp->rdir[0];
121 vtmp[1] = ldir[1] - np->rp->rdir[1];
122 vtmp[2] = ldir[2] - np->rp->rdir[2];
123 d2 = DOT(vtmp, np->pnorm);
124 d2 *= d2;
125 d2 = (DOT(vtmp,vtmp) - d2) / d2;
126 /* gaussian */
127 dtmp = exp(-d2/dtmp)/(4.*PI * np->pdot * dtmp);
128 /* worth using? */
129 if (dtmp > FTINY) {
130 copycolor(ctmp, np->scolor);
131 dtmp *= omega;
132 scalecolor(ctmp, dtmp);
133 addcolor(cval, ctmp);
134 }
135 }
136 if (ldot < -FTINY && np->tdiff > FTINY) {
137 /*
138 * Compute diffuse transmission.
139 */
140 copycolor(ctmp, np->mcolor);
141 dtmp = -ldot * omega * np->tdiff * (1.0/PI);
142 scalecolor(ctmp, dtmp);
143 addcolor(cval, ctmp);
144 }
145 if (ldot < -FTINY && (np->specfl&(SP_TRAN|SP_PURE)) == SP_TRAN) {
146 /*
147 * Compute specular transmission. Specular transmission
148 * is always modified by material color.
149 */
150 /* roughness + source */
151 dtmp = np->alpha2 + omega*(1.0/PI);
152 /* gaussian */
153 dtmp = exp((2.*DOT(np->prdir,ldir)-2.)/dtmp) /
154 (PI*np->pdot*dtmp);
155 /* worth using? */
156 if (dtmp > FTINY) {
157 copycolor(ctmp, np->mcolor);
158 dtmp *= np->tspec * omega;
159 scalecolor(ctmp, dtmp);
160 addcolor(cval, ctmp);
161 }
162 }
163 }
164
165
166 extern int
167 m_normal( /* color a ray that hit something normal */
168 register OBJREC *m,
169 register RAY *r
170 )
171 {
172 NORMDAT nd;
173 double fest;
174 double transtest, transdist;
175 double mirtest, mirdist;
176 int hastexture;
177 double d;
178 COLOR ctmp;
179 register int i;
180 /* easy shadow test */
181 if (r->crtype & SHADOW && m->otype != MAT_TRANS)
182 return(1);
183
184 if (m->oargs.nfargs != (m->otype == MAT_TRANS ? 7 : 5))
185 objerror(m, USER, "bad number of arguments");
186 /* check for back side */
187 if (r->rod < 0.0) {
188 if (!backvis && m->otype != MAT_TRANS) {
189 raytrans(r);
190 return(1);
191 }
192 raytexture(r, m->omod);
193 flipsurface(r); /* reorient if backvis */
194 } else
195 raytexture(r, m->omod);
196 nd.mp = m;
197 nd.rp = r;
198 /* get material color */
199 setcolor(nd.mcolor, m->oargs.farg[0],
200 m->oargs.farg[1],
201 m->oargs.farg[2]);
202 /* get roughness */
203 nd.specfl = 0;
204 nd.alpha2 = m->oargs.farg[4];
205 if ((nd.alpha2 *= nd.alpha2) <= FTINY)
206 nd.specfl |= SP_PURE;
207
208 if ( (hastexture = (DOT(r->pert,r->pert) > FTINY*FTINY)) ) {
209 nd.pdot = raynormal(nd.pnorm, r); /* perturb normal */
210 } else {
211 VCOPY(nd.pnorm, r->ron);
212 nd.pdot = r->rod;
213 }
214 if (r->ro != NULL && isflat(r->ro->otype))
215 nd.specfl |= SP_FLAT;
216 if (nd.pdot < .001)
217 nd.pdot = .001; /* non-zero for dirnorm() */
218 multcolor(nd.mcolor, r->pcol); /* modify material color */
219 mirtest = transtest = 0;
220 mirdist = transdist = r->rot;
221 nd.rspec = m->oargs.farg[3];
222 /* compute Fresnel approx. */
223 if (nd.specfl & SP_PURE && nd.rspec > FTINY) {
224 fest = FRESNE(r->rod);
225 nd.rspec += fest*(1. - nd.rspec);
226 } else
227 fest = 0.;
228 /* compute transmission */
229 if (m->otype == MAT_TRANS) {
230 nd.trans = m->oargs.farg[5]*(1.0 - nd.rspec);
231 nd.tspec = nd.trans * m->oargs.farg[6];
232 nd.tdiff = nd.trans - nd.tspec;
233 if (nd.tspec > FTINY) {
234 nd.specfl |= SP_TRAN;
235 /* check threshold */
236 if (!(nd.specfl & SP_PURE) &&
237 specthresh >= nd.tspec-FTINY)
238 nd.specfl |= SP_TBLT;
239 if (!hastexture || r->crtype & SHADOW) {
240 VCOPY(nd.prdir, r->rdir);
241 transtest = 2;
242 } else {
243 for (i = 0; i < 3; i++) /* perturb */
244 nd.prdir[i] = r->rdir[i] - r->pert[i];
245 if (DOT(nd.prdir, r->ron) < -FTINY)
246 normalize(nd.prdir); /* OK */
247 else
248 VCOPY(nd.prdir, r->rdir);
249 }
250 }
251 } else
252 nd.tdiff = nd.tspec = nd.trans = 0.0;
253 /* transmitted ray */
254 if ((nd.specfl&(SP_TRAN|SP_PURE|SP_TBLT)) == (SP_TRAN|SP_PURE)) {
255 RAY lr;
256 if (rayorigin(&lr, r, TRANS, nd.tspec) == 0) {
257 VCOPY(lr.rdir, nd.prdir);
258 rayvalue(&lr);
259 scalecolor(lr.rcol, nd.tspec);
260 multcolor(lr.rcol, nd.mcolor); /* modified by color */
261 addcolor(r->rcol, lr.rcol);
262 transtest *= bright(lr.rcol);
263 transdist = r->rot + lr.rt;
264 }
265 } else
266 transtest = 0;
267
268 if (r->crtype & SHADOW) { /* the rest is shadow */
269 r->rt = transdist;
270 return(1);
271 }
272 /* get specular reflection */
273 if (nd.rspec > FTINY) {
274 nd.specfl |= SP_REFL;
275 /* compute specular color */
276 if (m->otype != MAT_METAL) {
277 setcolor(nd.scolor, nd.rspec, nd.rspec, nd.rspec);
278 } else if (fest > FTINY) {
279 d = nd.rspec*(1. - fest);
280 for (i = 0; i < 3; i++)
281 nd.scolor[i] = fest + nd.mcolor[i]*d;
282 } else {
283 copycolor(nd.scolor, nd.mcolor);
284 scalecolor(nd.scolor, nd.rspec);
285 }
286 /* check threshold */
287 if (!(nd.specfl & SP_PURE) && specthresh >= nd.rspec-FTINY)
288 nd.specfl |= SP_RBLT;
289 /* compute reflected ray */
290 for (i = 0; i < 3; i++)
291 nd.vrefl[i] = r->rdir[i] + 2.*nd.pdot*nd.pnorm[i];
292 /* penetration? */
293 if (hastexture && DOT(nd.vrefl, r->ron) <= FTINY)
294 for (i = 0; i < 3; i++) /* safety measure */
295 nd.vrefl[i] = r->rdir[i] + 2.*r->rod*r->ron[i];
296 }
297 /* reflected ray */
298 if ((nd.specfl&(SP_REFL|SP_PURE|SP_RBLT)) == (SP_REFL|SP_PURE)) {
299 RAY lr;
300 if (rayorigin(&lr, r, REFLECTED, nd.rspec) == 0) {
301 VCOPY(lr.rdir, nd.vrefl);
302 rayvalue(&lr);
303 multcolor(lr.rcol, nd.scolor);
304 addcolor(r->rcol, lr.rcol);
305 if (!hastexture && nd.specfl & SP_FLAT) {
306 mirtest = 2.*bright(lr.rcol);
307 mirdist = r->rot + lr.rt;
308 }
309 }
310 }
311 /* diffuse reflection */
312 nd.rdiff = 1.0 - nd.trans - nd.rspec;
313
314 if (nd.specfl & SP_PURE && nd.rdiff <= FTINY && nd.tdiff <= FTINY)
315 return(1); /* 100% pure specular */
316
317 if (!(nd.specfl & SP_PURE))
318 gaussamp(r, &nd); /* checks *BLT flags */
319
320 if (nd.rdiff > FTINY) { /* ambient from this side */
321 ambient(ctmp, r, hastexture?nd.pnorm:r->ron);
322 if (nd.specfl & SP_RBLT)
323 scalecolor(ctmp, 1.0-nd.trans);
324 else
325 scalecolor(ctmp, nd.rdiff);
326 multcolor(ctmp, nd.mcolor); /* modified by material color */
327 addcolor(r->rcol, ctmp); /* add to returned color */
328 }
329 if (nd.tdiff > FTINY) { /* ambient from other side */
330 flipsurface(r);
331 if (hastexture) {
332 FVECT bnorm;
333 bnorm[0] = -nd.pnorm[0];
334 bnorm[1] = -nd.pnorm[1];
335 bnorm[2] = -nd.pnorm[2];
336 ambient(ctmp, r, bnorm);
337 } else
338 ambient(ctmp, r, r->ron);
339 if (nd.specfl & SP_TBLT)
340 scalecolor(ctmp, nd.trans);
341 else
342 scalecolor(ctmp, nd.tdiff);
343 multcolor(ctmp, nd.mcolor); /* modified by color */
344 addcolor(r->rcol, ctmp);
345 flipsurface(r);
346 }
347 /* add direct component */
348 direct(r, dirnorm, &nd);
349 /* check distance */
350 d = bright(r->rcol);
351 if (transtest > d)
352 r->rt = transdist;
353 else if (mirtest > d)
354 r->rt = mirdist;
355
356 return(1);
357 }
358
359
360 static void
361 gaussamp( /* sample gaussian specular */
362 RAY *r,
363 register NORMDAT *np
364 )
365 {
366 RAY sr;
367 FVECT u, v, h;
368 double rv[2];
369 double d, sinp, cosp;
370 int niter;
371 register int i;
372 /* quick test */
373 if ((np->specfl & (SP_REFL|SP_RBLT)) != SP_REFL &&
374 (np->specfl & (SP_TRAN|SP_TBLT)) != SP_TRAN)
375 return;
376 /* set up sample coordinates */
377 v[0] = v[1] = v[2] = 0.0;
378 for (i = 0; i < 3; i++)
379 if (np->pnorm[i] < 0.6 && np->pnorm[i] > -0.6)
380 break;
381 v[i] = 1.0;
382 fcross(u, v, np->pnorm);
383 normalize(u);
384 fcross(v, np->pnorm, u);
385 /* compute reflection */
386 if ((np->specfl & (SP_REFL|SP_RBLT)) == SP_REFL &&
387 rayorigin(&sr, r, SPECULAR, np->rspec) == 0) {
388 dimlist[ndims++] = (int)np->mp;
389 for (niter = 0; niter < MAXITER; niter++) {
390 if (niter)
391 d = frandom();
392 else
393 d = urand(ilhash(dimlist,ndims)+samplendx);
394 multisamp(rv, 2, d);
395 d = 2.0*PI * rv[0];
396 cosp = tcos(d);
397 sinp = tsin(d);
398 rv[1] = 1.0 - specjitter*rv[1];
399 if (rv[1] <= FTINY)
400 d = 1.0;
401 else
402 d = sqrt( np->alpha2 * -log(rv[1]) );
403 for (i = 0; i < 3; i++)
404 h[i] = np->pnorm[i] + d*(cosp*u[i] + sinp*v[i]);
405 d = -2.0 * DOT(h, r->rdir) / (1.0 + d*d);
406 for (i = 0; i < 3; i++)
407 sr.rdir[i] = r->rdir[i] + d*h[i];
408 if (DOT(sr.rdir, r->ron) > FTINY) {
409 rayvalue(&sr);
410 multcolor(sr.rcol, np->scolor);
411 addcolor(r->rcol, sr.rcol);
412 break;
413 }
414 }
415 ndims--;
416 }
417 /* compute transmission */
418 if ((np->specfl & (SP_TRAN|SP_TBLT)) == SP_TRAN &&
419 rayorigin(&sr, r, SPECULAR, np->tspec) == 0) {
420 dimlist[ndims++] = (int)np->mp;
421 for (niter = 0; niter < MAXITER; niter++) {
422 if (niter)
423 d = frandom();
424 else
425 d = urand(ilhash(dimlist,ndims)+1823+samplendx);
426 multisamp(rv, 2, d);
427 d = 2.0*PI * rv[0];
428 cosp = tcos(d);
429 sinp = tsin(d);
430 rv[1] = 1.0 - specjitter*rv[1];
431 if (rv[1] <= FTINY)
432 d = 1.0;
433 else
434 d = sqrt( np->alpha2 * -log(rv[1]) );
435 for (i = 0; i < 3; i++)
436 sr.rdir[i] = np->prdir[i] + d*(cosp*u[i] + sinp*v[i]);
437 if (DOT(sr.rdir, r->ron) < -FTINY) {
438 normalize(sr.rdir); /* OK, normalize */
439 rayvalue(&sr);
440 scalecolor(sr.rcol, np->tspec);
441 multcolor(sr.rcol, np->mcolor); /* modified */
442 addcolor(r->rcol, sr.rcol);
443 break;
444 }
445 }
446 ndims--;
447 }
448 }