| 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 |  |  | *  fvect.c - routines for float vector calculations | 
| 9 |  |  | * | 
| 10 |  |  | *     8/14/85 | 
| 11 |  |  | */ | 
| 12 |  |  |  | 
| 13 |  |  | #include  "fvect.h" | 
| 14 |  |  |  | 
| 15 |  |  |  | 
| 16 |  |  | double | 
| 17 |  |  | fdot(v1, v2)                    /* return the dot product of two vectors */ | 
| 18 |  |  | register FVECT  v1, v2; | 
| 19 |  |  | { | 
| 20 |  |  | return(DOT(v1,v2)); | 
| 21 |  |  | } | 
| 22 |  |  |  | 
| 23 |  |  |  | 
| 24 |  |  | double | 
| 25 |  |  | dist2(p1, p2)                   /* return square of distance between points */ | 
| 26 |  |  | register FVECT  p1, p2; | 
| 27 |  |  | { | 
| 28 |  |  | static FVECT  delta; | 
| 29 |  |  |  | 
| 30 |  |  | delta[0] = p2[0] - p1[0]; | 
| 31 |  |  | delta[1] = p2[1] - p1[1]; | 
| 32 |  |  | delta[2] = p2[2] - p1[2]; | 
| 33 |  |  | return(DOT(delta, delta)); | 
| 34 |  |  | } | 
| 35 |  |  |  | 
| 36 |  |  |  | 
| 37 |  |  | double | 
| 38 |  |  | dist2line(p, ep1, ep2)          /* return square of distance to line */ | 
| 39 |  |  | FVECT  p;               /* the point */ | 
| 40 |  |  | FVECT  ep1, ep2;        /* points on the line */ | 
| 41 |  |  | { | 
| 42 |  |  | static double  d, d1, d2; | 
| 43 |  |  |  | 
| 44 |  |  | d = dist2(ep1, ep2); | 
| 45 |  |  | d1 = dist2(ep1, p); | 
| 46 |  |  | d2 = dist2(ep2, p); | 
| 47 |  |  |  | 
| 48 |  |  | return(d1 - (d+d1-d2)*(d+d1-d2)/d/4); | 
| 49 |  |  | } | 
| 50 |  |  |  | 
| 51 |  |  |  | 
| 52 |  |  | double | 
| 53 |  |  | dist2lseg(p, ep1, ep2)          /* return square of distance to line segment */ | 
| 54 |  |  | FVECT  p;               /* the point */ | 
| 55 |  |  | FVECT  ep1, ep2;        /* the end points */ | 
| 56 |  |  | { | 
| 57 |  |  | static double  d, d1, d2; | 
| 58 |  |  |  | 
| 59 |  |  | d = dist2(ep1, ep2); | 
| 60 |  |  | d1 = dist2(ep1, p); | 
| 61 |  |  | d2 = dist2(ep2, p); | 
| 62 |  |  |  | 
| 63 |  |  | if (d2 > d1) {                  /* check if past endpoints */ | 
| 64 |  |  | if (d2 - d1 > d) | 
| 65 |  |  | return(d1); | 
| 66 |  |  | } else { | 
| 67 |  |  | if (d1 - d2 > d) | 
| 68 |  |  | return(d2); | 
| 69 |  |  | } | 
| 70 |  |  |  | 
| 71 |  |  | return(d1 - (d+d1-d2)*(d+d1-d2)/d/4);   /* distance to line */ | 
| 72 |  |  | } | 
| 73 |  |  |  | 
| 74 |  |  |  | 
| 75 |  |  | fcross(vres, v1, v2)            /* vres = v1 X v2 */ | 
| 76 |  |  | register FVECT  vres, v1, v2; | 
| 77 |  |  | { | 
| 78 |  |  | vres[0] = v1[1]*v2[2] - v1[2]*v2[1]; | 
| 79 |  |  | vres[1] = v1[2]*v2[0] - v1[0]*v2[2]; | 
| 80 |  |  | vres[2] = v1[0]*v2[1] - v1[1]*v2[0]; | 
| 81 |  |  | } | 
| 82 |  |  |  | 
| 83 |  |  |  | 
| 84 | greg | 1.4 | fvsum(vres, v0, v1, f)          /* vres = v0 + f*v1 */ | 
| 85 |  |  | FVECT  vres, v0, v1; | 
| 86 |  |  | double  f; | 
| 87 |  |  | { | 
| 88 |  |  | vres[0] = v0[0] + f*v1[0]; | 
| 89 |  |  | vres[1] = v0[1] + f*v1[1]; | 
| 90 |  |  | vres[2] = v0[2] + f*v1[2]; | 
| 91 |  |  | } | 
| 92 |  |  |  | 
| 93 |  |  |  | 
| 94 | greg | 1.1 | double | 
| 95 |  |  | normalize(v)                    /* normalize a vector, return old magnitude */ | 
| 96 |  |  | register FVECT  v; | 
| 97 |  |  | { | 
| 98 |  |  | static double  len; | 
| 99 |  |  |  | 
| 100 |  |  | len = DOT(v, v); | 
| 101 |  |  |  | 
| 102 | greg | 1.3 | if (len <= 0.0) | 
| 103 | greg | 1.1 | return(0.0); | 
| 104 |  |  |  | 
| 105 | greg | 1.2 | /****** problematic | 
| 106 | greg | 1.1 | if (len >= (1.0-FTINY)*(1.0-FTINY) && | 
| 107 |  |  | len <= (1.0+FTINY)*(1.0+FTINY)) | 
| 108 |  |  | return(1.0); | 
| 109 | greg | 1.2 | ******/ | 
| 110 | greg | 1.1 |  | 
| 111 |  |  | len = sqrt(len); | 
| 112 |  |  | v[0] /= len; | 
| 113 |  |  | v[1] /= len; | 
| 114 |  |  | v[2] /= len; | 
| 115 |  |  | return(len); | 
| 116 |  |  | } | 
| 117 | greg | 1.5 |  | 
| 118 |  |  |  | 
| 119 |  |  | spinvector(vres, vorig, vnorm, theta)   /* rotate vector around normal */ | 
| 120 |  |  | FVECT  vres, vorig, vnorm; | 
| 121 |  |  | double  theta; | 
| 122 |  |  | { | 
| 123 | greg | 1.6 | extern double  cos(), sin(); | 
| 124 |  |  | double  sint, cost, normprod; | 
| 125 | greg | 1.5 | FVECT  vperp; | 
| 126 |  |  | register int  i; | 
| 127 |  |  |  | 
| 128 |  |  | if (theta == 0.0) { | 
| 129 | greg | 1.6 | if (vres != vorig) | 
| 130 |  |  | VCOPY(vres, vorig); | 
| 131 | greg | 1.5 | return; | 
| 132 |  |  | } | 
| 133 | greg | 1.6 | cost = cos(theta); | 
| 134 | greg | 1.5 | sint = sin(theta); | 
| 135 | greg | 1.6 | normprod = DOT(vorig, vnorm)*(1.-cost); | 
| 136 | greg | 1.5 | fcross(vperp, vnorm, vorig); | 
| 137 |  |  | for (i = 0; i < 3; i++) | 
| 138 | greg | 1.6 | vres[i] = vorig[i]*cost + vnorm[i]*normprod + vperp[i]*sint; | 
| 139 | greg | 1.5 | } |