From 22d90711dab0a4a98c201051d0cee0d7eef68ca8 Mon Sep 17 00:00:00 2001 From: lawnjelly Date: Tue, 9 May 2023 13:27:09 +0100 Subject: [PATCH] Make acos and asin safe A common bug with using acos and asin is that input outside -1 to 1 range will result in Nan output. This can occur due to floating point error in the input. The standard solution is to provide safe_acos function with clamped input. For Godot it may make more sense to make the standard functions safe. --- core/math/basis.cpp | 1 + core/math/math_funcs.h | 10 ++++++---- core/math/quaternion.cpp | 4 +++- modules/gdscript/doc_classes/@GDScript.xml | 4 ++-- modules/skeleton_3d/nodes/skeleton_ik.cpp | 3 ++- scene/3d/path.cpp | 4 +++- scene/resources/sky.cpp | 6 ++++-- 7 files changed, 21 insertions(+), 11 deletions(-) diff --git a/core/math/basis.cpp b/core/math/basis.cpp index 4623ce029..5e79f27b0 100644 --- a/core/math/basis.cpp +++ b/core/math/basis.cpp @@ -974,6 +974,7 @@ void Basis::get_axis_angle(Vector3 &r_axis, real_t &r_angle) const { // as we have reached here there are no singularities so we can handle normally real_t s = Math::sqrt((rows[1][2] - rows[2][1]) * (rows[1][2] - rows[2][1]) + (rows[2][0] - rows[0][2]) * (rows[2][0] - rows[0][2]) + (rows[0][1] - rows[1][0]) * (rows[0][1] - rows[1][0])); // s=|axis||sin(angle)|, used to normalise + // acos does clamping. angle = Math::acos((rows[0][0] + rows[1][1] + rows[2][2] - 1) / 2); if (angle < 0) { s = -s; diff --git a/core/math/math_funcs.h b/core/math/math_funcs.h index 5dd1abb25..450b09bc4 100644 --- a/core/math/math_funcs.h +++ b/core/math/math_funcs.h @@ -73,11 +73,13 @@ public: static _ALWAYS_INLINE_ double tanh(double p_x) { return ::tanh(p_x); } static _ALWAYS_INLINE_ float tanh(float p_x) { return ::tanhf(p_x); } - static _ALWAYS_INLINE_ double asin(double p_x) { return ::asin(p_x); } - static _ALWAYS_INLINE_ float asin(float p_x) { return ::asinf(p_x); } + // Always does clamping so always safe to use. + static _ALWAYS_INLINE_ double asin(double p_x) { return p_x < -1 ? (-Math_PI / 2) : (p_x > 1 ? (Math_PI / 2) : ::asin(p_x)); } + static _ALWAYS_INLINE_ float asin(float p_x) { return p_x < -1 ? (-Math_PI / 2) : (p_x > 1 ? (Math_PI / 2) : ::asinf(p_x)); } - static _ALWAYS_INLINE_ double acos(double p_x) { return ::acos(p_x); } - static _ALWAYS_INLINE_ float acos(float p_x) { return ::acosf(p_x); } + // Always does clamping so always safe to use. + static _ALWAYS_INLINE_ double acos(double p_x) { return p_x < -1 ? Math_PI : (p_x > 1 ? 0 : ::acos(p_x)); } + static _ALWAYS_INLINE_ float acos(float p_x) { return p_x < -1 ? Math_PI : (p_x > 1 ? 0 : ::acosf(p_x)); } static _ALWAYS_INLINE_ double atan(double p_x) { return ::atan(p_x); } static _ALWAYS_INLINE_ float atan(float p_x) { return ::atanf(p_x); } diff --git a/core/math/quaternion.cpp b/core/math/quaternion.cpp index 1b221b159..dbb5ec8a6 100644 --- a/core/math/quaternion.cpp +++ b/core/math/quaternion.cpp @@ -35,7 +35,9 @@ real_t Quaternion::angle_to(const Quaternion &p_to) const { real_t d = dot(p_to); - return Math::acos(CLAMP(d * d * 2 - 1, -1, 1)); + + // acos does clamping. + return Math::acos(d * d * 2 - 1); } // set_euler_xyz expects a vector containing the Euler angles in the format diff --git a/modules/gdscript/doc_classes/@GDScript.xml b/modules/gdscript/doc_classes/@GDScript.xml index 65384d6e1..ba5259fa4 100644 --- a/modules/gdscript/doc_classes/@GDScript.xml +++ b/modules/gdscript/doc_classes/@GDScript.xml @@ -52,7 +52,7 @@ - Returns the arc cosine of [code]s[/code] in radians. Use to get the angle of cosine [code]s[/code]. [code]s[/code] must be between [code]-1.0[/code] and [code]1.0[/code] (inclusive), otherwise, [method acos] will return [constant NAN]. + Returns the arc cosine of [code]s[/code] in radians. Use to get the angle of cosine [code]s[/code]. [code]s[/code] will be clamped between [code]-1.0[/code] and [code]1.0[/code] (inclusive), in order to prevent [method acos] from returning [constant NAN]. [codeblock] # c is 0.523599 or 30 degrees if converted with rad2deg(s) c = acos(0.866025) @@ -63,7 +63,7 @@ - Returns the arc sine of [code]s[/code] in radians. Use to get the angle of sine [code]s[/code]. [code]s[/code] must be between [code]-1.0[/code] and [code]1.0[/code] (inclusive), otherwise, [method asin] will return [constant NAN]. + Returns the arc sine of [code]s[/code] in radians. Use to get the angle of sine [code]s[/code]. [code]s[/code] will be clamped between [code]-1.0[/code] and [code]1.0[/code] (inclusive), in order to prevent [method asin] from returning [constant NAN]. [codeblock] # s is 0.523599 or 30 degrees if converted with rad2deg(s) s = asin(0.5) diff --git a/modules/skeleton_3d/nodes/skeleton_ik.cpp b/modules/skeleton_3d/nodes/skeleton_ik.cpp index 2c8f81202..26f14f1e6 100644 --- a/modules/skeleton_3d/nodes/skeleton_ik.cpp +++ b/modules/skeleton_3d/nodes/skeleton_ik.cpp @@ -298,7 +298,8 @@ void FabrikInverseKinematic::solve(Task *p_task, real_t blending_delta, bool ove const Vector3 rot_axis(initial_ori.cross(ci->current_ori).normalized()); if (rot_axis[0] != 0 && rot_axis[1] != 0 && rot_axis[2] != 0) { - const real_t rot_angle(Math::acos(CLAMP(initial_ori.dot(ci->current_ori), -1, 1))); + // acos does clamping. + const real_t rot_angle(Math::acos(initial_ori.dot(ci->current_ori))); new_bone_pose.basis.rotate(rot_axis, rot_angle); } diff --git a/scene/3d/path.cpp b/scene/3d/path.cpp index 2ecd5c78c..004f7d520 100644 --- a/scene/3d/path.cpp +++ b/scene/3d/path.cpp @@ -176,7 +176,9 @@ void PathFollow::_update_transform(bool p_update_xyz_rot) { Vector3 axis = t_prev.cross(t_cur); float dot = t_prev.dot(t_cur); - float angle = Math::acos(CLAMP(dot, -1, 1)); + + // acos does clamping. + float angle = Math::acos(dot); if (likely(!Math::is_zero_approx(angle))) { if (rotation_mode == ROTATION_Y) { diff --git a/scene/resources/sky.cpp b/scene/resources/sky.cpp index a5d19552d..054dcdb62 100644 --- a/scene/resources/sky.cpp +++ b/scene/resources/sky.cpp @@ -173,7 +173,8 @@ Ref ProceduralSky::_generate_sky() { normal.normalize(); - float v_angle = Math::acos(CLAMP(normal.y, -1.0, 1.0)); + // acos does clamping. + float v_angle = Math::acos(normal.y); Color color; @@ -192,7 +193,8 @@ Ref ProceduralSky::_generate_sky() { color.g *= sky_energy; color.b *= sky_energy; - float sun_angle = Math::rad2deg(Math::acos(CLAMP(sun.dot(normal), -1.0, 1.0))); + // acos does clamping. + float sun_angle = Math::rad2deg(Math::acos(sun.dot(normal))); if (sun_angle < sun_angle_min) { color = color.blend(sun_linear);