ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/common/fvect.c
(Generate patch)

Comparing ray/src/common/fvect.c (file contents):
Revision 1.4 by greg, Tue Oct 17 13:35:30 1989 UTC vs.
Revision 2.5 by gwlarson, Wed Nov 25 16:56:33 1998 UTC

# Line 1 | Line 1
1 < /* Copyright (c) 1986 Regents of the University of California */
1 > /* Copyright (c) 1998 Silicon Graphics, Inc. */
2  
3   #ifndef lint
4 < static char SCCSid[] = "$SunId$ LBL";
4 > static char SCCSid[] = "$SunId$ SGI";
5   #endif
6  
7   /*
# Line 10 | Line 10 | static char SCCSid[] = "$SunId$ LBL";
10   *     8/14/85
11   */
12  
13 + #include  <math.h>
14   #include  "fvect.h"
15  
15 #define  FTINY          1e-7
16  
17
17   double
18   fdot(v1, v2)                    /* return the dot product of two vectors */
19   register FVECT  v1, v2;
# Line 27 | Line 26 | double
26   dist2(p1, p2)                   /* return square of distance between points */
27   register FVECT  p1, p2;
28   {
29 <        static FVECT  delta;
29 >        FVECT  delta;
30  
31          delta[0] = p2[0] - p1[0];
32          delta[1] = p2[1] - p1[1];
33          delta[2] = p2[2] - p1[2];
34 +
35          return(DOT(delta, delta));
36   }
37  
# Line 41 | Line 41 | dist2line(p, ep1, ep2)         /* return square of distance t
41   FVECT  p;               /* the point */
42   FVECT  ep1, ep2;        /* points on the line */
43   {
44 <        static double  d, d1, d2;
44 >        register double  d, d1, d2;
45  
46          d = dist2(ep1, ep2);
47          d1 = dist2(ep1, p);
48 <        d2 = dist2(ep2, p);
48 >        d2 = d + d1 - dist2(ep2, p);
49  
50 <        return(d1 - (d+d1-d2)*(d+d1-d2)/d/4);
50 >        return(d1 - 0.25*d2*d2/d);
51   }
52  
53  
# Line 56 | Line 56 | dist2lseg(p, ep1, ep2)         /* return square of distance t
56   FVECT  p;               /* the point */
57   FVECT  ep1, ep2;        /* the end points */
58   {
59 <        static double  d, d1, d2;
59 >        register double  d, d1, d2;
60  
61          d = dist2(ep1, ep2);
62          d1 = dist2(ep1, p);
# Line 69 | Line 69 | FVECT  ep1, ep2;       /* the end points */
69                  if (d1 - d2 > d)
70                          return(d2);
71          }
72 +        d2 = d + d1 - d2;
73  
74 <        return(d1 - (d+d1-d2)*(d+d1-d2)/d/4);   /* distance to line */
74 >        return(d1 - 0.25*d2*d2/d);      /* distance to line */
75   }
76  
77  
# Line 84 | Line 85 | register FVECT  vres, v1, v2;
85  
86  
87   fvsum(vres, v0, v1, f)          /* vres = v0 + f*v1 */
88 < FVECT  vres, v0, v1;
89 < double  f;
88 > register FVECT  vres, v0, v1;
89 > register double  f;
90   {
91          vres[0] = v0[0] + f*v1[0];
92          vres[1] = v0[1] + f*v1[1];
# Line 97 | Line 98 | double
98   normalize(v)                    /* normalize a vector, return old magnitude */
99   register FVECT  v;
100   {
101 <        static double  len;
101 >        register double  len, d;
102          
103 <        len = DOT(v, v);
103 >        d = DOT(v, v);
104          
105 <        if (len <= 0.0)
105 >        if (d <= 0.0)
106                  return(0.0);
107          
108 <        /****** problematic
109 <        if (len >= (1.0-FTINY)*(1.0-FTINY) &&
110 <                        len <= (1.0+FTINY)*(1.0+FTINY))
111 <                return(1.0);
111 <        ******/
108 >        if (d <= 1.0+FTINY && d >= 1.0-FTINY)
109 >                len = 0.5 + 0.5*d;      /* first order approximation */
110 >        else
111 >                len = sqrt(d);
112  
113 <        len = sqrt(len);
114 <        v[0] /= len;
115 <        v[1] /= len;
116 <        v[2] /= len;
113 >        v[0] *= d = 1.0/len;
114 >        v[1] *= d;
115 >        v[2] *= d;
116 >
117          return(len);
118 + }
119 +
120 +
121 + spinvector(vres, vorig, vnorm, theta)   /* rotate vector around normal */
122 + FVECT  vres, vorig, vnorm;
123 + double  theta;
124 + {
125 +        double  sint, cost, normprod;
126 +        FVECT  vperp;
127 +        register int  i;
128 +        
129 +        if (theta == 0.0) {
130 +                if (vres != vorig)
131 +                        VCOPY(vres, vorig);
132 +                return;
133 +        }
134 +        cost = cos(theta);
135 +        sint = sin(theta);
136 +        normprod = DOT(vorig, vnorm)*(1.-cost);
137 +        fcross(vperp, vnorm, vorig);
138 +        for (i = 0; i < 3; i++)
139 +                vres[i] = vorig[i]*cost + vnorm[i]*normprod + vperp[i]*sint;
140   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines