1 /////////////////////////////////////////////////////////////////////////// 2 // 3 // Copyright (c) 2002, Industrial Light & Magic, a division of Lucas 4 // Digital Ltd. LLC 5 // 6 // All rights reserved. 7 // 8 // Redistribution and use in source and binary forms, with or without 9 // modification, are permitted provided that the following conditions are 10 // met: 11 // * Redistributions of source code must retain the above copyright 12 // notice, this list of conditions and the following disclaimer. 13 // * Redistributions in binary form must reproduce the above 14 // copyright notice, this list of conditions and the following disclaimer 15 // in the documentation and/or other materials provided with the 16 // distribution. 17 // * Neither the name of Industrial Light & Magic nor the names of 18 // its contributors may be used to endorse or promote products derived 19 // from this software without specific prior written permission. 20 // 21 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 22 // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT 23 // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR 24 // A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT 25 // OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, 26 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT 27 // LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, 28 // DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY 29 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 30 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE 31 // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 32 // 33 /////////////////////////////////////////////////////////////////////////// 34 35 36 37 #ifndef INCLUDED_IMATHLINEALGO_H 38 #define INCLUDED_IMATHLINEALGO_H 39 40 //------------------------------------------------------------------ 41 // 42 // This file contains algorithms applied to or in conjunction 43 // with lines (Imath::Line). These algorithms may require 44 // more headers to compile. The assumption made is that these 45 // functions are called much less often than the basic line 46 // functions or these functions require more support classes 47 // 48 // Contains: 49 // 50 // bool closestPoints(const Line<T>& line1, 51 // const Line<T>& line2, 52 // Vec3<T>& point1, 53 // Vec3<T>& point2) 54 // 55 // bool intersect( const Line3<T> &line, 56 // const Vec3<T> &v0, 57 // const Vec3<T> &v1, 58 // const Vec3<T> &v2, 59 // Vec3<T> &pt, 60 // Vec3<T> &barycentric, 61 // bool &front) 62 // 63 // V3f 64 // closestVertex(const Vec3<T> &v0, 65 // const Vec3<T> &v1, 66 // const Vec3<T> &v2, 67 // const Line3<T> &l) 68 // 69 // V3f 70 // rotatePoint(const Vec3<T> p, Line3<T> l, float angle) 71 // 72 //------------------------------------------------------------------ 73 74 #include "ImathLine.h" 75 #include "ImathVecAlgo.h" 76 #include "ImathFun.h" 77 78 namespace Imath { 79 80 81 template <class T> 82 bool 83 closestPoints 84 (const Line3<T>& line1, 85 const Line3<T>& line2, 86 Vec3<T>& point1, 87 Vec3<T>& point2) 88 { 89 // 90 // Compute point1 and point2 such that point1 is on line1, point2 91 // is on line2 and the distance between point1 and point2 is minimal. 92 // This function returns true if point1 and point2 can be computed, 93 // or false if line1 and line2 are parallel or nearly parallel. 94 // This function assumes that line1.dir and line2.dir are normalized. 95 // 96 97 Vec3<T> w = line1.pos - line2.pos; 98 T d1w = line1.dir ^ w; 99 T d2w = line2.dir ^ w; 100 T d1d2 = line1.dir ^ line2.dir; 101 T n1 = d1d2 * d2w - d1w; 102 T n2 = d2w - d1d2 * d1w; 103 T d = 1 - d1d2 * d1d2; 104 T absD = abs (d); 105 106 if ((absD > 1) || 107 (abs (n1) < limits<T>::max() * absD && 108 abs (n2) < limits<T>::max() * absD)) 109 { 110 point1 = line1 (n1 / d); 111 point2 = line2 (n2 / d); 112 return true; 113 } 114 else 115 { 116 return false; 117 } 118 } 119 120 121 template <class T> 122 bool 123 intersect 124 (const Line3<T> &line, 125 const Vec3<T> &v0, 126 const Vec3<T> &v1, 127 const Vec3<T> &v2, 128 Vec3<T> &pt, 129 Vec3<T> &barycentric, 130 bool &front) 131 { 132 // 133 // Given a line and a triangle (v0, v1, v2), the intersect() function 134 // finds the intersection of the line and the plane that contains the 135 // triangle. 136 // 137 // If the intersection point cannot be computed, either because the 138 // line and the triangle's plane are nearly parallel or because the 139 // triangle's area is very small, intersect() returns false. 140 // 141 // If the intersection point is outside the triangle, intersect 142 // returns false. 143 // 144 // If the intersection point, pt, is inside the triangle, intersect() 145 // computes a front-facing flag and the barycentric coordinates of 146 // the intersection point, and returns true. 147 // 148 // The front-facing flag is true if the dot product of the triangle's 149 // normal, (v2-v1)%(v1-v0), and the line's direction is negative. 150 // 151 // The barycentric coordinates have the following property: 152 // 153 // pt = v0 * barycentric.x + v1 * barycentric.y + v2 * barycentric.z 154 // 155 156 Vec3<T> edge0 = v1 - v0; 157 Vec3<T> edge1 = v2 - v1; 158 Vec3<T> normal = edge1 % edge0; 159 160 T l = normal.length(); 161 162 if (l != 0) 163 normal /= l; 164 else 165 return false; // zero-area triangle 166 167 // 168 // d is the distance of line.pos from the plane that contains the triangle. 169 // The intersection point is at line.pos + (d/nd) * line.dir. 170 // 171 172 T d = normal ^ (v0 - line.pos); 173 T nd = normal ^ line.dir; 174 175 if (abs (nd) > 1 || abs (d) < limits<T>::max() * abs (nd)) 176 pt = line (d / nd); 177 else 178 return false; // line and plane are nearly parallel 179 180 // 181 // Compute the barycentric coordinates of the intersection point. 182 // The intersection is inside the triangle if all three barycentric 183 // coordinates are between zero and one. 184 // 185 186 { 187 Vec3<T> en = edge0.normalized(); 188 Vec3<T> a = pt - v0; 189 Vec3<T> b = v2 - v0; 190 Vec3<T> c = (a - en * (en ^ a)); 191 Vec3<T> d = (b - en * (en ^ b)); 192 T e = c ^ d; 193 T f = d ^ d; 194 195 if (e >= 0 && e <= f) 196 barycentric.z = e / f; 197 else 198 return false; // outside 199 } 200 201 { 202 Vec3<T> en = edge1.normalized(); 203 Vec3<T> a = pt - v1; 204 Vec3<T> b = v0 - v1; 205 Vec3<T> c = (a - en * (en ^ a)); 206 Vec3<T> d = (b - en * (en ^ b)); 207 T e = c ^ d; 208 T f = d ^ d; 209 210 if (e >= 0 && e <= f) 211 barycentric.x = e / f; 212 else 213 return false; // outside 214 } 215 216 barycentric.y = 1 - barycentric.x - barycentric.z; 217 218 if (barycentric.y < 0) 219 return false; // outside 220 221 front = ((line.dir ^ normal) < 0); 222 return true; 223 } 224 225 226 template <class T> 227 Vec3<T> 228 closestVertex 229 (const Vec3<T> &v0, 230 const Vec3<T> &v1, 231 const Vec3<T> &v2, 232 const Line3<T> &l) 233 { 234 Vec3<T> nearest = v0; 235 T neardot = (v0 - l.closestPointTo(v0)).length2(); 236 237 T tmp = (v1 - l.closestPointTo(v1)).length2(); 238 239 if (tmp < neardot) 240 { 241 neardot = tmp; 242 nearest = v1; 243 } 244 245 tmp = (v2 - l.closestPointTo(v2)).length2(); 246 if (tmp < neardot) 247 { 248 neardot = tmp; 249 nearest = v2; 250 } 251 252 return nearest; 253 } 254 255 256 template <class T> 257 Vec3<T> 258 rotatePoint (const Vec3<T> p, Line3<T> l, T angle) 259 { 260 // 261 // Rotate the point p around the line l by the given angle. 262 // 263 264 // 265 // Form a coordinate frame with <x,y,a>. The rotation is the in xy 266 // plane. 267 // 268 269 Vec3<T> q = l.closestPointTo(p); 270 Vec3<T> x = p - q; 271 T radius = x.length(); 272 273 x.normalize(); 274 Vec3<T> y = (x % l.dir).normalize(); 275 276 T cosangle = Math<T>::cos(angle); 277 T sinangle = Math<T>::sin(angle); 278 279 Vec3<T> r = q + x * radius * cosangle + y * radius * sinangle; 280 281 return r; 282 } 283 284 285 } // namespace Imath 286 287 #endif 288
