--- ray/src/rt/srcsamp.c	2003/02/22 02:07:29	2.7
+++ ray/src/rt/srcsamp.c	2009/06/06 02:11:44	2.16
@@ -1,5 +1,5 @@
 #ifndef lint
-static const char	RCSid[] = "$Id: srcsamp.c,v 2.7 2003/02/22 02:07:29 greg Exp $";
+static const char	RCSid[] = "$Id: srcsamp.c,v 2.16 2009/06/06 02:11:44 greg Exp $";
 #endif
 /*
  * Source sampling routines
@@ -7,62 +7,7 @@ static const char	RCSid[] = "$Id: srcsamp.c,v 2.7 2003
  *  External symbols declared in source.h
  */
 
-/* ====================================================================
- * The Radiance Software License, Version 1.0
- *
- * Copyright (c) 1990 - 2002 The Regents of the University of California,
- * through Lawrence Berkeley National Laboratory.   All rights reserved.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- *
- * 1. Redistributions of source code must retain the above copyright
- *         notice, this list of conditions and the following disclaimer.
- *
- * 2. Redistributions in binary form must reproduce the above copyright
- *       notice, this list of conditions and the following disclaimer in
- *       the documentation and/or other materials provided with the
- *       distribution.
- *
- * 3. The end-user documentation included with the redistribution,
- *           if any, must include the following acknowledgment:
- *             "This product includes Radiance software
- *                 (http://radsite.lbl.gov/)
- *                 developed by the Lawrence Berkeley National Laboratory
- *               (http://www.lbl.gov/)."
- *       Alternately, this acknowledgment may appear in the software itself,
- *       if and wherever such third-party acknowledgments normally appear.
- *
- * 4. The names "Radiance," "Lawrence Berkeley National Laboratory"
- *       and "The Regents of the University of California" must
- *       not be used to endorse or promote products derived from this
- *       software without prior written permission. For written
- *       permission, please contact radiance@radsite.lbl.gov.
- *
- * 5. Products derived from this software may not be called "Radiance",
- *       nor may "Radiance" appear in their name, without prior written
- *       permission of Lawrence Berkeley National Laboratory.
- *
- * THIS SOFTWARE IS PROVIDED ``AS IS'' AND ANY EXPRESSED OR IMPLIED
- * WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
- * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
- * DISCLAIMED.   IN NO EVENT SHALL Lawrence Berkeley National Laboratory OR
- * ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
- * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
- * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
- * USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
- * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
- * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
- * OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
- * SUCH DAMAGE.
- * ====================================================================
- *
- * This software consists of voluntary contributions made by many
- * individuals on behalf of Lawrence Berkeley National Laboratory.   For more
- * information on Lawrence Berkeley National Laboratory, please see
- * <http://www.lbl.gov/>.
- */
+#include "copyright.h"
 
 #include  "ray.h"
 
