The GtkRadiant sources as originally released under the GPL license.

This commit is contained in:
Travis Bradshaw
2012-01-31 15:20:35 -06:00
commit 0991a5ce8b
1590 changed files with 431941 additions and 0 deletions

23
libs/math/aabb.cpp Normal file
View File

@@ -0,0 +1,23 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#include "aabb.h"

301
libs/math/aabb.h Normal file
View File

@@ -0,0 +1,301 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#if !defined(INCLUDED_MATH_AABB_H)
#define INCLUDED_MATH_AABB_H
/// \file
/// \brief Axis-aligned bounding-box data types and related operations.
#include "math/matrix.h"
#include "math/plane.h"
class AABB
{
public:
Vector3 origin, extents;
AABB() : origin(0, 0, 0), extents(-1,-1,-1)
{
}
AABB(const Vector3& origin_, const Vector3& extents_) :
origin(origin_), extents(extents_)
{
}
};
const float c_aabb_max = FLT_MAX;
inline bool extents_valid(float f)
{
return f >= 0.0f && f <= c_aabb_max;
}
inline bool origin_valid(float f)
{
return f >= -c_aabb_max && f <= c_aabb_max;
}
inline bool aabb_valid(const AABB& aabb)
{
return origin_valid(aabb.origin[0])
&& origin_valid(aabb.origin[1])
&& origin_valid(aabb.origin[2])
&& extents_valid(aabb.extents[0])
&& extents_valid(aabb.extents[1])
&& extents_valid(aabb.extents[2]);
}
inline AABB aabb_for_minmax(const Vector3& min, const Vector3& max)
{
AABB aabb;
aabb.origin = vector3_mid(min, max);
aabb.extents = vector3_subtracted(max, aabb.origin);
return aabb;
}
template<typename Index>
class AABBExtend
{
public:
static void apply(AABB& aabb, const Vector3& point)
{
float displacement = point[Index::VALUE] - aabb.origin[Index::VALUE];
float half_difference = static_cast<float>(0.5 * (fabs(displacement) - aabb.extents[Index::VALUE]));
if(half_difference > 0.0f)
{
aabb.origin[Index::VALUE] += (displacement >= 0.0f) ? half_difference : -half_difference;
aabb.extents[Index::VALUE] += half_difference;
}
}
static void apply(AABB& aabb, const AABB& other)
{
float displacement = other.origin[Index::VALUE] - aabb.origin[Index::VALUE];
float difference = other.extents[Index::VALUE] - aabb.extents[Index::VALUE];
if(fabs(displacement) > fabs(difference))
{
float half_difference = static_cast<float>(0.5 * (fabs(displacement) + difference));
if(half_difference > 0.0f)
{
aabb.origin[Index::VALUE] += (displacement >= 0.0f) ? half_difference : -half_difference;
aabb.extents[Index::VALUE] += half_difference;
}
}
else if(difference > 0.0f)
{
aabb.origin[Index::VALUE] = other.origin[Index::VALUE];
aabb.extents[Index::VALUE] = other.extents[Index::VALUE];
}
}
};
inline void aabb_extend_by_point(AABB& aabb, const Vector3& point)
{
AABBExtend< IntegralConstant<0> >::apply(aabb, point);
AABBExtend< IntegralConstant<1> >::apply(aabb, point);
AABBExtend< IntegralConstant<2> >::apply(aabb, point);
}
inline void aabb_extend_by_point_safe(AABB& aabb, const Vector3& point)
{
if(aabb_valid(aabb))
{
aabb_extend_by_point(aabb, point);
}
else
{
aabb.origin = point;
aabb.extents = Vector3(0, 0, 0);
}
}
class AABBExtendByPoint
{
AABB& m_aabb;
public:
AABBExtendByPoint(AABB& aabb) : m_aabb(aabb)
{
}
void operator()(const Vector3& point) const
{
aabb_extend_by_point_safe(m_aabb, point);
}
};
inline void aabb_extend_by_aabb(AABB& aabb, const AABB& other)
{
AABBExtend< IntegralConstant<0> >::apply(aabb, other);
AABBExtend< IntegralConstant<1> >::apply(aabb, other);
AABBExtend< IntegralConstant<2> >::apply(aabb, other);
}
inline void aabb_extend_by_aabb_safe(AABB& aabb, const AABB& other)
{
if(aabb_valid(aabb) && aabb_valid(other))
{
aabb_extend_by_aabb(aabb, other);
}
else if(aabb_valid(other))
{
aabb = other;
}
}
inline void aabb_extend_by_vec3(AABB& aabb, const Vector3& extension)
{
vector3_add(aabb.extents, extension);
}
template<typename Index>
inline bool aabb_intersects_point_dimension(const AABB& aabb, const Vector3& point)
{
return fabs(point[Index::VALUE] - aabb.origin[Index::VALUE]) < aabb.extents[Index::VALUE];
}
inline bool aabb_intersects_point(const AABB& aabb, const Vector3& point)
{
return aabb_intersects_point_dimension< IntegralConstant<0> >(aabb, point)
&& aabb_intersects_point_dimension< IntegralConstant<1> >(aabb, point)
&& aabb_intersects_point_dimension< IntegralConstant<2> >(aabb, point);
}
template<typename Index>
inline bool aabb_intersects_aabb_dimension(const AABB& aabb, const AABB& other)
{
return fabs(other.origin[Index::VALUE] - aabb.origin[Index::VALUE]) < (aabb.extents[Index::VALUE] + other.extents[Index::VALUE]);
}
inline bool aabb_intersects_aabb(const AABB& aabb, const AABB& other)
{
return aabb_intersects_aabb_dimension< IntegralConstant<0> >(aabb, other)
&& aabb_intersects_aabb_dimension< IntegralConstant<1> >(aabb, other)
&& aabb_intersects_aabb_dimension< IntegralConstant<2> >(aabb, other);
}
inline unsigned int aabb_classify_plane(const AABB& aabb, const Plane3& plane)
{
double distance_origin = vector3_dot(plane.normal(), aabb.origin) + plane.dist();
if(fabs(distance_origin) < (fabs(plane.a * aabb.extents[0])
+ fabs(plane.b * aabb.extents[1])
+ fabs(plane.c * aabb.extents[2])))
{
return 1; // partially inside
}
else if (distance_origin < 0)
{
return 2; // totally inside
}
return 0; // totally outside
}
inline unsigned int aabb_oriented_classify_plane(const AABB& aabb, const Matrix4& transform, const Plane3& plane)
{
double distance_origin = vector3_dot(plane.normal(), aabb.origin) + plane.dist();
if(fabs(distance_origin) < (fabs(aabb.extents[0] * vector3_dot(plane.normal(), vector4_to_vector3(transform.x())))
+ fabs(aabb.extents[1] * vector3_dot(plane.normal(), vector4_to_vector3(transform.y())))
+ fabs(aabb.extents[2] * vector3_dot(plane.normal(), vector4_to_vector3(transform.z())))))
{
return 1; // partially inside
}
else if (distance_origin < 0)
{
return 2; // totally inside
}
return 0; // totally outside
}
inline void aabb_corners(const AABB& aabb, Vector3 corners[8])
{
Vector3 min(vector3_subtracted(aabb.origin, aabb.extents));
Vector3 max(vector3_added(aabb.origin, aabb.extents));
corners[0] = Vector3(min[0], max[1], max[2]);
corners[1] = Vector3(max[0], max[1], max[2]);
corners[2] = Vector3(max[0], min[1], max[2]);
corners[3] = Vector3(min[0], min[1], max[2]);
corners[4] = Vector3(min[0], max[1], min[2]);
corners[5] = Vector3(max[0], max[1], min[2]);
corners[6] = Vector3(max[0], min[1], min[2]);
corners[7] = Vector3(min[0], min[1], min[2]);
}
inline void aabb_planes(const AABB& aabb, Plane3 planes[6])
{
planes[0] = Plane3(g_vector3_axes[0], aabb.origin[0] + aabb.extents[0]);
planes[1] = Plane3(vector3_negated(g_vector3_axes[0]), -(aabb.origin[0] - aabb.extents[0]));
planes[2] = Plane3(g_vector3_axes[1], aabb.origin[1] + aabb.extents[1]);
planes[3] = Plane3(vector3_negated(g_vector3_axes[1]), -(aabb.origin[1] - aabb.extents[1]));
planes[4] = Plane3(g_vector3_axes[2], aabb.origin[2] + aabb.extents[2]);
planes[5] = Plane3(vector3_negated(g_vector3_axes[2]), -(aabb.origin[2] - aabb.extents[2]));
}
const Vector3 aabb_normals[6] = {
Vector3( 1, 0, 0 ),
Vector3( 0, 1, 0 ),
Vector3( 0, 0, 1 ),
Vector3(-1, 0, 0 ),
Vector3( 0,-1, 0 ),
Vector3( 0, 0,-1 ),
};
const float aabb_texcoord_topleft[2] = { 0, 0 };
const float aabb_texcoord_topright[2] = { 1, 0 };
const float aabb_texcoord_botleft[2] = { 0, 1 };
const float aabb_texcoord_botright[2] = { 1, 1 };
inline AABB aabb_for_oriented_aabb(const AABB& aabb, const Matrix4& transform)
{
return AABB(
matrix4_transformed_point(transform, aabb.origin),
Vector3(
static_cast<float>(fabs(transform[0] * aabb.extents[0])
+ fabs(transform[4] * aabb.extents[1])
+ fabs(transform[8] * aabb.extents[2])),
static_cast<float>(fabs(transform[1] * aabb.extents[0])
+ fabs(transform[5] * aabb.extents[1])
+ fabs(transform[9] * aabb.extents[2])),
static_cast<float>(fabs(transform[2] * aabb.extents[0])
+ fabs(transform[6] * aabb.extents[1])
+ fabs(transform[10] * aabb.extents[2]))
)
);
}
inline AABB aabb_for_oriented_aabb_safe(const AABB& aabb, const Matrix4& transform)
{
if(aabb_valid(aabb))
{
return aabb_for_oriented_aabb(aabb, transform);
}
return aabb;
}
inline AABB aabb_infinite()
{
return AABB(Vector3(0, 0, 0), Vector3(c_aabb_max, c_aabb_max, c_aabb_max));
}
#endif

