Index: trunk/VUT/GtpVisibilityPreprocessor/src/VspBspTree.h
===================================================================
--- trunk/VUT/GtpVisibilityPreprocessor/src/VspBspTree.h	(revision 437)
+++ trunk/VUT/GtpVisibilityPreprocessor/src/VspBspTree.h	(revision 437)
@@ -0,0 +1,787 @@
+#ifndef _VspBspTree_H__
+#define _VspBspTree_H__
+
+#include "Mesh.h"
+#include "Containers.h"
+#include "Polygon3.h"
+#include <stack>
+#include "Statistics.h"
+#include "VssRay.h"
+#include "RayInfo.h"
+
+class ViewCell;
+class BspViewCell;
+class Plane3;
+class VspBspTree;  
+class VspBspInterior;
+class VspBspNode;
+class AxisAlignedBox3;
+class Ray;
+
+class VspBspNodeGeometry
+{
+public:
+	VspBspNodeGeometry()
+	{};  
+
+	~VspBspNodeGeometry();
+
+	float GetArea() const;
+	
+	/** Computes new cell based on the old cell definition and a new split plane
+		@param side indicates which side of the halfspace 
+	*/
+	void SplitGeometry(VspBspNodeGeometry &front, 
+					   VspBspNodeGeometry &back,
+					   const VspBspTree &tree,
+					   const Plane3 &splitPlane) const;
+
+	Polygon3 *SplitPolygon(Polygon3 *poly, const VspBspTree &tree) const;
+
+	PolygonContainer mPolys;
+};
+
+/** Data structure used for optimized ray casting.
+*/
+struct VspBspRayTraversalData 
+{
+    VspBspNode *mNode;
+    Vector3 mExitPoint;
+    float mMaxT;
+    
+    VspBspRayTraversalData() {}
+
+    VspBspRayTraversalData(VspBspNode *n, const Vector3 &extp, const float maxt):
+    mNode(n), mExitPoint(extp), mMaxT(maxt) 
+	{}
+};
+
+class VspBspTreeStatistics: public StatisticsBase
+{
+public:
+	// total number of nodes
+	int nodes;
+	// number of splits
+	int splits;
+	// totals number of rays
+	int rays;
+	// maximal reached depth
+	int maxDepth;
+	// minimal depth
+	int minDepth;
+	
+	// max depth nodes
+	int maxDepthNodes;
+	// minimum depth nodes
+	int minDepthNodes;
+	// max depth nodes
+	int minPvsNodes;
+	// nodes with minimum PVS
+	int minRaysNodes;
+	// max ray contribution nodes
+	int maxRayContribNodes;
+	// minimum area nodes
+	int minAreaNodes;
+
+	// max number of rays per node
+	int maxObjectRefs;
+	// accumulated depth (used to compute average)
+	int accumDepth;
+	// number of initial polygons
+	int polys;
+	/// samples contributing to pvs
+	int contributingSamples;
+	/// sample contributions to pvs
+	int sampleContributions;
+	/// largest pvs
+	int largestPvs;
+
+	// Constructor
+	VspBspTreeStatistics() 
+	{
+		Reset();
+	}
+
+	int Nodes() const {return nodes;}
+	int Interior() const { return nodes / 2; }
+	int Leaves() const { return (nodes / 2) + 1; }
+	
+	// TODO: computation wrong
+	double AvgDepth() const { return accumDepth / (double)Leaves();}; 
+  
+	void Reset() 
+	{
+		nodes = 0;
+		splits = 0;
+		
+		maxDepth = 0;
+		minDepth = 99999;
+		polys = 0;
+		accumDepth = 0;
+        
+		maxDepthNodes = 0;
+		minPvsNodes = 0;
+		minRaysNodes = 0;
+		maxRayContribNodes = 0;
+		minAreaNodes = 0;
+
+		contributingSamples = 0;
+		sampleContributions = 0;
+	}
+
+	void Print(ostream &app) const;
+
+	friend ostream &operator<<(ostream &s, const VspBspTreeStatistics &stat) 
+	{
+		stat.Print(s);
+		return s;
+	} 
+};
+
+class VspBspViewCellsStatistics: public StatisticsBase
+{
+public:
+
+	/// number of view cells
+	int viewCells;
+
+	/// size of the PVS
+	int pvs;
+  
+	/// largest PVS of all view cells
+	int maxPvs;
+
+	/// smallest PVS of all view cells
+	int minPvs;
+
+	/// view cells with empty PVS
+	int emptyPvs;
+
+	/// number of bsp leaves covering the view space
+	int bspLeaves;
+
+	/// largest number of leaves covered by one view cell  
+	int maxVspBspLeaves;
+
+    // Constructor
+	VspBspViewCellsStatistics() 
+	{
+		Reset();
+	}
+
+	double AvgVspBspLeaves() const {return (double)bspLeaves / (double)viewCells;};
+	double AvgPvs() const {return (double)pvs / (double)viewCells;};
+  
+	void Reset() 
+	{
+		viewCells = 0;	
+		pvs = 0;  
+		maxPvs = 0;
+
+		minPvs = 999999;
+		emptyPvs = 0;
+		bspLeaves = 0;
+		maxVspBspLeaves = 0;
+	}
+
+	void Print(ostream &app) const;
+
+	friend ostream &operator<<(ostream &s, const VspBspViewCellsStatistics &stat) 
+	{
+		stat.Print(s);
+		return s;
+	}  
+};
+
+/**
+    VspBspNode abstract class serving for interior and leaf node implementation
+*/
+class VspBspNode 
+{
+	friend class VspBspTree;
+
+public:
+	VspBspNode();
+	virtual ~VspBspNode(){};
+	VspBspNode(VspBspInterior *parent);
+
+	/** Determines whether this node is a leaf or not
+	@return true if leaf
+	*/
+	virtual bool IsLeaf() const = 0;
+
+	/** Determines whether this node is a root
+	@return true if root
+	*/
+	virtual bool IsRoot() const;
+
+	/** Returns parent node.
+	*/
+	VspBspInterior *GetParent();
+
+	/** Sets parent node.
+	*/
+	void SetParent(VspBspInterior *parent);
+
+	
+	static int sMailId;
+	int mMailbox;
+  
+	void Mail() { mMailbox = sMailId; }
+	static void NewMail() { ++ sMailId; }
+	bool Mailed() const { return mMailbox == sMailId; }
+
+protected:
+
+	/// parent of this node
+	VspBspInterior *mParent;
+};
+
+/** BSP interior node implementation 
+*/
+class VspBspInterior : public VspBspNode 
+{
+	friend class VspBspTree;
+public:
+	/** Standard contructor taking split plane as argument.
+	*/
+	VspBspInterior(const Plane3 &plane);
+	~VspBspInterior();
+	/** @return false since it is an interior node 
+	*/
+	bool IsLeaf() const;
+
+	VspBspNode *GetBack();
+	VspBspNode *GetFront();
+
+	Plane3 *GetPlane();
+
+	void ReplaceChildLink(VspBspNode *oldChild, VspBspNode *newChild);
+	void SetupChildLinks(VspBspNode *b, VspBspNode *f);
+
+	/** Splits polygons with respect to the split plane.
+		@param polys the polygons to be split. the polygons are consumed and
+			   distributed to the containers frontPolys, backPolys, coincident.
+		@param frontPolys returns the polygons in the front of the split plane
+		@param backPolys returns the polygons in the back of the split plane
+		@param coincident returns the polygons coincident to the split plane
+
+		@returns the number of splits	
+	*/
+	int SplitPolygons(PolygonContainer &polys, 
+					  PolygonContainer &frontPolys, 
+					  PolygonContainer &backPolys, 
+					  PolygonContainer &coincident);
+
+	friend ostream &operator<<(ostream &s, const VspBspInterior &A)
+	{
+		return s << A.mPlane;
+	}
+
+protected:
+
+	/// Splitting plane corresponding to this node
+	Plane3 mPlane;
+	/// back node
+	VspBspNode *mBack;
+	/// front node
+	VspBspNode *mFront;
+};
+
+/** BSP leaf node implementation.
+*/
+class VspBspLeaf : public VspBspNode 
+{
+	friend class VspBspTree;
+
+public:
+	VspBspLeaf();
+	VspBspLeaf(BspViewCell *viewCell);
+	VspBspLeaf(VspBspInterior *parent);
+	VspBspLeaf(VspBspInterior *parent, BspViewCell *viewCell);
+
+	/** @return true since it is an interior node 
+	*/
+	bool IsLeaf() const;
+	
+	/** Returns pointer of view cell.
+	*/
+	BspViewCell *GetViewCell() const;
+
+	/** Sets pointer to view cell.
+	*/
+	void SetViewCell(BspViewCell *viewCell);
+
+	/** Adds ray sample contributions to the PVS.
+		@param sampleContributions the number contributions of the samples
+		@param contributingSampels the number of contributing rays
+		
+	*/
+	void AddToPvs(const RayInfoContainer &rays, 
+				  int &sampleContributions,
+				  int &contributingSamples, 
+				  bool storeLeavesWithRays = false);
+
+protected:
+
+	/// if NULL this does not correspond to feasible viewcell
+	BspViewCell *mViewCell;
+};
+
+/** Implementation of the view cell BSP tree. 
+*/
+class VspBspTree 
+{
+public:
+	
+	/** Additional data which is passed down the BSP tree during traversal.
+	*/
+	struct VspBspTraversalData
+	{  
+		/// the current node
+		VspBspNode *mNode;
+		/// polygonal data for splitting
+		PolygonContainer *mPolygons;
+		/// current depth
+		int mDepth;
+		/// the view cell associated with this subdivsion
+		ViewCell *mViewCell;
+		/// rays piercing this node
+		RayInfoContainer *mRays;
+		/// area of current node
+		float mArea;
+		/// geometry of current node induced by split planes
+		VspBspNodeGeometry *mGeometry;
+
+		/// pvs size
+		int mPvs;
+		
+		/** Returns average ray contribution.
+		*/
+		float GetAvgRayContribution() const
+		{
+			return (float)mPvs / ((float)mRays->size() + Limits::Small);
+		}
+
+
+		VspBspTraversalData():
+		mNode(NULL),
+		mPolygons(NULL),
+		mDepth(0),
+		mViewCell(NULL),
+		mRays(NULL),
+		mPvs(0),
+		mArea(0.0),
+		mGeometry(NULL)
+		{}
+		
+		VspBspTraversalData(VspBspNode *node, 
+						 PolygonContainer *polys, 
+						 const int depth, 
+						 ViewCell *viewCell,
+						 RayInfoContainer *rays,
+						 int pvs,
+						 float area,
+						 VspBspNodeGeometry *cell): 
+		mNode(node), 
+		mPolygons(polys), 
+		mDepth(depth), 
+		mViewCell(viewCell),
+		mRays(rays),
+		mPvs(pvs),
+		mArea(area),
+		mGeometry(cell)
+		{}
+    };
+	
+	typedef std::stack<VspBspTraversalData> VspBspTraversalStack;
+
+	/** Default constructor creating an empty tree.
+		@param viewCell view cell corresponding to unbounded space
+	*/ 
+	VspBspTree(BspViewCell *viewCell);
+
+	~VspBspTree();
+
+	const VspBspTreeStatistics &GetStatistics() const; 
+  
+
+	/** Constructs the tree from a given set of rays.
+		@param sampleRays the set of sample rays the construction is based on
+		@param viewCells if not NULL, new view cells are 
+		created in the leafs and stored in the conatainer
+	*/
+	void Construct(const VssRayContainer &sampleRays);
+
+	/** Returns list of BSP leaves.
+	*/
+	void CollectLeaves(vector<VspBspLeaf *> &leaves) const;
+
+	/** Returns box which bounds the whole tree.
+	*/
+	AxisAlignedBox3 GetBoundingBox()const;
+
+	/** Returns root of BSP tree.
+	*/
+	VspBspNode *GetRoot() const;
+
+	/** Exports VspBsp tree to file.
+	*/
+	bool Export(const string filename);
+
+	/** Collects the leaf view cells of the tree
+		@param viewCells returns the view cells 
+	*/
+	void CollectViewCells(ViewCellContainer &viewCells) const;
+
+	/** A ray is cast possible intersecting the tree.
+		@param the ray that is cast.
+		@returns the number of intersections with objects stored in the tree.
+	*/
+	int CastRay(Ray &ray);
+
+	/** Set to true if new view cells shall be generated in each leaf.
+	*/
+	void SetGenerateViewCells(int generateViewCells);
+
+	/// bsp tree construction types
+	enum {FROM_INPUT_VIEW_CELLS, FROM_SCENE_GEOMETRY, FROM_SAMPLES};
+
+	/** Returns statistics.
+	*/
+	VspBspTreeStatistics &GetStat();
+
+	/** finds neighbouring leaves of this tree node.
+	*/
+	int FindNeighbors(VspBspNode *n, vector<VspBspLeaf *> &neighbors, 
+					  const bool onlyUnmailed) const;
+
+	/** Constructs geometry associated with the half space intersections 
+		leading to this node.
+	*/
+	void ConstructGeometry(VspBspNode *n, PolygonContainer &cell) const;
+
+	/** Construct geometry and stores it in a geometry node container.
+	*/
+	void ConstructGeometry(VspBspNode *n, VspBspNodeGeometry &cell) const;
+
+	/** Returns random leaf of BSP tree.
+		@param halfspace defines the halfspace from which the leaf is taken.
+	*/
+	VspBspLeaf *GetRandomLeaf(const Plane3 &halfspace);
+
+	/** Returns random leaf of BSP tree.
+		@param onlyUnmailed if only unmailed leaves should be returned.
+	*/
+	VspBspLeaf *GetRandomLeaf(const bool onlyUnmailed = false);
+
+	/** Returns true if merge criteria are reached.
+	*/
+    bool ShouldMerge(VspBspLeaf *front, VspBspLeaf *back) const;
+
+	/** Merges view cells based on some criteria
+	    E.g., empty view cells can pe purged, view cells which have 
+		a very similar PVS can be merged to one larger view cell.
+
+		@returns true if merge was successful.
+	*/
+	bool MergeViewCells(VspBspLeaf *front, VspBspLeaf *back) const;
+
+	/** Traverses tree and counts all view cells as well as their PVS size.
+	*/
+	void EvaluateViewCellsStats(VspBspViewCellsStatistics &stat) const;
+
+
+protected:
+
+	// --------------------------------------------------------------
+	// For sorting objects
+	// --------------------------------------------------------------
+	struct SortableEntry
+	{
+		enum {POLY_MIN, POLY_MAX};
+    
+		int type;
+		float value;
+		Polygon3 *poly;
+		SortableEntry() {}
+		SortableEntry(const int t, const float v, Polygon3 *poly): 
+		type(t), value(v), poly(poly) {}
+		
+		bool operator<(const SortableEntry &b) const 
+		{
+			return value < b.value;
+		}  
+	};
+
+	/** Evaluates tree stats in the BSP tree leafs.
+	*/
+	void EvaluateLeafStats(const VspBspTraversalData &data);
+
+	/** Subdivides node with respect to the traversal data.
+	    @param tStack current traversal stack
+		@param tData traversal data also holding node to be subdivided
+		@returns new root of the subtree
+	*/
+	VspBspNode *Subdivide(VspBspTraversalStack &tStack, VspBspTraversalData &tData);
+
+	/** Constructs the tree from the given list of polygons and rays.
+		@param polys stores set of polygons on which subdivision may be based
+		@param rays storesset of rays on which subdivision may be based
+	*/
+	void Construct(PolygonContainer *polys, RayInfoContainer *rays);
+
+	/** Selects the best possible splitting plane. 
+		@param leaf the leaf to be split
+		@param polys the polygon list on which the split decition is based
+		@param rays ray container on which selection may be based
+		@note the polygons can be reordered in the process
+		@returns the split plane
+	*/
+	Plane3 SelectPlane(VspBspLeaf *leaf, 
+					   VspBspTraversalData &data);
+
+	
+	/** Strategies where the effect of the split plane is tested
+	    on all input rays.
+
+		@returns the cost of the candidate split plane
+	*/
+	float SplitPlaneCost(const Plane3 &candidatePlane, 
+						 const VspBspTraversalData &data);
+
+
+
+	/** Subdivide leaf.
+		@param leaf the leaf to be subdivided
+		
+		@param polys the polygons to be split
+		@param frontPolys returns the polygons in front of the split plane
+		@param backPolys returns the polygons in the back of the split plane
+		
+		@param rays the polygons to be filtered
+		@param frontRays returns the polygons in front of the split plane
+		@param backRays returns the polygons in the back of the split plane
+
+		@returns the root of the subdivision
+	*/
+
+	VspBspInterior *SubdivideNode(VspBspTraversalData &tData,
+							   VspBspTraversalData &frontData,
+							   VspBspTraversalData &backData,
+							   PolygonContainer &coincident);
+
+	/** Selects the split plane in order to construct a tree with
+		certain characteristics (e.g., balanced tree, least splits, 
+		2.5d aligned)
+		@param polygons container of polygons
+		@param rays bundle of rays on which the split can be based
+	*/
+	Plane3 SelectPlaneHeuristics(VspBspLeaf *leaf,
+								 VspBspTraversalData &data);
+
+	/** Extracts the meshes of the objects and adds them to polygons. 
+		Adds object aabb to the aabb of the tree.
+		@param maxPolys the maximal number of objects to be stored as polygons
+		@returns the number of polygons
+	*/
+	int AddToPolygonSoup(const ObjectContainer &objects, 
+						 PolygonContainer &polys, 
+						 int maxObjects = 0);
+
+	/** Extracts the meshes of the view cells and and adds them to polygons.
+		Adds view cell aabb to the aabb of the tree.
+		@param maxPolys the maximal number of objects to be stored as polygons
+		@returns the number of polygons
+	*/
+	int AddToPolygonSoup(const ViewCellContainer &viewCells, 
+						 PolygonContainer &polys, 
+						 int maxObjects = 0);
+
+	/** Extract polygons of this mesh and add to polygon container.
+		@param mesh the mesh that drives the polygon construction
+		@param parent the parent intersectable this polygon is constructed from
+		@returns number of polygons
+	*/
+	int AddMeshToPolygons(Mesh *mesh, PolygonContainer &polys, MeshInstance *parent);
+
+	/** returns next candidate index and reorders polygons so no candidate is chosen two times
+		@param the current candidate index
+		@param max the range of candidates
+	*/
+	int GetNextCandidateIdx(int currentIdx, PolygonContainer &polys);
+	
+	/** Computes best cost ratio for the suface area heuristics for axis aligned
+		splits. This heuristics minimizes the cost for ray traversal.
+		@param polys the polygons guiding the ratio computation
+		@param box the bounding box of the leaf
+		@param axis the current split axis
+		@param position returns the split position
+		@param objectsBack the number of objects in the back of the split plane
+		@param objectsFront the number of objects in the front of the split plane
+	*/
+	float BestCostRatio(const PolygonContainer &polys,
+						const AxisAlignedBox3 &box,
+						const int axis,
+						float &position,
+						int &objectsBack,
+						int &objectsFront) const;
+	
+	/** Sorts split candidates for surface area heuristics for axis aligned splits.
+		@param polys the input for choosing split candidates
+		@param axis the current split axis
+		@param splitCandidates returns sorted list of split candidates
+	*/
+	void SortSplitCandidates(const PolygonContainer &polys, 
+							 const int axis, 
+							 vector<SortableEntry> &splitCandidates) const;
+
+	/** Selects an axis aligned split plane.
+		Returns true if split is valied
+	*/
+	bool SelectAxisAlignedPlane(Plane3 &plane, const PolygonContainer &polys) const;
+
+	/** Bounds ray and returns minT and maxT.
+		@returns true if ray hits BSP tree bounding box
+	*/
+	bool BoundRay(const Ray &ray, float &minT, float &maxT) const;
+
+	/** Subdivides the rays into front and back rays according to the split plane.
+		
+		@param plane the split plane
+		@param rays contains the rays to be split. The rays are 
+			   distributed into front and back rays.
+		@param frontRays returns rays on the front side of the plane
+		@param backRays returns rays on the back side of the plane
+		
+		@returns the number of splits
+	*/
+	int SplitRays(const Plane3 &plane,
+				  RayInfoContainer &rays, 
+			      RayInfoContainer &frontRays, 
+				  RayInfoContainer &backRays);
+
+
+	/** Extracts the split planes representing the space bounded by node n.
+	*/
+	void ExtractHalfSpaces(VspBspNode *n, vector<Plane3> &halfSpaces) const;
+
+	/** Adds the object to the pvs of the front and back leaf with a given classification.
+
+		@param obj the object to be added
+		@param cf the ray classification regarding the split plane
+		@param frontPvs returns the PVS of the front partition
+		@param backPvs returns the PVS of the back partition
+	
+	*/
+	void AddObjToPvs(Intersectable *obj, const int cf, int &frontPvs, int &backPvs) const;
+	
+	/** Computes PVS size induced by the rays.
+	*/
+	int ComputePvsSize(const RayInfoContainer &rays) const;
+
+	/** Returns true if tree can be terminated.
+	*/
+	inline bool TerminationCriteriaMet(const VspBspTraversalData &data) const;
+
+	/** Computes accumulated ray lenght of this rays.
+	*/
+	float AccumulatedRayLength(const RayInfoContainer &rays) const;
+
+
+
+	/// Pointer to the root of the tree
+	VspBspNode *mRoot;
+		
+	VspBspTreeStatistics mStat;
+
+	/// Strategies for choosing next split plane.
+	enum {NO_STRATEGY = 0,
+		  RANDOM_POLYGON = 1, 
+		  AXIS_ALIGNED = 2,
+		  LEAST_SPLITS = 4, 
+		  BALANCED_POLYS = 8,
+		  BALANCED_VIEW_CELLS = 16,
+		  LARGEST_POLY_AREA = 32,
+		  VERTICAL_AXIS = 64,
+		  BLOCKED_RAYS = 128,
+		  LEAST_RAY_SPLITS = 256,
+		  BALANCED_RAYS = 512,
+		  PVS = 1024
+		};
+
+	/// box around the whole view domain
+	AxisAlignedBox3 mBox;
+
+	/// view cell corresponding to unbounded space
+	BspViewCell *mRootCell;
+
+	/// minimal number of rays before subdivision termination
+	int mTermMinRays;
+	/// maximal possible depth
+	int mTermMaxDepth;
+	/// mininum area
+	float mTermMinArea;
+	/// mininum PVS
+	int mTermMinPvs;
+
+	/// minimal number of rays for axis aligned split
+	int mTermMinRaysForAxisAligned;
+	/// minimal number of objects for axis aligned split
+	int mTermMinObjectsForAxisAligned;
+	/// maximal contribution per ray
+	float mTermMaxRayContribution;
+	/// minimal accumulated ray length
+	float mTermMinAccRayLength;
+
+
+	/// strategy to get the best split plane
+	int mSplitPlaneStrategy;
+	/// number of candidates evaluated for the next split plane
+	int mMaxPolyCandidates;
+	/// number of candidates for split planes evaluated using the rays
+	int mMaxRayCandidates;
+	
+	float mCtDivCi;
+
+	/// if intersected leaves should be stored with a sample
+	bool mStoreLeavesWithRays;
+
+	/// axis aligned split criteria
+	float mAaCtDivCi;
+	float mSplitBorder;
+	float mMaxCostRatio;
+
+	// factors guiding the split plane heuristics
+	float mBlockedRaysFactor;
+	float mLeastRaySplitsFactor;
+	float mBalancedRaysFactor;
+	float mPvsFactor;
+	float mLeastSplitsFactor;
+
+	//-- thresholds used for view cells merge
+	int mMinPvsDif;
+	int mMinPvs;
+	int mMaxPvs;
+
+	/// if area or accumulated ray lenght should be used for PVS heuristics
+	bool mPvsUseArea;
+
+private:
+	
+	static const float sLeastRaySplitsTable[5];
+	/** Evaluates split plane classification with respect to the plane's 
+		contribution for balanced rays.
+	*/
+	static const float sBalancedRaysTable[5];
+
+	/// Generates unique ids for PVS criterium
+	static void GenerateUniqueIdsForPvs();
+
+	//-- unique ids for PVS criterium
+	static int sFrontId; 
+	static int sBackId;
+	static int sFrontAndBackId;
+};
+
+#endif
Index: trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.h
===================================================================
--- trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.h	(revision 436)
+++ trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.h	(revision 437)
@@ -333,6 +333,8 @@
 		
 	*/
