diff --git a/intern/cycles/kernel/shaders/node_sky_texture.osl b/intern/cycles/kernel/shaders/node_sky_texture.osl index 08bc8f85120..20d379939ab 100644 --- a/intern/cycles/kernel/shaders/node_sky_texture.osl +++ b/intern/cycles/kernel/shaders/node_sky_texture.osl @@ -122,6 +122,11 @@ vector geographical_to_direction(float lat, float lon) return vector(cos(lat) * cos(lon), cos(lat) * sin(lon), sin(lat)); } +float precise_angle(vector a, vector b) +{ + return 2.0 * atan2(length(a - b), length(a + b)); +} + color sky_radiance_nishita(vector dir, float nishita_data[9], string filename) { /* definitions */ @@ -138,7 +143,7 @@ color sky_radiance_nishita(vector dir, float nishita_data[9], string filename) if (dir[2] >= 0.0) { /* definitions */ vector sun_dir = geographical_to_direction(sun_elevation, sun_rotation + M_PI_2); - float sun_dir_angle = acos(dot(dir, sun_dir)); + float sun_dir_angle = precise_angle(dir, sun_dir); float half_angular = angular_diameter / 2.0; float dir_elevation = M_PI_2 - direction[0]; diff --git a/intern/cycles/kernel/svm/svm_sky.h b/intern/cycles/kernel/svm/svm_sky.h index e877bd9a5c8..214c0cd1a9a 100644 --- a/intern/cycles/kernel/svm/svm_sky.h +++ b/intern/cycles/kernel/svm/svm_sky.h @@ -145,7 +145,7 @@ ccl_device float3 sky_radiance_nishita(KernelGlobals *kg, if (dir.z >= 0.0f) { /* definitions */ float3 sun_dir = geographical_to_direction(sun_elevation, sun_rotation + M_PI_2_F); - float sun_dir_angle = acos(dot(dir, sun_dir)); + float sun_dir_angle = precise_angle(dir, sun_dir); float half_angular = angular_diameter / 2.0f; float dir_elevation = M_PI_2_F - direction.x; diff --git a/intern/cycles/util/util_math.h b/intern/cycles/util/util_math.h index 737c834e073..8caabf6eac3 100644 --- a/intern/cycles/util/util_math.h +++ b/intern/cycles/util/util_math.h @@ -787,6 +787,16 @@ ccl_device_inline float compare_floats(float a, float b, float abs_diff, int ulp return (abs(__float_as_int(a) - __float_as_int(b)) < ulp_diff); } +/* Calculate the angle between the two vectors a and b. + * The usual approach acos(dot(a, b)) has severe precision issues for small angles, + * which are avoided by this method. + * Based on "Mangled Angles" from https://people.eecs.berkeley.edu/~wkahan/Mindless.pdf + */ +ccl_device_inline float precise_angle(float3 a, float3 b) +{ + return 2.0f * atan2f(len(a - b), len(a + b)); +} + CCL_NAMESPACE_END #endif /* __UTIL_MATH_H__ */