23
libs/math/curve.cpp Normal file
View File

@@ -0,0 +1,23 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#include "curve.h"

273
libs/math/curve.h Normal file
View File

@@ -0,0 +1,273 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#if !defined(INCLUDED_MATH_CURVE_H)
#define INCLUDED_MATH_CURVE_H
/// \file
/// \brief Curve data types and related operations.
#include "debugging/debugging.h"
#include "container/array.h"
#include <math/matrix.h>
template<typename I, typename Degree>
struct BernsteinPolynomial
{
static double apply(double t)
{
return 1; // general case not implemented
}
};
typedef IntegralConstant<0> Zero;
typedef IntegralConstant<1> One;
typedef IntegralConstant<2> Two;
typedef IntegralConstant<3> Three;
typedef IntegralConstant<4> Four;
template<>
struct BernsteinPolynomial<Zero, Zero>
{
static double apply(double t)
{
return 1;
}
};
template<>
struct BernsteinPolynomial<Zero, One>
{
static double apply(double t)
{
return 1 - t;
}
};
template<>
struct BernsteinPolynomial<One, One>
{
static double apply(double t)
{
return t;
}
};
template<>
struct BernsteinPolynomial<Zero, Two>
{
static double apply(double t)
{
return (1 - t) * (1 - t);
}
};
template<>
struct BernsteinPolynomial<One, Two>
{
static double apply(double t)
{
return 2 * (1 - t) * t;
}
};
template<>
struct BernsteinPolynomial<Two, Two>
{
static double apply(double t)
{
return t * t;
}
};
template<>
struct BernsteinPolynomial<Zero, Three>
{
static double apply(double t)
{
return (1 - t) * (1 - t) * (1 - t);
}
};
template<>
struct BernsteinPolynomial<One, Three>
{
static double apply(double t)
{
return 3 * (1 - t) * (1 - t) * t;
}
};
template<>
struct BernsteinPolynomial<Two, Three>
{
static double apply(double t)
{
return 3 * (1 - t) * t * t;
}
};
template<>
struct BernsteinPolynomial<Three, Three>
{
static double apply(double t)
{
return t * t * t;
}
};
typedef Array<Vector3> ControlPoints;
inline Vector3 CubicBezier_evaluate(const Vector3* firstPoint, double t)
{
Vector3 result(0, 0, 0);
double denominator = 0;
{
double weight = BernsteinPolynomial<Zero, Three>::apply(t);
result += vector3_scaled(*firstPoint++, weight);
denominator += weight;
}
{
double weight = BernsteinPolynomial<One, Three>::apply(t);
result += vector3_scaled(*firstPoint++, weight);
denominator += weight;
}
{
double weight = BernsteinPolynomial<Two, Three>::apply(t);
result += vector3_scaled(*firstPoint++, weight);
denominator += weight;
}
{
double weight = BernsteinPolynomial<Three, Three>::apply(t);
result += vector3_scaled(*firstPoint++, weight);
denominator += weight;
}
return result / denominator;
}
inline Vector3 CubicBezier_evaluateMid(const Vector3* firstPoint)
{
return vector3_scaled(firstPoint[0], 0.125)
+ vector3_scaled(firstPoint[1], 0.375)
+ vector3_scaled(firstPoint[2], 0.375)
+ vector3_scaled(firstPoint[3], 0.125);
}
inline Vector3 CatmullRom_evaluate(const ControlPoints& controlPoints, double t)
{
// scale t to be segment-relative
t *= double(controlPoints.size() - 1);
// subtract segment index;
std::size_t segment = 0;
for(std::size_t i = 0; i < controlPoints.size() - 1; ++i)
{
if(t <= double(i+1))
{
segment = i;
break;
}
}
t -= segment;
const double reciprocal_alpha = 1.0 / 3.0;
Vector3 bezierPoints[4];
bezierPoints[0] = controlPoints[segment];
bezierPoints[1] = (segment > 0)
? controlPoints[segment] + vector3_scaled(controlPoints[segment + 1] - controlPoints[segment - 1], reciprocal_alpha * 0.5)
: controlPoints[segment] + vector3_scaled(controlPoints[segment + 1] - controlPoints[segment], reciprocal_alpha);
bezierPoints[2] = (segment < controlPoints.size() - 2)
? controlPoints[segment + 1] + vector3_scaled(controlPoints[segment] - controlPoints[segment + 2], reciprocal_alpha * 0.5)
: controlPoints[segment + 1] + vector3_scaled(controlPoints[segment] - controlPoints[segment + 1], reciprocal_alpha);
bezierPoints[3] = controlPoints[segment + 1];
return CubicBezier_evaluate(bezierPoints, t);
}
typedef Array<float> Knots;
inline double BSpline_basis(const Knots& knots, std::size_t i, std::size_t degree, double t)
{
if(degree == 0)
{
if(knots[i] <= t
&& t < knots[i + 1]
&& knots[i] < knots[i + 1])
{
return 1;
}
return 0;
}
double leftDenom = knots[i + degree] - knots[i];
double left = (leftDenom == 0) ? 0 : ((t - knots[i]) / leftDenom) * BSpline_basis(knots, i, degree - 1, t);
double rightDenom = knots[i + degree + 1] - knots[i + 1];
double right = (rightDenom == 0) ? 0 : ((knots[i + degree + 1] - t) / rightDenom) * BSpline_basis(knots, i + 1, degree - 1, t);
return left + right;
}
inline Vector3 BSpline_evaluate(const ControlPoints& controlPoints, const Knots& knots, std::size_t degree, double t)
{
Vector3 result(0, 0, 0);
for(std::size_t i = 0; i < controlPoints.size(); ++i)
{
result += vector3_scaled(controlPoints[i], BSpline_basis(knots, i, degree, t));
}
return result;
}
typedef Array<float> NURBSWeights;
inline Vector3 NURBS_evaluate(const ControlPoints& controlPoints, const NURBSWeights& weights, const Knots& knots, std::size_t degree, double t)
{
Vector3 result(0, 0, 0);
double denominator = 0;
for(std::size_t i = 0; i < controlPoints.size(); ++i)
{
double weight = weights[i] * BSpline_basis(knots, i, degree, t);
result += vector3_scaled(controlPoints[i], weight);
denominator += weight;
}
return result / denominator;
}
inline void KnotVector_openUniform(Knots& knots, std::size_t count, std::size_t degree)
{
knots.resize(count + degree + 1);
std::size_t equalKnots = 1;
for(std::size_t i = 0; i < equalKnots; ++i)
{
knots[i] = 0;
knots[knots.size() - (i + 1)] = 1;
}
std::size_t difference = knots.size() - 2 * (equalKnots);
for(std::size_t i = 0; i < difference; ++i)
{
knots[i + equalKnots] = Knots::value_type(double(i + 1) * 1.0 / double(difference + 1));
}
}
#endif

243
libs/math/expression.cpp Normal file
View File

