--- ray/src/gen/mkillum2.c 2007/09/30 19:05:26 2.22 +++ ray/src/gen/mkillum2.c 2007/12/08 01:43:09 2.26 @@ -1,5 +1,5 @@ #ifndef lint -static const char RCSid[] = "$Id: mkillum2.c,v 2.22 2007/09/30 19:05:26 greg Exp $"; +static const char RCSid[] = "$Id: mkillum2.c,v 2.26 2007/12/08 01:43:09 greg Exp $"; #endif /* * Routines to do the actual calculation for mkillum @@ -10,7 +10,7 @@ static const char RCSid[] = "$Id: mkillum2.c,v 2.22 20 #include "mkillum.h" #include "face.h" #include "cone.h" -#include "random.h" +#include "source.h" COLORV * distarr = NULL; /* distribution array */ @@ -22,7 +22,7 @@ newdist( /* allocate & clear distribution array */ int siz ) { - if (siz == 0) { + if (siz <= 0) { if (distsiz > 0) free((void *)distarr); distarr = NULL; @@ -42,16 +42,20 @@ newdist( /* allocate & clear distribution array */ int -process_ray(RAY *r, int rv) +process_ray( /* process a ray result or report error */ + RAY *r, + int rv +) { COLORV *colp; - if (rv == 0) + if (rv == 0) /* no result ready */ return(0); if (rv < 0) error(USER, "ray tracing process died"); if (r->rno >= distsiz) error(INTERNAL, "bad returned index in process_ray"); + multcolor(r->rcol, r->rcoef); /* in case it's a source ray */ colp = &distarr[r->rno * 3]; addcolor(colp, r->rcol); return(1); @@ -59,7 +63,7 @@ process_ray(RAY *r, int rv) void -raysamp( /* queue a ray sample */ +raysamp( /* queue a ray sample */ int ndx, FVECT org, FVECT dir @@ -81,6 +85,81 @@ raysamp( /* queue a ray sample */ void +srcsamps( /* sample sources from this surface position */ + struct illum_args *il, + FVECT org, + FVECT nrm, + MAT4 ixfm +) +{ + int nalt, nazi; + SRCINDEX si; + RAY sr; + FVECT v; + double d; + int i, j; + /* get sampling density */ + if (il->sampdens <= 0) { + nalt = nazi = 1; + } else { + i = PI * il->sampdens; + nalt = sqrt(i/PI) + .5; + nazi = PI*nalt + .5; + } + 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]; + if (!srcray(&sr, NULL, &si)) + break; /* end of sources */ + /* index direction */ + if (ixfm != NULL) + 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; + } + d *= si.dom; /* solid angle correction */ + scalecolor(sr.rcoef, d); + process_ray(&sr, ray_pqueue(&sr)); + } +} + + +void rayclean() /* finish all pending rays */ { RAY myRay; @@ -127,9 +206,9 @@ rounddir( /* compute uniform spherical direction */ } -static void +void flatdir( /* compute uniform hemispherical direction */ - register FVECT dv, + FVECT dv, double alt, double azi ) @@ -166,17 +245,17 @@ my_face( /* make an illum face */ char *nm ) { -#define MAXMISS (5*n*il->nsamps) int dim[2]; - int n, nalt, nazi, h, alti; + int n, nalt, nazi, alti; double sp[2], r1, r2; + int h; FVECT dn, org, dir; FVECT u, v; double ur[2], vr[2]; MAT4 xfm; - int nmisses; + int nallow; FACE *fa; - register int i, j; + int i, j; /* get/check arguments */ fa = getface(ob); if (fa->area == 0.0) { @@ -229,10 +308,10 @@ my_face( /* make an illum face */ } dim[0] = random(); /* sample polygon */ - nmisses = 0; + nallow = 5*n*il->nsamps; for (dim[1] = 0; dim[1] < n; dim[1]++) for (i = 0; i < il->nsamps; i++) { - /* random direction */ + /* randomize direction */ h = ilhash(dim, 2) + i; if (il->sd != NULL) { r_BSDF_incvec(dir, il->sd, dim[1], urand(h), xfm); @@ -246,30 +325,68 @@ my_face( /* make an illum face */ dir[j] = -dn[0]*u[j] - dn[1]*v[j] - dn[2]*fa->norm[j]; } - /* random location */ + /* randomize location */ do { - multisamp(sp, 2, urand(h+4862+nmisses)); + multisamp(sp, 2, urand(h+4862+nallow)); r1 = ur[0] + (ur[1]-ur[0]) * sp[0]; r2 = vr[0] + (vr[1]-vr[0]) * sp[1]; for (j = 0; j < 3; j++) org[j] = r1*u[j] + r2*v[j] + fa->offset*fa->norm[j]; - } while (!inface(org, fa) && nmisses++ < MAXMISS); - if (nmisses > MAXMISS) { + } while (!inface(org, fa) && nallow-- > 0); + if (nallow < 0) { objerror(ob, WARNING, "bad aspect"); rayclean(); freeface(ob); return(my_default(ob, il, nm)); } if (il->sd != NULL && DOT(dir, fa->norm) < -FTINY) - r1 = -1.0001*il->thick - .0001; + r1 = -1.0001*il->thick - 5.*FTINY; else - r1 = .0001; + 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) { + 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"); + dim[0] = random(); + nallow = 10*il->nsamps; + for (i = 0; i < il->nsamps; i++) { + /* randomize location */ + h = dim[0] + samplendx++; + do { + multisamp(sp, 2, urand(h+nallow)); + r1 = ur[0] + (ur[1]-ur[0]) * sp[0]; + r2 = vr[0] + (vr[1]-vr[0]) * sp[1]; + for (j = 0; j < 3; j++) + org[j] = r1*u[j] + r2*v[j] + + fa->offset*fa->norm[j]; + } while (!inface(org, fa) && nallow-- > 0); + if (nallow < 0) { + objerror(ob, WARNING, "bad aspect"); + rayclean(); + freeface(ob); + return(my_default(ob, il, nm)); + } + /* sample source rays */ + srcsamps(il, org, fa->norm, ixfm); + } + } + /* wait for all rays to finish */ rayclean(); if (il->sd != NULL) { /* run distribution through BSDF */ nalt = sqrt(il->sd->nout/PI) + .5; @@ -286,7 +403,6 @@ my_face( /* make an illum face */ /* clean up */ freeface(ob); return(0); -#undef MAXMISS } @@ -344,6 +460,7 @@ my_sphere( /* make an illum sphere */ /* send sample */ raysamp(dim[1]*nazi+dim[2], org, dir); } + /* wait for all rays to finish */ rayclean(); /* write out the sphere and its distribution */ if (average(il, distarr, n)) { @@ -374,7 +491,7 @@ my_ring( /* make an illum ring */ FVECT u, v; MAT4 xfm; CONE *co; - register int i, j; + int i, j; /* get/check arguments */ co = getcone(ob, 0); /* set up sampling */ @@ -403,7 +520,7 @@ my_ring( /* make an illum ring */ for (i = 0; i < il->nsamps; i++) { /* next sample point */ h = ilhash(dim,2) + i; - /* random direction */ + /* randomize direction */ if (il->sd != NULL) { r_BSDF_incvec(dir, il->sd, dim[1], urand(h), xfm); } else { @@ -413,9 +530,9 @@ my_ring( /* make an illum ring */ 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]; + dir[j] = -dn[0]*u[j] - dn[1]*v[j] - dn[2]*co->ad[j]; } - /* random location */ + /* randomize location */ multisamp(sp, 2, urand(h+8371)); r3 = sqrt(CO_R0(co)*CO_R0(co) + sp[0]*(CO_R1(co)*CO_R1(co) - CO_R0(co)*CO_R0(co))); @@ -423,15 +540,46 @@ my_ring( /* make an illum ring */ r1 = r3*cos(r2); r2 = r3*sin(r2); if (il->sd != NULL && DOT(dir, co->ad) < -FTINY) - r3 = -1.0001*il->thick - .0001; + r3 = -1.0001*il->thick - 5.*FTINY; else - r3 = .0001; + 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]; /* send sample */ raysamp(dim[1], org, dir); } + /* add in direct component? */ + if (!directvis && 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"); + dim[0] = random(); + for (i = 0; i < il->nsamps; i++) { + /* randomize location */ + h = dim[0] + samplendx++; + multisamp(sp, 2, urand(h)); + r3 = sqrt(CO_R0(co)*CO_R0(co) + + sp[0]*(CO_R1(co)*CO_R1(co) - CO_R0(co)*CO_R0(co))); + r2 = 2.*PI*sp[1]; + r1 = r3*cos(r2); + r2 = r3*sin(r2); + for (j = 0; j < 3; j++) + org[j] = CO_P0(co)[j] + r1*u[j] + r2*v[j]; + /* sample source rays */ + srcsamps(il, org, co->ad, ixfm); + } + } + /* wait for all rays to finish */ rayclean(); if (il->sd != NULL) { /* run distribution through BSDF */ nalt = sqrt(il->sd->nout/PI) + .5;