ViewVC Help

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

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

Comparing ray/src/gen/mksource.c (file contents):
Revision 2.4 by greg, Fri Sep 23 19:22:37 2005 UTC vs.
Revision 2.5 by greg, Fri Aug 24 04:41:36 2007 UTC

#	Line 17 | Line 17 | static const char RCSid[] = "$Id$";
17		/* Data structure for geodesic samples */
18
19		typedef struct tritree {
20	<	FVECT           gdv[3];         /* spherical triangle vertex direc. */
21	<	FVECT           sd;             /* sample direction if leaf */
20	>	int32           gdv[3];         /* spherical triangle vertex direc. */
21	>	int32           sd;             /* sample direction if leaf */
22		struct tritree  *kid;           /* 4 children if branch node */
23		COLR            val;            /* sampled color value */
24		} TRITREE;
25
26		typedef struct lostlight {
27		struct lostlight        *next;  /* next in list */
28	<	FVECT           sd;             /* lost source direction */
28	>	int32           sd;             /* lost source direction */
29		COLOR           intens;         /* output times solid angle */
30		} LOSTLIGHT;
31
#	Line 56 | Line 56 | vol_sign(const FVECT v1, const FVECT v2, const FVECT v
56
57		/* Is the given direction contained within the specified spherical triangle? */
58		int
59	<	intriv(FVECT tri[3], const FVECT sdir)
59	>	intriv(const int32 trid[3], const FVECT sdir)
60		{
61		int     sv[3];
62	+	FVECT   tri[3];
63
64	+	decodedir(tri[0], trid[0]);
65	+	decodedir(tri[1], trid[1]);
66	+	decodedir(tri[2], trid[2]);
67		sv[0] = vol_sign(sdir, tri[0], tri[1]);
68		sv[1] = vol_sign(sdir, tri[1], tri[2]);
69		sv[2] = vol_sign(sdir, tri[2], tri[0]);
#	Line 91 | Line 95 | leafsample(TRITREE *leaf)
95		{
96		RAY     myray;
97		RREAL   wt[3];
98	+	FVECT   sdir;
99		int     i, j;
100		/* random point on triangle */
101		i = random() % 3;
#	Line 98 | Line 103 | leafsample(TRITREE *leaf)
103		j = random() & 1;
104		wt[(i+2-j)%3] = 1. - wt[i] -
105		(wt[(i+1+j)%3] = (1.-wt[i])*frandom());
106	<	leaf->sd[0] = leaf->sd[1] = leaf->sd[2] = .0;
107	<	for (i = 0; i < 3; i++)
108	<	VSUM(leaf->sd, leaf->sd, leaf->gdv[i], wt[i]);
109	<	normalize(leaf->sd);                    /* record sample direction */
106	>	sdir[0] = sdir[1] = sdir[2] = .0;
107	>	for (i = 0; i < 3; i++) {
108	>	FVECT   vt;
109	>	decodedir(vt, leaf->gdv[i]);
110	>	VSUM(sdir, sdir, vt, wt[i]);
111	>	}
112	>	normalize(sdir);                        /* record sample direction */
113	>	leaf->sd = encodedir(sdir);
114		/* evaluate at inf. */
115	<	VSUM(myray.rorg, scene_cent, leaf->sd, scene_rad);
116	<	VCOPY(myray.rdir, leaf->sd);
115	>	VSUM(myray.rorg, scene_cent, sdir, scene_rad);
116	>	VCOPY(myray.rdir, sdir);
117		myray.rmax = 0.;
118		ray_trace(&myray);
119		setcolr(leaf->val, colval(myray.rcol,RED),
#	Line 114 | Line 123 | leafsample(TRITREE *leaf)
123
124		/* Initialize a branch node contained in the given spherical triangle */
125		void
126	<	subdivide(TRITREE *branch, FVECT dv[3])
126	>	subdivide(TRITREE *branch, const int32 dv[3])
127		{
128	<	FVECT   sdv[3];
128	>	FVECT   dvv[3], sdv[3];
129	>	int32   sd[3];
130		int     i;
131
132	<	for (i = 0; i < 3; i++)         /* copy spherical triangle */
133	<	VCOPY(branch->gdv[i], dv[i]);
132	>	for (i = 0; i < 3; i++) {       /* copy spherical triangle */
133	>	branch->gdv[i] = dv[i];
134	>	decodedir(dvv[i], dv[i]);
135	>	}
136		for (i = 0; i < 3; i++) {       /* create new vertices */
137		int     j = (i+1)%3;
138	<	VADD(sdv[i], dv[i], dv[j]);
138	>	VADD(sdv[i], dvv[i], dvv[j]);
139		normalize(sdv[i]);
140	+	sd[i] = encodedir(sdv[i]);
141		}
142		/* allocate leaves */
143		branch->kid = (TRITREE *)calloc(4, sizeof(TRITREE));
144		if (branch->kid == NULL)
145		error(SYSTEM, "out of memory in subdivide()");
146		/* assign subtriangle directions */
147	<	VCOPY(branch->kid[0].gdv[0], dv[0]);
148	<	VCOPY(branch->kid[0].gdv[1], sdv[0]);
149	<	VCOPY(branch->kid[0].gdv[2], sdv[2]);
150	<	VCOPY(branch->kid[1].gdv[0], sdv[0]);
151	<	VCOPY(branch->kid[1].gdv[1], dv[1]);
152	<	VCOPY(branch->kid[1].gdv[2], sdv[1]);
153	<	VCOPY(branch->kid[2].gdv[0], sdv[1]);
154	<	VCOPY(branch->kid[2].gdv[1], dv[2]);
155	<	VCOPY(branch->kid[2].gdv[2], sdv[2]);
156	<	VCOPY(branch->kid[3].gdv[0], sdv[0]);
157	<	VCOPY(branch->kid[3].gdv[1], sdv[1]);
158	<	VCOPY(branch->kid[3].gdv[2], sdv[2]);
147	>	branch->kid[0].gdv[0] = dv[0];
148	>	branch->kid[0].gdv[1] = sd[0];
149	>	branch->kid[0].gdv[2] = sd[2];
150	>	branch->kid[1].gdv[0] = sd[0];
151	>	branch->kid[1].gdv[1] = dv[1];
152	>	branch->kid[1].gdv[2] = sd[1];
153	>	branch->kid[2].gdv[0] = sd[1];
154	>	branch->kid[2].gdv[1] = dv[2];
155	>	branch->kid[2].gdv[2] = sd[2];
156	>	branch->kid[3].gdv[0] = sd[0];
157	>	branch->kid[3].gdv[1] = sd[1];
158	>	branch->kid[3].gdv[2] = sd[2];
159		}
160
161		/* Recursively subdivide the given node to the specified quadtree depth */
#	Line 155 | Line 168 | branchsample(TRITREE *node, int depth)
168		return;
169		if (isleaf(node)) {                     /* subdivide leaf node */
170		TRITREE branch, *moved_leaf;
171	+	FVECT   sdir;
172		subdivide(&branch, node->gdv);
173	<	moved_leaf = findleaf(&branch, node->sd);
173	>	decodedir(sdir, node->sd);
174	>	moved_leaf = findleaf(&branch, sdir);
175		if (moved_leaf != NULL) {       /* bequeath old sample */
176	<	VCOPY(moved_leaf->sd, node->sd);
176	>	moved_leaf->sd = node->sd;
177		copycolr(moved_leaf->val, node->val);
178		}
179		for (i = 0; i < 4; i++)         /* compute new samples */
#	Line 193 | Line 208 | geosample(int nsamps)
208		trunk[1][1] = sqrt(1. - trunk[1][2]*trunk[1][2]);
209		spinvector(trunk[2], trunk[1], trunk[0], 2.*PI/3.);
210		spinvector(trunk[3], trunk[1], trunk[0], 4.*PI/3.);
211	<	VCOPY(tree->gdv[0], trunk[0]);
197	<	VCOPY(tree->gdv[1], trunk[0]);
198	<	VCOPY(tree->gdv[2], trunk[0]);
211	>	tree->gdv[0] = tree->gdv[1] = tree->gdv[2] = encodedir(trunk[0]);
212		tree->kid = (TRITREE *)calloc(NTRUNKBR, sizeof(TRITREE));
213		if (tree->kid == NULL) goto memerr;
214		/* grow our tree from trunk */
215		for (i = 0; i < NTRUNKBR; i++) {
216		for (j = 0; j < 3; j++)         /* XXX works for tetra only */
217	<	VCOPY(tree->kid[i].gdv[j], trunk[(i+j)%NTRUNKVERT]);
217	>	tree->kid[i].gdv[j] = encodedir(trunk[(i+j)%NTRUNKVERT]);
218		leafsample(&tree->kid[i]);
219		branchsample(&tree->kid[i], depth);
220		}
#	Line 269 | Line 282 | findemax(TRITREE *node, int *expp)
282
283		/* Compute solid angle of spherical triangle (approx.) */
284		double
285	<	tri_omegav(FVECT v[3])
285	>	tri_omegav(const int32 vd[3])
286		{
287	<	FVECT   e1, e2, vcross;
288	<
287	>	FVECT   v[3], e1, e2, vcross;
288	>
289	>	decodedir(v[0], vd[0]);
290	>	decodedir(v[1], vd[1]);
291	>	decodedir(v[2], vd[2]);
292		VSUB(e1, v[1], v[0]);
293		VSUB(e2, v[2], v[1]);
294		fcross(vcross, e1, e2);
295		return(.5*VLEN(vcross));
296		}
297
298	<	/* Sum intensity times direction for non-zero leaves */
298	>	/* Sum intensity times direction for above-threshold perimiter within radius */
299		void
300		vector_sum(FVECT vsum, TRITREE *node,
301	<	const FVECT cent, double mincos, int ethresh)
301	>	FVECT cent, double maxr2, int ethresh)
302		{
303		if (isleaf(node)) {
304		double  intens;
305	+	FVECT   sdir;
306		if (node->val[EXP] < ethresh)
307	<	return;
308	<	if (DOT(node->sd,cent) < mincos)
309	<	return;
307	>	return;                         /* below threshold */
308	>	if (fdir2diff(node->sd,cent) > maxr2)
309	>	return;                         /* too far away */
310		intens = colrval(node->val,GRN) * tri_omegav(node->gdv);
311	<	VSUM(vsum, vsum, node->sd, intens);
311	>	decodedir(sdir, node->sd);
312	>	VSUM(vsum, vsum, sdir, intens);
313		return;
314		}
315	<	if (DOT(node->gdv[0],node->gdv[1]) < mincos &&
316	<	DOT(node->gdv[0],cent) > mincos &&
317	<	DOT(node->gdv[1],cent) > mincos &&
318	<	DOT(node->gdv[2],cent) > mincos)
319	<	return;
320	<	vector_sum(vsum, &node->kid[0], cent, mincos, ethresh);
321	<	vector_sum(vsum, &node->kid[1], cent, mincos, ethresh);
322	<	vector_sum(vsum, &node->kid[2], cent, mincos, ethresh);
323	<	vector_sum(vsum, &node->kid[3], cent, mincos, ethresh);
315	>	if (dir2diff(node->gdv[0],node->gdv[1]) > maxr2 &&
316	>	fdir2diff(node->gdv[0],cent) < maxr2 &&
317	>	fdir2diff(node->gdv[1],cent) < maxr2 &&
318	>	fdir2diff(node->gdv[2],cent) < maxr2)
319	>	return;                                 /* containing node */
320	>	vector_sum(vsum, &node->kid[0], cent, maxr2, ethresh);
321	>	vector_sum(vsum, &node->kid[1], cent, maxr2, ethresh);
322	>	vector_sum(vsum, &node->kid[2], cent, maxr2, ethresh);
323	>	vector_sum(vsum, &node->kid[3], cent, maxr2, ethresh);
324		}
325
326		/* Claim source contributions within the given solid angle */
327		void
328	<	claimlight(COLOR intens, TRITREE *node, const FVECT cent, double mincos)
328	>	claimlight(COLOR intens, TRITREE *node, FVECT cent, double maxr2)
329		{
330		int     remaining;
331		int     i;
332		if (isleaf(node)) {             /* claim contribution */
333		COLOR   contrib;
334		if (node->val[EXP] <= 0)
335	<	return;
336	<	if (DOT(node->sd,cent) < mincos)
337	<	return;
335	>	return;         /* already claimed */
336	>	if (fdir2diff(node->sd,cent) > maxr2)
337	>	return;         /* too far away */
338		colr_color(contrib, node->val);
339		scalecolor(contrib, tri_omegav(node->gdv));
340		addcolor(intens, contrib);
341		copycolr(node->val, blkclr);
342		return;
343		}
344	<	if (DOT(node->gdv[0],node->gdv[1]) < mincos &&
345	<	DOT(node->gdv[0],cent) > mincos &&
346	<	DOT(node->gdv[1],cent) > mincos &&
347	<	DOT(node->gdv[2],cent) > mincos)
348	<	return;
344	>	if (dir2diff(node->gdv[0],node->gdv[1]) > maxr2 &&
345	>	fdir2diff(node->gdv[0],cent) < maxr2 &&
346	>	fdir2diff(node->gdv[1],cent) < maxr2 &&
347	>	fdir2diff(node->gdv[2],cent) < maxr2)
348	>	return;                 /* previously claimed node */
349		remaining = 0;                  /* recurse on children */
350		for (i = 0; i < 4; i++) {
351	<	claimlight(intens, &node->kid[i], cent, mincos);
351	>	claimlight(intens, &node->kid[i], cent, maxr2);
352		if (!isleaf(&node->kid[i]) || node->kid[i].val[EXP] != 0)
353		++remaining;
354		}
#	Line 339 | Line 357 | claimlight(COLOR intens, TRITREE *node, const FVECT ce
357		/* consolidate empties */
358		free((void *)node->kid); node->kid = NULL;
359		copycolr(node->val, blkclr);
360	<	VCOPY(node->sd, node->gdv[0]);  /* doesn't really matter */
360	>	node->sd = node->gdv[0];        /* doesn't really matter */
361		}
362
363		/* Add lost light contribution to the given list */
364		void
365	<	add2lost(LOSTLIGHT **llp, COLOR intens, const FVECT cent)
365	>	add2lost(LOSTLIGHT **llp, COLOR intens, FVECT cent)
366		{
367		LOSTLIGHT       *newll = (LOSTLIGHT *)malloc(sizeof(LOSTLIGHT));
368
369		if (newll == NULL)
370		return;
371		copycolor(newll->intens, intens);
372	<	VCOPY(newll->sd, cent);
372	>	newll->sd = encodedir(cent);
373		newll->next = *llp;
374		*llp = newll;
375		}
376
377		/* Check lost light list for contributions */
378		void
379	<	getlost(LOSTLIGHT **llp, COLOR intens, const FVECT cent, double omega)
379	>	getlost(LOSTLIGHT **llp, COLOR intens, FVECT cent, double omega)
380		{
381	<	const double    mincos = 1. - omega/(2.*PI);
381	>	const double    maxr2 = omega/PI;
382		LOSTLIGHT       lhead, *lastp, *thisp;
383
384		lhead.next = *llp;
385		lastp = &lhead;
386		while ((thisp = lastp->next) != NULL)
387	<	if (DOT(thisp->sd,cent) >= mincos) {
387	>	if (fdir2diff(thisp->sd,cent) <= maxr2) {
388		LOSTLIGHT       *mynext = thisp->next;
389		addcolor(intens, thisp->intens);
390		free((void *)thisp);
#	Line 419 | Line 437 | mksources(TRITREE *samptree, double thresh, double max
437		if (startleaf == NULL)
438		break;
439		/* claim it */
440	<	VCOPY(curcent, startleaf->sd);
440	>	decodedir(curcent, startleaf->sd);
441		curomega = tri_omegav(startleaf->gdv);
442		currad = sqrt(curomega/PI);
443		growstep = 3.*currad;
#	Line 435 | Line 453 | mksources(TRITREE *samptree, double thresh, double max
453		vsum[0] = vsum[1] = vsum[2] = .0;
454		for (i = 0; i < NTRUNKBR; i++)
455		vector_sum(vsum, &samptree->kid[i],
456	<	curcent, cos(currad+growstep),
456	>	curcent, 2.-2.*cos(currad+growstep),
457		thisethresh);
458		if (normalize(vsum) == .0)
459		break;
#	Line 451 | Line 469 | mksources(TRITREE *samptree, double thresh, double max
469		curomega = 2.*PI*(1. - cos(currad));
470		for (i = 0; i < NTRUNKBR; i++)
471		claimlight(curintens, &samptree->kid[i],
472	<	curcent, cos(currad));
472	>	curcent, 2.-2.*cos(currad));
473		} while (curomega < maxomega &&
474		bright(curintens)/curomega > thisthresh);
475		if (bright(curintens) < minintens) {

Diff Legend

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