Changeset 1370 for GTP/trunk/Lib/Vis/Preprocessing/src/BvHierarchy.cpp

Timestamp:

09/14/06 18:55:38 (18 years ago)

Author:

mattausch

Message:

debugged global sorting, worked on object-viewspace subdivision

File:

• GTP/trunk/Lib/Vis/Preprocessing/src/BvHierarchy.cpp (modified) (33 diffs)

GTP/trunk/Lib/Vis/Preprocessing/src/BvHierarchy.cpp

@@ -17 +17 @@
 #include "Beam.h"
 #include "VspTree.h"
+#include "HierarchyManager.h"
 
 
@@ -22 +23 @@
 
 
+#define PROBABILIY_IS_BV_VOLUME 1
 #define USE_FIXEDPOINT_T 0
 
@@ -197 +199 @@
                 Randomize(); // initialise random generator for heuristics
 
+
+        /////////////////////////////////////////////////////////////
         //-- termination criteria for autopartition
         Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.maxDepth", mTermMaxDepth);
         Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.maxLeaves", mTermMaxLeaves);
         Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.minObjects", mTermMinObjects);
+        Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.minRays", mTermMinRays);
         Environment::GetSingleton()->GetFloatValue("BvHierarchy.Termination.minProbability", mTermMinProbability);
         
         Environment::GetSingleton()->GetIntValue("BvHierarchy.Termination.missTolerance", mTermMissTolerance);
-        
+
+
+        ////////////////////////////////////////////////
         //-- max cost ratio for early tree termination
+
         Environment::GetSingleton()->GetFloatValue("BvHierarchy.Termination.maxCostRatio", mTermMaxCostRatio);
-
         Environment::GetSingleton()->GetFloatValue("BvHierarchy.Termination.minGlobalCostRatio",
                 mTermMinGlobalCostRatio);
@@ -213 +220 @@
                 mTermGlobalCostMissTolerance);
 
-        //-- factors for bsp tree split plane heuristics
-
-        // if only the driving axis is used for axis aligned split
+
+        /////////////////////////////////////////
+        //-- factors for subdivision heuristics
+
+        // if only the driving axis is used for splits
         Environment::GetSingleton()->GetBoolValue("BvHierarchy.splitUseOnlyDrivingAxis", mOnlyDrivingAxis);
         Environment::GetSingleton()->GetFloatValue("BvHierarchy.maxStaticMemory", mMaxMemory);
@@ -263 +272 @@
 }
 
