source: GTP/trunk/Lib/Vis/Preprocessing/src/GvsPreprocessor.cpp @ 2745

Revision 2745, 57.9 KB checked in by mattausch, 16 years ago (diff)
  • Property svn:executable set to *
Line 
1#include "Environment.h"
2#include "GvsPreprocessor.h"
3#include "GlRenderer.h"
4#include "VssRay.h"
5#include "ViewCellsManager.h"
6#include "Triangle3.h"
7#include "IntersectableWrapper.h"
8#include "Plane3.h"
9#include "RayCaster.h"
10#include "Exporter.h"
11#include "SamplingStrategy.h"
12#include "BvHierarchy.h"
13#include "Polygon3.h"
14#include "IntersectableWrapper.h"
15#include "Timer/PerfTimer.h"
16#include "Vector2.h"
17#include "RndGauss.h"
18
19
20
21namespace GtpVisibilityPreprocessor
22{
23 
24#define GVS_DEBUG 1
25#define NOT_ACCOUNTED_OBJECT 0
26#define ACCOUNTED_OBJECT 2
27#define USE_REVERSE_SAMPLING 1
28#define USE_SUBDIVISION 1
29
30static const float MIN_DIST = 1e-8f;
31
32static ObjectContainer myobjects;
33static int sInvalidSamples = 0;
34static VizStruct currentViz;
35
36static int successfullReverseRays = 0;
37
38
39///////////
40//-- timers
41
42
43static PerfTimer rayTimer;
44static PerfTimer kdPvsTimer;
45static PerfTimer gvsTimer;
46static PerfTimer initialTimer;
47static PerfTimer intersectionTimer;
48static PerfTimer preparationTimer;
49static PerfTimer mainLoopTimer;
50static PerfTimer contributionTimer;
51static PerfTimer castTimer;
52static PerfTimer generationTimer;
53
54
55/////////////////////
56
57
58
59
60
61GvsPreprocessor::GvsPreprocessor():
62Preprocessor(),
63mProcessedViewCells(0),
64mCurrentViewCell(NULL),
65mCurrentViewPoint(Vector3(0.0f, 0.0f, 0.0f)),
66mOnlyRandomSampling(false),
67//mOnlyRandomSampling(true),
68mUseProbablyVisibleSampling(false)
69{
70        Environment::GetSingleton()->GetIntValue("GvsPreprocessor.totalSamples", mTotalSamples);
71        Environment::GetSingleton()->GetIntValue("GvsPreprocessor.initialSamples", mInitialSamples);
72        Environment::GetSingleton()->GetIntValue("GvsPreprocessor.gvsSamplesPerPass", mGvsSamplesPerPass);
73        Environment::GetSingleton()->GetIntValue("GvsPreprocessor.minContribution", mMinContribution);
74        Environment::GetSingleton()->GetFloatValue("GvsPreprocessor.epsilon", mEps);
75        Environment::GetSingleton()->GetFloatValue("GvsPreprocessor.threshold", mThreshold);   
76        Environment::GetSingleton()->GetBoolValue("GvsPreprocessor.perViewCell", mPerViewCell);
77        Environment::GetSingleton()->GetIntValue("GvsPreprocessor.maxViewCells", mMaxViewCells);
78
79        Environment::GetSingleton()->GetBoolValue("Preprocessor.evaluatePixelError", mEvaluatePixelError);
80
81        Environment::GetSingleton()->GetBoolValue("ViewCells.useKdPvs", mUseKdPvs);
82
83        Environment::GetSingleton()->GetFloatValue("GvsPreprocessor.initialJitter", mInitialJitter);
84        Environment::GetSingleton()->GetBoolValue("GvsPreprocessor.useDeterministicGvs", mUseDeterministicGvs);
85
86        Environment::GetSingleton()->GetFloatValue("GvsPreprocessor.radiusOfInfluence", mRadiusOfInfluence);
87
88        char gvsStatsLog[100];
89        Environment::GetSingleton()->GetStringValue("GvsPreprocessor.stats", gvsStatsLog);
90        mGvsStatsStream.open(gvsStatsLog);
91
92        cout << "number of gvs samples per pass: " << mGvsSamplesPerPass << endl;
93        cout << "number of samples per pass: " << mSamplesPerPass << endl;
94        cout << "only random sampling: " << mOnlyRandomSampling << endl;
95
96        Debug << "Gvs preprocessor options" << endl;
97        Debug << "number of total samples: " << mTotalSamples << endl;
98        Debug << "number of initial samples: " << mInitialSamples << endl;
99        Debug << "threshold: " << mThreshold << endl;
100        Debug << "epsilon: " << mEps << endl;
101        Debug << "stats: " << gvsStatsLog << endl;
102        Debug << "per view cell: " << mPerViewCell << endl;
103        Debug << "max view cells: " << mMaxViewCells << endl;
104        Debug << "min contribution: " << mMinContribution << endl;
105
106        mDistribution = new ViewCellBasedDistribution(*this, NULL);
107
108        mGvsStats.Reset();
109
110        // hack: some generic statistical values that can be used
111        // by an application using the preprocessor
112        for (int i = 0; i < 2; ++ i)
113                mGenericStats[i] = 0;
114}
115
116
117GvsPreprocessor::~GvsPreprocessor()
118{
119        delete mDistribution;
120        ClearRayQueue();
121}
122
123
124void GvsPreprocessor::ClearRayQueue()
125{
126        // clean ray queue
127        while (!mRayQueue.empty())
128        {
129                // handle next ray
130                VssRay *ray = mRayQueue.top();
131                mRayQueue.pop();
132
133                // note: deletion of the ray is handled by ray pool
134                mRayCaster->SheduleRayForDeletion(ray);
135        }
136}
137
138
139ViewCell *GvsPreprocessor::NextViewCell()
140{
141        if (mProcessedViewCells == (int)mViewCells.size())
142                return NULL; // no more view cells
143
144        ViewCell *vc = mViewCells[mProcessedViewCells];
145
146   if (!vc->GetMesh())
147                mViewCellsManager->CreateMesh(vc);
148
149        mGvsStats.mViewCellId = vc->GetId();
150
151        Debug << "vc: " << vc->GetId() << endl;
152
153        ++ mProcessedViewCells;
154   
155        return vc;
156}
157
158
159int GvsPreprocessor::CheckDiscontinuity(const VssRay &currentRay,
160                                                                                const Triangle3 &hitTriangle,
161                                                                                const VssRay &oldRay)
162{
163        // the predicted hitpoint is the point where the ray would have intersected the hit triangle
164        Vector3 predictedHit = CalcPredictedHitPoint(currentRay, hitTriangle, oldRay);
165
166        float predictedLen = Magnitude(predictedHit - currentRay.mOrigin);
167        float len = Magnitude(currentRay.mTermination - currentRay.mOrigin);
168       
169        // discrepancy between predicted and actual hit point large => this is likely to be a discontinuity
170
171        // matt: problem: this only finds CLOSER discontinuities, but not discontinuities which are FARTHER away
172        if ((predictedLen - len) < mThreshold)
173        // q: rather use relative distance?
174        //if ((predictedLen / len) < mThreshold)
175                return 0;
176
177        SimpleRay simpleRay;
178
179        // apply reverse sampling to find the gap
180        if (!ComputeReverseRay(currentRay, hitTriangle, oldRay, simpleRay))
181                return 0;
182
183        static VssRayContainer reverseRays;
184        reverseRays.clear();
185
186        if (mUseDeterministicGvs)
187        {
188                VssRay *reverseRay =
189                        mRayCaster->CastRay(simpleRay, mViewCellsManager->GetViewSpaceBox(), !mPerViewCell);
190
191                reverseRays.push_back(reverseRay);
192        }
193        else
194        {
195                float scale = 0.01f;
196
197                if (0)
198                        CastRayBundle4(simpleRay, reverseRays, mViewCellsManager->GetViewSpaceBox());
199                else
200                        CastRayBundle16(simpleRay, reverseRays, scale);
201        }
202
203        mGvsStats.mReverseSamples += (int)reverseRays.size();
204        successfullReverseRays = EnqueueRays(reverseRays);
205
206        return (int)reverseRays.size();
207}
208
209
210bool GvsPreprocessor::CheckDiscontinuity2(const VssRay &currentRay,
211                                                                                  const Triangle3 &hitTriangle,
212                                                                                  const VssRay &oldRay,
213                                                                                  SimpleRay &reverseRay)
214{
215        // the predicted hitpoint is the point where the ray would have intersected the hit triangle
216        Vector3 predictedHit = CalcPredictedHitPoint(currentRay, hitTriangle, oldRay);
217
218        float predictedLen = Magnitude(predictedHit - currentRay.mOrigin);
219        float len = Magnitude(currentRay.mTermination - currentRay.mOrigin);
220       
221        // discrepancy between predicted and actual hit point large
222        // => this is likely to be a discontinuity
223        if (fabs(predictedLen - len) < mThreshold)
224                return false;
225
226        return CreateRandomizedReverseRay(currentRay, predictedHit, reverseRay, hitTriangle);
227}
228
229
230int GvsPreprocessor::CountObject(Intersectable *triObj)
231{
232        int numObjects = 0;
233
234        if ((triObj->mCounter != (NOT_ACCOUNTED_OBJECT + 1)) &&
235                (triObj->mCounter != (ACCOUNTED_OBJECT + 1)))
236        {
237                ++ triObj->mCounter;
238                ++ numObjects;
239        }
240
241        mGenericStats[1] += numObjects;
242
243        return numObjects;
244}
245
246
247void GvsPreprocessor::UpdateStatsForVisualization(KdIntersectable *kdInt)
248{
249        int numObj = 0;
250
251        // count new objects in pvs due to kd node conversion   
252        myobjects.clear();
253        // also gather duplicates to avoid mailing
254        mKdTree->CollectObjectsWithDublicates(kdInt->GetItem(), myobjects);
255
256        ObjectContainer::const_iterator oit, oit_end = myobjects.end();
257
258        for (oit = myobjects.begin(); oit != oit_end; ++ oit)
259                numObj += CountObject(*oit);
260
261        mViewCellsManager->UpdateStatsForViewCell(mCurrentViewCell, kdInt, numObj);
262}       
263
264
265void GvsPreprocessor::AddKdNodeToPvs(const Vector3 &termination)
266{
267        kdPvsTimer.Entry();
268
269        // exchange the triangle with the node in the pvs for kd pvs
270        KdNode *node = mKdTree->GetPvsNode(termination);
271
272        //KdNode *node = mKdTree->GetPvsNodeCheck(ray.mTermination, obj);
273        //if (!node) cerr << "e " << obj->GetBox() << " " << ray.mTermination << endl; else
274        if (!node->Mailed())
275        {
276                // add to pvs
277                node->Mail();
278                KdIntersectable *kdInt = mKdTree->GetOrCreateKdIntersectable(node);
279                mCurrentViewCell->GetPvs().AddSampleDirty(kdInt, 1.0f);
280
281                if (QT_VISUALIZATION_SHOWN) UpdateStatsForVisualization(kdInt);
282        }
283
284        kdPvsTimer.Exit();
285}
286
287
288bool GvsPreprocessor::HasContribution(VssRay &ray)
289{
290        if (!ray.mTerminationObject)
291                return false;
292
293        contributionTimer.Entry();
294
295        bool result;
296
297        if (!mPerViewCell)
298        {
299                // store the rays + the intersected view cells
300                const bool storeViewCells = false;
301
302                mViewCellsManager->ComputeSampleContribution(ray,
303                                                                                                         true,
304                                                                                                         storeViewCells,
305                                                                                                         true);
306
307                result = ray.mPvsContribution > 0;
308        }
309        else // original per view cell gvs
310        {
311                // test if triangle was accounted for yet
312                result = AddTriangleObject(ray.mTerminationObject);
313               
314                if (mUseKdPvs)
315                        AddKdNodeToPvs(ray.mTermination);
316        }
317
318        contributionTimer.Exit();
319
320        return result;
321}
322
323
324bool GvsPreprocessor::HandleRay(VssRay *vssRay)
325{
326        if (!vssRay)
327                return false;
328
329        if (!HasContribution(*vssRay))
330        {
331                mRayCaster->SheduleRayForDeletion(vssRay);
332                return false;
333        }
334       
335        if (GVS_DEBUG)
336        {
337                if (mVssRays.size() < 1000)
338                        mVssRays.push_back(new VssRay(*vssRay));
339        }
340
341        // add new ray to ray queue
342        mRayQueue.push(vssRay);
343
344        ++ mGvsStats.mTotalContribution;
345
346        if (mRayQueue.size() % 50000 == 0)
347                cout << "ray queue size: " << mRayQueue.size() << endl;
348
349        if (mGvsStats.mTotalContribution % 10000 == 0)
350                cout << "contri: " << mGvsStats.mTotalContribution << endl;
351
352        return true;
353}
354
355
356/** Creates 3 new vertices for triangle vertex with specified index.
357*/
358void GvsPreprocessor::CreateDisplacedVertices(VertexContainer &vertices,
359                                                                                          const Triangle3 &hitTriangle,
360                                                                                          const VssRay &ray,
361                                                                                          int index) const
362{
363        const int indexU = (index + 1) % 3;
364        const int indexL = (index - 1) % 3; //(index == 0) ? 2 : index - 1;
365
366        const Vector3 v  = hitTriangle.mVertices[index];
367        const Vector3 vu = hitTriangle.mVertices[indexU];
368        const Vector3 vl = hitTriangle.mVertices[indexL];
369
370        const Vector3 a = v - ray.GetOrigin();
371        const Vector3 b = vu - v;
372        const Vector3 c = v - vl;
373       
374        const float len = Magnitude(a);
375        //if (len > Limits::Small) a /= len;
376
377        const Vector3 dir1 = Normalize(CrossProd(a, b)); // N((pi-xp)×(pi+1 - pi));
378        const Vector3 dir2 = Normalize(CrossProd(a, c)); // N((pi-xp)×(pi - pi-1))
379        const Vector3 dir3 = DotProd(dir2, dir1) > 0 ? // N((pi-xp)×di,i-1+di,i+1×(pi-xp))
380                Normalize(dir2 + dir1) : Normalize(CrossProd(a, dir1) + CrossProd(dir2, a));
381#if 1
382
383        // compute the new three hit points
384        // pi, i + 1:  pi + e·|pi-xp|·di, j
385        const Vector3 pt1 = v + mEps * len * dir1;
386        // pi, i - 1:  pi + e·|pi-xp|·di, j
387    const Vector3 pt2 = v + mEps * len * dir2;
388        // pi, i:      pi + e·|pi-xp|·di, j
389        const Vector3 pt3 = v + mEps * len * dir3;
390
391#else
392
393        const Vector3 pt1 = v + mEps * dir1;
394        // pi, i - 1:  pi + e·|pi-xp|·di, j
395    const Vector3 pt2 = v + mEps * dir2;
396        // pi, i:      pi + e·|pi-xp|·di, j
397        const Vector3 pt3 = v + mEps * dir3;
398
399#endif
400
401        vertices.push_back(pt2);
402        vertices.push_back(pt3);
403        vertices.push_back(pt1);
404}
405
406
407void GvsPreprocessor::EnlargeTriangle(VertexContainer &vertices,
408                                                                          const Triangle3 &hitTriangle,
409                                                                          const VssRay &ray) const
410{
411        CreateDisplacedVertices(vertices, hitTriangle, ray, 0);
412        CreateDisplacedVertices(vertices, hitTriangle, ray, 1);
413        CreateDisplacedVertices(vertices, hitTriangle, ray, 2);
414}
415
416
417Vector3 GvsPreprocessor::CalcPredictedHitPoint(const VssRay &newRay,
418                                                                                           const Triangle3 &hitTriangle,
419                                                                                           const VssRay &oldRay) const
420{
421        // find the intersection of the plane induced by the
422        // hit triangle with the new ray
423        Plane3 plane(hitTriangle.GetNormal(), hitTriangle.mVertices[0]);
424
425        const Vector3 hitPt =
426                plane.FindIntersection(newRay.mTermination, newRay.mOrigin);
427       
428        return hitPt;
429}
430
431
432bool EqualVisibility(const VssRay &a, const VssRay &b)
433{
434        return a.mTerminationObject == b.mTerminationObject;
435}
436
437
438int GvsPreprocessor::SubdivideEdge(const Triangle3 &hitTriangle,
439                                                                   const Vector3 &p1,
440                                                                   const Vector3 &p2,
441                                                                   const VssRay &ray1,
442                                                                   const VssRay &ray2,
443                                                                   const VssRay &oldRay,
444                                                                   int level)
445{
446        int castRays = 0;
447       
448        if (EqualVisibility(ray1, ray2) || (SqrDistance(p1, p2) <= MIN_DIST)  || (level ++ > 10))
449        {
450                //if (SqrMagnitude(p1 - p2) <= MIN_DIST)
451                //      cerr << "min dist reached!! " << SqrMagnitude(p1 - p2) << " " << MIN_DIST << " " << ray1.mTerminationObject << " " << ray2.mTerminationObject << " level: " << level << " " << EqualVisibility(ray1, ray2) << endl;
452
453                return castRays;
454        }
455
456        // the new subdivision point
457        const Vector3 p = (p1 + p2) * 0.5f;
458       
459        SimpleRay sray(oldRay.mOrigin, p - oldRay.mOrigin, SamplingStrategy::GVS, 1.0f);
460
461        VssRay *newRay;
462
463        // cast single ray
464        castTimer.Entry();
465
466        // cast single ray to check if a triangle was missed. note: not efficient!!
467        newRay = mRayCaster->CastRay(sray, mViewCellsManager->GetViewSpaceBox(), !mPerViewCell);
468        ++ castRays;
469       
470        castTimer.Exit();
471
472        if (!newRay)
473        {
474                //cout<<"w";
475                return castRays;
476        }
477
478        //newRay->mFlags |= VssRay::BorderSample;
479
480        // cast reverse rays if necessary
481        if (USE_REVERSE_SAMPLING)
482        {
483                castRays += CheckDiscontinuity(*newRay, hitTriangle, oldRay);
484       
485                // also test reverse direction
486                if (newRay->mTerminationObject)
487                {
488                        TriangleIntersectable *tObj = static_cast<TriangleIntersectable *>(newRay->mTerminationObject);
489                        castRays += CheckDiscontinuity(oldRay, tObj->GetItem(), *newRay);
490                }
491        }
492
493        if (GVS_DEBUG)
494                currentViz.enlargedTriangle.push_back(p);
495
496        // subdivide further
497        castRays += SubdivideEdge(hitTriangle, p1, p, ray1, *newRay, oldRay, level);
498        castRays += SubdivideEdge(hitTriangle, p, p2, *newRay, ray2, oldRay, level);
499
500        rayTimer.Entry();
501
502        // new triangle discovered? => add new ray to queue
503        HandleRay(newRay);
504
505        rayTimer.Exit();
506
507        return castRays;
508}
509
510
511int GvsPreprocessor::AdaptiveBorderSampling(const VssRay &currentRay)
512{
513        Intersectable *tObj = currentRay.mTerminationObject;
514
515        // question matt: can this be possible?
516        if (!tObj) return 0;
517
518        Triangle3 hitTriangle;
519
520        // other types not implemented yet
521        if (tObj->Type() == Intersectable::TRIANGLE_INTERSECTABLE)
522                hitTriangle = static_cast<TriangleIntersectable *>(tObj)->GetItem();
523        else
524                cout << "border sampling: " << Intersectable::GetTypeName(tObj) << " not yet implemented" << endl;
525
526        //cout << "type: " << Intersectable::GetTypeName(tObj) << endl;
527        generationTimer.Entry();
528
529        static VertexContainer enlargedTriangle;
530        enlargedTriangle.clear();
531
532        /// create 3 new hit points for each vertex
533        EnlargeTriangle(enlargedTriangle, hitTriangle, currentRay);
534       
535        /// create rays from sample points and handle them
536        static SimpleRayContainer simpleRays;
537        simpleRays.clear();
538
539        VertexContainer::const_iterator vit, vit_end = enlargedTriangle.end();
540
541        for (vit = enlargedTriangle.begin(); vit != vit_end; ++ vit)
542        {
543                const Vector3 rayDir = (*vit) - currentRay.GetOrigin();
544
545                SimpleRay sr(currentRay.GetOrigin(), rayDir, SamplingStrategy::GVS, 1.0f);
546                simpleRays.AddRay(sr);
547        }
548
549        //if (rand() < RAND_MAX / 10)   cout << "o: " << currentRay.GetOrigin() << endl;
550
551        //////////
552        //-- fill up simple rays with random rays so we can cast 16 rays at a time
553
554        const int numRandomRays = 16 - (int)simpleRays.size();
555        GenerateRays(numRandomRays, *mDistribution, simpleRays, sInvalidSamples);
556
557        generationTimer.Exit();
558
559
560        /////////////////////
561        //-- cast rays to vertices and determine visibility
562       
563        static VssRayContainer vssRays;
564        vssRays.clear();
565
566        // for the original implementation we don't cast double rays
567        const bool castDoubleRays = !mPerViewCell;
568        // cannot prune invalid rays because we have to compare adjacent rays
569        const bool pruneInvalidRays = false;
570
571        //gvsTimer.Entry();
572        castTimer.Entry();
573
574        CastRays(simpleRays, vssRays, castDoubleRays, pruneInvalidRays);
575       
576        castTimer.Exit();
577        //gvsTimer.Exit();
578
579        const int numBorderSamples = (int)vssRays.size() - numRandomRays;
580        int castRays = (int)simpleRays.size();
581
582        if (GVS_DEBUG)
583        {
584                /*visTriangles.push_back(tObj);
585
586                for (int i = 0; i < 3; ++ i)
587                {
588                        Triangle3 t1(hitTriangle.mVertices[i], enlargedTriangle[3 * i + 0], enlargedTriangle[3 * i + 1]);
589                        Triangle3 t2(hitTriangle.mVertices[i], enlargedTriangle[3 * i + 1], enlargedTriangle[3 * i + 2]);
590
591                        TriangleIntersectable *to1 = new TriangleIntersectable(t1); to1->SetId(tObj->GetId());
592                        TriangleIntersectable *to2 = new TriangleIntersectable(t2); to2->SetId(tObj->GetId());
593
594                        visTriangles.push_back(to1);
595                        visTriangles.push_back(to2);
596                }
597               
598                // set flags
599                VssRayContainer::const_iterator rit, rit_end = vssRays.end();
600
601                for (rit = vssRays.begin(); rit != rit_end; ++ rit)
602                        (*rit)->mFlags |= VssRay::BorderSample;
603                */
604
605                // visualize enlarged triangles
606                currentViz.enlargedTriangle = enlargedTriangle;
607                currentViz.originalTriangle = static_cast<TriangleIntersectable *>(tObj);
608        }
609       
610        if (USE_SUBDIVISION)
611        {
612                //int oldsize = mTrianglePvs.size();
613
614                // recursivly subdivide each edge
615                for (int i = 0; i < numBorderSamples; ++ i)
616                {
617                        if (USE_REVERSE_SAMPLING)
618                        {
619                                castRays += CheckDiscontinuity(*vssRays[i], hitTriangle, currentRay);
620                               
621                                // also test reverse direction
622                                if (vssRays[i]->mTerminationObject)
623                                {
624                                        TriangleIntersectable *tObj = static_cast<TriangleIntersectable *>(vssRays[i]->mTerminationObject);
625                                        castRays += CheckDiscontinuity(currentRay, tObj->GetItem(), *vssRays[i]);
626                                }
627                        }
628
629                        castRays += SubdivideEdge(hitTriangle,
630                                                      enlargedTriangle[i],
631                                                                          enlargedTriangle[(i + 1) % numBorderSamples],
632                                                                          *vssRays[i],
633                                                                          *vssRays[(i + 1) % numBorderSamples],
634                                                                          currentRay,
635                                                                          0);
636                }
637
638                //int newsize = mTrianglePvs.size();if(newsize - oldsize) mcout<<"found with subdiv: " << newsize - oldsize << endl;
639        }
640
641        if (GVS_DEBUG)
642                visTriangles.push_back(currentViz);
643
644        mGvsStats.mBorderSamples += castRays;
645       
646        // handle cast rays
647        EnqueueRays(vssRays);
648
649        return castRays;
650}
651
652
653int GvsPreprocessor::AdaptiveBorderSamplingOpt(const VssRay &currentRay)
654{
655        Intersectable *tObj = currentRay.mTerminationObject;
656
657        // question matt: can this be possible?
658        if (!tObj) return 0;
659
660        Triangle3 hitTriangle;
661
662        // other types not implemented yet
663        if (tObj->Type() == Intersectable::TRIANGLE_INTERSECTABLE)
664                hitTriangle = static_cast<TriangleIntersectable *>(tObj)->GetItem();
665        else
666                cout << "error: gvs sampling for " << Intersectable::GetTypeName(tObj) << " not implemented" << endl;
667
668       
669        static SimpleRayContainer simpleRays;
670        simpleRays.clear();
671       
672
673        generationTimer.Entry();
674#if 0
675        static VertexContainer enlargedTriangle;
676        enlargedTriangle.clear();
677       
678        /// create 3 new hit points for each vertex
679        EnlargeTriangle(enlargedTriangle, hitTriangle, currentRay);
680
681        static VertexContainer subdivEnlargedTri;       
682        subdivEnlargedTri.clear();
683
684        for (size_t i = 0; i < enlargedTriangle.size(); ++ i)
685        {
686                const Vector3 p1 = enlargedTriangle[i];
687                const Vector3 p2 = enlargedTriangle[(i + 1) % enlargedTriangle.size()];
688
689                // the new subdivision point
690                const Vector3 p = (p1 + p2) * 0.5f;
691
692                subdivEnlargedTri.push_back(p1);
693                subdivEnlargedTri.push_back(p);
694        }
695
696        VertexContainer::const_iterator vit, vit_end = subdivEnlargedTri.end();
697
698        for (vit = subdivEnlargedTri.begin(); vit != vit_end; ++ vit)
699        {
700                const Vector3 rayDir = (*vit) - currentRay.GetOrigin();
701
702                SimpleRay sr(currentRay.GetOrigin(), rayDir, SamplingStrategy::GVS, 1.0f);
703                simpleRays.AddRay(sr); 
704        }
705
706#endif 
707
708        //////////
709        //-- fill up simple rays with random rays so we can cast 16 rays at a time
710
711        const int numRandomRays = 128;//16 - ((int)simpleRays.size() % 16);
712
713        GenerateImportanceSamples(currentRay, hitTriangle, numRandomRays, simpleRays);
714       
715        generationTimer.Exit();
716
717
718
719        /////////////////////
720        //-- cast rays to vertices and determine visibility
721       
722        static VssRayContainer vssRays;
723        vssRays.clear();
724
725        // for the original implementation we don't cast double rays
726        const bool castDoubleRays = !mPerViewCell;
727        // cannot prune invalid rays because we have to compare adjacent rays
728        const bool pruneInvalidRays = false;
729
730        //if ((simpleRays.size() % 16) != 0) cerr << "ray size not aligned" << endl;
731       
732        //gvsTimer.Entry();
733        castTimer.Entry();
734
735        CastRays(simpleRays, vssRays, castDoubleRays, pruneInvalidRays);
736       
737        castTimer.Exit();
738        //gvsTimer.Exit();
739
740        const int numBorderSamples = (int)vssRays.size() - numRandomRays;
741        int castRays = (int)simpleRays.size();
742
743        mGvsStats.mBorderSamples += castRays;
744
745        if (USE_REVERSE_SAMPLING)
746        {
747                static SimpleRayContainer simpleRays;
748                simpleRays.clear();
749
750                SimpleRay simpleRay;
751                for (size_t i = 0; i < vssRays.size(); ++ i)
752                {
753                        // check for discontinuity and cast reverse rays if necessary
754                        //castRays += CheckDiscontinuity2(*vssRays[i], hitTriangle, currentRay);
755                        if (CheckDiscontinuity2(*vssRays[i], hitTriangle, currentRay, simpleRay))
756                                simpleRays.push_back(simpleRay);
757                }
758               
759                static VssRayContainer reverseRays;
760                reverseRays.clear();
761                CastRays(simpleRays, reverseRays, false, false);
762
763                mGvsStats.mReverseSamples += (int)reverseRays.size();
764       
765                castRays +=(int)reverseRays.size();
766        }
767
768#if 0
769
770    // recursivly subdivide each edge
771        for (int i = 0; i < numBorderSamples; ++ i)
772        {
773                castRays += SubdivideEdge(hitTriangle,
774                                                                  subdivEnlargedTri[i],
775                                                                  subdivEnlargedTri[(i + 1) % numBorderSamples],
776                                                                  *vssRays[i],
777                                                                  *vssRays[(i + 1) % numBorderSamples],
778                                                                  currentRay, 0);
779        }
780       
781#endif
782
783        // handle cast rays
784        EnqueueRays(vssRays);
785
786        return castRays;
787}
788
789
790bool GvsPreprocessor::GetPassingPoint(const VssRay &currentRay,
791                                                                          const Triangle3 &occluder,
792                                                                          const VssRay &oldRay,
793                                                                          Vector3 &newPoint) const
794{
795        /////////////////////////
796        //-- We interserct the plane p = (xp, hit(x), hit(xold))
797        //-- with the newly found occluder (xold is the previous ray from
798        //-- which x was generated). On the intersecting line, we select a point
799        //-- pnew which lies just outside of the new triangle so the ray
800        //-- just passes through the gap
801
802        const Plane3 plane(currentRay.GetOrigin(),
803                                           currentRay.GetTermination(),
804                                           oldRay.GetTermination());
805       
806        Vector3 pt1, pt2;
807
808        if (!occluder.GetPlaneIntersection(plane, pt1, pt2))
809                return false;
810
811        // get the intersection point on the old ray
812        const Plane3 triPlane(occluder.GetNormal(), occluder.mVertices[0]);
813
814        const float t = triPlane.FindT(oldRay.mOrigin, oldRay.mTermination);
815        const Vector3 pt3 = oldRay.mOrigin + t * (oldRay.mTermination - oldRay.mOrigin);
816
817        // evaluate new hitpoint just outside the triangle
818        // the point is chosen to be on the side closer to the original ray
819        if (SqrDistance(pt1, pt3) < SqrDistance(pt2, pt3))
820                newPoint = pt1 + mEps * (pt1 - pt2);
821        else
822                newPoint = pt2 + mEps * (pt2 - pt1);
823
824        //cout << "passing point: " << newPoint << endl << endl;
825        return true;
826}
827
828
829bool GvsPreprocessor::ComputeReverseRay(const VssRay &currentRay,
830                                                                                const Triangle3 &hitTriangle,
831                                                                                const VssRay &oldRay,
832                                                                                SimpleRay &reverseRay)
833{
834        // optain triangle that occludes the path to the currently processed mesh
835        Triangle3 occluder;
836        Intersectable *tObj = currentRay.mTerminationObject;
837
838        // q: why can this happen?
839        if (!tObj)
840                return false;
841
842        // other types not implemented yet
843        if (tObj->Type() == Intersectable::TRIANGLE_INTERSECTABLE)
844                occluder = static_cast<TriangleIntersectable *>(tObj)->GetItem();
845        else
846                cout << "reverse sampling: " << tObj->Type() << " not yet implemented" << endl;
847       
848        // get a point which is passing just outside of the occluder
849    Vector3 newPoint;
850
851        // q: why is there sometimes no intersecton found?
852        if (!GetPassingPoint(currentRay, occluder, oldRay, newPoint))
853                return false;
854
855        const Vector3 predicted = CalcPredictedHitPoint(currentRay, hitTriangle, oldRay);
856
857        //////////////
858        //-- Construct the mutated ray with xnew,
859        //-- dir = predicted(x)- pnew as direction vector
860
861        Vector3 newDir = newDir = predicted - newPoint;
862
863        // take xnew, p = intersect(viewcell, line(pnew, predicted(x)) as origin ?
864        // difficult to say!!
865        const float offset = 0.5f;
866        Vector3 newOrigin = newPoint - newDir * offset;
867       
868
869        //////////////
870        //-- for per view cell sampling, we must intersect
871        //-- with the current view cell
872
873    if (mPerViewCell)
874        {
875                if (!IntersectsViewCell(newOrigin, -newDir))
876                        return false;
877        }
878
879        reverseRay = SimpleRay(newOrigin, newDir, SamplingStrategy::GVS, 1.0f);
880
881        return true;
882}
883
884
885bool GvsPreprocessor::IntersectsViewCell(Vector3 &origin, const Vector3 &dir) const
886{
887        AxisAlignedBox3 box = mCurrentViewCell->GetBox();
888
889        if (box.IsInside(origin))
890        {
891                //cout<<"z";
892                return true;
893        }
894       
895        // intersect ray with view cell and
896        SimpleRay ray(origin, dir, 0, 1.0f);
897
898        float tnear, tfar;
899       
900        // check if ray intersects view cell
901        if (!mCurrentViewCell->GetBox().Intersects(ray, tnear, tfar))
902                return false;
903
904        origin = origin + dir * tnear;
905        //cout<<"r";
906        return true;
907}
908
909
910int GvsPreprocessor::CastInitialSamples(int numSamples)
911{       
912        initialTimer.Entry();
913
914        // generate simple rays
915        static SimpleRayContainer simpleRays;
916        simpleRays.clear();
917
918        generationTimer.Entry();
919
920        if (mUseProbablyVisibleSampling)// && (RandomValue(0, 1) > 0.5))
921        {
922                ProbablyVisibleDistribution distr(*this, mCurrentViewCell, &mProbablyVisibleTriangles);
923                GenerateRays(numSamples, distr, simpleRays, sInvalidSamples);
924        }
925        else
926        {
927                GenerateRays(mUseDeterministicGvs ?
928                        numSamples : numSamples / 16, *mDistribution, simpleRays, sInvalidSamples);
929        }
930
931        generationTimer.Exit();
932
933        // cast generated samples
934        static VssRayContainer vssRays;
935        vssRays.clear();
936
937        castTimer.Entry();
938
939        if (mUseDeterministicGvs)
940        {
941                const bool castDoubleRays = !mPerViewCell;
942                const bool pruneInvalidRays = false;
943                //const bool pruneInvalidRays = true;
944
945                CastRays(simpleRays, vssRays, castDoubleRays, pruneInvalidRays);
946        }
947        else
948        {
949                if (0)
950                        // casting bundles of 4 rays generated from the simple rays
951                        CastRayBundles4(simpleRays, vssRays);
952                else
953                        CastRayBundles16(simpleRays, vssRays, mInitialJitter);
954        }
955
956        castTimer.Exit();
957//for (int i = 0; i< simpleRays.size();++ i)
958//      if (rand() < RAND_MAX / 10)     cout << "o: " << simpleRays[i].mOrigin << endl;
959
960        // add to ray queue
961        EnqueueRays(vssRays);
962
963        initialTimer.Exit();
964
965        return (int)vssRays.size();
966}
967
968
969int GvsPreprocessor::EnqueueRays(const VssRayContainer &samples)
970{
971        int success = 0;
972        rayTimer.Entry();
973
974        // add samples to ray queue
975        VssRayContainer::const_iterator vit, vit_end = samples.end();
976        for (vit = samples.begin(); vit != vit_end; ++ vit)
977        {
978                VssRay *ray = *vit;
979                if (HandleRay(ray)) ++ success;
980        }
981
982        rayTimer.Exit();
983
984        return success;
985}
986
987
988int GvsPreprocessor::ProcessQueue()
989{
990        gvsTimer.Entry();
991
992        ++ mGvsStats.mGvsRuns;
993
994        int castSamples = 0;
995
996        while (!mRayQueue.empty())
997        {
998                size_t oldContri = mTrianglePvs.size();
999                // handle next ray
1000                VssRay *ray = mRayQueue.top();
1001                mRayQueue.pop();
1002               
1003                if (mUseDeterministicGvs)
1004                {
1005                        castSamples += AdaptiveBorderSampling(*ray);
1006                        mRayCaster->SheduleRayForDeletion(ray);
1007                }
1008                else
1009                {
1010                        castSamples += AdaptiveBorderSamplingOpt(*ray);
1011
1012                        // ray still gives contribution => enqeue again
1013                        if ((mTrianglePvs.size() - oldContri) > 0)
1014                        {
1015                                mRayQueue.push(ray);
1016                        }
1017                        else
1018                        {
1019                                mRayCaster->SheduleRayForDeletion(ray);
1020                        }
1021                }
1022
1023
1024        }
1025
1026        gvsTimer.Exit();
1027
1028        return castSamples;
1029}
1030
1031
1032void ExportVssRays(Exporter *exporter, const VssRayContainer &vssRays)
1033{
1034        VssRayContainer vcRays, vcRays2, vcRays3;
1035
1036        VssRayContainer::const_iterator rit, rit_end = vssRays.end();
1037
1038        // prepare some rays for output
1039        for (rit = vssRays.begin(); rit != rit_end; ++ rit)
1040        {
1041                //const float p = RandomValue(0.0f, (float)vssRays.size());
1042
1043                if (1)//(p < raysOut)
1044                {
1045                        if ((*rit)->mFlags & VssRay::BorderSample)
1046                                vcRays.push_back(*rit);
1047                        else if ((*rit)->mFlags & VssRay::ReverseSample)
1048                                vcRays2.push_back(*rit);
1049                        else
1050                                vcRays3.push_back(*rit);
1051                }
1052        }
1053
1054        exporter->ExportRays(vcRays, RgbColor(1, 0, 0));
1055        exporter->ExportRays(vcRays2, RgbColor(0, 1, 0));
1056        exporter->ExportRays(vcRays3, RgbColor(1, 1, 1));
1057}
1058
1059
1060void GvsPreprocessor::VisualizeViewCell(const ObjectContainer &objects)
1061{
1062    Intersectable::NewMail();
1063        Material m;
1064       
1065        char str[64]; sprintf_s(str, "pass%06d.wrl", mProcessedViewCells);
1066
1067        Exporter *exporter = Exporter::GetExporter(str);
1068        if (!exporter)
1069                return;
1070
1071        ObjectContainer::const_iterator oit, oit_end = objects.end();
1072
1073        for (oit = objects.begin(); oit != oit_end; ++ oit)
1074        {
1075                Intersectable *intersect = *oit;
1076               
1077                m = RandomMaterial();
1078                exporter->SetForcedMaterial(m);
1079                exporter->ExportIntersectable(intersect);
1080        }
1081
1082        cout << "vssrays: " << (int)mVssRays.size() << endl;
1083        ExportVssRays(exporter, mVssRays);
1084
1085
1086        /////////////////
1087        //-- export view cell geometry
1088
1089        //exporter->SetWireframe();
1090
1091        m.mDiffuseColor = RgbColor(0, 1, 0);
1092        exporter->SetForcedMaterial(m);
1093
1094        AxisAlignedBox3 bbox = mCurrentViewCell->GetMesh()->mBox;
1095        exporter->ExportBox(bbox);
1096        //exporter->SetFilled();
1097
1098        delete exporter;
1099}
1100
1101
1102void GvsPreprocessor::VisualizeViewCells()
1103{
1104        char str[64]; sprintf_s(str, "tmp/pass%06d_%04d-", mProcessedViewCells, mPass);
1105                       
1106        // visualization
1107        if (mGvsStats.mPassContribution > 0)
1108        {
1109                const bool exportRays = true;
1110                const bool exportPvs = true;
1111
1112                mViewCellsManager->ExportSingleViewCells(mObjects,
1113                                                                                                 10,
1114                                                                                                 false,
1115                                                                                                 exportPvs,
1116                                                                                                 exportRays,
1117                                                                                                 1000,
1118                                                                                                 str);
1119        }
1120
1121        // remove pass samples
1122        ViewCellContainer::const_iterator vit, vit_end = mViewCellsManager->GetViewCells().end();
1123
1124        for (vit = mViewCellsManager->GetViewCells().begin(); vit != vit_end; ++ vit)
1125        {
1126                (*vit)->DelRayRefs();
1127        }
1128}
1129
1130
1131void GvsPreprocessor::CompileViewCellsFromPointList()
1132{
1133        ViewCell::NewMail();
1134
1135        // Receive list of view cells from view cells manager
1136        ViewCellPointsList *vcPoints = mViewCellsManager->GetViewCellPointsList();
1137
1138        vector<ViewCellPoints *>::const_iterator vit, vit_end = vcPoints->end();
1139
1140        for (vit = vcPoints->begin(); vit != vit_end; ++ vit)
1141        {
1142                ViewCellPoints *vp = *vit;
1143
1144                if (vp->first) // view cell  specified
1145                {
1146                        ViewCell *viewCell = vp->first;
1147
1148                        if (!viewCell->Mailed())
1149                        {
1150                                viewCell->Mail();
1151                                mViewCells.push_back(viewCell);
1152
1153                                if (mViewCells.size() >= mMaxViewCells)
1154                                        break;
1155                        }
1156                }
1157                else // no view cell specified => compute view cells using view point
1158                {
1159                        SimpleRayContainer::const_iterator rit, rit_end = vp->second.end();
1160
1161                        for (rit = vp->second.begin(); rit != rit_end; ++ rit)
1162                        {
1163                                SimpleRay ray = *rit;
1164
1165                                ViewCell *viewCell = mViewCellsManager->GetViewCell(ray.mOrigin);
1166
1167                                if (viewCell && !viewCell->Mailed())
1168                                {
1169                                        viewCell->Mail();
1170                                        mViewCells.push_back(viewCell);
1171
1172                                        if (mViewCells.size() >= mMaxViewCells)
1173                                                break;
1174                                }
1175                        }
1176                }
1177        }
1178}
1179
1180
1181void GvsPreprocessor::CompileViewCellsList()
1182{
1183        if (!mViewCellsManager->GetViewCellPointsList()->empty())
1184        {
1185                cout << "processing view point list" << endl;
1186                CompileViewCellsFromPointList();
1187        }
1188        else
1189        {
1190                ViewCell::NewMail();
1191
1192                while ((int)mViewCells.size() < mMaxViewCells)
1193                {
1194                        const int tries = 10000;
1195                        int i = 0;
1196
1197                        for (i = 0; i < tries; ++ i)
1198                        {
1199                                const int idx = (int)RandomValue(0.0f, (float)mViewCellsManager->GetNumViewCells() - 0.5f);
1200
1201                                ViewCell *viewCell = mViewCellsManager->GetViewCell(idx);
1202
1203                                if (!viewCell->Mailed())
1204                                {
1205                                        viewCell->Mail();
1206                                        break;
1207                                }
1208
1209                                mViewCells.push_back(viewCell);
1210                        }
1211
1212                        if (i == tries)
1213                        {
1214                                cerr << "big error! no view cell found" << endl;
1215                                break;
1216                        }
1217                }
1218        }
1219
1220        cout << "\ncomputing list of " << mViewCells.size() << " view cells" << endl;
1221}
1222
1223
1224void GvsPreprocessor::ComputeViewCellGeometryIntersection()
1225{
1226        intersectionTimer.Entry();
1227
1228        AxisAlignedBox3 box = mCurrentViewCell->GetBox();
1229
1230        int oldTrianglePvs = (int)mTrianglePvs.size();
1231        int newKdObj = 0;
1232
1233        // compute pvs kd nodes that intersect view cell
1234        ObjectContainer kdobjects;
1235
1236        // find intersecting objects not in pvs (i.e., only unmailed)
1237        mKdTree->CollectKdObjects(box, kdobjects);
1238
1239        ObjectContainer pvsKdObjects;
1240
1241        ObjectContainer::const_iterator oit, oit_end = kdobjects.end();
1242
1243        for (oit = kdobjects.begin(); oit != oit_end; ++ oit)
1244        {
1245        KdIntersectable *kdInt = static_cast<KdIntersectable *>(*oit);
1246       
1247                myobjects.clear();
1248                mKdTree->CollectObjectsWithDublicates(kdInt->GetItem(), myobjects);
1249
1250                // account for kd object pvs
1251                bool addkdobj = false;
1252
1253                ObjectContainer::const_iterator oit, oit_end = myobjects.end();
1254
1255                for (oit = myobjects.begin(); oit != oit_end; ++ oit)
1256                {
1257                        TriangleIntersectable *triObj = static_cast<TriangleIntersectable *>(*oit);
1258             
1259                        // test if the triangle itself intersects
1260                        if (box.Intersects(triObj->GetItem()))
1261                                addkdobj = AddTriangleObject(triObj);
1262                }
1263               
1264                // add node to kd pvs
1265                if (1 && addkdobj)
1266                {
1267                        ++ newKdObj;
1268                        mCurrentViewCell->GetPvs().AddSampleDirty(kdInt, 1.0f);
1269                        //mCurrentViewCell->GetPvs().AddSampleDirtyCheck(kdInt, 1.0f);
1270                        if (QT_VISUALIZATION_SHOWN) UpdateStatsForVisualization(kdInt);
1271                }
1272        }
1273
1274        cout << "added " << (int)mTrianglePvs.size() - oldTrianglePvs << " triangles (" << newKdObj << " objects) by intersection" << endl;
1275
1276        intersectionTimer.Exit();
1277}
1278
1279
1280void GvsPreprocessor::PerViewCellComputation()
1281{
1282        while (mCurrentViewCell = NextViewCell())
1283        {
1284                preparationTimer.Entry();
1285
1286                // hack: reset counters (could be done with a mailing-like approach
1287                ObjectContainer::const_iterator oit, oit_end = mObjects.end();
1288                for (oit = mObjects.begin(); oit != oit_end; ++ oit)
1289                        (*oit)->mCounter = NOT_ACCOUNTED_OBJECT;
1290
1291                preparationTimer.Exit();
1292
1293                mDistribution->SetViewCell(mCurrentViewCell);
1294
1295        ComputeViewCell();
1296        }
1297}
1298
1299
1300void GvsPreprocessor::PerViewCellComputation2()
1301{
1302        while (1)
1303        {
1304                if (!mRendererWidget)
1305                        continue;
1306
1307        mCurrentViewCell = mViewCellsManager->GetViewCell(mRendererWidget->GetViewPoint());
1308
1309                // no valid view cell or view cell already computed
1310                if (!mCurrentViewCell || !mCurrentViewCell->GetPvs().Empty() || !mRendererWidget->mComputeGVS)
1311                        continue;
1312
1313                mRendererWidget->mComputeGVS = false;
1314                // hack: reset counters
1315                ObjectContainer::const_iterator oit, oit_end = mObjects.end();
1316
1317                for (oit = mObjects.begin(); oit != oit_end; ++ oit)
1318                        (*oit)->mCounter = NOT_ACCOUNTED_OBJECT;
1319
1320                ComputeViewCell();
1321                ++ mProcessedViewCells;
1322        }
1323}
1324
1325
1326void GvsPreprocessor::StorePvs(const ObjectContainer &objectPvs)
1327{
1328        ObjectContainer::const_iterator oit, oit_end = objectPvs.end();
1329
1330        for (oit = objectPvs.begin(); oit != oit_end; ++ oit)
1331        {
1332                mCurrentViewCell->GetPvs().AddSample(*oit, 1);
1333        }
1334}
1335
1336
1337void GvsPreprocessor::UpdatePvs(ViewCell *currentViewCell)
1338{
1339        ObjectPvs newPvs;
1340        BvhLeaf::NewMail();
1341
1342        ObjectPvsIterator pit = currentViewCell->GetPvs().GetIterator();
1343
1344        // output PVS of view cell
1345        while (pit.HasMoreEntries())
1346        {               
1347                Intersectable *intersect = pit.Next();
1348
1349                BvhLeaf *bv = intersect->mBvhLeaf;
1350
1351                if (!bv || bv->Mailed())
1352                        continue;
1353               
1354                bv->Mail();
1355
1356                //m.mDiffuseColor = RgbColor(1, 0, 0);
1357                newPvs.AddSampleDirty(bv, 1.0f);
1358        }
1359
1360        newPvs.SimpleSort();
1361
1362        currentViewCell->SetPvs(newPvs);
1363}
1364
1365 
1366void GvsPreprocessor::GetObjectPvs(ObjectContainer &objectPvs) const
1367{
1368        BvhLeaf::NewMail();
1369
1370        objectPvs.reserve((int)mTrianglePvs.size());
1371
1372        ObjectContainer::const_iterator oit, oit_end = mTrianglePvs.end();
1373
1374        for (oit = mTrianglePvs.begin(); oit != oit_end; ++ oit)
1375        {
1376                Intersectable *intersect = *oit;
1377       
1378                BvhLeaf *bv = intersect->mBvhLeaf;
1379
1380                // hack: reset counter
1381                (*oit)->mCounter = NOT_ACCOUNTED_OBJECT;
1382
1383                if (!bv || bv->Mailed())
1384                        continue;
1385
1386                bv->Mail();
1387                objectPvs.push_back(bv);
1388        }
1389}
1390
1391
1392void GvsPreprocessor::GlobalComputation()
1393{
1394        int passSamples = 0;
1395        int randomSamples = 0;
1396        int gvsSamples = 0;
1397        int oldContribution = 0;
1398
1399        while (mGvsStats.mTotalSamples < mTotalSamples)
1400        {
1401                mRayCaster->InitPass();
1402
1403                int newRandomSamples, newGvsSamples;
1404
1405                // Ray queue empty =>
1406                // cast a number of uniform samples to fill ray queue
1407                newRandomSamples = CastInitialSamples(mInitialSamples);
1408               
1409                if (!mOnlyRandomSampling)
1410                        newGvsSamples = ProcessQueue();
1411
1412                passSamples += newRandomSamples + newGvsSamples;
1413                mGvsStats.mTotalSamples += newRandomSamples + newGvsSamples;
1414
1415                mGvsStats.mRandomSamples += newRandomSamples;
1416                mGvsStats.mGvsSamples += newGvsSamples;
1417
1418
1419                if (passSamples % (mGvsSamplesPerPass + 1) == mGvsSamplesPerPass)
1420                {
1421                        ++ mPass;
1422
1423                        mGvsStats.mPassContribution = mGvsStats.mTotalContribution - oldContribution;
1424                       
1425                        ////////
1426                        //-- stats
1427
1428                        //cout << "\nPass " << mPass << " #samples: " << mGvsStats.mTotalSamples << " of " << mTotalSamples << endl;
1429                        mGvsStats.mPass = mPass;
1430                        mGvsStats.Stop();
1431                        mGvsStats.Print(mGvsStatsStream);
1432
1433                        // reset
1434                        oldContribution = mGvsStats.mTotalContribution;
1435                        mGvsStats.mPassContribution = 0;
1436                        passSamples = 0;
1437
1438                        if (GVS_DEBUG)
1439                                VisualizeViewCells();
1440                }
1441        }
1442}
1443
1444
1445bool GvsPreprocessor::ComputeVisibility()
1446{
1447        cout << "Gvs Preprocessor started\n" << flush;
1448        const long startTime = GetTime();
1449
1450        //Randomize(0);
1451        mGvsStats.Reset();
1452        mGvsStats.Start();
1453
1454        if (!mLoadViewCells)
1455        {       
1456                /// construct the view cells from the scratch
1457                ConstructViewCells();
1458                // reset pvs already gathered during view cells construction
1459                mViewCellsManager->ResetPvs();
1460                cout << "finished view cell construction" << endl;
1461        }
1462
1463        if (mPerViewCell)
1464        {
1465#if 1
1466                // provide list of view cells to compute
1467                CompileViewCellsList();
1468                // start per view cell gvs
1469                PerViewCellComputation();
1470#else
1471                PerViewCellComputation2();
1472
1473#endif
1474        }
1475        else
1476        {
1477                GlobalComputation();
1478        }
1479
1480        cout << "cast " << mGvsStats.mTotalSamples / (1e3f * TimeDiff(startTime, GetTime())) << "M single rays/s" << endl;
1481
1482        if (GVS_DEBUG)
1483        {
1484                CLEAR_CONTAINER(mVssRays);
1485        }
1486       
1487        // export the preprocessed information to a file
1488        if (0 && mExportVisibility)
1489                ExportPreprocessedData(mVisibilityFileName);
1490       
1491        // compute the pixel error of this visibility solution
1492        if (0 && mEvaluatePixelError)
1493                ComputeRenderError();
1494
1495        return true;
1496}
1497
1498
1499void GvsPreprocessor::DeterminePvsObjects(VssRayContainer &rays)
1500{
1501        // store triangle directly
1502        mViewCellsManager->DeterminePvsObjects(rays, true);
1503}
1504
1505
1506void GvsStatistics::Print(ostream &app) const
1507{
1508        app << "#ViewCells\n" << mViewCells << endl;
1509        app << "#Id\n" << mViewCellId << endl;
1510        app << "#Time\n" << mTimePerViewCell << endl;;
1511        app << "#TriPvs\n" << mTrianglePvs << endl;
1512        app << "#ObjectPvs\n" << mPerViewCellPvs << endl;
1513        app << "#UpdatedPvs\n" << (int)mPvsCost << endl;
1514
1515        app << "#PerViewCellSamples\n" << mPerViewCellSamples << endl;
1516        app << "#RaysPerSec\n" << RaysPerSec() << endl;
1517       
1518        app << "#TotalPvs\n" << mTotalPvs << endl;
1519        app << "#TotalTime\n" << mTotalTime << endl;
1520        app << "#TotalSamples\n" << mTotalSamples << endl;
1521
1522        app << "#AvgPvs: " << mPvsCost / mViewCells << endl;
1523        app << "#TotalTrianglePvs\n" << mTotalTrianglePvs << endl;
1524       
1525        if (0)
1526        {
1527                app << "#ReverseSamples\n" << mReverseSamples << endl;
1528                app << "#BorderSamples\n" << mBorderSamples << endl;
1529
1530                app << "#Pass\n" << mPass << endl;
1531                app << "#ScDiff\n" << mPassContribution << endl;
1532                app << "#GvsRuns\n" << mGvsRuns << endl;       
1533
1534                app     << "#SamplesContri\n" << mTotalContribution << endl;
1535        }
1536
1537        app << endl;
1538}
1539
1540
1541void GvsPreprocessor::ComputeViewCell()
1542{
1543        const long startTime = GetTime();
1544
1545        // initialise
1546        mGvsStats.mPerViewCellSamples = 0;
1547        mGvsStats.mRandomSamples = 0;
1548        mGvsStats.mGvsSamples = 0;
1549        mUseProbablyVisibleSampling = false;
1550        //renderer->mPvsStat.currentPass = 0;
1551
1552        mPass = 0;
1553
1554        int oldContribution = 0;
1555        int passSamples = 0;
1556
1557        for (int i = 0; i < 2; ++ i)
1558                mGenericStats[i] = 0;
1559
1560        mTrianglePvs.clear();
1561
1562
1563        // hack: take bounding box of view cell as view cell bounds
1564        mCurrentViewCell->GetMesh()->ComputeBoundingBox();
1565
1566        // use mailing for kd node pvs
1567        KdNode::NewMail();
1568
1569        cout << "\n***********************\n"
1570                << "computing view cell " << mProcessedViewCells
1571                << " (id: " << mCurrentViewCell->GetId() << ")" << endl;
1572        cout << "bb: " << mCurrentViewCell->GetBox() << endl;
1573
1574        mainLoopTimer.Entry();
1575
1576        while (1)
1577        {
1578                sInvalidSamples = 0;
1579                int oldPvsSize = mCurrentViewCell->GetPvs().GetSize();
1580               
1581                mRayCaster->InitPass();
1582
1583                int newRandomSamples, newGvsSamples, newSamples;
1584
1585                // Ray queue empty =>
1586                // cast a number of uniform samples to fill ray queue
1587                newRandomSamples = CastInitialSamples(mInitialSamples);
1588               
1589                //if (!mUseProbablyVisibleSampling && !mOnlyRandomSampling)
1590                if (!mOnlyRandomSampling)
1591                        newGvsSamples = ProcessQueue();
1592
1593                newSamples = newRandomSamples + newGvsSamples;
1594
1595                passSamples += newSamples;
1596                mGvsStats.mPerViewCellSamples += newSamples;
1597
1598                mGvsStats.mRandomSamples += newRandomSamples;
1599                mGvsStats.mGvsSamples += newGvsSamples;
1600
1601
1602                if (passSamples >= mGvsSamplesPerPass)
1603                {
1604                        if (0 && GVS_DEBUG) VisualizeViewCell(mTrianglePvs);
1605
1606                        //const bool convertPerPass = true;
1607                        const bool convertPerPass = false;
1608
1609                        if (convertPerPass)
1610                        {       
1611                                int newObjects = ConvertObjectPvs();
1612
1613                                cout << "added " << newObjects << " triangles to triangle pvs" << endl;
1614                                mGvsStats.mTotalContribution += newObjects;
1615                        }
1616
1617                        ++ mPass;
1618                        mGvsStats.mPassContribution = mGvsStats.mTotalContribution - oldContribution;
1619
1620                        if (!mUseDeterministicGvs && (mGvsStats.mPassContribution < 2000))
1621                                mUseProbablyVisibleSampling = true;
1622                       
1623                        if (mUseProbablyVisibleSampling)
1624                                PrepareProbablyVisibleSampling();
1625
1626
1627                        ////////
1628                        //-- stats
1629
1630                        mGvsStats.mPass = mPass;
1631
1632                        cout << "\nPass " << mPass << " #samples: " << mGvsStats.mPerViewCellSamples << endl;
1633                        cout << "contribution=" << mGvsStats.mPassContribution << " (of " << mMinContribution << ")" << endl;
1634                        cout << "obj contri=" << mCurrentViewCell->GetPvs().GetSize() - oldPvsSize << " (of " << mMinContribution << ")" << endl;
1635                        cout << "successful reverse rays=" << successfullReverseRays << endl;
1636                        cout << "cast reverse rays=" << mGvsStats.mReverseSamples << endl;
1637
1638                        // termination criterium
1639                        if (mGvsStats.mPassContribution < mMinContribution)
1640                                break;
1641
1642                        // reset stats
1643                        oldContribution = mGvsStats.mTotalContribution;
1644                        mGvsStats.mPassContribution = 0;
1645                        passSamples = 0;
1646                        successfullReverseRays = 0;
1647                       
1648
1649                        // compute the pixel error of this visibility solution
1650                        if (1 && mEvaluatePixelError)
1651                                ComputeRenderError();
1652                }
1653        }
1654       
1655        mainLoopTimer.Exit();
1656
1657        // at last compute objects that directly intersect view cell
1658        ComputeViewCellGeometryIntersection();
1659               
1660        // compute the pixel error of this visibility solution
1661        if (mEvaluatePixelError)
1662                ComputeRenderError();
1663
1664        ////////
1665        //-- stats
1666
1667        // timing
1668        mGvsStats.mTimePerViewCell = TimeDiff(startTime, GetTime()) * 1e-3f;
1669
1670        ComputeStats();
1671}
1672
1673
1674void GvsPreprocessor::ComputeStats()
1675{
1676        // compute pvs size using larger (bvh objects)
1677        // note: for kd pvs we had to do this during gvs computation
1678        if (!mUseKdPvs)
1679        {
1680                ObjectContainer objectPvs;
1681
1682                // optain object pvs
1683                GetObjectPvs(objectPvs);
1684
1685                // add pvs
1686                ObjectContainer::const_iterator it, it_end = objectPvs.end();
1687
1688                for (it = objectPvs.begin(); it != it_end; ++ it)
1689                {
1690                        mCurrentViewCell->GetPvs().AddSampleDirty(*it, 1.0f);
1691                }
1692
1693                cout << "triangle pvs of " << (int)mTrianglePvs.size()
1694                        << " was converted to object pvs of " << (int)objectPvs.size() << endl;
1695
1696                mGvsStats.mPerViewCellPvs = (int)objectPvs.size();
1697                mGvsStats.mPvsCost = 0; // todo objectPvs.EvalPvsCost();
1698        }
1699        else
1700        {
1701                 if (0) // compute pvs kd nodes that intersect view cell
1702                 {
1703                         ObjectContainer mykdobjects;
1704
1705                         // extract kd pvs
1706                         ObjectContainer::const_iterator it, it_end = mTrianglePvs.end();
1707
1708                         for (it = mTrianglePvs.begin(); it != it_end; ++ it)
1709                         {
1710                                 Intersectable *obj = *it;
1711                                 // find intersecting objects not yet in pvs (i.e., only unmailed)
1712                                 mKdTree->CollectKdObjects(obj->GetBox(), mykdobjects);
1713                         }
1714
1715                         // add pvs
1716                         it_end = mykdobjects.end();
1717
1718                         for (it = mykdobjects.begin(); it != it_end; ++ it)
1719                         {
1720                                 mCurrentViewCell->GetPvs().AddSampleDirty(*it, 1.0f);
1721                         }
1722                 
1723                         cout << "added " << (int)mykdobjects.size() << " new objects " << endl;
1724                 }
1725
1726                 mGvsStats.mPerViewCellPvs = mCurrentViewCell->GetPvs().GetSize();
1727                 mGvsStats.mPvsCost = mCurrentViewCell->GetPvs().EvalPvsCost();
1728
1729        }
1730
1731        mGvsStats.mTrianglePvs = (int)mTrianglePvs.size();
1732
1733
1734        ////////////////
1735        // global values
1736
1737        mGvsStats.mViewCells = mProcessedViewCells;
1738        mGvsStats.mTotalTrianglePvs += mGvsStats.mTrianglePvs;
1739        mGvsStats.mTotalTime += mGvsStats.mTimePerViewCell;
1740    mGvsStats.mTotalPvs += mGvsStats.mPerViewCellPvs;
1741        mGvsStats.mTotalSamples += mGvsStats.mPerViewCellSamples;
1742
1743        cout << "id: " << mCurrentViewCell->GetId() << endl;
1744        cout << "pvs cost "  << mGvsStats.mPvsCost / 1000000 << " M" << endl;
1745        cout << "pvs tri: " << mTrianglePvs.size() << endl;
1746        cout << "samples: " << mGvsStats.mTotalSamples / 1000000 << " M" << endl;
1747        printf("contri: %f tri / K rays\n", 1000.0f * (float)mTrianglePvs.size() / mGvsStats.mTotalSamples);
1748
1749        //cout << "invalid samples ratio: " << (float)sInvalidSamples / mGvsStats.mPerViewCellSamples << endl;
1750
1751        double rayTime = rayTimer.TotalTime();
1752        double kdPvsTime = kdPvsTimer.TotalTime();
1753        double gvsTime = gvsTimer.TotalTime();
1754        double initialTime = initialTimer.TotalTime();
1755        double intersectionTime = intersectionTimer.TotalTime();
1756        double preparationTime = preparationTimer.TotalTime();
1757        double mainLoopTime = mainLoopTimer.TotalTime();
1758        double contributionTime = contributionTimer.TotalTime();
1759        double castTime = castTimer.TotalTime();
1760        double generationTime = generationTimer.TotalTime();
1761        double rawCastTime = mRayCaster->rawCastTimer.TotalTime();
1762
1763        cout << "handleRay              : " << rayTime << endl;
1764        cout << "kd pvs                 : " << kdPvsTime << endl;
1765        cout << "gvs sampling           : " << gvsTime << endl;
1766        cout << "initial sampling       : " << initialTime << endl;
1767        cout << "view cell intersection : " << intersectionTime << endl;
1768        cout << "preparation            : " << preparationTime << endl;
1769        cout << "main loop              : " << mainLoopTime << endl;
1770        cout << "has contribution       : " << contributionTime << endl;
1771        cout << "cast time              : " << castTime << endl;
1772        cout << "generation time        : " << generationTime << endl;
1773        cout << "raw cast time          : " << rawCastTime << endl;
1774
1775        double randomRaysPerSec = mGvsStats.mRandomSamples / (1e6 * initialTime);
1776        double gvsRaysPerSec = mGvsStats.mGvsSamples / (1e6 * gvsTime);
1777        double rawRaysPerSec = mGvsStats.mPerViewCellSamples / (1e6 * rawCastTime);
1778
1779        cout << "cast " << randomRaysPerSec << " M random rays/s: " << endl;
1780        cout << "cast " << gvsRaysPerSec << " M gvs rays/s: " << endl;
1781        cout << "cast " << rawRaysPerSec << " M raw rays/s: " << endl;
1782
1783        mGvsStats.Stop();
1784        mGvsStats.Print(mGvsStatsStream);
1785}
1786
1787
1788int GvsPreprocessor::ConvertObjectPvs()
1789{
1790        static ObjectContainer triangles;
1791
1792        int newObjects = 0;
1793
1794        ObjectPvsIterator pit = mCurrentViewCell->GetPvs().GetIterator();
1795
1796        triangles.clear();
1797
1798        // output PVS of view cell
1799        while (pit.HasMoreEntries())
1800        {               
1801                KdIntersectable *kdInt = static_cast<KdIntersectable *>(pit.Next());
1802
1803                mKdTree->CollectObjectsWithDublicates(kdInt->GetItem(), triangles);
1804
1805                ObjectContainer::const_iterator oit, oit_end = triangles.end();
1806
1807                for (oit = triangles.begin(); oit != oit_end; ++ oit)
1808                {
1809                        if (AddTriangleObject(*oit))
1810                                ++ newObjects;
1811                }
1812        }
1813
1814        return newObjects;
1815}
1816
1817
1818bool GvsPreprocessor::AddTriangleObject(Intersectable *triObj)
1819{
1820        if ((triObj->mCounter < ACCOUNTED_OBJECT))
1821        {
1822                triObj->mCounter += ACCOUNTED_OBJECT;
1823
1824                mTrianglePvs.push_back(triObj);
1825                mGenericStats[0] = (int)mTrianglePvs.size();
1826
1827                return true;
1828        }
1829
1830        return false;
1831}
1832
1833
1834void GvsPreprocessor::CastRayBundles4(const SimpleRayContainer &rays, VssRayContainer &vssRays)
1835{
1836        AxisAlignedBox3 box = mViewCellsManager->GetViewSpaceBox();
1837
1838        SimpleRayContainer::const_iterator it, it_end = rays.end();
1839
1840        for (it = rays.begin(); it != it_end; ++ it)
1841                CastRayBundle4(*it, vssRays, box);
1842}
1843
1844
1845void GvsPreprocessor::CastRayBundle4(const SimpleRay &ray,
1846                                                                         VssRayContainer &vssRays,
1847                                                                         const AxisAlignedBox3 &box)
1848{
1849        const float pertubDir = 0.01f;
1850        static Vector3 pertub;
1851
1852        static int hit_triangles[4];
1853        static float dist[4];
1854        static Vector3 dirs[4];
1855        static Vector3 origs[4];
1856
1857        origs[0] = ray.mOrigin;
1858        dirs[0] = ray.mDirection;
1859
1860        for (int i = 1; i < 4; ++ i)
1861        {
1862                origs[i] = ray.mOrigin;
1863
1864                pertub.x = RandomValue(-pertubDir, pertubDir);
1865                pertub.y = RandomValue(-pertubDir, pertubDir);
1866                pertub.z = RandomValue(-pertubDir, pertubDir);
1867
1868                dirs[i] = ray.mDirection + pertub;
1869                dirs[i] *= 1.0f / Magnitude(dirs[i]);
1870        }
1871
1872        Cast4Rays(dist, dirs, origs, vssRays, box);
1873}
1874
1875
1876void GvsPreprocessor::CastRays4(const SimpleRayContainer &rays, VssRayContainer &vssRays)
1877{
1878        SimpleRayContainer::const_iterator it, it_end = rays.end();
1879       
1880        static float dist[4];
1881        static Vector3 dirs[4];
1882        static Vector3 origs[4];
1883
1884        AxisAlignedBox3 box = mViewCellsManager->GetViewSpaceBox();
1885
1886        for (it = rays.begin(); it != it_end; ++ it)
1887        {
1888                for (int i = 1; i < 4; ++ i)
1889                {
1890                        origs[i] = rays[i].mOrigin;
1891                        dirs[i] = rays[i].mDirection;
1892                }
1893        }
1894
1895        Cast4Rays(dist, dirs, origs, vssRays, box);
1896}
1897
1898
1899void GvsPreprocessor::Cast4Rays(float *dist,
1900                                                                Vector3 *dirs,
1901                                                                Vector3 *origs,
1902                                                                VssRayContainer &vssRays,
1903                                                                const AxisAlignedBox3 &box)
1904{
1905        static int hit_triangles[4];
1906
1907        mRayCaster->CastRaysPacket4(box.Min(), box.Max(), origs, dirs, hit_triangles,   dist);       
1908
1909        VssRay *ray;
1910
1911        for (int i = 0; i < 4; ++ i)
1912        {
1913                if (hit_triangles[i] != -1)
1914                {
1915                        TriangleIntersectable *triObj = static_cast<TriangleIntersectable *>(mObjects[hit_triangles[i]]);
1916
1917                        if (DotProd(triObj->GetItem().GetNormal(), dirs[i] >= 0))
1918                                ray = NULL;
1919                        else
1920                                ray = mRayCaster->RequestRay(origs[i], origs[i] + dirs[i] * dist[i], NULL, triObj, mPass, 1.0f);
1921                }
1922                else
1923                {
1924                        ray = NULL;
1925                }
1926
1927                vssRays.push_back(ray);
1928        }
1929}
1930
1931
1932void GvsPreprocessor::CastRayBundle16(const SimpleRay &ray, VssRayContainer &vssRays, float scale)
1933{
1934        static SimpleRayContainer simpleRays;
1935        simpleRays.clear();
1936
1937        simpleRays.push_back(ray);
1938        GenerateJitteredRays(simpleRays, ray, 15, 0, scale);
1939
1940        CastRays(simpleRays, vssRays, false, false);
1941}
1942
1943
1944void GvsPreprocessor::CastRayBundles16(const SimpleRayContainer &rays, VssRayContainer &vssRays, float scale)
1945{
1946        SimpleRayContainer::const_iterator it, it_end = rays.end();
1947
1948        for (it = rays.begin(); it != it_end; ++ it)
1949                CastRayBundle16(*it, vssRays, scale);
1950}
1951
1952
1953void GvsPreprocessor::ComputeRenderError()
1954{
1955        // compute rendering error     
1956        if (renderer && renderer->mPvsStatFrames)
1957        {
1958                if (!mViewCellsManager->GetViewCellPointsList()->empty())
1959                {
1960                        ViewCellPointsList *vcPoints = mViewCellsManager->GetViewCellPointsList();
1961                        ViewCellPointsList::const_iterator vit = vcPoints->begin(),     vit_end = vcPoints->end();
1962
1963                        SimpleRayContainer viewPoints;
1964
1965                        for (; vit != vit_end; ++ vit)
1966                        {
1967                                ViewCellPoints *vp = *vit;
1968                               
1969                                SimpleRayContainer::const_iterator rit = vp->second.begin(), rit_end = vp->second.end();
1970                       
1971                                for (; rit!=rit_end; ++rit)
1972                                {
1973                                        ViewCell *vc = mViewCellsManager->GetViewCell((*rit).mOrigin);
1974                                        if (vc == mCurrentViewCell)
1975                                                viewPoints.push_back(*rit);
1976                                }
1977                        }
1978
1979                        renderer->mPvsErrorBuffer.clear();
1980
1981                        if (viewPoints.size() != renderer->mPvsErrorBuffer.size())
1982                        {
1983                                renderer->mPvsErrorBuffer.resize(viewPoints.size());
1984                                renderer->ClearErrorBuffer();
1985                        }
1986
1987                        cout << "evaluating list of " << viewPoints.size() << " pts" << endl;
1988                        renderer->EvalPvsStat(viewPoints);
1989                }
1990                else
1991                {
1992                        cout << "evaluating random points" << endl;
1993                        renderer->EvalPvsStat();
1994                }
1995
1996                mStats <<
1997                        "#AvgPvsRenderError\n" <<renderer->mPvsStat.GetAvgError()<<endl<<
1998                        "#AvgPixelError\n" <<renderer->GetAvgPixelError()<<endl<<
1999                        "#MaxPixelError\n" <<renderer->GetMaxPixelError()<<endl<<
2000                        "#MaxPvsRenderError\n" <<renderer->mPvsStat.GetMaxError()<<endl<<
2001                        "#ErrorFreeFrames\n" <<renderer->mPvsStat.GetErrorFreeFrames()<<endl<<
2002                        "#AvgRenderPvs\n" <<renderer->mPvsStat.GetAvgPvs()<<endl;
2003        }
2004}
2005
2006
2007void GvsPreprocessor::GenerateImportanceSamples(const VssRay &ray,
2008                                                                                                const Triangle3 &triangle,
2009                                                                                                int numSamples,
2010                                                                                                SimpleRayContainer &simpleRays)
2011{
2012        const size_t samplesSize = simpleRays.size();
2013
2014        while (simpleRays.size() < (samplesSize + numSamples))
2015        {
2016                SimpleRay sray;
2017                if (GenerateImportanceSample(ray, triangle, sray))
2018                        simpleRays.push_back(sray);
2019        }
2020}
2021
2022
2023bool GvsPreprocessor::GenerateImportanceSample(const VssRay &ray,
2024                                                                                           const Triangle3 &triangle,
2025                                                                                           SimpleRay &sray)
2026{
2027        Vector3 d = ray.GetDir();
2028 
2029        // Compute right handed coordinate system from direction
2030        Vector3 U, V;
2031
2032        Vector3 nd = Normalize(d);
2033        nd.RightHandedBase(U, V);
2034       
2035        Vector3 origin = ray.mOrigin;
2036        Vector3 termination = ray.mTermination;
2037
2038        float rr[2];
2039
2040        rr[0] = RandomValue(0, 1);
2041        rr[1] = RandomValue(0, 1);
2042
2043        Vector2 vr2(rr[0], rr[1]);
2044        const float sigma = triangle.GetBoundingBox().Radius() * mRadiusOfInfluence;
2045        Vector2 gaussvec2;
2046
2047        GaussianOn2D(vr2, sigma, // input
2048                         gaussvec2); // output
2049       
2050        const Vector3 shift = gaussvec2.xx * U + gaussvec2.yy * V;
2051
2052        //cout << "t: " << termination;
2053        termination += shift;
2054        //cout << " new t: " << termination << endl;
2055        Vector3 direction = termination - origin;
2056
2057        const float len = Magnitude(direction);
2058       
2059        if (len < Limits::Small)
2060                return false;
2061 
2062        direction /= len;
2063
2064        // $$ jb the pdf is yet not correct for all sampling methods!
2065        const float pdf = 1.0f;
2066        sray = SimpleRay(origin, direction, SamplingStrategy::GVS, pdf);
2067       
2068        return true;
2069}
2070
2071
2072void GvsPreprocessor::PrepareProbablyVisibleSampling()
2073{
2074        // warning: using mailing!
2075        Intersectable::NewMail();
2076
2077        mProbablyVisibleTriangles.clear();
2078        CollectProbablyVisibleTriangles(mProbablyVisibleTriangles);
2079}
2080
2081
2082void GvsPreprocessor::CollectProbablyVisibleTriangles(ObjectContainer &triangles)
2083{
2084        ObjectPvsIterator pit = mCurrentViewCell->GetPvs().GetIterator();
2085
2086        static ObjectContainer tmpTriangles;
2087
2088        while (pit.HasMoreEntries())
2089        {       
2090                tmpTriangles.clear();
2091
2092                KdIntersectable *kdObj = static_cast<KdIntersectable *>(pit.Next());
2093                mKdTree->CollectObjectsWithDublicates(kdObj->GetItem(), tmpTriangles);
2094
2095                ObjectContainer::const_iterator oit, oit_end = tmpTriangles.end();
2096
2097                for (oit = tmpTriangles.begin(); oit != oit_end; ++ oit)
2098                {
2099                        TriangleIntersectable *triObj = static_cast<TriangleIntersectable *>(*oit);
2100
2101                        // find objects which are not yet accounted for yet contained in kd pvs objects
2102                        if (!triObj->Mailed() && (triObj->mCounter < ACCOUNTED_OBJECT))
2103                        {
2104                                triObj->Mail();
2105                                triangles.push_back(triObj);
2106                        }
2107                }
2108        }
2109}
2110
2111
2112void GvsPreprocessor::CollectProbablyVisibleItems(ProbablyVisibleQueue &vqueue)
2113{
2114        ObjectPvsIterator pit = mCurrentViewCell->GetPvs().GetIterator();
2115        static ObjectContainer tmpTriangles;
2116
2117        while (pit.HasMoreEntries())
2118        {
2119                tmpTriangles.clear();
2120
2121                Intersectable::NewMail();
2122
2123                ProbablyVisibleItem *item = new ProbablyVisibleItem();
2124
2125                KdIntersectable *kdObj = static_cast<KdIntersectable *>(pit.Next());
2126                mKdTree->CollectObjects(kdObj->GetItem(), tmpTriangles);
2127                item->mNode = kdObj->GetItem();
2128
2129                ObjectContainer::const_iterator oit, oit_end = tmpTriangles.end();
2130
2131                for (oit = tmpTriangles.begin(); oit != oit_end; ++ oit)
2132                {
2133                        TriangleIntersectable *triObj = static_cast<TriangleIntersectable *>(*oit);
2134
2135                        // find objects which are not yet accounted for yet contained in kd pvs objects
2136                        if (triObj->mCounter < ACCOUNTED_OBJECT)
2137                        {
2138                                item->mTriangles.push_back(triObj);
2139                        }
2140                }
2141        }
2142}
2143
2144
2145void GvsPreprocessor::CreateRandomizedReverseRays(const VssRay &ray,
2146                                                                                                  const Vector3 &predictedPoint,
2147                                                                                                  SimpleRayContainer &rays,
2148                                                                                                  int number,
2149                                                                                                  const Triangle3 &triangle)
2150{
2151        size_t currentNumber = rays.size();
2152
2153        int tries = 0;
2154
2155        while (((rays.size() - currentNumber) < number) && (tries ++ < 20))
2156        {
2157                SimpleRay simpleRay;
2158                if (CreateRandomizedReverseRay(ray, predictedPoint, simpleRay, triangle))
2159                        rays.push_back(simpleRay);
2160        }
2161
2162        if (rays.size() != 16)
2163                cout << "x";
2164        //else cout <<"y";
2165        //cout << "generated " << rays.size() << " rays" << endl;
2166}
2167
2168
2169bool GvsPreprocessor::CreateRandomizedReverseRay(const VssRay &ray,
2170                                                                                                 const Vector3 &predictedPoint,
2171                                                                                                 SimpleRay &simpleRay,
2172                                                                                                 const Triangle3 &triangle)
2173{
2174        /*Vector3 nd = ray.GetNormalizedDir();
2175        const AxisAlignedBox3 box = mCurrentViewCell->GetBox();
2176        //Vector3 nd = box.GetNearestFace(origin).GetNormal();
2177        //cout << "pt: " << origin << " face: " << box.GetNearestFace(origin) << endl;
2178
2179        // Compute right handed coordinate system from direction
2180        Vector3 U, V;
2181        nd.RightHandedBase(U, V);
2182       
2183        float rr[2];
2184
2185        rr[0] = RandomValue(0, 1);
2186        rr[1] = RandomValue(0, 1);
2187
2188        Vector2 vr2(rr[0], rr[1]);
2189        const float sigma = triangle.GetBoundingBox().Radius() * 1.0f;
2190        Vector2 gaussvec2;
2191
2192        GaussianOn2D(vr2, sigma, // input
2193                         gaussvec2); // output
2194        */
2195        const Vector3 shift = 0;//gaussvec2.xx * U + gaussvec2.yy * V;
2196
2197        const Vector3 newTermination = ray.mTermination + shift;
2198       
2199        //cout << "o: " << origin << " new o: " << newOrigin << endl;
2200        Vector3 newOrigin = mCurrentViewCell->GetBox().GetRandomPoint();
2201        Vector3 direction = newTermination - newOrigin;
2202
2203        const float len = Magnitude(direction);
2204       
2205        if (len < Limits::Small)
2206                return false;
2207 
2208        direction /= len;
2209
2210        if (0 && !IntersectsViewCell(newOrigin, -direction))
2211                return false;
2212
2213        //cout << "y";
2214        // $$ jb the pdf is yet not correct for all sampling methods!
2215        const float pdf = 1.0f;
2216        simpleRay = SimpleRay(newOrigin, direction, SamplingStrategy::GVS, pdf);
2217       
2218        return true;
2219}
2220
2221}
Note: See TracBrowser for help on using the repository browser.