-	void AddToPvs(const BoundedRayContainer &rays, int &sampleContributions,
-				  int &contributingSamples, bool storeLeavesWithRays = false);
+	void AddToPvs(const BoundedRayContainer &rays, 
+				  int &sampleContributions,
+				  int &contributingSamples, 
+				  bool storeLeavesWithRays = false);
 
 protected:
@@ -471,8 +473,4 @@
 	int CastRay(Ray &ray);
 
-	/** Set to true if new view cells shall be generated in each leaf.
-	*/
-	void SetGenerateViewCells(int generateViewCells);
-
 	/// bsp tree construction types
 	enum {FROM_INPUT_VIEW_CELLS, FROM_SCENE_GEOMETRY, FROM_SAMPLES};
@@ -701,5 +699,5 @@
 						  BspTraversalData &backData, 
 						  const PolygonContainer &coincident,
-						  const Plane3 splitPlane) const;
+						  const Plane3 &splitPlane) const;
 	
 	/** Computes best cost ratio for the suface area heuristics for axis aligned
@@ -767,14 +765,21 @@
 	*/
 	void AddObjToPvs(Intersectable *obj, const int cf, int &frontPvs, int &backPvs) const;
-	
+
+	/** Computes PVS size induced by the rays.
+	*/
 	int ComputePvsSize(const BoundedRayContainer &rays) const;
 
+	/** Returns true if tree can be terminated.
+	*/
 	inline bool TerminationCriteriaMet(const BspTraversalData &data) const;
 
+	/** Computes accumulated ray lenght of this rays.
+	*/
 	float AccumulatedRayLength(BoundedRayContainer &rays) const;
 
-	/// Pointer to the root of the tree
+	/// Pointer to the root of the tree.
 	BspNode *mRoot;
-		
+
+	/// Stores statistics during traversal.
 	BspTreeStatistics mStat;
 
@@ -800,5 +805,5 @@
 	BspViewCell *mRootCell;
 
-	/// should view cells be stored or generated in the leaves?
+	/// if view cells should be generated or the given view cells should be used.
 	bool mGenerateViewCells;
 