+static int CountRays(const ObjectContainer &objects)
+{
+        int nRays = 0;
+
+        ObjectContainer::const_iterator oit, oit_end = objects.end();
+
+        for (oit = objects.begin(); oit != oit_end; ++ oit)
+        {
+                nRays += (int)(*oit)->mVssRays.size();
+        }
+
+        return nRays;
+}
 
 BvhInterior *BvHierarchy::SubdivideNode(const BvhSubdivisionCandidate &sc,
@@ -268 +290 @@
                                                                                 BvhTraversalData &backData)
 {
+        mBvhStats.nodes += 2; // we have two new leaves
+
         const BvhTraversalData &tData = sc.mParentData;
         BvhLeaf *leaf = tData.mNode;
-    mBvhStats.nodes += 2; // we have two new leaves
+        AxisAlignedBox3 parentBox = leaf->GetBoundingBox();
 
         // add the new nodes to the tree
-        BvhInterior *node = new BvhInterior(tData.mBoundingBox, leaf->GetParent());
-        
+        BvhInterior *node = new BvhInterior(parentBox, leaf->GetParent());
+        
+
+        ////////////////////////////////
+        //-- create front and back leaf
+
+        AxisAlignedBox3 fbox = ComputeBoundingBox(sc.mFrontObjects, &parentBox);
+        AxisAlignedBox3 bbox = ComputeBoundingBox(sc.mBackObjects, &parentBox);
+
+        BvhLeaf *back = 
+                new BvhLeaf(bbox, node, (int)sc.mBackObjects.size());
+        BvhLeaf *front = 
+                new BvhLeaf(fbox, node, (int)sc.mFrontObjects.size());
+
+        BvhInterior *parent = leaf->GetParent();
+
+        // replace a link from node's parent
+        if (parent)
+        {
+                parent->ReplaceChildLink(leaf, node);
+                node->SetParent(parent);
+        }
+        else // no parent => this node is the root
+        {
+                mRoot = node;
+        }
+
+        // and setup child links
+        node->SetupChildLinks(front, back);
+
+        ++ mBvhStats.splits;
+
 
         ///////////////////////////////////////////////////////////////////
@@ -281 +335 @@
         frontData.mDepth = backData.mDepth = tData.mDepth + 1;
 
-        frontData.mBoundingBox = ComputeBoundingBox(sc.mFrontObjects, &tData.mBoundingBox);
-        backData.mBoundingBox = ComputeBoundingBox(sc.mBackObjects, &tData.mBoundingBox);
-        
-
-        ////////////////////////////////
-        //-- create front and back leaf
-
-        BvhLeaf *back = 
-                new BvhLeaf(backData.mBoundingBox, node, (int)sc.mBackObjects.size());
-        BvhLeaf *front = 
-                new BvhLeaf(frontData.mBoundingBox, node, (int)sc.mFrontObjects.size());
-
-        BvhInterior *parent = leaf->GetParent();
-
-        // replace a link from node's parent
-        if (parent)
-        {
-                parent->ReplaceChildLink(leaf, node);
-                node->SetParent(parent);
-        }
-        else // no parent => this node is the root
-        {
-                mRoot = node;
-        }
-
-        // and setup child links
-        node->SetupChildLinks(front, back);
-
-        ++ mBvhStats.splits;
-
-        ////////////////////////////////////
-        //-- fill traversal data
-
         frontData.mNode = front;
         backData.mNode = back;
@@ -319 +340 @@
         back->mObjects = sc.mBackObjects;
         front->mObjects = sc.mFrontObjects;
+
+        // if the number of rays is too low, no assumptions can be made
+        // (=> switch to surface area heuristics?)
+        frontData.mNumRays = CountRays(sc.mFrontObjects);
+        backData.mNumRays = CountRays(sc.mBackObjects);
 
         AssociateObjectsWithLeaf(back);
         AssociateObjectsWithLeaf(front);
    
+#if PROBABILIY_IS_BV_VOLUME
+        // volume of bvh (= probability that this bvh can be seen)
+        frontData.mProbability = fbox.GetVolume();
+        backData.mProbability = bbox.GetVolume();
+#else
         // compute probability of this node being visible, 
         // i.e., volume of the view cells that can see this node
         frontData.mProbability = EvalViewCellsVolume(sc.mFrontObjects);
         backData.mProbability = EvalViewCellsVolume(sc.mBackObjects);
+#endif
 
     // how often was max cost ratio missed in this branch?
@@ -389 +421 @@
                 tBackData.mNode->SetSubdivisionCandidate(backCandidate);
 
-                Debug << "leaf: " << tFrontData.mNode << " setting f candidate: " 
-                          << tFrontData.mNode->GetSubdivisionCandidate() << " type: " 
-                          << tFrontData.mNode->GetSubdivisionCandidate()->Type() << endl;
-
-                Debug << "leaf: " << tBackData.mNode << " setting b candidate: " 
-                          << tBackData.mNode->GetSubdivisionCandidate() << " type: " 
-                          << tBackData.mNode->GetSubdivisionCandidate()->Type() << endl;
-                
                 tQueue.Push(frontCandidate);
                 tQueue.Push(backCandidate);
@@ -460 +484 @@
         const float renderCostDecr = oldRenderCost - newRenderCost;
 
-        Debug << "\nbvh render cost decr: " << renderCostDecr << endl;
+        //Debug << "\nbvh render cost decr: " << renderCostDecr << endl;
         splitCandidate.SetRenderCostDecrease(renderCostDecr);
 
@@ -481 +505 @@
         // matt: TODO
         return ( 0
-                //|| ((int)data.mNode->mObjects.size() < mTermMinObjects)
-                //|| (data.mProbability <= mTermMinProbability)
-                //|| (data.mDepth >= mTermMaxDepth)
+                || ((int)data.mNode->mObjects.size() < mTermMinObjects)
+                || (data.mProbability <= mTermMinProbability)
+                || (data.mDepth >= mTermMaxDepth)
+                || (data.mNumRays <= mTermMinRays)
                  );
 }
@@ -493 +518 @@
         return (0
                 || (mBvhStats.Leaves() >= mTermMaxLeaves)
-                //|| (mGlobalCostMisses >= mTermGlobalCostMissTolerance)
+                || (mGlobalCostMisses >= mTermGlobalCostMissTolerance)
                 //|| mOutOfMemory 
                 );
