Marlin: Updated bsp tree with ray intersection, faster generation

git-svn-id: http://www.openflipper.org/svnrepo/OpenFlipper/branches/Free@9435 383ad7c9-94d9-4d36-a494-682f7c89f535

common/UpdateType.hh

@@ -34,7 +34,7 @@ const UpdateType UPDATE_ALL(              UpdateTypeSet(1)       );
 /// This is the update identifier for global Object selections(source/target)
 const UpdateType UPDATE_OBJECT_SELECTION( UpdateTypeSet(1) << 1  );
 
-/// This is the update identifier for global Object visibiölity ( show/hide )
+/// This is the update identifier for global Object visibility ( show/hide )
 const UpdateType UPDATE_VISIBILITY(       UpdateTypeSet(1) << 2  );
 
 

common/bsp/BSPImplT.cc

@@ -98,8 +98,6 @@ _nearest(Node* _node, NearestNeighborData& _data) const
     }
   }
   
-
-  
   // non-terminal node
   else
   {
@@ -120,5 +118,66 @@ _nearest(Node* _node, NearestNeighborData& _data) const
   }
 }
 
+//-----------------------------------------------------------------------------
+
+template <class BSPCore>
+typename BSPImplT<BSPCore>::NearestNeighbor
+BSPImplT<BSPCore>::
+raycollision(const Point& _p, const Point& _r) const
+{
+  RayCollisionData  data;
+  data.ref  = _p;
+  data.dist = FLT_MAX;
+  data.ray  = _r;
+  
+  _raycollision(this->root_, data);
+  return NearestNeighbor(data.nearest, data.dist);
+}
+
+
+//-----------------------------------------------------------------------------
+
+
+template <class BSPCore>
+void
+BSPImplT<BSPCore>::
+_raycollision(Node* _node, RayCollisionData& _data) const
+{
+    // terminal node
+    if (!_node->left_child_)
+    {
+        Scalar dist;
+        Point v0, v1, v2;
+	Scalar u, v;
+
+        for (HandleIter it=_node->begin(); it!=_node->end(); ++it)
+        {
+            this->traits_.points(*it, v0, v1, v2);
+            if (ACG::Geometry::triangleIntersection(_data.ref, _data.ray, v0, v1, v2, dist, u, v)) {
+		
+                //face intersects with ray. But is it closer than any that we have found so far?
+                if (dist < _data.dist)
+                {
+                    _data.dist = dist;
+                    _data.nearest = *it;
+                }
+            }
+        }
+    }
+  
+  // non-terminal node
+  else
+    {
+        Scalar tmin, tmax;
+	bool used = false;
+        if ( _node->left_child_ && ACG::Geometry::axisAlignedBBIntersection( _data.ref, _data.ray, _node->left_child_->bb_min, _node->left_child_->bb_max, tmin, tmax)) {
+            _raycollision(_node->left_child_, _data);
+        }
+	if ( _node->right_child_ && ACG::Geometry::axisAlignedBBIntersection( _data.ref, _data.ray, _node->right_child_->bb_min, _node->right_child_->bb_max, tmin, tmax)) {
+            _raycollision(_node->right_child_, _data);
+        }
+    }
+}
+
 
 //=============================================================================
\ No newline at end of file

common/bsp/BSPImplT.hh

@@ -91,10 +91,11 @@ public: //---------------------------------------------------------------------
     Scalar  dist;
   };
 
-
-  /// Return handle of the nearest neighbor
+  /// Return handle of the nearest neighbor face
   NearestNeighbor nearest(const Point& _p) const;
   
+  /// Return handle of the nearest collided face
+  NearestNeighbor raycollision (const Point& _p, const Point& _r) const;
   
 
 private: //---------------------------------------------------------------------
@@ -108,9 +109,21 @@ private: //---------------------------------------------------------------------
     Handle  nearest;
   };
   
+  /// Store ray collide information
+  struct RayCollisionData
+  {
+    Point   ref;
+    Point   ray;
+    Scalar  dist;
+    Handle  nearest;
+  };
+
 
   // Recursive part of nearest()
   void _nearest(Node* _node, NearestNeighborData& _data) const;