@@ -0,0 +1,243 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#include "expression.h"
Vector3 testAdded1(const Vector3& a, const Vector3& b)
{
return vector3_added(a, vector3_added(a, b));
}
Vector3 testAdded2(const Vector3& a, const Vector3& b)
{
return vector3_for_expression( vector_added( vector3_identity(a), vector_added( vector3_identity(a), vector3_identity(b) ) ) );
}
Vector3 testMultiplied1(const Vector3& a, const Vector3& b)
{
return vector3_scaled(a, b);
}
Vector3 testMultiplied2(const Vector3& a, const Vector3& b)
{
return vector3_for_expression( vector_multiplied( vector3_identity(a), vector3_identity(b) ) );
}
Vector3 testCross1(const Vector3& a, const Vector3& b)
{
return vector3_cross(a, b);
}
Vector3 testCross2(const Vector3& a, const Vector3& b)
{
return vector3_for_expression( vector_cross( vector3_identity(a), vector3_identity(b) ) );
}
double testDot1(const Vector3& a, const Vector3& b)
{
return vector3_dot(a, b);
}
double testDot2(const Vector3& a, const Vector3& b)
{
return float_for_expression( vector_dot( vector3_identity(a), vector3_identity(b) ) );
}
double testLength1(const Vector3& a)
{
return vector3_length(a);
}
double testLength2(const Vector3& a)
{
return float_for_expression( vector_length( vector3_identity(a) ) );
}
Vector3 testNormalised1(const Vector3& a)
{
return vector3_normalised(a);
}
Vector3 testNormalised2(const Vector3& a)
{
return vector3_for_expression( vector_normalised( vector3_identity(a) ) );
}
Vector3 testNegated1(const Vector3& a)
{
return vector3_negated(a);
}
Vector3 testNegated2(const Vector3& a)
{
return vector3_for_expression( vector_negated( vector3_identity(a) ) );
}
Vector3 testScaled1(const Vector3& a, const double& b)
{
return vector3_scaled(a, b);
}
Vector3 testScaled2(const Vector3& a, const double& b)
{
return vector3_for_expression( vector_scaled( vector3_identity(a), float_literal(b) ) );
}
Vector3 testMatrixMultiplied1(const Vector3& a, const Matrix4& b)
{
return matrix4_transformed_point(b, vector3_added(a, Vector3(1, 0, 0)));
}
Vector3 testMatrixMultiplied2(const Vector3& a, const Matrix4& b)
{
return vector3_for_expression(
point_multiplied(
vector_added(
vector3_identity(a),
vector3_literal(Vector3(1, 0, 0))
),
matrix4_identity(b)
)
);
}
Matrix4 testMatrix4Multiplied1(const Matrix4& a, const Matrix4& b)
{
return matrix4_multiplied_by_matrix4(a, matrix4_multiplied_by_matrix4(a, b));
}
Matrix4 testMatrix4Multiplied2(const Matrix4& a, const Matrix4& b)
{
return matrix4_for_expression(
matrix4_multiplied(
matrix4_identity(a),
matrix4_identity(b)
)
);
}
Matrix4 testMatrix4AffineMultiplied1(const Matrix4& a, const Matrix4& b)
{
return matrix4_affine_multiplied_by_matrix4(a, b);
}
Matrix4 testMatrix4AffineMultiplied2(const Matrix4& a, const Matrix4& b)
{
return matrix4_affine_for_expression(
matrix4_multiplied(
matrix4_identity(a),
matrix4_identity(b)
)
);
}
Matrix4 testMatrix4MultipliedConstant1(const Matrix4& a)
{
return matrix4_multiplied_by_matrix4(a, g_matrix4_identity);
}
Matrix4 testMatrix4MultipliedConstant2(const Matrix4& a)
{
return matrix4_for_expression(
matrix4_multiplied(
matrix4_identity(a),
matrix4_identity(g_matrix4_identity)
)
);
}
Matrix4 testMatrix4Transposed1(const Matrix4& a)
{
return matrix4_transposed(a);
}
Matrix4 testMatrix4Transposed2(const Matrix4& a)
{
return matrix4_for_expression( matrix_transposed( matrix4_identity(a) ) );
}
Vector3 testMulti1(const Matrix4& a, const Vector3& b, const Vector3& c)
{
return vector3_added(matrix4_transformed_point(matrix4_transposed(a), b), c);
}
Vector3 testMulti2(const Matrix4& a, const Vector3& b, const Vector3& c)
{
return vector3_for_expression(
vector_added(
point_multiplied(
vector3_identity(b),
matrix_transposed(matrix4_identity(a))
),
vector3_identity(c)
)
);
}
template<typename Value, typename First, typename Second>
class TestBinaryFunction
{
typedef Value(*Function)(const First&, const Second&);
Function m_function;
public:
TestBinaryFunction(Function function) : m_function(function)
{
}
Value run(const First& first, const Second& second) const
{
return m_function(first, second);
}
};
template<typename Value, typename First>
class TestUnaryFunction
{
typedef Value(*Function)(const First&);
Function m_function;
public:
TestUnaryFunction(Function function) : m_function(function)
{
}
Value run(const First& first) const
{
return m_function(first);
}
};
class TestAll
{
public:
TestAll()
{
Vector3 result1 = TestBinaryFunction<Vector3, Vector3, Vector3>(testAdded1).run(Vector3(0, 0, 0), Vector3(1, 1, 1));
Vector3 result2 = TestBinaryFunction<Vector3, Vector3, Vector3>(testAdded2).run(Vector3(0, 0, 0), Vector3(1, 1, 1));
Vector3 result3 = TestBinaryFunction<Vector3, Vector3, Vector3>(testMultiplied1).run(Vector3(1, 2, 3), Vector3(2, 1, 0.5f));
Vector3 result4 = TestBinaryFunction<Vector3, Vector3, Vector3>(testMultiplied2).run(Vector3(1, 2, 3), Vector3(2, 1, 0.5f));
Vector3 result5 = TestBinaryFunction<Vector3, Vector3, double>(testScaled1).run(Vector3(1, 2, 3), 2.0);
Vector3 result6 = TestBinaryFunction<Vector3, Vector3, double>(testScaled2).run(Vector3(1, 2, 3), 2.0);
Vector3 result7 = TestBinaryFunction<Vector3, Vector3, Matrix4>(testMatrixMultiplied1).run(Vector3(1, 2, 3), matrix4_rotation_for_x_degrees(90));
Vector3 result8 = TestBinaryFunction<Vector3, Vector3, Matrix4>(testMatrixMultiplied2).run(Vector3(1, 2, 3), matrix4_rotation_for_x_degrees(90));
Vector3 result9 = TestUnaryFunction<Vector3, Vector3>(testNormalised1).run(Vector3(1, 2, 3));
Vector3 result10 = TestUnaryFunction<Vector3, Vector3>(testNormalised2).run(Vector3(1, 2, 3));
}
} g_testAll;

617
libs/math/expression.h Normal file
View File