Index: trunk/VUT/GtpVisibilityPreprocessor/src/VspBspTree.cpp
===================================================================
--- trunk/VUT/GtpVisibilityPreprocessor/src/VspBspTree.cpp	(revision 437)
+++ trunk/VUT/GtpVisibilityPreprocessor/src/VspBspTree.cpp	(revision 437)
@@ -0,0 +1,2006 @@
+#include "Plane3.h"
+#include "VspBspTree.h"
+#include "Mesh.h"
+#include "common.h"
+#include "ViewCell.h"
+#include "Environment.h"
+#include "Polygon3.h"
+#include "Ray.h"
+#include "AxisAlignedBox3.h"
+#include <stack>
+#include <time.h>
+#include <iomanip>
+#include "Exporter.h"
+#include "Plane3.h"
+
+//-- static members
+/** Evaluates split plane classification with respect to the plane's 
+	contribution for a minimum number of ray splits.
+*/
+const float VspBspTree::sLeastRaySplitsTable[] = {0, 0, 1, 1, 0};
+/** Evaluates split plane classification with respect to the plane's 
+	contribution for balanced rays.
+*/
+const float VspBspTree::sBalancedRaysTable[] = {1, -1, 0, 0, 0};
+
+int VspBspNode::sMailId = 1;
+
+int VspBspTree::sFrontId = 0; 
+int VspBspTree::sBackId = 0;
+int VspBspTree::sFrontAndBackId = 0;
+
+/****************************************************************/
+/*                  class VspBspNode implementation                */
+/****************************************************************/
+
+VspBspNode::VspBspNode(): 
+mParent(NULL)
+{}
+
+VspBspNode::VspBspNode(VspBspInterior *parent): 
+mParent(parent)
+{}
+
+
+bool VspBspNode::IsRoot() const 
+{
+	return mParent == NULL;
+}
+
+VspBspInterior *VspBspNode::GetParent() 
+{ 
+	return mParent; 
+}
+
+void VspBspNode::SetParent(VspBspInterior *parent)
+{
+	mParent = parent;
+}
+
+/****************************************************************/
+/*              class VspBspInterior implementation                */
+/****************************************************************/
+
+
+VspBspInterior::VspBspInterior(const Plane3 &plane): 
+mPlane(plane), mFront(NULL), mBack(NULL)
+{}
+
+VspBspInterior::~VspBspInterior() 
+{ 
+	DEL_PTR(mFront); 
+	DEL_PTR(mBack);
+}
+
+bool VspBspInterior::IsLeaf() const
+{ 
+	return false; 
+}
+
+VspBspNode *VspBspInterior::GetBack() 
+{
+	return mBack;
+}
+
+VspBspNode *VspBspInterior::GetFront() 
+{
+	return mFront;
+}
+
+Plane3 *VspBspInterior::GetPlane()
+{
+	return &mPlane;
+}
+
+void VspBspInterior::ReplaceChildLink(VspBspNode *oldChild, VspBspNode *newChild) 
+{
+	if (mBack == oldChild)
+		mBack = newChild;
+	else
+		mFront = newChild;
+}
+
+void VspBspInterior::SetupChildLinks(VspBspNode *b, VspBspNode *f) 
+{
+    mBack = b;
+    mFront = f;
+}
+
+int VspBspInterior::SplitPolygons(PolygonContainer &polys, 
+								  PolygonContainer &frontPolys, 
+								  PolygonContainer &backPolys, 
+								  PolygonContainer &coincident)
+{
+	int splits = 0;
+
+	while (!polys.empty())
+	{
+		Polygon3 *poly = polys.back();
+		polys.pop_back();
+
+		// classify polygon
+		const int cf = poly->ClassifyPlane(mPlane);
+
+		switch (cf)
+		{
+			case Polygon3::COINCIDENT:
+				coincident.push_back(poly);
+				break;			
+			case Polygon3::FRONT_SIDE:	
+				frontPolys.push_back(poly);
+				break;
+			case Polygon3::BACK_SIDE:
+				backPolys.push_back(poly);
+				break;
+			case Polygon3::SPLIT:
+				backPolys.push_back(poly);
+				frontPolys.push_back(poly);
+				++ splits;
+				break;
+			default:
+                Debug << "SHOULD NEVER COME HERE\n";
+				break;
+		}
+	}
+
+	return splits;
+}
+
+/****************************************************************/
+/*                  class VspBspLeaf implementation                */
+/****************************************************************/
+VspBspLeaf::VspBspLeaf(): mViewCell(NULL)
+{
+}
+
+VspBspLeaf::VspBspLeaf(BspViewCell *viewCell): 
+mViewCell(viewCell)
+{
+}
+
+VspBspLeaf::VspBspLeaf(VspBspInterior *parent): 
+VspBspNode(parent), mViewCell(NULL)
+{}
+
+VspBspLeaf::VspBspLeaf(VspBspInterior *parent, BspViewCell *viewCell): 
+VspBspNode(parent), mViewCell(viewCell)
+{
+}
+
+BspViewCell *VspBspLeaf::GetViewCell() const
+{
+	return mViewCell;
+}
+
+void VspBspLeaf::SetViewCell(BspViewCell *viewCell)
+{
+	mViewCell = viewCell;
+}
+
+bool VspBspLeaf::IsLeaf() const 
+{ 
+	return true; 
+}
+
+void VspBspLeaf::AddToPvs(const RayInfoContainer &rays, 
+						  int &sampleContributions,
+						  int &contributingSamples,
+						  bool storeLeavesWithRays)
+{
+	sampleContributions = 0;
+	contributingSamples = 0;
+
+    RayInfoContainer::const_iterator it, it_end = rays.end();
+
+	// add contributions from samples to the PVS
+	for (it = rays.begin(); it != it_end; ++ it)
+	{
+		int contribution = 0;
+		VssRay *ray = (*it).mRay;
+			
+		if (ray->mTerminationObject)
+			contribution += mViewCell->GetPvs().AddSample(ray->mTerminationObject);
+		
+		if (ray->mOriginObject)
+			contribution += mViewCell->GetPvs().AddSample(ray->mOriginObject);
+
+		if (contribution > 0)
+		{
+			sampleContributions += contribution;
+			++ contributingSamples;
+		}
+	}
+}
+
+
+/****************************************************************/
+/*                  class VspBspTree implementation             */
+/****************************************************************/
+
+VspBspTree::VspBspTree(BspViewCell *viewCell): 
+mRootCell(viewCell), 
+mRoot(NULL), 
+mStoreLeavesWithRays(false),
+mPvsUseArea(true)
+{
+	Randomize(); // initialise random generator for heuristics
+
+	//-- termination criteria for autopartition
+	environment->GetIntValue("VspBspTree.Termination.maxDepth", mTermMaxDepth);
+	environment->GetIntValue("VspBspTree.Termination.minPvs", mTermMinPvs);
+	environment->GetIntValue("VspBspTree.Termination.minRays", mTermMinRays);
+	environment->GetFloatValue("VspBspTree.Termination.minArea", mTermMinArea);	
+	environment->GetFloatValue("VspBspTree.Termination.maxRayContribution", mTermMaxRayContribution);
+	environment->GetFloatValue("VspBspTree.Termination.minAccRayLenght", mTermMinAccRayLength);
+
+	//-- factors for bsp tree split plane heuristics
+	environment->GetFloatValue("VspBspTree.Factor.balancedRays", mBalancedRaysFactor);
+	environment->GetFloatValue("VspBspTree.Factor.pvsFactor", mPvsFactor);
+	environment->GetFloatValue("VspBspTree.Factor.leastSplits" , mLeastSplitsFactor);
+	environment->GetFloatValue("VspBspTree.Termination.ct_div_ci", mCtDivCi);
+
+	//-- termination criteria for axis aligned split
+	environment->GetFloatValue("VspBspTree.Termination.AxisAligned.ct_div_ci", mAaCtDivCi);
+	environment->GetFloatValue("VspBspTree.Termination.AxisAligned.maxCostRatio", mMaxCostRatio);
+	
+	environment->GetIntValue("VspBspTree.Termination.AxisAligned.minRays", 
+							 mTermMinRaysForAxisAligned);
+	
+	//-- partition criteria
+	environment->GetIntValue("VspBspTree.maxPolyCandidates", mMaxPolyCandidates);
+	environment->GetIntValue("VspBspTree.maxRayCandidates", mMaxRayCandidates);
+	environment->GetIntValue("VspBspTree.splitPlaneStrategy", mSplitPlaneStrategy);
+	environment->GetFloatValue("VspBspTree.AxisAligned.splitBorder", mSplitBorder);
+
+	environment->GetFloatValue("VspBspTree.Construction.sideTolerance", Vector3::sDistTolerance);
+	Vector3::sDistToleranceSqrt = Vector3::sDistTolerance * Vector3::sDistTolerance;
+
+	// post processing stuff
+	environment->GetIntValue("ViewCells.PostProcessing.minPvsDif", mMinPvsDif);
+	environment->GetIntValue("ViewCells.PostProcessing.minPvs", mMinPvs);
+	environment->GetIntValue("ViewCells.PostProcessing.maxPvs", mMaxPvs);
+
+    Debug << "BSP max depth: " << mTermMaxDepth << endl;
+	Debug << "BSP min PVS: " << mTermMinPvs << endl;
+	Debug << "BSP min area: " << mTermMinArea << endl;
+	Debug << "BSP min rays: " << mTermMinRays << endl;
+	Debug << "BSP max polygon candidates: " << mMaxPolyCandidates << endl;
+	Debug << "BSP max plane candidates: " << mMaxRayCandidates << endl;
+
+	Debug << "Split plane strategy: ";
+	if (mSplitPlaneStrategy & RANDOM_POLYGON)
+		Debug << "random polygon ";
+	if (mSplitPlaneStrategy & AXIS_ALIGNED)
+		Debug << "axis aligned ";
+	if (mSplitPlaneStrategy & LEAST_SPLITS)     
+		Debug << "least splits ";
+	if (mSplitPlaneStrategy & VERTICAL_AXIS)
+		Debug << "vertical axis ";
+	if (mSplitPlaneStrategy & BLOCKED_RAYS)
+		Debug << "blocked rays ";
+	if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
+		Debug << "least ray splits ";
+	if (mSplitPlaneStrategy & BALANCED_RAYS)
+		Debug << "balanced rays ";
+	if (mSplitPlaneStrategy & PVS)
+		Debug << "pvs";
+
+	Debug << endl;
+}
+
+
+const VspBspTreeStatistics &VspBspTree::GetStatistics() const 
+{
+	return mStat;
+}
+
+void VspBspTreeStatistics::Print(ostream &app) const
+{
+	app << "===== VspBspTree statistics ===============\n";
+
+	app << setprecision(4);
+
+	app << "#N_CTIME  ( Construction time [s] )\n" << Time() << " \n";
+
+	app << "#N_NODES ( Number of nodes )\n" << nodes << "\n";
+
+	app << "#N_INTERIORS ( Number of interior nodes )\n" << Interior() << "\n";
+
+	app << "#N_LEAVES ( Number of leaves )\n" << Leaves() << "\n";
+
+	app << "#N_SPLITS ( Number of splits )\n" << splits << "\n";
+
+	app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maximum depth )\n" 
+		<<	maxDepthNodes * 100 / (double)Leaves() << endl;
+
+	app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maximum depth )\n" 
+		<<	maxDepthNodes * 100 / (double)Leaves() << endl;
+
+	app << "#N_PMINPVSLEAVES  ( Percentage of leaves with mininimal PVS )\n" 
+		<< minPvsNodes * 100 / (double)Leaves() << endl;
+
+	app << "#N_PMINRAYSLEAVES  ( Percentage of leaves with minimal number of rays)\n" 
+		<<	minRaysNodes * 100 / (double)Leaves() << endl;
+
+	app << "#N_PMINAREALEAVES  ( Percentage of leaves with mininum area )\n"
+		<< minAreaNodes * 100 / (double)Leaves() << endl;
+
+	app << "#N_PMAXRAYCONTRIBLEAVES  ( Percentage of leaves with maximal ray contribution )\n" 
+		<<	maxRayContribNodes * 100 / (double)Leaves() << endl;
+
+	app << "#N_PMAXDEPTH ( Maximal reached depth )\n" << maxDepth << endl;
+
+	app << "#N_PMINDEPTH ( Minimal reached depth )\n" << minDepth << endl;
+
+	app << "#AVGDEPTH ( average depth )\n" << AvgDepth() << endl;
+
+	app << "#N_INPUT_POLYGONS (number of input polygons )\n" <<	polys << endl;
+
+	//app << "#N_PVS: " << pvs << endl;
+
+	app << "#N_ROUTPUT_INPUT_POLYGONS ( ratio polygons after subdivision / input polygons )\n" <<
+		 (polys + splits) / (double)polys << endl;
+	
+	app << "===== END OF VspBspTree statistics ==========\n";
+}
+
+
+VspBspTree::~VspBspTree()
+{
+	DEL_PTR(mRoot);
+}
+
+int VspBspTree::AddMeshToPolygons(Mesh *mesh, PolygonContainer &polys, MeshInstance *parent)
+{
+	FaceContainer::const_iterator fi;
+	
+	// copy the face data to polygons
+	for (fi = mesh->mFaces.begin(); fi != mesh->mFaces.end(); ++ fi)
+	{
+		Polygon3 *poly = new Polygon3((*fi), mesh);
+		
+		if (poly->Valid())
+		{
+			poly->mParent = parent; // set parent intersectable
+			polys.push_back(poly);
+		}
+		else
+			DEL_PTR(poly);
+	}
+	return (int)mesh->mFaces.size();
+}
+
+int VspBspTree::AddToPolygonSoup(const ViewCellContainer &viewCells, 
+							  PolygonContainer &polys, 
+							  int maxObjects)
+{
+	int limit = (maxObjects > 0) ? 
+		Min((int)viewCells.size(), maxObjects) : (int)viewCells.size();
+  
+	int polysSize = 0;
+
+	for (int i = 0; i < limit; ++ i)
+	{
+		if (viewCells[i]->GetMesh()) // copy the mesh data to polygons
+		{
+			mBox.Include(viewCells[i]->GetBox()); // add to BSP tree aabb
+			polysSize += AddMeshToPolygons(viewCells[i]->GetMesh(), polys, viewCells[i]);
+		}
+	}
+
+	return polysSize;
+}
+
+int VspBspTree::AddToPolygonSoup(const ObjectContainer &objects, PolygonContainer &polys, int maxObjects)
+{
+	int limit = (maxObjects > 0) ? Min((int)objects.size(), maxObjects) : (int)objects.size();
+  
+	for (int i = 0; i < limit; ++i)
+	{
+		Intersectable *object = objects[i];//*it;
+		Mesh *mesh = NULL;
+
+		switch (object->Type()) // extract the meshes
+		{
+		case Intersectable::MESH_INSTANCE:
+			mesh = dynamic_cast<MeshInstance *>(object)->GetMesh();
+			break;
+		case Intersectable::VIEW_CELL:
+			mesh = dynamic_cast<ViewCell *>(object)->GetMesh();
+			break;
+			// TODO: handle transformed mesh instances
+		default:
+			Debug << "intersectable type not supported" << endl;
+			break;
+		}
+		
+        if (mesh) // copy the mesh data to polygons
+		{
+			mBox.Include(object->GetBox()); // add to BSP tree aabb
+			AddMeshToPolygons(mesh, polys, mRootCell);
+		}
+	}
+
+	return (int)polys.size();
+}
+
+void VspBspTree::Construct(const VssRayContainer &sampleRays)
+{
+    mStat.nodes = 1;
+	mBox.Initialize(); 	// initialise BSP tree bounding box
+	
+	PolygonContainer polys;
+	RayInfoContainer *rays = new RayInfoContainer();
+
+	VssRayContainer::const_iterator rit, rit_end = sampleRays.end();
+
+	long startTime = GetTime();
+
+	Debug << "**** Extracting polygons from rays ****\n";
+
+	//-- extract polygons intersected by the rays
+	for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
+	{
+		VssRay *ray = *rit;
+	
+		if (ray->mTerminationObject)
+		{
+			MeshInstance *obj = dynamic_cast<MeshInstance *>(ray->mTerminationObject);
+			AddMeshToPolygons(obj->GetMesh(), polys, obj);
+		}
+
+		if (ray->mOriginObject)
+		{
+			MeshInstance *obj = dynamic_cast<MeshInstance *>(ray->mOriginObject);
+			AddMeshToPolygons(obj->GetMesh(), polys, obj);
+		}
+	}
+
+	//-- store rays
+	for (rit = sampleRays.begin(); rit != rit_end; ++ rit)
+	{
+		rays->push_back(RayInfo(*rit));
+	}
+
+	mStat.polys = (int)polys.size();
+
+	Debug << "**** Finished polygon extraction ****" << endl;
+	Debug << (int)polys.size() << " polys extracted from " << (int)sampleRays.size() << " rays" << endl;
+	Debug << "extraction time: " << TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
+
+	Construct(new PolygonContainer(polys), rays);
+
+	// clean up polygons because they are not deleted in the tree
+	CLEAR_CONTAINER(polys);
+}
+
+void VspBspTree::Construct(PolygonContainer *polys, RayInfoContainer *rays)
+{
+	std::stack<VspBspTraversalData> tStack;
+
+	mRoot = new VspBspLeaf();
+
+	VspBspNodeGeometry *cell = new VspBspNodeGeometry();
+	ConstructGeometry(mRoot, *cell);
+
+	VspBspTraversalData tData(mRoot, polys, 0, mRootCell, rays, 
+						   ComputePvsSize(*rays), cell->GetArea(), cell);
+
+	tStack.push(tData);
+
+	mStat.Start();
+	cout << "Contructing bsp tree ... ";
+	long startTime = GetTime();
+	while (!tStack.empty()) 
+	{
+		tData = tStack.top();
+
+	    tStack.pop();
+
+		// subdivide leaf node
+		VspBspNode *r = Subdivide(tStack, tData);
+
+		if (r == mRoot)
+			Debug << "BSP tree construction time spent at root: " << TimeDiff(startTime, GetTime())*1e-3 << "s" << endl;
+	}
+
+	cout << "finished\n";
+
+	mStat.Stop();
+}
+
+bool VspBspTree::TerminationCriteriaMet(const VspBspTraversalData &data) const
+{
+	return 
+		(((int)data.mRays->size() <= mTermMinRays) ||
+		 (data.mPvs <= mTermMinPvs) ||
+		 (data.mArea <= mTermMinArea) ||
+		 (data.mDepth >= mTermMaxDepth) ||
+		 (data.GetAvgRayContribution() < mTermMaxRayContribution));
+}
+
+VspBspNode *VspBspTree::Subdivide(VspBspTraversalStack &tStack, VspBspTraversalData &tData)
+{
+	//-- terminate traversal  
+	if (TerminationCriteriaMet(tData))		
+	{
+		VspBspLeaf *leaf = dynamic_cast<VspBspLeaf *>(tData.mNode);
+	
+		BspViewCell *viewCell =
+			dynamic_cast<BspViewCell *>(ViewCell::Generate());
+		
+		leaf->SetViewCell(viewCell);
+		//TODO viewCell->mVspBspLeaves.push_back(leaf);
+
+		//-- add pvs
+		if (viewCell != mRootCell)
+		{
+			int conSamp = 0, sampCon = 0;
+			leaf->AddToPvs(*tData.mRays, conSamp, sampCon);
+			
+			mStat.contributingSamples += conSamp;
+			mStat.sampleContributions += sampCon;
+		}
+
+		EvaluateLeafStats(tData);
+		
+		//-- clean up
+
+		DEL_PTR(tData.mPolygons);
+		DEL_PTR(tData.mRays);
+		DEL_PTR(tData.mGeometry);
+		
+		return leaf;
+	}
+
+	//-- continue subdivision
+	PolygonContainer coincident;
+	
+	VspBspTraversalData tFrontData(NULL, new PolygonContainer(), tData.mDepth + 1, mRootCell, 
+								new RayInfoContainer(), 0, 0, new VspBspNodeGeometry());
+	VspBspTraversalData tBackData(NULL, new PolygonContainer(), tData.mDepth + 1, mRootCell, 
+							   new RayInfoContainer(), 0, 0, new VspBspNodeGeometry());
+
+	// create new interior node and two leaf nodes
+	VspBspInterior *interior = 
+		SubdivideNode(tData, tFrontData, tBackData, coincident);
+
+#ifdef _DEBUG	
+	if (frontPolys->empty() && backPolys->empty() && (coincident.size() > 2))
+	{	for (PolygonContainer::iterator it = coincident.begin(); it != coincident.end(); ++it)
+			Debug << (*it) << " " << (*it)->GetArea() << " " << (*it)->mParent << endl ;
+		Debug << endl;}
+#endif
+
+	// don't need coincident polygons anymory
+	CLEAR_CONTAINER(coincident);
+
+	// push the children on the stack
+	tStack.push(tFrontData);
+	tStack.push(tBackData);
+
+	// cleanup
+	DEL_PTR(tData.mNode);
+	DEL_PTR(tData.mPolygons);
+	DEL_PTR(tData.mRays);
+	DEL_PTR(tData.mGeometry);		
+
+	return interior;
+}
+
+VspBspInterior *VspBspTree::SubdivideNode(VspBspTraversalData &tData,
+									VspBspTraversalData &frontData,
+									VspBspTraversalData &backData,
+									PolygonContainer &coincident)
+{
+	mStat.nodes += 2;
+	
+	VspBspLeaf *leaf = dynamic_cast<VspBspLeaf *>(tData.mNode);
+	// select subdivision plane
+	VspBspInterior *interior = 
+		new VspBspInterior(SelectPlane(leaf, tData)); 
+
+#ifdef _DEBUG
+	Debug << interior << endl;
+#endif
+	
+	// subdivide rays into front and back rays
+	SplitRays(interior->mPlane, *tData.mRays, *frontData.mRays, *backData.mRays);
+	
+	// subdivide polygons with plane
+	mStat.splits += interior->SplitPolygons(*tData.mPolygons, 
+		                                    *frontData.mPolygons, 
+											*backData.mPolygons, 
+											coincident);
+
+    // compute pvs
+	frontData.mPvs = ComputePvsSize(*frontData.mRays);
+	backData.mPvs = ComputePvsSize(*backData.mRays);
+
+	// split geometry and compute area
+	if (1)
+	{
+		tData.mGeometry->SplitGeometry(*frontData.mGeometry, 
+									   *backData.mGeometry, 
+									   *this, 
+									   interior->mPlane);
+	
+		
+		frontData.mArea = frontData.mGeometry->GetArea();
+		backData.mArea = backData.mGeometry->GetArea();
+	}
+
+	// compute accumulated ray length
+	//frontData.mAccRayLength = AccumulatedRayLength(*frontData.mRays);
+	//backData.mAccRayLength = AccumulatedRayLength(*backData.mRays);
+
+	//-- create front and back leaf
+
+	VspBspInterior *parent = leaf->GetParent();
+
+	// replace a link from node's parent
+	if (!leaf->IsRoot())
+	{
+		parent->ReplaceChildLink(leaf, interior);
+		interior->SetParent(parent);
+	}
+	else // new root
+	{
+		mRoot = interior;
+	}
+
+	// and setup child links
+	interior->SetupChildLinks(new VspBspLeaf(interior), new VspBspLeaf(interior));
+	
+	frontData.mNode = interior->mFront;
+	backData.mNode = interior->mBack;
+	
+	//DEL_PTR(leaf);
+	return interior;
+}
+
+void VspBspTree::SortSplitCandidates(const PolygonContainer &polys, 
+								  const int axis, 
+								  vector<SortableEntry> &splitCandidates) const
+{
+  splitCandidates.clear();
+  
+  int requestedSize = 2 * (int)polys.size();
+  // creates a sorted split candidates array  
+  splitCandidates.reserve(requestedSize);
+  
+  PolygonContainer::const_iterator it, it_end = polys.end();
+
+  AxisAlignedBox3 box;
+
+  // insert all queries 
+  for(it = polys.begin(); it != it_end; ++ it)
+  {
+	  box.Initialize();
+	  box.Include(*(*it));
+	  
+	  splitCandidates.push_back(SortableEntry(SortableEntry::POLY_MIN, box.Min(axis), *it));
+      splitCandidates.push_back(SortableEntry(SortableEntry::POLY_MAX, box.Max(axis), *it));
+  }
+  
+  stable_sort(splitCandidates.begin(), splitCandidates.end());
+}
+
+
+float VspBspTree::BestCostRatio(const PolygonContainer &polys,
+							 const AxisAlignedBox3 &box,
+							 const int axis,
+							 float &position,
+							 int &objectsBack,
+							 int &objectsFront) const
+{
+	vector<SortableEntry> splitCandidates;
+
+	SortSplitCandidates(polys, axis, splitCandidates);
+	
+	// go through the lists, count the number of objects left and right
+	// and evaluate the following cost funcion:
+	// C = ct_div_ci  + (ol + or)/queries
+	
+	int objectsLeft = 0, objectsRight = (int)polys.size();
+	
+	float minBox = box.Min(axis);
+	float maxBox = box.Max(axis);
+	float boxArea = box.SurfaceArea();
+  
+	float minBand = minBox + mSplitBorder * (maxBox - minBox);
+	float maxBand = minBox + (1.0f - mSplitBorder) * (maxBox - minBox);
+	
+	float minSum = 1e20f;
+	vector<SortableEntry>::const_iterator ci, ci_end = splitCandidates.end();
+
+	for(ci = splitCandidates.begin(); ci != ci_end; ++ ci) 
+	{
+		switch ((*ci).type) 
+		{
+			case SortableEntry::POLY_MIN:
+				++ objectsLeft;
+				break;
+			case SortableEntry::POLY_MAX:
+			    -- objectsRight;
+				break;
+			default:
+				break;
+		}
+		
+		if ((*ci).value > minBand && (*ci).value < maxBand) 
+		{
+			AxisAlignedBox3 lbox = box;
+			AxisAlignedBox3 rbox = box;
+			lbox.SetMax(axis, (*ci).value);
+			rbox.SetMin(axis, (*ci).value);
+
+			float sum = objectsLeft * lbox.SurfaceArea() + 
+						objectsRight * rbox.SurfaceArea();
+      
+			if (sum < minSum) 
+			{
+				minSum = sum;
+				position = (*ci).value;
+
+				objectsBack = objectsLeft;
+				objectsFront = objectsRight;
+			}
+		}
+	}
+  
+	float oldCost = (float)polys.size();
+	float newCost = mAaCtDivCi + minSum / boxArea;
+	float ratio = newCost / oldCost;
+
+
+#if 0
+  Debug << "====================" << endl;
+  Debug << "costRatio=" << ratio << " pos=" << position<<" t=" << (position - minBox)/(maxBox - minBox)
+      << "\t o=(" << objectsBack << "," << objectsFront << ")" << endl;
+#endif
+  return ratio;
+}
+
+bool VspBspTree::SelectAxisAlignedPlane(Plane3 &plane, 
+									 const PolygonContainer &polys) const
+{
+	AxisAlignedBox3 box;
+	box.Initialize();
+	
+	// create bounding box of region
+	Polygon3::IncludeInBox(polys, box);
+        
+	int objectsBack = 0, objectsFront = 0;
+	int axis = 0;
+	float costRatio = MAX_FLOAT;
+	Vector3 position;
+
+	//-- area subdivision
+	for (int i = 0; i < 3; ++ i) 
+	{
+		float p = 0;
+		float r = BestCostRatio(polys, box, i, p, objectsBack, objectsFront);
+		
+		if (r < costRatio)
+		{
+			costRatio = r;
+			axis = i;
+			position = p;
+		}
+	}
+	
+	if (costRatio >= mMaxCostRatio)
+		return false;
+
+	Vector3 norm(0,0,0); norm[axis] = 1.0f;
+	plane = Plane3(norm, position);
+
+	return true;
+}
+
+Plane3 VspBspTree::SelectPlane(VspBspLeaf *leaf, VspBspTraversalData &data)
+{
+	if ((mSplitPlaneStrategy & AXIS_ALIGNED) &&
+		((int)data.mRays->size() > mTermMinRaysForAxisAligned))
+	{
+		Plane3 plane;
+		if (SelectAxisAlignedPlane(plane, *data.mPolygons)) // TODO: should be rays
+			return plane;
+	}
+
+	// simplest strategy: just take next polygon
+	if (mSplitPlaneStrategy & RANDOM_POLYGON)
+	{
+        if (!data.mPolygons->empty())
+		{
+			Polygon3 *nextPoly = (*data.mPolygons)[Random((int)data.mPolygons->size())];
+			return nextPoly->GetSupportingPlane();
+		}
+		else
+		{
+			//-- choose plane on midpoint of a ray
+			const int candidateIdx = Random((int)data.mRays->size());
+								
+			const Vector3 minPt = (*data.mRays)[candidateIdx].ExtrapOrigin();
+			const Vector3 maxPt = (*data.mRays)[candidateIdx].ExtrapTermination();
+
+			const Vector3 pt = (maxPt + minPt) * 0.5;
+
+			const Vector3 normal =  (*data.mRays)[candidateIdx].mRay->GetDir();
+			
+			return Plane3(normal, pt);
+		}
+
+		return Plane3();
+	}
+
+	// use heuristics to find appropriate plane
+	return SelectPlaneHeuristics(leaf, data); 
+}
+
+Plane3 VspBspTree::SelectPlaneHeuristics(VspBspLeaf *leaf, VspBspTraversalData &data)
+{
+	float lowestCost = MAX_FLOAT;
+	Plane3 bestPlane;
+	Plane3 plane;
+
+	int limit = Min((int)data.mPolygons->size(), mMaxPolyCandidates);
+	
+	int candidateIdx = limit;
+
+	for (int i = 0; i < limit; ++ i)
+	{
+		candidateIdx = GetNextCandidateIdx(candidateIdx, *data.mPolygons);
+		
+		Polygon3 *poly = (*data.mPolygons)[candidateIdx];
+		// evaluate current candidate
+		const float candidateCost = 
+			SplitPlaneCost(poly->GetSupportingPlane(), data);
+
+		if (candidateCost < lowestCost)
+		{
+			bestPlane = poly->GetSupportingPlane();
+			lowestCost = candidateCost;
+		}
+	}
+	
+	//Debug << "lowest: " << lowestCost << endl;
+
+	//-- choose candidate planes extracted from rays
+	// we currently use different two methods chosen with
+	// equal probability
+	
+	RayInfoContainer *rays = data.mRays;
+	// take 3 ray endpoints, where two are minimum and one a maximum
+	// point or the other way round
+	for (int i = 0; i < mMaxRayCandidates / 2; ++ i)
+	{
+		candidateIdx = Random((int)rays->size());
+	
+		RayInfo rayInf = (*rays)[candidateIdx];
+		const Vector3 minPt = rayInf.ExtrapOrigin();
+		const Vector3 maxPt = rayInf.ExtrapTermination();
+
+		const Vector3 pt = (maxPt + minPt) * 0.5;
+		const Vector3 normal = Normalize(rayInf.mRay->GetDir());
+
+		plane = Plane3(normal, pt);
+	
+		const float candidateCost = SplitPlaneCost(plane, data);
+
+		if (candidateCost < lowestCost)
+		{
+			bestPlane = plane;	
+			lowestCost = candidateCost;
+		}
+	}
+
+	// take plane normal as plane normal and the midpoint of the ray.
+	// PROBLEM: does not resemble any point where visibility is likely to change
+	//Debug << "lowest: " << lowestCost << endl;
+	for (int i = 0; i < mMaxRayCandidates / 2; ++ i)
+	{
+		Vector3 pt[3];
+		int idx[3];
+		int cmaxT = 0;
+		int cminT = 0;
+		bool chooseMin = false;
+
+		for (int j = 0; j < 3; j ++)
+		{
+			idx[j] = Random((int)rays->size() * 2);
+				
+			if (idx[j] >= (int)rays->size())
+			{
+				idx[j] -= (int)rays->size();
+				
+				chooseMin = (cminT < 2);
+			}
+			else
+				chooseMin = (cmaxT < 2);
+
+			RayInfo rayInf = (*rays)[idx[j]];
+			pt[j] = chooseMin ? rayInf.ExtrapOrigin() : rayInf.ExtrapTermination();
+		}	
+			
+		plane = Plane3(pt[0], pt[1], pt[2]);
+
+		const float candidateCost = SplitPlaneCost(plane, data);
+
+		if (candidateCost < lowestCost)
+		{
+			//Debug << "choose ray plane 2: " << candidateCost << endl;
+			bestPlane = plane;
+			
+			lowestCost = candidateCost;
+		}
+	}	
+
+#ifdef _DEBUG
+	Debug << "plane lowest cost: " << lowestCost << endl;
+#endif
+	return bestPlane;
+}
+
+int VspBspTree::GetNextCandidateIdx(int currentIdx, PolygonContainer &polys)
+{
+	const int candidateIdx = Random(currentIdx --);
+
+	// swap candidates to avoid testing same plane 2 times
+	std::swap(polys[currentIdx], polys[candidateIdx]);
+	
+	return currentIdx;
+	//return Random((int)polys.size());
+}
+
+
+bool VspBspTree::BoundRay(const Ray &ray, float &minT, float &maxT) const
+{
+	maxT = 1e6;
+	minT = 0;
+	
+	// test with tree bounding box
+	if (!mBox.GetMinMaxT(ray, &minT, &maxT))
+		return false;
+
+	if (minT < 0) // start ray from origin
+		minT = 0;
+
+	// bound ray or line segment
+	if ((ray.GetType() == Ray::LOCAL_RAY) && 
+	    !ray.intersections.empty() && 
+		(ray.intersections[0].mT <= maxT))
+	{
+		maxT = ray.intersections[0].mT;
+	}
+
+	return true;
+}
+
+inline void VspBspTree::GenerateUniqueIdsForPvs()
+{
+	Intersectable::NewMail(); sBackId = ViewCell::sMailId;
+	Intersectable::NewMail(); sFrontId = ViewCell::sMailId;
+	Intersectable::NewMail(); sFrontAndBackId = ViewCell::sMailId;
+}
+
+float VspBspTree::SplitPlaneCost(const Plane3 &candidatePlane, 
+								 const VspBspTraversalData &data)
+{
+	float val = 0;
+
+	float sumBalancedRays = 0;
+	float sumRaySplits = 0;
+	
+	int backId = 0;
+	int frontId = 0;
+	int frontAndBackId = 0;
+
+	int frontPvs = 0;
+	int backPvs = 0;
+
+	// probability that view point lies in child
+	float pOverall = 0;
+	float pFront = 0;
+	float pBack = 0;
+
+	if (mSplitPlaneStrategy & PVS)
+	{
+		// create unique ids for pvs heuristics
+		GenerateUniqueIdsForPvs();
+
+		if (mPvsUseArea) // use front and back cell areas to approximate volume
+		{	
+			// construct child geometry with regard to the candidate split plane
+			VspBspNodeGeometry frontCell;
+			VspBspNodeGeometry backCell;
+
+			data.mGeometry->SplitGeometry(frontCell, backCell, *this, candidatePlane);
+		
+			pFront = frontCell.GetArea();
+			pBack = backCell.GetArea();
+
+			pOverall = data.mArea;
+		}
+	}
+			
+	RayInfoContainer::const_iterator rit, rit_end = data.mRays->end();
+
+	for (rit = data.mRays->begin(); rit != data.mRays->end(); ++ rit)
+	{
+		VssRay *ray = (*rit).mRay;
+		const int cf = (*rit).ComputeRayIntersection(candidatePlane, (*rit).mRay->mT);
+
+		if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
+		{
+			sumBalancedRays += cf;
+		}
+		
+		if (mSplitPlaneStrategy & BALANCED_RAYS)
+		{
+			if (cf == 0)
+				++ sumRaySplits;
+		}
+
+		if (mSplitPlaneStrategy & PVS)
+		{
+			// in case the ray intersects an object
+			// assure that we only count the object 
+			// once for the front and once for the back side of the plane
+			
+			// add the termination object
+			AddObjToPvs(ray->mTerminationObject, cf, frontPvs, backPvs);
+			
+			// add the source object
+			AddObjToPvs(ray->mOriginObject, cf, frontPvs, backPvs);
+			
+			if (!mPvsUseArea) // use front and back cell areas to approximate volume
+			{	
+				const float len = (*rit).SegmentLength();
+				pOverall += len;
+
+				// use length of rays to approximate volume
+				if (cf == 1)
+					pFront += len;
+				if (cf == -1)
+					pBack += len;
+				if (cf == 0)
+				{
+					float newLen = len * 
+						((*rit).GetMaxT() - (*rit).mRay->mT) / 
+						((*rit).GetMaxT() - (*rit).GetMinT());
+
+					if (candidatePlane.Side((*rit).ExtrapOrigin()) <= 0)
+					{
+						pBack += newLen;
+						pFront += len - newLen;
+					}
+					else
+					{
+						pFront += newLen;
+						pBack += len - newLen;
+					}
+				}				
+			}
+		}
+	}
+
+	const float raysSize = (float)data.mRays->size() + Limits::Small;
+
+	if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
+		val += mLeastRaySplitsFactor * sumRaySplits / raysSize;
+
+	if (mSplitPlaneStrategy & BALANCED_RAYS)
+		val += mBalancedRaysFactor * fabs(sumBalancedRays) /  raysSize;
+
+	float denom = pOverall * (float)data.mPvs * 2.0f + Limits::Small;
+	if ((mSplitPlaneStrategy & PVS) && data.mArea && data.mPvs)
+	{
+		val += mPvsFactor * (frontPvs * pFront + (backPvs * pBack)) / denom;
+
+		// give penalty to unbalanced split
+		if (1)
+		if (((pFront * 0.2 + Limits::Small) > pBack) || (pFront < (pBack * 0.2 + Limits::Small)))
+			val += 0.5;
+	}
+
+#ifdef _DEBUG
+	Debug << "totalpvs: " << pvs << " ptotal: " << pOverall
+		  << " frontpvs: " << frontPvs << " pFront: " << pFront 
+		  << " backpvs: " << backPvs << " pBack: " << pBack << endl << endl;
+#endif
+	return val;
+}
+
+void VspBspTree::AddObjToPvs(Intersectable *obj,
+						  const int cf,
+						  int &frontPvs,
+						  int &backPvs) const
+{
+	if (!obj)
+		return;
+	// TODO: does this really belong to no pvs?
+	//if (cf == Ray::COINCIDENT) return;
+
+	// object belongs to both PVS
+	if (cf >= 0)
+	{
+		if ((obj->mMailbox != sFrontId) && 
+			(obj->mMailbox != sFrontAndBackId))
+		{
+			++ frontPvs;
+
+			if (obj->mMailbox == sBackId)
+				obj->mMailbox = sFrontAndBackId;	
+			else
+				obj->mMailbox = sFrontId;								
+		}
+	}
+	
+	if (cf <= 0)
+	{
+		if ((obj->mMailbox != sBackId) &&
+			(obj->mMailbox != sFrontAndBackId))
+		{
+			++ backPvs;
+
+			if (obj->mMailbox == sFrontId)
+				obj->mMailbox = sFrontAndBackId; 
+			else
+				obj->mMailbox = sBackId;				
+		}
+	}
+}
+
+void VspBspTree::CollectLeaves(vector<VspBspLeaf *> &leaves) const
+{
+	stack<VspBspNode *> nodeStack;
+	nodeStack.push(mRoot);
+ 
+	while (!nodeStack.empty()) 
+	{
+		VspBspNode *node = nodeStack.top();
+    
+		nodeStack.pop();
+    
+		if (node->IsLeaf()) 
+		{
+			VspBspLeaf *leaf = (VspBspLeaf *)node;		
+			leaves.push_back(leaf);
+		} 
+		else 
+		{
+			VspBspInterior *interior = dynamic_cast<VspBspInterior *>(node);
+
+			nodeStack.push(interior->GetBack());
+			nodeStack.push(interior->GetFront());
+		}
+	}
+}
+
+AxisAlignedBox3 VspBspTree::GetBoundingBox() const
+{
+	return mBox;
+}
+
+VspBspNode *VspBspTree::GetRoot() const
+{
+	return mRoot;
+}
+
+void VspBspTree::EvaluateLeafStats(const VspBspTraversalData &data)
+{
+	// the node became a leaf -> evaluate stats for leafs
+	VspBspLeaf *leaf = dynamic_cast<VspBspLeaf *>(data.mNode);
+
+	// store maximal and minimal depth
+	if (data.mDepth > mStat.maxDepth)
+		mStat.maxDepth = data.mDepth;
+
+	if (data.mDepth < mStat.minDepth)
+		mStat.minDepth = data.mDepth;
+
+	// accumulate depth to compute average depth
+	mStat.accumDepth += data.mDepth;
+	
+
+	if (data.mDepth >= mTermMaxDepth)
+		++ mStat.maxDepthNodes;
+
+	if (data.mPvs < mTermMinPvs)
+		++ mStat.minPvsNodes;
+
+	if ((int)data.mRays->size() < mTermMinRays)
+		++ mStat.minRaysNodes;
+
+	if (data.GetAvgRayContribution() > mTermMaxRayContribution)
+		++ mStat.maxRayContribNodes;
+	
+	if (data.mGeometry->GetArea() <= mTermMinArea) 
+		++ mStat.minAreaNodes;
+
+#ifdef _DEBUG
+	Debug << "BSP stats: "
+		  << "Depth: " << data.mDepth << " (max: " << mTermMaxDepth << "), "
+		  << "PVS: " << data.mPvs << " (min: " << mTermMinPvs << "), "
+		  << "Area: " << data.mArea << " (min: " << mTermMinArea << "), "
+		  << "#polygons: " << (int)data.mPolygons->size() << " (max: " << mTermMinPolys << "), "
+		  << "#rays: " << (int)data.mRays->size() << " (max: " << mTermMinRays << "), " 
+		  << "#pvs: " << leaf->GetViewCell()->GetPvs().GetSize() << "=, "
+		  << "#avg ray contrib (pvs): " << (float)data.mPvs / (float)data.mRays->size() << endl;
+#endif
+}
+
+int VspBspTree::CastRay(Ray &ray)
+{
+	int hits = 0;
+  
+	stack<VspBspRayTraversalData> tStack;
+  
+	float maxt, mint;
+
+	if (!BoundRay(ray, mint, maxt))
+		return 0;
+
+	Intersectable::NewMail();
+
+	Vector3 entp = ray.Extrap(mint);
+	Vector3 extp = ray.Extrap(maxt);
+  
+	VspBspNode *node = mRoot;
+	VspBspNode *farChild = NULL;
+	
+	while (1)
+	{
+		if (!node->IsLeaf()) 
+		{
+			VspBspInterior *in = (VspBspInterior *) node;
+			
+			Plane3 *splitPlane = in->GetPlane();
+
+			int entSide = splitPlane->Side(entp);
+			int extSide = splitPlane->Side(extp);
+
+			Vector3 intersection;
+
+			if (entSide < 0)
+			{
+				node = in->GetBack();
+
+				if(extSide <= 0) // plane does not split ray => no far child
+					continue;
+					
+				farChild = in->GetFront(); // plane splits ray
+
+			} else if (entSide > 0)
+			{
+				node = in->GetFront();
+
+				if (extSide >= 0) // plane does not split ray => no far child
+					continue;
+
+				farChild = in->GetBack(); // plane splits ray			
+			}
+			else // ray and plane are coincident
+			{
+				// WHAT TO DO IN THIS CASE ?
+				//break;
+				node = in->GetFront();
+				continue;
+			}
+
+			// push data for far child
+			tStack.push(VspBspRayTraversalData(farChild, extp, maxt));
+
+			// find intersection of ray segment with plane
+			float t;
+			extp = splitPlane->FindIntersection(ray.GetLoc(), extp, &t);
+			maxt *= t;
+			
+		} else // reached leaf => intersection with view cell
+		{
+			VspBspLeaf *leaf = dynamic_cast<VspBspLeaf *>(node);
+      
+			if (!leaf->mViewCell->Mailed())
+			{
+				//ray.bspIntersections.push_back(Ray::VspBspIntersection(maxt, leaf));
+				leaf->mViewCell->Mail();
+				++ hits;
+			}
+			
+   			//-- fetch the next far child from the stack
+			if (tStack.empty())
+				break;
+      
+			entp = extp;
+			mint = maxt; // NOTE: need this?
+
+			if (ray.GetType() == Ray::LINE_SEGMENT && mint > 1.0f)
+				break;
+
+			VspBspRayTraversalData &s = tStack.top();
+
+			node = s.mNode;
+			extp = s.mExitPoint;
+			maxt = s.mMaxT;
+
+			tStack.pop();
+		}
+	}
+
+	return hits;
+}
+
+bool VspBspTree::Export(const string filename)
+{
+	Exporter *exporter = Exporter::GetExporter(filename);
+
+	if (exporter) 
+	{
+		//exporter->ExportVspBspTree(*this);
+		return true;
+	}	
+
+	return false;
+}
+
+void VspBspTree::CollectViewCells(ViewCellContainer &viewCells) const
+{
+	stack<VspBspNode *> nodeStack;
+	nodeStack.push(mRoot);
+
+	ViewCell::NewMail();
+
+	while (!nodeStack.empty()) 
+	{
+		VspBspNode *node = nodeStack.top();
+		nodeStack.pop();
+
+		if (node->IsLeaf()) 
+		{
+			ViewCell *viewCell = dynamic_cast<VspBspLeaf *>(node)->mViewCell;
+
+			if (!viewCell->Mailed()) 
+			{
+				viewCell->Mail();
+				viewCells.push_back(viewCell);
+			}
+		}
+		else 
+		{
+			VspBspInterior *interior = dynamic_cast<VspBspInterior *>(node);
+
+			nodeStack.push(interior->mFront);
+			nodeStack.push(interior->mBack);
+		}
+	}
+}
+
+void VspBspTree::EvaluateViewCellsStats(VspBspViewCellsStatistics &stat) const
+{
+	stat.Reset();
+
+	stack<VspBspNode *> nodeStack;
+	nodeStack.push(mRoot);
+
+	ViewCell::NewMail();
+
+	// exclude root cell
+	mRootCell->Mail();
+
+	while (!nodeStack.empty()) 
+	{
+		VspBspNode *node = nodeStack.top();
+		nodeStack.pop();
+
+		if (node->IsLeaf()) 
+		{
+			++ stat.bspLeaves;
+
+			BspViewCell *viewCell = dynamic_cast<VspBspLeaf *>(node)->mViewCell;
+
+			if (!viewCell->Mailed()) 
+			{
+				viewCell->Mail();
+				
+				++ stat.viewCells;
+				const int pvsSize = viewCell->GetPvs().GetSize();
+
+                stat.pvs += pvsSize;
+
+				if (pvsSize < 1)
+					++ stat.emptyPvs;
+
+				if (pvsSize > stat.maxPvs)
+					stat.maxPvs = pvsSize;
+
+				if (pvsSize < stat.minPvs)
+					stat.minPvs = pvsSize;
+
+				//TODO
+				//if ((int)viewCell->mVspBspLeaves.size() > stat.maxVspBspLeaves)
+				//	stat.maxVspBspLeaves = (int)viewCell->mVspBspLeaves.size();		
+			}
+		}
+		else 
+		{
+			VspBspInterior *interior = dynamic_cast<VspBspInterior *>(node);
+
+			nodeStack.push(interior->mFront);
+			nodeStack.push(interior->mBack);
+		}
+	}
+}
+
+bool VspBspTree::MergeViewCells(VspBspLeaf *front, VspBspLeaf *back) const
+{
+	BspViewCell *viewCell = 
+		dynamic_cast<BspViewCell *>(ViewCell::Merge(*front->mViewCell, *back->mViewCell));
+	
+	if (!viewCell)
+		return false;
+
+	// change view cells of all leaves associated with the previous view cells
+	BspViewCell *fVc = front->mViewCell;
+	BspViewCell *bVc = back->mViewCell;
+
+	//TODO 
+	/*
+	vector<VspBspLeaf *> fLeaves = fVc->mVspBspLeaves;
+	vector<VspBspLeaf *> bLeaves = bVc->mVspBspLeaves;
+
+	vector<VspBspLeaf *>::const_iterator it;
+	
+	for (it = fLeaves.begin(); it != fLeaves.end(); ++ it)
+	{
+		(*it)->SetViewCell(viewCell);
+		viewCell->mVspBspLeaves.push_back(*it);
+	}
+	for (it = bLeaves.begin(); it != bLeaves.end(); ++ it)
+	{
+		(*it)->SetViewCell(viewCell);
+		viewCell->mVspBspLeaves.push_back(*it);
+	}
+	
+	DEL_PTR(fVc);
+	DEL_PTR(bVc);
+*/
+	return true;
+}
+
+bool VspBspTree::ShouldMerge(VspBspLeaf *front, VspBspLeaf *back) const
+{
+	ViewCell *fvc = front->mViewCell;
+	ViewCell *bvc = back->mViewCell;
+
+	if ((fvc == mRootCell) || (bvc == mRootCell) || (fvc == bvc))
+		return false;
+
+	int fdiff = fvc->GetPvs().Diff(bvc->GetPvs());
+
+	if (fvc->GetPvs().GetSize() + fdiff < mMaxPvs)
+	{
+		if ((fvc->GetPvs().GetSize() < mMinPvs) ||	
+			(bvc->GetPvs().GetSize() < mMinPvs) ||
+			((fdiff < mMinPvsDif) && (bvc->GetPvs().Diff(fvc->GetPvs()) < mMinPvsDif)))
+		{
+			return true;
+		}
+	}
+	
+	return false;
+}
+
+
+VspBspTreeStatistics &VspBspTree::GetStat()
+{
+	return mStat;
+}
+
+float VspBspTree::AccumulatedRayLength(const RayInfoContainer &rays) const
+{
+	float len = 0;
+
+	RayInfoContainer::const_iterator it, it_end = rays.end();
+
+	for (it = rays.begin(); it != it_end; ++ it)
+		len += (*it).SegmentLength();
+
+	return len;
+}
+
+int VspBspTree::SplitRays(const Plane3 &plane,
+					   RayInfoContainer &rays, 
+					   RayInfoContainer &frontRays, 
+					   RayInfoContainer &backRays)
+{
+	int splits = 0;
+
+	while (!rays.empty())
+	{
+		RayInfo bRay = rays.back();
+		VssRay *ray = bRay.mRay;
+
+		rays.pop_back();
+	
+		// get classification and receive new t
+		const int cf = 	bRay.ComputeRayIntersection(plane, ray->mT);
+
+		switch(cf)
+		{
+		case -1:
+			backRays.push_back(bRay);
+			break;
+		case 1:
+			break;
+			frontRays.push_back(bRay);
+		case 0: 
+			//-- split ray
+
+			// if start point behind plane
+			if (plane.Side(bRay.ExtrapOrigin()) <= 0)
+			{
+				backRays.push_back(RayInfo(ray, bRay.GetMinT(), ray->mT));
+				frontRays.push_back(RayInfo(ray, ray->mT, bRay.GetMaxT()));
+			}
+			else
+			{
+				frontRays.push_back(RayInfo(ray, bRay.GetMinT(), ray->mT));
+				backRays.push_back(RayInfo(ray, ray->mT, bRay.GetMaxT()));
+			}
+			break;
+		default:
+			Debug << "Should not come here 4" << endl;
+			break;
+		}
+	}
+
+	return splits;
+}
+
+void VspBspTree::ExtractHalfSpaces(VspBspNode *n, vector<Plane3> &halfSpaces) const
+{
+	VspBspNode *lastNode;
+	do
+	{
+		lastNode = n;
+
+		// want to get planes defining geometry of this node => don't take
+		// split plane of node itself
+		n = n->GetParent();
+		
+		if (n)
+		{
+			VspBspInterior *interior = dynamic_cast<VspBspInterior *>(n);
+			Plane3 halfSpace = *dynamic_cast<VspBspInterior *>(interior)->GetPlane();
+
+            if (interior->mFront != lastNode)
+				halfSpace.ReverseOrientation();
+
+			halfSpaces.push_back(halfSpace);
+		}
+	}
+	while (n);
+}
+
+void VspBspTree::ConstructGeometry(VspBspNode *n, VspBspNodeGeometry &cell) const
+{
+	PolygonContainer polys;
+	ConstructGeometry(n, polys);
+	cell.mPolys = polys;
+}
+
+void VspBspTree::ConstructGeometry(VspBspNode *n, PolygonContainer &cell) const
+{
+	vector<Plane3> halfSpaces;
+	ExtractHalfSpaces(n, halfSpaces);
+
+	PolygonContainer candidatePolys;
+
+	// bounded planes are added to the polygons (reverse polygons
+	// as they have to be outfacing
+	for (int i = 0; i < (int)halfSpaces.size(); ++ i)
+	{
+		Polygon3 *p = GetBoundingBox().CrossSection(halfSpaces[i]);
+		
+		if (p->Valid())
+		{
+			candidatePolys.push_back(p->CreateReversePolygon());
+			DEL_PTR(p);
+		}
+	}
+
+	// add faces of bounding box (also could be faces of the cell)
+	for (int i = 0; i < 6; ++ i)
+	{
+		VertexContainer vertices;
+	
+		for (int j = 0; j < 4; ++ j)
+			vertices.push_back(mBox.GetFace(i).mVertices[j]);
+
+		candidatePolys.push_back(new Polygon3(vertices));
+	}
+
+	for (int i = 0; i < (int)candidatePolys.size(); ++ i)
+	{
+		// polygon is split by all other planes
+		for (int j = 0; (j < (int)halfSpaces.size()) && candidatePolys[i]; ++ j)
+		{
+			if (i == j) // polygon and plane are coincident
+				continue;
+
+			VertexContainer splitPts;
+			Polygon3 *frontPoly, *backPoly;
+
+			const int cf = candidatePolys[i]->ClassifyPlane(halfSpaces[j]);
+			
+			switch (cf)
+			{
+				case Polygon3::SPLIT:
+					frontPoly = new Polygon3();
+					backPoly = new Polygon3();
+
+					candidatePolys[i]->Split(halfSpaces[j], 
+											 *frontPoly, 
+											 *backPoly);
+
+					DEL_PTR(candidatePolys[i]);
+
+					if (frontPoly->Valid())
+						candidatePolys[i] = frontPoly;
+					else
+						DEL_PTR(frontPoly);
+
+					DEL_PTR(backPoly);
+					break;
+				case Polygon3::BACK_SIDE:
+					DEL_PTR(candidatePolys[i]);
+					break;
+				// just take polygon as it is
+				case Polygon3::FRONT_SIDE:
+				case Polygon3::COINCIDENT:
+				default:
+					break;
+			}
+		}
+		
+		if (candidatePolys[i])
+			cell.push_back(candidatePolys[i]);
+	}
+}
+
+
+int VspBspTree::FindNeighbors(VspBspNode *n, vector<VspBspLeaf *> &neighbors, 
+						   const bool onlyUnmailed) const
+{
+	PolygonContainer cell;
+
+	ConstructGeometry(n, cell);
+
+	stack<VspBspNode *> nodeStack;
+	nodeStack.push(mRoot);
+		
+	// planes needed to verify that we found neighbor leaf.
+	vector<Plane3> halfSpaces;
+	ExtractHalfSpaces(n, halfSpaces);
+
+	while (!nodeStack.empty()) 
+	{
+		VspBspNode *node = nodeStack.top();
+		nodeStack.pop();
+
+		if (node->IsLeaf())
+		{
+            if (node != n && (!onlyUnmailed || !node->Mailed()))
+			{
+				// test all planes of current node if candidate really
+				// is neighbour
+				PolygonContainer neighborCandidate;
+				ConstructGeometry(node, neighborCandidate);
+				
+				bool isAdjacent = true;
+				for (int i = 0; (i < halfSpaces.size()) && isAdjacent; ++ i)
+				{
+					const int cf = 
+						Polygon3::ClassifyPlane(neighborCandidate, halfSpaces[i]);
+
+					if (cf == Polygon3::BACK_SIDE)
+						isAdjacent = false;
+				}
+
+				if (isAdjacent)
+					neighbors.push_back(dynamic_cast<VspBspLeaf *>(node));
+
+				CLEAR_CONTAINER(neighborCandidate);
+			}
+		}
+		else 
+		{
+			VspBspInterior *interior = dynamic_cast<VspBspInterior *>(node);
+	
+			const int cf = Polygon3::ClassifyPlane(cell, interior->mPlane);
+
+			if (cf == Polygon3::FRONT_SIDE)
+				nodeStack.push(interior->mFront);
+			else
+				if (cf == Polygon3::BACK_SIDE)
+					nodeStack.push(interior->mBack);
+				else 
+				{
+					// random decision
+					nodeStack.push(interior->mBack);
+					nodeStack.push(interior->mFront);
+				}
+		}
+	}
+	
+	CLEAR_CONTAINER(cell);
+	return (int)neighbors.size();
+}
+
+VspBspLeaf *VspBspTree::GetRandomLeaf(const Plane3 &halfspace)
+{
+    stack<VspBspNode *> nodeStack;
+  	nodeStack.push(mRoot);
+	
+	int mask = rand();
+  
+	while (!nodeStack.empty()) 
+	{
+		VspBspNode *node = nodeStack.top();
+		nodeStack.pop();
+	  
+		if (node->IsLeaf()) 
+		{
+			return dynamic_cast<VspBspLeaf *>(node);
+		} 
+		else 
+		{
+			VspBspInterior *interior = dynamic_cast<VspBspInterior *>(node);
+			
+			VspBspNode *next;
+	
+			PolygonContainer cell;
+
+			// todo: not very efficient: constructs full cell everytime
+			ConstructGeometry(interior, cell);
+
+			const int cf = Polygon3::ClassifyPlane(cell, halfspace);
+
+			if (cf == Polygon3::BACK_SIDE)
+				next = interior->mFront;
+			else
+				if (cf == Polygon3::FRONT_SIDE)
+					next = interior->mFront;
+			else 
+			{
+				// random decision
+				if (mask & 1)
+					next = interior->mBack;
+				else
+					next = interior->mFront;
+				mask = mask >> 1;
+			}
+
+			nodeStack.push(next);
+		}
+	}
+	
+	return NULL;
+}
+
+VspBspLeaf *VspBspTree::GetRandomLeaf(const bool onlyUnmailed)
+{
+	stack<VspBspNode *> nodeStack;
+	
+	nodeStack.push(mRoot);
+
+	int mask = rand();
+	
+	while (!nodeStack.empty()) 
+	{
+		VspBspNode *node = nodeStack.top();
+		nodeStack.pop();
+		
+		if (node->IsLeaf()) 
+		{
+			if ( (!onlyUnmailed || !node->Mailed()) )
+				return dynamic_cast<VspBspLeaf *>(node);
+		}
+		else 
+		{
+			VspBspInterior *interior = dynamic_cast<VspBspInterior *>(node);
+
+			// random decision
+			if (mask & 1)
+				nodeStack.push(interior->mBack);
+			else
+				nodeStack.push(interior->mFront);
+
+			mask = mask >> 1;
+		}
+	}
+	
+	return NULL;
+}
+
+int VspBspTree::ComputePvsSize(const RayInfoContainer &rays) const
+{
+	int pvsSize = 0;
+
+	RayInfoContainer::const_iterator rit, rit_end = rays.end();
+
+	Intersectable::NewMail();
+
+	for (rit = rays.begin(); rit != rays.end(); ++ rit)
+	{
+		VssRay *ray = (*rit).mRay;
+		
+		if (ray->mOriginObject)
+		{
+			if (!ray->mOriginObject->Mailed())
+			{
+				ray->mOriginObject->Mail();
+				++ pvsSize;
+			}
+		}
+		if (ray->mTerminationObject)
+		{
+			if (!ray->mTerminationObject->Mailed())
+			{
+				ray->mTerminationObject->Mail();
+				++ pvsSize;
+			}
+		}
+	}
+
+	return pvsSize;
+}
+
+/*************************************************************
+ *            VspBspNodeGeometry Implementation                 *
+ *************************************************************/
+
+VspBspNodeGeometry::~VspBspNodeGeometry()
+{
+	CLEAR_CONTAINER(mPolys);
+}
+
+float VspBspNodeGeometry::GetArea() const 
+{ 
+	return Polygon3::GetArea(mPolys);
+}
+
+void VspBspNodeGeometry::SplitGeometry(VspBspNodeGeometry &front,
+									VspBspNodeGeometry &back,
+									const VspBspTree &tree,						 
+									const Plane3 &splitPlane) const
+{	
+	// get cross section of new polygon
+	Polygon3 *planePoly = tree.GetBoundingBox().CrossSection(splitPlane);
+
+	planePoly = SplitPolygon(planePoly, tree);
+
+	//-- plane poly splits all other cell polygons
+	for (int i = 0; i < (int)mPolys.size(); ++ i)
+	{
+		const int cf = mPolys[i]->ClassifyPlane(splitPlane, 0.00001f);
+			
+		// split new polygon with all previous planes
+		switch (cf)
+		{
+			case Polygon3::SPLIT:
+				{
+					Polygon3 *poly = new Polygon3(mPolys[i]->mVertices);
+
+					Polygon3 *frontPoly = new Polygon3();
+					Polygon3 *backPoly = new Polygon3();
+				
+					poly->Split(splitPlane, *frontPoly, *backPoly);
+
+					DEL_PTR(poly);
+
+					if (frontPoly->Valid())
+						front.mPolys.push_back(frontPoly);
+					if (backPoly->Valid())
+						back.mPolys.push_back(backPoly);
+				}
+				
+				break;
+			case Polygon3::BACK_SIDE:
+				back.mPolys.push_back(new Polygon3(mPolys[i]->mVertices));			
+				break;
+			case Polygon3::FRONT_SIDE:
+				front.mPolys.push_back(new Polygon3(mPolys[i]->mVertices));	
+				break;
+			case Polygon3::COINCIDENT:
+				//front.mPolys.push_back(CreateReversePolygon(mPolys[i]));
+				back.mPolys.push_back(new Polygon3(mPolys[i]->mVertices));
+				break;
+			default:
+				break;
+		}
+	}
+
+	//-- finally add the new polygon to the child cells
+	if (planePoly)
+	{
+		// add polygon with normal pointing into positive half space to back cell
+		back.mPolys.push_back(planePoly);
+		// add polygon with reverse orientation to front cell
+		front.mPolys.push_back(planePoly->CreateReversePolygon());
+	}
+
+	//Debug << "returning new geometry " << mPolys.size() << " f: " << front.mPolys.size() << " b: " << back.mPolys.size() << endl;
+	//Debug << "old area " << GetArea() << " f: " << front.GetArea() << " b: " << back.GetArea() << endl;
+}
+
+Polygon3 *VspBspNodeGeometry::SplitPolygon(Polygon3 *planePoly,
+										const VspBspTree &tree) const
+{
+	// polygon is split by all other planes
+	for (int i = 0; (i < (int)mPolys.size()) && planePoly; ++ i)
+	{
+		Plane3 plane = mPolys[i]->GetSupportingPlane();
+
+		const int cf = 
+			planePoly->ClassifyPlane(plane, 0.00001f);
+			
+		// split new polygon with all previous planes
+		switch (cf)
+		{
+			case Polygon3::SPLIT:
+				{
+					Polygon3 *frontPoly = new Polygon3();
+					Polygon3 *backPoly = new Polygon3();
+
+					planePoly->Split(plane, *frontPoly, *backPoly);
+					DEL_PTR(planePoly);
+
+					if (backPoly->Valid())
+						planePoly = backPoly;
+					else
+						DEL_PTR(backPoly);
+				}
+				break;
+			case Polygon3::FRONT_SIDE:
+				DEL_PTR(planePoly);
+                break;
+			// polygon is taken as it is
+			case Polygon3::BACK_SIDE:
+			case Polygon3::COINCIDENT:
+			default:
+				break;
+		}
+	}
+
+	return planePoly;
+}
+
+void VspBspViewCellsStatistics::Print(ostream &app) const
+{
+	app << "===== VspBspViewCells statistics ===============\n";
+
+	app << setprecision(4);
+
+	//app << "#N_CTIME  ( Construction time [s] )\n" << Time() << " \n";
+
+	app << "#N_OVERALLPVS ( objects in PVS )\n" << pvs << endl;
+
+	app << "#N_PMAXPVS ( largest PVS )\n" << maxPvs << endl;
+
+	app << "#N_PMINPVS ( smallest PVS )\n" << minPvs << endl;
+
+	app << "#N_PAVGPVS ( average PVS )\n" << AvgPvs() << endl;
+
+	app << "#N_PEMPTYPVS ( view cells with PVS smaller 2 )\n" << emptyPvs << endl;
+
+	app << "#N_VIEWCELLS ( number of view cells)\n" << viewCells << endl;
+
+	app << "#N_AVGBSPLEAVES (average number of BSP leaves per view cell )\n" << AvgVspBspLeaves() << endl;
+
+	app << "#N_MAXBSPLEAVES ( maximal number of BSP leaves per view cell )\n" << maxVspBspLeaves << endl;
+	
+	app << "===== END OF VspBspViewCells statistics ==========\n";
+}
Index: trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.cpp
===================================================================
--- trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.cpp	(revision 436)
+++ trunk/VUT/GtpVisibilityPreprocessor/src/ViewCellBsp.cpp	(revision 437)
@@ -124,6 +124,4 @@
 							   PolygonContainer &coincident)
 {
-	Polygon3 *splitPoly = NULL;
-
 	int splits = 0;
 
@@ -140,7 +138,4 @@
 		// classify polygon
 		const int cf = poly->ClassifyPlane(mPlane);
-
-		Polygon3 *front_piece = NULL;
-		Polygon3 *back_piece = NULL;
 
 		switch (cf)
@@ -156,33 +151,33 @@
 				break;
 			case Polygon3::SPLIT:
-				front_piece = new Polygon3(poly->mParent);
-				back_piece = new Polygon3(poly->mParent);
-
-				//-- split polygon into front and back part
-				poly->Split(mPlane, 
-							*front_piece, 
-							*back_piece);
+				{
+					Polygon3 *front_piece = new Polygon3(poly->mParent);
+					Polygon3 *back_piece = new Polygon3(poly->mParent);
+
+					//-- split polygon into front and back part
+					poly->Split(mPlane, *front_piece, *back_piece);
 					
-				++ splits; // increase number of splits
-
-				//-- inherit rays from parent polygon for blocked ray criterium
-				poly->InheritRays(*front_piece, *back_piece);
-				//Debug << "p: " << poly->mPiercingRays.size() << " f: " << front_piece->mPiercingRays.size() << " b: " << back_piece->mPiercingRays.size() << endl;
+					++ splits; // increase number of splits
+
+					//-- inherit rays from parent polygon for blocked ray criterium
+					poly->InheritRays(*front_piece, *back_piece);
+					//Debug << "p: " << poly->mPiercingRays.size() << " f: " << front_piece->mPiercingRays.size() << " b: " << back_piece->mPiercingRays.size() << endl;
 				
-				// check if polygons still valid 
-				if (front_piece->Valid())
-					frontPolys.push_back(front_piece);
-				else
-					DEL_PTR(front_piece);
+					// check if polygons still valid 
+					if (front_piece->Valid())
+						frontPolys.push_back(front_piece);
+					else
+						DEL_PTR(front_piece);
 				
-				if (back_piece->Valid())
-					backPolys.push_back(back_piece);
-				else				
-					DEL_PTR(back_piece);
+					if (back_piece->Valid())
+						backPolys.push_back(back_piece);
+					else				
+						DEL_PTR(back_piece);
 				
 #ifdef _DEBUG
-				Debug << "split " << *poly << endl << *front_piece << endl << *back_piece << endl;
+					Debug << "split " << *poly << endl << *front_piece << endl << *back_piece << endl;
 #endif
-				DEL_PTR(poly);			
+					DEL_PTR(poly);
+				}
 				break;
 			default:
@@ -273,7 +268,7 @@
 mRootCell(viewCell), 
 mRoot(NULL), 
-mGenerateViewCells(false),
 mStoreLeavesWithRays(false),
-mPvsUseArea(true)
+mPvsUseArea(true),
+mGenerateViewCells(true)
 {
 	Randomize(); // initialise random generator for heuristics
@@ -447,4 +442,6 @@
 	PolygonContainer *polys = new PolygonContainer();
 
+	// don't generate new view cell, insert this one
+	mGenerateViewCells = false;
 	// extract polygons that guide the split process
 	mStat.polys += AddMeshToPolygons(viewCell->GetMesh(), *polys, viewCell);
@@ -457,5 +454,5 @@
 {	
 	std::stack<BspTraversalData> tStack;
-	
+
 	// traverse existing tree or create new tree
     if (!mRoot)
@@ -624,4 +621,6 @@
 	mStat.polys = AddToPolygonSoup(viewCells, *polys);
 
+	// view cells are given
+	mGenerateViewCells = false;
 	// construct tree from the view cell polygons
 	Construct(polys, new BoundedRayContainer());
@@ -635,5 +634,6 @@
 	
 	PolygonContainer *polys = new PolygonContainer();
-	
+
+	mGenerateViewCells = true;
 	// copy mesh instance polygons into one big polygon soup
 	mStat.polys = AddToPolygonSoup(objects, *polys);
@@ -652,4 +652,7 @@
 
 	RayContainer::const_iterator rit, rit_end = sampleRays.end();
+
+	// generate view cells
+	mGenerateViewCells = false;
 
 	long startTime = GetTime();
@@ -726,4 +729,6 @@
 	PolygonContainer *polys = new PolygonContainer();
 	
+	mGenerateViewCells = true;
+
 	// copy mesh instance polygons into one big polygon soup
 	mStat.polys = AddToPolygonSoup(objects, *polys);
@@ -844,10 +849,4 @@
 	PolygonContainer coincident;
 	
-	PolygonContainer *frontPolys = new PolygonContainer();
-	PolygonContainer *backPolys = new PolygonContainer();
-
-	BoundedRayContainer *frontRays = new BoundedRayContainer();
-	BoundedRayContainer *backRays = new BoundedRayContainer();
-	
 	BspTraversalData tFrontData(NULL, new PolygonContainer(), tData.mDepth + 1, mRootCell, 
 								new BoundedRayContainer(), 0, 0, new BspNodeGeometry());
@@ -856,8 +855,6 @@
 
 	// create new interior node and two leaf nodes
-	BspInterior *interior = SubdivideNode(tData,
-										  tFrontData,
-										  tBackData,
-										  coincident);
+	BspInterior *interior = 
+		SubdivideNode(tData, tFrontData, tBackData, coincident);
 
 #ifdef _DEBUG	
@@ -897,5 +894,5 @@
 							   BspTraversalData &backData, 
 							   const PolygonContainer &coincident,
-							   const Plane3 splitPlane) const
+							   const Plane3 &splitPlane) const
 {
 	// if not empty, tree is further subdivided => don't have to find view cell
@@ -1534,5 +1531,6 @@
 			
 			// add the termination object
-			AddObjToPvs(ray->intersections[0].mObject, cf, frontPvs, backPvs);
+			if (!ray->intersections.empty())
+				AddObjToPvs(ray->intersections[0].mObject, cf, frontPvs, backPvs);
 			
 			// add the source object
@@ -1592,4 +1590,5 @@
 
 	const float raysSize = (float)rays.size() + Limits::Small;
+
 	if (mSplitPlaneStrategy & LEAST_RAY_SPLITS)
 		val += mLeastRaySplitsFactor * sumRaySplits / raysSize;
@@ -1604,5 +1603,5 @@
 
 		// give penalty to unbalanced split
-		if (1)
+		if (0)
 		if (((pFront * 0.2 + Limits::Small) > pBack) || (pFront < (pBack * 0.2 + Limits::Small)))
 			val += 0.5;
@@ -1614,4 +1613,5 @@
 		  << " backpvs: " << backPvs << " pBack: " << pBack << endl << endl;
 #endif
+	
 	return val;
 }
@@ -2036,8 +2036,4 @@
 }
 
