glsl: Rewrite atan2 implementation to fix accuracy and handling of zero/infinity.
authorFrancisco Jerez <currojerez@riseup.net>
Sat, 21 Jan 2017 21:41:08 +0000 (13:41 -0800)
committerFrancisco Jerez <currojerez@riseup.net>
Tue, 31 Jan 2017 18:32:45 +0000 (10:32 -0800)
commite9ffd12827ac11a2d2002a42fa8eb1df847153ba
tree51670ec05f1d7d1f4903dba0524a7c2522998d70
parent69042a5be4664c7928a21bd23e4f6795bfb19f60
glsl: Rewrite atan2 implementation to fix accuracy and handling of zero/infinity.

This addresses several issues of the current atan2 implementation:

 - Negative zero (and negative denorms which end up getting flushed to
   zero) isn't handled correctly by the current implementation.  The
   reason is that it does 'y >= 0' and 'x < 0' comparisons to decide
   on which side of the branch cut the argument is, which causes us to
   return incorrect results (off by up to 2π) for very small negative
   values.

 - There is a serious precision problem for x values of large enough
   magnitude introduced by the floating point division operation being
   implemented as a mul+rcp sequence.  This can lead to the quotient
   getting flushed to zero in some cases introducing an error of over
   8e6 ULP in the result -- Or in the most catastrophic case will
   cause us to return NaN instead of the correct value ±π/2 for y=±∞
   and x very large.  We can fix this easily by scaling down both
   arguments when the absolute value of the denominator goes above
   certain threshold.  The error of this atan2 implementation remains
   below 25 ULP in most of its domain except for a neighborhood of y=0
   where it reaches a maximum error of about 180 ULP.

 - It emits a bunch of instructions including no less than three
   if-else branches per scalar component that don't seem to get
   optimized out later on.  This implementation uses about 13% less
   instructions on Intel SKL hardware and doesn't emit any control
   flow instructions.

v2: Fix up argument scaling to take into account the range and
    precision of exotic FP24 hardware.  Flip coordinate system for
    arguments along the vertical line as if they were on the left
    half-plane in order to avoid division by zero which may give
    unspecified results on non-GLSL 4.1-capable hardware.  Sprinkle in
    some more comments.

Reviewed-by: Ian Romanick <ian.d.romanick@intel.com>
src/compiler/glsl/builtin_functions.cpp