source: GTP/trunk/Lib/Vis/Preprocessing/src/AxisAlignedBox3.cpp @ 1581

Revision 1581, 60.0 KB checked in by bittner, 18 years ago (diff)

Qtglwidget changes - second view + trackball - ray casting seems not functional

RevLine 
[162]1
2// GOLEM library
3#include <assert.h>
4#include <iostream>
5using namespace std;
6#include "AxisAlignedBox3.h"
7#include "Ray.h"
[295]8#include "Polygon3.h"
[486]9#include "Mesh.h"
[162]10
[863]11namespace GtpVisibilityPreprocessor {
[860]12
[162]13#define FATAL Debug
14#define FATAL_ABORT exit(1)
15
16
17// AxisAlignedBox3 implementations
18
19// Overload << operator for C++-style output
20ostream&
21operator<< (ostream &s, const AxisAlignedBox3 &A)
22{
23  return s << '[' << A.mMin.x << ", " << A.mMin.y << ", " << A.mMin.z << "]["
24           << A.mMax.x << ", " << A.mMax.y << ", " << A.mMax.z << ']';
25}
26
27// Overload >> operator for C++-style input
28istream&
29operator>> (istream &s, AxisAlignedBox3 &A)
30{
31  char a;
32  // read "[min.x, min.y, min.z][mMax.x, mMax.y, mMax.z]"
33  return s >> a >> A.mMin.x >> a >> A.mMin.y >> a >> A.mMin.z >> a >> a
34    >> A.mMax.x >> a >> A.mMax.y >> a >> A.mMax.z >> a;
35}
36
37bool
38AxisAlignedBox3::Unbounded() const
39{
40  return (mMin == Vector3(-MAXFLOAT)) ||
41         (mMax == Vector3(-MAXFLOAT));
42}
43
44void
45AxisAlignedBox3::Include(const Vector3 &newpt)
46{
47  Minimize(mMin, newpt);
48  Maximize(mMax, newpt);
49}
50
51void
[295]52AxisAlignedBox3::Include(const Polygon3 &newpoly)
53{
54        VertexContainer::const_iterator it, it_end = newpoly.mVertices.end();
55
56        for (it = newpoly.mVertices.begin(); it != it_end; ++ it)
57                Include(*it);
58}
59
60void
[860]61AxisAlignedBox3::Include(const PolygonContainer &polys)
62{
63        PolygonContainer::const_iterator it, it_end = polys.end();
64
65        for (it = polys.begin(); it != it_end; ++ it)
66                Include(*(*it));
67}
68
69void
[538]70AxisAlignedBox3::Include(Mesh *mesh)
71{
72        VertexContainer::const_iterator it, it_end = mesh->mVertices.end();
73
74        for (it = mesh->mVertices.begin(); it != it_end; ++ it)
75                Include(*it);
76}
77
78
79void
[162]80AxisAlignedBox3::Include(const AxisAlignedBox3 &bbox)
81{
82  Minimize(mMin, bbox.mMin);
83  Maximize(mMax, bbox.mMax);
84}
85
86bool
87AxisAlignedBox3::IsCorrect()
88{
89  if ( (mMin.x > mMax.x) ||
90       (mMin.y > mMax.y) ||
91       (mMin.z > mMax.z) )
92    return false; // box is not formed
93  return true;
94}
95
96void
97AxisAlignedBox3::GetEdge(const int edge, Vector3 *a, Vector3 *b) const
98{
99  switch(edge) {
100  case 0:
101    a->SetValue(mMin.x, mMin.y, mMin.z);
102    b->SetValue(mMin.x, mMin.y, mMax.z);
103    break;
104  case 1:
105    a->SetValue(mMin.x, mMin.y, mMin.z);
106    b->SetValue(mMin.x, mMax.y, mMin.z);
107    break;
108  case 2:
109    a->SetValue(mMin.x, mMin.y, mMin.z);
110    b->SetValue(mMax.x, mMin.y, mMin.z);
111    break;
112  case 3:
113    a->SetValue(mMax.x, mMax.y, mMax.z);
114    b->SetValue(mMax.x, mMax.y, mMin.z);
115    break;
116  case 4:
117    a->SetValue(mMax.x, mMax.y, mMax.z);
118    b->SetValue(mMax.x, mMin.y, mMax.z);
119    break;
120  case 5:
121    a->SetValue(mMax.x, mMax.y, mMax.z);
122    b->SetValue(mMin.x, mMax.y, mMax.z);
123    break;
124   
125  case 6:
126    a->SetValue(mMin.x, mMin.y, mMax.z);
127    b->SetValue(mMin.x, mMax.y, mMax.z);
128    break;
129  case 7:
130    a->SetValue(mMin.x, mMin.y, mMax.z);
131    b->SetValue(mMax.x, mMin.y, mMax.z);
132    break;
133  case 8:
134    a->SetValue(mMin.x, mMax.y, mMin.z);
135    b->SetValue(mMin.x, mMax.y, mMax.z);
136    break;
137  case 9:
138    a->SetValue(mMin.x, mMax.y, mMin.z);
139    b->SetValue(mMax.x, mMax.y, mMin.z);
140    break;
141  case 10:
142    a->SetValue(mMax.x, mMin.y, mMin.z);
143    b->SetValue(mMax.x, mMax.y, mMin.z);
144    break;
145  case 11:
146    a->SetValue(mMax.x, mMin.y, mMin.z);
147    b->SetValue(mMax.x, mMin.y, mMax.z);
148    break;
149  }
150}
151
[312]152// returns the vertex indices in the range <0..7>, v = 4.x + 2.y + z, where
153// x,y,z are either 0 or 1; (0 .. min coordinate, 1 .. max coordinate)
[162]154void
[312]155AxisAlignedBox3::GetEdge(const int edge, int  &aIdx, int &bIdx) const
156{
157  switch(edge) {
158  case 0: aIdx = 0; bIdx = 1; break;
159  case 1: aIdx = 0; bIdx = 2; break;
160  case 2: aIdx = 0; bIdx = 4; break;
161  case 3: aIdx = 7; bIdx = 6; break;
162  case 4: aIdx = 7; bIdx = 5; break;
163  case 5: aIdx = 7; bIdx = 3; break;
164  case 6: aIdx = 1; bIdx = 3; break;
165  case 7: aIdx = 1; bIdx = 5; break;
166  case 8: aIdx = 2; bIdx = 3; break;
167  case 9: aIdx = 2; bIdx = 6; break;
168  case 10: aIdx = 4; bIdx = 6; break;
169  case 11: aIdx = 4; bIdx = 5; break;
170  }
171}
172
173void
[162]174AxisAlignedBox3::Include(const int &axis, const float &newBound)
175{
176  switch (axis) {
177    case 0: { // x-axis
178      if (mMin.x > newBound)
179        mMin.x = newBound;
180      if (mMax.x < newBound)
181        mMax.x = newBound;
182      break;
183    }
184    case 1: { // y-axis
185      if (mMin.y > newBound)
186        mMin.y = newBound;
187      if (mMax.y < newBound)
188        mMax.y = newBound;
189      break;
190    }
191    case 2: { // z-axis
192      if (mMin.z > newBound)
193        mMin.z = newBound;
194      if (mMax.z < newBound)
195        mMax.z = newBound;
196      break;
197    }
198  }
199}
200
201#if 0
202// ComputeMinMaxT computes the minimum and maximum signed distances
203// of intersection with the ray; it returns 1 if the ray hits
204// the bounding box and 0 if it does not.
205int
206AxisAlignedBox3::ComputeMinMaxT(const Ray &ray,
207                                float *tmin,
208                                float *tmax) const
209{
210  float minx, maxx, miny, maxy, minz, maxz;
211
212  if (fabs(ray.dir.x) < 0.001) {
213    if (mMin.x < ray.loc.x && mMax.x > ray.loc.x) {
214      minx = -MAXFLOAT;
215      maxx = MAXFLOAT;
216    }
217    else
218      return 0;
219  }
220  else {
221    float t1 = (mMin.x - ray.loc.x) / ray.dir.x;
222    float t2 = (mMax.x - ray.loc.x) / ray.dir.x;
223    if (t1 < t2) {
224      minx = t1;
225      maxx = t2;
226    }
227    else {
228      minx = t2;
229      maxx = t1;
230    }
231    if (maxx < 0)
232      return 0;
233  }
234
235  if (fabs(ray.dir.y) < 0.001) {
236    if (mMin.y < ray.loc.y && mMax.y > ray.loc.y) {
237      miny = -MAXFLOAT;
238      maxy = MAXFLOAT;
239    }
240    else
241      return 0;
242  }
243  else {
244    float t1 = (mMin.y - ray.loc.y) / ray.dir.y;
245    float t2 = (mMax.y - ray.loc.y) / ray.dir.y;
246    if (t1 < t2) {
247      miny = t1;
248      maxy = t2;
249    }
250    else {
251      miny = t2;
252      maxy = t1;
253    }
254    if (maxy < 0.0)
255      return 0;
256  }
257
258  if (fabs(ray.dir.z) < 0.001) {
259    if (mMin.z < ray.loc.z && mMax.z > ray.loc.z) {
260      minz = -MAXFLOAT;
261      maxz = MAXFLOAT;
262    }
263    else
264      return 0;
265  }
266  else {
267    float t1 = (mMin.z - ray.loc.z) / ray.dir.z;
268    float t2 = (mMax.z - ray.loc.z) / ray.dir.z;
269    if (t1 < t2) {
270      minz = t1;
271      maxz = t2;
272    }
273    else {
274      minz = t2;
275      maxz = t1;
276    }
277    if (maxz < 0.0)
278      return 0;
279  }
280
281  *tmin = minx;
282  if (miny > *tmin)
283    *tmin = miny;
284  if (minz > *tmin)
285    *tmin = minz;
286
287  *tmax = maxx;
288  if (maxy < *tmax)
289    *tmax = maxy;
290  if (maxz < *tmax)
291    *tmax = maxz;
292
293  return 1; // yes, intersection was found
294}
295#else
296// another variant of the same, with less variables
297int
298AxisAlignedBox3::ComputeMinMaxT(const Ray &ray,
299                                float *tmin,
300                                float *tmax) const
301{
[534]302  const float dirEps = 1e-8f;
[162]303  register float minx, maxx;
304 
[534]305  if (fabs(ray.dir.x) < dirEps) {
[162]306    if (mMin.x < ray.loc.x && mMax.x > ray.loc.x) {
307      minx = -MAXFLOAT;
308      maxx = MAXFLOAT;
309    }
310    else
311      return 0;
312  }
313  else {
314    float t1 = (mMin.x - ray.loc.x) * ray.invDir.x;
315    float t2 = (mMax.x - ray.loc.x) * ray.invDir.x;
316    if (t1 < t2) {
317      minx = t1;
318      maxx = t2;
319    }
320    else {
321      minx = t2;
322      maxx = t1;
323    }
[492]324        //    if (maxx < 0.0)
325        //      return 0;
[162]326  }
327
328  *tmin = minx;
329  *tmax = maxx;
330 
[534]331  if (fabs(ray.dir.y) < dirEps) {
[162]332    if (mMin.y < ray.loc.y && mMax.y > ray.loc.y) {
333      minx = -MAXFLOAT;
334      maxx = MAXFLOAT;
335    }
336    else
337      return 0;
338  }
339  else {
340    float t1 = (mMin.y - ray.loc.y) * ray.invDir.y;
341    float t2 = (mMax.y - ray.loc.y) * ray.invDir.y;
342    if (t1 < t2) {
343      minx = t1;
344      maxx = t2;
345    }
346    else {
347      minx = t2;
348      maxx = t1;
349    }
[492]350        //    if (maxx < 0.0)
351        //      return 0;
[162]352  }
353
354  if (minx > *tmin)
355    *tmin = minx;
356  if (maxx < *tmax)
357    *tmax = maxx;
358 
[534]359  if (fabs(ray.dir.z) < dirEps) {
[162]360    if (mMin.z < ray.loc.z && mMax.z > ray.loc.z) {
361      minx = -MAXFLOAT;
362      maxx = MAXFLOAT;
363    }
364    else
365      return 0;
366  }
367  else {
368    float t1 = (mMin.z - ray.loc.z) * ray.invDir.z;
369    float t2 = (mMax.z - ray.loc.z) * ray.invDir.z;
370    if (t1 < t2) {
371      minx = t1;
372      maxx = t2;
373    }
374    else {
375      minx = t2;
376      maxx = t1;
377    }
[492]378        //    if (maxx < 0.0)
379        //      return 0;
[162]380  }
381
382  if (minx > *tmin)
383    *tmin = minx;
384  if (maxx < *tmax)
385    *tmax = maxx;
386
387  return 1; // yes, intersection was found
388}
389#endif
390
391// ComputeMinMaxT computes the minimum and maximum parameters
392// of intersection with the ray; it returns 1 if the ray hits
393// the bounding box and 0 if it does not.
394int
395AxisAlignedBox3::ComputeMinMaxT(const Ray &ray, float *tmin, float *tmax,
396                                EFaces &entryFace, EFaces &exitFace) const
397{
398  float minx, maxx, miny, maxy, minz, maxz;
399  int swapped[3];
400 
401  if (fabs(ray.dir.x) < 0.001) {
402    if (mMin.x < ray.loc.x && mMax.x > ray.loc.x) {
403      minx = -MAXFLOAT;
404      maxx = MAXFLOAT;
405    }
406    else
407      return 0;
408  }
409  else {
410    float t1 = (mMin.x - ray.loc.x) / ray.dir.x;
411    float t2 = (mMax.x - ray.loc.x) / ray.dir.x;
412    if (t1 < t2) {
413      minx = t1;
414      maxx = t2;
415      swapped[0] = 0;
416    }
417    else {
418      minx = t2;
419      maxx = t1;
420      swapped[0] = 1;
421    }
422    if (maxx < 0)
423      return 0;
424  }
425
426  if (fabs(ray.dir.y) < 0.001) {
427    if (mMin.y < ray.loc.y && mMax.y > ray.loc.y) {
428      miny = -MAXFLOAT;
429      maxy = MAXFLOAT;
430    }
431    else
432      return 0;
433  }
434  else {
435    float t1 = (mMin.y - ray.loc.y) / ray.dir.y;
436    float t2 = (mMax.y - ray.loc.y) / ray.dir.y;
437    if (t1 < t2) {
438      miny = t1;
439      maxy = t2;
440      swapped[1] = 0;
441    }
442    else {
443      miny = t2;
444      maxy = t1;
445      swapped[1] = 1;
446    }
447    if (maxy < 0.0)
448      return 0;
449  }
450
451  if (fabs(ray.dir.z) < 0.001) {
452    if (mMin.z < ray.loc.z && mMax.z > ray.loc.z) {
453      minz = -MAXFLOAT;
454      maxz = MAXFLOAT;
455    }
456    else
457      return 0;
458  }
459  else {
460    float t1 = (mMin.z - ray.loc.z) / ray.dir.z;
461    float t2 = (mMax.z - ray.loc.z) / ray.dir.z;
462    if (t1 < t2) {
463      minz = t1;
464      maxz = t2;
465      swapped[2] = 0;
466    }
467    else {
468      minz = t2;
469      maxz = t1;
470      swapped[2] = 1;
471    }
472    if (maxz < 0.0)
473      return 0;
474  }
475
476  *tmin = minx;
477  entryFace = ID_Back;
478  if (miny > *tmin) {
479    *tmin = miny;
480    entryFace = ID_Left;
481  }
482  if (minz > *tmin) {
483    *tmin = minz;
484    entryFace = ID_Bottom;
485  }
486 
487  *tmax = maxx;
488  exitFace = ID_Back;
489  if (maxy < *tmax) {
490    *tmax = maxy;
491    exitFace = ID_Left;
492  }
493  if (maxz < *tmax) {
494    *tmax = maxz;
495    exitFace = ID_Bottom;
496  }
497
498  if (swapped[entryFace])
499    entryFace = (EFaces)(entryFace + 3);
500 
501  if (!swapped[exitFace])
502    exitFace = (EFaces)(exitFace + 3);
503 
504  return 1; // yes, intersection was found
505}
506
507void
508AxisAlignedBox3::Describe(ostream &app, int ind) const
509{
510  indent(app, ind);
511  app << "AxisAlignedBox3: min at(" << mMin << "), max at(" << mMax << ")\n";
512}
513
514// computes the passing through parameters for case tmin<tmax and tmax>0
515int
516AxisAlignedBox3::GetMinMaxT(const Ray &ray, float *tmin, float *tmax) const
517{
518  if (!ComputeMinMaxT(ray, tmin, tmax))
519    return 0;
520  if ( *tmax < *tmin)
521    return 0; // the ray passes outside the box
522
523  if ( *tmax < 0.0)
524    return 0; // the intersection is not on the positive halfline
525
526  return 1; // ray hits the box .. origin can be outside or inside the box
527}
528
529// computes the signed distances for case tmin<tmax and tmax>0
530int
531AxisAlignedBox3::GetMinMaxT(const Ray &ray, float *tmin, float *tmax,
532                            EFaces &entryFace, EFaces &exitFace) const
533{
534  if (!ComputeMinMaxT(ray, tmin, tmax, entryFace, exitFace))
535    return 0;
536  if ( *tmax < *tmin)
537    return 0; // the ray passes outside the box
538 
539  if ( *tmax < 0.0)
540    return 0; // the intersection is not on the positive halfline
541 
542  return 1; // ray hits the box .. origin can be outside or inside the box
543}
544
545#if 0
546int
547AxisAlignedBox3::IsInside(const Vector3 &point) const
548{
549  return (point.x >= mMin.x) && (point.x <= mMax.x) &&
550         (point.y >= mMin.y) && (point.y <= mMax.y) &&
551         (point.z >= mMin.z) && (point.z <= mMax.z);
552}
553#else
554int
555AxisAlignedBox3::IsInside(const Vector3 &v) const
556{
557  return ! (v.x < mMin.x ||
558            v.x > mMax.x ||
559            v.y < mMin.y ||
560            v.y > mMax.y ||
561            v.z < mMin.z ||
562            v.z > mMax.z);
563}
564#endif
565
566bool
567AxisAlignedBox3::Includes(const AxisAlignedBox3 &b) const
568{
569  return (b.mMin.x >= mMin.x &&
570      b.mMin.y >= mMin.y &&
571      b.mMin.z >= mMin.z &&
572      b.mMax.x <= mMax.x &&
573      b.mMax.y <= mMax.y &&
574      b.mMax.z <= mMax.z);
575
576}
577
578
579// compute the coordinates of one vertex of the box
580Vector3
581AxisAlignedBox3::GetVertex(int xAxis, int yAxis, int zAxis) const
582{
583  Vector3 p;
584  if (xAxis)
585    p.x = mMax.x;
586  else
587    p.x = mMin.x;
588
589  if (yAxis)
590    p.y = mMax.y;
591  else
592    p.y = mMin.y;
593
594  if (zAxis)
595    p.z = mMax.z;
596  else
597    p.z = mMin.z;
598  return p;
599}
600
601// compute the vertex for number N = <0..7>, N = 4.x + 2.y + z, where
602// x,y,z are either 0 or 1; (0 .. min coordinate, 1 .. max coordinate)
603void
[191]604AxisAlignedBox3::GetVertex(const int N, Vector3 &vertex) const
[162]605{
606  switch (N) {
607    case 0: vertex = mMin; break;
608    case 1: vertex.SetValue(mMin.x, mMin.y, mMax.z); break;
609    case 2: vertex.SetValue(mMin.x, mMax.y, mMin.z); break;
610    case 3: vertex.SetValue(mMin.x, mMax.y, mMax.z); break;
611    case 4: vertex.SetValue(mMax.x, mMin.y, mMin.z); break;
612    case 5: vertex.SetValue(mMax.x, mMin.y, mMax.z); break;
613    case 6: vertex.SetValue(mMax.x, mMax.y, mMin.z); break;
614    case 7: vertex = mMax; break;
615    default: {
616      FATAL << "ERROR in AxisAlignedBox3::GetVertex N=" << N <<  "\n";
617      FATAL_ABORT;
618    }
619  }
620}
621
622// Returns region 0 .. 26 ; R = 9*x + 3*y + z ; (x,y,z) \in {0,1,2}
623int
624AxisAlignedBox3::GetRegionID(const Vector3 &point) const
625{
626  int ID = 0;
627
628  if (point.z >= mMin.z) {
629    if (point.z <= mMax.z)
630      ID += 1; // inside the two boundary planes
631    else
632      ID += 2; // outside
633  }
634
635  if (point.y >= mMin.y) {
636    if (point.y <= mMax.y)
637      ID += 3; // inside the two boundary planes
638    else
639      ID += 6; // outside
640  }
641
642  if (point.x >= mMin.x) {
643    if (point.x <= mMax.x)
644      ID += 9; // inside the two boundary planes
645    else
646      ID += 18; // outside
647  }
648  return ID;
649}
650
651
652
653// computes if a given box (smaller) lies at least in one
654// projection whole in the box (larger) = this encompasses given
655// ;-----;
656// | ;-; |
657// | '-' |
658// '-----'
659int
660AxisAlignedBox3::IsPiercedByBox(const AxisAlignedBox3 &box, int &axis) const
661{
662  // test on x-axis
663  if ( (mMax.x < box.mMin.x) ||
664       (mMin.x > box.mMax.x) )
665    return 0; // the boxes do not overlap at all at x-axis
666  if ( (box.mMin.y > mMin.y) &&
667       (box.mMax.y < mMax.y) &&
668       (box.mMin.z > mMin.z) &&
669       (box.mMax.z < mMax.z) ) {
670    axis = 0;
671    return 1; // the boxes overlap in x-axis
672  }
673  // test on y-axis
674  if ( (mMax.y < box.mMin.y) ||
675       (mMin.y > box.mMax.y) )
676    return 0; // the boxes do not overlap at all at y-axis
677  if ( (box.mMin.x > mMin.x) &&
678       (box.mMax.x < mMax.x) &&
679       (box.mMin.z > mMin.z) &&
680       (box.mMax.z < mMax.z) ) {
681    axis = 1;
682    return 1; // the boxes overlap in y-axis
683  }
684  // test on z-axis
685  if ( (mMax.z < box.mMin.z) ||
686       (mMin.z > box.mMax.z) )
687    return 0; // the boxes do not overlap at all at y-axis
688  if ( (box.mMin.x > mMin.x) &&
689       (box.mMax.x < mMax.x) &&
690       (box.mMin.y > mMin.y) &&
691       (box.mMax.y < mMax.y) ) {
692    axis = 2;
693    return 1; // the boxes overlap in z-axis
694  }
695  return 0;
696}
697
698float
699AxisAlignedBox3::SurfaceArea() const
700{
701  Vector3 ext = mMax - mMin;
[302]702 
[469]703  return 2.0f * (ext.x * ext.y +
[162]704                ext.x * ext.z +
705                ext.y * ext.z);
706}
707
708const int AxisAlignedBox3::bvertices[27][9] =
709{                   // region number.. position
710  {5,1,3,2,6,4,-1,-1,-1},      // 0 .. x=0 y=0 z=0
711  {4,5,1,3,2,0,-1,-1,-1},      // 1 .. x=0 y=0 z=1
712  {4,5,7,3,2,0,-1,-1,-1},      // 2 .. x=0 y=0 z=2
713
714  {0,1,3,2,6,4,-1,-1,-1},      // 3 .. x=0 y=1 z=0
715  {0,1,3,2,-1,-1,-1,-1,-1},    // 4 .. x=0 y=1 z=1
716  {1,5,7,3,2,0,-1,-1,-1},      // 5 .. x=0 y=1 z=2
717
718  {0,1,3,7,6,4,-1,-1,-1},      // 6 .. x=0 y=2 z=0
719  {0,1,3,7,6,2,-1,-1,-1},      // 7 .. x=0 y=2 z=1
720  {1,5,7,6,2,0,-1,-1,-1},      // 8 .. x=0 y=2 z=2
721
722  // the regions number <9,17>
723  {5,1,0,2,6,4,-1,-1,-1},      // 9  .. x=1 y=0 z=0
724  {5,1,0,4,-1,-1,-1,-1,-1},    // 10 .. x=1 y=0 z=1
725  {7,3,1,0,4,5,-1,-1,-1},      // 11 .. x=1 y=0 z=2
726
727  {4,0,2,6,-1,-1,-1,-1,-1},    // 12 .. x=1 y=1 z=0
728  {0,2,3,1,5,4,6,7,-1},        // 13 .. x=1 y=1 z=1 .. inside the box
729  {1,5,7,3,-1,-1,-1,-1,-1},    // 14 .. x=1 y=1 z=2
730
731  {4,0,2,3,7,6,-1,-1,-1},      // 15 .. x=1 y=2 z=0
732  {6,2,3,7,-1,-1,-1,-1,-1},    // 16 .. x=1 y=2 z=1
733  {1,5,7,6,2,3,-1,-1,-1},      // 17 .. x=1 y=2 z=2
734
735  // the regions number <18,26>
736  {1,0,2,6,7,5,-1,-1,-1},      // 18 .. x=2 y=0 z=0
737  {1,0,4,6,7,5,-1,-1,-1},      // 19 .. x=2 y=0 z=1
738  {0,4,6,7,3,1,-1,-1,-1},      // 20 .. x=2 y=0 z=2
739
740  {4,0,2,6,7,5,-1,-1,-1},      // 21 .. x=2 y=1 z=0
741  {5,4,6,7,-1,-1,-1,-1,-1},    // 22 .. x=2 y=1 z=1
742  {5,4,6,7,3,1,-1,-1,-1},      // 23 .. x=2 y=1 z=2
743
744  {4,0,2,3,7,5,-1,-1,-1},      // 24 .. x=2 y=2 z=0
745  {5,4,6,2,3,7,-1,-1,-1},      // 25 .. x=2 y=2 z=1
746  {5,4,6,2,3,1,-1,-1,-1},      // 26 .. x=2 y=2 z=2
747};
748
749// the visibility of boundary faces from a given region
750// one to three triples: (axis, min_vertex, max_vertex), axis==-1(terminator)
751const int AxisAlignedBox3::bfaces[27][10] =
752{                              // region number .. position
753  {0,0,3,1,0,5,2,0,6,-1},      // 0 .. x=0 y=0 z=0
754  {0,0,3,1,0,5,-1,-1,-1,-1},   // 1 .. x=0 y=0 z=1
755  {0,0,3,1,0,5,2,1,7,-1},      // 2 .. x=0 y=0 z=2
756
757  {0,0,3,2,0,6,-1,-1,-1,-1},   // 3 .. x=0 y=1 z=0
758  {0,0,3,-1,-1,-1,-1,-1,-1,-1},// 4 .. x=0 y=1 z=1
759  {0,0,3,2,1,7,-1,-1,-1,-1},   // 5 .. x=0 y=1 z=2
760
761  {0,0,3,1,2,7,2,0,6,-1},      // 6 .. x=0 y=2 z=0
762  {0,0,3,1,2,7,-1,-1,-1,-1},   // 7 .. x=0 y=2 z=1
763  {0,0,3,1,2,7,2,1,7,-1},      // 8 .. x=0 y=2 z=2
764
765  // the regions number <9,17>
766  {1,0,5,2,0,6,-1,-1,-1,-1},   // 9  .. x=1 y=0 z=0
767  {1,0,5,-1,-1,-1,-1,-1,-1,-1},// 10 .. x=1 y=0 z=1
768  {1,0,5,2,1,7,-1,-1,-1,-1},   // 11 .. x=1 y=0 z=2
769
770  {2,0,6,-1,-1,-1,-1,-1,-1,-1},// 12 .. x=1 y=1 z=0
771  {-1,-1,-1,-1,-1,-1,-1,-1,-1,-1},// 13 .. x=1 y=1 z=1 .. inside the box
772  {2,1,7,-1,-1,-1,-1,-1,-1,-1},// 14 .. x=1 y=1 z=2
773
774  {1,2,7,2,0,6,-1,-1,-1,-1},   // 15 .. x=1 y=2 z=0
775  {1,2,7,-1,-1,-1,-1,-1,-1,-1},// 16 .. x=1 y=2 z=1
776  {1,2,7,2,1,7,-1,-1,-1,-1},   // 17 .. x=1 y=2 z=2
777
778  // the region number <18,26>
779  {0,4,7,1,0,5,2,0,6,-1},      // 18 .. x=2 y=0 z=0
780  {0,4,7,1,0,5,-1,-1,-1,-1},   // 19 .. x=2 y=0 z=1
781  {0,4,7,1,0,5,2,1,7,-1},      // 20 .. x=2 y=0 z=2
782
783  {0,4,7,2,0,6,-1,-1,-1,-1},   // 21 .. x=2 y=1 z=0
784  {0,4,7,-1,-1,-1,-1,-1,-1,-1},// 22 .. x=2 y=1 z=1
785  {0,4,7,2,1,7,-1,-1,-1,-1},   // 23 .. x=2 y=1 z=2
786
787  {0,4,7,1,2,7,2,0,6,-1},      // 24 .. x=2 y=2 z=0
788  {0,4,7,1,2,7,-1,-1,-1,-1},   // 25 .. x=2 y=2 z=1
789  {0,4,7,1,2,7,2,1,7,-1},      // 26 .. x=2 y=2 z=2
790};
791
792// the correct corners indexing from entry face to exit face
793// first index determines entry face, second index exit face, and
794// the two numbers (indx, inc) determines: ind = the index on the exit
795// face, when starting from the vertex 0 on entry face, 'inc' is
796// the increment when we go on entry face in order 0,1,2,3 to create
797// convex shaft with the rectangle on exit face. That is, inc = -1 or 1.
798const int AxisAlignedBox3::pairFaceRects[6][6][2] = {
799  { // entry face = 0
800    {-1,0}, // exit face 0 .. no meaning
801    {0,-1}, // 1
802    {0,-1}, // 2
803    {0,1}, // 3 .. opposite face
804    {3,1}, // 4
805    {1,1} // 5
806  },
807  { // entry face = 1
808    {0,-1}, // exit face 0
809    {-1,0}, // 1 .. no meaning
810    {0,-1}, // 2
811    {1,1}, // 3
812    {0,1}, // 4 .. opposite face
813    {3,1} // 5
814  },
815  { // entry face = 2
816    {0,-1}, // 0
817    {0,-1}, // 1
818    {-1,0}, // 2 .. no meaning
819    {3,1}, // 3
820    {1,1}, // 4
821    {0,1} // 5 .. opposite face
822  },
823  { // entry face = 3
824    {0,1}, // 0 .. opposite face
825    {3,-1}, // 1
826    {1,1}, // 2
827    {-1,0}, // 3  .. no meaning
828    {0,-1}, // 4
829    {0,-1} // 5
830  },
831  { // entry face = 4
832    {1,1}, // 0
833    {0,1}, // 1 .. opposite face
834    {3,1}, // 2
835    {0,-1}, // 3
836    {-1,0}, // 4 .. no meaning
837    {0,-1} // 5
838  },
839  { // entry face = 5
840    {3,-1}, // 0
841    {1,1}, // 1
842    {0,1}, // 2  .. opposite face
843    {0,-1}, // 3
844    {0,-1}, // 4
845    {-1,0} // 5 .. no meaning
846  }
847};
848
849
850// ------------------------------------------------------------
851// The vertices that form CLOSEST points with respect to the region
852// for all the regions possible, number of regions is 3^3 = 27,
853// since two parallel sides of bbox forms three disjoint spaces.
854// The vertices are given in anti-clockwise order, stopped by value 15,
855// at most 8 points, at least 1 point.
856// The table includes the closest 1/2/4/8 points, followed possibly
857// by the set of coordinates that should be used for testing for
858// the proximity queries. The coordinates to be tested are described by
859// the pair (a,b), when a=0, we want to test min vector of the box,
860//                 when a=1, we want to test max vector of the box
861//                 b=0,1,2 corresponds to the axis (0=x,1=y,2=z)
862// The sequence is ended by 15, number -1 is used as the separator
863// between the vertices and coordinates.
864const int
865AxisAlignedBox3::cvertices[27][9] =
866{                   // region number.. position
867  {0,15,15,15,15,15,15,15,15}, // 0 .. x=0 y=0 z=0 D one vertex
868  {0,1,-1,0,0,0,1,15,15},      // 1 .. x=0 y=0 z=1 D two vertices foll. by 2
869  {1,15,15,15,15,15,15,15,15}, // 2 .. x=0 y=0 z=2 D one vertex
870
871  {0,2,-1,0,0,0,2,15,15},      // 3 .. x=0 y=1 z=0 D
872  {0,1,3,2,-1,0,0,15,15},      // 4 .. x=0 y=1 z=1 D
873  {1,3,-1,0,0,1,2,15,15},      // 5 .. x=0 y=1 z=2 D
874
875  {2,15,15,15,15,15,15,15,15}, // 6 .. x=0 y=2 z=0 D
876  {2,3,-1,0,0,1,1,15,15},      // 7 .. x=0 y=2 z=1 D
877  {3,15,15,15,15,15,15,15,15}, // 8 .. x=0 y=2 z=2 D
878
879  // the regions number <9,17>
880  {0,4,-1,0,1,0,2,15,15},      // 9  .. x=1 y=0 z=0 D
881  {5,1,0,4,-1,0,1,15,15},      // 10 .. x=1 y=0 z=1 D
882  {1,5,-1,0,1,1,2,15,15},      // 11 .. x=1 y=0 z=2 D
883
884  {4,0,2,6,-1,0,2,15,15},      // 12 .. x=1 y=1 z=0 D
885  {0,2,3,1,5,4,6,7,15},        // 13 .. x=1 y=1 z=1 .. inside the box
886  {1,5,7,3,-1,1,2,15,15},      // 14 .. x=1 y=1 z=2 D
887
888  {6,2,-1,0,2,1,1,15,15},      // 15 .. x=1 y=2 z=0 D
889  {6,2,3,7,-1,1,1,15,15},      // 16 .. x=1 y=2 z=1 D
890  {3,7,-1,1,1,1,2,15,15},      // 17 .. x=1 y=2 z=2 D
891
892  // the regions number <18,26>
893  {4,15,15,15,15,15,15,15,15}, // 18 .. x=2 y=0 z=0 D
894  {4,5,-1,0,1,1,0,15,15},      // 19 .. x=2 y=0 z=1 D
895  {5,15,15,15,15,15,15,15,15}, // 20 .. x=2 y=0 z=2 D
896
897  {4,6,-1,0,2,1,0,15,15},      // 21 .. x=2 y=1 z=0 D
898  {5,4,6,7,-1,1,0,15,15},      // 22 .. x=2 y=1 z=1 D
899  {7,5,-1,1,0,1,2,15,15},      // 23 .. x=2 y=1 z=2 D
900
901  {6,15,15,15,15,15,15,15,15}, // 24 .. x=2 y=2 z=0 D
902  {6,7,-1,1,0,1,1,15,15},      // 25 .. x=2 y=2 z=1 D
903  {7,15,15,15,15,15,15,15,15}, // 26 .. x=2 y=2 z=2 D
904};
905
906// Table for Sphere-AABB intersection based on the region knowledge
907// Similar array to previous cvertices, but we omit the surfaces
908// which are not necessary for testing. First are vertices,
909// they are finished with -1. Second, there are indexes in
910// the pair (a,b), when a=0, we want to test min vector of the box,
911//                 when a=1, we want to test max vector of the box
912//                 b=0,1,2 corresponds to the axis (0=x,1=y,2=z)
913//
914// So either we check the vertices or only the distance in specified
915// dimensions. There are at all four possible cases:
916//
917//   1) we check one vertex - then sequence start with non-negative index
918//      and is finished with 15
919//   2) we check two coordinates of min/max vector describe by the pair
920//         (a,b) .. a=min/max(0/1) b=x/y/z (0/1/2), sequence starts with 8
921//      and finishes with 15
922//   3) we check only one coordinate of min/max, as for 2), sequence start
923//      with 9 and ends with 15
924//   4) Position 13 - sphere is inside the box, intersection always exist
925//        the sequence start with 15 .. no further testing is necessary
926//        in this case
927const int
928AxisAlignedBox3::csvertices[27][6] =
929{                   // region number.. position
930  {0,15,15,15,15,15},  // 0 .. x=0 y=0 z=0 D vertex only
931  {8,0,0,0,1,15},      // 1 .. x=0 y=0 z=1 D two coords.
932  {1,15,15,15,15,15},  // 2 .. x=0 y=0 z=2 D vertex only
933
934  {8,0,0,0,2,15},      // 3 .. x=0 y=1 z=0 D two coords
935  {9,0,0,15,15,15},    // 4 .. x=0 y=1 z=1 D one coord
936  {8,0,0,1,2,15},      // 5 .. x=0 y=1 z=2 D two coords.
937
938  {2,15,15,15,15,15},  // 6 .. x=0 y=2 z=0 D one vertex
939  {8,0,0,1,1,15},      // 7 .. x=0 y=2 z=1 D two coords
940  {3,15,15,15,15,15},  // 8 .. x=0 y=2 z=2 D one vertex
941
942  // the regions number <9,17>
943  {8,0,1,0,2,15},      // 9  .. x=1 y=0 z=0 D two coords
944  {9,0,1,15,15,15},    // 10 .. x=1 y=0 z=1 D one coord
945  {8,0,1,1,2,15},      // 11 .. x=1 y=0 z=2 D two coords
946
947  {9,0,2,15,15,15},    // 12 .. x=1 y=1 z=0 D one coord
948  {15,15,15,15,15,15}, // 13 .. x=1 y=1 z=1 inside the box, special case/value
949  {9,1,2,15,15,15},    // 14 .. x=1 y=1 z=2 D one corrd
950
951  {8,0,2,1,1,15},      // 15 .. x=1 y=2 z=0 D two coords
952  {9,1,1,15,15},       // 16 .. x=1 y=2 z=1 D one coord
953  {8,1,1,1,2,15},      // 17 .. x=1 y=2 z=2 D two coords
954
955  // the regions number <18,26>
956  {4,15,15,15,15,15},  // 18 .. x=2 y=0 z=0 D one vertex
957  {8,0,1,1,0,15},      // 19 .. x=2 y=0 z=1 D two coords
958  {5,15,15,15,15,15},  // 20 .. x=2 y=0 z=2 D one vertex
959
960  {8,0,2,1,0,15},      // 21 .. x=2 y=1 z=0 D two coords
961  {9,1,0,15,15,15},    // 22 .. x=2 y=1 z=1 D one coord
962  {8,1,0,1,2,15},      // 23 .. x=2 y=1 z=2 D two coords
963
964  {6,15,15,15,15,15},  // 24 .. x=2 y=2 z=0 D one vertex
965  {8,1,0,1,1,15},      // 25 .. x=2 y=2 z=1 D two coords
966  {7,15,15,15,15,15},  // 26 .. x=2 y=2 z=2 D one vertex
967};
968
969
970// The vertices that form all FARTHEST points with respect to the region
971// for all the regions possible, number of regions is 3^3 = 27,
972// since two parallel sides of bbox forms three disjoint spaces.
973// The vertices are given in anti-clockwise order, stopped by value 15,
974// at most 8 points, at least 1 point.
975// For testing, if the AABB is whole in the sphere, it is enough
976// to test only vertices, either 1,2,4, or 8.
977const int
978AxisAlignedBox3::fvertices[27][9] =
979{                   // region number.. position
980  {7,15,15,15,15,15,15,15,15}, // 0 .. x=0 y=0 z=0 D
981  {6,7,15,15,15,15,15,15,15},  // 1 .. x=0 y=0 z=1 D
982  {6,15,15,15,15,15,15,15,15}, // 2 .. x=0 y=0 z=2 D
983
984  {5,7,15,15,15,15,15,15,15},  // 3 .. x=0 y=1 z=0 D
985  {4,5,7,6,15,15,15,15,15},    // 4 .. x=0 y=1 z=1 D
986  {4,6,15,15,15,15,15,15,15},  // 5 .. x=0 y=1 z=2 D
987
988  {5,15,15,15,15,15,15,15,15}, // 6 .. x=0 y=2 z=0 D
989  {4,5,15,15,15,15,15,15,15},  // 7 .. x=0 y=2 z=1 D
990  {4,15,15,15,15,15,15,15,15}, // 8 .. x=0 y=2 z=2 D
991
992  // the regions number <9,17>
993  {3,7,15,15,15,15,15,15,15},  // 9  .. x=1 y=0 z=0 D
994  {7,3,2,6,15,15,15,15,15},    // 10 .. x=1 y=0 z=1 D
995  {2,6,15,15,15,15,15,15,15},  // 11 .. x=1 y=0 z=2 D
996
997  {5,1,3,7,15,15,15,15,15},    // 12 .. x=1 y=1 z=0 D
998  {0,7,1,6,3,4,5,2,15},        // 13 .. x=1 y=1 z=1 .. inside the box
999  {0,4,6,2,15,15,15,15,15},    // 14 .. x=1 y=1 z=2 D
1000
1001  {5,1,15,15,15,15,15,15,15},  // 15 .. x=1 y=2 z=0 D
1002  {4,0,1,5,15,15,15,15,15},    // 16 .. x=1 y=2 z=1 D
1003  {4,0,15,15,15,15,15,15,15},  // 17 .. x=1 y=2 z=2 D
1004
1005  // the regions number <18,26>
1006  {3,15,15,15,15,15,15,15,15}, // 18 .. x=2 y=0 z=0 D
1007  {2,3,15,15,15,15,15,15,15},  // 19 .. x=2 y=0 z=1 D
1008  {2,15,15,15,15,15,15,15,15}, // 20 .. x=2 y=0 z=2 D
1009
1010  {1,3,15,15,15,15,15,15,15},  // 21 .. x=2 y=1 z=0 D
1011  {1,0,2,3,15,15,15,15,15},    // 22 .. x=2 y=1 z=1 D
1012  {2,0,15,15,15,15,15,15,15},  // 23 .. x=2 y=1 z=2 D
1013
1014  {1,15,15,15,15,15,15,15,15}, // 24 .. x=2 y=2 z=0 D
1015  {0,1,15,15,15,15,15,15,15},  // 25 .. x=2 y=2 z=1 D
1016  {0,15,15,15,15,15,15,15,15}, // 26 .. x=2 y=2 z=2 D
1017};
1018
1019// Similar table as above, farthest points, but only the ones
1020// necessary for testing the intersection problem. If we do
1021// not consider the case 13, center of the sphere is inside the
1022// box, then we can always test at most 2 box vertices to say whether
1023// the whole box is inside the sphere.
1024// The number of vertices is minimized using some assumptions
1025// about the ortogonality of vertices and sphere properties.
1026const int
1027AxisAlignedBox3::fsvertices[27][9] =
1028{                   // region number.. position
1029  {7,15,15,15,15,15,15,15,15}, // 0 .. x=0 y=0 z=0 D 1 vertex
1030  {6,7,15,15,15,15,15,15,15},  // 1 .. x=0 y=0 z=1 D 2 vertices
1031  {6,15,15,15,15,15,15,15,15}, // 2 .. x=0 y=0 z=2 D 1 vertex
1032
1033  {5,7,15,15,15,15,15,15,15},  // 3 .. x=0 y=1 z=0 D 2 vertices
1034  {4,7,15,5,6,15,15,15,15},    // 4 .. x=0 y=1 z=1 D 4/2 vertices
1035  {4,6,15,15,15,15,15,15,15},  // 5 .. x=0 y=1 z=2 D 2 vertices
1036
1037  {5,15,15,15,15,15,15,15,15}, // 6 .. x=0 y=2 z=0 D 1 vertex
1038  {4,5,15,15,15,15,15,15,15},  // 7 .. x=0 y=2 z=1 D 2 vertices
1039  {4,15,15,15,15,15,15,15,15}, // 8 .. x=0 y=2 z=2 D 1 vertex
1040
1041  // the regions number <9,17>
1042  {3,7,15,15,15,15,15,15,15},  // 9  .. x=1 y=0 z=0 D 2 vertices
1043  {7,2,15,3,6,15,15,15,15},    // 10 .. x=1 y=0 z=1 D 4/2 vertices
1044  {2,6,15,15,15,15,15,15,15},  // 11 .. x=1 y=0 z=2 D 2 vertices
1045
1046  {5,3,15,1,7,15,15,15,15},    // 12 .. x=1 y=1 z=0 D 4/2 vertices
1047  {0,7,1,6,3,4,5,2,15},        // 13 .. x=1 y=1 z=1 .. inside the box
1048  {0,6,15,4,2,15,15,15,15},    // 14 .. x=1 y=1 z=2 D 4/2 vertices
1049
1050  {5,1,15,15,15,15,15,15,15},  // 15 .. x=1 y=2 z=0 D 2 vertices
1051  {4,1,15,0,5,15,15,15,15},    // 16 .. x=1 y=2 z=1 D 4/2 vertices
1052  {4,0,15,15,15,15,15,15,15},  // 17 .. x=1 y=2 z=2 D 2 vertices
1053
1054  // the regions number <18,26>
1055  {3,15,15,15,15,15,15,15,15}, // 18 .. x=2 y=0 z=0 D 1 vertex
1056  {2,3,15,15,15,15,15,15,15},  // 19 .. x=2 y=0 z=1 D 2 vertices
1057  {2,15,15,15,15,15,15,15,15}, // 20 .. x=2 y=0 z=2 D 1 vertex
1058
1059  {1,3,15,15,15,15,15,15,15},  // 21 .. x=2 y=1 z=0 D 2 vertices
1060  {1,2,15,0,3,15,15,15,15},    // 22 .. x=2 y=1 z=1 D 4/2 vertices
1061  {2,0,15,15,15,15,15,15,15},  // 23 .. x=2 y=1 z=2 D 2 vertices
1062
1063  {1,15,15,15,15,15,15,15,15}, // 24 .. x=2 y=2 z=0 D 1 vertex
1064  {0,1,15,15,15,15,15,15,15},  // 25 .. x=2 y=2 z=1 D 2 vertices
1065  {0,15,15,15,15,15,15,15,15}, // 26 .. x=2 y=2 z=2 D 1 vertex
1066};
1067
[312]1068
[162]1069// The fast computation of arctangent .. the maximal error is less
1070// than 4.1 degrees, according to Graphics GEMSII, 1991, pages 389--391
1071// Ron Capelli: "Fast approximation to the arctangent"
1072float
1073atan22(const float& y)
1074{
1075  const float x = 1.0;
1076  const float c = (float)(M_PI * 0.25);
1077 
1078  if (y < 0.0) {
1079    if (y < -1.0)
[469]1080      return c * (-2.0f + x / y); // for angle in <-PI/2, -PI/4)
[162]1081    else
1082      return c * (y / x); // for angle in <-PI/4 , 0>
1083  }
1084  else {
1085    if (y > 1.0)
[469]1086      return c * (2.0f - x / y); // for angle in <PI/4, PI/2>
[162]1087    else
1088      return c * (y / x); // for angle in <0, PI/2>
1089  }
1090}
1091
1092
1093float
1094AxisAlignedBox3::ProjectToSphereSA(const Vector3 &viewpoint, int *tcase) const
1095{
1096  int id = GetRegionID(viewpoint);
1097  *tcase = id;
1098
1099  // spherical projection .. SA represents solid angle
1100  if (id == 13) // .. inside the box
1101    return (float)(4.0*M_PI); // the whole sphere
1102  float SA = 0.0; // inital value
1103   
1104  int i = 0; // the pointer in the array of vertices
1105  while (bfaces[id][i] >= 0) {
1106    int axisO = bfaces[id][i++];
1107    int minvIdx = bfaces[id][i++];
1108    int maxvIdx = bfaces[id][i++];
1109    Vector3 vmin, vmax;
1110    GetVertex(minvIdx, vmin);
1111    GetVertex(maxvIdx, vmax);
1112    float h = fabs(vmin[axisO] - viewpoint[axisO]);
1113    int axis = (axisO + 1) % 3; // next axis
1114    float a = (vmin[axis] - viewpoint[axis]) / h; // minimum for v-range
1115    float b = (vmax[axis] - viewpoint[axis]) / h; // maximum for v-range
1116    //if (a > b) {
1117    //  FATAL << "ProjectToSphereSA::Error a > b\n";
1118    //  FATAL_ABORT;
1119    //}
1120    //if (vmin[axisO] != vmax[axisO]) {
1121    //  FATAL << "ProjectToSphereSA::Error a-axis != b-axis\n";
1122    //  FATAL_ABORT;     
1123    //}
1124    axis = (axisO + 2) % 3; // next second axis       
1125    float c = (vmin[axis] - viewpoint[axis]) / h; // minimum for u-range
1126    float d = (vmax[axis] - viewpoint[axis]) / h; // maximum for u-range
1127    //if (c > d) {
1128    //  FATAL << "ProjectToSphereSA::Error c > d\n";
1129    //  FATAL_ABORT;
1130    //}
[469]1131    SA +=atan22(d*b/sqrt(b*b + d*d + 1.0f)) - atan22(b*c/sqrt(b*b + c*c + 1.0f))
1132      - atan22(d*a/sqrt(a*a + d*d + 1.0f)) + atan22(a*c/sqrt(a*a + c*c + 1.0f));
[162]1133  }
1134
1135#if 0
1136  if ((SA > 2.0*M_PI) ||
1137      (SA < 0.0)) {
1138    FATAL << "The solid angle has strange value: ";
1139    FATAL << "SA = "<< SA << endl;
1140    FATAL_ABORT;
1141  }
1142#endif
1143 
1144  return SA;
1145}
1146
1147// Projects the box to a plane given a normal vector only and
1148// computes the surface area of the projected silhouette
1149// no clipping of the box is performed.
1150float
1151AxisAlignedBox3::ProjectToPlaneSA(const Vector3 &normal) const
1152{
1153  Vector3 size = Size();
1154
1155  // the surface area of the box to a yz-plane - perpendicular to x-axis
1156  float sax = size.y * size.z;
1157
1158  // the surface area of the box to a zx-plane - perpendicular to y-axis
1159  float say = size.z * size.x;
1160
1161  // the surface area of the box to a xy-plane - perpendicular to z-axis
1162  float saz = size.x * size.y;
1163 
1164  return sax * fabs(normal.x) + say * fabs(normal.y) + saz * fabs(normal.z);
1165}
1166
1167
1168
1169// This definition allows to be a point when answering true
1170bool
1171AxisAlignedBox3::IsCorrectAndNotPoint() const
1172{
1173  if ( (mMin.x > mMax.x) ||
1174       (mMin.y > mMax.y) ||
1175       (mMin.z > mMax.z) )
1176    return false; // box is not formed
1177
1178  if ( (mMin.x == mMax.x) &&
1179       (mMin.y == mMax.y) &&
1180       (mMin.z == mMax.z) )
1181    return false; // degenerates to a point
1182 
1183  return true;
1184}
1185
1186// This definition allows to be a point when answering true
1187bool
1188AxisAlignedBox3::IsPoint() const
1189{
1190  if ( (mMin.x == mMax.x) &&
1191       (mMin.y == mMax.y) &&
1192       (mMin.z == mMax.z) )
1193    return true; // degenerates to a point
1194 
1195  return false;
1196}
1197
1198// This definition requires shape of non-zero volume
1199bool
1200AxisAlignedBox3::IsSingularOrIncorrect() const
1201{
1202  if ( (mMin.x >= mMax.x) ||
1203       (mMin.y >= mMax.y) ||
1204       (mMin.z >= mMax.z) )
1205    return true; // box is not formed
1206
1207  return false; // has non-zero volume
1208}
1209
1210// returns true, when the sphere specified by the origin and radius
1211// fully contains the box
1212bool
1213AxisAlignedBox3::IsFullyContainedInSphere(const Vector3 &center, float rad) const
1214{
1215  int region = GetRegionID(center);
1216  float rad2 = rad*rad;
1217
1218  // vertex of the box
1219  Vector3 vertex;
1220
1221  int i = 0;
1222  for (i = 0 ; ; i++) {
1223    int a = fsvertices[region][i];
1224    if (a == 15)
1225      return true; // if was not false untill now, it must be contained
1226
1227    assert( (a>=0) && (a<8) );
1228   
1229    // normal vertex
1230    GetVertex(a, vertex);
1231
1232    if (SqrMagnitude(vertex - center) > rad2)
1233      return false;   
1234  } // for
1235 
1236}
1237
1238// returns true, when the volume of the sphere and volume of the
1239// axis aligned box has no intersection
1240bool
1241AxisAlignedBox3::HasNoIntersectionWithSphere(const Vector3 &center, float rad) const
1242{
1243  int region = GetRegionID(center);
1244  float rad2 = rad*rad;
1245
1246  // vertex of the box
1247  Vector3 vertex;
1248
1249  switch (csvertices[region][0]) {
1250  case 8: {
1251    // test two coordinates described within the field
1252    int face = 3*csvertices[region][1] + csvertices[region][2];
1253    float dist = GetExtent(face) - center[csvertices[region][2]];
1254    dist *= dist;
1255    face = 3 * (csvertices[region][3]) + csvertices[region][4];
1256    float dist2 = GetExtent(face) - center[csvertices[region][4]];
1257    dist += (dist2 * dist2);
1258    if (dist > rad2)
1259      return true; // no intersection is possible
1260  }
1261  case 9: {
1262    // test one coordinate described within the field
1263    int face = 3*csvertices[region][1] + csvertices[region][2];
1264    float dist = fabs(GetExtent(face) - center[csvertices[region][2]]);
1265    if (dist > rad)
1266      return true; // no intersection is possible
1267  }
1268  case 15:
1269    return false; // box and sphere surely has intersection
1270  default: {
1271    // test using normal vertices
1272    assert( (csvertices[region][0]>=0) && (csvertices[region][0]<8) );
1273
1274    // normal vertex
1275    GetVertex(csvertices[region][0], vertex);
1276
1277    if (SqrMagnitude(vertex - center) > rad2)
1278      return true; // no intersectino is possible   
1279    }
1280  } // switch
1281
1282  return false; // partial or full containtment
1283}
1284
1285#if 0
1286// Given the sphere, determine the mutual position between the
1287// sphere and box
1288
1289// SOME BUG IS INSIDE !!!! V.H. 25/4/2001
1290int
1291AxisAlignedBox3::MutualPositionWithSphere(const Vector3 &center, float rad) const
1292{
1293  int region = GetRegionID(center);
1294  float rad2 = rad*rad;
1295
1296  // vertex of the box
1297  Vector3 vertex;
1298
1299  // first testing for  full containtment - whether sphere fully
1300  // contains the box
1301  int countInside = 0; // how many points were found inside
1302 
1303  int i = 0;
1304  for (i = 0 ; ; i++) {
1305    int a = fsvertices[region][i];
1306    if (a == 15)
1307      return 1; // the sphere fully contain the box
1308
1309    assert( (a>=0) && (a<8) );
1310   
1311    // normal vertex
1312    GetVertex(a, vertex);
1313
1314    if (SqrMagnitude(vertex - center) <= rad2)
1315      countInside++; // the number of vertices inside the sphere
1316    else {
1317      if (countInside)
1318        return 0; // partiall overlap has been found
1319      // the sphere does not fully contain the box .. only way to break
1320      // this loop and go for other testing
1321      break;
1322    }
1323  } // for
1324
1325  // now only box and sphere can partially overlap or no intersection
1326  switch (csvertices[region][0]) {
1327  case 8: {
1328    // test two coordinates described within the field
1329    int face = 3*csvertices[region][1] + csvertices[region][2];
1330    float dist = GetExtent(face) - center[csvertices[region][2]];
1331    dist *= dist;
1332    face = 3 * (csvertices[region][3]) + csvertices[region][4];
1333    float dist2 = GetExtent(face) - center[csvertices[region][4]];
1334    dist += (dist2 * dist2);
1335    if (dist > rad2 )
1336      return -1; // no intersection is possible
1337  }
1338  case 9: {
1339    // test one coordinate described within the field
1340    int face = 3*csvertices[region][1] + csvertices[region][2];
1341    float dist = fabs(GetExtent(face) - center[csvertices[region][2]]);
1342    if (dist > rad)
1343      return -1; // no intersection is possible
1344  }
1345  case 15:
1346    return 0 ; // partial overlap is now guaranteed
1347  default: {
1348    // test using normal vertices
1349    assert( (csvertices[region][0]>=0) && (csvertices[region][0]<8) );
1350   
1351    // normal vertex
1352    GetVertex(csvertices[region][0], vertex);
1353
1354    if (SqrMagnitude(vertex - center) > rad2)
1355      return -1; // no intersection is possible   
1356    }
1357  } // switch
1358
1359  return 0; // partial intersection is guaranteed
1360}
1361#else
1362
1363// Some maybe smarter version, extendible easily to d-dimensional
1364// space!
1365// Given a sphere described by the center and radius,
1366// the fullowing function returns:
1367//   -1 ... the sphere and the box are completely separate
1368//    0 ... the sphere and the box only partially overlap
1369//    1 ... the sphere contains fully the box
1370//  Note: the case when box fully contains the sphere is not reported
1371//        since it was not required.
1372int
1373AxisAlignedBox3::MutualPositionWithSphere(const Vector3 &center, float rad) const
1374{
1375//#define SPEED_UP
1376
1377#ifndef SPEED_UP
1378  // slow version, instructively written 
1379#if 0
1380  // does it make sense to test
1381  // checking the sides of the box for possible non-intersection
1382  if ( ((center.x + rad) < mMin.x) ||
1383       ((center.x - rad) > mMax.x) ||
1384       ((center.y + rad) < mMin.y) ||
1385       ((center.y - rad) > mMax.y) ||
1386       ((center.z + rad) < mMin.z) ||
1387       ((center.z - rad) > mMax.z) ) {
1388    // cout << "r ";
1389    return -1;  // no overlap is possible
1390  }
1391#endif
1392 
1393  // someoverlap is possible, check the distance of vertices
1394  rad = rad*rad;
1395  float sumMin = 0;
1396  // Try to minimize the function of a distance
1397  // from the sphere center
1398
1399  // for x-axis
1400  float minSqrX = sqr(mMin.x - center.x);
1401  float maxSqrX = sqr(mMax.x - center.x);
1402  if (center.x < mMin.x)
1403    sumMin = minSqrX;
1404  else
1405    if (center.x > mMax.x)
1406      sumMin = maxSqrX;
1407
1408  // for y-axis
1409  float minSqrY = sqr(mMin.y - center.y);
1410  float maxSqrY = sqr(mMax.y - center.y);
1411  if (center.y < mMin.y)
1412    sumMin += minSqrY;
1413  else
1414    if (center.y > mMax.y)
1415      sumMin += maxSqrY;
1416
1417  // for z-axis
1418  float minSqrZ = sqr(mMin.z - center.z);
1419  float maxSqrZ = sqr(mMax.z - center.z);
1420  if (center.z < mMin.z)
1421    sumMin += minSqrZ;
1422  else
1423    if (center.z > mMax.z)
1424      sumMin += maxSqrZ;
1425
1426  if (sumMin > rad)
1427    return -1; // no intersection between sphere and box
1428
1429  // try to find out the maximum distance between the
1430  // sphere center and vertices
1431  float sumMax = 0;
1432 
1433  if (minSqrX > maxSqrX)
1434    sumMax = minSqrX;
1435  else
1436    sumMax = maxSqrX;
1437 
1438  if (minSqrY > maxSqrY)
1439    sumMax += minSqrY;
1440  else
1441    sumMax += maxSqrY;
1442 
1443  if (minSqrZ > maxSqrZ)
1444    sumMax += minSqrZ;
1445  else
1446    sumMax += maxSqrZ;
1447 
1448  // sumMin < rad
1449  if (sumMax < rad)
1450    return 1; // the sphere contains the box completely
1451
1452  // partial intersection, part of the box is outside the sphere
1453  return 0;
1454#else
1455
1456  // Optimized version of the test
1457 
1458#ifndef __VECTOR_HACK
1459#error "__VECTOR_HACK for Vector3 was not defined"
1460#endif
1461
1462  // some overlap is possible, check the distance of vertices
1463  rad = rad*rad;
1464  float sumMin = 0;
1465  float sumMax = 0;
1466  // Try to minimize the function of a distance
1467  // from the sphere center
1468
1469  const float *minp = &(min[0]);
1470  const float *maxp = &(max[0]);
1471  const float *pcenter = &(center[0]);
1472 
1473  // for x-axis
1474  for (int i = 0; i < 3; i++, minp++, maxp++, pcenter++) {
1475    float minsqr = sqr(*minp - *pcenter);
1476    float maxsqr = sqr(*maxp - *pcenter);
1477    if (*pcenter < *minp)
1478      sumMin += minsqr;
1479    else
1480      if (*pcenter > *maxp)
1481        sumMin += maxsqr;
1482    sumMax += (minsqr > maxsqr) ? minsqr : maxsqr;
1483  }
1484
1485  if (sumMin > rad)
1486    return -1; // no intersection between sphere and box
1487
1488  // sumMin < rad
1489  if (sumMax < rad)
1490    return 1; // the sphere contains the box completely
1491
1492  // partial intersection, part of the box is outside the sphere
1493  return 0;
1494#endif
1495}
1496#endif
1497
1498// Given the cube describe by the center and the half-sie,
1499// determine the mutual position between the cube and the box
1500int
1501AxisAlignedBox3::MutualPositionWithCube(const Vector3 &center, float radius) const
1502{
1503  // the cube is described by the center and the distance to the any face
1504  // along the axes
1505
1506  // Note on efficiency!
1507  // Can be quite optimized using tables, but I do not have time
1508  // V.H. 18/11/2001
1509
1510  AxisAlignedBox3 a =
1511    AxisAlignedBox3(Vector3(center.x - radius, center.y - radius, center.z - radius),
1512           Vector3(center.x + radius, center.y + radius, center.z + radius));
1513
1514  if (a.Includes(*this))
1515    return 1; // cube contains the box
1516
1517  if (OverlapS(a,*this))
1518    return 0; // cube partially overlap the box
1519
1520  return -1; // completely separate 
1521}
1522
1523void
1524AxisAlignedBox3::GetSqrDistances(const Vector3 &point,
1525                                 float &minDistance,
1526                                 float &maxDistance
[209]1527                                 ) const
[162]1528{
1529
1530 
1531#ifndef __VECTOR_HACK
1532#error "__VECTOR_HACK for Vector3 was not defined"
1533#endif
1534
1535  // some overlap is possible, check the distance of vertices
1536  float sumMin = 0;
1537  float sumMax = 0;
1538
1539  // Try to minimize the function of a distance
1540  // from the sphere center
1541
1542  const float *minp = &(mMin[0]);
1543  const float *maxp = &(mMax[0]);
1544  const float *pcenter = &(point[0]);
1545 
1546  // for x-axis
1547  for (int i = 0; i < 3; i++, minp++, maxp++, pcenter++) {
1548    float minsqr = sqr(*minp - *pcenter);
1549    float maxsqr = sqr(*maxp - *pcenter);
1550    if (*pcenter < *minp)
1551      sumMin += minsqr;
1552    else
1553      if (*pcenter > *maxp)
1554        sumMin += maxsqr;
1555    sumMax += (minsqr > maxsqr) ? minsqr : maxsqr;
1556  }
1557
1558  minDistance = sumMin;
1559  maxDistance = sumMax;
1560}
1561
1562
[177]1563int
1564AxisAlignedBox3::Side(const Plane3 &plane) const
1565{
1566  Vector3 v;
1567  int i, m=3, M=-3, s;
1568 
1569  for (i=0;i<8;i++) {
1570    GetVertex(i, v);
1571    if((s = plane.Side(v)) < m)
1572      m=s;
1573    if(s > M)
1574      M=s;
1575    if (m && m==-M)
1576      return 0;
1577  }
1578 
1579  return (m == M) ? m : m + M;
1580}
[191]1581
[697]1582
1583Plane3 AxisAlignedBox3::GetPlane(const int face) const
1584{
1585        switch (face)
1586        {
1587       
1588        case 0:
1589                return Plane3(Vector3(-1, 0, 0), mMin);
1590        case 1:
1591                return Plane3(Vector3(1, 0, 0), mMax);
1592        case 2:
1593                return Plane3(Vector3(0, -1, 0), mMin);
1594        case 3:
1595                return Plane3(Vector3(0, 1, 0), mMax);
1596        case 4:
1597                return Plane3(Vector3(0, 0, -1), mMin);
1598        case 5:
1599                return Plane3(Vector3(0, 0, 1), mMax);
1600        }
1601
1602        // should not come here
1603        return Plane3();
1604}
1605
1606
[191]1607int
1608AxisAlignedBox3::GetFaceVisibilityMask(const Rectangle3 &rectangle) const
1609{
1610  int mask = 0;
1611  for (int i=0; i < 4; i++)
1612    mask |= GetFaceVisibilityMask(rectangle.mVertices[i]);
1613  return mask;
1614}
1615
[697]1616
[191]1617int
1618AxisAlignedBox3::GetFaceVisibilityMask(const Vector3 &position) const {
1619 
1620  // assume that we are not inside the box
1621  int c=0;
1622 
1623  if (position.x<(mMin.x-Limits::Small))
1624    c|=1;
1625  else
1626    if (position.x>(mMax.x+Limits::Small))
1627      c|=2;
1628 
1629  if (position.y<(mMin.y-Limits::Small))
1630    c|=4;
1631  else
1632    if (position.y>(mMax.y+Limits::Small))
1633      c|=8;
1634 
1635  if (position.z<(mMin.z-Limits::Small))
1636    c|=16;
1637  else
1638    if (position.z>(mMax.z+Limits::Small))
1639      c|=32;
1640 
1641  return c;
1642}
1643
1644
1645Rectangle3
1646AxisAlignedBox3::GetFace(const int face) const
1647{
1648  Vector3 v[4];
1649  switch (face) {
1650   
1651  case 0:
1652    v[3].SetValue(mMin.x,mMin.y,mMin.z);
1653    v[2].SetValue(mMin.x,mMax.y,mMin.z);
1654    v[1].SetValue(mMin.x,mMax.y,mMax.z);
1655    v[0].SetValue(mMin.x,mMin.y,mMax.z);
1656    break;
1657   
1658  case 1:
1659    v[0].SetValue(mMax.x,mMin.y,mMin.z);
1660    v[1].SetValue(mMax.x,mMax.y,mMin.z);
1661    v[2].SetValue(mMax.x,mMax.y,mMax.z);
1662    v[3].SetValue(mMax.x,mMin.y,mMax.z);
1663    break;
1664   
1665  case 2:
1666    v[3].SetValue(mMin.x,mMin.y,mMin.z);
1667    v[2].SetValue(mMin.x,mMin.y,mMax.z);
1668    v[1].SetValue(mMax.x,mMin.y,mMax.z);
1669    v[0].SetValue(mMax.x,mMin.y,mMin.z);
1670    break;
1671 
1672  case 3:
1673    v[0].SetValue(mMin.x,mMax.y,mMin.z);
1674    v[1].SetValue(mMin.x,mMax.y,mMax.z);
1675    v[2].SetValue(mMax.x,mMax.y,mMax.z);
1676    v[3].SetValue(mMax.x,mMax.y,mMin.z);
1677    break;
1678
1679  case 4:
1680    v[3].SetValue(mMin.x,mMin.y,mMin.z);
1681    v[2].SetValue(mMax.x,mMin.y,mMin.z);
1682    v[1].SetValue(mMax.x,mMax.y,mMin.z);
1683    v[0].SetValue(mMin.x,mMax.y,mMin.z);
1684    break;
1685 
1686  case 5:
1687    v[0].SetValue(mMin.x,mMin.y,mMax.z);
1688    v[1].SetValue(mMax.x,mMin.y,mMax.z);
1689    v[2].SetValue(mMax.x,mMax.y,mMax.z);
1690    v[3].SetValue(mMin.x,mMax.y,mMax.z);
1691    break;
1692  }
1693 
1694  return Rectangle3(v[0], v[1], v[2], v[3]);
1695}
[313]1696
[318]1697struct VertexData
1698{
1699        Vector3 mVertex;
1700        float mAngle;
1701       
1702        VertexData(Vector3 vtx, float angle): mVertex(vtx), mAngle(angle)
1703        {}
1704
1705        bool operator<(const VertexData &b) const
[372]1706        {
1707                return mAngle > b.mAngle;
[318]1708        }
1709};
1710
[312]1711// TODO: use a table to avoid normal and distance computations
[448]1712Polygon3 *AxisAlignedBox3::CrossSection(const Plane3 &plane) const
[372]1713{
1714        Polygon3 *planePoly = new Polygon3();
1715       
1716        int side[8];
1717        bool onFrontSide = false, onBackSide = false;
1718       
1719        Vector3 vtx;
1720        //-- compute classification of vertices
1721        for (int i = 0; i < 8; ++i)
1722        {
1723                GetVertex(i, vtx);
1724                side[i] = plane.Side(vtx);
1725                if (side[i] > 0)
1726                        onFrontSide = true;
1727                else if (side[i] < 0)
1728                        onBackSide = true;
1729                else // vertex coincident => push_back
1730                        planePoly->mVertices.push_back(vtx);
1731        }
1732
1733        //-- find intersections
1734        if (onFrontSide && onBackSide)
1735        {
1736                Vector3 ptA, ptB;
1737                for (int i = 0; i < 12; ++ i)
1738                {
1739                        int aIdx, bIdx;
1740                        GetEdge(i, aIdx, bIdx);
1741
1742                        ptA = GetVertex(aIdx);
1743                        ptB = GetVertex(bIdx);
1744
1745                        int sideA = side[aIdx];
1746                        int sideB = side[bIdx];
1747
1748                        if (((sideA > 0) && (sideB < 0)) || (sideA < 0) && (sideB > 0))
1749                                planePoly->mVertices.push_back(plane.FindIntersection(ptA, ptB));       
1750                }
1751        }
1752
[540]1753        // order intersections 
[372]1754        if (planePoly->mVertices.size() > 3)
1755        {
1756                Vector3 centerOfMass(0);
1757                int i;
1758                // compute center of mass
1759                for (i = 0; i < (int)planePoly->mVertices.size(); ++ i)
1760                        centerOfMass += planePoly->mVertices[i];
1761               
1762                centerOfMass /= (float)planePoly->mVertices.size();
1763
1764                vector<VertexData> vertexData;
1765               
1766                Vector3 refVec = Normalize(centerOfMass - planePoly->mVertices[0]);
1767
1768                // compute angle to reference point
1769                for (i = 1; i < (int)planePoly->mVertices.size(); ++ i)
1770                  {
1771                    float angle =
1772                      Angle(refVec, centerOfMass - planePoly->mVertices[i], plane.mNormal);
1773                   
1774                    vertexData.push_back(VertexData(planePoly->mVertices[i], angle));
1775                  }
1776               
1777                std::stable_sort(vertexData.begin(), vertexData.end());
1778
1779                // update vertices
1780                for (i = 1; i < (int)planePoly->mVertices.size(); ++ i)
1781                        planePoly->mVertices[i] = vertexData[i - 1].mVertex;
1782        }
1783        else if (planePoly->mVertices.size() == 3)
1784        {
1785                // fix orientation if needed
[318]1786                if (DotProd(planePoly->GetNormal(), plane.mNormal) < 0)
1787                {
1788                        Vector3 v = planePoly->mVertices[1];
1789                        planePoly->mVertices[1] = planePoly->mVertices[2];
1790                        planePoly->mVertices[2] = v;
[372]1791                }
1792        }
1793       
1794        return planePoly;
1795}
[448]1796
[557]1797
[463]1798bool AxisAlignedBox3::GetRaySegment(const Ray &ray,
1799                                                                        float &minT,
1800                                                                        float &maxT) const
1801{
[704]1802        maxT = 1e15f;
[463]1803        minT = 0;
1804       
[557]1805        // test with bounding box
[463]1806        if (!GetMinMaxT(ray, &minT, &maxT))
1807                return false;
1808
[557]1809        if (minT < 0)
[463]1810                minT = 0;
1811
1812        // bound ray or line segment
1813        if (//(ray.GetType() == Ray::LOCAL_RAY) &&
1814            !ray.intersections.empty() &&
1815                (ray.intersections[0].mT <= maxT))
1816        {
1817                maxT = ray.intersections[0].mT;
1818        }
1819
1820        return true;
1821}
[482]1822
1823
[492]1824  // Compute tmin and tmax for a ray, whenever required .. need not pierce box
1825int
1826AxisAlignedBox3::ComputeMinMaxT(const Vector3 &origin,
1827                                                                const Vector3 &direction,
1828                                                                float *tmin,
1829                                                                float *tmax) const
1830{
1831
1832  register float minx, maxx;
1833
1834 
1835  Vector3 invDirection;
[503]1836  const float eps = 1e-6f;
1837  const float invEps = 1e6f;
[492]1838 
1839  // it does change the ray direction very slightly,
1840  // but the size direction vector is not practically changed
1841 
1842  if (fabs(direction.x) < eps) {
[503]1843    if (direction.x < 0.0f)
[492]1844      invDirection.x = -invEps;
1845    else
1846      invDirection.x = invEps;
1847  }
1848  else
[503]1849    invDirection.x = 1.0f / direction.x;
[492]1850 
1851  if (fabs(direction.y) < eps) {
1852    if (direction.y < 0.0)
1853      invDirection.y = -invEps;
1854    else
1855      invDirection.y = invEps;
1856  }
1857  else
[503]1858    invDirection.y = 1.0f / direction.y;
[492]1859 
1860  if (fabs(direction.z) < eps) {
[503]1861    if (direction.z < 0.0f)
[492]1862      invDirection.z = -invEps;
1863    else
1864      invDirection.z = invEps;
1865  }
1866  else
[503]1867    invDirection.z = 1.0f / direction.z;
[492]1868
1869
1870 
[503]1871  if (fabs(direction.x) < 0.001f) {
[492]1872    if (mMin.x < origin.x && mMax.x > origin.x) {
1873      minx = -MAXFLOAT;
1874      maxx = MAXFLOAT;
1875    }
1876    else
1877      return 0;
1878  }
1879  else {
1880    float t1 = (mMin.x - origin.x) * invDirection.x;
1881    float t2 = (mMax.x - origin.x) * invDirection.x;
1882    if (t1 < t2) {
1883      minx = t1;
1884      maxx = t2;
1885    }
1886    else {
1887      minx = t2;
1888      maxx = t1;
1889    }
1890        //    if (maxx < 0.0)
1891        //      return 0;
1892  }
1893
1894  *tmin = minx;
1895  *tmax = maxx;
1896 
1897  if (fabs(direction.y) < 0.001) {
1898    if (mMin.y < origin.y && mMax.y > origin.y) {
1899      minx = -MAXFLOAT;
1900      maxx = MAXFLOAT;
1901    }
1902    else
1903      return 0;
1904  }
1905  else {
1906    float t1 = (mMin.y - origin.y) * invDirection.y;
1907    float t2 = (mMax.y - origin.y) * invDirection.y;
1908    if (t1 < t2) {
1909      minx = t1;
1910      maxx = t2;
1911    }
1912    else {
1913      minx = t2;
1914      maxx = t1;
1915    }
1916        //    if (maxx < 0.0)
1917        //      return 0;
1918  }
1919
1920  if (minx > *tmin)
1921    *tmin = minx;
1922  if (maxx < *tmax)
1923    *tmax = maxx;
1924 
1925  if (fabs(direction.z) < 0.001) {
1926    if (mMin.z < origin.z && mMax.z > origin.z) {
1927      minx = -MAXFLOAT;
1928      maxx = MAXFLOAT;
1929    }
1930    else
1931      return 0;
1932  }
1933  else {
1934    float t1 = (mMin.z - origin.z) * invDirection.z;
1935    float t2 = (mMax.z - origin.z) * invDirection.z;
1936    if (t1 < t2) {
1937      minx = t1;
1938      maxx = t2;
1939    }
1940    else {
1941      minx = t2;
1942      maxx = t1;
1943    }
1944        //    if (maxx < 0.0)
1945        //      return 0;
1946  }
1947
1948  if (minx > *tmin)
1949    *tmin = minx;
1950  if (maxx < *tmax)
1951    *tmax = maxx;
1952
1953  return 1; // yes, intersection was found
1954}
1955
1956
[482]1957bool AxisAlignedBox3::GetIntersectionFace(Rectangle3 &face,
1958                                                                                  const AxisAlignedBox3 &neighbour) const
1959                                                                                   
1960{
1961        if (EpsilonEqual(mMin[0], neighbour.Max(0)))
1962        {
1963                float maxy = min(mMax.y, neighbour.mMax.y);
1964                float maxz = min(mMax.z, neighbour.mMax.z);
1965                float miny = max(mMin.y, neighbour.mMin.y);
1966                float minz = max(mMin.z, neighbour.mMin.z);
1967
1968                face.mVertices[3].SetValue(mMin.x, miny, minz);
1969                face.mVertices[2].SetValue(mMin.x, maxy, minz);
1970                face.mVertices[1].SetValue(mMin.x, maxy, maxz);
1971                face.mVertices[0].SetValue(mMin.x, miny, maxz);
1972
1973        return true;
1974        }
1975    if (EpsilonEqual(mMax[0], neighbour.Min(0)))
1976        {
1977                float maxy = min(mMax.y, neighbour.mMax.y);
1978                float maxz = min(mMax.z, neighbour.mMax.z);
1979                float miny = max(mMin.y, neighbour.mMin.y);
1980                float minz = max(mMin.z, neighbour.mMin.z);
1981
1982                face.mVertices[0].SetValue(mMax.x, miny, minz);
1983                face.mVertices[1].SetValue(mMax.x, maxy, minz);
1984                face.mVertices[2].SetValue(mMax.x, maxy, maxz);
1985                face.mVertices[3].SetValue(mMax.x, miny, maxz);
1986
1987                return true;
1988        }
1989    if (EpsilonEqual(mMin[1], neighbour.Max(1)))
1990        {
1991                float maxx = min(mMax.x, neighbour.mMax.x);
1992                float maxz = min(mMax.z, neighbour.mMax.z);
1993                float minx = max(mMin.x, neighbour.mMin.x);
1994                float minz = max(mMin.z, neighbour.mMin.z);
1995
1996                face.mVertices[3].SetValue(minx, mMin.y, minz);
1997                face.mVertices[2].SetValue(minx, mMin.y, maxz);
1998                face.mVertices[1].SetValue(maxx, mMin.y, maxz);
1999                face.mVertices[0].SetValue(maxx, mMin.y, minz);
2000               
2001                return true;
2002        }
2003        if (EpsilonEqual(mMax[1], neighbour.Min(1)))
2004        {               
2005                float maxx = min(mMax.x, neighbour.mMax.x);
2006                float maxz = min(mMax.z, neighbour.mMax.z);
2007                float minx = max(mMin.x, neighbour.mMin.x);
2008                float minz = max(mMin.z, neighbour.mMin.z);
2009
2010                face.mVertices[0].SetValue(minx, mMax.y, minz);
2011                face.mVertices[1].SetValue(minx, mMax.y, maxz);
2012                face.mVertices[2].SetValue(maxx, mMax.y, maxz);
2013                face.mVertices[3].SetValue(maxx, mMax.y, minz);
2014
2015                return true;
2016        }
2017        if (EpsilonEqual(mMin[2], neighbour.Max(2)))
2018        {
2019                float maxx = min(mMax.x, neighbour.mMax.x);
2020                float maxy = min(mMax.y, neighbour.mMax.y);
2021                float minx = max(mMin.x, neighbour.mMin.x);
2022                float miny = max(mMin.y, neighbour.mMin.y);
2023
2024                face.mVertices[3].SetValue(minx, miny, mMin.z);
2025                face.mVertices[2].SetValue(maxx, miny, mMin.z);
2026                face.mVertices[1].SetValue(maxx, maxy, mMin.z);
2027                face.mVertices[0].SetValue(minx, maxy, mMin.z);
2028       
2029                return true;
2030        }
2031        if (EpsilonEqual(mMax[2], neighbour.Min(2)))
2032        {
2033                float maxx = min(mMax.x, neighbour.mMax.x);
2034                float maxy = min(mMax.y, neighbour.mMax.y);
2035                float minx = max(mMin.x, neighbour.mMin.x);
2036                float miny = max(mMin.y, neighbour.mMin.y);
2037
2038                face.mVertices[0].SetValue(minx, miny, mMax.z);
2039                face.mVertices[1].SetValue(maxx, miny, mMax.z);
2040                face.mVertices[2].SetValue(maxx, maxy, mMax.z);
2041                face.mVertices[3].SetValue(minx, maxy, mMax.z);
2042
2043                return true;
2044        }
2045
2046        return false;
2047}
2048
[486]2049
[860]2050void IncludeBoxInMesh(const AxisAlignedBox3 &box, Mesh &mesh)
[486]2051{
2052        // add 6 vertices of the box
[860]2053        int index = (int)mesh.mVertices.size();
[486]2054       
[860]2055        for (int i=0; i < 8; ++ i)
[486]2056        {
2057                Vector3 v;
[860]2058                box.GetVertex(i, v);
2059                mesh.mVertices.push_back(v);
[486]2060        }
2061       
[860]2062        mesh.AddFace(new Face(index + 0, index + 1, index + 3, index + 2) );
2063        mesh.AddFace(new Face(index + 0, index + 2, index + 6, index + 4) );
2064        mesh.AddFace(new Face(index + 4, index + 6, index + 7, index + 5) );
[486]2065       
[860]2066        mesh.AddFace(new Face(index + 3, index + 1, index + 5, index + 7) );
2067        mesh.AddFace(new Face(index + 0, index + 4, index + 5, index + 1) );
2068        mesh.AddFace(new Face(index + 2, index + 3, index + 7, index + 6) );
[542]2069}
[486]2070
[542]2071
2072void AxisAlignedBox3::ExtractPolys(PolygonContainer &polys) const
2073{
2074        Polygon3 *face1 = new Polygon3();
2075        polys.push_back(face1);   
2076
[545]2077    face1->mVertices.push_back(Vector3(mMin.x, mMin.y, mMax.z));
2078        face1->mVertices.push_back(Vector3(mMin.x, mMax.y, mMax.z));
2079        face1->mVertices.push_back(Vector3(mMin.x, mMax.y ,mMin.z));
2080        face1->mVertices.push_back(Vector3(mMin.x, mMin.y, mMin.z));
[542]2081
2082        Polygon3 *face2 = new Polygon3(); 
2083        polys.push_back(face2);
[545]2084   
2085    face2->mVertices.push_back(Vector3(mMax.x, mMin.y, mMin.z));
2086    face2->mVertices.push_back(Vector3(mMax.x, mMax.y, mMin.z));
2087    face2->mVertices.push_back(Vector3(mMax.x, mMax.y, mMax.z));
2088    face2->mVertices.push_back(Vector3(mMax.x, mMin.y, mMax.z));
[542]2089 
2090        Polygon3 *face3 = new Polygon3(); 
2091        polys.push_back(face3);
2092
[545]2093    face3->mVertices.push_back(Vector3(mMax.x, mMin.y ,mMin.z));
2094        face3->mVertices.push_back(Vector3(mMax.x, mMin.y, mMax.z));
2095        face3->mVertices.push_back(Vector3(mMin.x, mMin.y, mMax.z));
2096        face3->mVertices.push_back(Vector3(mMin.x, mMin.y, mMin.z));
[542]2097
2098        Polygon3 *face4 = new Polygon3(); 
2099        polys.push_back(face4);
2100
[545]2101        face4->mVertices.push_back(Vector3(mMin.x, mMax.y, mMin.z));
2102        face4->mVertices.push_back(Vector3(mMin.x, mMax.y, mMax.z));
2103        face4->mVertices.push_back(Vector3(mMax.x, mMax.y, mMax.z));
2104        face4->mVertices.push_back(Vector3(mMax.x, mMax.y, mMin.z));
[542]2105   
2106        Polygon3 *face5 = new Polygon3(); 
2107        polys.push_back(face5);
2108
[545]2109        face5->mVertices.push_back(Vector3(mMin.x, mMax.y, mMin.z));
2110    face5->mVertices.push_back(Vector3(mMax.x, mMax.y, mMin.z));
2111    face5->mVertices.push_back(Vector3(mMax.x, mMin.y, mMin.z));
2112        face5->mVertices.push_back(Vector3(mMin.x, mMin.y, mMin.z));
[542]2113
2114        Polygon3 *face6 = new Polygon3(); 
2115        polys.push_back(face6); 
2116 
[545]2117    face6->mVertices.push_back(Vector3(mMin.x, mMin.y, mMax.z));
2118    face6->mVertices.push_back(Vector3(mMax.x, mMin.y, mMax.z));
2119    face6->mVertices.push_back(Vector3(mMax.x, mMax.y, mMax.z));
2120    face6->mVertices.push_back(Vector3(mMin.x, mMax.y, mMax.z));
[547]2121
[860]2122}
2123
[1047]2124
[1006]2125void AxisAlignedBox3::Split(const int axis,
2126                                                        const float value,
[1047]2127                                                        AxisAlignedBox3 &front,
2128                                                        AxisAlignedBox3 &back) const
[1006]2129{
[1415]2130        if ((value >= mMin[axis]) && (value <= mMax[axis]))
[1006]2131        {
[1047]2132                front.mMin = mMin; front.mMax = mMax;
2133                back.mMin = mMin; back.mMax = mMax;
[1006]2134
[1047]2135                back.mMax[axis] = value;
2136                front.mMin[axis] = value;
[1006]2137        }
2138}
[1047]2139
[1415]2140 
2141void AxisAlignedBox3::Scale(const float scale)
2142{
2143        Vector3 newSize = Size()*(scale*0.5f);
2144        Vector3 center = Center();
2145        mMin = center - newSize;
2146        mMax = center + newSize;
2147}
[1047]2148
[1415]2149
2150void AxisAlignedBox3::Scale(const Vector3 &scale)
2151{
2152        Vector3 newSize = Size()*(scale*0.5f);
2153        Vector3 center = Center();
2154        mMin = center - newSize;
2155        mMax = center + newSize;
2156}
2157
2158
2159void AxisAlignedBox3::Translate(const Vector3 &shift)
2160{
2161         mMin += shift;
2162         mMax += shift;
2163}
2164
2165
[1486]2166void AxisAlignedBox3::Initialize() {
2167    mMin = Vector3(MAXFLOAT);
2168    mMax = Vector3(-MAXFLOAT);
2169  }
2170
2171
2172Vector3 AxisAlignedBox3::Center() const
2173{
2174        return 0.5 * (mMin + mMax);
2175}
2176 
2177Vector3 AxisAlignedBox3::Diagonal() const
2178{
2179        return (mMax -mMin);
2180}
2181
2182
2183float AxisAlignedBox3::Center(const int axis) const
2184{
2185        return  0.5f * (mMin[axis] + mMax[axis]);
2186}
2187
2188
2189float AxisAlignedBox3::Min(const int axis) const
2190{
2191        return mMin[axis];
2192}
2193
2194float AxisAlignedBox3::Max(const int axis) const {
2195        return  mMax[axis];
2196}
2197
2198float AxisAlignedBox3::Size(const int axis) const {
2199        return  Max(axis) - Min(axis);
2200}
2201
2202// Read-only const access tomMin and max vectors using references
2203const Vector3& AxisAlignedBox3::Min() const
2204{
2205        return mMin;
2206}
2207
2208
2209const Vector3& AxisAlignedBox3::Max() const
2210{
2211        return mMax;
2212}
2213
2214
2215void AxisAlignedBox3::Enlarge (const Vector3 &v)
2216{
2217        mMax += v;
2218        mMin -= v;
2219}
2220
2221
2222void AxisAlignedBox3::SetMin(const Vector3 &v)
2223{
2224        mMin = v;
2225}
2226
2227
2228void AxisAlignedBox3::SetMax(const Vector3 &v)
2229{
2230        mMax = v;
2231}
2232
2233
2234void AxisAlignedBox3::SetMin(int axis, const float value)
2235{
2236        mMin[axis] = value;
2237}
2238
2239
2240void AxisAlignedBox3::SetMax(int axis, const float value)
2241{
2242        mMax[axis] = value;
2243}
2244
2245// Decrease box by given splitting plane
2246void AxisAlignedBox3::Reduce(int axis, int right, float value)
2247{
2248        if ( (value >=mMin[axis]) && (value <= mMax[axis]) )
2249                if (right)
2250                        mMin[axis] = value;
2251                else
2252                        mMax[axis] = value;
2253}
2254
2255
2256Vector3 AxisAlignedBox3::Size() const
2257{
2258        return mMax - mMin;
2259}
2260
[1581]2261}
Note: See TracBrowser for help on using the repository browser.