@@ -0,0 +1,617 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#if !defined (INCLUDED_EXPRESSION_H)
#define INCLUDED_EXPRESSION_H
#include <math/matrix.h>
template<typename Value>
class Literal
{
Value m_value;
public:
typedef Value value_type;
Literal(const Value& value)
: m_value(value)
{
}
const value_type& eval() const
{
return m_value;
}
};
template<typename Value>
inline Literal<Value> float_literal(const Value& value)
{
return Literal<Value>(value);
}
template<typename Expression>
inline float float_for_expression(const Expression& expression)
{
return expression.eval();
}
template<typename First, typename Second>
class ScalarDivided
{
First first;
Second second;
public:
typedef typename First::value_type value_type;
ScalarDivided(const First& first_, const Second& second_) : first(first_), second(second_)
{
}
value_type eval() const
{
return static_cast<value_type>(first.eval() / second.eval());
}
};
template<typename First, typename Second>
inline ScalarDivided<First, Second> float_divided(const First& first, const Second& second)
{
return ScalarDivided<First, Second>(first, second);
}
template<typename First>
inline ScalarDivided<Literal<typename First::value_type>, First> float_reciprocal(const First& first)
{
typedef typename First::value_type first_value_type;
return ScalarDivided<Literal<first_value_type>, First>(float_literal(first_value_type(1.0)), first);
}
template<typename First>
class SquareRoot
{
First first;
public:
typedef typename First::value_type value_type;
SquareRoot(const First& first_) : first(first_)
{
}
value_type eval() const
{
return static_cast<value_type>(sqrt(first.eval()));
}
};
template<typename First>
inline SquareRoot<First> float_square_root(const First& first)
{
return SquareRoot<First>(first);
}
template<typename Element>
class BasicVector3Literal
{
const BasicVector3<Element> m_value;
public:
typedef Element value_type;
typedef IntegralConstant<3> dimension;
BasicVector3Literal(const BasicVector3<Element>& value)
: m_value(value)
{
}
const value_type& eval(unsigned int i) const
{
return m_value[i];
}
};
template<typename Element>
inline BasicVector3Literal<Element> vector3_literal(const BasicVector3<Element>& value)
{
return BasicVector3Literal<Element>(value);
}
typedef BasicVector3Literal<float> Vector3Literal;
template<typename Element>
class BasicVector3Identity
{
const BasicVector3<Element>& m_value;
public:
typedef Element value_type;
typedef IntegralConstant<3> dimension;
BasicVector3Identity(const BasicVector3<Element>& value)
: m_value(value)
{
}
const value_type& eval(unsigned int i) const
{
return m_value[i];
}
};
template<typename Element>
inline BasicVector3Identity<Element> vector3_identity(const BasicVector3<Element>& value)
{
return BasicVector3Identity<Element>(value);
}
typedef BasicVector3Identity<float> Vector3Identity;
template<typename Expression>
inline BasicVector3<typename Expression::value_type> vector3_for_expression(const Expression& expression)
{
return Vector3(expression.eval(0), expression.eval(1), expression.eval(2));
}
template<typename Operation, typename First, typename Second>
class VectorScalar
{
First first;
Literal<typename Second::value_type> second;
public:
typedef typename First::value_type value_type;
typedef typename First::dimension dimension;
VectorScalar(const First& first_, const Second& second_)
: first(first_), second(second_.eval())
{
}
value_type eval(unsigned int i) const
{
return Operation::apply( first.eval(i), second.eval() );
}
};
template<typename Operation, typename First, typename Second>
class VectorVector
{
First first;
Second second;
public:
typedef typename First::value_type value_type;
typedef typename First::dimension dimension;
VectorVector(const First& first_, const Second& second_)
: first(first_), second(second_)
{
}
value_type eval(unsigned int i) const
{
return Operation::apply(first.eval(i), second.eval(i));
}
};
template<typename First, typename Second>
class Added
{
public:
typedef First value_type;
static value_type apply(const First& first, const Second& second)
{
return static_cast<value_type>(first + second);
}
};
template<typename First, typename Second>
inline VectorVector<Added<typename First::value_type, typename Second::value_type>, First, Second>
vector_added(const First& first, const Second& second)
{
typedef typename First::value_type first_value_type;
typedef typename Second::value_type second_value_type;
return VectorVector<Added<first_value_type, second_value_type>, First, Second>(first, second);
}
template<typename First, typename Second>
class Multiplied
{
public:
typedef First value_type;
static value_type apply(const First& first, const Second& second)
{
return static_cast<value_type>(first * second);
}
};
template<typename First, typename Second>
inline VectorVector<Multiplied<typename First::value_type, typename Second::value_type>, First, Second>
vector_multiplied(const First& first, const Second& second)
{
typedef typename First::value_type first_value_type;
typedef typename Second::value_type second_value_type;
return VectorVector<Multiplied<first_value_type, second_value_type>, First, Second>(first, second);
}
template<typename First, typename Second>
inline VectorScalar<Multiplied<typename First::value_type, typename Second::value_type>, First, Second>
vector_scaled(const First& first, const Second& second)
{
typedef typename First::value_type first_value_type;
typedef typename Second::value_type second_value_type;
return VectorScalar<Multiplied<first_value_type, second_value_type>, First, Second>(first, second);
}
template<typename First>
class Negated
{
public:
typedef First value_type;
static value_type apply(const First& first)
{
return -first;
}
};
template<typename First, typename Operation>
class VectorUnary
{
First first;
public:
typedef typename First::value_type value_type;
typedef typename First::dimension dimension;
VectorUnary(const First& first_) : first(first_)
{
}
value_type eval(unsigned int i) const
{
return Operation::apply(first.eval(i));
}
};
template<typename First>
inline VectorUnary<First, Negated<typename First::value_type> >
vector_negated(const First& first)
{
typedef typename First::value_type first_value_type;
return VectorUnary<First, Negated<first_value_type> >(first);
}
template<typename First, typename Second>
class VectorCross
{
First first;
Second second;
public:
typedef typename First::value_type value_type;
typedef typename First::dimension dimension;
VectorCross(const First& first_, const Second& second_)
: first(first_), second(second_)
{
}
value_type eval(unsigned int i) const
{
return first.eval((i+1)%3) * second.eval((i+2)%3) - first.eval((i+2)%3) * second.eval((i+1)%3);
}
};
template<typename First, typename Second>
inline VectorCross<First, Second>
vector_cross(const First& first, const Second& second)
{
return VectorCross<First, Second>(first, second);
}
template<typename First, typename Second>
class VectorDot
{
First first;
Second second;
public:
typedef typename First::value_type value_type;
typedef typename First::dimension dimension;
VectorDot(const First& first_, const Second& second_)
: first(first_), second(second_)
{
}
template<typename Index>
struct eval_dot
{
static value_type apply(const First& first, const Second& second)
{
return static_cast<value_type>(
first.eval(Index::VALUE) * second.eval(Index::VALUE)
+ eval_dot< IntegralConstant<Index::VALUE-1> >::apply(first, second)
);
}
};
template<>
struct eval_dot< IntegralConstant<0> >
{
static value_type apply(const First& first, const Second& second)
{
return first.eval(0) * second.eval(0);
}
};
value_type eval() const
{
return eval_dot< IntegralConstant<dimension::VALUE - 1> >::apply(first, second);
}
};
template<typename First, typename Second>
inline VectorDot<First, Second> vector_dot(const First& first, const Second& second)
{
return VectorDot<First, Second>(first, second);
}
template<typename First>
class VectorLengthSquared
{
First first;
public:
typedef typename First::value_type value_type;
typedef typename First::dimension dimension;
VectorLengthSquared(const First& first_)
: first(first_)
{
}
static value_type squared(const value_type& value)
{
return value * value;
}
template<typename Index>
struct eval_squared
{
static value_type apply(const First& first)
{
return static_cast<value_type>(
squared(first.eval(Index::VALUE))
+ eval_squared< IntegralConstant<Index::VALUE - 1> >::apply(first)
);
}
};
template<>
struct eval_squared< IntegralConstant<0> >
{
static value_type apply(const First& first)
{
return squared(first.eval(0));
}
};
value_type eval() const
{
return eval_squared< IntegralConstant<dimension::VALUE - 1> >::apply(first);
}
};
template<typename First>
inline VectorLengthSquared<First> vector_length_squared(const First& first)
{
return VectorLengthSquared<First>(first);
}
template<typename First>
inline SquareRoot< VectorLengthSquared<First> > vector_length(const First& first)
{
return float_square_root(vector_length_squared(first));
}
#if 1
template<typename First>
inline VectorScalar<
Multiplied<typename First::value_type, typename First::value_type>,
First,
// multiple evaulations of subexpression
ScalarDivided<
Literal<typename First::value_type>,
SquareRoot<
VectorLengthSquared<First>
>
>
> vector_normalised(const First& first)
{
typedef typename First::value_type first_value_type;
return vector_scaled(first, float_reciprocal(vector_length(first)));
}
#else
template<typename First>
inline VectorScalar<
Multiplied<typename First::value_type, typename First::value_type>,
First,
// single evaluation of subexpression
Literal<typename First::value_type>
>
vector_normalised(const First& first)
{
typedef typename First::value_type first_value_type;
return vector_scaled(first, float_literal(static_cast<first_value_type>(first_value_type(1.0) / vector_length(first).eval())));
}
#endif
class Matrix4Literal
{
const Matrix4 m_value;
public:
typedef float value_type;
typedef IntegralConstant<4> dimension0;
typedef IntegralConstant<4> dimension1;
Matrix4Literal(const Matrix4& value)
: m_value(value)
{
}
const value_type& eval(unsigned int r, unsigned int c) const
{
return m_value[r*4+c];
}
};
inline Matrix4Literal matrix4_literal(const Matrix4& value)
{
return Matrix4Literal(value);
}
class Matrix4Identity
{
const Matrix4& m_value;
public:
typedef float value_type;
typedef IntegralConstant<4> dimension0;
typedef IntegralConstant<4> dimension1;
Matrix4Identity(const Matrix4& value)
: m_value(value)
{
}
const value_type& eval(unsigned int r, unsigned int c) const
{
return m_value[r*4+c];
}
};
inline Matrix4Identity matrix4_identity(const Matrix4& value)
{
return Matrix4Identity(value);
}
template<typename Expression>
inline Matrix4 matrix4_for_expression(const Expression& expression)
{
return Matrix4(
expression.eval(0, 0), expression.eval(0, 1), expression.eval(0, 2), expression.eval(0, 3),
expression.eval(1, 0), expression.eval(1, 1), expression.eval(1, 2), expression.eval(1, 3),
expression.eval(2, 0), expression.eval(2, 1), expression.eval(2, 2), expression.eval(2, 3),
expression.eval(3, 0), expression.eval(3, 1), expression.eval(3, 2), expression.eval(3, 3)
);
}
template<typename Expression>
inline Matrix4 matrix4_affine_for_expression(const Expression& expression)
{
return Matrix4(
expression.eval(0, 0), expression.eval(0, 1), expression.eval(0, 2), 0,
expression.eval(1, 0), expression.eval(1, 1), expression.eval(1, 2), 0,
expression.eval(2, 0), expression.eval(2, 1), expression.eval(2, 2), 0,
expression.eval(3, 0), expression.eval(3, 1), expression.eval(3, 2), 1
);
}
template<typename First, typename Second>
class PointMultiplied
{
const First& first;
const Second& second;
public:
typedef typename First::value_type value_type;
typedef typename First::dimension dimension;
PointMultiplied(const First& first_, const Second& second_)
: first(first_), second(second_)
{
}
value_type eval(unsigned int i) const
{
return static_cast<value_type>(second.eval(0, i) * first.eval(0)
+ second.eval(1, i) * first.eval(1)
+ second.eval(2, i) * first.eval(2)
+ second.eval(3, i));
}
};
template<typename First, typename Second>
inline PointMultiplied<First, Second> point_multiplied(const First& point, const Second& matrix)
{
return PointMultiplied<First, Second>(point, matrix);
}
template<typename First, typename Second>
class Matrix4Multiplied
{
const First& first;
const Second& second;
public:
typedef typename First::value_type value_type;
typedef typename First::dimension0 dimension0;
typedef typename First::dimension1 dimension1;
Matrix4Multiplied(const First& first_, const Second& second_)
: first(first_), second(second_)
{
}
value_type eval(unsigned int r, unsigned int c) const
{
return static_cast<value_type>(
second.eval(r, 0) * first.eval(0, c)
+ second.eval(r, 1) * first.eval(1, c)
+ second.eval(r, 2) * first.eval(2, c)
+ second.eval(r, 3) * first.eval(3, c)
);
}
};
template<typename First, typename Second>
inline Matrix4Multiplied<First, Second> matrix4_multiplied(const First& first, const Second& second)
{
return Matrix4Multiplied<First, Second>(first, second);
}
template<typename First>
class MatrixTransposed
{
const First& first;
public:
typedef typename First::value_type value_type;
typedef typename First::dimension0 dimension0;
typedef typename First::dimension1 dimension1;
MatrixTransposed(const First& first_)
: first(first_)
{
}
value_type eval(unsigned int r, unsigned int c) const
{
return first.eval(c, r);
}
};
template<typename First>
inline MatrixTransposed<First> matrix_transposed(const First& first)
{
return MatrixTransposed<First>(first);
}
#endif

23
libs/math/frustum.cpp Normal file
View File

@@ -0,0 +1,23 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#include "frustum.h"

629
libs/math/frustum.h Normal file
View File

