aabb.cc 11.5 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
/*
 * File:   aabb.cc
 * Author: Jordan P. Lefebvre, lefebvre.jordan@gmail.com
 *
 * Created on August 13, 2012, 8:23 PM
 */

#include "radixgeometry/aabb.hh"
#include "radixgeometry/interval.hh"

#include <sstream>
#include <iomanip>
#include <algorithm>
namespace radix
{
AABB::AABB(void)
17
18
    : x0(-1), x1(1), y0(-1), y1(1), z0(-1), z1(1)
{
19
20
21
22
23
}

AABB::AABB(const Real _x0, const Real _x1,
           const Real _y0, const Real _y1,
           const Real _z0, const Real _z1)
24
25
    : x0(_x0), x1(_x1), y0(_y0), y1(_y1), z0(_z0), z1(_z1)
{
26
27
28
}

AABB::AABB(const Point3D p0, const Point3D p1)
29
30
    : x0(p0.x), x1(p1.x), y0(p0.y), y1(p1.y), z0(p0.z), z1(p1.z)
{
31
32
33
}

AABB::AABB(const AABB& aabb)
34
35
    : x0(aabb.x0), x1(aabb.x1), y0(aabb.y0), y1(aabb.y1), z0(aabb.z0), z1(aabb.z1)
{
36
37
}

38
39
AABB& AABB::operator=(const AABB& rhs)
{
40
41
42
43
44
45
46
47
48
49
50
51
52
    if (this == &rhs)
        return (*this);

    x0 = rhs.x0;
    x1 = rhs.x1;
    y0 = rhs.y0;
    y1 = rhs.y1;
    z0 = rhs.z0;
    z1 = rhs.z1;

    return (*this);
}//operation=

53
54
AABB::~AABB(void)
{
55
56
}

57
58
bool AABB::hit(const Ray& ray, Real& t) const
{
59
60
61
62

    Real tx_min, ty_min, tz_min;
    Real tx_max, ty_max, tz_max;
    Real t_near, t_far;
63
64
    if( ray.d.x == 1 || ray.d.x == -1)
    {
65
66
        if( ray.o.y < y0 || ray.o.y > y1 ) return false;
        if( ray.o.z < z0 || ray.o.z > z1 ) return false;
67
68
        if( ray.d.x > 0 )
        {
69
70
            tx_min = (x0-ray.o.x) * ray.d.x;
            tx_max = (x1-ray.o.x) * ray.d.x;
71
72
73
        }
        else
        {
74
75
76
77
            tx_min = (x1-ray.o.x) * ray.d.x;
            tx_max = (x0-ray.o.x) * ray.d.x;
        }
        if( tx_max < kEpsilon ) return false;
78
79
        if( tx_min <= kEpsilon )
        {
80
            t = tx_max;
81
82
83
        }
        else
        {
84
85
86
87
            t = tx_min;
        }
        return true;
    }
88
89
    if(  ray.d.y == 1 || ray.d.y == -1)
    {
90
91
        if( ray.o.x < x0 || ray.o.x > x1 ) return false;
        if( ray.o.z < z0 || ray.o.z > z1 ) return false;
92
93
        if( ray.d.y > 0 )
        {
94
95
            ty_min = (y0-ray.o.y) * ray.d.y;
            ty_max = (y1-ray.o.y) * ray.d.y;
96
97
98
        }
        else
        {
99
100
101
102
            ty_min = (y1-ray.o.y) * ray.d.y;
            ty_max = (y0-ray.o.y) * ray.d.y;
        }
        if( ty_max < kEpsilon ) return false;
103
104
        if( ty_min <= kEpsilon )
        {
105
            t = ty_max;
106
107
108
        }
        else
        {
109
110
111
112
            t = ty_min;
        }
        return true;
    }
113
114
    if(  ray.d.z == 1 || ray.d.z == -1 )
    {
115
116
        if( ray.o.y < y0 || ray.o.y > y1 ) return false;
        if( ray.o.x < x0 || ray.o.x > x1 ) return false;
117
118
        if( ray.d.z > 0 )
        {
119
120
            tz_min = (z0-ray.o.z) * ray.d.z;
            tz_max = (z1-ray.o.z) * ray.d.z;
121
122
123
        }
        else
        {
124
125
126
127
            tz_min = (z1-ray.o.z) * ray.d.z;
            tz_max = (z0-ray.o.z) * ray.d.z;
        }
        if( tz_max < kEpsilon ) return false;
128
129
        if( tz_min <= kEpsilon )
        {
130
            t = tz_max;
131
132
133
        }
        else
        {
134
135
136
137
138
139
            t = tz_min;
        }

        return true;
    }
    Real a = 1.0 / ray.d.x;
140
141
    if (a >= 0)
    {
142
143
144
145
        Real a1 = x0 - ray.o.x;
        Real a2 = x1 - ray.o.x;
        tx_min = (a1) * a;
        tx_max = (a2) * a;
146
147
148
    }
    else
    {
149
150
151
152
153
154
155
        Real a1 = x0 - ray.o.x;
        Real a2 = x1 - ray.o.x;
        tx_min = (a2) * a;
        tx_max = (a1) * a;
    }

    Real b = 1.0 / ray.d.y;
156
157
    if (b >= 0)
    {
158
159
160
161
        Real b1 = y0 - ray.o.y;
        Real b2 = y1 - ray.o.y;
        ty_min = (b1) * b;
        ty_max = (b2) * b;
162
163
164
    }
    else
    {
165
166
167
168
169
170
171
        Real b1 = y0 - ray.o.y;
        Real b2 = y1 - ray.o.y;
        ty_min = (b2) * b;
        ty_max = (b1) * b;
    }

    Real c = 1.0 / ray.d.z;
172
173
    if (c >= 0)
    {
174
175
176
177
        Real c1 = z0 - ray.o.z;
        Real c2 = z1 - ray.o.z;
        tz_min = (c1) * c;
        tz_max = (c2) * c;
178
179
180
    }
    else
    {
181
182
183
184
185
186
187
188
189
        Real c1 = z0 - ray.o.z;
        Real c2 = z1 - ray.o.z;
        tz_min = (c2) * c;
        tz_max = (c1) * c;
    }


    //int face_in, face_out;
    // find the largest entering t value
190
191
    if (tx_min > ty_min)
    {
192
193
        t_near = tx_min;
        //face_in = (a >= 0.0) ? 0 : 3;
194
195
196
    }
    else
    {
197
198
199
200
        t_near = ty_min;
        //face_in = (b >= 0.0) ? 1 : 4;
    }

201
202
    if (tz_min > t_near)
    {
203
204
205
206
207
        t_near = tz_min;
        //face_in = (c >= 0.0) ? 2 : 5;
    }

    //find smallest exiting t value
208
209
    if (tx_max < ty_max)
    {
210
211
        t_far = tx_max;
        //face_out = (a >= 0.0) ? 3 : 0;
212
213
214
    }
    else
    {
215
216
217
        t_far = ty_max;
        //face_out = (b >= 0.0) ? 4 : 1;
    }
218
219
    if (tz_max < t_far)
    {
220
221
222
        t_far = tz_max;
        //face_out = (c >= 0.0) ? 5 : 2;
    }
223
224
225
226
    if (t_near < t_far && t_far > kEpsilon)
    {
        if (t_near > kEpsilon)
        {
227
            t = t_near;
228
229
230
        }
        else
        {
231
232
233
234
235
236
237
            t = t_far;
        }
        return true;
    }
    return false;
}//hit

238
239
bool AABB::inside(const Point3D& p) const
{
240
241
242
243
    return ((p.x > x0 && p.x < x1)
            && (p.y > y0 && p.y < y1)
            && (p.z > z0 && p.z < z1));
}//inside
244
245
bool AABB::on(const Point3D& p) const
{
246
247
248
249
250
251
252
253
254
255
    return ( ((p.x == x0 || p.x == x1)
              && (p.y >= y0 && p.y <= y1)
              && (p.z >= z0 && p.z <= z1))
             || ((p.y == y0 || p.y == y1)
                 && (p.x >= x0 && p.x <= x1)
                 && (p.z >= z0 && p.z <= z1))
             || ((p.z == z0 || p.z == z1)
                 && (p.x >= x0 && p.x <= x1)
                 && (p.y >= y0 && p.y <= y1)));
}
256
257
std::string AABB::toString()const
{
258
259
260
261
262
263
264
    std::stringstream stream;
    stream << std::setprecision(15)
           << "aabb x0="<<x0<<",x1="<<x1
           << ",y0="<<y0<<",y1="<<y1
           << ",z0="<<z0<<",z1="<<z1;
    return stream.str();
}
265
266
AABB AABB::operator+(const AABB& rhs)
{
267
268
269
270
271
272
    // union with a null AABB results in this
    if( isNull(rhs) ) return *this;
    return AABB(std::min(rhs.x0,this->x0), std::max(rhs.x1,this->x1)
                ,std::min(rhs.y0,this->y0), std::max(rhs.y1,this->y1)
                ,std::min(rhs.z0,this->z0), std::max(rhs.z1,this->z1));
}
273
274
AABB AABB::operator&(const AABB& rhs)
{
275
276
277
278
279
280
281
282
283
284
285
    // intersection with a null AABB results in null
    if( isNull(rhs) ) return AABB(0,0,0,0,0,0);
    Interval x,y,z;
    bool xb = Interval(x0,x1,0).intersect(Interval(rhs.x0,rhs.x1,1),x);
    bool yb = Interval(y0,y1,0).intersect(Interval(rhs.y0,rhs.y1,1),y);
    bool zb = Interval(z0,z1,0).intersect(Interval(rhs.z0,rhs.z1,1),z);
    if( !xb || !yb || !zb ) return AABB(0,0,0,0,0,0); // no bounding box, it is a point
    return AABB(x.a,x.b
                ,y.a,y.b
                ,z.a,z.b);
}
286
287
bool AABB::isNull(const AABB &aabb)
{
288
289
    return !aabb.hasVolume();
}
290
291
Point3D AABB::origin() const
{
292
293
294
    return Point3D((x0+x1)/2.0
                   ,(y0+y1)/2.0
                   ,(z0+z1)/2.0
295
                  );
296
297
298
299
}
bool AABB::intersects(const AABB& o)const
{
    bool xIntersects = (x0 >= o.x0 && x0 < o.x1)
300
301
302
                       || (x1 > o.x0 && x1 <= o.x1)
                       || (o.x0 >= x0 && o.x0 < x1)
                       || (o.x1 > x0 && o.x1 <= x1);
303
304

    bool yIntersects = (y0 >= o.y0 && y0 < o.y1)
305
306
307
                       || (y1 > o.y0 && y1 <= o.y1)
                       || (o.y0 >= y0 && o.y0 < y1)
                       || (o.y1 > y0 && o.y1 <= y1);
308
    bool zIntersects = (z0 >= o.z0 && z0 < o.z1)
309
310
311
                       || (z1 > o.z0 && z1 <= o.z1)
                       || (o.z0 >= z0 && o.z0 < z1)
                       || (o.z1 > z0 && o.z1 <= z1);
312
313
314
    return (xIntersects && yIntersects && zIntersects);
}

315
316
AABB AABB::operator*(const Matrix & matrix)const
{
317
318
319
320
321
322
323
324
325
326
327
328
329
    Point3D x0y0z0= matrix * Point3D(x0,y0,z0);
    Point3D x1y0z0= matrix * Point3D(x1,y0,z0);
    Point3D x0y1z0= matrix * Point3D(x0,y1,z0);
    Point3D x0y1z1= matrix * Point3D(x0,y1,z1);
    Point3D x1y0z1= matrix * Point3D(x1,y0,z1);
    Point3D x0y0z1= matrix * Point3D(x0,y0,z1);
    Point3D x1y1z1= matrix * Point3D(x1,y1,z1);
    Point3D x1y1z0= matrix * Point3D(x1,y1,z0);

    AABB result;
    result.x0 = std::min(x0y0z0.x
                         ,std::min(x1y0z0.x
                                   ,std::min(x0y1z0.x
330
331
332
333
                                           ,std::min(x0y1z1.x
                                                   ,std::min(x1y0z1.x
                                                           ,std::min(x0y0z1.x
                                                                   ,std::min(x1y1z1.x,x1y1z0.x)))))));
334
335
336
337
338
    if (std::abs(result.x0) == infinity) result.x0 = (result.x0 < 0) ? -maxReal : maxReal;

    result.x1 = std::max(x0y0z0.x
                         ,std::max(x1y0z0.x
                                   ,std::max(x0y1z0.x
339
340
341
342
                                           ,std::max(x0y1z1.x
                                                   ,std::max(x1y0z1.x
                                                           ,std::max(x0y0z1.x
                                                                   ,std::max(x1y1z1.x,x1y1z0.x)))))));
343
344
345
346
347
    if (std::abs(result.x1) == infinity) result.x1 = (result.x1 < 0) ? -maxReal : maxReal;

    result.y0 = std::min(x0y0z0.y
                         ,std::min(x1y0z0.y
                                   ,std::min(x0y1z0.y
348
349
350
351
                                           ,std::min(x0y1z1.y
                                                   ,std::min(x1y0z1.y
                                                           ,std::min(x0y0z1.y
                                                                   ,std::min(x1y1z1.y,x1y1z0.y)))))));
352
353
354
355
356
    if (std::abs(result.y0) == infinity) result.y0 = (result.y0 < 0) ? -maxReal : maxReal;

    result.y1 = std::max(x0y0z0.y
                         ,std::max(x1y0z0.y
                                   ,std::max(x0y1z0.y
357
358
359
360
                                           ,std::max(x0y1z1.y
                                                   ,std::max(x1y0z1.y
                                                           ,std::max(x0y0z1.y
                                                                   ,std::max(x1y1z1.y,x1y1z0.y)))))));
361
362
363
364
365
    if (std::abs(result.y1) == infinity) result.y1 = (result.y1 < 0) ? -maxReal : maxReal;

    result.z0 = std::min(x0y0z0.z
                         ,std::min(x1y0z0.z
                                   ,std::min(x0y1z0.z
366
367
368
369
                                           ,std::min(x0y1z1.z
                                                   ,std::min(x1y0z1.z
                                                           ,std::min(x0y0z1.z
                                                                   ,std::min(x1y1z1.z,x1y1z0.z)))))));
370
371
372
373
374
    if (std::abs(result.z0) == infinity) result.z0 = (result.z0 < 0) ? -maxReal : maxReal;

    result.z1 = std::max(x0y0z0.z
                         ,std::max(x1y0z0.z
                                   ,std::max(x0y1z0.z
375
376
377
378
                                           ,std::max(x0y1z1.z
                                                   ,std::max(x1y0z1.z
                                                           ,std::max(x0y0z1.z
                                                                   ,std::max(x1y1z1.z,x1y1z0.z)))))));
379
380
381
382
383
    if (std::abs(result.z1) == infinity) result.z1 = (result.z1 < 0) ? -maxReal : maxReal;

    return result;
}

384
385
std::vector<Point3D> AABB::vertices() const
{
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
    std::vector<Point3D> vertices;
    vertices.push_back(Point3D(x0,y0,z0));
    vertices.push_back(Point3D(x0,y0,z1));
    vertices.push_back(Point3D(x0,y1,z0));
    vertices.push_back(Point3D(x1,y0,z0));
    vertices.push_back(Point3D(x1,y1,z0));
    vertices.push_back(Point3D(x1,y0,z1));
    vertices.push_back(Point3D(x1,y1,z1));
    vertices.push_back(Point3D(x0,y1,z1));
    return vertices;
}
}//namespace radix