Index: /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCell.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCell.cpp	(revision 366)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCell.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Polygon3.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Polygon3.h	(revision 362)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Polygon3.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/SceneGraph.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/SceneGraph.h	(revision 339)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/SceneGraph.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.h	(revision 367)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/AxisAlignedBox3.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/AxisAlignedBox3.h	(revision 318)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/AxisAlignedBox3.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Vector3.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Vector3.h	(revision 349)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Vector3.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Makefile
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Makefile	(revision 369)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Makefile	(revision 372)
@@ -1,5 +1,5 @@
 #############################################################################
 # Makefile for building: preprocessor
-# Generated by qmake (1.07a) (Qt 3.3.2) on: Wed Nov 02 15:58:09 2005
+# Generated by qmake (1.07a) (Qt 3.3.2) on: Thu Nov 03 11:34:08 2005
 # Project:  preprocessor.pro
 # Template: app
@@ -39,5 +39,7 @@
 ####### Files
 
-HEADERS =	Halton.h
+HEADERS =	Halton.h \
+		VssRay.h \
+		VssPreprocessor.h
 SOURCES =	Preprocessor.cpp \
 		SamplingPreprocessor.cpp \
@@ -69,5 +71,8 @@
 		ViewCell.cpp \
 		ViewCellBsp.cpp \
-		Halton.cpp
+		Halton.cpp \
+		VssRay.cpp \
+		VssTree.cpp \
+		VssPreprocessor.cpp
 OBJECTS =	Preprocessor.obj \
 		SamplingPreprocessor.obj \
@@ -99,5 +104,8 @@
 		ViewCell.obj \
 		ViewCellBsp.obj \
-		Halton.obj
+		Halton.obj \
+		VssRay.obj \
+		VssTree.obj \
+		VssPreprocessor.obj
 FORMS =	
 UICDECLS =	
@@ -191,4 +199,7 @@
 	-$(DEL_FILE) ViewCellBsp.obj
 	-$(DEL_FILE) Halton.obj