@@ -0,0 +1,629 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#if !defined(INCLUDED_MATH_FRUSTUM_H)
#define INCLUDED_MATH_FRUSTUM_H
/// \file
/// \brief View-frustum data types and related operations.
#include "generic/enumeration.h"
#include "math/matrix.h"
#include "math/plane.h"
#include "math/aabb.h"
#include "math/line.h"
inline Matrix4 matrix4_frustum(float left, float right, float bottom, float top, float nearval, float farval)
{
return Matrix4(
static_cast<float>( (2*nearval) / (right-left) ),
0,
0,
0,
0,
static_cast<float>( (2*nearval) / (top-bottom) ),
0,
0,
static_cast<float>( (right+left) / (right-left) ),
static_cast<float>( (top+bottom) / (top-bottom) ),
static_cast<float>( -(farval+nearval) / (farval-nearval) ),
-1,
0,
0,
static_cast<float>( -(2*farval*nearval) / (farval-nearval) ),
0
);
}
typedef unsigned char ClipResult;
const ClipResult c_CLIP_PASS = 0x00; // 000000
const ClipResult c_CLIP_LT_X = 0x01; // 000001
const ClipResult c_CLIP_GT_X = 0x02; // 000010
const ClipResult c_CLIP_LT_Y = 0x04; // 000100
const ClipResult c_CLIP_GT_Y = 0x08; // 001000
const ClipResult c_CLIP_LT_Z = 0x10; // 010000
const ClipResult c_CLIP_GT_Z = 0x20; // 100000
const ClipResult c_CLIP_FAIL = 0x3F; // 111111
template<typename Index>
class Vector4ClipLT
{
public:
static bool compare(const Vector4& self)
{
return self[Index::VALUE] < self[3];
}
static double scale(const Vector4& self, const Vector4& other)
{
return (self[Index::VALUE] - self[3]) / (other[3] - other[Index::VALUE]);
}
};
template<typename Index>
class Vector4ClipGT
{
public:
static bool compare(const Vector4& self)
{
return self[Index::VALUE] > -self[3];
}
static double scale(const Vector4& self, const Vector4& other)
{
return (self[Index::VALUE] + self[3]) / (-other[3] - other[Index::VALUE]);
}
};
template<typename ClipPlane>
class Vector4ClipPolygon
{
public:
typedef Vector4* iterator;
typedef const Vector4* const_iterator;
static std::size_t apply(const_iterator first, const_iterator last, iterator out)
{
const_iterator next = first, i = last - 1;
iterator tmp(out);
bool b0 = ClipPlane::compare(*i);
while(next != last)
{
bool b1 = ClipPlane::compare(*next);
if(b0 ^ b1)
{
*out = vector4_subtracted(*next, *i);
double scale = ClipPlane::scale(*i, *out);
(*out)[0] = static_cast<float>((*i)[0] + scale*((*out)[0]));
(*out)[1] = static_cast<float>((*i)[1] + scale*((*out)[1]));
(*out)[2] = static_cast<float>((*i)[2] + scale*((*out)[2]));
(*out)[3] = static_cast<float>((*i)[3] + scale*((*out)[3]));
++out;
}
if(b1)
{
*out = *next;
++out;
}
i = next;
++next;
b0 = b1;
}
return out - tmp;
}
};
#define CLIP_X_LT_W(p) (Vector4ClipLT< IntegralConstant<0> >::compare(p))
#define CLIP_X_GT_W(p) (Vector4ClipGT< IntegralConstant<0> >::compare(p))
#define CLIP_Y_LT_W(p) (Vector4ClipLT< IntegralConstant<1> >::compare(p))
#define CLIP_Y_GT_W(p) (Vector4ClipGT< IntegralConstant<1> >::compare(p))
#define CLIP_Z_LT_W(p) (Vector4ClipLT< IntegralConstant<2> >::compare(p))
#define CLIP_Z_GT_W(p) (Vector4ClipGT< IntegralConstant<2> >::compare(p))
inline ClipResult homogenous_clip_point(const Vector4& clipped)
{
ClipResult result = c_CLIP_FAIL;
if(CLIP_X_LT_W(clipped)) result &= ~c_CLIP_LT_X; // X < W
if(CLIP_X_GT_W(clipped)) result &= ~c_CLIP_GT_X; // X > -W
if(CLIP_Y_LT_W(clipped)) result &= ~c_CLIP_LT_Y; // Y < W
if(CLIP_Y_GT_W(clipped)) result &= ~c_CLIP_GT_Y; // Y > -W
if(CLIP_Z_LT_W(clipped)) result &= ~c_CLIP_LT_Z; // Z < W
if(CLIP_Z_GT_W(clipped)) result &= ~c_CLIP_GT_Z; // Z > -W
return result;
}
/// \brief Clips \p point by canonical matrix \p self.
/// Stores the result in \p clipped.
/// Returns a bitmask indicating which clip-planes the point was outside.
inline ClipResult matrix4_clip_point(const Matrix4& self, const Vector3& point, Vector4& clipped)
{
clipped[0] = point[0];
clipped[1] = point[1];
clipped[2] = point[2];
clipped[3] = 1;
matrix4_transform_vector4(self, clipped);
return homogenous_clip_point(clipped);
}
inline std::size_t homogenous_clip_triangle(Vector4 clipped[9])
{
Vector4 buffer[9];
std::size_t count = 3;
count = Vector4ClipPolygon< Vector4ClipLT< IntegralConstant<0> > >::apply(clipped, clipped + count, buffer);
count = Vector4ClipPolygon< Vector4ClipGT< IntegralConstant<0> > >::apply(buffer, buffer + count, clipped);
count = Vector4ClipPolygon< Vector4ClipLT< IntegralConstant<1> > >::apply(clipped, clipped + count, buffer);
count = Vector4ClipPolygon< Vector4ClipGT< IntegralConstant<1> > >::apply(buffer, buffer + count, clipped);
count = Vector4ClipPolygon< Vector4ClipLT< IntegralConstant<2> > >::apply(clipped, clipped + count, buffer);
return Vector4ClipPolygon< Vector4ClipGT< IntegralConstant<2> > >::apply(buffer, buffer + count, clipped);
}
/// \brief Transforms and clips the triangle formed by \p p0, \p p1, \p p2 by the canonical matrix \p self.
/// Stores the resulting polygon in \p clipped.
/// Returns the number of points in the resulting polygon.
inline std::size_t matrix4_clip_triangle(const Matrix4& self, const Vector3& p0, const Vector3& p1, const Vector3& p2, Vector4 clipped[9])
{
clipped[0][0] = p0[0];
clipped[0][1] = p0[1];
clipped[0][2] = p0[2];
clipped[0][3] = 1;
clipped[1][0] = p1[0];
clipped[1][1] = p1[1];
clipped[1][2] = p1[2];
clipped[1][3] = 1;
clipped[2][0] = p2[0];
clipped[2][1] = p2[1];
clipped[2][2] = p2[2];
clipped[2][3] = 1;
matrix4_transform_vector4(self, clipped[0]);
matrix4_transform_vector4(self, clipped[1]);
matrix4_transform_vector4(self, clipped[2]);
return homogenous_clip_triangle(clipped);
}
inline std::size_t homogenous_clip_line(Vector4 clipped[2])
{
const Vector4& p0 = clipped[0];
const Vector4& p1 = clipped[1];
// early out
{
ClipResult mask0 = homogenous_clip_point(clipped[0]);
ClipResult mask1 = homogenous_clip_point(clipped[1]);
if((mask0 | mask1) == c_CLIP_PASS) // both points passed all planes
return 2;
if(mask0 & mask1) // both points failed any one plane
return 0;
}
{
const bool index = CLIP_X_LT_W(p0);
if(index ^ CLIP_X_LT_W(p1))
{
Vector4 clip(vector4_subtracted(p1, p0));
double scale = (p0[0] - p0[3]) / (clip[3] - clip[0]);
clip[0] = static_cast<float>(p0[0] + scale * clip[0]);
clip[1] = static_cast<float>(p0[1] + scale * clip[1]);
clip[2] = static_cast<float>(p0[2] + scale * clip[2]);
clip[3] = static_cast<float>(p0[3] + scale * clip[3]);
clipped[index] = clip;
}
else if(index == 0)
return 0;
}
{
const bool index = CLIP_X_GT_W(p0);
if(index ^ CLIP_X_GT_W(p1))
{
Vector4 clip(vector4_subtracted(p1, p0));
double scale = (p0[0] + p0[3]) / (-clip[3] - clip[0]);
clip[0] = static_cast<float>(p0[0] + scale * clip[0]);
clip[1] = static_cast<float>(p0[1] + scale * clip[1]);
clip[2] = static_cast<float>(p0[2] + scale * clip[2]);
clip[3] = static_cast<float>(p0[3] + scale * clip[3]);
clipped[index] = clip;
}
else if(index == 0)
return 0;
}
{
const bool index = CLIP_Y_LT_W(p0);
if(index ^ CLIP_Y_LT_W(p1))
{
Vector4 clip(vector4_subtracted(p1, p0));
double scale = (p0[1] - p0[3]) / (clip[3] - clip[1]);
clip[0] = static_cast<float>(p0[0] + scale * clip[0]);
clip[1] = static_cast<float>(p0[1] + scale * clip[1]);
clip[2] = static_cast<float>(p0[2] + scale * clip[2]);
clip[3] = static_cast<float>(p0[3] + scale * clip[3]);
clipped[index] = clip;
}
else if(index == 0)
return 0;
}
{
const bool index = CLIP_Y_GT_W(p0);
if(index ^ CLIP_Y_GT_W(p1))
{
Vector4 clip(vector4_subtracted(p1, p0));
double scale = (p0[1] + p0[3]) / (-clip[3] - clip[1]);
clip[0] = static_cast<float>(p0[0] + scale * clip[0]);
clip[1] = static_cast<float>(p0[1] + scale * clip[1]);
clip[2] = static_cast<float>(p0[2] + scale * clip[2]);
clip[3] = static_cast<float>(p0[3] + scale * clip[3]);
clipped[index] = clip;
}
else if(index == 0)
return 0;
}
{
const bool index = CLIP_Z_LT_W(p0);
if(index ^ CLIP_Z_LT_W(p1))
{
Vector4 clip(vector4_subtracted(p1, p0));
double scale = (p0[2] - p0[3]) / (clip[3] - clip[2]);
clip[0] = static_cast<float>(p0[0] + scale * clip[0]);
clip[1] = static_cast<float>(p0[1] + scale * clip[1]);
clip[2] = static_cast<float>(p0[2] + scale * clip[2]);
clip[3] = static_cast<float>(p0[3] + scale * clip[3]);
clipped[index] = clip;
}
else if(index == 0)
return 0;
}
{
const bool index = CLIP_Z_GT_W(p0);
if(index ^ CLIP_Z_GT_W(p1))
{
Vector4 clip(vector4_subtracted(p1, p0));
double scale = (p0[2] + p0[3]) / (-clip[3] - clip[2]);
clip[0] = static_cast<float>(p0[0] + scale * clip[0]);
clip[1] = static_cast<float>(p0[1] + scale * clip[1]);
clip[2] = static_cast<float>(p0[2] + scale * clip[2]);
clip[3] = static_cast<float>(p0[3] + scale * clip[3]);
clipped[index] = clip;
}
else if(index == 0)
return 0;
}
return 2;
}
/// \brief Transforms and clips the line formed by \p p0, \p p1 by the canonical matrix \p self.
/// Stores the resulting line in \p clipped.
/// Returns the number of points in the resulting line.
inline std::size_t matrix4_clip_line(const Matrix4& self, const Vector3& p0, const Vector3& p1, Vector4 clipped[2])
{
clipped[0][0] = p0[0];
clipped[0][1] = p0[1];
clipped[0][2] = p0[2];
clipped[0][3] = 1;
clipped[1][0] = p1[0];
clipped[1][1] = p1[1];
clipped[1][2] = p1[2];
clipped[1][3] = 1;
matrix4_transform_vector4(self, clipped[0]);
matrix4_transform_vector4(self, clipped[1]);
return homogenous_clip_line(clipped);
}
struct Frustum
{
Plane3 right, left, bottom, top, back, front;
Frustum()
{
}
Frustum(const Plane3& _right,
const Plane3& _left,
const Plane3& _bottom,
const Plane3& _top,
const Plane3& _back,
const Plane3& _front)
: right(_right), left(_left), bottom(_bottom), top(_top), back(_back), front(_front)
{
}
};
inline Frustum frustum_transformed(const Frustum& frustum, const Matrix4& transform)
{
return Frustum(
plane3_transformed(frustum.right, transform),
plane3_transformed(frustum.left, transform),
plane3_transformed(frustum.bottom, transform),
plane3_transformed(frustum.top, transform),
plane3_transformed(frustum.back, transform),
plane3_transformed(frustum.front, transform)
);
}
inline Frustum frustum_inverse_transformed(const Frustum& frustum, const Matrix4& transform)
{
return Frustum(
plane3_inverse_transformed(frustum.right, transform),
plane3_inverse_transformed(frustum.left, transform),
plane3_inverse_transformed(frustum.bottom, transform),
plane3_inverse_transformed(frustum.top, transform),
plane3_inverse_transformed(frustum.back, transform),
plane3_inverse_transformed(frustum.front, transform)
);
}
inline bool viewproj_test_point(const Matrix4& viewproj, const Vector3& point)
{
Vector4 hpoint(matrix4_transformed_vector4(viewproj, Vector4(point, 1.0f)));
if(fabs(hpoint[0]) < fabs(hpoint[3])
&& fabs(hpoint[1]) < fabs(hpoint[3])
&& fabs(hpoint[2]) < fabs(hpoint[3]))
return true;
return false;
}
inline bool viewproj_test_transformed_point(const Matrix4& viewproj, const Vector3& point, const Matrix4& localToWorld)
{
return viewproj_test_point(viewproj, matrix4_transformed_point(localToWorld, point));
}
inline Frustum frustum_from_viewproj(const Matrix4& viewproj)
{
return Frustum
(
plane3_normalised(Plane3(viewproj[ 3] - viewproj[ 0], viewproj[ 7] - viewproj[ 4], viewproj[11] - viewproj[ 8], viewproj[15] - viewproj[12])),
plane3_normalised(Plane3(viewproj[ 3] + viewproj[ 0], viewproj[ 7] + viewproj[ 4], viewproj[11] + viewproj[ 8], viewproj[15] + viewproj[12])),
plane3_normalised(Plane3(viewproj[ 3] + viewproj[ 1], viewproj[ 7] + viewproj[ 5], viewproj[11] + viewproj[ 9], viewproj[15] + viewproj[13])),
plane3_normalised(Plane3(viewproj[ 3] - viewproj[ 1], viewproj[ 7] - viewproj[ 5], viewproj[11] - viewproj[ 9], viewproj[15] - viewproj[13])),
plane3_normalised(Plane3(viewproj[ 3] - viewproj[ 2], viewproj[ 7] - viewproj[ 6], viewproj[11] - viewproj[10], viewproj[15] - viewproj[14])),
plane3_normalised(Plane3(viewproj[ 3] + viewproj[ 2], viewproj[ 7] + viewproj[ 6], viewproj[11] + viewproj[10], viewproj[15] + viewproj[14]))
);
}
struct VolumeIntersection
{
enum Value
{
OUTSIDE,
INSIDE,
PARTIAL
};
};
typedef EnumeratedValue<VolumeIntersection> VolumeIntersectionValue;
const VolumeIntersectionValue c_volumeOutside(VolumeIntersectionValue::OUTSIDE);
const VolumeIntersectionValue c_volumeInside(VolumeIntersectionValue::INSIDE);
const VolumeIntersectionValue c_volumePartial(VolumeIntersectionValue::PARTIAL);
inline VolumeIntersectionValue frustum_test_aabb(const Frustum& frustum, const AABB& aabb)
{
VolumeIntersectionValue result = c_volumeInside;
switch(aabb_classify_plane(aabb, frustum.right))
{
case 2:
return c_volumeOutside;
case 1:
result = c_volumePartial;
}
switch(aabb_classify_plane(aabb, frustum.left))
{
case 2:
return c_volumeOutside;
case 1:
result = c_volumePartial;
}
switch(aabb_classify_plane(aabb, frustum.bottom))
{
case 2:
return c_volumeOutside;
case 1:
result = c_volumePartial;
}
switch(aabb_classify_plane(aabb, frustum.top))
{
case 2:
return c_volumeOutside;
case 1:
result = c_volumePartial;
}
switch(aabb_classify_plane(aabb, frustum.back))
{
case 2:
return c_volumeOutside;
case 1:
result = c_volumePartial;
}
switch(aabb_classify_plane(aabb, frustum.front))
{
case 2:
return c_volumeOutside;
case 1:
result = c_volumePartial;
}
return result;
}
inline double plane_distance_to_point(const Plane3& plane, const Vector3& point)
{
return vector3_dot(plane.normal(), point) + plane.d;
}
inline double plane_distance_to_oriented_extents(const Plane3& plane, const Vector3& extents, const Matrix4& orientation)
{
return fabs(extents[0] * vector3_dot(plane.normal(), vector4_to_vector3(orientation.x())))
+ fabs(extents[1] * vector3_dot(plane.normal(), vector4_to_vector3(orientation.y())))
+ fabs(extents[2] * vector3_dot(plane.normal(), vector4_to_vector3(orientation.z())));
}
/// \brief Return false if \p aabb with \p orientation is partially or completely outside \p plane.
inline bool plane_contains_oriented_aabb(const Plane3& plane, const AABB& aabb, const Matrix4& orientation)
{
double dot = plane_distance_to_point(plane, aabb.origin);
return !(dot > 0 || -dot < plane_distance_to_oriented_extents(plane, aabb.extents, orientation));
}
inline VolumeIntersectionValue frustum_intersects_transformed_aabb(const Frustum& frustum, const AABB& aabb, const Matrix4& localToWorld)
{
AABB aabb_world(aabb);
matrix4_transform_point(localToWorld, aabb_world.origin);
if(plane_contains_oriented_aabb(frustum.right, aabb_world, localToWorld)
|| plane_contains_oriented_aabb(frustum.left, aabb_world, localToWorld)
|| plane_contains_oriented_aabb(frustum.bottom, aabb_world, localToWorld)
|| plane_contains_oriented_aabb(frustum.top, aabb_world, localToWorld)
|| plane_contains_oriented_aabb(frustum.back, aabb_world, localToWorld)
|| plane_contains_oriented_aabb(frustum.front, aabb_world, localToWorld))
return c_volumeOutside;
return c_volumeInside;
}
inline bool plane3_test_point(const Plane3& plane, const Vector3& point)
{
return vector3_dot(point, plane.normal()) + plane.dist() <= 0;
}
inline bool plane3_test_line(const Plane3& plane, const Segment& segment)
{
return segment_classify_plane(segment, plane) == 2;
}
inline bool frustum_test_point(const Frustum& frustum, const Vector3& point)
{
return !plane3_test_point(frustum.right, point)
&& !plane3_test_point(frustum.left, point)
&& !plane3_test_point(frustum.bottom, point)
&& !plane3_test_point(frustum.top, point)
&& !plane3_test_point(frustum.back, point)
&& !plane3_test_point(frustum.front, point);
}
inline bool frustum_test_line(const Frustum& frustum, const Segment& segment)
{
return !plane3_test_line(frustum.right, segment)
&& !plane3_test_line(frustum.left, segment)
&& !plane3_test_line(frustum.bottom, segment)
&& !plane3_test_line(frustum.top, segment)
&& !plane3_test_line(frustum.back, segment)
&& !plane3_test_line(frustum.front, segment);
}
inline bool viewer_test_plane(const Vector4& viewer, const Plane3& plane)
{
return ((plane.a * viewer[0])
+ (plane.b * viewer[1])
+ (plane.c * viewer[2])
+ (plane.d * viewer[3])) > 0;
}
inline Vector3 triangle_cross(const Vector3& p0, const Vector3& p1, const Vector3& p2)
{
return vector3_cross(vector3_subtracted(p1, p0), vector3_subtracted(p1, p2));
}
inline bool viewer_test_triangle(const Vector4& viewer, const Vector3& p0, const Vector3& p1, const Vector3& p2)
{
Vector3 cross(triangle_cross(p0, p1, p2));
return ((viewer[0] * cross[0])
+ (viewer[1] * cross[1])
+ (viewer[2] * cross[2])
+ (viewer[3] * 0)) > 0;
}
inline Vector4 viewer_from_transformed_viewer(const Vector4& viewer, const Matrix4& transform)
{
if(viewer[3] == 0)
{
return Vector4(matrix4_transformed_direction(transform, vector4_to_vector3(viewer)), 0);
}
else
{
return Vector4(matrix4_transformed_point(transform, vector4_to_vector3(viewer)), viewer[3]);
}
}
inline bool viewer_test_transformed_plane(const Vector4& viewer, const Plane3& plane, const Matrix4& localToWorld)
{
#if 0
return viewer_test_plane(viewer_from_transformed_viewer(viewer, matrix4_affine_inverse(localToWorld)), plane);
#else
return viewer_test_plane(viewer, plane3_transformed(plane, localToWorld));
#endif
}
inline Vector4 viewer_from_viewproj(const Matrix4& viewproj)
{
// get viewer pos in object coords
Vector4 viewer(matrix4_transformed_vector4(matrix4_full_inverse(viewproj), Vector4(0, 0, -1, 0)));
if(viewer[3] != 0) // non-affine matrix
{
viewer[0] /= viewer[3];
viewer[1] /= viewer[3];
viewer[2] /= viewer[3];
viewer[3] /= viewer[3];
}
return viewer;
}
#endif