-void BspTree::SetGenerateViewCells(int generateViewCells)
-{
-	mGenerateViewCells = generateViewCells;
-}
 
 BspTreeStatistics &BspTree::GetStat()
Index: trunk/VUT/GtpVisibilityPreprocessor/src/Plane3.h
===================================================================
--- trunk/VUT/GtpVisibilityPreprocessor/src/Plane3.h	(revision 396)
+++ trunk/VUT/GtpVisibilityPreprocessor/src/Plane3.h	(revision 437)
@@ -15,6 +15,6 @@
   
   Plane3(const Vector3 &a,
-	 const Vector3 &b,
-	 const Vector3 &c
+		 const Vector3 &b,
+		 const Vector3 &c
 	 ) {
     Vector3 v1=a-b, v2=c-b;
@@ -44,4 +44,10 @@
   }
 
+  /** Finds intersection of line segment between points a and b with plane.
+	  @param a start point
+	  @param b end point
+	  @param t if not NULL, returns parameter value of intersections
+	  @param coplanar if not NULL, returns true if plane and line segments are coplanar.
+  */
   Vector3 FindIntersection(const Vector3 &a,
 			   const Vector3 &b,
@@ -56,5 +62,5 @@
 	{
 		if (coplanar) (*coplanar) = true;	
-		if (t) (*t) = 1;
+		if (t) (*t) = 0;
 		return a;
 	}
@@ -68,4 +74,7 @@
   }
   