@@ -506 +531 @@
         ++ mCreatedLeaves;
 
+        
+        if (data.mProbability <= mTermMinProbability)
+        {
+                ++ mBvhStats.minProbabilityNodes;
+        }
+
+        ////////////////////////////////////////////
+        // depth related stuff
+
+        if (data.mDepth > mBvhStats.maxDepth)
+        {
+                mBvhStats.maxDepth = data.mDepth;
+        }
+
+        if (data.mDepth < mBvhStats.minDepth)
+        {
+                mBvhStats.minDepth = data.mDepth;
+        }
+
         if (data.mDepth >= mTermMaxDepth)
         {
         ++ mBvhStats.maxDepthNodes;
-                //Debug << "new max depth: " << mVspStats.maxDepthNodes << endl;
-        }
-
-        if (data.mDepth < mTermMaxDepth)
-        {
-        ++ mBvhStats.minDepthNodes;
-        }
-
-        if (data.mProbability <= mTermMinProbability)
-                ++ mBvhStats.minProbabilityNodes;
-        
+        }
+
         // accumulate depth to compute average depth
         mBvhStats.accumDepth += data.mDepth;
- 
-        if ((int)(leaf->mObjects.size()) < mTermMinObjects)
+
+
+        ////////////////////////////////////////////
+        // objects related stuff
+
+        // note: this number should always accumulate to the total number of objects
+        mBvhStats.objectRefs += (int)leaf->mObjects.size();
+
+        if ((int)leaf->mObjects.size() <= mTermMinObjects)
+        {
              ++ mBvhStats.minObjectsNodes;
- 
-        if ((int)(leaf->mObjects.size()) > mBvhStats.maxObjectRefs)
+        }
+
+        if ((int)leaf->mObjects.size() > mBvhStats.maxObjectRefs)
+        {
                 mBvhStats.maxObjectRefs = (int)leaf->mObjects.size();
+        }
+
+        if ((int)leaf->mObjects.size() < mBvhStats.minObjectRefs)
+        {
+                mBvhStats.minObjectRefs = (int)leaf->mObjects.size();
+        }
+
+        ////////////////////////////////////////////
+        // ray related stuff
+
+        // note: this number should always accumulate to the total number of rays
+        mBvhStats.rayRefs += data.mNumRays;
+        
+        if (data.mNumRays <= mTermMinRays)
+        {
+             ++ mBvhStats.minRaysNodes;
+        }
+
+        if (data.mNumRays > mBvhStats.maxRayRefs)
+        {
+                mBvhStats.maxRayRefs = data.mNumRays;
+        }
+
+        if (data.mNumRays < mBvhStats.minRayRefs)
+        {
+                mBvhStats.minRayRefs = data.mNumRays;
+        }
+
+        cout << "depth: " << data.mDepth << " objects: " << (int)leaf->mObjects.size() 
+                 << " rays: " << data.mNumRays << " rays / objects " 
+                 << (float)data.mNumRays / (float)leaf->mObjects.size() << endl;
 }
 
 
 #if 0
+
+/// compute object boundaries using spatial mid split
 float BvHierarchy::EvalLocalObjectPartition(const BvhTraversalData &tData,
                                                                                         const int axis,
@@ -553 +631 @@
                 // object mailed => belongs to back objects
                 if (objMid < midPoint) 
+                {
                         objectsBack.push_back(obj);
+                }
                 else
+                {
                         objectsFront.push_back(obj);
-        }
-
-        const float oldRenderCost = tData.mProbability * (float)tData.mNode->mObjects.size();
+                }
+        }
+
+        const float oldProp = EvalViewCellsVolume(tData.mNode->mObjects);
+        //const float oldProp = tData.mProbability;
+
+        const float oldRenderCost = oldProb * (float)tData.mNode->mObjects.size();
         const float newRenderCost = 
                 EvalRenderCost(tData, objectsFront, objectsBack);
@@ -568 +653 @@
 #else
 