23
libs/math/line.cpp Normal file
View File

@@ -0,0 +1,23 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#include "line.h"

151
libs/math/line.h Normal file
View File

@@ -0,0 +1,151 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#if !defined(INCLUDED_MATH_LINE_H)
#define INCLUDED_MATH_LINE_H
/// \file
/// \brief Line data types and related operations.
#include "math/vector.h"
#include "math/plane.h"
/// \brief A line segment defined by a start point and and end point.
class Line
{
public:
Vector3 start, end;
Line()
{
}
Line(const Vector3& start_, const Vector3& end_) : start(start_), end(end_)
{
}
};
inline Vector3 line_closest_point(const Line& line, const Vector3& point)
{
Vector3 v = line.end - line.start;
Vector3 w = point - line.start;
double c1 = vector3_dot(w,v);
if ( c1 <= 0 )
return line.start;
double c2 = vector3_dot(v,v);
if ( c2 <= c1 )
return line.end;
return Vector3(line.start + v * (c1 / c2));
}
class Segment
{
public:
Vector3 origin, extents;
Segment()
{
}
Segment(const Vector3& origin_, const Vector3& extents_) :
origin(origin_), extents(extents_)
{
}
};
inline Segment segment_for_startend(const Vector3& start, const Vector3& end)
{
Segment segment;
segment.origin = vector3_mid(start, end);
segment.extents = vector3_subtracted(end, segment.origin);
return segment;
}
inline unsigned int segment_classify_plane(const Segment& segment, const Plane3& plane)
{
double distance_origin = vector3_dot(plane.normal(), segment.origin) + plane.dist();
if (fabs(distance_origin) < fabs(vector3_dot(plane.normal(), segment.extents)))
{
return 1; // partially inside
}
else if (distance_origin < 0)
{
return 2; // totally inside
}
return 0; // totally outside
}
class Ray
{
public:
Vector3 origin, direction;
Ray()
{
}
Ray(const Vector3& origin_, const Vector3& direction_) :
origin(origin_), direction(direction_)
{
}
};
inline Ray ray_for_points(const Vector3& origin, const Vector3& p2)
{
return Ray(origin, vector3_normalised(vector3_subtracted(p2, origin)));
}
inline void ray_transform(Ray& ray, const Matrix4& matrix)
{
matrix4_transform_point(matrix, ray.origin);
matrix4_transform_direction(matrix, ray.direction);
}
// closest-point-on-line
inline double ray_squared_distance_to_point(const Ray& ray, const Vector3& point)
{
return vector3_length_squared(
vector3_subtracted(
point,
vector3_added(
ray.origin,
vector3_scaled(
ray.direction,
vector3_dot(
vector3_subtracted(point, ray.origin),
ray.direction
)
)
)
)
);
}
inline double ray_distance_to_plane(const Ray& ray, const Plane3& plane)
{
return -(vector3_dot(plane.normal(), ray.origin) - plane.dist()) / vector3_dot(ray.direction, plane.normal());
}
#endif