+  /** Finds value of intersection parameter t on line segment from a to b.
+      @returns -1 if coplanar, else parameter t
+  */
   float FindT(const Vector3 &a,
 			  const Vector3 &b) const 
@@ -74,6 +83,7 @@
 	  const float dv = DotProd(v, mNormal);
     
+	  // does not intersect
 	  if (signum(dv) == 0)
-		  return 1;
+		  return -1;
 	
 	  return - Distance(a) / dv; // TODO: could be done more efficiently
Index: trunk/VUT/GtpVisibilityPreprocessor/src/RayInfo.h
===================================================================
--- trunk/VUT/GtpVisibilityPreprocessor/src/RayInfo.h	(revision 420)
+++ trunk/VUT/GtpVisibilityPreprocessor/src/RayInfo.h	(revision 437)
@@ -7,4 +7,6 @@
 class VssRay;
 class RayInfo;
+class Plane3;
+class Vector3;
 
 typedef vector<RayInfo> RayInfoContainer;
@@ -40,8 +42,20 @@
 	}
 	
+	/** Extracts the scalar of the starting point of the ray segment
+		that lies in the axis.
+	*/
+	float ExtrapOrigin(const int axis) const;
+	/** Extracts the scalar of the termination point of the ray segment
+		that lies in the axis.
+	*/
+	float ExtrapTermination(const int axis) const;
 	