+/// compute object partition by getting balanced objects on the left and right side
 float BvHierarchy::EvalLocalObjectPartition(const BvhTraversalData &tData,
                                                                                         const int axis,
@@ -592 +678 @@
 
         const float oldProp = EvalViewCellsVolume(tData.mNode->mObjects);
-        //const float oldProp2 = tData.mProbability;
+        //const float oldProp = tData.mProbability;
 
         const float oldRenderCost = oldProp * (float)tData.mNode->mObjects.size();
@@ -649 +735 @@
 
         vector<float>::const_iterator bit = bordersRight.begin();
-
+        cout << "here42" << endl;
         SortableEntryContainer::const_iterator cit, cit_end = mSubdivisionCandidates->end();
         for (cit = mSubdivisionCandidates->begin(); cit != cit_end; ++ cit, ++ bit) 
@@ -671 +757 @@
                 const float sum = objectsLeft * lbox.SurfaceArea() + objectsRight * rbox.SurfaceArea();
       
-            // cout<<"pos="<<(*ci).value<<"\t q=("<<ql<<","<<qr<<")\t r=("<<rl<<","<<rr<<")"<<endl;
-            // cout<<"cost= "<<sum<<endl;
+                cout << "pos=" << (*cit).mPos << "\t q=(" << objectsLeft << "," << objectsRight <<")\t r=(" 
+                         << lbox.SurfaceArea() << "," << rbox.SurfaceArea() << ")" << endl;
+            cout << "cost= " << sum << endl;
       
                 if (sum < minSum) 
                 {
                         minSum = sum;   
-                        // objects belongs to left side now
+                        // objects belong to left side now
                         for (; currentPos != (cit + 1); ++ currentPos);
                 }
         }
 
+        ////////////////////////////////////////////
         //-- assign object to front and back volume
 
@@ -796 +884 @@
         {
                 Intersectable *object = (*cit).mObject;
-                
+        
                 // evaluate change in l and r volume
                 // voll = view cells that see only left node (i.e., left pvs)
@@ -864 +952 @@
 {
         //-- insert object queries
-        ObjectContainer *objects;
-
-        if (!mUseGlobalSorting)
-                objects = &tData.mNode->mObjects;
-        else
-                objects = tData.mSortedObjects[axis];
-        
-        CreateLocalSubdivisionCandidates(*objects, &mSubdivisionCandidates, mUseGlobalSorting, axis);
+        ObjectContainer *objects = mUseGlobalSorting ? tData.mSortedObjects[axis] : &tData.mNode->mObjects;
+
+        CreateLocalSubdivisionCandidates(*objects, &mSubdivisionCandidates, !mUseGlobalSorting, axis);
 }
 
@@ -937 +1020 @@
         }
 
-        // mail the objects on the left side
-        Intersectable::NewMail();
-        // mail view cells on the left side
+        // we will mail view cells switching to the back side
         ViewCell::NewMail();
         
@@ -1009 +1090 @@
         ObjectContainer nFrontObjects[3];
         ObjectContainer nBackObjects[3];
-
         float nCostRatio[3];
 
-        // create bounding box of node geometry
-        AxisAlignedBox3 box = tData.mBoundingBox;
-                
         int sAxis = 0;
         int bestAxis = -1;
@@ -1020 +1097 @@
         if (mOnlyDrivingAxis)
         {
+                const AxisAlignedBox3 box = tData.mNode->GetBoundingBox();
                 sAxis = box.Size().DrivingAxis();
         }
