22 res = q * Quaternion(p) * conjugate(q); |
22 res = q * Quaternion(p) * conjugate(q); |
23 p = res.toVector(); |
23 p = res.toVector(); |
24 */ |
24 */ |
25 |
25 |
26 // optimized |
26 // optimized |
27 float t2 = q.a*q.b; |
27 Float t2 = q.a*q.b; |
28 float t3 = q.a*q.c; |
28 Float t3 = q.a*q.c; |
29 float t4 = q.a*q.d; |
29 Float t4 = q.a*q.d; |
30 float t5 = -q.b*q.b; |
30 Float t5 = -q.b*q.b; |
31 float t6 = q.b*q.c; |
31 Float t6 = q.b*q.c; |
32 float t7 = q.b*q.d; |
32 Float t7 = q.b*q.d; |
33 float t8 = -q.c*q.c; |
33 Float t8 = -q.c*q.c; |
34 float t9 = q.c*q.d; |
34 Float t9 = q.c*q.d; |
35 float t10 = -q.d*q.d; |
35 Float t10 = -q.d*q.d; |
36 float x,y,z; |
36 Float x,y,z; |
37 x = 2*( (t8 + t10)*p.x + (t6 - t4)*p.y + (t3 + t7)*p.z ) + p.x; |
37 x = 2*( (t8 + t10)*p.x + (t6 - t4)*p.y + (t3 + t7)*p.z ) + p.x; |
38 y = 2*( (t4 + t6)*p.x + (t5 + t10)*p.y + (t9 - t2)*p.z ) + p.y; |
38 y = 2*( (t4 + t6)*p.x + (t5 + t10)*p.y + (t9 - t2)*p.z ) + p.y; |
39 z = 2*( (t7 - t3)*p.x + (t2 + t9)*p.y + (t5 + t8)*p.z ) + p.z; |
39 z = 2*( (t7 - t3)*p.x + (t2 + t9)*p.y + (t5 + t8)*p.z ) + p.z; |
40 p = Vector3(x,y,z); |
40 p = Vector3(x,y,z); |
41 x = 2*( (t8 + t10)*u.x + (t6 - t4)*u.y + (t3 + t7)*u.z ) + u.x; |
41 x = 2*( (t8 + t10)*u.x + (t6 - t4)*u.y + (t3 + t7)*u.z ) + u.x; |
49 p.normalize(); |
49 p.normalize(); |
50 u.normalize(); |
50 u.normalize(); |
51 v.normalize(); |
51 v.normalize(); |
52 } |
52 } |
53 |
53 |
54 void Camera::move(const float fw, const float left, const float up) |
54 void Camera::move(const Float fw, const Float left, const Float up) |
55 { |
55 { |
56 eye = eye + fw*p + left*u + up*v; |
56 eye = eye + fw*p + left*u + up*v; |
57 } |
57 } |
58 |
58 |
59 /* http://www.siggraph.org/education/materials/HyperGraph/raytrace/rtinter3.htm */ |
59 /* http://www.siggraph.org/education/materials/HyperGraph/raytrace/rtinter3.htm */ |
60 bool BBox::intersect(const Ray &ray, float &a, float &b) |
60 bool BBox::intersect(const Ray &ray, Float &a, Float &b) |
61 { |
61 { |
62 float tnear = -FLT_MAX; |
62 Float tnear = -FLT_MAX; |
63 float tfar = FLT_MAX; |
63 Float tfar = FLT_MAX; |
64 float t1, t2; |
64 Float t1, t2; |
65 |
65 |
66 for (int i = 0; i < 3; i++) |
66 for (int i = 0; i < 3; i++) |
67 { |
67 { |
68 if (ray.dir.cell[i] == 0) { |
68 if (ray.dir.cell[i] == 0) { |
69 /* ray is parallel to these planes */ |
69 /* ray is parallel to these planes */ |
92 a = tnear; |
92 a = tnear; |
93 b = tfar; |
93 b = tfar; |
94 return true; |
94 return true; |
95 } |
95 } |
96 |
96 |
97 bool Sphere::intersect(const Ray &ray, float &dist) |
97 bool Sphere::intersect(const Ray &ray, Float &dist) |
98 { |
98 { |
99 Vector3 V = ray.o - center; |
99 Vector3 V = ray.o - center; |
100 register float d = -dot(V, ray.dir); |
100 register Float d = -dot(V, ray.dir); |
101 register float Det = d * d - (dot(V,V) - sqr_radius); |
101 register Float Det = d * d - (dot(V,V) - sqr_radius); |
102 if (Det > 0) { |
102 if (Det > 0) { |
103 d -= sqrtf(Det); |
103 d -= sqrtf(Det); |
104 if (d > 0 && d < dist) |
104 if (d > 0 && d < dist) |
105 { |
105 { |
106 dist = d; |
106 dist = d; |
108 } |
108 } |
109 } |
109 } |
110 return false; |
110 return false; |
111 } |
111 } |
112 |
112 |
113 bool Sphere::intersect_all(const Ray &ray, float dist, vector<float> &allts) |
113 bool Sphere::intersect_all(const Ray &ray, Float dist, vector<Float> &allts) |
114 { |
114 { |
115 //allts = new vector<float>(); |
115 //allts = new vector<Float>(); |
116 |
116 |
117 Vector3 V = ((Ray)ray).o - center; |
117 Vector3 V = ((Ray)ray).o - center; |
118 float Vd = - dot(V, ray.dir); |
118 Float Vd = - dot(V, ray.dir); |
119 float Det = Vd * Vd - (dot(V,V) - sqr_radius); |
119 Float Det = Vd * Vd - (dot(V,V) - sqr_radius); |
120 |
120 |
121 if (Det > 0) { |
121 if (Det > 0) { |
122 Det = sqrtf(Det); |
122 Det = sqrtf(Det); |
123 float t1 = Vd - Det; |
123 Float t1 = Vd - Det; |
124 float t2 = Vd + Det; |
124 Float t2 = Vd + Det; |
125 if (t1 < 0) |
125 if (t1 < 0) |
126 { |
126 { |
127 if (t2 > 0) |
127 if (t2 > 0) |
128 { |
128 { |
129 // ray from inside of the sphere |
129 // ray from inside of the sphere |
150 bbox.L = center - radius; |
150 bbox.L = center - radius; |
151 bbox.H = center + radius; |
151 bbox.H = center + radius; |
152 return bbox; |
152 return bbox; |
153 } |
153 } |
154 |
154 |
155 bool Box::intersect(const Ray &ray, float &dist) |
155 bool Box::intersect(const Ray &ray, Float &dist) |
156 { |
156 { |
157 float a,b; |
157 Float a,b; |
158 bool res = get_bbox().intersect(ray, a, b); |
158 bool res = get_bbox().intersect(ray, a, b); |
159 if (res && a < dist) |
159 if (res && a < dist) |
160 { |
160 { |
161 dist = a; |
161 dist = a; |
162 return true; |
162 return true; |
214 } |
214 } |
215 |
215 |
216 int u = (k + 1) % 3; |
216 int u = (k + 1) % 3; |
217 int v = (k + 2) % 3; |
217 int v = (k + 2) % 3; |
218 |
218 |
219 float krec = 1.0f / N[k]; |
219 Float krec = 1.0f / N[k]; |
220 nu = N[u] * krec; |
220 nu = N[u] * krec; |
221 nv = N[v] * krec; |
221 nv = N[v] * krec; |
222 nd = dot(N, A) * krec; |
222 nd = dot(N, A) * krec; |
223 |
223 |
224 // first line equation |
224 // first line equation |
225 float reci = 1.0f / (b[u] * c[v] - b[v] * c[u]); |
225 Float reci = 1.0f / (b[u] * c[v] - b[v] * c[u]); |
226 bnu = b[u] * reci; |
226 bnu = b[u] * reci; |
227 bnv = -b[v] * reci; |
227 bnv = -b[v] * reci; |
228 |
228 |
229 // second line equation |
229 // second line equation |
230 cnu = c[v] * reci; |
230 cnu = c[v] * reci; |
233 // finalize normal |
233 // finalize normal |
234 N.normalize(); |
234 N.normalize(); |
235 } |
235 } |
236 |
236 |
237 // see comment for previous method |
237 // see comment for previous method |
238 bool Triangle::intersect(const Ray &ray, float &dist) |
238 bool Triangle::intersect(const Ray &ray, Float &dist) |
239 { |
239 { |
240 Vector3 O = ray.o; |
240 Vector3 O = ray.o; |
241 Vector3 D = ray.dir; |
241 Vector3 D = ray.dir; |
242 |
242 |
243 const int modulo3[5] = {0,1,2,0,1}; |
243 const int modulo3[5] = {0,1,2,0,1}; |
244 const int ku = modulo3[k+1]; |
244 const int ku = modulo3[k+1]; |
245 const int kv = modulo3[k+2]; |
245 const int kv = modulo3[k+2]; |
246 const float lnd = 1.0f / (D[k] + nu * D[ku] + nv * D[kv]); |
246 const Float lnd = 1.0f / (D[k] + nu * D[ku] + nv * D[kv]); |
247 const float t = (nd - O[k] - nu * O[ku] - nv * O[kv]) * lnd; |
247 const Float t = (nd - O[k] - nu * O[ku] - nv * O[kv]) * lnd; |
248 |
248 |
249 if (!(t < dist && t > 0)) |
249 if (!(t < dist && t > 0)) |
250 return false; |
250 return false; |
251 |
251 |
252 float hu = O[ku] + t * D[ku] - A[ku]; |
252 Float hu = O[ku] + t * D[ku] - A[ku]; |
253 float hv = O[kv] + t * D[kv] - A[kv]; |
253 Float hv = O[kv] + t * D[kv] - A[kv]; |
254 float beta = hv * bnu + hu * bnv; |
254 Float beta = hv * bnu + hu * bnv; |
255 |
255 |
256 if (beta < 0) |
256 if (beta < 0) |
257 return false; |
257 return false; |
258 |
258 |
259 float gamma = hu * cnu + hv * cnv; |
259 Float gamma = hu * cnu + hv * cnv; |
260 if (gamma < 0) |
260 if (gamma < 0) |
261 return false; |
261 return false; |
262 |
262 |
263 if ((beta + gamma) > 1) |
263 if ((beta + gamma) > 1) |
264 return false; |
264 return false; |