Commit 79774a67 authored by Lukas Werner's avatar Lukas Werner
Browse files

WIP commit of AABB intersection check

parent d305e418
Pipeline #6424 failed with stage
in 65 minutes and 14 seconds
......@@ -19,7 +19,9 @@ namespace walberla {
namespace pe {
namespace raytracing {
inline bool intersects(SphereID sphere, Ray* ray, real_t &t);
inline bool intersects(PlaneID plane, Ray* ray, real_t &t);
inline bool intersects(PlaneID plane, Ray* ray, real_t* t);
inline bool intersects(BoxID box, Ray* ray, real_t* t);
inline bool intersects(GeomPrimitive* body, Ray* ray, real_t* t);
struct IntersectsFunctor
{
......@@ -55,7 +57,7 @@ namespace walberla {
real_t t_;
t_ = (t0 < t1) ? t0 : t1; // assign the closest hit point to t
if (t_ < 0) {
// at least one of the calculated hit point is behind the rays origin
// at least one of the calculated hit points is behind the rays origin
if (t1 < 0) {
// both of the points are behind the origin (ray does not hit sphere)
*t = real_t(INFINITY);
......@@ -87,6 +89,74 @@ namespace walberla {
*t = INFINITY;
return false;
}
// does not work yet
inline bool intersects(BoxID box, Ray* ray, real_t* t) {
return intersects((GeomPrimitive*)box, ray, t);
}
// does not work yet
inline bool intersects(GeomPrimitive* body, Ray* ray, real_t* t) {
// An Efficient and Robust Ray–Box Intersection Algorithm: http://people.csail.mit.edu/amy/papers/box-jgt.pdf
AABB aabb = body->getAABB();
Vec3 bounds[2] = {
aabb.min(),
aabb.max()
};
real_t tmincoord = 0, tmaxcoord = 0; // x
real_t txmin, txmax;
real_t tmin = txmin = (bounds[ray->sign[0]][0] - ray->origin[0]) * ray->inv_direction[0];
real_t tmax = txmax = (bounds[1-ray->sign[0]][0] - ray->origin[0]) * ray->inv_direction[0];
real_t tymin = (bounds[ray->sign[1]][1] - ray->origin[1]) * ray->inv_direction[1];
real_t tymax = (bounds[1-ray->sign[1]][1] - ray->origin[1]) * ray->inv_direction[1];
if (tmin > tymax || tymin > tmax) {
*t = INFINITY;
return false;
}
if (tymin > tmin) {
tmincoord = 1; //y
tmin = tymin;
}
if (tymax < tmax) {
tmaxcoord = 1; //y
tmax = tymax;
}
real_t tzmin = (bounds[ray->sign[2]][2] - ray->origin[2]) * ray->inv_direction[2];
real_t tzmax = (bounds[1-ray->sign[2]][2] - ray->origin[2]) * ray->inv_direction[2];
if (tmin > tzmax || tzmin > tmax) {
*t = INFINITY;
return false;
}
if (tzmin > tmin) {
tmincoord = 2; //z
tmin = tzmin;
}
if (tzmax < tmax) {
tmaxcoord = 2; //z
tmax = tzmax;
}
WALBERLA_LOG_INFO("tmax: " << tmax << "; tmin: " << tmin);
WALBERLA_LOG_INFO("tmincoord: " << tmincoord << "; tmaxcoord: " << tmaxcoord);
real_t t_;
if (tmin > 0) {
// ray hit box from outside
t_ = tmin;
} else if (tmax < 0) {
// tmin and tmax are smaller than 0 -> box is in rays negative direction
*t = INFINITY;
return false;
} else {
// ray origin within box
t_ = tmax;
}
*t = t_;
return true;
}
}
}
}
......@@ -7,11 +7,16 @@ namespace pe {
namespace raytracing {
class Ray {
public:
Vec3 origin; // coordinates of screen pixel
Vec3 origin;
Vec3 direction; // e.g. (1,0,0) for ray in x direction
Vec3 inv_direction;
int sign[3];
Ray (Vec3 origin_, Vec3 direction_) : origin(origin_), direction(direction_) {
inv_direction = Vec3(1/direction[0], 1/direction[1], 1/direction[2]);
sign[0] = (inv_direction[0] < 0);
sign[1] = (inv_direction[1] < 0);
sign[2] = (inv_direction[2] < 0);
}
};
}
......
......@@ -48,7 +48,6 @@ void SphereIntersectsTest()
// sphere around ray origin
Sphere sp3(123, 1, Vec3(-5,3,3), Vec3(0,0,0), Quat(), 2, iron, false, true, false);
WALBERLA_CHECK(intersects(&sp3, &ray1, &t));
WALBERLA_LOG_INFO(t);
WALBERLA_CHECK(realIsEqual(t, real_t(2)));
}
......@@ -77,6 +76,41 @@ void PlaneIntersectsTest() {
WALBERLA_CHECK(!intersects(&pl4, &ray1, &t), "ray hit plane behind origin");
}
void BoxIntersectsTest() {
MaterialID iron = Material::find("iron");
Box box1(127, 5, Vec3(-15, 0, 0), Vec3(0, 0, 0), Quat(), Vec3(10, 10, 10), iron, false, true, false);
Ray ray1(Vec3(-5,3,3), Vec3(1,0,0));
real_t t;
WALBERLA_LOG_INFO("RAY -> BOX");
WALBERLA_CHECK(!intersects(&box1, &ray1, &t));
Box box2(128, 5, Vec3(-2, 0, 0), Vec3(0, 0, 0), Quat(), Vec3(10, 10, 10), iron, false, true, false);
WALBERLA_CHECK(intersects(&box2, &ray1, &t));
WALBERLA_CHECK(realIsEqual(t, real_t(8), 1e-7));
Box box3(128, 5, Vec3(5, 0, 0), Vec3(0, 0, 0), Quat(), Vec3(10, 10, 10), iron, false, true, false);
WALBERLA_CHECK(intersects(&box3, &ray1, &t));
WALBERLA_CHECK(realIsEqual(t, real_t(5)));
Ray ray2(Vec3(-2,0,0), Vec3(1,0,1));
WALBERLA_LOG_INFO(intersects(&box3, &ray2, &t));
WALBERLA_LOG_INFO(t);
Ray ray3(Vec3(-5,3,3), Vec3(2, -1.5, 0.5));
Box box4(128, 5, Vec3(8, 0, 0), Vec3(0, 0, 0), Quat(), Vec3(10, 10, 10), iron, false, true, false);
WALBERLA_LOG_INFO("-- -- -- --");
WALBERLA_LOG_INFO(intersects(&box4, &ray3, &t));
WALBERLA_LOG_INFO("t: " << t);
Ray ray4(Vec3(-5,3,3), Vec3(2, -1, 0.5));
WALBERLA_LOG_INFO("-- -- -- --");
WALBERLA_LOG_INFO(intersects(&box4, &ray4, &t));
WALBERLA_LOG_INFO("t: " << t);
}
int main( int argc, char** argv )
{
walberla::debug::enterTestMode();
......@@ -86,6 +120,7 @@ int main( int argc, char** argv )
SphereIntersectsTest();
PlaneIntersectsTest();
//BoxIntersectsTest(); // does not work yet
return EXIT_SUCCESS;
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment