1 #ifndef CYLINDER_HXX
2 #define CYLINDER_HXX
3
4 #include "Primitive.hxx"
5 #include "InfinitePlane.hxx"
6
7 /**
8 * axis aligned cylinder
9 */
10 typedef enum align
11 {
12 xAxis,
13 yAxis,
14 zAxis
15 };
16
17 class Cylinder : public Primitive
18 {
19 Vec3f _down;
20 Vec3f _up;
21 Vec3f _center;
22
23 enum align _dim;
24 float radius;
25 float height;
26
27 /**
28 * compute the box of a cylinder
29 * therefore take the center of the sphere upwards and do the same as for the sphere
30 * and then continue for the center of the sphere downwards
31 */
32 Box CalcBounds()
33 {
34 Box BoundingBox;
35
36 Vec3f circle = Vec3f( radius );
37 circle [ _dim ] = 0;
38
39 BoundingBox.Extend( _down + circle);
40 BoundingBox.Extend( _down - circle);
41
42 BoundingBox.Extend( _up + circle );
43 BoundingBox.Extend( _up - circle );
44
45 return BoundingBox;
46 }
47
48 public:
49
50 Cylinder(Vec3f point, enum align _dim, float radius, float height)
51 : _dim(_dim), radius(radius), height(height)
52 {
53 // unit vector in _dim direction
54 Vec3f add = Vec3f(0);
55 add[_dim] = 1;
56
57 // compute center of the cylinder
58 _center = point;
59 // center of the sphere downwards
60 _down = _center - height * add;
61 // center of the sphere upwards
62 _up = _center + height * add;
63 };
64
65 virtual ~Cylinder()
66 {}
67
68 virtual bool Intersect(Ray &ray)
69 {
70 /**
71 * first: infinite long cylinder
72 * computation as for sphere, just one dimension missing
73 */
74
75 /**
76 * x^2 + y^2 = radius^2
77 */
78 Vec3f diff = ray.org - _down;
79 Vec3f dir = ray.dir;
80
81 diff[_dim] = 0.0;
82 dir [_dim] = 0.0;
83
84 float a = Dot(dir,dir);
85 float b = 2 * Dot(dir,diff);
86 float c = Dot(diff,diff) - radius * radius;
87
88 float inRoot = b*b - 4*a*c;
89 if (inRoot < 0)
90 return false;
91 float root = sqrt(inRoot);
92
93 float dist = (-b - root)/(2*a);
94 if (dist > ray.t)
95 return false;
96
97 if (dist < Epsilon)
98 {
99 dist = (-b + root)/(2*a);
100 if (dist < Epsilon || dist > ray.t)
101 return false;
102 }
103
104 /**
105 * side
106 * cut off, so that a 2*height cylinder stays
107 * which is open to both sides
108 */
109 Vec3f hit = ray.org + dist * ray.dir;
110 float t1 = (hit - _center)[_dim];
111
112 if ( - height < t1 && t1 < height )
113 {
114 ray.t = dist;
115 ray.hit = this;
116
117 return true;
118 }
119
120 /**
121 * compute top and bottom surfaces => close the cylinder
122 * check if ray hits planes located at _down and _up within the radius of the cylinder
123 */
124 Ray ray1 = ray;
125 Ray ray2 = ray;
126
127 Vec3f normal = Vec3f(0);
128 normal[_dim] = 1.0;
129
130 Primitive* plane1 = new InfinitePlane( _down, normal );
131 bool hit_plane1 = plane1->Intersect(ray1);
132 delete plane1;
133
134 Primitive* plane2 = new InfinitePlane( _up, normal );
135 bool hit_plane2 = plane2->Intersect(ray2);
136 delete plane2;
137
138 // bottom surface
139 if (hit_plane1)
140 {
141 Vec3f hit1 = (ray1.org + ray1.t * ray1.dir) - _down;
142 hit1[_dim] = 0.0;
143
144 if ( Length(hit1) < radius * radius )
145 {
146 ray.t = ray1.t;
147 ray.hit = this;
148
149 return true;
150 }
151 }
152
153 // top surface
154 if (hit_plane2)
155 {
156 Vec3f hit2 = (ray2.org + ray2.t * ray2.dir) - _up;
157 hit2[_dim] = 0.0;
158
159 if ( Length(hit2) < radius * radius )
160 {
161 ray.t = ray2.t;
162 ray.hit = this;
163
164 return true;
165 }
166 }
167
168 return false;
169 };
170
171 virtual Vec3f GetNormal(Ray &ray)
172 {
173 Vec3f hit = (ray.org + ray.t * ray.dir) - _center;
174 Vec3f normal = Vec3f(0);
175
176 // normal vector of bottom surface
177 if ( fabs(hit[_dim] - height) < Epsilon)
178 {
179 normal[_dim] = +1.0;
180 }
181 // normal vector of top surface
182 else if ( fabs(hit[_dim] + height) < Epsilon)
183 {
184 normal[_dim] = -1.0;
185 }
186 // normal vector of side
187 else
188 {
189 normal[_dim] = hit[_dim];
190 normal = hit - normal;
191 Normalize(normal);
192 }
193
194 return normal;
195 };
196
197 };
198
199 #endif
syntax highlighted by Code2HTML, v. 0.9.1