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

Revision 386, 18.6 KB checked in by bittner, 19 years ago (diff)

VssPreprocessor? updates - directional rays, x3dexporter updates - removed duplicated code for ray exports

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