ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/aniso.c
Revision: 2.25
Committed: Thu May 27 15:27:57 1993 UTC (30 years, 11 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.24: +2 -14 lines
Log Message:
removed approximation to Fresnel reflection, since the direct
component was not being calculated accordingly (BRDF not bidirectional)
changed BRDTfunc arguments and operation, fixed one or two bugs there also

File Contents

# Content
1 /* Copyright (c) 1992 Regents of the University of California */
2
3 #ifndef lint
4 static char SCCSid[] = "$SunId$ LBL";
5 #endif
6
7 /*
8 * Shading functions for anisotropic materials.
9 */
10
11 #include "ray.h"
12
13 #include "otypes.h"
14
15 #include "func.h"
16
17 #include "random.h"
18
19 extern double specthresh; /* specular sampling threshold */
20 extern double specjitter; /* specular sampling jitter */
21
22 static agaussamp();
23
24 /*
25 * This routine implements the anisotropic Gaussian
26 * model described by Ward in Siggraph `92 article.
27 * We orient the surface towards the incoming ray, so a single
28 * surface can be used to represent an infinitely thin object.
29 *
30 * Arguments for MAT_PLASTIC2 and MAT_METAL2 are:
31 * 4+ ux uy uz funcfile [transform...]
32 * 0
33 * 6 red grn blu specular-frac. u-facet-slope v-facet-slope
34 *
35 * Real arguments for MAT_TRANS2 are:
36 * 8 red grn blu rspec u-rough v-rough trans tspec
37 */
38
39 /* specularity flags */
40 #define SP_REFL 01 /* has reflected specular component */
41 #define SP_TRAN 02 /* has transmitted specular */
42 #define SP_FLAT 04 /* reflecting surface is flat */
43 #define SP_RBLT 010 /* reflection below sample threshold */
44 #define SP_TBLT 020 /* transmission below threshold */
45 #define SP_BADU 040 /* bad u direction calculation */
46
47 typedef struct {
48 OBJREC *mp; /* material pointer */
49 RAY *rp; /* ray pointer */
50 short specfl; /* specularity flags, defined above */
51 COLOR mcolor; /* color of this material */
52 COLOR scolor; /* color of specular component */
53 FVECT vrefl; /* vector in reflected direction */
54 FVECT prdir; /* vector in transmitted direction */
55 FVECT u, v; /* u and v vectors orienting anisotropy */
56 double u_alpha; /* u roughness */
57 double v_alpha; /* v roughness */
58 double rdiff, rspec; /* reflected specular, diffuse */
59 double trans; /* transmissivity */
60 double tdiff, tspec; /* transmitted specular, diffuse */
61 FVECT pnorm; /* perturbed surface normal */
62 double pdot; /* perturbed dot product */
63 } ANISODAT; /* anisotropic material data */
64
65
66 diraniso(cval, np, ldir, omega) /* compute source contribution */
67 COLOR cval; /* returned coefficient */
68 register ANISODAT *np; /* material data */
69 FVECT ldir; /* light source direction */
70 double omega; /* light source size */
71 {
72 double ldot;
73 double dtmp, dtmp1, dtmp2;
74 FVECT h;
75 double au2, av2;
76 COLOR ctmp;
77
78 setcolor(cval, 0.0, 0.0, 0.0);
79
80 ldot = DOT(np->pnorm, ldir);
81
82 if (ldot < 0.0 ? np->trans <= FTINY : np->trans >= 1.0-FTINY)
83 return; /* wrong side */
84
85 if (ldot > FTINY && np->rdiff > FTINY) {
86 /*
87 * Compute and add diffuse reflected component to returned
88 * color. The diffuse reflected component will always be
89 * modified by the color of the material.
90 */
91 copycolor(ctmp, np->mcolor);
92 dtmp = ldot * omega * np->rdiff / PI;
93 scalecolor(ctmp, dtmp);
94 addcolor(cval, ctmp);
95 }
96 if (ldot > FTINY && (np->specfl&(SP_REFL|SP_BADU)) == SP_REFL) {
97 /*
98 * Compute specular reflection coefficient using
99 * anisotropic gaussian distribution model.
100 */
101 /* add source width if flat */
102 if (np->specfl & SP_FLAT)
103 au2 = av2 = omega/(4.0*PI);
104 else
105 au2 = av2 = 0.0;
106 au2 += np->u_alpha*np->u_alpha;
107 av2 += np->v_alpha*np->v_alpha;
108 /* half vector */
109 h[0] = ldir[0] - np->rp->rdir[0];
110 h[1] = ldir[1] - np->rp->rdir[1];
111 h[2] = ldir[2] - np->rp->rdir[2];
112 /* ellipse */
113 dtmp1 = DOT(np->u, h);
114 dtmp1 *= dtmp1 / au2;
115 dtmp2 = DOT(np->v, h);
116 dtmp2 *= dtmp2 / av2;
117 /* gaussian */
118 dtmp = DOT(np->pnorm, h);
119 dtmp = (dtmp1 + dtmp2) / (dtmp*dtmp);
120 dtmp = exp(-dtmp) * (0.25/PI)
121 * sqrt(ldot/(np->pdot*au2*av2));
122 /* worth using? */
123 if (dtmp > FTINY) {
124 copycolor(ctmp, np->scolor);
125 dtmp *= omega;
126 scalecolor(ctmp, dtmp);
127 addcolor(cval, ctmp);
128 }
129 }
130 if (ldot < -FTINY && np->tdiff > FTINY) {
131 /*
132 * Compute diffuse transmission.
133 */
134 copycolor(ctmp, np->mcolor);
135 dtmp = -ldot * omega * np->tdiff / PI;
136 scalecolor(ctmp, dtmp);
137 addcolor(cval, ctmp);
138 }
139 if (ldot < -FTINY && (np->specfl&(SP_TRAN|SP_BADU)) == SP_TRAN) {
140 /*
141 * Compute specular transmission. Specular transmission
142 * is always modified by material color.
143 */
144 /* roughness + source */
145 au2 = av2 = omega / PI;
146 au2 += np->u_alpha*np->u_alpha;
147 av2 += np->v_alpha*np->v_alpha;
148 /* "half vector" */
149 h[0] = ldir[0] - np->prdir[0];
150 h[1] = ldir[1] - np->prdir[1];
151 h[2] = ldir[2] - np->prdir[2];
152 dtmp = DOT(h,h);
153 if (dtmp > FTINY*FTINY) {
154 dtmp1 = DOT(h,np->pnorm);
155 dtmp = 1.0 - dtmp1*dtmp1/dtmp;
156 if (dtmp > FTINY*FTINY) {
157 dtmp1 = DOT(h,np->u);
158 dtmp1 *= dtmp1 / au2;
159 dtmp2 = DOT(h,np->v);
160 dtmp2 *= dtmp2 / av2;
161 dtmp = (dtmp1 + dtmp2) / dtmp;
162 }
163 } else
164 dtmp = 0.0;
165 /* gaussian */
166 dtmp = exp(-dtmp) * (1.0/PI)
167 * sqrt(-ldot/(np->pdot*au2*av2));
168 /* worth using? */
169 if (dtmp > FTINY) {
170 copycolor(ctmp, np->mcolor);
171 dtmp *= np->tspec * omega;
172 scalecolor(ctmp, dtmp);
173 addcolor(cval, ctmp);
174 }
175 }
176 }
177
178
179 m_aniso(m, r) /* shade ray that hit something anisotropic */
180 register OBJREC *m;
181 register RAY *r;
182 {
183 ANISODAT nd;
184 COLOR ctmp;
185 register int i;
186 /* easy shadow test */
187 if (r->crtype & SHADOW)
188 return;
189
190 if (m->oargs.nfargs != (m->otype == MAT_TRANS2 ? 8 : 6))
191 objerror(m, USER, "bad number of real arguments");
192 nd.mp = m;
193 nd.rp = r;
194 /* get material color */
195 setcolor(nd.mcolor, m->oargs.farg[0],
196 m->oargs.farg[1],
197 m->oargs.farg[2]);
198 /* get roughness */
199 nd.specfl = 0;
200 nd.u_alpha = m->oargs.farg[4];
201 nd.v_alpha = m->oargs.farg[5];
202 if (nd.u_alpha < FTINY || nd.v_alpha <= FTINY)
203 objerror(m, USER, "roughness too small");
204 /* reorient if necessary */
205 if (r->rod < 0.0)
206 flipsurface(r);
207 /* get modifiers */
208 raytexture(r, m->omod);
209 nd.pdot = raynormal(nd.pnorm, r); /* perturb normal */
210 if (nd.pdot < .001)
211 nd.pdot = .001; /* non-zero for diraniso() */
212 multcolor(nd.mcolor, r->pcol); /* modify material color */
213 /* get specular component */
214 if ((nd.rspec = m->oargs.farg[3]) > FTINY) {
215 nd.specfl |= SP_REFL;
216 /* compute specular color */
217 if (m->otype == MAT_METAL2)
218 copycolor(nd.scolor, nd.mcolor);
219 else
220 setcolor(nd.scolor, 1.0, 1.0, 1.0);
221 scalecolor(nd.scolor, nd.rspec);
222 /* check threshold */
223 if (specthresh >= nd.rspec-FTINY)
224 nd.specfl |= SP_RBLT;
225 /* compute refl. direction */
226 for (i = 0; i < 3; i++)
227 nd.vrefl[i] = r->rdir[i] + 2.0*nd.pdot*nd.pnorm[i];
228 if (DOT(nd.vrefl, r->ron) <= FTINY) /* penetration? */
229 for (i = 0; i < 3; i++) /* safety measure */
230 nd.vrefl[i] = r->rdir[i] + 2.*r->rod*r->ron[i];
231 }
232 /* compute transmission */
233 if (m->otype == MAT_TRANS2) {
234 nd.trans = m->oargs.farg[6]*(1.0 - nd.rspec);
235 nd.tspec = nd.trans * m->oargs.farg[7];
236 nd.tdiff = nd.trans - nd.tspec;
237 if (nd.tspec > FTINY) {
238 nd.specfl |= SP_TRAN;
239 /* check threshold */
240 if (specthresh >= nd.tspec-FTINY)
241 nd.specfl |= SP_TBLT;
242 if (DOT(r->pert,r->pert) <= FTINY*FTINY) {
243 VCOPY(nd.prdir, r->rdir);
244 } else {
245 for (i = 0; i < 3; i++) /* perturb */
246 nd.prdir[i] = r->rdir[i] - r->pert[i];
247 if (DOT(nd.prdir, r->ron) < -FTINY)
248 normalize(nd.prdir); /* OK */
249 else
250 VCOPY(nd.prdir, r->rdir);
251 }
252 }
253 } else
254 nd.tdiff = nd.tspec = nd.trans = 0.0;
255
256 /* diffuse reflection */
257 nd.rdiff = 1.0 - nd.trans - nd.rspec;
258
259 if (r->ro != NULL && (r->ro->otype == OBJ_FACE ||
260 r->ro->otype == OBJ_RING))
261 nd.specfl |= SP_FLAT;
262
263 getacoords(r, &nd); /* set up coordinates */
264
265 if (nd.specfl & (SP_REFL|SP_TRAN) && !(nd.specfl & SP_BADU))
266 agaussamp(r, &nd);
267
268 if (nd.rdiff > FTINY) { /* ambient from this side */
269 ambient(ctmp, r);
270 if (nd.specfl & SP_RBLT)
271 scalecolor(ctmp, 1.0-nd.trans);
272 else
273 scalecolor(ctmp, nd.rdiff);
274 multcolor(ctmp, nd.mcolor); /* modified by material color */
275 addcolor(r->rcol, ctmp); /* add to returned color */
276 }
277 if (nd.tdiff > FTINY) { /* ambient from other side */
278 flipsurface(r);
279 ambient(ctmp, r);
280 if (nd.specfl & SP_TBLT)
281 scalecolor(ctmp, nd.trans);
282 else
283 scalecolor(ctmp, nd.tdiff);
284 multcolor(ctmp, nd.mcolor); /* modified by color */
285 addcolor(r->rcol, ctmp);
286 flipsurface(r);
287 }
288 /* add direct component */
289 direct(r, diraniso, &nd);
290 }
291
292
293 static
294 getacoords(r, np) /* set up coordinate system */
295 RAY *r;
296 register ANISODAT *np;
297 {
298 register MFUNC *mf;
299 register int i;
300
301 mf = getfunc(np->mp, 3, 0x7, 1);
302 setfunc(np->mp, r);
303 errno = 0;
304 for (i = 0; i < 3; i++)
305 np->u[i] = evalue(mf->ep[i]);
306 if (errno) {
307 objerror(np->mp, WARNING, "compute error");
308 np->specfl |= SP_BADU;
309 return;
310 }
311 if (mf->f != &unitxf)
312 multv3(np->u, np->u, mf->f->xfm);
313 fcross(np->v, np->pnorm, np->u);
314 if (normalize(np->v) == 0.0) {
315 objerror(np->mp, WARNING, "illegal orientation vector");
316 np->specfl |= SP_BADU;
317 return;
318 }
319 fcross(np->u, np->v, np->pnorm);
320 }
321
322
323 static
324 agaussamp(r, np) /* sample anisotropic gaussian specular */
325 RAY *r;
326 register ANISODAT *np;
327 {
328 RAY sr;
329 FVECT h;
330 double rv[2];
331 double d, sinp, cosp;
332 register int i;
333 /* compute reflection */
334 if ((np->specfl & (SP_REFL|SP_RBLT)) == SP_REFL &&
335 rayorigin(&sr, r, SPECULAR, np->rspec) == 0) {
336 dimlist[ndims++] = (int)np->mp;
337 d = urand(ilhash(dimlist,ndims)+samplendx);
338 multisamp(rv, 2, d);
339 d = 2.0*PI * rv[0];
340 cosp = cos(d) * np->u_alpha;
341 sinp = sin(d) * np->v_alpha;
342 d = sqrt(cosp*cosp + sinp*sinp);
343 cosp /= d;
344 sinp /= d;
345 rv[1] = 1.0 - specjitter*rv[1];
346 if (rv[1] <= FTINY)
347 d = 1.0;
348 else
349 d = sqrt(-log(rv[1]) /
350 (cosp*cosp/(np->u_alpha*np->u_alpha) +
351 sinp*sinp/(np->v_alpha*np->v_alpha)));
352 for (i = 0; i < 3; i++)
353 h[i] = np->pnorm[i] +
354 d*(cosp*np->u[i] + sinp*np->v[i]);
355 d = -2.0 * DOT(h, r->rdir) / (1.0 + d*d);
356 for (i = 0; i < 3; i++)
357 sr.rdir[i] = r->rdir[i] + d*h[i];
358 if (DOT(sr.rdir, r->ron) <= FTINY) /* penetration? */
359 VCOPY(sr.rdir, np->vrefl); /* jitter no good */
360 rayvalue(&sr);
361 multcolor(sr.rcol, np->scolor);
362 addcolor(r->rcol, sr.rcol);
363 ndims--;
364 }
365 /* compute transmission */
366 if ((np->specfl & (SP_TRAN|SP_TBLT)) == SP_TRAN &&
367 rayorigin(&sr, r, SPECULAR, np->tspec) == 0) {
368 dimlist[ndims++] = (int)np->mp;
369 d = urand(ilhash(dimlist,ndims)+1823+samplendx);
370 multisamp(rv, 2, d);
371 d = 2.0*PI * rv[0];
372 cosp = cos(d) * np->u_alpha;
373 sinp = sin(d) * np->v_alpha;
374 d = sqrt(cosp*cosp + sinp*sinp);
375 cosp /= d;
376 sinp /= d;
377 rv[1] = 1.0 - specjitter*rv[1];
378 if (rv[1] <= FTINY)
379 d = 1.0;
380 else
381 d = sqrt(-log(rv[1]) /
382 (cosp*cosp/(np->u_alpha*np->u_alpha) +
383 sinp*sinp/(np->v_alpha*np->u_alpha)));
384 for (i = 0; i < 3; i++)
385 sr.rdir[i] = np->prdir[i] +
386 d*(cosp*np->u[i] + sinp*np->v[i]);
387 if (DOT(sr.rdir, r->ron) < -FTINY)
388 normalize(sr.rdir); /* OK, normalize */
389 else
390 VCOPY(sr.rdir, np->prdir); /* else no jitter */
391 rayvalue(&sr);
392 scalecolor(sr.rcol, np->tspec);
393 multcolor(sr.rcol, np->mcolor); /* modify by color */
394 addcolor(r->rcol, sr.rcol);
395 ndims--;
396 }
397 }