@@ -80,6 +25,7 @@ register RAY  *r;		/* origin is read, direction is set
 register SRCINDEX  *si;		/* source index (modified to current) */
 {
 	int  cent[3], size[3], parr[2];
+	SRCREC  *srcp;
 	FVECT  vpos;
 	double  d;
 	register int  i;
@@ -106,13 +52,14 @@ nextsample:
 	if (!skipparts(cent, size, parr, si->spt))
 		error(CONSISTENCY, "bad source partition in nextssamp");
 					/* compute sample */
+	srcp = source + si->sn;
 	if (dstrsrc > FTINY) {			/* jitter sample */
 		dimlist[ndims] = si->sn + 8831;
 		dimlist[ndims+1] = si->sp + 3109;
 		d = urand(ilhash(dimlist,ndims+2)+samplendx);
-		if (source[si->sn].sflags & SFLAT) {
+		if (srcp->sflags & SFLAT) {
 			multisamp(vpos, 2, d);
-			vpos[2] = 0.5;
+			vpos[SW] = 0.5;
 		} else
 			multisamp(vpos, 3, d);
 		for (i = 0; i < 3; i++)
@@ -123,14 +70,40 @@ nextsample:
 
 	for (i = 0; i < 3; i++)
 		vpos[i] += (double)cent[i]/MAXSPART;
+					/* avoid circular aiming failures */
+	if ((srcp->sflags & SCIR) && (si->np > 1 || dstrsrc > 0.7)) {
+		FVECT	trim;
+		if (srcp->sflags & (SFLAT|SDISTANT)) {
+			d = 1.12837917;		/* correct setflatss() */
+			trim[SU] = d*sqrt(1.0 - 0.5*vpos[SV]*vpos[SV]);
+			trim[SV] = d*sqrt(1.0 - 0.5*vpos[SU]*vpos[SU]);
+			trim[SW] = 0.0;
+		} else {
+			trim[SW] = trim[SU] = vpos[SU]*vpos[SU];
+			d = vpos[SV]*vpos[SV];
+			if (d > trim[SW]) trim[SW] = d;
+			trim[SU] += d;
+			d = vpos[SW]*vpos[SW];
+			if (d > trim[SW]) trim[SW] = d;
+			trim[SU] += d;
+			if (trim[SU] > FTINY*FTINY) {
+				d = 1.0/0.7236;	/* correct sphsetsrc() */
+				trim[SW] = trim[SV] = trim[SU] =
+						d*sqrt(trim[SW]/trim[SU]);
+			} else
+				trim[SW] = trim[SV] = trim[SU] = 0.0;
+		}
+		for (i = 0; i < 3; i++)
+			vpos[i] *= trim[i];
+	}
 					/* compute direction */
 	for (i = 0; i < 3; i++)
-		r->rdir[i] = source[si->sn].sloc[i] +
-				vpos[SU]*source[si->sn].ss[SU][i] +
-				vpos[SV]*source[si->sn].ss[SV][i] +
-				vpos[SW]*source[si->sn].ss[SW][i];
+		r->rdir[i] = srcp->sloc[i] +
+				vpos[SU]*srcp->ss[SU][i] +
+				vpos[SV]*srcp->ss[SV][i] +
+				vpos[SW]*srcp->ss[SW][i];
 
-	if (!(source[si->sn].sflags & SDISTANT))
+	if (!(srcp->sflags & SDISTANT))
 		for (i = 0; i < 3; i++)
 			r->rdir[i] -= r->rorg[i];
 					/* compute distance */
@@ -138,23 +111,23 @@ nextsample:
 		goto nextsample;		/* at source! */
 
 					/* compute sample size */
-	if (source[si->sn].sflags & SFLAT) {
+	if (srcp->sflags & SFLAT) {
 		si->dom = sflatform(si->sn, r->rdir);
 		si->dom *= size[SU]*size[SV]/(MAXSPART*(double)MAXSPART);
-	} else if (source[si->sn].sflags & SCYL) {
+	} else if (srcp->sflags & SCYL) {
 		si->dom = scylform(si->sn, r->rdir);
 		si->dom *= size[SU]/(double)MAXSPART;
 	} else {
 		si->dom = size[SU]*size[SV]*(double)size[SW] /
 				(MAXSPART*MAXSPART*(double)MAXSPART) ;
 	}
-	if (source[si->sn].sflags & SDISTANT) {
-		si->dom *= source[si->sn].ss2;
+	if (srcp->sflags & SDISTANT) {
+		si->dom *= srcp->ss2;
 		return(FHUGE);
 	}
 	if (si->dom <= 1e-4)
 		goto nextsample;		/* behind source? */
-	si->dom *= source[si->sn].ss2/(d*d);
+	si->dom *= srcp->ss2/(d*d);
 	return(d);		/* sample OK, return distance */
 }
 
@@ -169,12 +142,13 @@ unsigned char  *pt;		/* partition array */
 					/* check this partition */
 	p = spart(pt, pp[0]);
 	pp[0]++;
-	if (p == S0)			/* leaf partition */
+	if (p == S0) {			/* leaf partition */
 		if (pp[1]) {
 			pp[1]--;
 			return(0);	/* not there yet */
 		} else
 			return(1);	/* we've arrived */
+	}
 				/* else check lower */
 	sz[p] >>= 1;
 	ct[p] -= sz[p];
@@ -284,7 +258,7 @@ flatpart(si, r)				/* partition a flat source */
 register SRCINDEX  *si;
 register RAY  *r;
 {
-	register FLOAT  *vp;
+	register RREAL  *vp;
 	FVECT  v;
 	double  du2, dv2;
 	int  pi;
@@ -295,7 +269,7 @@ register RAY  *r;
 	v[1] = r->rorg[1] - vp[1];
 	v[2] = r->rorg[2] - vp[2];
 	vp = source[si->sn].snorm;
-	if (DOT(v,vp) <= FTINY) {	/* behind source */
+	if (DOT(v,vp) <= 0.) {		/* behind source */
 		si->np = 0;
 		return;
 	}
@@ -368,7 +342,7 @@ scylform(sn, dir)		/* compute cosine for cylinder's pr
 int  sn;
 register FVECT  dir;		/* assume normalized */
 {
-	register FLOAT  *dv;
+	register RREAL  *dv;
 	double  d;
 
 	dv = source[sn].ss[SU];