23
libs/math/matrix.cpp Normal file
View File

@@ -0,0 +1,23 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#include "matrix.h"

1322
libs/math/matrix.h Normal file

File diff suppressed because it is too large Load Diff

23
libs/math/pi.cpp Normal file
View File

@@ -0,0 +1,23 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#include "pi.h"

45
libs/math/pi.h Normal file
View File

@@ -0,0 +1,45 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#if !defined(INCLUDED_MATH_PI_H)
#define INCLUDED_MATH_PI_H
/// \file
/// \brief Pi constants and degrees/radians conversion.
const double c_pi = 3.1415926535897932384626433832795;
const double c_half_pi = c_pi / 2;
const double c_2pi = 2 * c_pi;
const double c_inv_2pi = 1 / c_2pi;
const double c_DEG2RADMULT = c_pi / 180.0;
const double c_RAD2DEGMULT = 180.0 / c_pi;
inline double radians_to_degrees(double radians)
{
return radians * c_RAD2DEGMULT;
}
inline double degrees_to_radians(double degrees)
{
return degrees * c_DEG2RADMULT;
}
#endif

23
libs/math/plane.cpp Normal file
View File

@@ -0,0 +1,23 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#include "plane.h"

153
libs/math/plane.h Normal file
View File

@@ -0,0 +1,153 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#if !defined(INCLUDED_MATH_PLANE_H)
#define INCLUDED_MATH_PLANE_H
/// \file
/// \brief Plane data types and related operations.
#include "math/matrix.h"
/// \brief A plane equation stored in double-precision floating-point.
class Plane3
{
public:
double a, b, c, d;
Plane3()
{
}
Plane3(double _a, double _b, double _c, double _d)
: a(_a), b(_b), c(_c), d(_d)
{
}
template<typename Element>
Plane3(const BasicVector3<Element>& normal, double dist)
: a(normal.x()), b(normal.y()), c(normal.z()), d(dist)
{
}
BasicVector3<double>& normal()
{
return reinterpret_cast<BasicVector3<double>&>(*this);
}
const BasicVector3<double>& normal() const
{
return reinterpret_cast<const BasicVector3<double>&>(*this);
}
double& dist()
{
return d;
}
const double& dist() const
{
return d;
}
};
inline Plane3 plane3_normalised(const Plane3& plane)
{
double rmagnitude = 1.0 / sqrt(plane.a * plane.a + plane.b * plane.b + plane.c * plane.c);
return Plane3(
plane.a * rmagnitude,
plane.b * rmagnitude,
plane.c * rmagnitude,
plane.d * rmagnitude
);
}
inline Plane3 plane3_translated(const Plane3& plane, const Vector3& translation)
{
Plane3 transformed;
transformed.a = plane.a;
transformed.b = plane.b;
transformed.c = plane.c;
transformed.d = -((-plane.d * transformed.a + translation.x()) * transformed.a +
(-plane.d * transformed.b + translation.y()) * transformed.b +
(-plane.d * transformed.c + translation.z()) * transformed.c);
return transformed;
}
inline Plane3 plane3_transformed(const Plane3& plane, const Matrix4& transform)
{
Plane3 transformed;
transformed.a = transform[0] * plane.a + transform[4] * plane.b + transform[8] * plane.c;
transformed.b = transform[1] * plane.a + transform[5] * plane.b + transform[9] * plane.c;
transformed.c = transform[2] * plane.a + transform[6] * plane.b + transform[10] * plane.c;
transformed.d = -((-plane.d * transformed.a + transform[12]) * transformed.a +
(-plane.d * transformed.b + transform[13]) * transformed.b +
(-plane.d * transformed.c + transform[14]) * transformed.c);
return transformed;
}
inline Plane3 plane3_inverse_transformed(const Plane3& plane, const Matrix4& transform)
{
return Plane3
(
transform[ 0] * plane.a + transform[ 1] * plane.b + transform[ 2] * plane.c + transform[ 3] * plane.d,
transform[ 4] * plane.a + transform[ 5] * plane.b + transform[ 6] * plane.c + transform[ 7] * plane.d,
transform[ 8] * plane.a + transform[ 9] * plane.b + transform[10] * plane.c + transform[11] * plane.d,
transform[12] * plane.a + transform[13] * plane.b + transform[14] * plane.c + transform[15] * plane.d
);
}
inline Plane3 plane3_flipped(const Plane3& plane)
{
return Plane3(vector3_negated(plane.normal()), -plane.dist());
}
const double c_PLANE_NORMAL_EPSILON = 0.0001f;
const double c_PLANE_DIST_EPSILON = 0.02;
inline bool plane3_equal(const Plane3& self, const Plane3& other)
{
return vector3_equal_epsilon(self.normal(), other.normal(), c_PLANE_NORMAL_EPSILON)
&& float_equal_epsilon(self.dist(), other.dist(), c_PLANE_DIST_EPSILON);
}
inline bool plane3_opposing(const Plane3& self, const Plane3& other)
{
return plane3_equal(self, plane3_flipped(other));
}
inline bool plane3_valid(const Plane3& self)
{
return float_equal_epsilon(vector3_dot(self.normal(), self.normal()), 1.0, 0.01);
}
template<typename Element>
inline Plane3 plane3_for_points(const BasicVector3<Element>& p0, const BasicVector3<Element>& p1, const BasicVector3<Element>& p2)
{
Plane3 self;
self.normal() = vector3_normalised(vector3_cross(vector3_subtracted(p1, p0), vector3_subtracted(p2, p0)));
self.dist() = vector3_dot(p0, self.normal());
return self;
}
template<typename Element>
inline Plane3 plane3_for_points(const BasicVector3<Element> planepts[3])
{
return plane3_for_points(planepts[2], planepts[1], planepts[0]);
}
#endif

23
libs/math/quaternion.cpp Normal file
View File

@@ -0,0 +1,23 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#include "quaternion.h"

309
libs/math/quaternion.h Normal file
View File

