--- ray/src/rt/sphere.c 1989/04/19 22:24:30 1.2 +++ ray/src/rt/sphere.c 2003/03/12 04:59:05 2.5 @@ -1,15 +1,12 @@ -/* Copyright (c) 1986 Regents of the University of California */ - #ifndef lint -static char SCCSid[] = "$SunId$ LBL"; +static const char RCSid[] = "$Id: sphere.c,v 2.5 2003/03/12 04:59:05 greg Exp $"; #endif - /* * sphere.c - compute ray intersection with spheres. - * - * 8/19/85 */ +#include "copyright.h" + #include "ray.h" #include "otypes.h" @@ -23,13 +20,19 @@ register RAY *r; double root[2]; /* quadratic roots */ int nroots; double t; - register double *ap; + register FLOAT *ap; register int i; - if (so->oargs.nfargs != 4 || so->oargs.farg[3] <= FTINY) - objerror(so, USER, "bad arguments"); - + if (so->oargs.nfargs != 4) + objerror(so, USER, "bad # arguments"); ap = so->oargs.farg; + if (ap[3] < -FTINY) { + objerror(so, WARNING, "negative radius"); + so->otype = so->otype == OBJ_SPHERE ? + OBJ_BUBBLE : OBJ_SPHERE; + ap[3] = -ap[3]; + } else if (ap[3] <= FTINY) + objerror(so, USER, "zero radius"); /* * We compute the intersection by substituting into @@ -37,12 +40,14 @@ register RAY *r; * quadratic equation in t is then solved for the * smallest positive root, which is our point of * intersection. - * Because the ray direction is normalized, a is always 1. + * Since the ray is normalized, a should always be + * one. We compute it here to prevent instability in the + * intersection calculation. */ - - a = 1.0; /* compute quadratic coefficients */ - b = c = 0.0; + /* compute quadratic coefficients */ + a = b = c = 0.0; for (i = 0; i < 3; i++) { + a += r->rdir[i]*r->rdir[i]; t = r->rorg[i] - ap[i]; b += 2.0*r->rdir[i]*t; c += t*t; @@ -71,8 +76,9 @@ register RAY *r; r->ron[i] = (r->rop[i] - ap[i]) / a; } r->rod = -DOT(r->rdir, r->ron); - r->rofs = 1.0; setident4(r->rofx); - r->robs = 1.0; setident4(r->robx); + r->rox = NULL; + r->pert[0] = r->pert[1] = r->pert[2] = 0.0; + r->uv[0] = r->uv[1] = 0.0; return(1); /* hit */ }