+  
+  //resursive part of raycollide()
+  void _raycollision(Node* _node, RayCollisionData& _data) const;
 };
 
 

common/bsp/TriangleBSPCoreT.cc

@@ -59,38 +59,6 @@
 
 //== CLASS DEFINITION =========================================================
 
-
-template <class BSPTraits>
-TriangleBSPCoreT<BSPTraits>::Node::
-Node(const Handles& _handles, Node* _parent)
-  : handles_(_handles), 
-    parent_(_parent), left_child_(0), right_child_(0)
-{}
-
-
-//-----------------------------------------------------------------------------
-
-
-template <class BSPTraits>
-TriangleBSPCoreT<BSPTraits>::Node::
-~Node() 
-{ 
-  delete left_child_; 
-  delete right_child_;
-
-  if (parent_) 
-  {
-    if (this == parent_->left_child_) 
-      parent_->left_child_ = 0;
-    else 
-      parent_->right_child_ = 0;
-  }
-}
-
-
-//-----------------------------------------------------------------------------
-
-
 template <class BSPTraits>
 void
 TriangleBSPCoreT<BSPTraits>::
@@ -103,8 +71,12 @@ build(unsigned int _max_handles, unsigned int _max_depth)
   // delete own handles (don't store them twice)
   handles_ = Handles();
 
+  nodes=1;
+  traits_.calculateBoundingBoxRoot (root_);
   // call recursive helper
   _build(root_, _max_handles, _max_depth);
+  
+  std::cout << "Number of nodes:" << nodes << std::endl;
 }
 
 