@@ -1032 +1110 @@
                         if (mUseCostHeuristics)
                         {
-                                //-- partition objects using heuristics
-                                nCostRatio[axis] =
-                                        EvalLocalCostHeuristics(
-                                                                                        tData,
-                                                                                        axis,
-                                                                                        nFrontObjects[axis],
-                                                                                        nBackObjects[axis]);
+                                //////////////////////////////////
+                //-- split objects using heuristics
+                                
+                                if (mHierarchyManager->GetViewSpaceSubdivisionType() == 
+                                        HierarchyManager::KD_BASED_VIEWSPACE_SUBDIV)
+                                {
+                                        //-- heuristics using objects weighted by view cells volume
+                                        nCostRatio[axis] =
+                                                EvalLocalCostHeuristics(
+                                                                                                tData,
+                                                                                                axis,
+                                                                                                nFrontObjects[axis],
+                                                                                                nBackObjects[axis]);
+                                }
+                                else
+                                {
+                                        //-- use surface area heuristic because view cells not constructed yet                          
+                                        nCostRatio[axis] =
+                                                EvalLocalCostHeuristics(
+                                                                                                tData,
+                                                                                                axis,
+                                                                                                nFrontObjects[axis],
+                                                                                                nBackObjects[axis]);
+                                }
                         }
                         else
                         {
+                                //-- split objects using some simple criteria
+
                                 nCostRatio[axis] =
                                         EvalLocalObjectPartition(
@@ -1072 +1169 @@
 
 
-void BvHierarchy::AssociateObjectsWithRays(const VssRayContainer &rays)
-{
-
+int BvHierarchy::AssociateObjectsWithRays(const VssRayContainer &rays) const
+{
+        int nRays = 0;
         VssRayContainer::const_iterator rit, rit_end = rays.end();
+
+        VssRay::NewMail();
 
     for (rit = rays.begin(); rit != rays.end(); ++ rit)
@@ -1084 +1183 @@
                 {
                         ray->mTerminationObject->mVssRays.push_back(ray);
-                }
-
-                if (0 && ray->mOriginObject)
+                        if (!ray->Mailed())
+                        {
+                                ray->Mail();
+                                ++ nRays;
+                        }
+                }
+
+                if (1 && ray->mOriginObject)
                 {
                         ray->mOriginObject->mVssRays.push_back(ray);
-                }
-        }
+
+                        if (!ray->Mailed())
+                        {
+                                ray->Mail();
+                                ++ nRays;
+                        }
+                }
+        }
+
+        return nRays;
 }
 
@@ -1137 +1249 @@
                                                                   const ObjectContainer &objectsBack) const
 {
-        BvhLeaf *leaf = tData.mNode;
-
         // probability that view point lies in a view cell which sees this node
         const float pFront = EvalViewCellsVolume(objectsFront);
         const float pBack = EvalViewCellsVolume(objectsBack);
                 
-        const int totalObjects = (int)leaf->mObjects.size();
-        const int nObjectsFront = (int)objectsFront.size();
-        const int nObjectsBack = (int)objectsBack.size();
-
         //-- pvs rendering heuristics
-        const float newRenderCost = nObjectsFront * pFront + nObjectsBack * pBack;
+        const float newRenderCost = (int)objectsFront.size() * pFront + (int)objectsBack.size() * pBack;
 
 #ifdef _DEBUG
@@ -1171 +1277 @@
         ObjectContainer::const_iterator oit, oit_end = objects.end();
 
-        //-- compute bounding box
         for (oit = objects.begin(); oit != oit_end; ++ oit)
         {
                 Intersectable *obj = *oit;
 
-                // compute bounding box of view space
+                // grow bounding box to include all objects
                 box.Include(obj->GetBox());
         }
@@ -1496 +1601 @@
         BvhSubdivisionCandidate::sBvHierarchy = this;
         mBvhStats.nodes = 1;
+        mGlobalCostMisses = 0;
 
         // compute bounding box from objects
@@ -1502 +1608 @@
         BvhLeaf *bvhleaf = dynamic_cast<BvhLeaf *>(mRoot);
 
+        // multiply termination criterium for comparison,
+        // so it can be set between zero and one and
+        // no division is necessary during traversal
+
+#if PROBABILIY_IS_BV_VOLUME
         mTermMinProbability *= mBoundingBox.GetVolume();
-        mGlobalCostMisses = 0;
-        
+        // probability that bounding volume is seen
+        const float prop = GetBoundingBox().GetVolume();
+#else
+        mTermMinProbability *= mVspTree->GetBoundingBox().GetVolume();
+        // probability that volume is "seen" from the view cells
+        const float prop = EvalViewCellsVolume(objects);
+#endif
+
         // only rays intersecting objects in node are interesting
-        AssociateObjectsWithRays(sampleRays);
-        
-        // probabilty is voume of all "seen" view cells
-#if 1
-        const float prop = EvalViewCellsVolume(objects);
-#else
-        const float prop = GetBoundingBox().GetVolume();
-#endif
+        const int nRays = AssociateObjectsWithRays(sampleRays);
+        //Debug << "using " << nRays << " of " << (int)sampleRays.size() << " rays" << endl;
 
         // create bvh traversal data
-        BvhTraversalData oData(bvhleaf, 0, mBoundingBox, prop);
+        BvhTraversalData oData(bvhleaf, 0, prop, nRays);
 
         // create sorted object lists for the first data
@@ -1540 +1651 @@
 
         return oSubdivisionCandidate;
-}
-
-
-bool BvHierarchy::AddLeafToPvs(BvhLeaf *leaf, 
-                                                           ViewCell *vc, 
-                                                           const float pdf, 
-                                                           float &contribution)
-{
-        // add kd intersecable to pvs
-        BvhIntersectable *bvhObj = GetOrCreateBvhIntersectable(leaf);
-        
-        return vc->AddPvsSample(bvhObj, pdf, contribution);
 }
 