-	float ExtrapOrigin(const int axis) const;
+	/** Extracts the starting point of the ray segment.
+	*/
+	Vector3 ExtrapOrigin() const;
 	
-	float ExtrapTermination(const int axis) const;
+	/** Extracts the end point of the ray segment.
+	*/
+	Vector3 ExtrapTermination() const;
 	
 	float GetMinT () const;
@@ -52,12 +66,25 @@
 	void SetMaxT (const float t);
 
+	float SegmentLength() const;
+	float SqrSegmentLength() const;
+
 	/** Computes intersection of this ray with the axis aligned split plane.
+	
 		@param axis axis of the split plane
-		@param position position of the split plane
-		@param t returns the t value of the ray intersection
+		@param position scalar position of the split plane for the chosen axis
+		@param t returns the t parameter value of the ray intersection
 
 		@returns 0 if ray intersects plane, -1 if on back side of plane, 1 if on front side
 	*/
 	int ComputeRayIntersection(const int axis, const float position, float &t) const;
+
+	/** Computes intersection of this ray with the split plane.
+
+		@param splitPlane the split plane
+		@param t returns the t parameter value of the ray intersection
+
+		@returns 0 if ray intersects plane, -1 if on back side of plane, 1 if on front side
+	*/
+	int ComputeRayIntersection(const Plane3 &splitPlane, float &t) const;
 };
 