@@ -122,43 +94,22 @@ _build(Node*         _node,
   if ((_depth == 0) || ((_node->end()-_node->begin()) <= (int)_max_handles))
     return;
   
-
-
-  // compute bounding box
-  HandleIter it;
-  Point p0, p1, p2;
-  Point bb_min; bb_min.vectorize(FLT_MAX);
-  Point bb_max; bb_max.vectorize(-FLT_MAX);
-  for (it=_node->begin(); it!=_node->end(); ++it)
-  {
-    traits_.points(*it, p0, p1, p2);
-    bb_min.minimize(p0);
-    bb_min.minimize(p1);
-    bb_min.minimize(p2);
-    bb_max.maximize(p0);
-    bb_max.maximize(p1);
-    bb_max.maximize(p2);
-  }
-
-
-  // split longest side of bounding box
-  Point   bb     = bb_max - bb_min;
-  Scalar  length = bb[0];
-  int     axis   = 0;
-  if (bb[1] > length) length = bb[ (axis=1) ];
-  if (bb[2] > length) length = bb[ (axis=2) ];
-
+  Point median;
+  int axis;
+  // compute bounding boxes for children
+  traits_.calculateBoundingBox (_node, median, axis);
   
   // construct splitting plane
   const Point XYZ[3] = { Point(1,0,0), Point(0,1,0), Point(0,0,1) };
-  _node->plane_ = Plane((bb_min+bb_max)*0.5, XYZ[axis]);
-
+  _node->plane_ = Plane(median, XYZ[axis]);
 
   // partition for left and right child
   Handles lhandles, rhandles;
   lhandles.reserve(_node->handles_.size()/2);
   rhandles.reserve(_node->handles_.size()/2);
 
+  HandleIter it;
+  Point p0, p1, p2;
   bool left, right;
   for (it=_node->begin(); it!=_node->end(); ++it)
   {
@@ -176,7 +127,6 @@ _build(Node*         _node,
     if (right) rhandles.push_back(*it);
   }
 
-
   // check it
   if (lhandles.size() == _node->handles_.size() ||
       rhandles.size() == _node->handles_.size())
@@ -188,8 +138,25 @@ _build(Node*         _node,
   // create children
   _node->left_child_  = new Node(lhandles, _node);  lhandles = Handles();
   _node->right_child_ = new Node(rhandles, _node);  rhandles = Handles();
-
-
+  nodes+=2;
+  
+  //save bounding boxes for children
+  /*
+  _node->left_child_->bb_min  = _node->bb_min;
+  _node->left_child_->bb_max  = _node->bb_max;
+  _node->left_child_->bb_max[axis] = median [axis];
+  
+  _node->right_child_->bb_min = _node->bb_min;
+  _node->right_child_->bb_min[axis] = median [axis];
+  _node->right_child_->bb_max = _node->bb_max;
+  */
+  _node->right_child_->bb_min  = _node->bb_min;
+  _node->right_child_->bb_max  = _node->bb_max;
+  _node->right_child_->bb_max[axis] = median [axis];
+  
+  _node->left_child_->bb_min = _node->bb_min;
+  _node->left_child_->bb_min[axis] = median [axis];
+  _node->left_child_->bb_max = _node->bb_max;
 
   // recurse to childen
   _build(_node->left_child_,  _max_handles, _depth-1);

common/bsp/TriangleBSPCoreT.hh

@@ -60,6 +60,8 @@
 #include <OpenMesh/Core/Geometry/VectorT.hh>
 #include <vector>
 
+#include "TriangleBSPT.hh"
+
 
 //== CLASS DEFINITION =========================================================
 
@@ -72,13 +74,13 @@ public: //---------------------------------------------------------------------
   typedef BSPTraits                      Traits;
   typedef typename BSPTraits::Point      Point;
   typedef typename BSPTraits::Handle     Handle;
+  typedef typename BSPTraits::Node	 Node;
   typedef typename Point::value_type     Scalar;
   typedef ACG::Geometry::PlaneT<Scalar>  Plane;
   typedef std::vector<Handle>            Handles;
   typedef typename Handles::iterator     HandleIter;
 
 
-
 public: //---------------------------------------------------------------------
 
 
@@ -98,26 +100,6 @@ public: //---------------------------------------------------------------------
   void build(unsigned int _max_handles, unsigned int _max_depth);
 
 
-
-protected: //-------------------------------------------------------------------
-
-
-  // Node of the tree: contains parent, children and splitting plane
-  struct Node
-  {
-    Node(const Handles& _handles, Node* _parent);
-    ~Node();
-
-    HandleIter begin() { return handles_.begin(); }
-    HandleIter end()   { return handles_.end();   }
-
-    Handles     handles_;
-    Node        *parent_, *left_child_, *right_child_;
-    Plane       plane_;
-  };
-
-  
-
 private: //---------------------------------------------------------------------
 
 
@@ -135,6 +117,8 @@ protected: //-------------------------------------------------------------------
   BSPTraits  traits_;
   Handles    handles_;
   Node*      root_;
+  int	nodes;
+  
 };
 
 

common/bsp/TriangleBSPT.hh

@@ -55,6 +55,7 @@
 
 //== INCLUDES =================================================================
 
+#include "BSPTreeNode.hh"
 #include "TriangleBSPCoreT.hh"
 #include "BSPImplT.hh"
 //#include "Distance.hh"
@@ -62,7 +63,7 @@
 
 
 //== CLASS DEFINITION =========================================================
-
+#include <list>
 
 template <class BSPTraits>
 class TriangleBSPT : public BSPImplT< TriangleBSPCoreT<BSPTraits> >
@@ -72,10 +73,8 @@ public:
   TriangleBSPT(const BSPTraits& _traits) : Base(_traits) {}
 };
 
-
 //== CLASS DEFINITION =========================================================
 
-
 template <class Mesh>
 class OpenMeshTriangleBSPTraits
 {
@@ -84,27 +83,105 @@ public:
     typedef typename Mesh::Scalar	Scalar;
     typedef typename Mesh::Point	Point;
     typedef typename Mesh::FaceHandle	Handle;
+    typedef std::vector<Handle>		Handles;
+    typedef typename Handles::iterator	HandleIter;
+    typedef TreeNode<Mesh>		Node;    
      
     OpenMeshTriangleBSPTraits(const Mesh& _mesh) : mesh_(_mesh) {}
 
-  void points(Handle _h, Point& _p0, Point& _p1, Point& _p2) const 
+    /// Returns the points belonging to the face handle _h
+    inline void points(const Handle _h, Point& _p0, Point& _p1, Point& _p2) const
     {
         typename Mesh::CFVIter fv_it(mesh_.cfv_iter(_h));
-    _p0 = mesh_.point(fv_it); ++fv_it;
-    _p1 = mesh_.point(fv_it); ++fv_it;
+        _p0 = mesh_.point(fv_it);
+        ++fv_it;
+        _p1 = mesh_.point(fv_it);
+        ++fv_it;
         _p2 = mesh_.point(fv_it);
     }
 
-  Scalar sqrdist(Handle _h, const Point& _p) const 
+    Scalar sqrdist(const Handle _h, const Point& _p) const
     {
         Point p0, p1, p2, q;
         points(_h, p0, p1, p2);
         return ACG::Geometry::distPointTriangleSquared(_p, p0, p1, p2, q);
     }
     
+    void calculateBoundingBox(
+        Node* _node,
+        Point& median,
+	int& axis)
+    {	
+	//determine splitting axis
+	HandleIter it_h;
+        Point p0, p1, p2, bb_min, bb_max;
+        bb_min.vectorize(FLT_MAX);
+        bb_max.vectorize(-FLT_MAX);
+	std::list<Point>* vertices = new std::list<Point>;
+	
+        for (it_h=_node->begin(); it_h!=_node->end(); ++it_h)
+        {
+            points(*it_h, p0, p1, p2);
+            bb_min.minimize(p0);
+            bb_min.minimize(p1);
+            bb_min.minimize(p2);
+            bb_max.maximize(p0);
+            bb_max.maximize(p1);
+            bb_max.maximize(p2);
+	    vertices->push_back (p0);
+	    vertices->push_back (p1);
+	    vertices->push_back (p2);
+        }
+	
+        // split longest side of bounding box
+        Point   bb     = bb_max - bb_min;
+        Scalar  length = bb[0];
+	axis = 0;
+        if (bb[1] > length) length = bb[ (axis=1) ];
+        if (bb[2] > length) length = bb[ (axis=2) ];
+	
+	//calculate the median value in axis-direction
+	switch (axis) {
+	  case 0: vertices->sort (x_sort() ); break;
+	  case 1: vertices->sort (y_sort() ); break;
+	  case 2: vertices->sort (z_sort() ); break; }
+	vertices->unique(); ///todo: does this work with Points?!
+	
+	int size = vertices->size();
+	typename std::list<Point>::iterator it_v;
+	it_v = vertices->begin();
+	std::advance(it_v, size/2);
+	median = *it_v;
+    }
+
+    void calculateBoundingBoxRoot(
+        Node* _node)
+    {
+        HandleIter it;
+        Point p0, p1, p2, bb_min, bb_max;
+        bb_min.vectorize(FLT_MAX);
+        bb_max.vectorize(-FLT_MAX);
+        for (it=_node->begin(); it!=_node->end(); ++it)
+        {
+            points(*it, p0, p1, p2);
+            bb_min.minimize(p0);
+            bb_min.minimize(p1);
+            bb_min.minimize(p2);
+            bb_max.maximize(p0);
+            bb_max.maximize(p1);
+            bb_max.maximize(p2);
+        }
+        _node->bb_min = bb_min;
+        _node->bb_max = bb_max;
+    }
+
 private:
 
     const Mesh& mesh_;
+    //functors for sorting in different directions
+    struct x_sort { bool operator()(const Point& first, const Point& second) { return (first[0] < second[0]); }  };
+    struct y_sort { bool operator()(const Point& first, const Point& second) { return (first[1] < second[1]); }  };
+    struct z_sort { bool operator()(const Point& first, const Point& second) { return (first[2] < second[2]); }  };
 };
 
 
@@ -121,7 +198,6 @@ public:
   OpenMeshTriangleBSPT(const Mesh& _mesh) : Base(Traits(_mesh)) {}
 };
 
-
 //=============================================================================
 #endif // MB_TRIANGLEBSP_HH defined
 //=============================================================================