1 | |
2 | // |
3 | // This source file is part of appleseed. |
4 | // Visit http://appleseedhq.net/ for additional information and resources. |
5 | // |
6 | // This software is released under the MIT license. |
7 | // |
8 | // Copyright (c) 2010-2013 Francois Beaune, Jupiter Jazz Limited |
9 | // Copyright (c) 2014-2017 Francois Beaune, The appleseedhq Organization |
10 | // |
11 | // Permission is hereby granted, free of charge, to any person obtaining a copy |
12 | // of this software and associated documentation files (the "Software"), to deal |
13 | // in the Software without restriction, including without limitation the rights |
14 | // to use, copy, modify, merge, publish, distribute, sublicense, and/or sell |
15 | // copies of the Software, and to permit persons to whom the Software is |
16 | // furnished to do so, subject to the following conditions: |
17 | // |
18 | // The above copyright notice and this permission notice shall be included in |
19 | // all copies or substantial portions of the Software. |
20 | // |
21 | // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR |
22 | // IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, |
23 | // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE |
24 | // AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER |
25 | // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, |
26 | // OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN |
27 | // THE SOFTWARE. |
28 | // |
29 | |
30 | #ifndef APPLESEED_FOUNDATION_MATH_VECTOR_H |
31 | #define APPLESEED_FOUNDATION_MATH_VECTOR_H |
32 | |
33 | // appleseed.foundation headers. |
34 | #include "foundation/math/scalar.h" |
35 | #include "foundation/utility/poison.h" |
36 | |
37 | // Imath headers. |
38 | #ifdef APPLESEED_ENABLE_IMATH_INTEROP |
39 | #include "foundation/platform/exrheaderguards.h" |
40 | BEGIN_EXR_INCLUDES |
41 | #include "OpenEXR/ImathVec.h" |
42 | END_EXR_INCLUDES |
43 | #endif |
44 | |
45 | // Standard headers. |
46 | #include <algorithm> |
47 | #include <cassert> |
48 | #include <cmath> |
49 | #include <cstddef> |
50 | |
51 | namespace foundation |
52 | { |
53 | |
54 | // |
55 | // N-dimensional vector class of arbitrary type. |
56 | // |
57 | |
58 | template <typename T, size_t N> |
59 | class Vector |
60 | { |
61 | public: |
62 | // Types. |
63 | typedef T ValueType; |
64 | typedef Vector<T, N> VectorType; |
65 | |
66 | // Dimension. |
67 | static const size_t Dimension = N; |
68 | |
69 | // Constructors. |
70 | Vector(); // leave all components uninitialized |
71 | explicit Vector(const ValueType* rhs); // initialize with array of N scalars |
72 | explicit Vector(const ValueType val); // set all components to 'val' |
73 | |
74 | // Construct a vector from another vector of a different type. |
75 | template <typename U> |
76 | explicit Vector(const Vector<U, N>& rhs); |
77 | |
78 | // Unchecked array subscripting. |
79 | ValueType& operator[](const size_t i); |
80 | const ValueType& operator[](const size_t i) const; |
81 | |
82 | private: |
83 | // Vector components. |
84 | ValueType m_comp[N]; |
85 | }; |
86 | |
87 | // Poisoning. |
88 | template <typename T, size_t N> |
89 | class PoisonImpl<Vector<T, N> > |
90 | { |
91 | public: |
92 | static void do_poison(Vector<T, N>& v); |
93 | }; |
94 | |
95 | // Exact inequality and equality tests. |
96 | template <typename T, size_t N> bool operator!=(const Vector<T, N>& lhs, const Vector<T, N>& rhs); |
97 | template <typename T, size_t N> bool operator==(const Vector<T, N>& lhs, const Vector<T, N>& rhs); |
98 | |
99 | // Approximate equality tests. |
100 | template <typename T, size_t N> bool feq(const Vector<T, N>& lhs, const Vector<T, N>& rhs); |
101 | template <typename T, size_t N> bool feq(const Vector<T, N>& lhs, const Vector<T, N>& rhs, const T eps); |
102 | |
103 | // Approximate zero tests. |
104 | template <typename T, size_t N> bool fz(const Vector<T, N>& v); |
105 | template <typename T, size_t N> bool fz(const Vector<T, N>& v, const T eps); |
106 | |
107 | // Vector arithmetic. |
108 | template <typename T, size_t N> Vector<T, N> operator+ (const Vector<T, N>& lhs, const Vector<T, N>& rhs); |
109 | template <typename T, size_t N> Vector<T, N> operator- (const Vector<T, N>& lhs, const Vector<T, N>& rhs); |
110 | template <typename T, size_t N> Vector<T, N> operator- (const Vector<T, N>& lhs); |
111 | template <typename T, size_t N> Vector<T, N> operator* (const Vector<T, N>& lhs, const T rhs); |
112 | template <typename T, size_t N> Vector<T, N> operator* (const T lhs, const Vector<T, N>& rhs); |
113 | template <typename T, size_t N> Vector<T, N> operator/ (const Vector<T, N>& lhs, const T rhs); |
114 | template <typename T, size_t N> Vector<T, N> operator* (const Vector<T, N>& lhs, const Vector<T, N>& rhs); |
115 | template <typename T, size_t N> Vector<T, N> operator/ (const Vector<T, N>& lhs, const Vector<T, N>& rhs); |
116 | template <typename T, size_t N> Vector<T, N>& operator+=(Vector<T, N>& lhs, const Vector<T, N>& rhs); |
117 | template <typename T, size_t N> Vector<T, N>& operator-=(Vector<T, N>& lhs, const Vector<T, N>& rhs); |
118 | template <typename T, size_t N> Vector<T, N>& operator*=(Vector<T, N>& lhs, const T rhs); |
119 | template <typename T, size_t N> Vector<T, N>& operator/=(Vector<T, N>& lhs, const T rhs); |
120 | template <typename T, size_t N> Vector<T, N>& operator*=(Vector<T, N>& lhs, const Vector<T, N>& rhs); |
121 | template <typename T, size_t N> Vector<T, N>& operator/=(Vector<T, N>& lhs, const Vector<T, N>& rhs); |
122 | |
123 | // Dot product. |
124 | template <typename T, size_t N> T dot(const Vector<T, N>& lhs, const Vector<T, N>& rhs); |
125 | |
126 | // Vector square norm, norm and normalization. |
127 | template <typename T, size_t N> T square_norm(const Vector<T, N>& v); |
128 | template <typename T, size_t N> T norm(const Vector<T, N>& v); |
129 | template <typename T, size_t N> Vector<T, N> normalize(const Vector<T, N>& v); |
130 | template <typename T, size_t N> Vector<T, N> safe_normalize(const Vector<T, N>& v, const Vector<T, N>& fallback); |
131 | template <typename T, size_t N> Vector<T, N> safe_normalize(const Vector<T, N>& v); |
132 | |
133 | // Bring the norm of a nearly-unit vector closer to 1 by performing a single Newton-Raphson step. |
134 | template <typename T, size_t N> |
135 | Vector<T, N> improve_normalization(const Vector<T, N>& v); |
136 | |
137 | // Bring the norm of a nearly-unit vector closer to 1 by performing a set number of Newton-Raphson steps. |
138 | template <size_t Steps, typename T, size_t N> |
139 | Vector<T, N> improve_normalization(const Vector<T, N>& v); |
140 | |
141 | // Return true if a vector is normalized (unit-length), false otherwise. |
142 | template <typename T, size_t N> bool is_normalized(const Vector<T, N>& v); |
143 | template <typename T, size_t N> bool is_normalized(const Vector<T, N>& v, const T eps); |
144 | |
145 | // Return n if dot(n, i) < 0, -n otherwise. |
146 | template <typename T, size_t N> |
147 | Vector<T, N> faceforward( |
148 | const Vector<T, N>& n, |
149 | const Vector<T, N>& i); |
150 | |
151 | // Return v if it is in the same hemisphere as ref, -v otherwise. |
152 | template <typename T, size_t N> |
153 | Vector<T, N> flip_to_same_hemisphere( |
154 | const Vector<T, N>& v, |
155 | const Vector<T, N>& ref); |
156 | |
157 | // Return v with the component along n zeroed. |
158 | template <typename T, size_t N> |
159 | Vector<T, N> project( |
160 | const Vector<T, N>& v, |
161 | const Vector<T, N>& n); |
162 | |
163 | // Return the reflection vector given an incoming vector i and a unit-length normal vector n. |
164 | template <typename T, size_t N> |
165 | Vector<T, N> reflect( |
166 | const Vector<T, N>& i, |
167 | const Vector<T, N>& n); |
168 | |
169 | // |
170 | // Sets t to the refracted vector given a unit-length incoming vector i, a unit-length normal |
171 | // vector n and a relative index of refraction eta. eta is the ratio of the index of |
172 | // refraction in the volume containing the incoming vector to that of the volume being entered. |
173 | // Returns false in the case of total internal reflection. |
174 | // |
175 | // Geometry: |
176 | // |
177 | // i n eta_i |
178 | // eta = ----- |
179 | // ^ ^ eta_t |
180 | // \ | |
181 | // \ | |
182 | // \ | |
183 | // \ | IOR = eta_i |
184 | // \ | |
185 | // \| |
186 | // ---------------+--------------- interface |
187 | // | |
188 | // | |
189 | // | IOR = eta_t |
190 | // | |
191 | // | |
192 | // | |
193 | // v |
194 | // |
195 | // t = refract(i, n, eta_i / eta_t) |
196 | // |
197 | |
198 | template <typename T, size_t N> |
199 | bool refract( |
200 | const Vector<T, N>& i, |
201 | const Vector<T, N>& n, |
202 | const T eta, |
203 | Vector<T, N>& t); |
204 | |
205 | // Clamp the argument to [min, max]. |
206 | template <typename T, size_t N> |
207 | Vector<T, N> clamp( |
208 | const Vector<T, N>& v, |
209 | const T min, |
210 | const T max); |
211 | |
212 | // Clamp the argument to [0,1]. |
213 | template <typename T, size_t N> Vector<T, N> saturate(const Vector<T, N>& v); |
214 | |
215 | // Return the smallest or largest signed component of a vector. |
216 | template <typename T, size_t N> T min_value(const Vector<T, N>& v); |
217 | template <typename T, size_t N> T max_value(const Vector<T, N>& v); |
218 | |
219 | // Return the index of the smallest or largest signed component of a vector. |
220 | template <typename T, size_t N> size_t min_index(const Vector<T, N>& v); |
221 | template <typename T, size_t N> size_t max_index(const Vector<T, N>& v); |
222 | |
223 | // Return the index of the smallest or largest component of a vector, in absolute value. |
224 | template <typename T, size_t N> size_t min_abs_index(const Vector<T, N>& v); |
225 | template <typename T, size_t N> size_t max_abs_index(const Vector<T, N>& v); |
226 | |
227 | // Component-wise min/max of two vectors. |
228 | template <typename T, size_t N> Vector<T, N> component_wise_min(const Vector<T, N>& lhs, const Vector<T, N>& rhs); |
229 | template <typename T, size_t N> Vector<T, N> component_wise_max(const Vector<T, N>& lhs, const Vector<T, N>& rhs); |
230 | |
231 | |
232 | // |
233 | // 2-dimensional vector class of arbitrary type. |
234 | // |
235 | |
236 | template <typename T> |
237 | class Vector<T, 2> |
238 | { |
239 | public: |
240 | // Types. |
241 | typedef T ValueType; |
242 | typedef Vector<T, 2> VectorType; |
243 | |
244 | // Dimension. |
245 | static const size_t Dimension = 2; |
246 | |
247 | // Vector components. |
248 | ValueType x, y; |
249 | |
250 | // Constructors. |
251 | Vector(); // leave all components uninitialized |
252 | explicit Vector(const ValueType* rhs); // initialize with array of 2 scalars |
253 | explicit Vector(const ValueType val); // set all components to 'val' |
254 | Vector( // set individual components |
255 | const ValueType x, |
256 | const ValueType y); |
257 | |
258 | // Construct a vector from another vector of a different type. |
259 | template <typename U> |
260 | explicit Vector(const Vector<U, 2>& rhs); |
261 | |
262 | #ifdef APPLESEED_ENABLE_IMATH_INTEROP |
263 | |
264 | // Implicit construction from an Imath::Vec2. |
265 | Vector(const Imath::Vec2<T>& rhs); |
266 | |
267 | // Reinterpret this vector as an Imath::Vec2. |
268 | operator Imath::Vec2<T>&(); |
269 | operator const Imath::Vec2<T>&() const; |
270 | |
271 | #endif |
272 | |
273 | // Unchecked array subscripting. |
274 | ValueType& operator[](const size_t i); |
275 | const ValueType& operator[](const size_t i) const; |
276 | }; |
277 | |
278 | // Determinant of the 2D matrix whose first column is lhs and second column is rhs. |
279 | template <typename T> T det(const Vector<T, 2>& lhs, const Vector<T, 2>& rhs); |
280 | |
281 | |
282 | // |
283 | // 3-dimensional vector class of arbitrary type. |
284 | // |
285 | |
286 | template <typename T> |
287 | class Vector<T, 3> |
288 | { |
289 | public: |
290 | // Types. |
291 | typedef T ValueType; |
292 | typedef Vector<T, 3> VectorType; |
293 | |
294 | // Dimension. |
295 | static const size_t Dimension = 3; |
296 | |
297 | // Vector components. |
298 | ValueType x, y, z; |
299 | |
300 | // Constructors. |
301 | Vector(); // leave all components uninitialized |
302 | explicit Vector(const ValueType* rhs); // initialize with array of 3 scalars |
303 | explicit Vector(const ValueType val); // set all components to 'val' |
304 | Vector( // set individual components |
305 | const ValueType x, |
306 | const ValueType y, |
307 | const ValueType z); |
308 | |
309 | // Construct a vector from another vector of a different type. |
310 | template <typename U> |
311 | explicit Vector(const Vector<U, 3>& rhs); |
312 | |
313 | #ifdef APPLESEED_ENABLE_IMATH_INTEROP |
314 | |
315 | // Implicit construction from an Imath::Vec3. |
316 | Vector(const Imath::Vec3<T>& rhs); |
317 | |
318 | // Reinterpret this vector as an Imath::Vec3. |
319 | operator Imath::Vec3<T>&(); |
320 | operator const Imath::Vec3<T>&() const; |
321 | |
322 | #endif |
323 | |
324 | // Build a unit vector from two angles. |
325 | static VectorType make_unit_vector( |
326 | const ValueType theta, // angle with Y basis vector, in radians |
327 | const ValueType phi); // angle with X basis vector, in radians |
328 | static VectorType make_unit_vector( |
329 | const ValueType cos_theta, // cosine of angle with Y basis vector |
330 | const ValueType sin_theta, // sine of angle with Y basis vector |
331 | const ValueType cos_phi, // cosine of angle with X basis vector |
332 | const ValueType sin_phi); // sine of angle with X basis vector |
333 | |
334 | // Unchecked array subscripting. |
335 | ValueType& operator[](const size_t i); |
336 | const ValueType& operator[](const size_t i) const; |
337 | }; |
338 | |
339 | // Cross product. |
340 | template <typename T> Vector<T, 3> cross(const Vector<T, 3>& lhs, const Vector<T, 3>& rhs); |
341 | |
342 | |
343 | // |
344 | // 4-dimensional vector class of arbitrary type. |
345 | // |
346 | |
347 | template <typename T> |
348 | class Vector<T, 4> |
349 | { |
350 | public: |
351 | // Types. |
352 | typedef T ValueType; |
353 | typedef Vector<T, 4> VectorType; |
354 | |
355 | // Dimension. |
356 | static const size_t Dimension = 4; |
357 | |
358 | // Vector components. |
359 | ValueType x, y, z, w; |
360 | |
361 | // Constructors. |
362 | Vector(); // leave all components uninitialized |
363 | explicit Vector(const ValueType* rhs); // initialize with array of 4 scalars |
364 | explicit Vector(const ValueType val); // set all components to 'val' |
365 | Vector( // set individual components |
366 | const ValueType x, |
367 | const ValueType y, |
368 | const ValueType z, |
369 | const ValueType w); |
370 | |
371 | // Construct a vector from another vector of a different type. |
372 | template <typename U> |
373 | explicit Vector(const Vector<U, 4>& rhs); |
374 | |
375 | // Unchecked array subscripting. |
376 | ValueType& operator[](const size_t i); |
377 | const ValueType& operator[](const size_t i) const; |
378 | }; |
379 | |
380 | |
381 | // |
382 | // Full specializations for 1D, 2D, 3D and 4D vectors of type int, float and double. |
383 | // |
384 | |
385 | typedef Vector<int, 1> Vector1i; |
386 | typedef Vector<size_t, 1> Vector1u; |
387 | typedef Vector<float, 1> Vector1f; |
388 | typedef Vector<double, 1> Vector1d; |
389 | |
390 | typedef Vector<int, 2> Vector2i; |
391 | typedef Vector<size_t, 2> Vector2u; |
392 | typedef Vector<float, 2> Vector2f; |
393 | typedef Vector<double, 2> Vector2d; |
394 | |
395 | typedef Vector<int, 3> Vector3i; |
396 | typedef Vector<size_t, 3> Vector3u; |
397 | typedef Vector<float, 3> Vector3f; |
398 | typedef Vector<double, 3> Vector3d; |
399 | |
400 | typedef Vector<int, 4> Vector4i; |
401 | typedef Vector<size_t, 4> Vector4u; |
402 | typedef Vector<float, 4> Vector4f; |
403 | typedef Vector<double, 4> Vector4d; |
404 | |
405 | |
406 | // |
407 | // Utility class to find the Imath equivalent of a given foundation::Vector type. |
408 | // |
409 | |
410 | #ifdef APPLESEED_ENABLE_IMATH_INTEROP |
411 | |
412 | template <typename T, size_t N> |
413 | struct ImathVecEquivalent; // no equivalent in the general case |
414 | |
415 | template <typename T> |
416 | struct ImathVecEquivalent<T, 2> |
417 | { |
418 | typedef Imath::Vec2<T> Type; |
419 | }; |
420 | |
421 | template <typename T> |
422 | struct ImathVecEquivalent<T, 3> |
423 | { |
424 | typedef Imath::Vec3<T> Type; |
425 | }; |
426 | |
427 | #endif |
428 | |
429 | |
430 | // |
431 | // N-dimensional vector implementation. |
432 | // |
433 | |
434 | template <typename T, size_t N> |
435 | inline Vector<T, N>::Vector() |
436 | { |
437 | } |
438 | |
439 | template <typename T, size_t N> |
440 | inline Vector<T, N>::Vector(const ValueType* rhs) |
441 | { |
442 | assert(rhs); |
443 | |
444 | for (size_t i = 0; i < N; ++i) |
445 | m_comp[i] = rhs[i]; |
446 | } |
447 | |
448 | template <typename T, size_t N> |
449 | inline Vector<T, N>::Vector(const ValueType val) |
450 | { |
451 | for (size_t i = 0; i < N; ++i) |
452 | m_comp[i] = val; |
453 | } |
454 | |
455 | template <typename T, size_t N> |
456 | template <typename U> |
457 | inline Vector<T, N>::Vector(const Vector<U, N>& rhs) |
458 | { |
459 | for (size_t i = 0; i < N; ++i) |
460 | m_comp[i] = static_cast<ValueType>(rhs[i]); |
461 | } |
462 | |
463 | template <typename T, size_t N> |
464 | inline T& Vector<T, N>::operator[](const size_t i) |
465 | { |
466 | assert(i < Dimension); |
467 | return m_comp[i]; |
468 | } |
469 | |
470 | template <typename T, size_t N> |
471 | inline const T& Vector<T, N>::operator[](const size_t i) const |
472 | { |
473 | assert(i < Dimension); |
474 | return m_comp[i]; |
475 | } |
476 | |
477 | template <typename T, size_t N> |
478 | void PoisonImpl<Vector<T, N> >::do_poison(Vector<T, N>& v) |
479 | { |
480 | for (size_t i = 0; i < N; ++i) |
481 | poison(v[i]); |
482 | } |
483 | |
484 | template <typename T, size_t N> |
485 | inline bool operator!=(const Vector<T, N>& lhs, const Vector<T, N>& rhs) |
486 | { |
487 | for (size_t i = 0; i < N; ++i) |
488 | { |
489 | if (lhs[i] != rhs[i]) |
490 | return true; |
491 | } |
492 | |
493 | return false; |
494 | } |
495 | |
496 | template <typename T, size_t N> |
497 | inline bool operator==(const Vector<T, N>& lhs, const Vector<T, N>& rhs) |
498 | { |
499 | return !(lhs != rhs); |
500 | } |
501 | |
502 | template <typename T, size_t N> |
503 | inline bool feq(const Vector<T, N>& lhs, const Vector<T, N>& rhs) |
504 | { |
505 | for (size_t i = 0; i < N; ++i) |
506 | { |
507 | if (!feq(lhs[i], rhs[i])) |
508 | return false; |
509 | } |
510 | |
511 | return true; |
512 | } |
513 | |
514 | template <typename T, size_t N> |
515 | inline bool feq(const Vector<T, N>& lhs, const Vector<T, N>& rhs, const T eps) |
516 | { |
517 | for (size_t i = 0; i < N; ++i) |
518 | { |
519 | if (!feq(lhs[i], rhs[i], eps)) |
520 | return false; |
521 | } |
522 | |
523 | return true; |
524 | } |
525 | |
526 | template <typename T, size_t N> |
527 | inline bool fz(const Vector<T, N>& v) |
528 | { |
529 | for (size_t i = 0; i < N; ++i) |
530 | { |
531 | if (!fz(v[i])) |
532 | return false; |
533 | } |
534 | |
535 | return true; |
536 | } |
537 | |
538 | template <typename T, size_t N> |
539 | inline bool fz(const Vector<T, N>& v, const T eps) |
540 | { |
541 | for (size_t i = 0; i < N; ++i) |
542 | { |
543 | if (!fz(v[i], eps)) |
544 | return false; |
545 | } |
546 | |
547 | return true; |
548 | } |
549 | |
550 | template <typename T, size_t N> |
551 | inline Vector<T, N> operator+(const Vector<T, N>& lhs, const Vector<T, N>& rhs) |
552 | { |
553 | Vector<T, N> result; |
554 | |
555 | for (size_t i = 0; i < N; ++i) |
556 | result[i] = lhs[i] + rhs[i]; |
557 | |
558 | return result; |
559 | } |
560 | |
561 | template <typename T, size_t N> |
562 | inline Vector<T, N> operator-(const Vector<T, N>& lhs, const Vector<T, N>& rhs) |
563 | { |
564 | Vector<T, N> result; |
565 | |
566 | for (size_t i = 0; i < N; ++i) |
567 | result[i] = lhs[i] - rhs[i]; |
568 | |
569 | return result; |
570 | } |
571 | |
572 | template <typename T, size_t N> |
573 | inline Vector<T, N> operator-(const Vector<T, N>& lhs) |
574 | { |
575 | Vector<T, N> result; |
576 | |
577 | for (size_t i = 0; i < N; ++i) |
578 | result[i] = -lhs[i]; |
579 | |
580 | return result; |
581 | } |
582 | |
583 | template <typename T, size_t N> |
584 | inline Vector<T, N> operator*(const Vector<T, N>& lhs, const T rhs) |
585 | { |
586 | Vector<T, N> result; |
587 | |
588 | for (size_t i = 0; i < N; ++i) |
589 | result[i] = lhs[i] * rhs; |
590 | |
591 | return result; |
592 | } |
593 | |
594 | template <typename T, size_t N> |
595 | inline Vector<T, N> operator*(const T lhs, const Vector<T, N>& rhs) |
596 | { |
597 | return rhs * lhs; |
598 | } |
599 | |
600 | template <typename T, size_t N> |
601 | inline Vector<T, N> operator/(const Vector<T, N>& lhs, const T rhs) |
602 | { |
603 | Vector<T, N> result; |
604 | |
605 | for (size_t i = 0; i < N; ++i) |
606 | result[i] = lhs[i] / rhs; |
607 | |
608 | return result; |
609 | } |
610 | |
611 | template <size_t N> |
612 | inline Vector<float, N> operator/(const Vector<float, N>& lhs, const float rhs) |
613 | { |
614 | return lhs * (1.0f / rhs); |
615 | } |
616 | |
617 | template <size_t N> |
618 | inline Vector<double, N> operator/(const Vector<double, N>& lhs, const double rhs) |
619 | { |
620 | return lhs * (1.0 / rhs); |
621 | } |
622 | |
623 | template <size_t N> |
624 | inline Vector<long double, N> operator/(const Vector<long double, N>& lhs, const long double rhs) |
625 | { |
626 | return lhs * (1.0L / rhs); |
627 | } |
628 | |
629 | template <typename T, size_t N> |
630 | inline Vector<T, N> operator*(const Vector<T, N>& lhs, const Vector<T, N>& rhs) |
631 | { |
632 | Vector<T, N> result; |
633 | |
634 | for (size_t i = 0; i < N; ++i) |
635 | result[i] = lhs[i] * rhs[i]; |
636 | |
637 | return result; |
638 | } |
639 | |
640 | template <typename T, size_t N> |
641 | inline Vector<T, N> operator/(const Vector<T, N>& lhs, const Vector<T, N>& rhs) |
642 | { |
643 | Vector<T, N> result; |
644 | |
645 | for (size_t i = 0; i < N; ++i) |
646 | result[i] = lhs[i] / rhs[i]; |
647 | |
648 | return result; |
649 | } |
650 | |
651 | template <typename T, size_t N> |
652 | inline Vector<T, N>& operator+=(Vector<T, N>& lhs, const Vector<T, N>& rhs) |
653 | { |
654 | for (size_t i = 0; i < N; ++i) |
655 | lhs[i] += rhs[i]; |
656 | |
657 | return lhs; |
658 | } |
659 | |
660 | template <typename T, size_t N> |
661 | inline Vector<T, N>& operator-=(Vector<T, N>& lhs, const Vector<T, N>& rhs) |
662 | { |
663 | for (size_t i = 0; i < N; ++i) |
664 | lhs[i] -= rhs[i]; |
665 | |
666 | return lhs; |
667 | } |
668 | |
669 | template <typename T, size_t N> |
670 | inline Vector<T, N>& operator*=(Vector<T, N>& lhs, const T rhs) |
671 | { |
672 | for (size_t i = 0; i < N; ++i) |
673 | lhs[i] *= rhs; |
674 | |
675 | return lhs; |
676 | } |
677 | |
678 | template <typename T, size_t N> |
679 | inline Vector<T, N>& operator/=(Vector<T, N>& lhs, const T rhs) |
680 | { |
681 | for (size_t i = 0; i < N; ++i) |
682 | lhs[i] /= rhs; |
683 | |
684 | return lhs; |
685 | } |
686 | |
687 | template <size_t N> |
688 | inline Vector<float, N>& operator/=(Vector<float, N>& lhs, const float rhs) |
689 | { |
690 | return lhs *= 1.0f / rhs; |
691 | } |
692 | |
693 | template <size_t N> |
694 | inline Vector<double, N>& operator/=(Vector<double, N>& lhs, const double rhs) |
695 | { |
696 | return lhs *= 1.0 / rhs; |
697 | } |
698 | |
699 | template <size_t N> |
700 | inline Vector<long double, N>& operator/=(Vector<long double, N>& lhs, const long double rhs) |
701 | { |
702 | return lhs *= 1.0L / rhs; |
703 | } |
704 | |
705 | template <typename T, size_t N> |
706 | inline Vector<T, N>& operator*=(Vector<T, N>& lhs, const Vector<T, N>& rhs) |
707 | { |
708 | for (size_t i = 0; i < N; ++i) |
709 | lhs[i] *= rhs[i]; |
710 | |
711 | return lhs; |
712 | } |
713 | |
714 | template <typename T, size_t N> |
715 | inline Vector<T, N>& operator/=(Vector<T, N>& lhs, const Vector<T, N>& rhs) |
716 | { |
717 | for (size_t i = 0; i < N; ++i) |
718 | lhs[i] /= rhs[i]; |
719 | |
720 | return lhs; |
721 | } |
722 | |
723 | template <typename T, size_t N> |
724 | inline T dot(const Vector<T, N>& lhs, const Vector<T, N>& rhs) |
725 | { |
726 | T result = T(0.0); |
727 | |
728 | for (size_t i = 0; i < N; ++i) |
729 | result += lhs[i] * rhs[i]; |
730 | |
731 | return result; |
732 | } |
733 | |
734 | template <typename T, size_t N> |
735 | inline T square_norm(const Vector<T, N>& v) |
736 | { |
737 | return dot(v, v); |
738 | } |
739 | |
740 | template <typename T, size_t N> |
741 | inline T norm(const Vector<T, N>& v) |
742 | { |
743 | return std::sqrt(square_norm(v)); |
744 | } |
745 | |
746 | template <typename T, size_t N> |
747 | inline Vector<T, N> normalize(const Vector<T, N>& v) |
748 | { |
749 | const T n = norm(v); |
750 | assert(n > T(0.0)); |
751 | return v / n; |
752 | } |
753 | |
754 | template <typename T, size_t N> |
755 | inline Vector<T, N> safe_normalize(const Vector<T, N>& v, const Vector<T, N>& fallback) |
756 | { |
757 | assert(is_normalized(fallback)); |
758 | const T n = norm(v); |
759 | return n > 0.0 ? v / n : fallback; |
760 | } |
761 | |
762 | template <typename T, size_t N> |
763 | inline Vector<T, N> safe_normalize(const Vector<T, N>& v) |
764 | { |
765 | Vector<T, N> result = v; |
766 | |
767 | const T n = norm(result); |
768 | |
769 | if (n > 0.0) |
770 | result /= n; |
771 | else result[0] = T(1.0); |
772 | |
773 | assert(is_normalized(result)); |
774 | |
775 | return result; |
776 | } |
777 | |
778 | template <typename T, size_t N> |
779 | inline Vector<T, N> improve_normalization(const Vector<T, N>& v) |
780 | { |
781 | return improve_normalization<1, T, N>(v); |
782 | } |
783 | |
784 | template <size_t Steps, typename T, size_t N> |
785 | inline Vector<T, N> improve_normalization(const Vector<T, N>& v) |
786 | { |
787 | Vector<T, N> result = v; |
788 | |
789 | for (size_t i = 0; i < Steps; ++i) |
790 | result *= (T(3.0) - square_norm(result)) * T(0.5); |
791 | |
792 | return result; |
793 | } |
794 | |
795 | template <typename T, size_t N> |
796 | inline bool is_normalized(const Vector<T, N>& v) |
797 | { |
798 | return feq(square_norm(v), T(1.0), make_eps<T>(1.0e-4f, 1.0e-5)); |
799 | } |
800 | |
801 | template <typename T, size_t N> |
802 | inline bool is_normalized(const Vector<T, N>& v, const T eps) |
803 | { |
804 | return feq(square_norm(v), T(1.0), eps); |
805 | } |
806 | |
807 | template <typename T, size_t N> |
808 | inline Vector<T, N> faceforward( |
809 | const Vector<T, N>& n, |
810 | const Vector<T, N>& i) |
811 | { |
812 | return dot(n, i) < T(0.0) ? n : -n; |
813 | } |
814 | |
815 | template <typename T, size_t N> |
816 | inline Vector<T, N> flip_to_same_hemisphere( |
817 | const Vector<T, N>& v, |
818 | const Vector<T, N>& ref) |
819 | { |
820 | return dot(v, ref) < T(0.0) ? -v : v; |
821 | } |
822 | |
823 | template <typename T, size_t N> |
824 | inline Vector<T, N> project( |
825 | const Vector<T, N>& v, |
826 | const Vector<T, N>& n) |
827 | { |
828 | return v - dot(v, n) * n; |
829 | } |
830 | |
831 | template <typename T, size_t N> |
832 | inline Vector<T, N> reflect( |
833 | const Vector<T, N>& i, |
834 | const Vector<T, N>& n) |
835 | { |
836 | assert(is_normalized(n)); |
837 | |
838 | const T dot_in = dot(i, n); |
839 | return (dot_in + dot_in) * n - i; |
840 | } |
841 | |
842 | template <typename T, size_t N> |
843 | inline bool refract( |
844 | const Vector<T, N>& i, |
845 | const Vector<T, N>& n, |
846 | const T eta, |
847 | Vector<T, N>& t) |
848 | { |
849 | // |
850 | // Reference: http://en.wikipedia.org/wiki/Snell's_law |
851 | // |
852 | |
853 | assert(is_normalized(i)); |
854 | assert(is_normalized(n)); |
855 | |
856 | const T cos_theta_i = dot(i, n); |
857 | assert(cos_theta_i >= T(0.0)); |
858 | |
859 | const T sin_theta_i2 = T(1.0) - square(cos_theta_i); |
860 | const T sin_theta_t2 = sin_theta_i2 * square(eta); |
861 | const T cos_theta_t2 = T(1.0) - sin_theta_t2; |
862 | |
863 | if (cos_theta_t2 < T(0.0)) |
864 | { |
865 | // Total internal reflection. |
866 | return false; |
867 | } |
868 | |
869 | const T cos_theta_t = std::sqrt(cos_theta_t2); |
870 | |
871 | t = (eta * cos_theta_i - cos_theta_t) * n - eta * i; |
872 | assert(is_normalized(t)); |
873 | |
874 | return true; |
875 | } |
876 | |
877 | template <typename T, size_t N> |
878 | inline Vector<T, N> clamp( |
879 | const Vector<T, N>& v, |
880 | const T min, |
881 | const T max) |
882 | { |
883 | Vector<T, N> result; |
884 | |
885 | for (size_t i = 0; i < N; ++i) |
886 | result[i] = clamp(v[i], min, max); |
887 | |
888 | return result; |
889 | } |
890 | |
891 | template <typename T, size_t N> |
892 | inline Vector<T, N> saturate(const Vector<T, N>& v) |
893 | { |
894 | Vector<T, N> result; |
895 | |
896 | for (size_t i = 0; i < N; ++i) |
897 | result[i] = saturate(v[i]); |
898 | |
899 | return result; |
900 | } |
901 | |
902 | template <typename T, size_t N> |
903 | inline T min_value(const Vector<T, N>& v) |
904 | { |
905 | T value = v[0]; |
906 | |
907 | for (size_t i = 1; i < N; ++i) |
908 | { |
909 | if (value > v[i]) |
910 | value = v[i]; |
911 | } |
912 | |
913 | return value; |
914 | } |
915 | |
916 | template <typename T, size_t N> |
917 | inline T max_value(const Vector<T, N>& v) |
918 | { |
919 | T value = v[0]; |
920 | |
921 | for (size_t i = 1; i < N; ++i) |
922 | { |
923 | if (value < v[i]) |
924 | value = v[i]; |
925 | } |
926 | |
927 | return value; |
928 | } |
929 | |
930 | template <typename T, size_t N> |
931 | inline size_t min_index(const Vector<T, N>& v) |
932 | { |
933 | size_t index = 0; |
934 | T value = v[0]; |
935 | |
936 | for (size_t i = 1; i < N; ++i) |
937 | { |
938 | const T x = v[i]; |
939 | |
940 | if (value > x) |
941 | { |
942 | value = x; |
943 | index = i; |
944 | } |
945 | } |
946 | |
947 | return index; |
948 | } |
949 | |
950 | template <typename T, size_t N> |
951 | inline size_t max_index(const Vector<T, N>& v) |
952 | { |
953 | size_t index = 0; |
954 | T value = v[0]; |
955 | |
956 | for (size_t i = 1; i < N; ++i) |
957 | { |
958 | const T x = v[i]; |
959 | |
960 | if (value < x) |
961 | { |
962 | value = x; |
963 | index = i; |
964 | } |
965 | } |
966 | |
967 | return index; |
968 | } |
969 | |
970 | template <typename T, size_t N> |
971 | inline size_t min_abs_index(const Vector<T, N>& v) |
972 | { |
973 | size_t index = 0; |
974 | T value = std::abs(v[0]); |
975 | |
976 | for (size_t i = 1; i < N; ++i) |
977 | { |
978 | const T x = std::abs(v[i]); |
979 | |
980 | if (value > x) |
981 | { |
982 | value = x; |
983 | index = i; |
984 | } |
985 | } |
986 | |
987 | return index; |
988 | } |
989 | |
990 | template <typename T, size_t N> |
991 | inline size_t max_abs_index(const Vector<T, N>& v) |
992 | { |
993 | size_t index = 0; |
994 | T value = std::abs(v[0]); |
995 | |
996 | for (size_t i = 1; i < N; ++i) |
997 | { |
998 | const T x = std::abs(v[i]); |
999 | |
1000 | if (value < x) |
1001 | { |
1002 | value = x; |
1003 | index = i; |
1004 | } |
1005 | } |
1006 | |
1007 | return index; |
1008 | } |
1009 | |
1010 | template <typename T, size_t N> |
1011 | inline Vector<T, N> component_wise_min(const Vector<T, N>& lhs, const Vector<T, N>& rhs) |
1012 | { |
1013 | Vector<T, N> result; |
1014 | |
1015 | for (size_t i = 0; i < N; ++i) |
1016 | result[i] = std::min(lhs[i], rhs[i]); |
1017 | |
1018 | return result; |
1019 | } |
1020 | |
1021 | template <typename T, size_t N> |
1022 | inline Vector<T, N> component_wise_max(const Vector<T, N>& lhs, const Vector<T, N>& rhs) |
1023 | { |
1024 | Vector<T, N> result; |
1025 | |
1026 | for (size_t i = 0; i < N; ++i) |
1027 | result[i] = std::max(lhs[i], rhs[i]); |
1028 | |
1029 | return result; |
1030 | } |
1031 | |
1032 | |
1033 | // |
1034 | // 2D vector implementation. |
1035 | // |
1036 | |
1037 | template <typename T> |
1038 | inline Vector<T, 2>::Vector() |
1039 | { |
1040 | } |
1041 | |
1042 | template <typename T> |
1043 | inline Vector<T, 2>::Vector(const ValueType* rhs) |
1044 | { |
1045 | assert(rhs); |
1046 | x = rhs[0]; |
1047 | y = rhs[1]; |
1048 | } |
1049 | |
1050 | template <typename T> |
1051 | inline Vector<T, 2>::Vector(const ValueType val) |
1052 | : x(val) |
1053 | , y(val) |
1054 | { |
1055 | } |
1056 | |
1057 | template <typename T> |
1058 | inline Vector<T, 2>::Vector( |
1059 | const ValueType x_, |
1060 | const ValueType y_) |
1061 | : x(x_) |
1062 | , y(y_) |
1063 | { |
1064 | } |
1065 | |
1066 | template <typename T> |
1067 | template <typename U> |
1068 | inline Vector<T, 2>::Vector(const Vector<U, 2>& rhs) |
1069 | : x(static_cast<ValueType>(rhs.x)) |
1070 | , y(static_cast<ValueType>(rhs.y)) |
1071 | { |
1072 | } |
1073 | |
1074 | #ifdef APPLESEED_ENABLE_IMATH_INTEROP |
1075 | |
1076 | template <typename T> |
1077 | inline Vector<T, 2>::Vector(const Imath::Vec2<T>& rhs) |
1078 | : x(rhs.x) |
1079 | , y(rhs.y) |
1080 | { |
1081 | } |
1082 | |
1083 | template <typename T> |
1084 | inline Vector<T, 2>::operator Imath::Vec2<T>&() |
1085 | { |
1086 | return reinterpret_cast<Imath::Vec2<T>&>(*this); |
1087 | } |
1088 | |
1089 | template <typename T> |
1090 | inline Vector<T, 2>::operator const Imath::Vec2<T>&() const |
1091 | { |
1092 | return reinterpret_cast<const Imath::Vec2<T>&>(*this); |
1093 | } |
1094 | |
1095 | #endif |
1096 | |
1097 | template <typename T> |
1098 | inline T& Vector<T, 2>::operator[](const size_t i) |
1099 | { |
1100 | assert(i < Dimension); |
1101 | return (&x)[i]; |
1102 | } |
1103 | |
1104 | template <typename T> |
1105 | inline const T& Vector<T, 2>::operator[](const size_t i) const |
1106 | { |
1107 | assert(i < Dimension); |
1108 | return (&x)[i]; |
1109 | } |
1110 | |
1111 | template <typename T> |
1112 | inline T det(const Vector<T, 2>& lhs, const Vector<T, 2>& rhs) |
1113 | { |
1114 | return lhs.x * rhs.y - lhs.y * rhs.x; |
1115 | } |
1116 | |
1117 | |
1118 | // |
1119 | // 3D vector implementation. |
1120 | // |
1121 | |
1122 | template <typename T> |
1123 | inline Vector<T, 3>::Vector() |
1124 | { |
1125 | } |
1126 | |
1127 | template <typename T> |
1128 | inline Vector<T, 3>::Vector(const ValueType* rhs) |
1129 | { |
1130 | assert(rhs); |
1131 | x = rhs[0]; |
1132 | y = rhs[1]; |
1133 | z = rhs[2]; |
1134 | } |
1135 | |
1136 | template <typename T> |
1137 | inline Vector<T, 3>::Vector(const ValueType val) |
1138 | : x(val) |
1139 | , y(val) |
1140 | , z(val) |
1141 | { |
1142 | } |
1143 | |
1144 | template <typename T> |
1145 | inline Vector<T, 3>::Vector( |
1146 | const ValueType x_, |
1147 | const ValueType y_, |
1148 | const ValueType z_) |
1149 | : x(x_) |
1150 | , y(y_) |
1151 | , z(z_) |
1152 | { |
1153 | } |
1154 | |
1155 | template <typename T> |
1156 | template <typename U> |
1157 | inline Vector<T, 3>::Vector(const Vector<U, 3>& rhs) |
1158 | : x(static_cast<ValueType>(rhs.x)) |
1159 | , y(static_cast<ValueType>(rhs.y)) |
1160 | , z(static_cast<ValueType>(rhs.z)) |
1161 | { |
1162 | } |
1163 | |
1164 | #ifdef APPLESEED_ENABLE_IMATH_INTEROP |
1165 | |
1166 | template <typename T> |
1167 | inline Vector<T, 3>::Vector(const Imath::Vec3<T>& rhs) |
1168 | : x(rhs.x) |
1169 | , y(rhs.y) |
1170 | , z(rhs.z) |
1171 | { |
1172 | } |
1173 | |
1174 | template <typename T> |
1175 | inline Vector<T, 3>::operator Imath::Vec3<T>&() |
1176 | { |
1177 | return reinterpret_cast<Imath::Vec3<T>&>(*this); |
1178 | } |
1179 | |
1180 | template <typename T> |
1181 | inline Vector<T, 3>::operator const Imath::Vec3<T>&() const |
1182 | { |
1183 | return reinterpret_cast<const Imath::Vec3<T>&>(*this); |
1184 | } |
1185 | |
1186 | #endif |
1187 | |
1188 | template <typename T> |
1189 | inline Vector<T, 3> Vector<T, 3>::make_unit_vector( |
1190 | const ValueType theta, |
1191 | const ValueType phi) |
1192 | { |
1193 | const ValueType sin_theta = std::sin(theta); |
1194 | |
1195 | return |
1196 | Vector<T, 3>( |
1197 | std::cos(phi) * sin_theta, |
1198 | std::cos(theta), |
1199 | std::sin(phi) * sin_theta); |
1200 | } |
1201 | |
1202 | template <typename T> |
1203 | inline Vector<T, 3> Vector<T, 3>::make_unit_vector( |
1204 | const ValueType cos_theta, |
1205 | const ValueType sin_theta, |
1206 | const ValueType cos_phi, |
1207 | const ValueType sin_phi) |
1208 | { |
1209 | return |
1210 | Vector<T, 3>( |
1211 | cos_phi * sin_theta, |
1212 | cos_theta, |
1213 | sin_phi * sin_theta); |
1214 | } |
1215 | |
1216 | template <typename T> |
1217 | inline T& Vector<T, 3>::operator[](const size_t i) |
1218 | { |
1219 | assert(i < Dimension); |
1220 | return (&x)[i]; |
1221 | } |
1222 | |
1223 | template <typename T> |
1224 | inline const T& Vector<T, 3>::operator[](const size_t i) const |
1225 | { |
1226 | assert(i < Dimension); |
1227 | return (&x)[i]; |
1228 | } |
1229 | |
1230 | template <typename T> |
1231 | inline Vector<T, 3> cross(const Vector<T, 3>& lhs, const Vector<T, 3>& rhs) |
1232 | { |
1233 | Vector<T, 3> result; |
1234 | result.x = lhs.y * rhs.z - rhs.y * lhs.z; |
1235 | result.y = lhs.z * rhs.x - rhs.z * lhs.x; |
1236 | result.z = lhs.x * rhs.y - rhs.x * lhs.y; |
1237 | return result; |
1238 | } |
1239 | |
1240 | |
1241 | // |
1242 | // 4D vector implementation. |
1243 | // |
1244 | |
1245 | template <typename T> |
1246 | inline Vector<T, 4>::Vector() |
1247 | { |
1248 | } |
1249 | |
1250 | template <typename T> |
1251 | inline Vector<T, 4>::Vector(const ValueType* rhs) |
1252 | { |
1253 | assert(rhs); |
1254 | x = rhs[0]; |
1255 | y = rhs[1]; |
1256 | z = rhs[2]; |
1257 | w = rhs[3]; |
1258 | } |
1259 | |
1260 | template <typename T> |
1261 | inline Vector<T, 4>::Vector(const ValueType val) |
1262 | : x(val) |
1263 | , y(val) |
1264 | , z(val) |
1265 | , w(val) |
1266 | { |
1267 | } |
1268 | |
1269 | template <typename T> |
1270 | inline Vector<T, 4>::Vector( |
1271 | const ValueType x_, |
1272 | const ValueType y_, |
1273 | const ValueType z_, |
1274 | const ValueType w_) |
1275 | : x(x_) |
1276 | , y(y_) |
1277 | , z(z_) |
1278 | , w(w_) |
1279 | { |
1280 | } |
1281 | |
1282 | template <typename T> |
1283 | template <typename U> |
1284 | inline Vector<T, 4>::Vector(const Vector<U, 4>& rhs) |
1285 | : x(static_cast<ValueType>(rhs.x)) |
1286 | , y(static_cast<ValueType>(rhs.y)) |
1287 | , z(static_cast<ValueType>(rhs.z)) |
1288 | , w(static_cast<ValueType>(rhs.w)) |
1289 | { |
1290 | } |
1291 | |
1292 | template <typename T> |
1293 | inline T& Vector<T, 4>::operator[](const size_t i) |
1294 | { |
1295 | assert(i < Dimension); |
1296 | return (&x)[i]; |
1297 | } |
1298 | |
1299 | template <typename T> |
1300 | inline const T& Vector<T, 4>::operator[](const size_t i) const |
1301 | { |
1302 | assert(i < Dimension); |
1303 | return (&x)[i]; |
1304 | } |
1305 | |
1306 | } // namespace foundation |
1307 | |
1308 | #endif // !APPLESEED_FOUNDATION_MATH_VECTOR_H |
1309 | |