@@ -0,0 +1,309 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#if !defined(INCLUDED_MATH_QUATERNION_H)
#define INCLUDED_MATH_QUATERNION_H
/// \file
/// \brief Quaternion data types and related operations.
#include "math/matrix.h"
/// \brief A quaternion stored in single-precision floating-point.
typedef Vector4 Quaternion;
const Quaternion c_quaternion_identity(0, 0, 0, 1);
inline Quaternion quaternion_multiplied_by_quaternion(const Quaternion& quaternion, const Quaternion& other)
{
return Quaternion(
quaternion[3]*other[0] + quaternion[0]*other[3] + quaternion[1]*other[2] - quaternion[2]*other[1],
quaternion[3]*other[1] + quaternion[1]*other[3] + quaternion[2]*other[0] - quaternion[0]*other[2],
quaternion[3]*other[2] + quaternion[2]*other[3] + quaternion[0]*other[1] - quaternion[1]*other[0],
quaternion[3]*other[3] - quaternion[0]*other[0] - quaternion[1]*other[1] - quaternion[2]*other[2]
);
}
inline void quaternion_multiply_by_quaternion(Quaternion& quaternion, const Quaternion& other)
{
quaternion = quaternion_multiplied_by_quaternion(quaternion, other);
}
/// \brief Constructs a quaternion which rotates between two points on the unit-sphere, \p from and \p to.
inline Quaternion quaternion_for_unit_vectors(const Vector3& from, const Vector3& to)
{
return Quaternion(vector3_cross(from, to), static_cast<float>(vector3_dot(from, to)));
}
inline Quaternion quaternion_for_axisangle(const Vector3& axis, double angle)
{
angle *= 0.5;
float sa = static_cast<float>(sin(angle));
return Quaternion(axis[0] * sa, axis[1] * sa, axis[2] * sa, static_cast<float>(cos(angle)));
}
inline Quaternion quaternion_inverse(const Quaternion& quaternion)
{
return Quaternion(vector3_negated(vector4_to_vector3(quaternion)), quaternion[3]);
}
inline void quaternion_conjugate(Quaternion& quaternion)
{
quaternion = quaternion_inverse(quaternion);
}
inline Quaternion quaternion_normalised(const Quaternion& quaternion)
{
const double n = (1.0 / (quaternion[0] * quaternion[0] + quaternion[1] * quaternion[1] + quaternion[2] * quaternion[2] + quaternion[3] * quaternion[3]));
return Quaternion(
static_cast<float>(quaternion[0] * n),
static_cast<float>(quaternion[1] * n),
static_cast<float>(quaternion[2] * n),
static_cast<float>(quaternion[3] * n)
);
}
inline void quaternion_normalise(Quaternion& quaternion)
{
quaternion = quaternion_normalised(quaternion);
}
/// \brief Constructs a pure-rotation matrix from \p quaternion.
inline Matrix4 matrix4_rotation_for_quaternion(const Quaternion& quaternion)
{
#if 0
const double xx = quaternion[0] * quaternion[0];
const double xy = quaternion[0] * quaternion[1];
const double xz = quaternion[0] * quaternion[2];
const double xw = quaternion[0] * quaternion[3];
const double yy = quaternion[1] * quaternion[1];
const double yz = quaternion[1] * quaternion[2];
const double yw = quaternion[1] * quaternion[3];
const double zz = quaternion[2] * quaternion[2];
const double zw = quaternion[2] * quaternion[3];
return Matrix4(
static_cast<float>( 1 - 2 * ( yy + zz ) ),
static_cast<float>( 2 * ( xy + zw ) ),
static_cast<float>( 2 * ( xz - yw ) ),
0,
static_cast<float>( 2 * ( xy - zw ) ),
static_cast<float>( 1 - 2 * ( xx + zz ) ),
static_cast<float>( 2 * ( yz + xw ) ),
0,
static_cast<float>( 2 * ( xz + yw ) ),
static_cast<float>( 2 * ( yz - xw ) ),
static_cast<float>( 1 - 2 * ( xx + yy ) ),
0,
0,
0,
0,
1
);
#else
const double x2 = quaternion[0] + quaternion[0];
const double y2 = quaternion[1] + quaternion[1];
const double z2 = quaternion[2] + quaternion[2];
const double xx = quaternion[0] * x2;
const double xy = quaternion[0] * y2;
const double xz = quaternion[0] * z2;
const double yy = quaternion[1] * y2;
const double yz = quaternion[1] * z2;
const double zz = quaternion[2] * z2;
const double wx = quaternion[3] * x2;
const double wy = quaternion[3] * y2;
const double wz = quaternion[3] * z2;
return Matrix4(
static_cast<float>( 1.0 - (yy + zz) ),
static_cast<float>(xy + wz),
static_cast<float>(xz - wy),
0,
static_cast<float>(xy - wz),
static_cast<float>( 1.0 - (xx + zz) ),
static_cast<float>(yz + wx),
0,
static_cast<float>(xz + wy),
static_cast<float>(yz - wx),
static_cast<float>( 1.0 - (xx + yy) ),
0,
0,
0,
0,
1
);
#endif
}
const double c_half_sqrt2 = 0.70710678118654752440084436210485;
const float c_half_sqrt2f = static_cast<float>(c_half_sqrt2);
inline bool quaternion_component_is_90(float component)
{
return (fabs(component) - c_half_sqrt2) < 0.001;
}
inline Matrix4 matrix4_rotation_for_quaternion_quantised(const Quaternion& quaternion)
{
if(quaternion.y() == 0
&& quaternion.z() == 0
&& quaternion_component_is_90(quaternion.x())
&& quaternion_component_is_90(quaternion.w()))
{
return matrix4_rotation_for_sincos_x((quaternion.x() > 0) ? 1.f : -1.f, 0);
}
if(quaternion.x() == 0
&& quaternion.z() == 0
&& quaternion_component_is_90(quaternion.y())
&& quaternion_component_is_90(quaternion.w()))
{
return matrix4_rotation_for_sincos_y((quaternion.y() > 0) ? 1.f : -1.f, 0);
}
if(quaternion.x() == 0
&& quaternion.y() == 0
&& quaternion_component_is_90(quaternion.z())
&& quaternion_component_is_90(quaternion.w()))
{
return matrix4_rotation_for_sincos_z((quaternion.z() > 0) ? 1.f : -1.f, 0);
}
return matrix4_rotation_for_quaternion(quaternion);
}
inline Quaternion quaternion_for_matrix4_rotation(const Matrix4& matrix4)
{
Matrix4 transposed = matrix4_transposed(matrix4);
double trace = transposed[0] + transposed[5] + transposed[10] + 1.0;
if(trace > 0.0001)
{
double S = 0.5 / sqrt(trace);
return Quaternion(
static_cast<float>((transposed[9] - transposed[6]) * S),
static_cast<float>((transposed[2] - transposed[8]) * S),
static_cast<float>((transposed[4] - transposed[1]) * S),
static_cast<float>(0.25 / S)
);
}
if(transposed[0] >= transposed[5] && transposed[0] >= transposed[10])
{
double S = 2.0 * sqrt(1.0 + transposed[0] - transposed[5] - transposed[10]);
return Quaternion(
static_cast<float>(0.25 / S),
static_cast<float>((transposed[1] + transposed[4]) / S),
static_cast<float>((transposed[2] + transposed[8]) / S),
static_cast<float>((transposed[6] + transposed[9]) / S)
);
}
if(transposed[5] >= transposed[0] && transposed[5] >= transposed[10])
{
double S = 2.0 * sqrt(1.0 + transposed[5] - transposed[0] - transposed[10]);
return Quaternion(
static_cast<float>((transposed[1] + transposed[4]) / S),
static_cast<float>(0.25 / S),
static_cast<float>((transposed[6] + transposed[9]) / S),
static_cast<float>((transposed[2] + transposed[8]) / S)
);
}
double S = 2.0 * sqrt(1.0 + transposed[10] - transposed[0] - transposed[5]);
return Quaternion(
static_cast<float>((transposed[2] + transposed[8]) / S),
static_cast<float>((transposed[6] + transposed[9]) / S),
static_cast<float>(0.25 / S),
static_cast<float>((transposed[1] + transposed[4]) / S)
);
}
/// \brief Returns \p self concatenated with the rotation transform produced by \p rotation.
/// The concatenated rotation occurs before \p self.
inline Matrix4 matrix4_rotated_by_quaternion(const Matrix4& self, const Quaternion& rotation)
{
return matrix4_multiplied_by_matrix4(self, matrix4_rotation_for_quaternion(rotation));
}
/// \brief Concatenates \p self with the rotation transform produced by \p rotation.
/// The concatenated rotation occurs before \p self.
inline void matrix4_rotate_by_quaternion(Matrix4& self, const Quaternion& rotation)
{
self = matrix4_rotated_by_quaternion(self, rotation);
}
/// \brief Rotates \p self by \p rotation, using \p pivotpoint.
inline void matrix4_pivoted_rotate_by_quaternion(Matrix4& self, const Quaternion& rotation, const Vector3& pivotpoint)
{
matrix4_translate_by_vec3(self, pivotpoint);
matrix4_rotate_by_quaternion(self, rotation);
matrix4_translate_by_vec3(self, vector3_negated(pivotpoint));
}
inline Vector3 quaternion_transformed_point(const Quaternion& quaternion, const Vector3& point)
{
double xx = quaternion.x() * quaternion.x();
double yy = quaternion.y() * quaternion.y();
double zz = quaternion.z() * quaternion.z();
double ww = quaternion.w() * quaternion.w();
double xy2 = quaternion.x() * quaternion.y() * 2;
double xz2 = quaternion.x() * quaternion.z() * 2;
double xw2 = quaternion.x() * quaternion.w() * 2;
double yz2 = quaternion.y() * quaternion.z() * 2;
double yw2 = quaternion.y() * quaternion.w() * 2;
double zw2 = quaternion.z() * quaternion.w() * 2;
return Vector3(
static_cast<float>(ww * point.x() + yw2 * point.z() - zw2 * point.y() + xx * point.x() + xy2 * point.y() + xz2 * point.z() - zz * point.x() - yy * point.x()),
static_cast<float>(xy2 * point.x() + yy * point.y() + yz2 * point.z() + zw2 * point.x() - zz * point.y() + ww * point.y() - xw2 * point.z() - xx * point.y()),
static_cast<float>(xz2 * point.x() + yz2 * point.y() + zz * point.z() - yw2 * point.x() - yy * point.z() + xw2 * point.y() - xx * point.z() + ww * point.z())
);
}
/// \brief Constructs a pure-rotation transform from \p axis and \p angle (radians).
inline Matrix4 matrix4_rotation_for_axisangle(const Vector3& axis, double angle)
{
return matrix4_rotation_for_quaternion(quaternion_for_axisangle(axis, angle));
}
/// \brief Rotates \p self about \p axis by \p angle.
inline void matrix4_rotate_by_axisangle(Matrix4& self, const Vector3& axis, double angle)
{
matrix4_multiply_by_matrix4(self, matrix4_rotation_for_axisangle(axis, angle));
}
/// \brief Rotates \p self about \p axis by \p angle using \p pivotpoint.
inline void matrix4_pivoted_rotate_by_axisangle(Matrix4& self, const Vector3& axis, double angle, const Vector3& pivotpoint)
{
matrix4_translate_by_vec3(self, pivotpoint);
matrix4_rotate_by_axisangle(self, axis, angle);
matrix4_translate_by_vec3(self, vector3_negated(pivotpoint));
}
#endif

23
libs/math/vector.cpp Normal file
View File

@@ -0,0 +1,23 @@
/*
Copyright (C) 2001-2006, William Joseph.
All Rights Reserved.
This file is part of GtkRadiant.
GtkRadiant is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
GtkRadiant is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with GtkRadiant; if not, write to the Free Software
Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
#include "vector.h"

1030
libs/math/vector.h Normal file

File diff suppressed because it is too large Load Diff