--- ray/src/gen/mkillum2.c 2007/12/08 01:43:09 2.26 +++ ray/src/gen/mkillum2.c 2012/10/13 20:15:43 2.38 @@ -1,5 +1,5 @@ #ifndef lint -static const char RCSid[] = "$Id: mkillum2.c,v 2.26 2007/12/08 01:43:09 greg Exp $"; +static const char RCSid[] = "$Id: mkillum2.c,v 2.38 2012/10/13 20:15:43 greg Exp $"; #endif /* * Routines to do the actual calculation for mkillum @@ -11,7 +11,11 @@ static const char RCSid[] = "$Id: mkillum2.c,v 2.26 20 #include "face.h" #include "cone.h" #include "source.h" +#include "paths.h" +#ifndef NBSDFSAMPS +#define NBSDFSAMPS 256 /* BSDF resampling count */ +#endif COLORV * distarr = NULL; /* distribution array */ int distsiz = 0; @@ -24,14 +28,14 @@ newdist( /* allocate & clear distribution array */ { if (siz <= 0) { if (distsiz > 0) - free((void *)distarr); + free(distarr); distarr = NULL; distsiz = 0; return; } if (distsiz < siz) { if (distsiz > 0) - free((void *)distarr); + free(distarr); distarr = (COLORV *)malloc(sizeof(COLOR)*siz); if (distarr == NULL) error(SYSTEM, "out of memory in newdist"); @@ -77,7 +81,7 @@ raysamp( /* queue a ray sample */ VCOPY(myRay.rorg, org); VCOPY(myRay.rdir, dir); myRay.rmax = .0; - rayorigin(&myRay, PRIMARY, NULL, NULL); + rayorigin(&myRay, PRIMARY|SPECULAR, NULL, NULL); myRay.rno = ndx; /* queue ray, check result */ process_ray(&myRay, ray_pqueue(&myRay)); @@ -92,16 +96,14 @@ srcsamps( /* sample sources from this surface positi MAT4 ixfm ) { - int nalt, nazi; + int nalt=1, nazi=1; SRCINDEX si; RAY sr; FVECT v; double d; int i, j; /* get sampling density */ - if (il->sampdens <= 0) { - nalt = nazi = 1; - } else { + if (il->sampdens > 0) { i = PI * il->sampdens; nalt = sqrt(i/PI) + .5; nazi = PI*nalt + .5; @@ -109,24 +111,9 @@ srcsamps( /* sample sources from this surface positi initsrcindex(&si); /* loop over (sub)sources */ for ( ; ; ) { VCOPY(sr.rorg, org); /* pick side to shoot from */ - if (il->sd != NULL) { - int sn = si.sn; - if (si.sp+1 >= si.np) ++sn; - if (sn >= nsources) break; - if (source[sn].sflags & SDISTANT) - d = DOT(source[sn].sloc, nrm); - else { - VSUB(v, source[sn].sloc, org); - d = DOT(v, nrm); - } - } else - d = 1.0; /* only transmission */ - if (d < 0.0) - d = -1.0001*il->thick - 5.*FTINY; - else - d = 5.*FTINY; - for (i = 3; i--; ) - sr.rorg[i] += d*nrm[i]; + d = 5.*FTINY; + VSUM(sr.rorg, sr.rorg, nrm, d); + samplendx++; /* increment sample counter */ if (!srcray(&sr, NULL, &si)) break; /* end of sources */ /* index direction */ @@ -134,24 +121,11 @@ srcsamps( /* sample sources from this surface positi multv3(v, sr.rdir, ixfm); else VCOPY(v, sr.rdir); - if (il->sd != NULL) { - i = getBSDF_incndx(il->sd, v); - if (i < 0) - continue; /* must not be important */ - sr.rno = i; - d = 1.0/getBSDF_incohm(il->sd, i); - } else { - if (v[2] >= -FTINY) - continue; /* only sample transmission */ - d = 1.0 - v[2]*v[2]; - i = d*nalt; - d = atan2(-v[1], -v[0])/(2.*PI); - if (d < 0.0) d += 1.0; - j = d*nazi + 0.5; - if (j >= nazi) j = 0; - sr.rno = i*nazi + j; - d = nalt*nazi/PI; - } + if (v[2] >= -FTINY) + continue; /* only sample transmission */ + v[0] = -v[0]; v[1] = -v[1]; v[2] = -v[2]; + sr.rno = flatindex(v, nalt, nazi); + d = nalt*nazi*(1./PI) * v[2]; d *= si.dom; /* solid angle correction */ scalecolor(sr.rcoef, d); process_ray(&sr, ray_pqueue(&sr)); @@ -224,6 +198,26 @@ flatdir( /* compute uniform hemispherical direction * int +flatindex( /* compute index for hemispherical direction */ + FVECT dv, + int nalt, + int nazi +) +{ + double d; + int i, j; + + d = 1.0 - dv[2]*dv[2]; + i = d*nalt; + d = atan2(dv[1], dv[0]) * (0.5/PI); + if (d < 0.0) d += 1.0; + j = d*nazi + 0.5; + if (j >= nazi) j = 0; + return(i*nazi + j); +} + + +int my_default( /* default illum action */ OBJREC *ob, struct illum_args *il, @@ -253,6 +247,7 @@ my_face( /* make an illum face */ FVECT u, v; double ur[2], vr[2]; MAT4 xfm; + char xfrot[64]; int nallow; FACE *fa; int i, j; @@ -263,23 +258,14 @@ my_face( /* make an illum face */ return(my_default(ob, il, nm)); } /* set up sampling */ - if (il->sd != NULL) { - if (!getBSDF_xfm(xfm, fa->norm, il->udir)) { - objerror(ob, WARNING, "illegal up direction"); - freeface(ob); - return(my_default(ob, il, nm)); - } - n = il->sd->ninc; + if (il->sampdens <= 0) { + nalt = nazi = 1; /* diffuse assumption */ } else { - if (il->sampdens <= 0) { - nalt = nazi = 1; /* diffuse assumption */ - } else { - n = PI * il->sampdens; - nalt = sqrt(n/PI) + .5; - nazi = PI*nalt + .5; - } - n = nazi*nalt; + n = PI * il->sampdens; + nalt = sqrt(n/PI) + .5; + nazi = PI*nalt + .5; } + n = nazi*nalt; newdist(n); /* take first edge >= sqrt(area) */ for (j = fa->nv-1, i = 0; i < fa->nv; j = i++) { @@ -313,18 +299,14 @@ my_face( /* make an illum face */ for (i = 0; i < il->nsamps; i++) { /* randomize direction */ h = ilhash(dim, 2) + i; - if (il->sd != NULL) { - r_BSDF_incvec(dir, il->sd, dim[1], urand(h), xfm); - } else { - multisamp(sp, 2, urand(h)); - alti = dim[1]/nazi; - r1 = (alti + sp[0])/nalt; - r2 = (dim[1] - alti*nazi + sp[1] - .5)/nazi; - flatdir(dn, r1, r2); - for (j = 0; j < 3; j++) + multisamp(sp, 2, urand(h)); + alti = dim[1]/nazi; + r1 = (alti + sp[0])/nalt; + r2 = (dim[1] - alti*nazi + sp[1] - .5)/nazi; + flatdir(dn, r1, r2); + for (j = 0; j < 3; j++) dir[j] = -dn[0]*u[j] - dn[1]*v[j] - dn[2]*fa->norm[j]; - } /* randomize location */ do { multisamp(sp, 2, urand(h+4862+nallow)); @@ -340,29 +322,23 @@ my_face( /* make an illum face */ freeface(ob); return(my_default(ob, il, nm)); } - if (il->sd != NULL && DOT(dir, fa->norm) < -FTINY) - r1 = -1.0001*il->thick - 5.*FTINY; - else - r1 = 5.*FTINY; + r1 = 5.*FTINY; for (j = 0; j < 3; j++) org[j] += r1*fa->norm[j]; /* send sample */ raysamp(dim[1], org, dir); } /* add in direct component? */ - if (!directvis && il->flags & IL_LIGHT) { + if (il->flags & IL_LIGHT) { MAT4 ixfm; - if (il->sd == NULL) { - for (i = 3; i--; ) { - ixfm[i][0] = u[i]; - ixfm[i][1] = v[i]; - ixfm[i][2] = fa->norm[i]; - ixfm[i][3] = 0.; - } - ixfm[3][0] = ixfm[3][1] = ixfm[3][2] = 0.; - ixfm[3][3] = 1.; - } else if (!invmat4(ixfm, xfm)) - objerror(ob, INTERNAL, "cannot invert BSDF transform"); + for (i = 3; i--; ) { + ixfm[i][0] = u[i]; + ixfm[i][1] = v[i]; + ixfm[i][2] = fa->norm[i]; + ixfm[i][3] = 0.; + } + ixfm[3][0] = ixfm[3][1] = ixfm[3][2] = 0.; + ixfm[3][3] = 1.; dim[0] = random(); nallow = 10*il->nsamps; for (i = 0; i < il->nsamps; i++) { @@ -388,11 +364,6 @@ my_face( /* make an illum face */ } /* wait for all rays to finish */ rayclean(); - if (il->sd != NULL) { /* run distribution through BSDF */ - nalt = sqrt(il->sd->nout/PI) + .5; - nazi = PI*nalt + .5; - redistribute(il->sd, nalt, nazi, u, v, fa->norm, xfm); - } /* write out the face and its distribution */ if (average(il, distarr, n)) { if (il->sampdens > 0) @@ -430,8 +401,6 @@ my_sphere( /* make an illum sphere */ nalt = sqrt(2./PI*n) + .5; nazi = PI/2.*nalt + .5; } - if (il->sd != NULL) - objerror(ob, WARNING, "BSDF ignored"); n = nalt*nazi; newdist(n); dim[0] = random(); @@ -495,23 +464,14 @@ my_ring( /* make an illum ring */ /* get/check arguments */ co = getcone(ob, 0); /* set up sampling */ - if (il->sd != NULL) { - if (!getBSDF_xfm(xfm, co->ad, il->udir)) { - objerror(ob, WARNING, "illegal up direction"); - freecone(ob); - return(my_default(ob, il, nm)); - } - n = il->sd->ninc; + if (il->sampdens <= 0) { + nalt = nazi = 1; /* diffuse assumption */ } else { - if (il->sampdens <= 0) { - nalt = nazi = 1; /* diffuse assumption */ - } else { - n = PI * il->sampdens; - nalt = sqrt(n/PI) + .5; - nazi = PI*nalt + .5; - } - n = nazi*nalt; + n = PI * il->sampdens; + nalt = sqrt(n/PI) + .5; + nazi = PI*nalt + .5; } + n = nazi*nalt; newdist(n); mkaxes(u, v, co->ad); dim[0] = random(); @@ -521,17 +481,13 @@ my_ring( /* make an illum ring */ /* next sample point */ h = ilhash(dim,2) + i; /* randomize direction */ - if (il->sd != NULL) { - r_BSDF_incvec(dir, il->sd, dim[1], urand(h), xfm); - } else { - multisamp(sp, 2, urand(h)); - alti = dim[1]/nazi; - r1 = (alti + sp[0])/nalt; - r2 = (dim[1] - alti*nazi + sp[1] - .5)/nazi; - flatdir(dn, r1, r2); - for (j = 0; j < 3; j++) - dir[j] = -dn[0]*u[j] - dn[1]*v[j] - dn[2]*co->ad[j]; - } + multisamp(sp, 2, urand(h)); + alti = dim[1]/nazi; + r1 = (alti + sp[0])/nalt; + r2 = (dim[1] - alti*nazi + sp[1] - .5)/nazi; + flatdir(dn, r1, r2); + for (j = 0; j < 3; j++) + dir[j] = -dn[0]*u[j] - dn[1]*v[j] - dn[2]*co->ad[j]; /* randomize location */ multisamp(sp, 2, urand(h+8371)); r3 = sqrt(CO_R0(co)*CO_R0(co) + @@ -539,10 +495,7 @@ my_ring( /* make an illum ring */ r2 = 2.*PI*sp[1]; r1 = r3*cos(r2); r2 = r3*sin(r2); - if (il->sd != NULL && DOT(dir, co->ad) < -FTINY) - r3 = -1.0001*il->thick - 5.*FTINY; - else - r3 = 5.*FTINY; + r3 = 5.*FTINY; for (j = 0; j < 3; j++) org[j] = CO_P0(co)[j] + r1*u[j] + r2*v[j] + r3*co->ad[j]; @@ -550,19 +503,16 @@ my_ring( /* make an illum ring */ raysamp(dim[1], org, dir); } /* add in direct component? */ - if (!directvis && il->flags & IL_LIGHT) { + if (il->flags & IL_LIGHT) { MAT4 ixfm; - if (il->sd == NULL) { - for (i = 3; i--; ) { - ixfm[i][0] = u[i]; - ixfm[i][1] = v[i]; - ixfm[i][2] = co->ad[i]; - ixfm[i][3] = 0.; - } - ixfm[3][0] = ixfm[3][1] = ixfm[3][2] = 0.; - ixfm[3][3] = 1.; - } else if (!invmat4(ixfm, xfm)) - objerror(ob, INTERNAL, "cannot invert BSDF transform"); + for (i = 3; i--; ) { + ixfm[i][0] = u[i]; + ixfm[i][1] = v[i]; + ixfm[i][2] = co->ad[i]; + ixfm[i][3] = 0.; + } + ixfm[3][0] = ixfm[3][1] = ixfm[3][2] = 0.; + ixfm[3][3] = 1.; dim[0] = random(); for (i = 0; i < il->nsamps; i++) { /* randomize location */ @@ -581,11 +531,6 @@ my_ring( /* make an illum ring */ } /* wait for all rays to finish */ rayclean(); - if (il->sd != NULL) { /* run distribution through BSDF */ - nalt = sqrt(il->sd->nout/PI) + .5; - nazi = PI*nalt + .5; - redistribute(il->sd, nalt, nazi, u, v, co->ad, xfm); - } /* write out the ring and its distribution */ if (average(il, distarr, n)) { if (il->sampdens > 0)