Reliable atan (y, x) on GLSL to convert XY coordinate to angle

Question

There are two versions of atan() in GLSL (specifically 3.00, which I use): atan(y_over_x) can only return angles between -PI / 2, PI / 2, while atan(y/x) can take into account all 4 quadrants, so the range of angles covers everything: -PI, PI, which is very similar to atan2() in C ++.

I would like to use the second atan to convert the XY coordinates to an angle. However, atan() in GLSL, also unable to handle when x = 0 , is not very stable. Especially where x is close to zero, the division can overflow, which leads to the opposite resulting angle (you get something close to -PI / 2, where you expect to get approximately PI / 2).

What a nice, simple implementation that we can build on top of the GLSL atan(y,x) to make it more reliable?

Answer 1

I am going to answer my question to share my knowledge. First of all, we note that instability occurs when x is near zero. However, we can also translate this as abs(x) << abs(y) . So, first we divide the plane (assuming that we are on the unit circle) into two areas: one, where |x| <= |y| |x| <= |y| and another, where |x| > |y| |x| > |y| as below:

We know that atan(x,y) is much more stable in the green area - when x is close to zero, we just have something close to atan (0.0), which is very stable numerically, while the usual atan(y,x) more stable in the orange area. You can also convince yourself that this relationship:

 atan(x,y) = PI/2 - atan(y,x) 

is executed for all non-origin (x, y), where it is undefined, and we are talking about atan(y,x) , which can return the value of the angle in the entire range -PI, PI, and not atan(y_over_x) which only returns the angle between -PI / 2, PI / 2. Therefore, our robust atan2() procedure for GLSL is pretty simple:

 float atan2(in float y, in float x) { bool s = (abs(x) > abs(y)); return mix(PI/2.0 - atan(x,y), atan(y,x), s); } 

As a side note, the identity for the mathematical function atan(x) is actually:

 atan(x) + atan(1/x) = sgn(x) * PI/2 

which is true because its range is (-PI / 2, PI / 2).

Answer 2

Your proposed solution still does not work in case x=y=0 . Here, both atan() functions return NaN.

Further, I would not rely on the mix to switch between the two cases. I'm not sure how this is implemented / compiled, but the IEEE float rules for x * NaN and x + NaN are again cast to NaN. Therefore, if your compiler really used mix / interpolation, the result should be NaN for x=0 or y=0 .

Here is another solution that solved the problem for me:

 float atan2(in float y, in float x) { return x == 0.0 ? sign(y)*PI/2 : atan(y, x); } 

For x=0 angle can be ± π / 2. Which of the two depends only on y . If y=0 too, the angle can be arbitrary (the vector has a length of 0). sign(y) returns 0 in this case, this is normal.

Answer 3

Depending on your target platform, this may be a problem. The OpenGL specification for atan (y, x) indicates that it should work in all quadrants, leaving undefined behavior only when x and y are 0.

So, it would be expected that any decent implementation would be stable near all axes, as this is the whole goal for the atan (or atan2 ) 2 argument.

The responder / responder is correct in that some implementations make keyboard shortcuts. However, the decision made suggests that a poor implementation will always be unstable when x is close to zero: on some hardware (for example, my Galaxy S4), the value is stable when x is close to zero, but unstable when y is about zero .

To test the GLSL rendering implementation of atan(y,x) , here's a WebGL test pattern. Follow the link below and as long as your OpenGL implementation is decent, you should see something like this:

Test pattern using native atan(y,x) : http://glslsandbox.com/e#26563.2

If all is well, you should see 8 different colors (ignoring the center).

Related atan(y,x) demos for multiple x and y values, including 0, are very large and very small. The central block atan(0.,0.) Is mathematically undefined, and the implementations differ. I saw 0 (red), PI / 2 (green), and NaN (black) on the equipment I tested.

Here is a test page for the decision made. Note: there is no mix(float,float,bool) in the WebGL version of the host, so I added an implementation that meets the specification.

Test pattern using atan2(y,x) from the accepted answer: http://glslsandbox.com/e#26666.0

Answer 4

Sometimes the best way to improve the performance of a piece of code is not to name it first. For example, one of the reasons you can determine the angle of a vector is that you can use this angle to construct a rotation matrix using combinations of the sine and cosine angles. However, the sine and cosine of the vector (relative to the origin) are already hidden in the forward direction inside the vector itself. All you have to do is create a normalized version of the vector by dividing each vector coordinate by the total length of the vector. Here is a two-dimensional example of calculating the sine and cosine of the angle of the vector [xy]:

 double length = sqrt(x*x + y*y); double cos = x / length; double sin = y / length; 

After you have the sine and cosine values, you can now directly fill the rotation matrix with these values ​​to rotate the arbitrary vectors clockwise or counterclockwise by the same angle, or you can combine the second rotation matrix to rotate it non-zero angle. In this case, you can imagine the rotation matrix as a “normalizing” angle to zero for an arbitrary vector. This approach can be expanded in the three-dimensional (or N-dimensional) case, although, for example, you will have three angles and six sin / cos pairs to calculate (one angle per plane) for 3D rotation.

In situations where you can use this approach, you get a big gain by completely bypassing the atan calculations, which is possible because the only reason you wanted to determine the angle was to calculate the sine and cosine values. By skipping the conversion to angular space and vice versa, you not only avoid worrying about division by zero, but also improve the accuracy of angles that are near the poles, otherwise they would suffer from multiplication / division by large numbers. I have successfully used this approach in the GLSL program, which rotates the scene to zero degrees to simplify the calculation.

It's easy to get into an immediate problem so that you can lose sight of why you need this information in the first place. Not that it works in every case, but sometimes it helps to think out of the box ...

Answer 5

A formula that gives an angle in four quadrants for any value

x and y coordinates. For x = y = 0, the result is not defined.

f (x, y) = pi () -pi () / 2 * (1 + sign (x)) * (1-sign (y ^ 2)) -pi () / 4 * (2 + sign (x) ) * sign (y)

  -sign(x*y)*atan((abs(x)-abs(y))/(abs(x)+abs(y)))