@@ -77,30 +77,34 @@ Vec trianglePoint( const Vec &v0, const Vec &v1, const Vec &v2, const Imath::Vec
7777
7878// / Implementation derived from Wild Magic (Version 2) Software Library, available
7979// / from http://www.geometrictools.com/Downloads/WildMagic2p5.zip under free license
80+ // / This link is now offline, but it was presumably similar to the code now found here:
81+ // / https://www.geometrictools.com/GTE/Mathematics/DistPointTriangle.h
82+ // / with explanation here:
83+ // / https://www.geometrictools.com/Documentation/DistancePoint3Triangle3.pdf
8084template <class Vec >
8185typename VectorTraits<Vec>::BaseType triangleClosestBarycentric ( const Vec &v0, const Vec &v1, const Vec &v2, const Vec &p, Imath::Vec3<typename VectorTraits<Vec>::BaseType> &barycentric )
8286{
8387 typedef typename VectorTraits<Vec>::BaseType Real;
8488
85- Vec triOrigin = v0;
8689 Vec triEdge0, triEdge1;
8790
91+ // \todo - do we really need to support vector types that don't have a subtraction operator?
92+ // This would be clearer if rewritten to use the operators on the types instead of these VectorOps,
93+ // and more consistent with everything else we do.
8894 vecSub ( v1, v0, triEdge0 );
8995 vecSub ( v2, v0, triEdge1 );
9096
9197 Vec kDiff ;
92- vecSub ( triOrigin , p, kDiff );
98+ vecSub ( v0 , p, kDiff );
9399
94100 Real a00 = vecDot ( triEdge0, triEdge0 );
95101 Real a01 = vecDot ( triEdge0, triEdge1 );
96102 Real a11 = vecDot ( triEdge1, triEdge1 );
97103 Real b0 = vecDot ( kDiff , triEdge0 );
98104 Real b1 = vecDot ( kDiff , triEdge1 );
99- Real c = vecDot ( kDiff , kDiff );
100105 Real det = Real (fabs (a00*a11-a01*a01));
101106 Real s = a01*b1-a11*b0;
102107 Real t = a01*b0-a00*b1;
103- Real distSqrd;
104108
105109 if ( s + t <= det )
106110 {
@@ -114,12 +118,10 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
114118 if ( -b0 >= a00 )
115119 {
116120 s = Real (1.0 );
117- distSqrd = a00+Real (2.0 )*b0+c;
118121 }
119122 else
120123 {
121124 s = -b0/a00;
122- distSqrd = b0*s+c;
123125 }
124126 }
125127 else
@@ -128,17 +130,14 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
128130 if ( b1 >= Real (0.0 ) )
129131 {
130132 t = Real (0.0 );
131- distSqrd = c;
132133 }
133134 else if ( -b1 >= a11 )
134135 {
135136 t = Real (1.0 );
136- distSqrd = a11+Real (2.0 )*b1+c;
137137 }
138138 else
139139 {
140140 t = -b1/a11;
141- distSqrd = b1*t+c;
142141 }
143142 }
144143 }
@@ -148,17 +147,14 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
148147 if ( b1 >= Real (0.0 ) )
149148 {
150149 t = Real (0.0 );
151- distSqrd = c;
152150 }
153151 else if ( -b1 >= a11 )
154152 {
155153 t = Real (1.0 );
156- distSqrd = a11+Real (2.0 )*b1+c;
157154 }
158155 else
159156 {
160157 t = -b1/a11;
161- distSqrd = b1*t+c;
162158 }
163159 }
164160 }
@@ -168,17 +164,14 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
168164 if ( b0 >= Real (0.0 ) )
169165 {
170166 s = Real (0.0 );
171- distSqrd = c;
172167 }
173168 else if ( -b0 >= a00 )
174169 {
175170 s = Real (1.0 );
176- distSqrd = a00+Real (2.0 )*b0+c;
177171 }
178172 else
179173 {
180174 s = -b0/a00;
181- distSqrd = b0*s+c;
182175 }
183176 }
184177 else // region 0
@@ -188,15 +181,12 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
188181 {
189182 s = Real (0.0 );
190183 t = Real (0.0 );
191- distSqrd = std::numeric_limits<Real>::max ();
192184 }
193185 else
194186 {
195187 Real invDet = Real (1.0 )/det;
196188 s *= invDet;
197189 t *= invDet;
198- distSqrd = s*(a00*s+a01*t+Real (2.0 )*b0) +
199- t*(a01*s+a11*t+Real (2.0 )*b1)+c;
200190 }
201191 }
202192 }
@@ -216,14 +206,11 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
216206 {
217207 s = Real (1.0 );
218208 t = Real (0.0 );
219- distSqrd = a00+Real (2.0 )*b0+c;
220209 }
221210 else
222211 {
223212 s = numer/denom;
224213 t = Real (1.0 ) - s;
225- distSqrd = s*(a00*s+a01*t+Real (2.0 )*b0) +
226- t*(a01*s+a11*t+Real (2.0 )*b1)+c;
227214 }
228215 }
229216 else
@@ -232,17 +219,14 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
232219 if ( tmp1 <= Real (0.0 ) )
233220 {
234221 t = Real (1.0 );
235- distSqrd = a11+Real (2.0 )*b1+c;
236222 }
237223 else if ( b1 >= Real (0.0 ) )
238224 {
239225 t = Real (0.0 );
240- distSqrd = c;
241226 }
242227 else
243228 {
244229 t = -b1/a11;
245- distSqrd = b1*t+c;
246230 }
247231 }
248232 }
@@ -258,14 +242,11 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
258242 {
259243 t = Real (1.0 );
260244 s = Real (0.0 );
261- distSqrd = a11+Real (2.0 )*b1+c;
262245 }
263246 else
264247 {
265248 t = numer/denom;
266249 s = Real (1.0 ) - t;
267- distSqrd = s*(a00*s+a01*t+Real (2.0 )*b0) +
268- t*(a01*s+a11*t+Real (2.0 )*b1)+c;
269250 }
270251 }
271252 else
@@ -274,17 +255,14 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
274255 if ( tmp1 <= Real (0.0 ) )
275256 {
276257 s = Real (1.0 );
277- distSqrd = a00+Real (2.0 )*b0+c;
278258 }
279259 else if ( b0 >= Real (0.0 ) )
280260 {
281261 s = Real (0.0 );
282- distSqrd = c;
283262 }
284263 else
285264 {
286265 s = -b0/a00;
287- distSqrd = b0*s+c;
288266 }
289267 }
290268 }
@@ -295,7 +273,6 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
295273 {
296274 s = Real (0.0 );
297275 t = Real (1.0 );
298- distSqrd = a11+Real (2.0 )*b1+c;
299276 }
300277 else
301278 {
@@ -304,14 +281,11 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
304281 {
305282 s = Real (1.0 );
306283 t = Real (0.0 );
307- distSqrd = a00+Real (2.0 )*b0+c;
308284 }
309285 else
310286 {
311287 s = numer/denom;
312288 t = Real (1.0 ) - s;
313- distSqrd = s*(a00*s+a01*t+Real (2.0 )*b0) +
314- t*(a01*s+a11*t+Real (2.0 )*b1)+c;
315289 }
316290 }
317291 }
@@ -321,7 +295,9 @@ typename VectorTraits<Vec>::BaseType triangleClosestBarycentric( const Vec &v0,
321295 barycentric.z = t;
322296 barycentric.x = 1 - s - t;
323297
324- return Real (fabs (distSqrd));
298+ Vec closest = vecAdd ( vecAdd ( v0, vecMul ( triEdge0, s ) ), vecMul (triEdge1, t ) );
299+ Vec diff = vecSub ( closest, p );
300+ return vecDot ( diff, diff );
325301}
326302
327303template <class Vec >
@@ -399,13 +375,13 @@ bool triangleContainsPoint( const Vec &v0, const Vec &v1, const Vec &v2, const V
399375 Real d = Real ( 1 ) / ( dotAA * dotBB - dotAB * dotAB );
400376
401377 barycentric.z = (dotBB * dotAC - dotAB * dotBC) * d;
402- if ( barycentric.z < Real ( 0 ) || barycentric.z > Real ( 1 ) )
378+ if ( !( barycentric.z >= Real ( 0 ) && barycentric.z <= Real ( 1 ) ) )
403379 {
404380 return false ;
405381 }
406382
407383 barycentric.y = (dotAA * dotBC - dotAB * dotAC) * d;
408- if ( barycentric.y < Real ( 0 ) || barycentric.z + barycentric.y > Real ( 1 ) )
384+ if ( !( barycentric.y >= Real ( 0 ) && barycentric.z + barycentric.y <= Real ( 1 ) ) )
409385 {
410386 return false ;
411387 }
0 commit comments