1 /**************************************************************************
3 * Copyright 2013-2014 RAD Game Tools and Valve Software
4 * Copyright 2010-2014 Rich Geldreich and Tenacious Software LLC
7 * Permission is hereby granted, free of charge, to any person obtaining a copy
8 * of this software and associated documentation files (the "Software"), to deal
9 * in the Software without restriction, including without limitation the rights
10 * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
11 * copies of the Software, and to permit persons to whom the Software is
12 * furnished to do so, subject to the following conditions:
14 * The above copyright notice and this permission notice shall be included in
15 * all copies or substantial portions of the Software.
17 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
18 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
19 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
20 * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
21 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
22 * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
25 **************************************************************************/
30 #include "vogl_core.h"
31 #include "vogl_matrix.h"
47 inline plane(const vec3F &norm, float dist)
48 : m_dist(dist), m_normal(norm)
52 inline plane(const vec4F &norm, float dist)
53 : m_dist(dist), m_normal(norm)
57 inline plane(float nx, float ny, float nz, float dist)
58 : m_dist(dist), m_normal(nx, ny, nz)
62 inline plane(const vec3F &norm, const vec3F &org)
63 : m_dist(norm.dot(org)), m_normal(norm)
67 inline plane(const vec3F &p, const vec3F &q, const vec3F &r)
72 inline plane(const vec4F &equation)
73 : m_dist(-equation[3]), m_normal(equation)
77 inline plane(const oriented_plane &op)
82 inline bool is_valid() const
84 return (m_normal[0] != 0.0f) || (m_normal[1] != 0.0f) || (m_normal[2] != 0.0f);
94 inline const vec3F &get_normal() const
99 inline vec3F &get_normal()
104 inline const float &get_distance() const
109 inline float &get_distance()
114 inline vec4F get_equation() const
116 return vec4F(m_normal[0], m_normal[1], m_normal[2], -m_dist);
119 inline vec3F get_point_on_plane() const
121 return m_normal * m_dist;
124 inline const plane &init(float nx, float ny, float nz, float dist)
126 m_normal.set(nx, ny, nz);
131 inline const plane &init(const vec4F &eq)
138 inline bool init(const vec3F &p, const vec3F &q, const vec3F &r, float epsilon = 1e-6f)
140 m_normal = vec3F::cross3((r - q), (p - q));
142 double normal_len2 = m_normal.norm();
143 if (normal_len2 <= epsilon)
149 m_normal *= static_cast<float>(1.0f / sqrt(normal_len2));
150 m_dist = m_normal.dot(q);
154 inline void init_denormalized(const vec3F &p, const vec3F &q, const vec3F &r)
156 m_normal = vec3F::cross3(r - q, p - q);
157 m_dist = m_normal.dot(q);
160 inline plane &init(const vec3F &norm, const vec3F &point_on_plane)
163 m_dist = norm.dot(point_on_plane);
167 inline plane &init(const oriented_plane &op);
169 inline plane &normalize_in_place()
171 double normal_len = m_normal.length();
172 if (normal_len == 0.0f)
176 double one_over_normal_len = 1.0f / normal_len;
177 m_normal *= static_cast<float>(one_over_normal_len);
178 m_dist *= static_cast<float>(one_over_normal_len);
183 inline bool operator==(const plane &rhs) const
185 return (m_normal == rhs.m_normal) && (m_dist == rhs.m_dist);
188 inline bool operator!=(const plane &rhs) const
190 return !(*this == rhs);
193 inline float get_distance(const vec3F &p) const
195 return p.dot(m_normal) - m_dist;
198 inline plane get_flipped() const
200 return plane(-m_normal, -m_dist);
203 // computes how much solve_axis changes with movement along exclusively movement_axis
204 inline float get_gradient(uint solve_axis, uint movement_axis, float epsilon = 1e-6f) const
206 const float dA = -m_normal[movement_axis];
207 const float dB = m_normal[solve_axis];
208 if (fabs(dB) <= epsilon)
216 inline plane get_transformed(const matrix44F &m)
218 matrix44F m_inverted;
219 if (!m.invert(m_inverted))
225 result.m_normal = matrix44F::transform_vector_transposed(m_normal, m_inverted).normalize_in_place();
226 result.m_dist = matrix44F::transform_point(get_point_on_plane(), m).dot3(result.m_normal);
230 inline plane get_transformed_ortho(const matrix44F &m)
233 result.m_normal = matrix44F::transform_vector(m_normal, m);
234 result.m_dist = m_dist + vec3F(m[3]).dot(result.m_normal);
238 inline vec3F project_to_plane(const vec3F &p) const
240 return p - m_normal * get_distance(p);
243 // computes z given xy
244 inline float solve_for_z(float x, float y) const
246 return m_normal[2] ? ((m_dist - x * m_normal[0] - y * m_normal[1]) / m_normal[2]) : 0.0f;
249 // generalized for any axis
250 inline float solve_for_axis(const vec3F &p, uint axis) const
252 VOGL_ASSERT(axis < 3);
253 static const uint8 axis2[] = { 1, 0, 0 };
254 static const uint8 axis3[] = { 2, 2, 1 };
255 return m_normal[axis] ? ((m_dist - p[axis2[axis]] * m_normal[axis2[axis]] - p[axis3[axis]] * m_normal[axis3[axis]]) / m_normal[axis]) : 0.0f;
265 inline oriented_plane()
269 inline oriented_plane(const plane &p)
274 inline oriented_plane(const vec3F &origon, const vec3F &u_axis, const vec3F &v_axis)
275 : m_u(u_axis), m_v(v_axis), m_origin(origon)
279 inline oriented_plane(const vec3F &origon, const vec3F &norm)
291 inline const vec3F &get_origin() const
295 inline vec3F &get_origin()
300 inline const vec3F &get_u() const
304 inline vec3F &get_u()
309 inline const vec3F &get_v() const
313 inline vec3F &get_v()
318 inline oriented_plane &init(const vec3F &origon, const vec3F &u_axis, const vec3F &v_axis)
326 inline oriented_plane &init(const plane &p)
328 m_origin = p.get_point_on_plane();
329 m_u = vec3F::make_axis(p.m_normal.get_minor_axis()).remove_unit_direction(p.m_normal).normalize_in_place();
330 m_v = vec3F::cross3(p.m_normal, m_u).normalize_in_place();
334 inline oriented_plane &init(const vec3F &origon, const vec3F &norm)
337 m_u = vec3F::make_axis(norm.get_minor_axis()).remove_unit_direction(norm).normalize_in_place();
338 m_v = vec3F::cross3(norm, m_u).normalize_in_place();
342 inline vec2F get_projected(const vec3F &p) const
344 vec3F rel(p - m_origin);
345 return vec2F(m_u.dot(rel), m_v.dot(rel));
348 inline vec3F get_point(const vec2F &p) const
350 return m_u * p[0] + m_v * p[1] + m_origin;
353 inline vec3F get_normal() const
355 return vec3F::cross3(m_u, m_v).normalize_in_place();
358 inline oriented_plane &normalize_in_place()
365 inline oriented_plane &init_planar_projection(const plane &p)
367 const uint major_axis = p.get_normal().get_major_axis();
368 uint axis0 = math::next_wrap<uint>(major_axis, 3U);
369 uint axis1 = math::next_wrap<uint>(axis0, 3U);
370 if (p.get_normal()[major_axis] > 0.0f)
371 std::swap(axis0, axis1);
372 return init(vec3F(0), vec3F::make_axis(axis0), vec3F::make_axis(axis1));
376 inline plane &plane::init(const oriented_plane &op)
378 m_normal = op.get_normal();
379 m_dist = op.get_origin().dot(op.get_normal());
383 inline void clip_convex_polygon_against_plane(vogl::vector<vec3F> &result, const vogl::vector<vec3F> &src_poly, const plane &p)
386 if (!src_poly.size())
389 uint prev_index = src_poly.size() - 1;
390 float prev_dist = p.get_distance(src_poly[prev_index]);
391 for (uint cur_index = 0; cur_index < src_poly.size(); ++cur_index)
393 if (prev_dist >= 0.0f)
394 result.push_back(src_poly[prev_index]);
396 float cur_dist = p.get_distance(src_poly[cur_index]);
397 if (((prev_dist < 0.0f) && (cur_dist > 0.0f)) || ((prev_dist > 0.0f) && (cur_dist < 0.0f)))
398 result.push_back(vec3F::lerp(src_poly[prev_index], src_poly[cur_index], prev_dist / (prev_dist - cur_dist)));
400 prev_index = cur_index;
401 prev_dist = cur_dist;
405 // This stuff was derived from the plane equation for a 3D triangle, with XY in screenspace and Z the component to interpolate:
406 //compute denormalized plane normal/distance:
407 //m_normal = vec3F::cross3(r - q, p - q);
408 //m_dist = m_normal * q;
409 //return vec<3, T>(a[1] * b[2] - a[2] * b[1], a[2] * b[0] - a[0] * b[2], a[0] * b[1] - a[1] * b[0]);
410 struct tri_gradient_common
413 double m_ax, m_ay, m_bx, m_by;
416 double m_one_over_norm_z;
420 utils::zero_object(*this);
423 // q is treated as the origin
424 // vector_type must be at least 2 dimension, xy can be in screenspace (for example), component to be interpolated is treated as z axis.
425 template <typename vector_type>
426 inline void init(const vector_type &p, const vector_type &q, const vector_type &r)
437 m_norm_z = m_ax * m_by - m_ay * m_bx;
438 m_one_over_norm_z = m_norm_z ? (1.0 / m_norm_z) : 0.0;
442 template <typename scalar_type>
443 struct tri_gradient_component
446 scalar_type m_norm_x, m_norm_y;
451 utils::zero_object(*this);
454 void init(const tri_gradient_common &common, float pz, float qz, float rz)
461 m_norm_x = static_cast<scalar_type>(common.m_ay * bz - az * common.m_by);
462 m_norm_y = static_cast<scalar_type>(az * common.m_bx - common.m_ax * bz);
464 m_dist = static_cast<scalar_type>(m_norm_x * common.m_qx + m_norm_y * common.m_qy + common.m_norm_z * m_qz);
467 // returns denormalized plane
468 inline plane get_plane(const tri_gradient_common &common) const
470 return plane(static_cast<float>(m_norm_x), static_cast<float>(m_norm_y), static_cast<float>(common.m_norm_z), static_cast<float>(m_dist));
473 // not accurate at float precision
474 inline double solve_for_z(const tri_gradient_common &common, float x, float y) const
476 return (m_dist - x * m_norm_x - y * m_norm_y) * common.m_one_over_norm_z;
479 // more accurate at float precision
480 inline double solve_for_z_alt(const tri_gradient_common &common, float x, float y) const
482 double dz_dx = get_dz_over_dx_gradient(common);
483 double dz_dy = get_dz_over_dy_gradient(common);
484 return m_qz + (x - common.m_qx) * dz_dx + (y - common.m_qy) * dz_dy;
487 inline double get_dz_over_dx_gradient(const tri_gradient_common &common) const
489 return -m_norm_x * common.m_one_over_norm_z;
492 inline double get_dz_over_dy_gradient(const tri_gradient_common &common) const
494 return -m_norm_y * common.m_one_over_norm_z;
508 inline plane2D(const vec2F &norm, float dist)
509 : m_dist(dist), m_normal(norm)
513 inline plane2D(float nx, float ny, float dist)
514 : m_dist(dist), m_normal(nx, ny)
518 inline plane2D(const vec2F &p, const vec2F &q, float epsilon = 1e-6f)
523 inline plane2D(const vec3F &equation)
524 : m_dist(-equation[2]), m_normal(equation)
528 inline bool is_valid() const
530 return (m_normal[0] != 0.0f) || (m_normal[1] != 0.0f);
533 inline plane2D &clear()
540 inline const vec2F &get_normal() const
545 inline vec2F &get_normal()
550 inline const float &get_distance() const
555 inline float &get_distance()
560 inline vec3F get_equation() const
562 return vec3F(m_normal[0], m_normal[1], -m_dist);
565 inline vec2F get_point_on_plane() const
567 return m_normal * m_dist;
570 inline const plane2D &init(float nx, float ny, float dist)
572 m_normal.set(nx, ny);
577 inline const plane2D &init(const vec3F &eq)
584 inline bool init(const vec2F &p, const vec2F &q, float epsilon = 1e-6f)
586 m_normal = (p - q).perp();
588 double normal_len2 = m_normal.norm();
589 if (normal_len2 <= epsilon)
595 m_normal *= static_cast<float>(1.0f / sqrt(normal_len2));
596 m_dist = m_normal.dot(q);
600 inline void init_denormalized(const vec2F &p, const vec2F &q)
602 m_normal = (p - q).perp();
603 m_dist = m_normal.dot(q);
606 inline plane2D &init_from_normal_point(const vec2F &norm, const vec2F &point_on_plane)
609 m_dist = norm.dot(point_on_plane);
613 inline plane2D &normalize_in_place()
615 double normal_len = m_normal.length();
616 if (normal_len == 0.0f)
620 double one_over_normal_len = 1.0f / normal_len;
621 m_normal *= static_cast<float>(one_over_normal_len);
622 m_dist *= static_cast<float>(one_over_normal_len);
627 inline bool operator==(const plane2D &rhs) const
629 return (m_normal == rhs.m_normal) && (m_dist == rhs.m_dist);
632 inline bool operator!=(const plane2D &rhs) const
634 return !(*this == rhs);
637 inline float get_distance(const vec2F &p) const
639 return p.dot(m_normal) - m_dist;
642 inline plane2D get_flipped() const
644 return plane2D(-m_normal, -m_dist);
647 // computes how much solve_axis changes with movement along exclusively movement_axis
648 inline float get_gradient(uint solve_axis, uint movement_axis, float epsilon = 1e-6f) const
650 const float dA = -m_normal[movement_axis];
651 const float dB = m_normal[solve_axis];
652 if (fabs(dB) <= epsilon)
660 inline plane2D get_transformed(const matrix33F &m)
662 matrix33F m_inverted;
663 if (!m.invert(m_inverted))
669 result.m_normal = matrix33F::transform_vector_transposed(m_normal, m_inverted).normalize_in_place();
670 result.m_dist = matrix33F::transform_point(get_point_on_plane(), m).dot2(result.m_normal);
674 inline plane2D get_transformed_ortho(const matrix33F &m)
677 result.m_normal = matrix33F::transform_vector(m_normal, m);
678 result.m_dist = m_dist + vec2F(m[2]).dot(result.m_normal);
682 inline vec2F project_to_plane(const vec2F &p) const
684 return p - m_normal * get_distance(p);
687 inline float solve_for_axis(const vec2F &p, uint axis) const
689 return m_normal[axis] ? ((m_dist - p[axis ^ 1] * m_normal[axis ^ 1]) / m_normal[axis]) : 0.0f;
692 inline float solve_for_x(float y) const
694 return m_normal[0] ? ((m_dist - y * m_normal[1]) / m_normal[0]) : 0.0f;
697 inline float solve_for_y(float x) const
699 return m_normal[1] ? ((m_dist - x * m_normal[0]) / m_normal[1]) : 0.0f;