+	-$(DEL_FILE) VssRay.obj
+	-$(DEL_FILE) VssTree.obj
+	-$(DEL_FILE) VssPreprocessor.obj
 
 
@@ -614,4 +625,53 @@
 		
 
+VssRay.obj: VssRay.cpp  \
+		VssRay.h \
+		AxisAlignedBox3.h \
+		Vector3.h \
+		common.h \
+		Rectangle3.h \
+		Matrix4x4.h \
+		Plane3.h \
+		
+
+VssTree.obj: VssTree.cpp  \
+		VssTree.h \
+		Environment.h \
+		VssRay.h \
+		AxisAlignedBox3.h \
+		Statistics.h \
+		Vector3.h \
+		common.h \
+		Rectangle3.h \
+		Matrix4x4.h \
+		Plane3.h \
+		
+
+VssPreprocessor.obj: VssPreprocessor.cpp  \
+		SceneGraph.h \
+		KdTree.h \
+		VssPreprocessor.h \
+		X3dExporter.h \
+		Environment.h \
+		MutualVisibility.h \
+		Polygon3.h \
+		ViewCell.h \
+		Containers.h \
+		AxisAlignedBox3.h \
+		Rectangle3.h \
+		Matrix4x4.h \
+		Vector3.h \
+		Plane3.h \
+		common.h \
+		Ray.h \
+		Pvs.h \
+		Preprocessor.h \
+		Mesh.h \
+		ViewCellBsp.h \
+		Intersectable.h \
+		Material.h \
+		Exporter.h \
+		
+
 ####### Install
 
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/VssTree.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/VssTree.h	(revision 372)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/VssTree.h	(revision 372)
@@ -0,0 +1,670 @@
+// ================================================================
+// $Id$
+// ****************************************************************
+// 
+// Initial coding by
+/**
+   @author Jiri Bittner
+*/
+
+#ifndef __VSSTREE_H__
+#define __VSSTREE_H__
+
+// Standard headers
+#include <iomanip>
+#include <vector>
+#include <functional>
+#include <stack>
+
+
+// User headers
+#include "VssRay.h"
+#include "AxisAlignedBox3.h"
+
+
+#define USE_KDNODE_VECTORS 1
+#define _RSS_STATISTICS
+#define _RSS_TRAVERSAL_STATISTICS
+
+
+#include "Statistics.h"
+
+// --------------------------------------------------------------
+// Static statistics for kd-tree search
+// --------------------------------------------------------------
+class VssStatistics :
+  public StatisticsBase
+{
+public:  
+  // total number of nodes
+  int nodes;
+  // number of splits along each of the axes
+  int splits[7];
+  // totals number of rays
+  int rays;
+  // total number of query domains
+  int queryDomains;
+  // total number of ray references
+  int rayRefs;
+  // refs in non empty leafs
+  int rayRefsNonZeroQuery;
+  // total number of query references
+  int queryDomainRefs;
+  // nodes with zero queries
+  int zeroQueryNodes;
+  // max depth nodes
+  int maxDepthNodes;
+  // max depth nodes
+  int minCostNodes;
+  // max number of rays per node
+  int maxRayRefs;
+  // number of dynamically added ray refs
+  int addedRayRefs;
+  // number of dynamically removed ray refs
+  int removedRayRefs;
+  
+  // Constructor
+  VssStatistics() {
+    Reset();
+  }
+
+  int Nodes() const {return nodes;}
+  int Interior() const { return nodes/2; }
+  int Leaves() const { return (nodes/2) + 1; }
+
+  void Reset() {
+    nodes = 0;
+    for (int i=0; i<7; i++)
+      splits[i] = 0;
+    rays = queryDomains = 0;
+    rayRefs = rayRefsNonZeroQuery = queryDomainRefs = 0;
+    zeroQueryNodes = 0;
+    maxDepthNodes = 0;
+    minCostNodes = 0;
+    maxRayRefs = 0;
+    addedRayRefs = removedRayRefs = 0;
+  }
+
+  
+  void
+  Print(ostream &app) const;
+
+  friend ostream &operator<<(ostream &s, const VssStatistics &stat) {
+    stat.Print(s);
+    return s;
+  }
+  
+};
+
+
+// --------------------------------------------------------------
+// For sorting rays
+// --------------------------------------------------------------
+struct  SSortableEntry
+{
+  enum EType {
+    ERayMin,
+    ERayMax
+  };
+
+  int type;
+  float value;
+  void *data;
+  
+  SSortableEntry() {}
+  SSortableEntry(const int t, const float v, void *d):type(t),
+																											value(v),
+						      data(d) {}
+	
+  friend bool operator<(const SSortableEntry &a, const SSortableEntry &b) {
+    return a.value < b.value;
+  }
+};
+
+
+class VssTreeInterior;
+
+
+// --------------------------------------------------------------
+// KD-tree node - abstract class
+// --------------------------------------------------------------
+class VssTreeNode
+{
+public:
+
+#define USE_FIXEDPOINT_T 1
+
+struct RayInfo
+{
+	// pointer to the actual ray
+	VssRay *mRay;
+	
+	// endpoints  - do we need them?
+#if USE_FIXEDPOINT_T
+	short mMinT, mMaxT;
+#else
+	float mMinT, mMaxT;
+#endif
+	
+	RayInfo():mRay(NULL) {}
+	
+	RayInfo(VssRay *r):mRay(r), mMinT(0),
+#if USE_FIXEDPOINT_T
+#define FIXEDPOINT_ONE 0x7FFE
+												//			  mMaxT(0xFFFF)
+																			mMaxT(FIXEDPOINT_ONE)
+#else
+		mMaxT(1.0f)
+#endif
+	{}
+	
+	RayInfo(VssRay *r,
+										const float _min,
+										const float _max
+										):mRay(r) {
+		SetMinT(_min);
+		SetMaxT(_max);
+	}
+	
+	RayInfo(VssRay *r,
+										const short _min,
+										const float _max
+										):mRay(r), mMinT(_min) {
+		SetMaxT(_max);
+	}
+	
+	RayInfo(VssRay *r,
+										const float _min,
+										const short _max
+										):mRay(r), mMaxT(_max) {
+		SetMinT(_min);
+	}
+	
+	friend bool operator<(const RayInfo &a, const RayInfo &b) {
+		return a.mRay < b.mRay;
+	}
+  
+	
+	float ExtrapOrigin(const int axis) const {
+		return mRay->GetOrigin(axis) + GetMinT()*mRay->GetDir(axis);
+	}
+	
+	float ExtrapTermination(const int axis) const {
+		return mRay->GetOrigin(axis) + GetMaxT()*mRay->GetDir(axis);
+	}
+	
+#if USE_FIXEDPOINT_T
+	float GetMinT () const { return mMinT/(float)(FIXEDPOINT_ONE); }
+	float GetMaxT () const { return mMaxT/(float)(FIXEDPOINT_ONE); }
+	
+	void SetMinT (const float t) {
+		mMinT = (short) (t*(float)(FIXEDPOINT_ONE));
+	}
+	
+	void SetMaxT (const float t) {
+		mMaxT = (short) (t*(float)(FIXEDPOINT_ONE));
+		mMaxT++;
+		//      if (mMaxT!=0xFFFF)
+		//	mMaxT++;
+	}
+#else
+	float GetMinT () const { return mMinT; }
+	float GetMaxT () const { return mMaxT; }
+	
+	void SetMinT (const float t) { mMinT = t; }
+	void SetMaxT (const float t) { mMaxT = t; }
+#endif 
+};
+
+
+	typedef vector<RayInfo> RayInfoContainer;
+	
+  enum { EInterior, ELeaf };
+
+  /////////////////////////////////
+  // The actual data goes here
+  
+  // link to the parent
+  VssTreeInterior *parent;
+
+  enum {SPLIT_X=0, SPLIT_Y, SPLIT_Z, SPLIT_DIRX, SPLIT_DIRY, SPLIT_DIRZ};
+  
+  // splitting axis
+  char axis;
+	
+  // depth
+  unsigned char depth;
+  
+  //  short depth;
+  //
+  /////////////////////////////////
+  
+  inline VssTreeNode(VssTreeInterior *p);
+
+  
+  virtual ~VssTreeNode() {};
+  virtual int Type() const  = 0;
+  
+
+  bool IsLeaf() const { return axis == -1; }
+  
+  virtual void Print(ostream &s) const = 0;
+
+  virtual int GetAccessTime() {
+    return 0x7FFFFFF;
+  }
+
+  
+	
+};
+
+// --------------------------------------------------------------
+// KD-tree node - interior node
+// --------------------------------------------------------------
+class VssTreeInterior :
+  public VssTreeNode
+{
+public:
+  // plane in local modelling coordinates
+  float position;
+
+  // pointers to children
+  VssTreeNode *back, *front;
+
+  // the bbox of the node
+  AxisAlignedBox3 bbox;
+
+  // the bbox of the node
+  AxisAlignedBox3 dirBBox;
+  
+  // data for caching
+  long accesses;
+  long lastAccessTime;
+  
+  VssTreeInterior(VssTreeInterior *p):VssTreeNode(p),
+					back(NULL),
+					front(NULL),
+					accesses(0),
+					lastAccessTime(-1)
+  { }
+
+  virtual int GetAccessTime() {
+    return lastAccessTime;
+  }
+
+  void SetupChildLinks(VssTreeNode *b, VssTreeNode *f) {
+    back = b;
+    front = f;
+    b->parent = f->parent = this;
+  }
+
+  void ReplaceChildLink(VssTreeNode *oldChild, VssTreeNode *newChild) {
+    if (back == oldChild)
+      back = newChild;
+    else
+      front = newChild;
+  }
+
+  virtual int Type() const  { return EInterior; }
+  
+  virtual ~VssTreeInterior() {
+    if (back)
+      delete back;
+    if (front)
+      delete front;
+  }
+  
+  virtual void Print(ostream &s) const {
+    if (axis == 0)
+      s<<"x ";
+    else
+      if (axis == 1)
+	s<<"y ";
+      else
+	s<<"z ";
+    s<<position<<" ";
+    back->Print(s);
+    front->Print(s);
+  }
+
+  
+  int ComputeRayIntersection(const RayInfo &rayData,
+														 float &t
+														 ) {
+		
+    // intersect the ray with the plane
+    float denom = rayData.mRay->GetDir(axis);
+    
+    if (fabs(denom) < 1e-20)
+      //if (denom == 0.0f)
+      return (rayData.mRay->GetOrigin(axis) > position) ? 1 : -1;
+    
+    t = (position - rayData.mRay->GetOrigin(axis))/denom;
+
+    if (t < rayData.GetMinT())
+      return (denom > 0) ? 1 : -1;
+
+    if (t > rayData.GetMaxT())
+      return (denom > 0) ? -1 : 1;
+    
+    return 0;
+  }
+
+};
+
+
+// --------------------------------------------------------------
+// KD-tree node - leaf node
+// --------------------------------------------------------------
+class VssTreeLeaf :
+  public VssTreeNode
+{
+public:
+  static int mailID;
+  int mailbox;
+  
+  RayInfoContainer rays;
+
+  
+  VssTreeLeaf(VssTreeInterior *p, const int n):VssTreeNode(p), rays() {
+    rays.reserve(n);
+  }
+  
+  virtual ~VssTreeLeaf() { }
+
+  virtual int Type() const  { return ELeaf; }
+
+  virtual void Print(ostream &s) const {
+    s<<endl<<"L: r="<<rays.size()<<endl;
+  };
+  
+  void AddRay(const RayInfo &data) {
+    rays.push_back(data);
+    data.mRay->Ref();
+  }
+
+  void Mail() { mailbox = mailID; }
+  static void NewMail() { mailID++; }
+  bool Mailed() const { return mailbox == mailID; }
+  
+  bool Mailed(const int mail) {
+    return mailbox >= mailID + mail;
+  }
+
+};
+
+// Inline functions
+inline
+VssTreeNode::VssTreeNode(VssTreeInterior *p):
+  parent(p), axis(-1), depth(p ? p->depth + 1 : 0) {}
+
+
+
+// ---------------------------------------------------------------
+// Main LSDS search class
+// ---------------------------------------------------------------
+class VssTree
+{
+  struct TraversalData
+  {  
+    VssTreeNode *node;
+    AxisAlignedBox3 bbox;
+    int depth;
+    float priority;
+    
+    TraversalData() {}
+
+    TraversalData(VssTreeNode *n, const float p):
+      node(n), priority(p)
+    {}
+
+    TraversalData(VssTreeNode *n,
+		   const AxisAlignedBox3 &b,
+		   const int d):
+      node(n), bbox(b), depth(d) {}
+    
+		
+    // comparator for the 
+    struct less_priority : public
+    binary_function<const TraversalData, const TraversalData, bool> {
+			
+      bool operator()(const TraversalData a, const TraversalData b) {
+				return a.priority < b.priority;
+      }
+      
+    };
+
+    //    ~TraversalData() {}
+    //    TraversalData(const TraversalData &s):node(s.node), bbox(s.bbox), depth(s.depth) {}
+    
+    friend bool operator<(const TraversalData &a,
+													const TraversalData &b) {
+      //      return a.node->queries.size() < b.node->queries.size();
+      VssTreeLeaf *leafa = (VssTreeLeaf *) a.node;
+      VssTreeLeaf *leafb = (VssTreeLeaf *) b.node;
+      //      return leafa->rays.size()*a.bbox.GetVolume() < leafb->rays.size()*b.bbox.GetVolume();
+      return 
+				leafa->rays.size()*a.bbox.GetVolume()
+				<
+				leafb->rays.size()*b.bbox.GetVolume();
+    }
+  };
+  
+  // simplified data for ray traversal only...
+
+  struct RayTraversalData {
+    
+    VssTreeNode::RayInfo rayData;
+    VssTreeNode *node;
+    
+    RayTraversalData() {}
+    RayTraversalData(VssTreeNode *n,
+										 const VssTreeNode::RayInfo &data):
+      rayData(data), node(n) {}
+  };
+	
+public:
+  /////////////////////////////
+  // The core pointer
+  VssTreeNode *root;
+  
+  /////////////////////////////
+  // Basic properties
+
+  // total number of nodes of the tree
+  int nodes;
+  // axis aligned bounding box of the scene
+  AxisAlignedBox3 bbox;
+
+  // axis aligned bounding box of directions
+  AxisAlignedBox3 dirBBox;
+  
+  /////////////////////////////
+  // Construction parameters
+
+  // epsilon used for the construction
+  float epsilon;
+
+  // ratio between traversal and intersection costs
+  float ct_div_ci;
+  // max depth of the tree
+  int termMaxDepth;
+  // minimal cost of the node to still get subdivided
+  int termMinCost;
+  // minimal ratio of the volume of the cell and the query volume
+  float termMinSize;
+
+  // minimal cost imporvement to subdivide a node
+  float maxCostRatio;
+  
+  // randomized construction
+  bool randomize;
+
+  // type of the splitting to use fo rthe tree construction
+  enum {ESplitRegular, ESplitVolume, ESplitQueries };
+  int splitType;
+	
+  // maximal size of the box on which the refdir splitting can be performed
+  // (relative to the scene bbox
+  float refDirBoxMaxSize;
+  
+  // maximum alovable memory in MB
+  float maxTotalMemory;
+
+  // maximum alovable memory for static kd tree in MB
+  float maxStaticMemory;
+
+  // this is used during the construction depending
+  // on the type of the tree and queries...
+  float maxMemory;
+
+
+  // minimal acess time for collapse
+  int accessTimeThreshold;
+
+	// minimal depth at which to perform collapse
+  int minCollapseDepth;
+
+  
+  // reusable array of split candidates
+  vector<SSortableEntry> *splitCandidates;
+  /////////////////////////////
+
+	VssStatistics stat;
+	
+  
+  VssTree();
+  virtual ~VssTree();
+
+  virtual void
+  Construct(
+						VssRayContainer &rays,
+						const bool onlyStaticPart
+						);
+	
+  // incemental construction
+  virtual void UpdateRays(VssRayContainer &remove,
+													VssRayContainer &add
+													);
+	
+  
+	
+  VssTreeNode *
+  Locate(const Vector3 &v);
+	
+  VssTreeNode *
+  SubdivideNode(VssTreeLeaf *leaf,
+								const AxisAlignedBox3 &box,
+								AxisAlignedBox3 &backBox,
+								AxisAlignedBox3 &frontBox
+								);
+	
+  VssTreeNode *
+  Subdivide(const TraversalData &tdata);
+	
+  int
+  SelectPlane(VssTreeLeaf *leaf,
+	      const AxisAlignedBox3 &box,
+	      float &position
+	      );
+
+  void
+  SortSplitCandidates(
+											VssTreeLeaf *node,
+											const int axis
+											);
+	
+  
+  // return memory usage in MB
+  float GetMemUsage() const {
+    return 
+      (sizeof(VssTree) +
+       stat.Leaves()*sizeof(VssTreeLeaf) +
+       stat.Interior()*sizeof(VssTreeInterior) +
+       stat.rayRefs*sizeof(VssTreeNode::RayInfo))/(1024.0f*1024.0f);
+  }
+	
+  float GetRayMemUsage() const {
+    return 
+      stat.rays*(sizeof(VssRay))/(1024.0f*1024.0f);
+  }
+  
+
+  AxisAlignedBox3 GetBBox(const VssTreeNode *node) {
+    if (node->parent == NULL)
+      return bbox;
+
+    if (!node->IsLeaf())
+      return ((VssTreeInterior *)node)->bbox;
+
+    if (node->parent->axis >= 3)
+      return node->parent->bbox;
+      
+    AxisAlignedBox3 box(node->parent->bbox);
+    if (node->parent->front == node)
+      box.SetMin(node->parent->axis, node->parent->position);
+    else
+      box.SetMax(node->parent->axis, node->parent->position);
+    return box;
+  }
+
+  AxisAlignedBox3 GetDirBBox(const VssTreeNode *node) {
+
+    if (node->parent == NULL)
+      return dirBBox;
+    
+    if (!node->IsLeaf() )
+      return ((VssTreeInterior *)node)->dirBBox;
+
+    if (node->parent->axis < 3)
+      return node->parent->dirBBox;
+    
+    AxisAlignedBox3 dBBox(node->parent->dirBBox);
+
+    if (node->parent->front == node)
+      dBBox.SetMin(node->parent->axis - 3, node->parent->position);
+    else
+      dBBox.SetMax(node->parent->axis - 3, node->parent->position);
+    return dBBox;
+  }
+  
+  int
+  ReleaseMemory(const int time);
+
+  int
+  CollapseSubtree(VssTreeNode *node, const int time);
+
+  void
+  CountAccess(VssTreeInterior *node, const long time) {
+    node->accesses++;
+    node->lastAccessTime = time;
+  }
+
+  VssTreeNode *
+  SubdivideLeaf(
+								VssTreeLeaf *leaf,
+								const float SAThreshold
+								);
+
+  void
+  RemoveRay(VssRay *ray,
+						vector<VssTreeLeaf *> *affectedLeaves,
+						const bool removeAllScheduledRays
+						);
+
+  void
+  AddRay(VssRay *ray);
+	
+  void
+  TraverseInternalNode(
+											 RayTraversalData &data,
+											 stack<RayTraversalData> &tstack);
+
+	void
+  EvaluateLeafStats(const TraversalData &data);
+
+};
+
+
+#endif // __LSDS_KDTREE_H__
+
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Mesh.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Mesh.h	(revision 354)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Mesh.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/gzstream.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/gzstream.h	(revision 162)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/gzstream.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Environment.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Environment.cpp	(revision 367)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Environment.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/KdTree.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/KdTree.cpp	(revision 362)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/KdTree.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.cpp	(revision 370)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCell.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCell.h	(revision 369)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/ViewCell.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Polytope.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Polytope.h	(revision 260)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Polytope.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/ExactPreprocessor.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/ExactPreprocessor.h	(revision 162)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/ExactPreprocessor.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Vector3.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Vector3.cpp	(revision 354)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Vector3.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Preprocessor.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Preprocessor.cpp	(revision 366)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Preprocessor.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Plane3.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Plane3.h	(revision 362)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Plane3.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/VssTree.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/VssTree.cpp	(revision 372)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/VssTree.cpp	(revision 372)
@@ -0,0 +1,925 @@
+
+// ================================================================
+// $Id: lsds_kdtree.cpp,v 1.18 2005/04/16 09:34:21 bittner Exp $
+// ****************************************************************
+/**
+   The KD tree based LSDS
+*/
+// Initial coding by
+/**
+   @author Jiri Bittner
+*/
+
+// Standard headers
+#include <stack>
+#include <queue>
+#include <algorithm>
+#include <fstream>
+#include <string>
+
+#include "VssTree.h"
+
+#include "Environment.h"
+#include "VssRay.h"
+
+
+#define DEBUG_DIR_SPLIT 0
+
+
+
+// Static variables
+int
+VssTreeLeaf::mailID = 0;
+
+
+// Constructor
+VssTree::VssTree()
+{
+  environment->GetIntValue("VssTree.maxDepth", termMaxDepth);
+
+	environment->GetIntValue("VssTree.minCost", termMinCost);
+
+  environment->GetFloatValue("VssTree.minSize", termMinSize);
+  termMinSize = sqr(termMinSize);
+
+  environment->GetFloatValue("VssTree.refDirBoxMaxSize", refDirBoxMaxSize);
+  refDirBoxMaxSize = sqr(refDirBoxMaxSize);
+  
+  environment->GetFloatValue("VssTree.epsilon", epsilon);
+  environment->GetFloatValue("VssTree.ct_div_ci", ct_div_ci);
+
+  environment->GetFloatValue("VssTree.maxTotalMemory", maxTotalMemory);
+  environment->GetFloatValue("VssTree.maxStaticMemory", maxStaticMemory);
+  
+  environment->GetFloatValue("VssTree.maxCostRatio", maxCostRatio);
+  
+
+
+  float refDirAngle;
+  environment->GetFloatValue("VssTree.refDirAngle", refDirAngle);
+
+  environment->GetIntValue("VssTree.accessTimeThreshold", accessTimeThreshold);
+  //= 1000;
+  environment->GetIntValue("VssTree.minCollapseDepth", minCollapseDepth);
+  //  int minCollapseDepth = 4;
+
+  //  pRefDirThresh = cos(0.5*M_PI - M_PI*refDirAngle/180.0);
+	//  cosRefDir = cos(M_PI*refDirAngle/180.0);
+	//  sinRefDir = sin(M_PI*refDirAngle/180.0);
+  
+  
+  // split type
+	char sname[128];
+	environment->GetStringValue("VssTree.splitType", sname);
+  string name(sname);
+	
+  if (name.compare("regular") == 0)
+    splitType = ESplitRegular;
+  else
+    if (name.compare("volume") == 0)
+      splitType = ESplitVolume;
+    else
+      if (name.compare("queries") == 0)
+				splitType = ESplitQueries;
+      else {
+				cerr<<"Invalid VssTree split type "<<name<<endl;
+				exit(1);
+      }
+	
+  environment->GetBoolValue("VssTree.randomize", randomize);
+
+  root = NULL;
+  
+  splitCandidates = new vector<SSortableEntry>;
+}
+
+
+VssTree::~VssTree()
+{
+  if (root)
+    delete root;
+}
+
+
+
+
+void
+VssStatistics::Print(ostream &app) const
+{
+  app << "===== VssTree statistics ===============\n";
+
+  app << "#N_RAYS Number of rays )\n"
+      << rays <<endl;
+  app << "#N_DOMAINS  ( Number of query domains )\n"
+      << queryDomains <<endl;
+  
+  app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
+
+  app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
+
+  app << "#N_SPLITS ( Number of splits in axes x y z dx dy dz \n";
+  for (int i=0; i<7; i++)
+    app << splits[i] <<" ";
+  app <<endl;
+
+  app << "#N_RAYREFS ( Number of rayRefs )\n" <<
+    rayRefs << "\n";
+
+  app << "#N_RAYRAYREFS  ( Number of rayRefs / ray )\n" <<
+    rayRefs/(double)rays << "\n";
+
+  app << "#N_LEAFRAYREFS  ( Number of rayRefs / leaf )\n" <<
+    rayRefs/(double)Leaves() << "\n";
+
+  app << "#N_MAXRAYREFS  ( Max number of rayRefs / leaf )\n" <<
+    maxRayRefs << "\n";
+
+  app << "#N_NONEMPTYRAYREFS  ( Number of rayRefs in nonEmpty leaves / non empty leaf )\n" <<
+    rayRefsNonZeroQuery/(double)(Leaves() - zeroQueryNodes) << "\n";
+
+  app << "#N_LEAFDOMAINREFS  ( Number of query domain Refs / leaf )\n" <<
+    queryDomainRefs/(double)Leaves() << "\n";
+
+  //  app << setprecision(4);
+
+  app << "#N_PEMPTYLEAVES  ( Percentage of leaves with zero query domains )\n"<<
+    zeroQueryNodes*100/(double)Leaves()<<endl;
+
+  app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maxdepth )\n"<<
+    maxDepthNodes*100/(double)Leaves()<<endl;
+
+  app << "#N_PMINCOSTLEAVES  ( Percentage of leaves with minCost )\n"<<
+    minCostNodes*100/(double)Leaves()<<endl;
+
+  app << "#N_ADDED_RAYREFS  (Number of dynamically added ray references )\n"<<
+    addedRayRefs<<endl;
+
+  app << "#N_REMOVED_RAYREFS  (Number of dynamically removed ray references )\n"<<
+    removedRayRefs<<endl;
+
+  //  app << setprecision(4);
+
+  app << "#N_CTIME  ( Construction time [s] )\n"
+      << Time() << " \n";
+
+  app << "===== END OF VssTree statistics ==========\n";
+
+}
+
+
+
+
+void
+VssTree::Construct(
+									 VssRayContainer &rays,
+									 const bool onlyStaticPart
+									 )
+{
+  stat.Start();
+  
+  if (onlyStaticPart)
+    maxMemory = maxStaticMemory;
+  else
+    maxMemory = maxTotalMemory;
+
+
+  if (root)
+    delete root;
+
+  root = new VssTreeLeaf(NULL, rays.size());
+  // first construct a leaf that will get subdivide
+  VssTreeLeaf *leaf = (VssTreeLeaf *) root;
+
+  stat.nodes = 1;
+  
+  bbox.Initialize();
+
+  dirBBox.Initialize();
+  
+  for(VssRayContainer::const_iterator ri = rays.begin();
+      ri != rays.end();
+      ri++) {
+    leaf->AddRay(VssTreeNode::RayInfo(*ri));
+    bbox.Include((*ri)->GetOrigin());
+    bbox.Include((*ri)->GetTermination());
+    
+    dirBBox.Include( (*ri)->GetNormalizedDir() );
+  }
+  
+  
+  stat.rays = leaf->rays.size();
+  
+  // Subdivide();
+  root = Subdivide(TraversalData(leaf, bbox, 0));
+
+  if (splitCandidates) {
+    // force realease of this vector
+    delete splitCandidates;
+    splitCandidates = new vector<SSortableEntry>;
+  }
+  
+  stat.Stop();
+
+  cout<<"##################################"<<endl;
+  
+}
+
+
+VssTreeNode *
+VssTree::Subdivide(const TraversalData &tdata)
+{
+  VssTreeNode *result = NULL;
+
+  priority_queue<TraversalData> tStack;
+  //  stack<TraversalData> tStack;
+  
+  tStack.push(tdata);
+
+  AxisAlignedBox3 backBox;
+  AxisAlignedBox3 frontBox;
+
+  
+  while (!tStack.empty()) {
+
+    if ( GetMemUsage() > maxMemory ) {
+      // count statistics on unprocessed leafs
+      while (!tStack.empty()) {
+				EvaluateLeafStats(tStack.top());
+				tStack.pop();
+      }
+      break;
+    }
+    
+    TraversalData data = tStack.top();
+    tStack.pop();
+    
+    VssTreeNode *node = SubdivideNode((VssTreeLeaf *) data.node,
+																			data.bbox,
+																			backBox,
+																			frontBox
+																			);
+    if (result == NULL)
+      result = node;
+    
+    if (!node->IsLeaf()) {
+			
+      VssTreeInterior *interior = (VssTreeInterior *) node;
+      // push the children on the stack
+      tStack.push(TraversalData(interior->back, backBox, data.depth+1));
+      tStack.push(TraversalData(interior->front, frontBox, data.depth+1));
+			
+    } else {
+      EvaluateLeafStats(data);
+    }
+  }
+
+  return result;
+}
+
+
+// returns selected plane for subdivision
+int
+VssTree::SelectPlane(
+										 VssTreeLeaf *leaf,
+										 const AxisAlignedBox3 &box,
+										 float &position
+										 )
+{
+
+  int minDirDepth = 6;
+
+  Vector3 size = box.Size();
+  
+  float pDirSubdivision = 0.0f;
+  if (leaf->depth >= minDirDepth && RandomValue(0, 1) < pDirSubdivision) {
+    AxisAlignedBox3 dirbox = GetDirBBox(leaf);
+    Vector3 size = dirbox.Size();
+    int axis = size.DrivingAxis();
+    position = (dirbox.Min()[axis] + dirbox.Max()[axis])*0.5f;
+    return axis + 3;
+  }
+  
+  int axis = size.DrivingAxis();
+
+  
+  if (splitType == ESplitRegular) {
+		position = (box.Min()[axis] + box.Max()[axis])*0.5f;
+  }
+  
+  return axis;
+}
+
+
+void
+VssTree::EvaluateLeafStats(const TraversalData &data)
+{
+
+  // the node became a leaf -> evaluate stats for leafs
+  VssTreeLeaf *leaf = (VssTreeLeaf *)data.node;
+
+  if (data.depth > termMaxDepth)
+    stat.maxDepthNodes++;
+  
+  if ( (int)(leaf->rays.size()) < termMinCost)
+    stat.minCostNodes++;
+  
+
+  if ( (int)(leaf->rays.size()) > stat.maxRayRefs)
+    stat.maxRayRefs = leaf->rays.size();
+  
+}
+
+
+
+VssTreeNode *
+VssTree::SubdivideNode(
+											 VssTreeLeaf *leaf,
+											 const AxisAlignedBox3 &box,
+											 AxisAlignedBox3 &backBBox,
+											 AxisAlignedBox3 &frontBBox
+											 )
+{
+  
+  if ( ((int)leaf->rays.size() < termMinCost) ||
+       (leaf->depth >= termMaxDepth) ) {
+    return leaf;
+  }
+  
+  float position;
+  
+  // select subdivision axis
+  int axis = SelectPlane( leaf, box, position);
+
+  if (axis == -1) {
+    return leaf;
+  }
+  
+  stat.nodes+=2;
+  stat.splits[axis]++;
+
+  // add the new nodes to the tree
+  VssTreeInterior *node = new VssTreeInterior(leaf->parent);
+
+  node->axis = axis;
+  node->position = position;
+  node->bbox = box;
+
+  node->dirBBox = GetDirBBox(leaf);
+  
+  backBBox = box;
+  frontBBox = box;
+
+  // first count ray sides
+  int raysBack = 0;
+  int raysFront = 0;
+  
+  if (axis <= VssTreeNode::SPLIT_Z) {
+    backBBox.SetMax(axis, position);
+    frontBBox.SetMin(axis, position);
+    
+    // this is the main ray classification loop!
+    for(VssTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
+				ri != leaf->rays.end();
+				ri++)
+      if ((*ri).mRay->IsActive()) {
+				
+				// determine the side of this ray with respect to the plane
+				int side = node->ComputeRayIntersection(*ri, (*ri).mRay->mT);
+				
+				(*ri).mRay->mSide = side;
+				
+				if (side <= 0)
+					raysBack++;
+				
+				if (side >= 0)
+					raysFront++;
+				
+#if 0
+	cout<<"-------------------"<<endl;
+	cout<<"axis = "<<(int)node->axis<<"\t position="<<node->position<<endl;
+	cout<<"origin="<<(*ri).ray->GetOrigin()<<"\t  term="<<(*ri).ray->GetTermination()<<endl;
+	cout<<"side = "<<side<<"\t t="<<t<<endl;
+	cout<<"-------------------"<<endl;
+#endif
+      }
+  } else {
+    // directional split
+    for(VssTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
+				ri != leaf->rays.end();
+				ri++)
+      if ((*ri).mRay->IsActive()) {
+				
+				// determine the side of this ray with respect to the plane
+				int side;
+				if ((*ri).mRay->GetNormalizedDir(axis - 3) > position)
+					side = 1;
+				else
+					side = -1;
+				
+				
+				(*ri).mRay->mSide = side;
+				
+				if (side <= 0)
+					raysBack++;
+				
+				if (side >= 0)
+					raysFront++;
+      }
+  }
+  
+  VssTreeLeaf *back = new VssTreeLeaf(node, raysBack);
+  VssTreeLeaf *front = new VssTreeLeaf(node, raysFront);
+
+  // replace a link from node's parent
+  if (  leaf->parent )
+    leaf->parent->ReplaceChildLink(leaf, node);
+  // and setup child links
+  node->SetupChildLinks(back, front);
+
+  
+  if (axis <= VssTreeNode::SPLIT_Z) {
+    for(VssTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
+				ri != leaf->rays.end();
+				ri++) {
+      if ((*ri).mRay->IsActive()) {
+				
+				// first unref ray from the former leaf
+				(*ri).mRay->Unref();
+				
+				if ((*ri).mRay->mSide == 0) {
+					if ((*ri).mRay->HasPosDir(axis)) {
+						back->AddRay(VssTreeNode::RayInfo((*ri).mRay,
+																							(*ri).mMinT,
+																							(*ri).mRay->mT)
+												 );
+						front->AddRay(VssTreeNode::RayInfo((*ri).mRay,
+																							 (*ri).mRay->mT,
+																							 (*ri).mMaxT));
+					} else {
+						back->AddRay(VssTreeNode::RayInfo((*ri).mRay,
+																							(*ri).mRay->mT,
+																							(*ri).mMaxT));
+						front->AddRay(VssTreeNode::RayInfo((*ri).mRay,
+																							 (*ri).mMinT,
+																							 (*ri).mRay->mT));
+					}
+				} else
+					if ((*ri).mRay->mSide == 1) 
+						front->AddRay(*ri);
+					else
+						back->AddRay(*ri);
+      } else
+				(*ri).mRay->Unref();
+    }
+  } else {
+    // rays front/back
+    
+#if DEBUG_DIR_SPLIT
+    cout<<"dir split, depth="<<(int)leaf->depth<<" front= "<<raysFront<<" back="<<raysBack<<endl;
+#endif
+    
+    for(VssTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
+				ri != leaf->rays.end();
+				ri++) {
+      if ((*ri).mRay->IsActive()) {
+				// first unref ray from the former leaf
+				(*ri).mRay->Unref();
+				
+				if ((*ri).mRay->mSide == 1)
+					front->AddRay(*ri);
+				else
+					back->AddRay(*ri);
+				
+      } else
+				(*ri).mRay->Unref();
+    }
+  }
+	
+  // update stats
+  stat.rayRefs -= leaf->rays.size();
+  stat.rayRefs += raysBack + raysFront;
+
+  
+  delete leaf;
+  return node;
+}
+
+
+
+
+
+
+int
+VssTree::ReleaseMemory(const int time)
+{
+  stack<VssTreeNode *> tstack;
+  
+  // find a node in the tree which subtree will be collapsed
+  int maxAccessTime = time - accessTimeThreshold;
+  int released;
+
+  tstack.push(root);
+
+  while (!tstack.empty()) {
+    VssTreeNode *node = tstack.top();
+    tstack.pop();
+    
+ 
+    if (!node->IsLeaf()) {
+      VssTreeInterior *in = (VssTreeInterior *)node;
+      //      cout<<"depth="<<(int)in->depth<<" time="<<in->lastAccessTime<<endl;
+      if (in->depth >= minCollapseDepth &&
+	  in->lastAccessTime <= maxAccessTime) {
+	released = CollapseSubtree(node, time);
+	break;
+      }
+      
+      if (in->back->GetAccessTime() <  
+	  in->front->GetAccessTime()) {
+	tstack.push(in->front);
+	tstack.push(in->back);
+      } else {
+	tstack.push(in->back);
+	tstack.push(in->front);
+      }
+    }
+  }
+
+  while (tstack.empty()) {
+    // could find node to collaps...
+    //    cout<<"Could not find a node to release "<<endl;
+    break;
+  }
+  
+  return released;
+}
+
+
+
+
+VssTreeNode *
+VssTree::SubdivideLeaf(
+											 VssTreeLeaf *leaf,
+											 const float sizeThreshold
+											 )
+{
+  VssTreeNode *node = leaf;
+
+  AxisAlignedBox3 leafBBox = GetBBox(leaf);
+
+	static int pass = 0;
+	pass ++;
+	
+  // check if we should perform a dynamic subdivision of the leaf
+  if (
+      leaf->rays.size() > (unsigned)termMinCost &&
+      SqrMagnitude(leafBBox.Size()) > sizeThreshold) {
+    
+		// memory check and realese...
+    if (GetMemUsage() > maxTotalMemory) {
+      ReleaseMemory( pass );
+    }
+    
+    AxisAlignedBox3 backBBox, frontBBox;
+
+    // subdivide the node
+    node =
+      SubdivideNode(leaf,
+										leafBBox,
+										backBBox,
+										frontBBox
+										);
+  }
+  return node;
+}
+
+
+
+void
+VssTree::UpdateRays(VssRayContainer &remove,
+										VssRayContainer &add
+										)
+{
+  VssTreeLeaf::NewMail();
+
+  // schedule rays for removal
+  for(VssRayContainer::const_iterator ri = remove.begin();
+      ri != remove.end();
+      ri++) {
+    (*ri)->ScheduleForRemoval();
+  }
+
+  int inactive=0;
+
+  for(VssRayContainer::const_iterator ri = remove.begin();
+      ri != remove.end();
+      ri++) {
+    if ((*ri)->ScheduledForRemoval())
+//      RemoveRay(*ri, NULL, false);
+// !!! BUG - with true it does not work correctly - aggreated delete
+      RemoveRay(*ri, NULL, true);
+    else
+      inactive++;
+  }
+
+
+  //  cout<<"all/inactive"<<remove.size()<<"/"<<inactive<<endl;
+  
+  for(VssRayContainer::const_iterator ri = add.begin();
+      ri != add.end();
+      ri++) {
+    AddRay(*ri);
+  }
+  
+
+}
+
+
+void
+VssTree::RemoveRay(VssRay *ray,
+									 vector<VssTreeLeaf *> *affectedLeaves,
+									 const bool removeAllScheduledRays
+									 )
+{
+	
+  stack<RayTraversalData> tstack;
+
+  tstack.push(RayTraversalData(root, VssTreeNode::RayInfo(ray)));
+  
+  RayTraversalData data;
+
+  // cout<<"Number of ray refs = "<<ray->RefCount()<<endl;
+
+  while (!tstack.empty()) {
+    data = tstack.top();
+    tstack.pop();
+
+    if (!data.node->IsLeaf()) {
+      // split the set of rays in two groups intersecting the
+      // two subtrees
+
+      TraverseInternalNode(data, tstack);
+      
+    } else {
+      // remove the ray from the leaf
+      // find the ray in the leaf and swap it with the last ray...
+      VssTreeLeaf *leaf = (VssTreeLeaf *)data.node;
+      
+      if (!leaf->Mailed()) {
+	leaf->Mail();
+	if (affectedLeaves)
+	  affectedLeaves->push_back(leaf);
+	
+	if (removeAllScheduledRays) {
+	  int tail = leaf->rays.size()-1;
+
+	  for (int i=0; i < (int)(leaf->rays.size()); i++) {
+	    if (leaf->rays[i].mRay->ScheduledForRemoval()) {
+	      // find a ray to replace it with
+	      while (tail >= i && leaf->rays[tail].mRay->ScheduledForRemoval()) {
+		stat.removedRayRefs++;
+		leaf->rays[tail].mRay->Unref();
+		leaf->rays.pop_back();
+		tail--;
+	      }
+
+	      if (tail < i)
+		break;
+	      
+	      stat.removedRayRefs++;
+	      leaf->rays[i].mRay->Unref();
+	      leaf->rays[i] = leaf->rays[tail];
+	      leaf->rays.pop_back();
+	      tail--;
+	    }
+	  }
+	}
+      }
+      
+      if (!removeAllScheduledRays)
+	for (int i=0; i < (int)leaf->rays.size(); i++) {
+	  if (leaf->rays[i].mRay == ray) {
+	    stat.removedRayRefs++;
+	    ray->Unref();
+	    leaf->rays[i] = leaf->rays[leaf->rays.size()-1];
+	    leaf->rays.pop_back();
+	    // check this ray again
+	    break;
+	  }
+	}
+      
+    }
+  }
+
+  if (ray->RefCount() != 0) {
+    cerr<<"Error: Number of remaining refs = "<<ray->RefCount()<<endl;
+    exit(1);
+  }
+  
+}
+
+
+void
+VssTree::AddRay(VssRay *ray)
+{
+
+  stack<RayTraversalData> tstack;
+  
+  tstack.push(RayTraversalData(root, VssTreeNode::RayInfo(ray)));
+  
+  RayTraversalData data;
+  
+  while (!tstack.empty()) {
+    data = tstack.top();
+    tstack.pop();
+
+    if (!data.node->IsLeaf()) {
+
+      TraverseInternalNode(data, tstack);
+
+    } else {
+      // remove the ray from the leaf
+      // find the ray in the leaf and swap it with the last ray...
+      VssTreeLeaf *leaf = (VssTreeLeaf *)data.node;
+      leaf->AddRay(data.rayData);
+      stat.addedRayRefs++;
+    }
+  }
+}
+
+void
+VssTree::TraverseInternalNode(
+															RayTraversalData &data,
+															stack<RayTraversalData> &tstack)
+{
+  VssTreeInterior *in =  (VssTreeInterior *) data.node;
+  
+  if (in->axis <= VssTreeNode::SPLIT_Z) {
+
+    // determine the side of this ray with respect to the plane
+    int side = in->ComputeRayIntersection(data.rayData,
+																					data.rayData.mRay->mT);
+    
+  
+    if (side == 0) {
+      if (data.rayData.mRay->HasPosDir(in->axis)) {
+				tstack.push(RayTraversalData(in->back,
+																		 VssTreeNode::RayInfo(data.rayData.mRay,
+																													data.rayData.mMinT,
+																													data.rayData.mRay->mT))
+										);
+				
+				tstack.push(RayTraversalData(in->front,
+																		 VssTreeNode::RayInfo(data.rayData.mRay,
+																													data.rayData.mRay->mT,
+																													data.rayData.mMaxT
+																													))
+										);
+	
+      } else {
+				tstack.push(RayTraversalData(in->back,
+																		 VssTreeNode::RayInfo(data.rayData.mRay,
+																													data.rayData.mRay->mT,
+																													data.rayData.mMaxT
+																													))
+										);
+				
+				tstack.push(RayTraversalData(in->front,
+																		 VssTreeNode::RayInfo(data.rayData.mRay,
+																													data.rayData.mMinT,
+																													data.rayData.mRay->mT))
+										);
+				
+				
+      }
+    } else
+      if (side == 1)
+				tstack.push(RayTraversalData(in->front, data.rayData));
+      else
+				tstack.push(RayTraversalData(in->back, data.rayData));
+  } 
+  else {
+    // directional split
+//      if (DotProd(data.rayData.mRay->GetNormalizedDir(), in->refDir) > sinRefDir )
+//        tstack.push(RayTraversalData(in->front, data.rayData));
+//      else
+//        tstack.push(RayTraversalData(in->back, data.rayData));
+  }
+  
+}
+
+
+int
+VssTree::CollapseSubtree(VssTreeNode *sroot, const int time)
+{
+  // first count all rays in the subtree
+  // use mail 1 for this purpose
+  stack<VssTreeNode *> tstack;
+  int rayCount = 0;
+  int totalRayCount = 0;
+  int collapsedNodes = 0;
+
+#if DEBUG_COLLAPSE
+  cout<<"Collapsing subtree"<<endl;
+  cout<<"acessTime="<<sroot->GetAccessTime()<<endl;
+  cout<<"depth="<<(int)sroot->depth<<endl;
+#endif
+
+//    tstat.collapsedSubtrees++;
+//    tstat.collapseDepths += (int)sroot->depth;
+//    tstat.collapseAccessTimes += time - sroot->GetAccessTime();
+  
+  tstack.push(sroot);
+  VssRay::NewMail();
+	
+  while (!tstack.empty()) {
+    collapsedNodes++;
+    VssTreeNode *node = tstack.top();
+    tstack.pop();
+
+    if (node->IsLeaf()) {
+      VssTreeLeaf *leaf = (VssTreeLeaf *) node;
+      for(VssTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
+					ri != leaf->rays.end();
+					ri++) {
+				
+				totalRayCount++;
+				if ((*ri).mRay->IsActive() && !(*ri).mRay->Mailed()) {
+					(*ri).mRay->Mail();
+					rayCount++;
+				}
+      }
+    } else {
+      tstack.push(((VssTreeInterior *)node)->back);
+      tstack.push(((VssTreeInterior *)node)->front);
+    }
+  }
+  
+  VssRay::NewMail();
+
+  // create a new node that will hold the rays
+  VssTreeLeaf *newLeaf = new VssTreeLeaf( sroot->parent, rayCount );
+  if (  newLeaf->parent )
+    newLeaf->parent->ReplaceChildLink(sroot, newLeaf);
+  
+
+  tstack.push( sroot );
+
+  while (!tstack.empty()) {
+
+    VssTreeNode *node = tstack.top();
+    tstack.pop();
+
+    if (node->IsLeaf()) {
+      VssTreeLeaf *leaf = (VssTreeLeaf *) node;
+      
+      for(VssTreeNode::RayInfoContainer::iterator ri = leaf->rays.begin();
+					ri != leaf->rays.end();
+					ri++) {
+				
+				// unref this ray from the old node
+				
+				if ((*ri).mRay->IsActive()) {
+					(*ri).mRay->Unref();
+					if (!(*ri).mRay->Mailed()) {
+						(*ri).mRay->Mail();
+						newLeaf->AddRay(*ri);
+					}
+				} else
+					(*ri).mRay->Unref();
+				
+      }
+    } else {
+      tstack.push(((VssTreeInterior *)node)->back);
+      tstack.push(((VssTreeInterior *)node)->front);
+    }
+  }
+  
+  // delete the node and all its children
+  delete sroot;
+  
+//   for(VssTreeNode::SRayContainer::iterator ri = newLeaf->rays.begin();
+//       ri != newLeaf->rays.end();
+//       ri++)
+//     (*ri).ray->UnMail(2);
+
+
+#if DEBUG_COLLAPSE
+  cout<<"Total memory before="<<GetMemUsage()<<endl;
+#endif
+
+  stat.nodes -= collapsedNodes - 1;
+  stat.rayRefs -= totalRayCount - rayCount;
+  
+#if DEBUG_COLLAPSE
+  cout<<"collapsed nodes"<<collapsedNodes<<endl;
+  cout<<"collapsed rays"<<totalRayCount - rayCount<<endl;
+  cout<<"Total memory after="<<GetMemUsage()<<endl;
+  cout<<"================================"<<endl;
+#endif
+
+	//  tstat.collapsedNodes += collapsedNodes;
+	//  tstat.collapsedRays += totalRayCount - rayCount;
+    
+  return totalRayCount - rayCount;
+}
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/VssRay.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/VssRay.h	(revision 372)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/VssRay.h	(revision 372)
@@ -0,0 +1,146 @@
+#ifndef __SIMPLE_RAY_H
+#define __SIMPLE_RAY_H
+
+#include <vector>
+using namespace std;
+#include "Vector3.h"
+
+class AxisAlignedBox3;
+class Intersectable;
+
+class VssRay {
+public:
+	// various flags
+  enum {
+    FPosDirX = 1,  // the direction of ray in X-axis is positive
+    FPosDirY = 2,  // the direction of ray in Y-axis is positive
+    FPosDirZ = 4   // the direction of ray in Z-axis is positive
+  };
+
+	static int mailID;
+  int mMailbox;
+	
+	// side of the ray - used for the ray classification
+  char mSide;
+	
+	// computed t
+  float mT;
+
+	// inverse of the ray size
+	float mInvSize;
+	
+	// counter of references to this ray
+  short mRefCount;
+
+	// various flags
+  char mFlags;
+	
+  Vector3 mOrigin;
+  Vector3 mTermination;
+	
+	Intersectable *mOriginObject;
+	Intersectable *mTerminationObject;
+	
+	//////////////////////////////
+	VssRay(const Vector3 &origin,
+								const Vector3 &termination,
+								Intersectable *originObject,
+								Intersectable *terminationObject):
+		mMailbox(-1),
+    mOrigin(origin),
+    mTermination(termination),
+		mOriginObject(originObject),
+		mTerminationObject(terminationObject),
+		mRefCount(0),
+    mFlags(0)
+	{
+    Precompute();
+  }
+	
+  void Precompute() {
+    mFlags = 0;
+    if (GetDir().x>0.0f) mFlags |= FPosDirX;
+    if (GetDir().y>0.0f) mFlags |= FPosDirY;
+    if (GetDir().z>0.0f) mFlags |= FPosDirZ;
+    mInvSize = 1.0f/Magnitude(GetDir());
+  }
+
+  void Mail() { mMailbox = mailID; }
+  static void NewMail() { mailID++; }
+  bool Mailed() const { return mMailbox == mailID; }
+
+  bool Mailed(const int mail) {
+    return mMailbox >= mailID + mail;
+  }
+
+	
+	
+	Vector3 GetOrigin() const { return mOrigin; }
+  Vector3 GetTermination() const { return mTermination; }
+  Vector3 GetDir() const { return mTermination - mOrigin; }
+  //  Vector3 GetNormalizedDir() const { return Normalize(termination - mOrigin); }
+  Vector3 GetNormalizedDir() const { return (mTermination - mOrigin)*mInvSize; }
+
+	float GetInvSize() const { return  mInvSize; }
+  float GetOrigin(const int axis) const { return mOrigin[axis]; }
+  float GetTermination(const int axis) const { return mTermination[axis]; }
+  float GetDir(const int axis) const { return mTermination[axis] - mOrigin[axis]; }
+  float GetNormalizedDir(const int axis) const {
+    return (mTermination[axis] - mOrigin[axis])*mInvSize;
+  }
+	
+  bool
+  ComputeMinMaxT(const AxisAlignedBox3 &box,
+								 float &tmin,
+								 float &tmax) const;
+	
+  bool
+  Intersects(const AxisAlignedBox3 &box,
+						 float &tmin,
+						 float &tmax) const;
+	
+  bool
+  IntersectsSphere(const Vector3 &center,
+									 const float sqrRadius,
+									 Vector3 &point,
+									 float &t) const;
+	
+  void
+  Translate(const Vector3 &translation) {
+    mOrigin += translation;
+    mTermination += translation;
+  }
+	
+  bool HasPosDir(const int axis) const { return mFlags & (1<<axis); }
+
+  char Flags() const { return mFlags;}
+  void SetFlags(char orFlag) { mFlags |= orFlag;}
+  
+  bool IsActive() const { return mRefCount>0; }
+	
+  // reference counting for leaf nodes
+  int RefCount() const { return mRefCount; }
+  int Ref() { return mRefCount++; }
+	
+  void ScheduleForRemoval() { if (mRefCount>0) mRefCount = -mRefCount; }
+  bool ScheduledForRemoval() const { return mRefCount<0; }
+  void Unref() {
+    if (mRefCount > 0)
+      mRefCount--;
+    else
+      if (mRefCount < 0)
+				mRefCount++;
+      else {
+				cerr<<"Trying to unref already deleted ray!"<<endl;
+				exit(1);
+      }
+  }
+};
+
+
+
+typedef  vector<VssRay *> VssRayContainer;
+
+
+
+#endif
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Environment.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Environment.h	(revision 162)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Environment.h	(revision 372)
@@ -311,6 +311,6 @@
   /// returns named option as a character string.
   bool GetStringValue(const char *name,
-		      char *value,
-		      const bool isFatal = false) const;
+											char *value,
+											const bool isFatal = false) const;
   
   //@}
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/ExactPreprocessor.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/ExactPreprocessor.cpp	(revision 162)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/ExactPreprocessor.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Polytope.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Polytope.cpp	(revision 260)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Polytope.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/VssPreprocessor.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/VssPreprocessor.cpp	(revision 372)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/VssPreprocessor.cpp	(revision 372)
@@ -0,0 +1,369 @@
+#include "SceneGraph.h"
+#include "KdTree.h"
+#include "VssPreprocessor.h"
+#include "X3dExporter.h"
+#include "Environment.h"
+#include "MutualVisibility.h"
+#include "Polygon3.h"
+#include "ViewCell.h"
+
+VssPreprocessor::VssPreprocessor():
+	mPass(0),
+	mSampleRays(NULL)
+{
+  // this should increase coherence of the samples
+  environment->GetIntValue("Sampling.samplesPerPass", mSamplesPerPass);
+  environment->GetIntValue("Sampling.totalSamples", mTotalSamples);
+  mKdPvsDepth = 100;
+  mStats.open("stats.log");
+
+}
+
+VssPreprocessor::~VssPreprocessor()
+{
+	CLEAR_CONTAINER(mSampleRays);
+}
+
+void
+VssPreprocessor::SetupRay(Ray &ray, 
+													const Vector3 &point, 
+													const Vector3 &direction,
+													const int type)
+{
+  ray.intersections.clear();
+  ray.kdLeaves.clear();
+  ray.meshes.clear();
+
+  //  cout<<point<<" "<<direction<<endl;
+  ray.Init(point, direction, type);
+}
+
+KdNode *
+VssPreprocessor::GetNodeForPvs(KdLeaf *leaf)
+{
+  KdNode *node = leaf;
+  while (node->mParent && node->mDepth > mKdPvsDepth)
+    node = node->mParent;
+  return node;
+}
+
+
+int
+VssPreprocessor::AddNodeSamples(const Ray &ray,
+																		 Intersectable *sObject,
+																		 Intersectable *tObject
+																		 )
+{
+  int contributingSamples = 0;
+  int j;
+	int objects = 0;
+	if (sObject)
+		objects++;
+	if (tObject)
+		objects++;
+
+	if (objects) {
+		for (j=0; j < ray.kdLeaves.size(); j++) {
+			KdNode *node = GetNodeForPvs( ray.kdLeaves[j] );
+			if (sObject)
+				contributingSamples += sObject->mKdPvs.AddSample(node);
+			if (tObject)
+				contributingSamples += tObject->mKdPvs.AddSample(node);
+		}
+	}
+	
+	for (j=1; j < ((int)ray.kdLeaves.size() - 1); j++) {
+		ray.kdLeaves[j]->AddPassingRay2(ray,
+																		objects,
+																		ray.kdLeaves.size()
+																		);
+  }
+  
+  return contributingSamples;
+}
+
+
+
+
+int
+VssPreprocessor::CastRay(Intersectable *object,
+												 Ray &ray
+												 )
+{
+	int sampleContributions = 0;
+
+	long t1 = GetRealTime();
+	// cast ray to KD tree to find intersection with other objects
+	mKdTree->CastRay(ray);
+	long t2 = GetRealTime();
+
+	if (0 && object && object->GetId() > 2197) {
+		object->Describe(cout)<<endl;
+		cout<<ray<<endl;
+	}
+
+	if (ray.kdLeaves.size()) {
+		Intersectable *terminator =
+			ray.intersections.size() ? ray.intersections[0].mObject: NULL;
+		
+		sampleContributions += AddNodeSamples(ray,
+																					object,
+																					terminator);
+	}
+	
+	
+	return sampleContributions;
+}
+
+
+KdNode *
+VssPreprocessor::GetNodeToSample(Intersectable *object)
+{
+	int pvsSize = object->mKdPvs.GetSize();
+	KdNode *nodeToSample = NULL;
+	
+	bool samplePvsBoundary = false;
+	if (pvsSize && samplePvsBoundary) {
+		// this samples the nodes from the boundary of the current PVS
+		// mail all nodes from the pvs
+		Intersectable::NewMail();
+		KdPvsMap::iterator i = object->mKdPvs.mEntries.begin();
+		
+		for (; i != object->mKdPvs.mEntries.end(); i++) {
+			KdNode *node = (*i).first;
+			node->Mail();
+		}
+		
+		int maxTries = 2*pvsSize;
+		Debug << "Finding random neighbour" << endl;	
+		for (int tries = 0; tries < 10; tries++) {
+			int index = RandomValue(0, pvsSize - 1);
+			KdPvsData data;
+			KdNode *node;
+			object->mKdPvs.GetData(index, node, data);
+			nodeToSample = mKdTree->FindRandomNeighbor(node, true);
+			if (nodeToSample)
+				break;
+		}
+	} else {
+		// just pickup a random node
+		//		nodeToSample = mKdTree->GetRandomLeaf(Plane3(normal, point));
+		nodeToSample = mKdTree->GetRandomLeaf();
+	}
+	return nodeToSample;
+}
+
+bool
+VssPreprocessor::ComputeVisibility()
+{
+  
+  // pickup an object
+  ObjectContainer objects;
+  
+  mSceneGraph->CollectObjects(&objects);
+
+  Vector3 point, normal, direction;
+  Ray ray;
+
+  long startTime = GetTime();
+  
+  int i;
+  int totalSamples = 0;
+
+  int pvsOut = Min((int)objects.size(), 10);
+
+  vector<Ray> rays[10];
+
+  while (totalSamples < mTotalSamples) {
+		int passContributingSamples = 0;
+		int passSampleContributions = 0;
+		int passSamples = 0;
+		int index = 0;
+		
+		int reverseSamples = 0;
+		
+			
+		//cout << "totalSamples: "  << totalSamples << endl;
+
+		for (i = 0; i < objects.size(); i++) {
+						
+			KdNode *nodeToSample = NULL;
+			Intersectable *object = objects[i];
+		
+			int pvsSize = 0;
+			if (ViewCell::sHierarchy == ViewCell::KD)
+				pvsSize = object->mKdPvs.GetSize();
+						
+			
+			
+			int faceIndex = object->GetRandomSurfacePoint(point, normal);
+			
+			bool viewcellSample = true;
+			int sampleContributions;
+			bool debug = false; //(object->GetId() >= 2199);
+			if (viewcellSample) {
+				//mKdTree->GetRandomLeaf(Plane3(normal, point));
+
+				nodeToSample = GetNodeToSample(object);
+
+				for (int k=0; k < mSamplesPerPass; k++) {
+					bool reverseSample = false;
+					
+
+					if (nodeToSample) {
+						AxisAlignedBox3 box = mKdTree->GetBox(nodeToSample);
+						Vector3 pointToSample = box.GetRandomPoint();
+						//						pointToSample.y = 0.9*box.Min().y + 0.1*box.Max().y;
+						if (object->GetRandomVisibleSurfacePoint( point, normal, pointToSample, 3 )) {
+							direction = pointToSample - point;
+						} else {
+							reverseSamples++;
+							reverseSample = true;
+							direction = point - pointToSample;
+							point = pointToSample;
+						}
+					}
+					else {
+						direction = UniformRandomVector(normal);
+					}
+					
+					// construct a ray
+					SetupRay(ray, point, direction, Ray::LOCAL_RAY);
+					
+					sampleContributions = CastRay(reverseSample ? NULL : object, ray);
+					
+					//-- CORR matt: put block inside loop
+					if (sampleContributions) {
+						passContributingSamples ++;
+						passSampleContributions += sampleContributions;
+					}
+
+					if ( i < pvsOut )
+						rays[i].push_back(ray);
+		
+					if (!ray.intersections.empty()) {
+						// check whether we can add this to the rays
+						for (int j = 0; j < pvsOut; j++) {
+							if (objects[j] == ray.intersections[0].mObject) {
+								rays[j].push_back(ray);
+							}
+						}
+					}
+					//-------------------
+				}
+			} else {
+				// edge samples
+				// get random visible mesh
+				//				object->GetRandomVisibleMesh(Plane3(normal, point));
+			}
+	
+			// CORR matt: must add all samples
+			passSamples += mSamplesPerPass;
+		}
+	
+		totalSamples += passSamples; 
+		
+		//    if (pass>10)
+		//      HoleSamplingPass();
+    
+		mPass++;
+
+		int pvsSize = 0;
+        
+		for (i=0; i < objects.size(); i++) {
+			Intersectable *object = objects[i];
+			pvsSize += object->mKdPvs.GetSize();
+		}
+		
+		
+		float avgRayContrib = (passContributingSamples > 0) ? 
+			passSampleContributions/(float)passContributingSamples : 0;
+
+		cout << "#Pass " << mPass << " : t = " << TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
+		cout << "#TotalSamples=" << totalSamples/1000 
+				 << "k   #SampleContributions=" << passSampleContributions << " (" 
+				 << 100*passContributingSamples/(float)passSamples<<"%)" << " avgPVS="
+				 << pvsSize/(float)objects.size() << endl 
+				 << "avg ray contrib=" << avgRayContrib << endl
+				 << "reverse samples [%]" << reverseSamples/(float)passSamples*100.0f << endl;
+
+		mStats <<
+			"#Pass\n" <<mPass<<endl<<
+			"#Time\n" << TimeDiff(startTime, GetTime())*1e-3 << endl<<
+			"#TotalSamples\n" << totalSamples<< endl<<
+			"#SampleContributions\n" << passSampleContributions << endl << 
+			"#PContributingSamples\n"<<100*passContributingSamples/(float)passSamples<<endl <<
+			"#AvgPVS\n"<< pvsSize/(float)objects.size() << endl <<
+			"#AvgRayContrib\n" << avgRayContrib << endl;
+	}
+	
+	if (ViewCell::sHierarchy == ViewCell::KD)	
+		cout << "#totalPvsSize=" << mKdTree->CollectLeafPvs() << endl;
+  
+	
+  //  HoleSamplingPass();
+	if (0) {
+		Exporter *exporter = Exporter::GetExporter("ray-density.x3d"); 
+		exporter->SetExportRayDensity(true);
+		exporter->ExportKdTree(*mKdTree);
+
+		delete exporter;
+	}
+  
+	bool exportRays = false;
+	if (exportRays) {
+		Exporter *exporter = NULL;
+		exporter = Exporter::GetExporter("sample-rays.x3d");
+		exporter->SetWireframe();
+		exporter->ExportKdTree(*mKdTree);
+
+		for (i=0; i < pvsOut; i++) 
+		exporter->ExportRays(rays[i], 1000, RgbColor(1, 0, 0));
+		exporter->SetFilled();
+		  
+		delete exporter;
+	}
+
+	//-- several visualizations and statistics
+	if (1) {
+		for (int k=0; k < pvsOut; k++) {
+			Intersectable *object = objects[k];
+			char s[64];
+			sprintf(s, "sample-pvs%04d.x3d", k);
+			Exporter *exporter = Exporter::GetExporter(s);
+			exporter->SetWireframe();
+			
+			
+			KdPvsMap::iterator i = object->mKdPvs.mEntries.begin();
+			Intersectable::NewMail();
+			
+			// avoid adding the object to the list
+			object->Mail();
+			ObjectContainer visibleObjects;
+			
+			for (; i != object->mKdPvs.mEntries.end(); i++) 
+				{
+					KdNode *node = (*i).first;
+					exporter->ExportBox(mKdTree->GetBox(node));
+					mKdTree->CollectObjects(node, visibleObjects);
+				}
+			
+			exporter->ExportRays(rays[k], 1000, RgbColor(0, 1, 0));
+			exporter->SetFilled();
+			
+			for (int j = 0; j < visibleObjects.size(); j++)
+				exporter->ExportIntersectable(visibleObjects[j]);
+			
+			
+			Material m;
+			m.mDiffuseColor = RgbColor(1, 0, 0);
+			exporter->SetForcedMaterial(m);
+			exporter->ExportIntersectable(object);
+			
+			delete exporter;
+		}
+  }
+  
+  return true;
+}
+
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/KdTree.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/KdTree.h	(revision 369)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/KdTree.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/main.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/main.cpp	(revision 366)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/main.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/SamplingPreprocessor.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/SamplingPreprocessor.h	(revision 371)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/SamplingPreprocessor.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/VssRay.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/VssRay.cpp	(revision 372)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/VssRay.cpp	(revision 372)
@@ -0,0 +1,142 @@
+#include "VssRay.h"
+#include "AxisAlignedBox3.h"
+
+// Static variables
+int
+VssRay::mailID = 0;
+
+
+bool
+VssRay::ComputeMinMaxT(const AxisAlignedBox3 &box,
+															float &tmin,
+															float &tmax) const
+{
+	
+  float minx, maxx, miny, maxy, minz, maxz;
+  Vector3 dir = GetDir();
+  
+  if (fabs(dir.x) < 0.001) {
+    if (box.Min().x < GetOrigin().x && box.Max().x > GetOrigin().x) {
+      minx = -MAXFLOAT;
+      maxx = MAXFLOAT;
+    }
+    else
+      return false;
+  }
+  else {
+    float t1 = (box.Min().x - GetOrigin().x) / dir.x;
+    float t2 = (box.Max().x - GetOrigin().x) / dir.x;
+    if (t1 < t2) {
+      minx = t1;
+      maxx = t2;
+    }
+    else {
+      minx = t2;
+      maxx = t1;
+    }
+    if (maxx < 0)
+      return false;
+  }
+
+  if (fabs(dir.y) < 0.001) {
+    if (box.Min().y < GetOrigin().y && box.Max().y > GetOrigin().y) {
+      miny = -MAXFLOAT;
+      maxy = MAXFLOAT;
+    }
+    else
+      return false;
+  }
+  else {
+    float t1 = (box.Min().y - GetOrigin().y) / dir.y;
+    float t2 = (box.Max().y - GetOrigin().y) / dir.y;
+    if (t1 < t2) {
+      miny = t1;
+      maxy = t2;
+    }
+    else {
+      miny = t2;
+      maxy = t1;
+    }
+    if (maxy < 0.0)
+      return false;
+  }
+	
+  if (fabs(dir.z) < 0.001) {
+    if (box.Min().z < GetOrigin().z && box.Max().z > GetOrigin().z) {
+      minz = -MAXFLOAT;
+      maxz = MAXFLOAT;
+    }
+    else
+      return false;
+  }
+  else {
+    float t1 = (box.Min().z - GetOrigin().z) / dir.z;
+    float t2 = (box.Max().z - GetOrigin().z) / dir.z;
+    if (t1 < t2) {
+      minz = t1;
+      maxz = t2;
+    }
+    else {
+      minz = t2;
+      maxz = t1;
+    }
+    if (maxz < 0.0)
+      return false;
+  }
+
+  tmin = minx;
+  if (miny > tmin)
+    tmin = miny;
+  if (minz > tmin)
+    tmin = minz;
+
+  tmax = maxx;
+  if (maxy < tmax)
+    tmax = maxy;
+  if (maxz < tmax)
+    tmax = maxz;
+
+  return 1; // yes, intersection was found
+}
+
+
+
+bool
+VssRay::Intersects(const AxisAlignedBox3 &box,
+													float &tmin,
+													float &tmax) const
+{
+  if (!ComputeMinMaxT(box, tmin, tmax))
+    return false;
+	
+  if ( tmax < tmin)
+    return false; // the ray passes outside the box
+  
+  if ( tmax < 0.0f)
+    return false; // the intersection is not on the positive halfline
+  
+  return true; // ray hits the box .. origin can be outside or inside the box
+}
+
+
+bool
+VssRay::IntersectsSphere(const Vector3 &center,
+																const float sqrRadius,
+																Vector3 &point,
+																float &t) const
+{
+  // compute ray/plane intersection
+  Vector3 dir = GetDir();
+	
+  t = -DotProd(GetOrigin() - center, dir)*sqr(GetInvSize());
+	
+  if (t < 0.0f)
+    t = 0.0f;
+  else
+    if (t > 1.0f)
+      t = 1.0f;
+  
+  point = GetOrigin() + t*dir;
+  
+  return SqrDistance(point, center) < sqrRadius;
+}
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/X3dExporter.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/X3dExporter.cpp	(revision 360)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/X3dExporter.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/common.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/common.cpp	(revision 162)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/common.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Preprocessor.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Preprocessor.h	(revision 366)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Preprocessor.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/preprocessor.pro
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/preprocessor.pro	(revision 364)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/preprocessor.pro	(revision 372)
@@ -28,5 +28,5 @@
 
 # Input
