ViewVC Help

View File | Revision Log | Show Annotations | Download File | Root Listing

root/Development/ray/src/gen/mkillum2.c

Comparing ray/src/gen/mkillum2.c (file contents):
Revision 2.25 by greg, Thu Dec 6 19:40:43 2007 UTC vs.
Revision 2.33 by greg, Tue Sep 8 23:05:47 2009 UTC

#	Line 12 | Line 12 | static const char      RCSid[] = "$Id$";
12		#include  "cone.h"
13		#include  "source.h"
14
15	+	#ifndef NBSDFSAMPS
16	+	#define NBSDFSAMPS      256             /* BSDF resampling count */
17	+	#endif
18
19		COLORV *        distarr = NULL;         /* distribution array */
20		int             distsiz = 0;
21	+	COLORV *        direct_discount = NULL; /* amount to take off direct */
22
23
24		void
#	Line 41 | Line 45 | newdist(                       /* allocate & clear distribution array */
45		}
46
47
48	+	static void
49	+	new_discount()                  /* allocate space for direct contrib. record */
50	+	{
51	+	if (distsiz <= 0)
52	+	return;
53	+	direct_discount = (COLORV *)calloc(distsiz, sizeof(COLOR));
54	+	if (direct_discount == NULL)
55	+	error(SYSTEM, "out of memory in new_discount");
56	+	}
57	+
58	+
59	+	static void
60	+	done_discount()                 /* clear off direct contrib. record */
61	+	{
62	+	if (direct_discount == NULL)
63	+	return;
64	+	free((void *)direct_discount);
65	+	direct_discount = NULL;
66	+	}
67	+
68	+
69		int
70		process_ray(                    /* process a ray result or report error */
71		RAY *r,
#	Line 58 | Line 83 | process_ray(                   /* process a ray result or report error
83		multcolor(r->rcol, r->rcoef);   /* in case it's a source ray */
84		colp = &distarr[r->rno * 3];
85		addcolor(colp, r->rcol);
86	+	if (r->rsrc >= 0 &&             /* remember source contrib. */
87	+	direct_discount != NULL) {
88	+	colp = &direct_discount[r->rno * 3];
89	+	addcolor(colp, r->rcol);
90	+	}
91		return(1);
92		}
93
#	Line 127 | Line 157 | srcsamps(                      /* sample sources from this surface positi
157		d = 5.*FTINY;
158		for (i = 3; i--; )
159		sr.rorg[i] += d*nrm[i];
160	+	samplendx++;                    /* increment sample counter */
161		if (!srcray(&sr, NULL, &si))
162		break;                  /* end of sources */
163		/* index direction */
#	Line 139 | Line 170 | srcsamps(                      /* sample sources from this surface positi
170		if (i < 0)
171		continue;       /* must not be important */
172		sr.rno = i;
173	<	d = getBSDF_incrad(il->sd, i);
143	<	d = 1.0/(PI*d*d);
173	>	d = 1.0/getBSDF_incohm(il->sd, i);
174		} else {
175		if (v[2] >= -FTINY)
176		continue;       /* only sample transmission */
177	<	d = 1.0 - v[2]*v[2];
178	<	i = d*nalt;
179	<	d = atan2(-v[1], -v[0])/(2.*PI);
150	<	if (d < 0.0) d += 1.0;
151	<	j = d*nazi + 0.5;
152	<	if (j >= nazi) j = 0;
153	<	sr.rno = i*nazi + j;
154	<	d = nalt*nazi/PI;
177	>	v[0] = -v[0]; v[1] = -v[1]; v[2] = -v[2];
178	>	sr.rno = flatindex(v, nalt, nazi);
179	>	d = nalt*nazi*(1./PI) * v[2];
180		}
181		d *= si.dom;                    /* solid angle correction */
182		scalecolor(sr.rcoef, d);
#	Line 223 | Line 248 | flatdir(               /* compute uniform hemispherical direction *
248		dv[2] = sqrt(1. - alt);
249		}
250
251	+	int
252	+	flatindex(              /* compute index for hemispherical direction */
253	+	FVECT   dv,
254	+	int     nalt,
255	+	int     nazi
256	+	)
257	+	{
258	+	double  d;
259	+	int     i, j;
260	+
261	+	d = 1.0 - dv[2]*dv[2];
262	+	i = d*nalt;
263	+	d = atan2(dv[1], dv[0]) * (0.5/PI);
264	+	if (d < 0.0) d += 1.0;
265	+	j = d*nazi + 0.5;
266	+	if (j >= nazi) j = 0;
267	+	return(i*nazi + j);
268	+	}
269
270	+
271		int
272		my_default(     /* default illum action */
273		OBJREC  *ob,
#	Line 351 | Line 395 | my_face(               /* make an illum face */
395		raysamp(dim[1], org, dir);
396		}
397		/* add in direct component? */
398	<	if (!directvis && il->flags & IL_LIGHT) {
398	>	if (!directvis && (il->flags & IL_LIGHT || il->sd != NULL)) {
399		MAT4    ixfm;
400		if (il->sd == NULL) {
401		for (i = 3; i--; ) {
#	Line 362 | Line 406 | my_face(               /* make an illum face */
406		}
407		ixfm[3][0] = ixfm[3][1] = ixfm[3][2] = 0.;
408		ixfm[3][3] = 1.;
409	<	} else if (!invmat4(ixfm, xfm))
410	<	objerror(ob, INTERNAL, "cannot invert BSDF transform");
409	>	} else {
410	>	if (!invmat4(ixfm, xfm))
411	>	objerror(ob, INTERNAL,
412	>	"cannot invert BSDF transform");
413	>	if (!(il->flags & IL_LIGHT))
414	>	new_discount();
415	>	}
416		dim[0] = random();
417		nallow = 10*il->nsamps;
418		for (i = 0; i < il->nsamps; i++) {
#	Line 393 | Line 442 | my_face(               /* make an illum face */
442		nalt = sqrt(il->sd->nout/PI) + .5;
443		nazi = PI*nalt + .5;
444		redistribute(il->sd, nalt, nazi, u, v, fa->norm, xfm);
445	+	done_discount();
446	+	if (!il->sampdens)
447	+	il->sampdens = nalt*nazi/PI + .999;
448		}
449		/* write out the face and its distribution */
450		if (average(il, distarr, n)) {
#	Line 551 | Line 603 | my_ring(               /* make an illum ring */
603		raysamp(dim[1], org, dir);
604		}
605		/* add in direct component? */
606	<	if (!directvis && il->flags & IL_LIGHT) {
606	>	if (!directvis && (il->flags & IL_LIGHT || il->sd != NULL)) {
607		MAT4    ixfm;
608		if (il->sd == NULL) {
609		for (i = 3; i--; ) {
#	Line 562 | Line 614 | my_ring(               /* make an illum ring */
614		}
615		ixfm[3][0] = ixfm[3][1] = ixfm[3][2] = 0.;
616		ixfm[3][3] = 1.;
617	<	} else if (!invmat4(ixfm, xfm))
618	<	objerror(ob, INTERNAL, "cannot invert BSDF transform");
617	>	} else {
618	>	if (!invmat4(ixfm, xfm))
619	>	objerror(ob, INTERNAL,
620	>	"cannot invert BSDF transform");
621	>	if (!(il->flags & IL_LIGHT))
622	>	new_discount();
623	>	}
624		dim[0] = random();
625		for (i = 0; i < il->nsamps; i++) {
626		/* randomize location */
#	Line 586 | Line 643 | my_ring(               /* make an illum ring */
643		nalt = sqrt(il->sd->nout/PI) + .5;
644		nazi = PI*nalt + .5;
645		redistribute(il->sd, nalt, nazi, u, v, co->ad, xfm);
646	+	done_discount();
647	+	if (!il->sampdens)
648	+	il->sampdens = nalt*nazi/PI + .999;
649		}
650		/* write out the ring and its distribution */
651		if (average(il, distarr, n)) {
#	Line 597 | Line 657 | my_ring(               /* make an illum ring */
657		/* clean up */
658		freecone(ob);
659		return(1);
660	+	}
661	+
662	+
663	+	void
664	+	redistribute(           /* pass distarr ray sums through BSDF */
665	+	struct BSDF_data *b,
666	+	int nalt,
667	+	int nazi,
668	+	FVECT u,
669	+	FVECT v,
670	+	FVECT w,
671	+	MAT4 xm
672	+	)
673	+	{
674	+	int     nout = 0;
675	+	MAT4    mymat, inmat;
676	+	COLORV  *idist;
677	+	COLORV  *cp;
678	+	FVECT   dv;
679	+	double  wt;
680	+	int     i, j, k, c, o;
681	+	COLOR   col, cinc;
682	+	/* copy incoming distribution */
683	+	if (b->ninc > distsiz)
684	+	error(INTERNAL, "error 1 in redistribute");
685	+	idist = (COLORV *)malloc(sizeof(COLOR)*b->ninc);
686	+	if (idist == NULL)
687	+	error(SYSTEM, "out of memory in redistribute");
688	+	memcpy(idist, distarr, sizeof(COLOR)*b->ninc);
689	+	/* compose direction transform */
690	+	for (i = 3; i--; ) {
691	+	mymat[i][0] = u[i];
692	+	mymat[i][1] = v[i];
693	+	mymat[i][2] = w[i];
694	+	mymat[i][3] = 0.;
695	+	}
696	+	mymat[3][0] = mymat[3][1] = mymat[3][2] = 0.;
697	+	mymat[3][3] = 1.;
698	+	if (xm != NULL)
699	+	multmat4(mymat, xm, mymat);
700	+	for (i = 3; i--; ) {            /* make sure it's normalized */
701	+	wt = 1./sqrt(   mymat[0][i]*mymat[0][i] +
702	+	mymat[1][i]*mymat[1][i] +
703	+	mymat[2][i]*mymat[2][i] );
704	+	for (j = 3; j--; )
705	+	mymat[j][i] *= wt;
706	+	}
707	+	if (!invmat4(inmat, mymat))     /* need inverse as well */
708	+	error(INTERNAL, "cannot invert BSDF transform");
709	+	newdist(nalt*nazi);             /* resample distribution */
710	+	for (i = b->ninc; i--; ) {
711	+	int     direct_out = -1;
712	+	COLOR   cdir;
713	+	getBSDF_incvec(dv, b, i);       /* compute incident irrad. */
714	+	multv3(dv, dv, mymat);
715	+	if (dv[2] < 0.0) {
716	+	dv[0] = -dv[0]; dv[1] = -dv[1]; dv[2] = -dv[2];
717	+	direct_out += (direct_discount != NULL);
718	+	}
719	+	wt = getBSDF_incohm(b, i);
720	+	wt *= dv[2];                    /* solid_angle*cosine(theta) */
721	+	cp = &idist[3*i];
722	+	copycolor(cinc, cp);
723	+	scalecolor(cinc, wt);
724	+	if (!direct_out) {              /* discount direct contr. */
725	+	cp = &direct_discount[3*i];
726	+	copycolor(cdir, cp);
727	+	scalecolor(cdir, -wt);
728	+	if (b->nout != b->ninc)
729	+	direct_out = flatindex(dv, nalt, nazi);
730	+	else
731	+	direct_out = i; /* assumes dist. mirroring */
732	+	}
733	+	for (k = nalt; k--; )           /* loop over distribution */
734	+	for (j = nazi; j--; ) {
735	+	int rstart = random();
736	+	for (c = NBSDFSAMPS; c--; ) {
737	+	double  sp[2];
738	+	multisamp(sp, 2, urand(rstart+c));
739	+	flatdir(dv, (k + sp[0])/nalt,
740	+	(j + .5 - sp[1])/nazi);
741	+	multv3(dv, dv, inmat);
742	+	/* evaluate BSDF @ outgoing */
743	+	o = getBSDF_outndx(b, dv);
744	+	if (o < 0) {
745	+	nout++;
746	+	continue;
747	+	}
748	+	wt = BSDF_value(b, i, o) * (1./NBSDFSAMPS);
749	+	copycolor(col, cinc);
750	+	if (b->nout != b->ninc)
751	+	o = k*nazi + j;
752	+	if (o == direct_out)
753	+	addcolor(col, cdir);    /* minus direct */
754	+	scalecolor(col, wt);
755	+	cp = &distarr[3*o];
756	+	addcolor(cp, col);      /* sum into distribution */
757	+	}
758	+	}
759	+	}
760	+	free(idist);                    /* free temp space */
761	+	if (nout) {
762	+	sprintf(errmsg, "missing %.1f%% of BSDF directions",
763	+	100.*nout/(b->ninc*nalt*nazi*NBSDFSAMPS));
764	+	error(WARNING, errmsg);
765	+	}
766		}

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines (old)
>	Changed lines (new)