Ellipsoids finally added to the ray-tracing engine

It took me a little while to correctly implement the ray-ellipsoid intersection function but there it is :-)

float4 O_C = ray.origin-ellipsoid.center;

float4 dir = ray.direction;

normalizeVector( dir );

float a =

        ((dir.x*dir.x)/(ellipsoid.size.x*ellipsoid.size.x))

      + ((dir.y*dir.y)/(ellipsoid.size.y*ellipsoid.size.y))

      + ((dir.z*dir.z)/(ellipsoid.size.z*ellipsoid.size.z));

float b =

        ((2.f*O_C.x*dir.x)/(ellipsoid.size.x*ellipsoid.size.x))

      + ((2.f*O_C.y*dir.y)/(ellipsoid.size.y*ellipsoid.size.y))

      + ((2.f*O_C.z*dir.z)/(ellipsoid.size.z*ellipsoid.size.z));

float c =

        ((O_C.x*O_C.x)/(ellipsoid.size.x*ellipsoid.size.x))

      + ((O_C.y*O_C.y)/(ellipsoid.size.y*ellipsoid.size.y))

      + ((O_C.z*O_C.z)/(ellipsoid.size.z*ellipsoid.size.z))

      - 1.f;

float d = ((b*b)-(4.f*a*c));

if ( d<0.f || a==0.f || b==0.f || c==0.f )

   return false;

d = sqrt(d);

float t1 = (-b+d)/(2.f*a);

float t2 = (-b-d)/(2.f*a);

if( t1<=EPSILON && t2<=EPSILON ) return false; // both intersections are behind the ray origin

back = (t1<=EPSILON || t2<=EPSILON); // If only one intersection (t>0) then we are inside the ellipsoid and the intersection is at the back of the ellipsoid

float t=0.f;

if( t1<=EPSILON )

   t = t2;

else

   if( t2<=EPSILON )

      t = t1;

   else

      t=(t1<t2) ? t1 : t2;

if( t<EPSILON ) return false; // Too close to intersection

intersection = ray.origin + t*dir;

normal = intersection-ellipsoid.center;

normal.x = 2.f*normal.x/(ellipsoid.size.x*ellipsoid.size.x);

normal.y = 2.f*normal.y/(ellipsoid.size.y*ellipsoid.size.y);

normal.z = 2.f*normal.z/(ellipsoid.size.z*ellipsoid.size.z);

normal.w = 0.f;

normal *= (back) ? -1.f : 1.f;

normalizeVector(normal);

return true;