Index: trunk/VUT/GtpVisibilityPreprocessor/src/VssPreprocessor.cpp
===================================================================
--- trunk/VUT/GtpVisibilityPreprocessor/src/VssPreprocessor.cpp	(revision 436)
+++ trunk/VUT/GtpVisibilityPreprocessor/src/VssPreprocessor.cpp	(revision 437)
@@ -463,6 +463,4 @@
 	  ObjectContainer objects;
 	  mSceneGraph->CollectObjects(&objects);
-	  
-	  mBspTree->SetGenerateViewCells(true);
 	  mBspTree->Construct(objects, bspRays);
 	}
Index: trunk/VUT/GtpVisibilityPreprocessor/src/RayInfo.cpp
===================================================================
--- trunk/VUT/GtpVisibilityPreprocessor/src/RayInfo.cpp	(revision 432)
+++ trunk/VUT/GtpVisibilityPreprocessor/src/RayInfo.cpp	(revision 437)
@@ -2,4 +2,5 @@
 #include "Ray.h"
 #include "VssRay.h"
+#include "Plane3.h"
 
 RayInfo::RayInfo(): mRay(NULL), mMinT(0),
@@ -51,5 +52,15 @@
 	return mRay->GetOrigin(axis) + GetMaxT()*mRay->GetDir(axis);
 }
+
+Vector3 RayInfo::ExtrapOrigin() const 
+{
+	return mRay->GetOrigin() * mRay->GetDir();
+}
 		
+Vector3 RayInfo::ExtrapTermination() const 
+{
+	return mRay->GetOrigin() + GetMaxT() * mRay->GetDir();
+}
+
 #if USE_FIXEDPOINT_T
 	float RayInfo::GetMinT () const 
@@ -117,2 +128,23 @@
 	return 0;
 }
+
+
+int RayInfo::ComputeRayIntersection(const Plane3 &splitPlane, float &t) const
+{	
+	t = splitPlane.FindT(mRay->GetOrigin(), mRay->GetTermination());
+
+	if ((t < GetMinT()) || (t > GetMaxT()))
+		return splitPlane.Side(ExtrapOrigin());
+
+	return 0;
+}
+
+float RayInfo::SegmentLength() const
+{
+	return Distance(ExtrapOrigin(), ExtrapTermination());
+}
+
+float RayInfo::SqrSegmentLength() const
+{
+	return SqrDistance(ExtrapOrigin(), ExtrapTermination());
+}
Index: trunk/VUT/GtpVisibilityPreprocessor/src/SamplingPreprocessor.cpp
===================================================================
--- trunk/VUT/GtpVisibilityPreprocessor/src/SamplingPreprocessor.cpp	(revision 433)
+++ trunk/VUT/GtpVisibilityPreprocessor/src/SamplingPreprocessor.cpp	(revision 437)
@@ -64,10 +64,8 @@
 	{
 	case BspTree::FROM_INPUT_VIEW_CELLS:
-		mBspTree->SetGenerateViewCells(false);
 		mBspTree->Construct(mViewCells);
 		break;
 	case BspTree::FROM_SCENE_GEOMETRY:
 		DeleteViewCells(); // we generate new view cells
-		mBspTree->SetGenerateViewCells(true);
 		mSceneGraph->CollectObjects(&objects);
 		mBspTree->Construct(objects);
@@ -75,5 +73,4 @@
 	case BspTree::FROM_SAMPLES:
 		DeleteViewCells(); // we generate new view cells
-		mBspTree->SetGenerateViewCells(true);
 		mBspTree->Construct(mSampleRays);
 		break;