-HEADERS += Halton.h
+HEADERS += Halton.h VssRay.h VssPreprocessor.h
 
 
@@ -38,5 +38,5 @@
 Exporter.cpp Camera.cpp X3dParser.cpp MeshKdTree.cpp Pvs.cpp \
 MutualVisibility.cpp Triangle3.cpp Rectangle3.cpp Plane3.cpp Polygon3.cpp \
-ViewCell.cpp ViewCellBsp.cpp Halton.cpp
+ViewCell.cpp ViewCellBsp.cpp Halton.cpp VssRay.cpp VssTree.cpp VssPreprocessor.cpp
 
 
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Polygon3.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Polygon3.cpp	(revision 367)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Polygon3.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/Containers.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/Containers.h	(revision 318)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/Containers.h	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/SamplingPreprocessor.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/SamplingPreprocessor.cpp	(revision 371)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/SamplingPreprocessor.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/AxisAlignedBox3.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/AxisAlignedBox3.cpp	(revision 359)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/AxisAlignedBox3.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/VssPreprocessor.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/VssPreprocessor.h	(revision 372)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/VssPreprocessor.h	(revision 372)
@@ -0,0 +1,55 @@
+#ifndef _VssPreprocessor_H__
+#define _VssPreprocessor_H__
+
+#include <fstream>
+using namespace std;
+
+#include "Preprocessor.h"
+
+
+/** Sampling based visibility preprocessing. The implementation is based on heuristical
+    sampling of view space */
+class VssPreprocessor : public Preprocessor {
+public:
+	int mPass;
+
+  int mSamplesPerPass;
+  int mTotalSamples;
+  int mKdPvsDepth;
+  ofstream mStats;
+  ObjectContainer mObjects;
+  RayContainer mSampleRays;
+	
+  VssPreprocessor();
+
+  ~VssPreprocessor();
+
+  virtual bool ComputeVisibility();
+
+  void
+  SetupRay(Ray &ray, 
+		   const Vector3 &point, 
+		   const Vector3 &direction,
+		   const int type);
+
+  KdNode *
+  GetNodeForPvs(KdLeaf *leaf);
+
+	int
+	AddNodeSamples(const Ray &ray,
+								 Intersectable *sObject,
+								 Intersectable *tObject
+								 );
+
+	int
+	CastRay(Intersectable *object,
+					Ray &ray);
+	
+	KdNode *
+	GetNodeToSample(Intersectable *object);
+
+};
+
+
+
+#endif
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/gzstream.cpp
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/gzstream.cpp	(revision 162)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/gzstream.cpp	(revision 372)
Index: /trunk/VUT/GtpVisibilityPreprocessor/src/common.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/common.h	(revision 329)
+++ /trunk/VUT/GtpVisibilityPreprocessor/src/common.h	(revision 372)
Index: unk/VUT/GtpVisibilityPreprocessor/src/VisibilityVector3.h
===================================================================
--- /trunk/VUT/GtpVisibilityPreprocessor/src/VisibilityVector3.h	(revision 162)
+++ 	(revision )
@@ -1,13 +1,0 @@
-#ifndef _VisibilityVector3_H__
-#define _VisibilityVector3_H__
-
-#include "OgreVector3.h"
-
-namespace GtpVisibility {
-  /** Vector3 class currently uses the native Ogre vector when compiled with the Ogre platform */
-  typedef Ogre::Vector3 Vector3;
-}
-
-
-
-#endif
