mesh: use faster acos() variant in normals calculation #114501

Merged
Aras Pranckevicius merged 7 commits from aras_p/blender:acos_approx into main 2023-11-07 18:22:25 +01:00
6 changed files with 40 additions and 15 deletions

View File

@ -23,6 +23,7 @@
#include "BLI_array_utils.hh"
#include "BLI_bit_vector.hh"
#include "BLI_linklist.h"
#include "BLI_math_base.hh"
#include "BLI_math_vector.hh"
#include "BLI_memarena.h"
#include "BLI_span.hh"
@ -243,7 +244,7 @@ static void accumulate_face_normal_to_vert(const Span<float3> positions,
}
/* Calculate angle between the two face edges incident on this vertex. */
const float fac = saacos(-dot_v3v3(edvec_prev, edvec_next));
const float fac = math::safe_acos_approx(-dot_v3v3(edvec_prev, edvec_next));
const float vnor_add[3] = {face_normal[0] * fac, face_normal[1] * fac, face_normal[2] * fac};
float *vnor = vert_normals[face_verts[i_curr]];
@ -535,7 +536,7 @@ static CornerNormalSpace lnor_space_define(const float lnor[3],
if (!edge_vectors.is_empty()) {
float alpha = 0.0f;
for (const float3 &vec : edge_vectors) {
alpha += saacosf(dot_v3v3(vec, lnor));
alpha += math::safe_acos_approx(dot_v3v3(vec, lnor));
}
/* This piece of code shall only be called for more than one loop. */
/* NOTE: In theory, this could be `count > 2`,
@ -546,8 +547,8 @@ static CornerNormalSpace lnor_space_define(const float lnor[3],
lnor_space.ref_alpha = alpha / float(edge_vectors.size());
}
else {
lnor_space.ref_alpha = (saacosf(dot_v3v3(vec_ref, lnor)) +
saacosf(dot_v3v3(vec_other, lnor))) /
lnor_space.ref_alpha = (math::safe_acos_approx(dot_v3v3(vec_ref, lnor)) +
math::safe_acos_approx(dot_v3v3(vec_other, lnor))) /
2.0f;
}
@ -567,7 +568,7 @@ static CornerNormalSpace lnor_space_define(const float lnor[3],
/* Beta is angle between ref_vec and other_vec, around lnor. */
dtp = dot_v3v3(lnor_space.vec_ref, vec_other);
if (LIKELY(dtp < LNOR_SPACE_TRIGO_THRESHOLD)) {
const float beta = saacos(dtp);
const float beta = math::safe_acos_approx(dtp);
lnor_space.ref_beta = (dot_v3v3(lnor_space.vec_ortho, vec_other) < 0.0f) ? pi2 - beta : beta;
}
else {
@ -698,7 +699,7 @@ short2 lnor_space_custom_normal_to_data(const CornerNormalSpace &lnor_space,
float vec[3], cos_beta;
float alpha;
alpha = saacosf(cos_alpha);
alpha = math::safe_acos_approx(cos_alpha);
if (alpha > lnor_space.ref_alpha) {
/* Note we could stick to [0, pi] range here,
* but makes decoding more complex, not worth it. */
@ -716,7 +717,7 @@ short2 lnor_space_custom_normal_to_data(const CornerNormalSpace &lnor_space,
cos_beta = dot_v3v3(lnor_space.vec_ref, vec);
if (cos_beta < LNOR_SPACE_TRIGO_THRESHOLD) {
float beta = saacosf(cos_beta);
float beta = math::safe_acos_approx(cos_beta);
if (dot_v3v3(lnor_space.vec_ortho, vec) < 0.0f) {
beta = pi2 - beta;
}
@ -1097,7 +1098,8 @@ static void split_loop_nor_fan_do(LoopSplitTaskDataCommon *common_data,
/* Code similar to accumulate_vertex_normals_poly_v3. */
/* Calculate angle between the two face edges incident on this vertex. */
lnor += face_normals[loop_to_face[mlfan_curr_index]] * saacos(math::dot(vec_curr, vec_prev));
lnor += face_normals[loop_to_face[mlfan_curr_index]] *
math::safe_acos_approx(math::dot(vec_curr, vec_prev));
processed_corners.append(mlfan_vert_index);

View File

@ -192,6 +192,24 @@ template<typename T> inline T safe_acos(const T &a)
return math::acos((a));
}
/* Faster/approximate version of acosf. Max error 4.51803e-5 (0.00258 degrees).*/
inline float safe_acos_approx(float x)
{
const float f = fabsf(x);
/* clamp and crush denormals. */
const float m = (f < 1.0f) ? 1.0f - (1.0f - f) : 1.0f;
/* Based on http://www.pouet.net/topic.php?which=9132&page=2
* 85% accurate (ULP 0)
* Examined 2130706434 values of acos:
* 15.2000597 avg ULP diff, 4492 max ULP, 4.51803e-05 max error // without "denormal crush"
* Examined 2130706434 values of acos:
* 15.2007108 avg ULP diff, 4492 max ULP, 4.51803e-05 max error // with "denormal crush"
*/
const float a = sqrtf(1.0f - m) *
(1.5707963267f + m * (-0.213300989f + m * (0.077980478f + m * -0.02164095f)));
return x < 0 ? (float)M_PI - a : a;
}
template<typename T> inline T asin(const T &a)
{
return std::asin(a);

View File

@ -8,6 +8,7 @@
#include "BLI_array.hh"
#include "BLI_math_base.h"
#include "BLI_math_base.hh"
#include "BLI_math_geom.h"
#include "BLI_math_bits.h"
@ -5015,7 +5016,7 @@ void accumulate_vertex_normals_tri_v3(float n1[3],
for (i = 0; i < nverts; i++) {
const float *cur_edge = vdiffs[i];
const float fac = saacos(-dot_v3v3(cur_edge, prev_edge));
const float fac = blender::math::safe_acos_approx(-dot_v3v3(cur_edge, prev_edge));
/* accumulate */
madd_v3_v3fl(vn[i], f_no, fac);
@ -5062,7 +5063,7 @@ void accumulate_vertex_normals_v3(float n1[3],
for (i = 0; i < nverts; i++) {
const float *cur_edge = vdiffs[i];
const float fac = saacos(-dot_v3v3(cur_edge, prev_edge));
const float fac = blender::math::safe_acos_approx(-dot_v3v3(cur_edge, prev_edge));
/* accumulate */
madd_v3_v3fl(vn[i], f_no, fac);
@ -5094,7 +5095,7 @@ void accumulate_vertex_normals_poly_v3(float **vertnos,
/* calculate angle between the two poly edges incident on
* this vertex */
const float fac = saacos(-dot_v3v3(cur_edge, prev_edge));
const float fac = blender::math::safe_acos_approx(-dot_v3v3(cur_edge, prev_edge));
/* accumulate */
madd_v3_v3fl(vertnos[i], polyno, fac);

View File

@ -16,6 +16,7 @@
#include "BLI_bitmap.h"
#include "BLI_linklist_stack.h"
#include "BLI_math_base.hh"
#include "BLI_math_vector.h"
#include "BLI_task.h"
#include "BLI_utildefines.h"
@ -72,7 +73,7 @@ BLI_INLINE void bm_vert_calc_normals_accum_loop(const BMLoop *l_iter,
if ((l_iter->prev->e->v1 == l_iter->prev->v) ^ (l_iter->e->v1 == l_iter->v)) {
dotprod = -dotprod;
}
const float fac = saacos(-dotprod);
const float fac = blender::math::safe_acos_approx(-dotprod);
/* Shouldn't happen as normalizing edge-vectors cause degenerate values to be zeroed out. */
BLI_assert(!isnan(fac));
madd_v3_v3fl(v_no, f_no, fac);
@ -642,7 +643,7 @@ static int bm_mesh_loops_calc_normals_for_loop(BMesh *bm,
/* Code similar to accumulate_vertex_normals_poly_v3. */
/* Calculate angle between the two face edges incident on this vertex. */
const BMFace *f = lfan_pivot->f;
const float fac = saacos(dot_v3v3(vec_next, vec_curr));
const float fac = blender::math::safe_acos_approx(dot_v3v3(vec_next, vec_curr));
const float *no = fnos ? fnos[BM_elem_index_get(f)] : f->no;
/* Accumulate */
madd_v3_v3fl(lnor, no, fac);

View File

@ -18,6 +18,7 @@
#include "BLI_alloca.h"
#include "BLI_heap.h"
#include "BLI_linklist.h"
#include "BLI_math_base.hh"
#include "BLI_math_geom.h"
#include "BLI_math_matrix.h"
#include "BLI_math_vector.h"
@ -615,7 +616,7 @@ static void bm_loop_normal_accum(const BMLoop *l, float no[3])
normalize_v3(vec1);
normalize_v3(vec2);
fac = saacos(-dot_v3v3(vec1, vec2));
fac = blender::math::safe_acos_approx(-dot_v3v3(vec1, vec2));
madd_v3_v3fl(no, l->f->no, fac);
}

View File

@ -8,6 +8,7 @@
* Method of smoothing deformation, also known as 'delta-mush'.
*/
#include "BLI_math_base.hh"
#include "BLI_math_matrix.h"
#include "BLI_math_vector.h"
#include "BLI_utildefines.h"
@ -470,7 +471,8 @@ static void calc_tangent_spaces(const Mesh *mesh,
if (calc_tangent_loop(v_dir_prev, v_dir_next, ts)) {
if (r_tangent_weights != nullptr) {
const float weight = fabsf(acosf(dot_v3v3(v_dir_next, v_dir_prev)));
const float weight = fabsf(
blender::math::safe_acos_approx(dot_v3v3(v_dir_next, v_dir_prev)));
r_tangent_weights[curr_corner] = weight;
r_tangent_weights_per_vertex[corner_verts[curr_corner]] += weight;
}