source: GTP/trunk/Lib/Vis/Preprocessing/src/Triangle3.cpp @ 2575

Revision 2575, 7.5 KB checked in by bittner, 17 years ago (diff)

big merge: preparation for havran ray caster, check if everything works

RevLine 
[191]1#include "Triangle3.h"
[1328]2#include "Ray.h"
[2575]3#include "SimpleRay.h"
[1328]4#include "AxisAlignedBox3.h"
[1344]5#include "Containers.h"
6#include "Polygon3.h"
[191]7
[1328]8
[2176]9using namespace std;
10
[863]11namespace GtpVisibilityPreprocessor {
[860]12
[1328]13       
14Triangle3::Triangle3(const Vector3 &a, const Vector3 &b, const Vector3 &c)
[191]15{
[1328]16        Init(a, b, c);
[191]17}
[860]18
[1328]19
20void Triangle3::Init(const Vector3 &a, const Vector3 &b, const Vector3 &c)
21{
22        mVertices[0] = a;
23        mVertices[1] = b;
24        mVertices[2] = c;
25}
26
27
28float Triangle3::GetSpatialAngle(const Vector3 &point) const
29{
30        return 0.0f;
31}
32
33
[1867]34int
35Triangle3::CastRay(const Ray &ray,
[2575]36                   float &t,
37                   const float nearestT,
38                   Vector3 &normal) const
[1328]39{
[2575]40  //////////////
41  // specialised triangle ray casting version
42  // using ray-plane intersection
43 
44  // get triangle edge vectors and plane normal
45  const Vector3 u = mVertices[0] - mVertices[1];
46  const Vector3 v = mVertices[2] - mVertices[1];
47 
48  // cross product
49  normal = Normalize(CrossProd(v, u));
50 
51  // ray direction vector
52  const Vector3 dir = ray.GetDir();
53  const Vector3 w0 = ray.GetLoc() - mVertices[1];
54 
55  // params to calc ray-plane intersect
56  const float a = -DotProd(normal, w0);
57  const float b = DotProd(normal, dir);
58 
59  // check for division by zero
60  if (fabs(b) < Limits::Small)
61    {   
62      // ray is parallel to triangle plane
63      if (a == 0)
[1328]64        {
[2575]65          // ray lies in triangle plane
[1867]66          return Ray::INTERSECTION_OUT_OF_LIMITS;
[1328]67        }
[2575]68      else
[1344]69        {
[2575]70          // ray disjoint from plane
71          return Ray::NO_INTERSECTION;
[1344]72        }
[2575]73    }
[1328]74
[2575]75  // distance from origin of ray to plane
76  t = a / b;
77 
78  if (t <= Limits::Small) // ray goes away from triangle
79    {
80      return Ray::INTERSECTION_OUT_OF_LIMITS;
81    }
82  // already found nearer intersection
83  else if ((ray.GetType() == Ray::LOCAL_RAY) && (t >= nearestT))
84    {
85      return Ray::NO_INTERSECTION;
86    }
87 
88  /////////////////
[1344]89    //-- found intersection point
[2575]90  //-- check if it is inside triangle
[1344]91 
[2575]92  const Vector3 pt = ray.GetLoc() + t * dir;
[1715]93#if GTP_DEBUG
[2575]94  if (!pt.CheckValidity())
95    {
96      cout << "tr: " << *this << endl;
97      cout << "v: " << pt << " t: " << t << " a: " << a << " b: " << b << " n: " << normal << endl;
98    }
[1421]99#endif
[2575]100  const Vector3 w = pt - mVertices[1];
101 
102  const float uu = DotProd(u, u);
103  const float uv = DotProd(u, v);
104  const float vv = DotProd(v, v);
105 
106  const float wu = DotProd(w, u);
107  const float wv = DotProd(w, v);
108 
109 
110  const float D = uv * uv - uu * vv;
111 
112  // get and test parametric coords
113  const float s = (uv * wv - vv * wu) / D;
114 
115  if ((s < -Limits::Small) || (s > 1.0f + Limits::Small)) // pt is outside triangle
116    {   
117      return Ray::NO_INTERSECTION;
118    }
119 
120  const float s2 = (uv * wu - uu * wv) / D;
121 
122  if ((s2 < -Limits::Small) || ((s + s2) > 1.0f + Limits::Small)) // pt is outside triangle
123    {   
124      return Ray::NO_INTERSECTION;
125    }
126 
127  return Ray::INTERSECTION; // I is in T
128}
[1328]129
130
[2575]131int
132Triangle3::CastSimpleRay(const SimpleRay &ray,
133                         float &t,
134                         const float nearestT) const
135{
136  //////////////
137  // specialised triangle ray casting version
138  // using ray-plane intersection
139 
140  // get triangle edge vectors and plane normal
141  const Vector3 u = mVertices[0] - mVertices[1];
142  const Vector3 v = mVertices[2] - mVertices[1];
[1867]143
[2575]144  // cross product
145  const Vector3 normal = Normalize(CrossProd(v, u));
[1867]146
[2575]147  // ray direction vector
148  const Vector3 dir = ray.mDirection;
149  const Vector3 w0 = ray.mOrigin - mVertices[1];
[1328]150
[2575]151  // params to calc ray-plane intersect
152  const float a = -DotProd(normal, w0);
153  const float b = DotProd(normal, dir);
[1328]154
[2575]155  // check for division by zero
156  if (fabs(b) < Limits::Small)
157  {   
158    // ray is parallel to triangle plane
159    if (a == 0)
160    {
161      // ray lies in triangle plane
162      return Ray::INTERSECTION_OUT_OF_LIMITS;
163    }
164    else {
165      // ray disjoint from plane
166      return Ray::NO_INTERSECTION;
167    }
168  }
[1328]169
[2575]170  // distance from origin of ray to plane
171  t = a / b;
[1328]172
[2575]173  if (t <= Limits::Small) // ray goes away from triangle
174  {
175    return Ray::INTERSECTION_OUT_OF_LIMITS;
176  }
177
178  // already found nearer intersection
179  if (t > nearestT + 1e-5)
180  {
181    return Ray::NO_INTERSECTION;
182  }
183
184  /////////////////
185  //-- found intersection point
186  //-- check if it is inside triangle
187 
188  const Vector3 pt = ray.mOrigin + t * dir;
189#if GTP_DEBUG
190  if (!pt.CheckValidity())
191  {
192    cout << "tr: " << *this << endl;
193    cout << "v: " << pt << " t: " << t << " a: " << a << " b: " << b << " n: " << normal << endl;
194  }
195#endif
196  const Vector3 w = pt - mVertices[1];
197
198  const float uu = DotProd(u, u);
199  const float uv = DotProd(u, v);
200  const float vv = DotProd(v, v);
201   
202  const float wu = DotProd(w, u);
203  const float wv = DotProd(w, v);
204 
205 
206  const float D = uv * uv - uu * vv;
207 
208  // get and test parametric coords
209  const float s = (uv * wv - vv * wu) / D;
210 
211  if ((s < -Limits::Small) ||
212      (s > 1.0f + Limits::Small)) // pt is outside triangle
213  {     
214    return Ray::NO_INTERSECTION;
215  }
216
217  const float s2 = (uv * wu - uu * wv) / D;
218 
219  if ((s2 < -Limits::Small) ||
220      ((s + s2) > 1.0f + Limits::Small)) // pt is outside triangle
221  {     
222    return Ray::NO_INTERSECTION;
223  }
224 
225  return Ray::INTERSECTION; // I is in T
[1328]226}
[2575]227 
[1328]228
229AxisAlignedBox3 Triangle3::GetBoundingBox() const
230{
231        AxisAlignedBox3 box;
232        box.Initialize();
233
234        box.Include(mVertices[0]);
235        box.Include(mVertices[1]);
236        box.Include(mVertices[2]);
237
[2575]238        box.EnlargeToMinSize();
239       
[1328]240        return box;
241}
242
243
244Vector3 Triangle3::GetNormal() const
245{
246        const Vector3 v1 = mVertices[0] - mVertices[1];
[1344]247        const Vector3 v2 = mVertices[2] - mVertices[1];
[1328]248
249        return Normalize(CrossProd(v2, v1));
250}
251
252
253Vector3 Triangle3::GetCenter() const
254{
255        return (mVertices[0] + mVertices[1] + mVertices[2]) / 3.0f;
256}
257
258
259float Triangle3::GetArea() const
260{
[1344]261        Vector3 v1 = mVertices[0] - mVertices[1], v2=mVertices[2] - mVertices[1];
[1328]262        return 0.5f * Magnitude(CrossProd(v2, v1));
263}
264
265
[1420]266bool Triangle3::CheckValidity() const
267{
268        return !(
269                EpsilonEqualV3(mVertices[0], mVertices[1]) ||
270                EpsilonEqualV3(mVertices[0], mVertices[2]) ||
271                EpsilonEqualV3(mVertices[1], mVertices[2])
272                );
273}
274
[1932]275
[1933]276bool Triangle3::GetPlaneIntersection(const Plane3 &plane,
277                                                                         Vector3 &intersectA,
278                                                                         Vector3 &intersectB) const
[1932]279{
280        int side[3];
281
282        // compute distance from plane
283        for (int i = 0; i < 3; ++ i)
284        {
285                side[i] = plane.Side(mVertices[i], Limits::Small);
286        }
287
288        /////
289        // no intersection => early exit
290        if (((side[0] > 0) && (side[1] > 0) && (side[2] > 0)) ||
291                ((side[0] < 0) && (side[1] < 0) && (side[2] < 0)))
292        {
293                return false;
294        }
295
296        /////////////
297        // at least 2 triangle vertices lie in plane => early exit
298
299        for (int i = 0; i < 3; ++ i)
300        {
301                if (!side[i] && !side[(i + 1) % 3])
302                {
303                        intersectA = mVertices[i];
304                        intersectB = mVertices[(i + 1) % 3];
[1933]305                       
[1932]306                        return true;
307                }
308        }
309
310        bool foundA = false;
311
312        // compute intersection points
313        for (int i = 0; i < 3; ++ i)
314        {
315                const int i_2 = (i + 1) % 3;
316
317                // intersection found
[1933]318                if ((side[i] >= 0) && (side[i_2] <= 0) ||
319                        (side[i] <= 0) && (side[i_2] >= 0))
[1932]320                {       
321                        const float t = plane.FindT(mVertices[i], mVertices[i_2]);
322                       
323                        if (!foundA)
324                        {
325                                intersectA = mVertices[i] + t * (mVertices[i_2] - mVertices[i]);
[1933]326                       
[1932]327                                foundA = true;
328                        }
329                        else
330                        {
331                                intersectB = mVertices[i] + t * (mVertices[i_2] - mVertices[i]);
[1933]332                       
[1932]333                                return true;
334                        }
335                }
336        }
337
338        cout << "warning! wrong triangle - plane intersection" << endl;
339        return false; // something went wrong!
[1715]340}
[1932]341
342
343}
Note: See TracBrowser for help on using the repository browser.