source: GTP/trunk/Lib/Vis/Preprocessing/src/MutualVisibility.cpp @ 863

Revision 863, 18.7 KB checked in by mattausch, 19 years ago (diff)

working on preprocessor integration
added iv stuff

Line 
1#include <assert.h>
2#include <stack>
3using namespace std;
4#include "KdTree.h"
5#include "AxisAlignedBox3.h"
6#include "Ray.h"
7#include "MutualVisibility.h"
8#include "Exporter.h"
9#include "Mesh.h"
10#include "Triangle3.h"
11#include "SceneGraph.h"
12
13namespace GtpVisibilityPreprocessor {
14
15void
16RayShaft::Init(
17               const Rectangle3 &source,
18               const Rectangle3 &target)
19{
20  mDepth = 0;
21  mSource = source;
22  mTarget = target;
23}
24
25
26Vector3
27RayShaft::GetIntersectionPoint(const int rayIndex,
28                               const int depth) const
29{
30  Vector3 origin, direction;
31  Ray ray;
32  GetRaySegment(rayIndex, ray);
33  if (depth >= mSamples[rayIndex].mIntersections.size()) {
34    cerr<<"depth of sample out of limits"<<endl;
35    exit(1);
36  }
37  return ray.Extrap(mSamples[rayIndex].mIntersections[depth].mT);
38}
39
40void
41RayShaft::GetRaySegment(const int i, Ray &ray) const
42{
43  Vector3 origin, direction;
44  GetRay(i, origin, direction);
45  ray.Init(origin, direction, Ray::LINE_SEGMENT);
46  if ( mSamples[i].IsValid() ) {
47    origin = ray.Extrap(mSamples[i].mMinT);
48    direction = ray.Extrap(mSamples[i].mMaxT) - origin;
49    ray.Init(origin, direction, Ray::LINE_SEGMENT);
50  }
51}
52
53void
54RayShaft::GetRay(const int rayIndex,
55                 Vector3 &origin,
56                 Vector3 &direction) const
57{
58 
59  assert(rayIndex < 4);
60 
61  origin = mSource.mVertices[rayIndex];
62  direction = mTarget.mVertices[rayIndex] - origin;
63}
64
65void
66MutualVisibilitySampler::PerformSplit(
67                                      const RayShaft &sample,
68                                      const bool splitSource,
69                                      const int axis,
70                                      RayShaft &sample1,
71                                      RayShaft &sample2
72                                      )
73{
74 
75
76  // split the triangles
77
78  if (splitSource) {
79    sample1.mTarget = sample.mTarget;
80    sample2.mTarget = sample.mTarget;
81    sample.mSource.Split(
82                         axis,
83                         sample1.mSource,
84                         sample2.mSource);
85  } else {
86   
87    sample1.mSource = sample.mSource;
88    sample2.mSource = sample.mSource;
89   
90    sample.mTarget.Split(
91                         axis,
92                         sample1.mTarget,
93                         sample2.mTarget);
94  }
95 
96  // split the intersections
97  switch (axis) {
98  case 0:
99    sample1.mSamples[0].mIntersections = sample.mSamples[0].mIntersections;
100    sample1.mSamples[3].mIntersections = sample.mSamples[3].mIntersections;
101   
102    sample2.mSamples[1].mIntersections = sample.mSamples[1].mIntersections;
103    sample2.mSamples[2].mIntersections = sample.mSamples[2].mIntersections;
104    break;
105
106  case 1:
107    sample1.mSamples[0].mIntersections = sample.mSamples[0].mIntersections;
108    sample1.mSamples[1].mIntersections = sample.mSamples[1].mIntersections;
109
110    sample2.mSamples[2].mIntersections = sample.mSamples[2].mIntersections;
111    sample2.mSamples[3].mIntersections = sample.mSamples[3].mIntersections;
112    break;
113  }
114
115  // the intersections for the new shaft rays will be established
116  // later
117  sample1.mDepth = sample2.mDepth = sample.mDepth+1;
118}
119
120float
121MutualVisibilitySampler::GetSpatialAngle(const RayShaft &shaft,
122                                         const Vector3 &point
123                                         )
124{
125  const int sampleIndices[]={
126    0,1,2,
127    0,2,3
128  };
129 
130  float sum = 0.0f;
131  int i, j=0;
132
133  if (!shaft.IsValid())
134    return 2.0f*mSolidAngleThreshold;
135
136  for (i=0; i < 2; i++) {
137    Triangle3 triangle(shaft.GetIntersectionPoint(sampleIndices[j++], 0),
138                       shaft.GetIntersectionPoint(sampleIndices[j++], 0),
139                       shaft.GetIntersectionPoint(sampleIndices[j++], 0));
140    sum += triangle.GetSpatialAngle(point);
141  }
142
143  return sum;
144}
145
146
147void
148MutualVisibilitySampler::ComputeError(RayShaft &shaft)
149{
150  // evaluate minimal error which can be achieved by more precise evaluation
151  // if this is above the threshold do not proceed further
152
153  // check whether the samples pierce the same mesh
154  Intersectable::NewMail();
155  int i, j;
156  for (i=0; i < 4; i++) {
157    RaySample *sample = &shaft.mSamples[i];
158    for (j=0; j < sample->mIntersections.size(); j++)
159      sample->mIntersections[j].mObject->Mail();
160  }
161
162  for (i=0; i < 4; i++) {
163    RaySample *sample = &shaft.mSamples[i];
164    for (j=0; j < sample->mIntersections.size(); j++) {
165      if (sample->mIntersections[j].mObject->IncMail() == 4) {
166        cerr<<"T";
167        shaft.mError = 0.0f;
168        return;
169      }
170    }
171  }
172   
173 
174  float maxAngle = GetSpatialAngle(shaft,
175                                   shaft.mTarget.GetCenter());
176 
177  for (i=0; i < 3; i++) {
178    float angle = GetSpatialAngle(shaft,
179                                  shaft.mTarget.mVertices[i]);
180    if (angle > maxAngle)
181      maxAngle = angle;
182  }
183
184  maxAngle = MAX_FLOAT;
185  shaft.mError = maxAngle;
186}
187
188
189void
190MutualVisibilitySampler::ConstructInitialSamples2(
191                                                  const AxisAlignedBox3 &source,
192                                                  const AxisAlignedBox3 &target,
193                                                  vector<RayShaft *> &samples
194                                                  )
195{
196  // get all rectangles potentially visible from the source box
197  int i;
198  int sourceMask = 0;
199  for (i=0; i < 8; i++)
200    sourceMask |= source.GetFaceVisibilityMask(target.GetVertex(i));
201
202  // now for each visble source face find all visible target faces
203  for (i=0; i < 6; i++) {
204    Rectangle3 sourceFace = source.GetFace(i);
205    if ( sourceMask &(1<<i) ) {
206      int mask = target.GetFaceVisibilityMask(sourceFace);
207      // construct triangle samples for all visible rectangles
208      for (int j=0; j < 6; j++)
209        if (mask & (1<<j)) {
210          AddInitialSamples2(sourceFace, target.GetFace(j), samples);
211        }
212    }
213  }
214}
215
216void
217MutualVisibilitySampler::ConstructInitialSamples3(
218                                                  const AxisAlignedBox3 &source,
219                                                  const AxisAlignedBox3 &target,
220                                                  vector<RayShaft *> &samples
221                                                  )
222{
223  // get all rectangles potentially visible from the source box
224  int i;
225  int sourceMask = 0;
226  for (i=0; i < 8; i++)
227    sourceMask |= source.GetFaceVisibilityMask(target.GetVertex(i));
228
229  // now for each visble source face find all visible target faces
230  int face;
231  for (face=0; face < 6; face++) {
232    if ( sourceMask & (1<<face) ) {
233      Rectangle3 sourceRect = source.GetFace(face);
234
235      Vector3 targetCenter = target.Center();
236     
237      Vector3 normal = sourceRect.GetNormal();
238     
239      Plane3 sourcePlane(normal, sourceRect.GetVertex(0));
240     
241      Plane3 targetPlane(normal, target.GetVertex(0));
242
243      int i;
244      for (i=1; i < 8; i++) {
245        Vector3 v = target.GetVertex(i);
246        if (targetPlane.Distance(v) < 0)
247          targetPlane = Plane3(normal, v);
248      }
249 
250      Vector3 xBasis = Normalize(sourceRect.GetVertex(1) - sourceRect.GetVertex(0));
251      Vector3 yBasis = Normalize(sourceRect.GetVertex(3) - sourceRect.GetVertex(0));
252
253      // cast rays between the centers of the boxes
254      Vector3 targetRCenter = targetPlane.FindIntersection(sourceRect.GetCenter(),
255                                                           targetCenter);
256     
257      Rectangle3 targetRect;
258     
259      float targetDist[4];
260      targetDist[0] = targetDist[1] = targetDist[2] = targetDist[3] = 0.0f;
261     
262      // cast rays between corresponding vertices of the boxes
263      int j;
264      int intersections=0;
265      for (i=0; i < 4; i++)
266        for (j=0; j < 8; j++) {
267          Vector3 v;
268          Vector3 diff;
269          bool coplanar;
270          float dist;
271         
272          v = targetPlane.FindIntersection(sourceRect.GetVertex(i),
273                                           target.GetVertex(j),
274                                           NULL,
275                                           &coplanar);
276          if (!coplanar) {
277            // evaluate target and
278            diff = targetRCenter - v;
279            dist = DotProd(diff, xBasis);
280            if (dist < targetDist[0])
281              targetDist[0] = dist;
282            if (dist > targetDist[1])
283              targetDist[1] = dist;
284
285            dist = DotProd(diff, yBasis);
286
287            if (dist < targetDist[2])
288              targetDist[2] = dist;
289            if (dist > targetDist[3])
290              targetDist[3] = dist;
291            intersections++;
292          }
293        }
294     
295      if (intersections>=4) {
296        targetRect.mVertices[0] = targetRCenter + targetDist[0]*xBasis + targetDist[2]*yBasis;
297        targetRect.mVertices[1] = targetRCenter + targetDist[1]*xBasis + targetDist[2]*yBasis;
298        targetRect.mVertices[2] = targetRCenter + targetDist[1]*xBasis + targetDist[3]*yBasis;
299        targetRect.mVertices[3] = targetRCenter + targetDist[0]*xBasis + targetDist[3]*yBasis;
300       
301        //      cout<<sourceRect<<targetRect<<endl;
302        AddInitialSamples(sourceRect, targetRect, samples);
303      }
304    }
305  }
306}
307
308void
309MutualVisibilitySampler::ConstructInitialSamples(
310                                                 const AxisAlignedBox3 &source,
311                                                 const AxisAlignedBox3 &target,
312                                                 vector<RayShaft *> &samples
313                                                 )
314{
315  Vector3 sourceCenter = source.Center();
316  Vector3 targetCenter = target.Center();
317  Vector3 normal = Normalize(sourceCenter - targetCenter );
318 
319  Plane3 sourcePlane(normal, source.GetVertex(0));
320  int i;
321  for (i=1; i < 8; i++) {
322    Vector3 v = source.GetVertex(i);
323    if (sourcePlane.Distance(v) > 0)
324      sourcePlane = Plane3(normal, v);
325  }
326 
327  Plane3 targetPlane(normal, target.GetVertex(0));
328  for (i=1; i < 8; i++) {
329    Vector3 v = target.GetVertex(i);
330    if (targetPlane.Distance(v) < 0)
331      targetPlane = Plane3(normal, v);
332  }
333
334 
335  Vector3 xBasis = CrossProd(Vector3(0,1,0), normal);
336
337  if (Magnitude(xBasis) > 1e-6)
338    xBasis.Normalize();
339  else {
340    xBasis = Normalize(CrossProd(Vector3(0,0,1), normal));
341  }
342
343  Vector3 yBasis = Normalize( CrossProd(normal, xBasis) );
344  // cast rays between the centers of the boxes
345  Vector3 targetRCenter = targetPlane.FindIntersection(sourceCenter,
346                                                       targetCenter);
347 
348  Vector3 sourceRCenter = sourcePlane.FindIntersection(sourceCenter,
349                                                       targetCenter);
350 
351
352  Rectangle3 sourceRect;
353  Rectangle3 targetRect;
354
355
356  if (0) {
357  float scale = Magnitude(source.Size())*0.7f;
358  sourceRect.mVertices[0] = sourceRCenter - scale*xBasis - scale*yBasis;
359  sourceRect.mVertices[1] = sourceRCenter + scale*xBasis - scale*yBasis;
360  sourceRect.mVertices[2] = sourceRCenter + scale*xBasis + scale*yBasis;
361  sourceRect.mVertices[3] = sourceRCenter - scale*xBasis + scale*yBasis;
362
363  scale = Magnitude(target.Size())*0.7f;
364  targetRect.mVertices[0] = targetRCenter - scale*xBasis - scale*yBasis;
365  targetRect.mVertices[1] = targetRCenter + scale*xBasis - scale*yBasis;
366  targetRect.mVertices[2] = targetRCenter + scale*xBasis + scale*yBasis;
367  targetRect.mVertices[3] = targetRCenter - scale*xBasis + scale*yBasis;
368  }
369
370
371  float sourceDist[4];
372  float targetDist[4];
373  sourceDist[0] = sourceDist[1] = sourceDist[2] = sourceDist[3] = 0.0f;
374  targetDist[0] = targetDist[1] = targetDist[2] = targetDist[3] = 0.0f;
375
376  // cast rays between corresponding vertices of the boxes
377  int j;
378  for (i=0; i < 8; i++)
379    for (j=0; j < 8; j++) {
380      Vector3 v;
381      Vector3 diff;
382      bool coplanar;
383      float dist;
384      v = sourcePlane.FindIntersection(source.GetVertex(i),
385                                       target.GetVertex(j),
386                                       NULL,
387                                       &coplanar);
388      if (!coplanar) {
389        // evaluate source and
390        diff = sourceRCenter - v;
391        dist = DotProd(diff, xBasis);
392        if (dist < sourceDist[0])
393          sourceDist[0] = dist;
394        if (dist > sourceDist[1])
395          sourceDist[1] = dist;
396
397        dist = DotProd(diff, yBasis);
398        if (dist < sourceDist[2])
399          sourceDist[2] = dist;
400        if (dist > sourceDist[3])
401          sourceDist[3] = dist;
402      }
403     
404      v = targetPlane.FindIntersection(source.GetVertex(i),
405                                       target.GetVertex(j),
406                                       NULL,
407                                       &coplanar);
408      if (!coplanar) {
409        // evaluate target and
410        diff = targetRCenter - v;
411        dist = DotProd(diff, xBasis);
412        if (dist < targetDist[0])
413          targetDist[0] = dist;
414        if (dist > targetDist[1])
415          targetDist[1] = dist;
416        dist = DotProd(diff, yBasis);
417        if (dist < targetDist[2])
418          targetDist[2] = dist;
419        if (dist > targetDist[3])
420          targetDist[3] = dist;
421      }
422    }
423
424  sourceRect.mVertices[0] = sourceRCenter + sourceDist[0]*xBasis + sourceDist[2]*yBasis;
425  sourceRect.mVertices[1] = sourceRCenter + sourceDist[1]*xBasis + sourceDist[2]*yBasis;
426  sourceRect.mVertices[2] = sourceRCenter + sourceDist[1]*xBasis + sourceDist[3]*yBasis;
427  sourceRect.mVertices[3] = sourceRCenter + sourceDist[0]*xBasis + sourceDist[3]*yBasis;
428
429  targetRect.mVertices[0] = targetRCenter + targetDist[0]*xBasis + targetDist[2]*yBasis;
430  targetRect.mVertices[1] = targetRCenter + targetDist[1]*xBasis + targetDist[2]*yBasis;
431  targetRect.mVertices[2] = targetRCenter + targetDist[1]*xBasis + targetDist[3]*yBasis;
432  targetRect.mVertices[3] = targetRCenter + targetDist[0]*xBasis + targetDist[3]*yBasis;
433
434
435  cout<<sourceRect<<targetRect<<endl;
436  AddInitialSamples(sourceRect, targetRect, samples);
437}
438
439void
440MutualVisibilitySampler::AddInitialSamples(
441                                           const Rectangle3 &sourceRect,
442                                           const Rectangle3 &targetRect,
443                                           vector<RayShaft *> &samples
444                                           )
445{
446
447  RayShaft *sample = new RayShaft(sourceRect,
448                                  targetRect);
449  samples.push_back(sample);
450}
451
452void
453MutualVisibilitySampler::AddInitialSamples2(
454                                            const Rectangle3 &sourceRect,
455                                            const Rectangle3 &targetRect,
456                                            vector<RayShaft *> &samples
457                                            )
458{
459  // align the rectangles properly
460  Rectangle3 sRect, tRect;
461  if (DotProd(sourceRect.GetNormal(), targetRect.GetNormal()) < -0.99f) {
462    int i;
463    for (i=0; i < 4; i++) {
464      tRect.mVertices[i] = targetRect.mVertices[( 7 - i )%4];
465    }
466   
467    RayShaft *sample = new RayShaft(sourceRect,
468                                    tRect);
469    samples.push_back(sample);
470  }
471 
472  return;
473 
474  float minDist = MAX_FLOAT;
475  int startI;
476  int startJ;
477  int i, j;
478
479 
480  for (i=0; i < 4; i++) {
481    for (j=0; j < 4; j++) {
482      float dist = Distance(sourceRect.mVertices[i], targetRect.mVertices[j]);
483      if (dist < minDist) {
484        minDist = dist;
485        startI = i;
486        startJ = j;
487      }
488    }
489  }
490  for (i=0; i < 4; i++) {
491    sRect.mVertices[i] = sourceRect.mVertices[(startI + i )%4];
492    tRect.mVertices[i] = targetRect.mVertices[(4 + startJ - i )%4];
493  }
494 
495  RayShaft *sample = new RayShaft(sRect,
496                                  tRect);
497  samples.push_back(sample);
498}
499
500
501void
502MutualVisibilitySampler::ExportShafts(vector<RayShaft *> &shafts,
503                                      const bool singleFile)
504{
505  static int id = 0;
506  char filename[64];
507  if (id > 20)
508    return;
509 
510  Exporter *exporter = NULL;
511  for (int i=0; i < shafts.size(); i++) {
512    if (!exporter) {
513      if (singleFile)
514        sprintf(filename, "shafts-single-%04d.x3d", id++);
515      else
516        sprintf(filename, "shafts%04d-%02d.x3d", id++, i);
517      exporter = Exporter::GetExporter(filename);
518    }
519   
520    exporter->SetWireframe();
521    //    exporter->ExportScene(mSceneGraph->mRoot);
522
523    exporter->ExportBox(mSource);
524    exporter->ExportBox(mTarget);
525
526    exporter->SetFilled();
527
528
529    RayShaft *shaft = shafts[i];
530    Mesh *mesh = new Mesh;
531    mesh->AddRectangle(shaft->mSource);
532    mesh->AddRectangle(shaft->mTarget);
533    vector<Ray *> rays;
534    for (int j=0; j < 4; j++) {
535      Ray *ray = new Ray;
536      shaft->GetRaySegment(j, *ray);
537      rays.push_back(ray);
538    }
539   
540    Material m = RandomMaterial();
541    exporter->SetForcedMaterial(m);
542    MeshInstance mi(mesh);
543    exporter->ExportIntersectable(&mi);
544    exporter->ExportRays(rays, -1.0f, m.mDiffuseColor);
545    if (!singleFile) {
546      delete exporter;
547      exporter = NULL;
548    }
549  }
550  if (exporter)
551    delete exporter;
552 
553}
554
555int
556MutualVisibilitySampler::CastRays(RayShaft &shaft)
557{
558  Ray ray;
559  int i;
560
561  for (i=0; i < 4; i++)
562    if (!shaft.mSamples[i].IsProcessed()) {
563      Vector3 origin, direction;
564      shaft.GetRay(i, origin, direction);
565      // determine intersections with the boxes
566      ray.Init(origin, direction, Ray::LINE_SEGMENT);
567      float stmin, stmax = 0.0f, ttmin=1.0f, ttmax;
568      bool valid = true;
569     
570      if (mUseBoxes) {
571        if (mSource.GetMinMaxT(ray, &stmin, &stmax) &&
572            mTarget.GetMinMaxT(ray, &ttmin, &ttmax)) {
573          shaft.mSamples[i].mMinT = stmax;
574          shaft.mSamples[i].mMaxT = ttmin;
575          origin = ray.Extrap(stmax);
576          direction = ray.Extrap(ttmin) - origin;
577          // reinit the ray
578          ray.Init(origin, direction, Ray::LINE_SEGMENT);
579        } else
580          valid = false;
581      } else {
582        shaft.mSamples[i].mMinT = 0.0f;
583        shaft.mSamples[i].mMaxT = 1.0f;
584      }
585      if (valid) {
586        if (!mKdTree->CastRay(ray)) {
587          cerr<<"V"<<endl;
588          return VISIBLE;
589        }
590        shaft.mSamples[i].mIntersections = ray.intersections;
591        cerr<<"I";
592      } else {
593        cerr<<"X";
594        shaft.mSamples[i].SetInvalid();
595      }
596    }
597  return INVISIBLE;
598}
599
600int
601MutualVisibilitySampler::ComputeVisibility()
602{
603  int result = INVISIBLE;
604
605  vector<RayShaft *> shafts;
606  ConstructInitialSamples3(mSource, mTarget, shafts);
607       
608  if (1)
609    ExportShafts(shafts, false);
610
611  stack<RayShaft *> shaftStack;
612  int i;
613  for (i=0; i < shafts.size(); i++)
614    shaftStack.push(shafts[i]);
615
616  shafts.clear();
617  Ray ray;
618 
619// now process the shafts as long as we have something to do
620  while (!shaftStack.empty()) {
621    RayShaft *shaft = shaftStack.top();
622    shaftStack.pop();
623   
624//      // cast a new ray
625//      int triangleSplitEdge = SetupExtremalRay(sample, source, ray);
626   
627    if (CastRays(*shaft) == VISIBLE) {
628      result = VISIBLE;
629      break;
630    }
631               
632    // compute error ....
633    ComputeError(*shaft);
634    cout<<shaft->mDepth<<"|";
635    if (shaft->IsValid())
636      shafts.push_back(shaft);
637
638    if (shaft->mDepth < 10 &&
639                                shaft->mError > mSolidAngleThreshold) {
640     
641      // generate 2 new samples
642      RayShaft *newSamples[2];
643      newSamples[0] = new RayShaft;
644      newSamples[1] = new RayShaft;
645     
646      // chose what to split
647      bool splitSource = shaft->mSource.GetArea() > shaft->mTarget.GetArea();
648      int axis;
649      if (splitSource) {
650                                axis = shaft->mSource.DominantAxis();
651      } else {
652                                axis = shaft->mTarget.DominantAxis();
653      }
654     
655      PerformSplit(*shaft, splitSource, axis, *newSamples[0], *newSamples[1]);
656      //      delete shaft;
657      shaftStack.push(newSamples[0]);
658      shaftStack.push(newSamples[1]);
659    } else {
660      // store a terminal shaft
661    }
662
663  }
664 
665  while (!shaftStack.empty()) {
666    RayShaft *shaft = shaftStack.top();
667    shaftStack.pop();
668    delete shaft;
669  }
670
671  if (0)
672    ExportShafts(shafts, true);
673
674  for (i=0; i < shafts.size(); i++) {
675    delete shafts[i];
676  }
677
678  return result;
679}
680
681MutualVisibilitySampler::MutualVisibilitySampler(SceneGraph *sceneGraph,
682                                                                                                                                                                                                 KdTree *kdTree,
683                                                                                                                                                                                                 const AxisAlignedBox3 &source,
684                                                                                                                                                                                                 const AxisAlignedBox3 &target,
685                                                                                                                                                                                                 const float solidAngleThreshold)
686{
687  mSceneGraph = sceneGraph;
688  mKdTree = kdTree;
689  mSource = source;
690  mTarget = target;
691  mUseBoxes = true;
692  mSolidAngleThreshold = solidAngleThreshold;
693}
694 
695
696int
697ComputeBoxVisibility(SceneGraph *sceneGraph,
698                                                                                 KdTree *kdTree,
699                                                                                 const AxisAlignedBox3 &source,
700                                                                                 const AxisAlignedBox3 &target,
701                                                                                 const float solidAngleThreshold)
702{
703  MutualVisibilitySampler
704                sampler(sceneGraph, kdTree, source, target, solidAngleThreshold);
705
706 
707  int visibility = sampler.ComputeVisibility();
708
709  return visibility;
710 
711}
712 
713}
Note: See TracBrowser for help on using the repository browser.