ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/sphere.c
Revision: 1.2
Committed: Wed Apr 19 22:24:30 1989 UTC (35 years ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 1.1: +17 -13 lines
Log Message:
moved setting of ray transformation to intersection routines (bug)

File Contents

# User Rev Content
1 greg 1.1 /* Copyright (c) 1986 Regents of the University of California */
2    
3     #ifndef lint
4     static char SCCSid[] = "$SunId$ LBL";
5     #endif
6    
7     /*
8     * sphere.c - compute ray intersection with spheres.
9     *
10     * 8/19/85
11     */
12    
13     #include "ray.h"
14    
15     #include "otypes.h"
16    
17    
18     o_sphere(so, r) /* compute intersection with sphere */
19     OBJREC *so;
20     register RAY *r;
21     {
22     double a, b, c; /* coefficients for quadratic equation */
23     double root[2]; /* quadratic roots */
24     int nroots;
25     double t;
26     register double *ap;
27     register int i;
28    
29     if (so->oargs.nfargs != 4 || so->oargs.farg[3] <= FTINY)
30     objerror(so, USER, "bad arguments");
31    
32     ap = so->oargs.farg;
33    
34     /*
35     * We compute the intersection by substituting into
36     * the surface equation for the sphere. The resulting
37     * quadratic equation in t is then solved for the
38     * smallest positive root, which is our point of
39     * intersection.
40     * Because the ray direction is normalized, a is always 1.
41     */
42    
43     a = 1.0; /* compute quadratic coefficients */
44     b = c = 0.0;
45     for (i = 0; i < 3; i++) {
46     t = r->rorg[i] - ap[i];
47     b += 2.0*r->rdir[i]*t;
48     c += t*t;
49     }
50     c -= ap[3] * ap[3];
51    
52     nroots = quadratic(root, a, b, c); /* solve quadratic */
53    
54     for (i = 0; i < nroots; i++) /* get smallest positive */
55     if ((t = root[i]) > FTINY)
56     break;
57     if (i >= nroots)
58     return(0); /* no positive root */
59    
60 greg 1.2 if (t >= r->rot)
61     return(0); /* other is closer */
62    
63     r->ro = so;
64     r->rot = t;
65     /* compute normal */
66     a = ap[3];
67     if (so->otype == OBJ_BUBBLE)
68     a = -a; /* reverse */
69     for (i = 0; i < 3; i++) {
70     r->rop[i] = r->rorg[i] + r->rdir[i]*t;
71     r->ron[i] = (r->rop[i] - ap[i]) / a;
72 greg 1.1 }
73 greg 1.2 r->rod = -DOT(r->rdir, r->ron);
74     r->rofs = 1.0; setident4(r->rofx);
75     r->robs = 1.0; setident4(r->robx);
76    
77     return(1); /* hit */
78 greg 1.1 }