1 #ifndef SPHERE_HXX
 2 #define SPHERE_HXX
 3 
 4 #include "Primitive.hxx"
 5 #include "Shader.hxx"
 6 #include "PhotonDistributor.hxx"
 7 
 8 class Sphere : public Primitive
 9 {
10   Vec3f center;
11   float radius;
12   Box   box;
13 public:
14   Sphere(Vec3f center, float radius, Shader *shader, PhotonDistributor *photon_distributor)
15     : Primitive(shader, photon_distributor),center(center),radius(radius)
16   {
17     box.Extend(center + Vec3f(radius)); 
18     box.Extend(center - Vec3f(radius)); 
19   };
20 
21   virtual bool Intersect(Ray &ray)
22   {
23     // mathematical derivation, numerically not very stable, but simple
24 
25     // --> find roots of f(t) = ((R+tD)-C)^2 - r^2
26     // f(t) = (R-C)^2 + 2(R-C)(tD) + (tD)^2 -r^2
27     // --> f(t) = [D^2] t^2 + [2D(R-C)] t + [(R-C)^2 - r^2]
28     Vec3f diff = ray.org - center;
29     float a = Dot(ray.dir,ray.dir);
30     float b = 2 * Dot(ray.dir,diff);
31     float c = Dot(diff,diff) - radius * radius;
32 
33     // use 'abc'-formula for finding root t_1,2 = (-b +/- sqrt(b^2-4ac))/(2a)
34     float inRoot = b*b - 4*a*c;
35     if (inRoot < 0) return false;
36     float root = sqrt(inRoot);
37 
38     float dist = (-b - root)/(2*a);
39     if (dist > ray.t)
40       return false;
41 
42     if (dist < Epsilon) {
43       dist = (-b + root)/(2*a);
44       if (dist < Epsilon || dist > ray.t)
45 	return false;
46     }
47 
48     ray.t = dist;
49     ray.hit = this;
50     return true;
51   };
52 
53   virtual Vec3f GetNormal(Ray &ray)
54   {
55     Vec3f hit = ray.org + ray.t * ray.dir;
56     Vec3f normal = hit - center;
57     Normalize(normal);
58     return normal;
59   };
60 
61   virtual Box CalcBounds()
62   {
63 	return box;
64   }
65 };
66 
67 #endif


syntax highlighted by Code2HTML, v. 0.9.1