@@ -1602 +1701 @@
                 backData.mSortedObjects[i]->reserve((int)sc.mFrontObjects.size());
 
-                ObjectContainer::const_iterator oit, oit_end = sc.mParentData.mNode->mObjects.end();
-
-                for (oit = sc.mParentData.mNode->mObjects.begin(); oit != oit_end; ++ oit)
+                ObjectContainer::const_iterator oit, oit_end = sc.mParentData.mSortedObjects[i]->end();
+
+                for (oit = sc.mParentData.mSortedObjects[i]->begin(); oit != oit_end; ++ oit)
                 {
                         if ((*oit)->Mailed())
@@ -1635 +1734 @@
         app << "#AXIS_ALIGNED_SPLITS (number of axis aligned splits)\n" << splits << endl;
 
-        app << "#N_PMINDEPTHLEAVES ( Percentage of leaves at minimum depth )\n" 
-                <<      minDepthNodes * 100 / (double)Leaves() << endl;
-
+        app << "#N_MAXCOSTNODES  ( Percentage of leaves with terminated because of max cost ratio )\n"
+                << maxCostNodes * 100 / (double)Leaves() << endl;
+
+        app << "#N_PMINPROBABILITYLEAVES  ( Percentage of leaves with mininum probability )\n"
+                << minProbabilityNodes * 100 / (double)Leaves() << endl;
+
+
+        //////////////////////////////////////////////////
+        
         app << "#N_PMAXDEPTHLEAVES ( Percentage of leaves at maximum depth )\n" 
                 <<      maxDepthNodes * 100 / (double)Leaves() << endl;
-
-        app << "#N_MAXCOSTNODES  ( Percentage of leaves with terminated because of max cost ratio )\n"
-                << maxCostNodes * 100 / (double)Leaves() << endl;
-
-        app << "#N_PMINPROBABILITYLEAVES  ( Percentage of leaves with mininum probability )\n"
-                << minProbabilityNodes * 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_PMINOBJECTSLEAVES  ( Percentage of leaves with mininum objects )\n"
                 << minObjectsNodes * 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_INVALIDLEAVES (number of invalid leaves )\n" << invalidLeaves << endl;
-
         app << "#N_MAXOBJECTREFS  ( Max number of object refs / leaf )\n" << maxObjectRefs << "\n";
 
-        //app << "#N_RAYS (number of rays / leaf)\n" << AvgRays() << endl;
-        
-        app << "========== END OF VspTree statistics ==========\n";
-}
-
-
-}
+        app << "#N_MINOBJECTREFS  ( Min number of object refs / leaf )\n" << minObjectRefs << "\n";
+        
+        app << "#N_PAVGOBJECTSLEAVES  ( average object refs / leaf)\n" << AvgObjectRefs() << endl;
+
+
+        ////////////////////////////////////////////////////////
+        
+        app << "#N_PMINRAYSLEAVES  ( Percentage of leaves with mininum rays )\n"
+                << minRaysNodes * 100 / (double)Leaves() << endl;
+
+        app << "#N_MAXRAYREFS  ( Max number of ray refs / leaf )\n" << maxRayRefs << "\n";
+
+        app << "#N_MINRAYREFS  ( Min number of ray refs / leaf )\n" << minRayRefs << "\n";
+        
+        app << "#N_PAVGRAYLEAVES  ( average ray refs / leaf )\n" << AvgRayRefs() << endl;
+        
+        app << "#N_PAVGRAYCONTRIBLEAVES  ( Average ray contribution)\n" <<
+                rayRefs / (double)objectRefs << endl;
+
+        app << "#N_PMAXRAYCONTRIBLEAVES  ( Percentage of leaves with maximal ray contribution )\n"<<
+                maxRayContriNodes * 100 / (double)Leaves() << endl;
+
+        app << "========== END OF BvHierarchy statistics ==========\n";
+}
+
+
+}