source: GTP/trunk/Lib/Vis/Preprocessing/src/Mesh.cpp @ 878

Revision 878, 13.7 KB checked in by bittner, 19 years ago (diff)

mesh inntersection bug fixed

Line 
1#include "Ray.h"
2#include "Mesh.h"
3#include "MeshKdTree.h"
4#include "Triangle3.h"
5
6namespace GtpVisibilityPreprocessor {
7
8bool MeshDebug = false;
9
10int Intersectable::sMailId = 21843194198;
11int Intersectable::sReservedMailboxes = 1;
12
13struct SortableVertex {
14
15  Vector3 vertex;
16
17  int originalId;
18  int newId;
19  int finalPos;
20
21  SortableVertex() {}
22
23  SortableVertex(const Vector3 &v,
24                                 const int id):
25        vertex(v),
26        originalId(id),
27        newId(id)
28  {}
29 
30  friend bool operator<(const SortableVertex &a,
31                                                const SortableVertex &b)
32  {
33        if (a.vertex.x < b.vertex.x)
34          return true;
35        else
36          if (a.vertex.x > b.vertex.x)
37                return false;
38       
39        if (a.vertex.y < b.vertex.y)
40          return true;
41        else
42          if (a.vertex.y > b.vertex.y)
43                return false;
44       
45        if (a.vertex.z < b.vertex.z)
46          return true;
47        else
48          //      if (a.z > b.z)
49          return false;
50       
51        //      return false;
52  }
53
54};
55
56
57void
58Mesh::ComputeBoundingBox()
59{
60
61  mBox.Initialize();
62  VertexContainer::const_iterator vi = mVertices.begin();
63  for (; vi != mVertices.end(); vi++) {
64    mBox.Include(*vi);
65  }
66//mBox.Enlarge(1e-4f);
67}
68
69void
70Mesh::Preprocess()
71{
72        Cleanup();
73   
74        ComputeBoundingBox();
75   
76        /** true if it is a watertight convex mesh
77        */
78        mIsConvex = false;
79
80        if (mFaces.size() > MeshKdTree::mTermMinCost)
81        {
82                mKdTree = new MeshKdTree(this);
83                MeshKdLeaf *root = (MeshKdLeaf *)mKdTree->GetRoot();
84               
85                for (int i = 0; i < mFaces.size(); i++)
86                        root->mFaces.push_back(i);
87               
88                cout<<"KD";
89                mKdTree->Construct();
90
91                if (mKdTree->GetRoot()->IsLeaf())
92                {
93                        cout<<"d";
94                        delete mKdTree;
95                        mKdTree = NULL;
96                }
97        }
98}
99
100
101void
102Mesh::IndexVertices()
103{
104  int i;
105  // check whether the vertices can be simplfied and reindexed
106  vector<SortableVertex> svertices(mVertices.size());
107
108  for (i=0; i < mVertices.size(); i++)
109        svertices[i] = SortableVertex(mVertices[i], i);
110
111  sort(svertices.begin(), svertices.end());
112
113  for (i=0; i < svertices.size() - 1; i++)
114        if (svertices[i].vertex == svertices[i+1].vertex)
115          svertices[i+1].newId = svertices[i].newId;
116
117  // remove the same vertices
118  int k = 0;
119  mVertices[0] = svertices[0].vertex;
120  svertices[0].finalPos = 0;
121 
122  for (i=1; i < svertices.size(); i++) {
123        if (svertices[i].newId != svertices[i-1].newId)
124          k++;
125       
126        mVertices[k] = svertices[i].vertex;
127        svertices[i].finalPos = k;
128  }
129
130  mVertices.resize(k + 1);
131 
132  vector<int> remapBuffer(svertices.size());
133  for (i = 0; i < svertices.size(); i++)
134        remapBuffer[svertices[i].originalId] = svertices[i].finalPos;
135 
136  // remap all faces
137 
138  for (int faceIndex = 0; faceIndex < mFaces.size(); faceIndex++) {
139        Face *face = mFaces[faceIndex];
140        for (int i = 0; i < face->mVertexIndices.size(); i++) {
141          face->mVertexIndices[i] = remapBuffer[face->mVertexIndices[i]];
142        }
143  }
144}
145
146AxisAlignedBox3
147Mesh::GetFaceBox(const int faceIndex)
148{
149  Face *face = mFaces[faceIndex];
150  AxisAlignedBox3 box;
151  box.SetMin( mVertices[face->mVertexIndices[0]] );
152  box.SetMax(box.Min());
153  for (int i = 1; i < face->mVertexIndices.size(); i++) {
154    box.Include(mVertices[face->mVertexIndices[i]]);
155  }
156  return box;
157}
158
159int
160Mesh::CastRayToFace(
161                                        const int faceIndex,
162                                        Ray &ray,
163                                        float &nearestT,
164                                        int &nearestFace,
165                                        Intersectable *instance
166                                        )
167{
168  float t;
169  int hit = 0;
170  if (RayFaceIntersection(faceIndex, ray, t, nearestT) == Ray::INTERSECTION) {
171    switch (ray.GetType()) {
172    case Ray::GLOBAL_RAY:
173      ray.intersections.push_back(Ray::Intersection(t, instance, faceIndex));
174      hit++;
175      break;
176    case Ray::LOCAL_RAY:
177      nearestT = t;
178      nearestFace = faceIndex;
179      hit++;
180      break;
181    case Ray::LINE_SEGMENT:
182      if (t <= 1.0f) {
183                ray.intersections.push_back(Ray::Intersection(t, instance, faceIndex));
184                hit++;
185      }
186      break;
187    }
188  }
189  return hit;
190}
191
192int
193Mesh::CastRay(
194                          Ray &ray,
195                          MeshInstance *instance
196                          )
197{
198  if (mKdTree) {
199    return mKdTree->CastRay(ray, instance);
200  }
201 
202  int faceIndex = 0;
203  int hits = 0;
204  float nearestT = MAX_FLOAT;
205  int nearestFace = -1;
206 
207  if (ray.GetType() == Ray::LOCAL_RAY && ray.intersections.size())
208    nearestT = ray.intersections[0].mT;
209
210       
211  for ( ;
212                faceIndex < mFaces.size();
213                faceIndex++) {
214    hits += CastRayToFace(faceIndex, ray, nearestT, nearestFace, instance);
215    if (mIsConvex && nearestFace != -1)
216      break;
217  }
218 
219  if ( hits && ray.GetType() == Ray::LOCAL_RAY ) {
220    if (ray.intersections.size())
221      ray.intersections[0] = Ray::Intersection(nearestT, instance, nearestFace);
222    else
223      ray.intersections.push_back(Ray::Intersection(nearestT, instance, nearestFace));
224  }
225 
226  return hits;
227}
228
229int
230Mesh::CastRayToSelectedFaces(
231                                                         Ray &ray,
232                                                         const vector<int> &faces,
233                                                         Intersectable *instance
234                                                         )
235{
236  vector<int>::const_iterator fi;
237  int faceIndex = 0;
238  int hits = 0;
239  float nearestT = MAX_FLOAT;
240  int nearestFace = -1;
241
242  if (ray.GetType() == Ray::LOCAL_RAY && ray.intersections.size())
243    nearestT = ray.intersections[0].mT;
244
245  for ( fi = faces.begin();
246                fi != faces.end();
247                fi++) {
248    hits += CastRayToFace(*fi, ray, nearestT, nearestFace, instance);
249    if (mIsConvex && nearestFace != -1)
250      break;
251  }
252 
253  if ( hits && ray.GetType() == Ray::LOCAL_RAY ) {
254    if (ray.intersections.size())
255      ray.intersections[0] = Ray::Intersection(nearestT, instance, nearestFace);
256    else
257      ray.intersections.push_back(Ray::Intersection(nearestT, instance, nearestFace));
258  }
259 
260  return hits;
261}
262
263
264// int_lineseg returns 1 if the given line segment intersects a 2D
265// ray travelling in the positive X direction.  This is used in the
266// Jordan curve computation for polygon intersection.
267inline int
268int_lineseg(float px,
269                        float py,
270                        float u1,
271                        float v1,
272                        float u2,
273                        float v2)
274{
275  float ydiff;
276
277  u1 -= px; u2 -= px;     // translate line
278  v1 -= py; v2 -= py;
279
280  if ((v1 > 0 && v2 > 0) ||
281      (v1 < 0 && v2 < 0) ||
282      (u1 < 0 && u2 < 0))
283    return 0;
284
285  if (u1 > 0 && u2 > 0)
286    return 1;
287
288  ydiff = v2 - v1;
289  if (fabs(ydiff) < Limits::Small) {      // denominator near 0
290    if (((fabs(v1) > Limits::Small) ||
291                 (u1 > 0) || (u2 > 0)))
292      return 0;
293    return 1;
294  }
295 
296  double t = -v1 / ydiff;                 // Compute parameter
297
298  double thresh;
299
300  if (ydiff < 0.0f)
301        thresh = -1e-20;
302  else
303        thresh = 1e-20;
304
305 
306  return (u1 + t * (u2 - u1)) > thresh; //-Limits::Small;
307}
308
309
310
311// intersection with the polygonal face of the mesh
312int
313Mesh::RayFaceIntersection(const int faceIndex,
314                                                  const Ray &ray,
315                                                  float &t,
316                                                  const float nearestT
317                                                  )
318{
319  Face *face  = mFaces[faceIndex];
320
321  Plane3 plane = GetFacePlane(faceIndex);
322  float dot = DotProd(plane.mNormal, ray.GetDir());
323
324  if (MeshDebug) {
325        cout<<endl<<endl;
326        cout<<"normal="<<plane.mNormal<<endl;
327        cout<<"dot="<<dot<<endl;
328  }
329 
330       
331  // Watch for near-zero denominator
332  // ONLY single sided polygons!!!!!
333  if (ray.mFlags & Ray::CULL_BACKFACES) {
334        if (dot > -Limits::Small)
335          //  if (fabs(dot) < Limits::Small)
336          return Ray::NO_INTERSECTION;
337  } else {
338        if (fabs(dot) < Limits::Small)
339          return Ray::NO_INTERSECTION;
340  }
341 
342  t = (-plane.mD - DotProd(plane.mNormal, ray.GetLoc())) / dot;
343
344  if (MeshDebug)
345        cout<<"t="<<t<<endl;
346
347  if (t <= Limits::Small)
348    return Ray::INTERSECTION_OUT_OF_LIMITS;
349 
350  if (t >= nearestT) {
351    return Ray::INTERSECTION_OUT_OF_LIMITS; // no intersection was found
352  }
353 
354  int count = 0;
355  float u, v, u1, v1, u2, v2;
356  int i;
357
358  int paxis = plane.mNormal.DrivingAxis();
359
360 
361  // Project the intersection point onto the coordinate plane
362  // specified by which.
363  ray.Extrap(t).ExtractVerts(&u, &v, paxis);
364
365
366  int size = (int)face->mVertexIndices.size();
367
368  if (MeshDebug)
369        cout<<"size="<<size<<endl;
370 
371  mVertices[face->mVertexIndices[size - 1]].
372    ExtractVerts(&u1, &v1, paxis );
373 
374  //$$JB changed 12.4.2006 from 0 ^^
375  if (0 && size <= 4) {
376    // assume a convex face
377    for (i = 0; i < size; i++) {
378      mVertices[face->mVertexIndices[i]].ExtractVerts(&u2, &v2, paxis);
379      // line u1, v1, u2, v2
380      if ((v1 - v2)*(u - u1) + (u2 - u1)*(v - v1) > 0) {
381                if (MeshDebug)
382                  cout<<"exit on "<<i<<endl;
383                return Ray::NO_INTERSECTION;
384          }
385      u1 = u2;
386      v1 = v2;
387    }
388   
389    return Ray::INTERSECTION;
390  }
391 
392  // We're stuck with the Jordan curve computation.  Count number
393  // of intersections between the line segments the polygon comprises
394  // with a ray originating at the point of intersection and
395  // travelling in the positive X direction.
396  for (i = 0; i < size; i++) {
397    mVertices[face->mVertexIndices[i]].ExtractVerts(&u2, &v2, paxis);
398    count += (int_lineseg(u, v, u1, v1, u2, v2) != 0);
399    u1 = u2;
400    v1 = v2;
401  }
402
403  if (MeshDebug)
404        cout<<"count="<<count<<endl;
405
406  // We hit polygon if number of intersections is odd.
407  return (count & 1) ? Ray::INTERSECTION : Ray::NO_INTERSECTION;
408}
409
410int
411Mesh::GetRandomSurfacePoint(Vector3 &point, Vector3 &normal)
412{
413  int faceIndex = (int)RandomValue(0, (Real)((int)mFaces.size()-1));
414 
415  // assume the face is convex and generate a convex combination
416  //
417  Face *face = mFaces[faceIndex];
418  point = Vector3(0,0,0);
419  float sum = 0.0f;
420  for (int i = 0; i < face->mVertexIndices.size(); i++) {
421    float r = RandomValue(0,1);
422    sum += r;
423    point += mVertices[face->mVertexIndices[i]]*r;
424  }
425  point *= 1.0f/sum;
426       
427        normal = GetFacePlane(faceIndex).mNormal;
428
429        return faceIndex;
430}
431
432int
433Mesh::GetRandomVisibleSurfacePoint(Vector3 &point,
434                                                                                                                                         Vector3 &normal,
435                                                                                                                                         const Vector3 &viewpoint,
436                                                                                                                                         const int maxTries
437                                                                                                                                         )
438{
439  Plane3 plane;
440  int faceIndex = (int)RandomValue(0, (Real)((int)mFaces.size()-1));
441        int tries;
442  for (tries = 0; tries < maxTries; tries++) {
443    Face *face = mFaces[faceIndex];
444    plane = GetFacePlane(faceIndex);
445   
446    if (plane.Side(viewpoint) > 0) {
447      point = Vector3(0,0,0);
448      float sum = 0.0f;
449      // pickup a point inside this triangle
450      for (int i = 0; i < face->mVertexIndices.size(); i++) {
451                                float r = RandomValue(0,1);
452                                sum += r;
453                                point += mVertices[face->mVertexIndices[i]]*r;
454      }
455      point *= 1.0f/sum;
456      break;
457    }
458  }
459 
460  normal = plane.mNormal;
461  return (tries < maxTries) ? faceIndex + 1 : 0;
462}
463
464
465int
466MeshInstance::CastRay(
467                                          Ray &ray
468                                          )
469{
470  int res = mMesh->CastRay(ray, this);
471  return res;
472}
473
474int
475MeshInstance::CastRay(
476                                          Ray &ray,
477                                          const vector<int> &faces
478                                          )
479{
480  return mMesh->CastRayToSelectedFaces(ray, faces, this);
481}
482
483
484
485int
486MeshInstance::GetRandomSurfacePoint(Vector3 &point, Vector3 &normal)
487{
488  return mMesh->GetRandomSurfacePoint(point, normal);
489}
490
491int
492MeshInstance::GetRandomVisibleSurfacePoint(Vector3 &point,
493                                                                                                                                                                         Vector3 &normal,
494                                                                                                                                                                         const Vector3 &viewpoint,
495                                                                                                                                                                         const int maxTries
496                                                                                                                                                                         )
497{
498        return mMesh->GetRandomVisibleSurfacePoint(point, normal, viewpoint, maxTries);
499}
500
501
502int
503TransformedMeshInstance::GetRandomSurfacePoint(Vector3 &point, Vector3 &normal)
504{
505  int index = mMesh->GetRandomSurfacePoint(point, normal);
506  point = mWorldTransform*point;
507  normal = TransformNormal(mWorldTransform, normal);
508  return index;
509}
510
511Plane3
512Mesh::GetFacePlane(const int faceIndex)
513{
514  Face *face = mFaces[faceIndex];
515        return Plane3(mVertices[face->mVertexIndices[0]],
516                                                                mVertices[face->mVertexIndices[1]],
517                                                                mVertices[face->mVertexIndices[2]]);
518}
519
520bool
521Mesh::ValidateFace(const int i)
522{
523        Face *face = mFaces[i];
524
525        Plane3 plane = Plane3(mVertices[face->mVertexIndices[0]],
526                                                                                                mVertices[face->mVertexIndices[1]],
527                                                                                                mVertices[face->mVertexIndices[2]]);
528       
529        if (!eq(Magnitude(plane.mNormal), 1.0f))
530                return false;
531
532        return true;
533}
534
535void
536Mesh::Cleanup()
537{
538        int toRemove = 0;
539        FaceContainer newFaces;
540        for (int i=0; i < mFaces.size(); i++)
541                if (ValidateFace(i)) {
542                        newFaces.push_back(mFaces[i]);
543                } else {
544                        cout<<"d";
545                        delete mFaces[i];
546                        toRemove++;
547                }
548       
549        if (toRemove) {
550                mFaces = newFaces;
551        }
552       
553        // cleanup vertices??
554}
555
556int
557TransformedMeshInstance::CastRay(
558                                                                 Ray &ray
559                                                                 )
560{
561  ray.ApplyTransform(Invert(mWorldTransform));
562  int res = mMesh->CastRay(ray, this);
563  ray.ApplyTransform(mWorldTransform);
564 
565  return res;
566}
567
568
569void
570Mesh::AddTriangle(const Triangle3 &triangle)
571{
572  int index = (int)mVertices.size();
573
574  for (int i=0; i < 3; i++) {
575    mVertices.push_back(triangle.mVertices[i]);
576  }
577 
578  AddFace(new Face(index + 0, index + 1, index + 2) );
579}
580
581void
582Mesh::AddRectangle(const Rectangle3 &rect)
583{
584  int index = (int)mVertices.size();
585
586  for (int i=0; i < 4; i++) {
587    mVertices.push_back(rect.mVertices[i]);
588  }
589 
590  AddFace(new Face(index + 0, index + 1, index + 2, index + 3) );
591}
592
593void
594Mesh::AssignRandomMaterial()
595{
596  if (!mMaterial)
597        mMaterial = new Material;
598  *mMaterial = RandomMaterial();
599
600}
601
602
603Mesh *CreateBox(const AxisAlignedBox3 &box)
604{
605        Mesh *mesh = new Mesh;
606  // add 8 vertices of the box
607  int index = (int)mesh->mVertices.size();
608  for (int i=0; i < 8; i++) {
609    Vector3 v;
610    box.GetVertex(i, v);
611    mesh->mVertices.push_back(v);
612  }
613 
614  mesh->AddFace(new Face(index + 0, index + 1, index + 3, index + 2) );
615  mesh->AddFace(new Face(index + 0, index + 2, index + 6, index + 4) );
616  mesh->AddFace(new Face(index + 4, index + 6, index + 7, index + 5) );
617 
618  mesh->AddFace(new Face(index + 3, index + 1, index + 5, index + 7) );
619  mesh->AddFace(new Face(index + 0, index + 4, index + 5, index + 1) );
620  mesh->AddFace(new Face(index + 2, index + 3, index + 7, index + 6) );
621 
622  return mesh;
623}
624
625}
Note: See TracBrowser for help on using the repository browser.