Snap/round 3D direction vector based on angle

I am trying to make a 3D direction adjustment algorithm that enters a unit vector and returns the closest (angle-wise) vector whose angle on each axis is a multiple of a given angle.

For example, if the given angle is 45º, I would snap the angle to 8 steps in each axis for a total of 26 possibilities. I hope that was clear enough. If it's not please let me know.

This is my current attempt:

Code (CSharp):

1. float3 SnapDirection(float3 input) {
2.  
3.    float3 angle = math.acos(input);
4.  
5.    float3 rounded = math.round(angle / radians(45)) * radians(45);
6.  
7.    float3 cosine = math.cos(rounded);
8.  
9.    return math.normalize(cosine);
10.  
11. }
12.  

What I expect it to output is the option with the smallest angle to the input. I measure it by math.acos(math.dot(input, output)). It mostly correct results but there are a good amount of cases where it doesn't get the best option.

This is because I'm rounding each axis independently, and so I get wrong results in the cases where the angle to one or 2 of the axes are longer than the optimal solution but the overall angle is shorter.

Consider the following case:

• Input => (0.3546994, 0.3626332, 0.861792)
• Current output => (0, 0, 1)
• Best output => (0.5773503, 0.5773503, 0.5773503)

Those vector correspond to the following angles:

• Input angles => ~(69.2º, 68.7º, 30.5º)
• Current output angles => (90º, 90º 45º)
• Best output angles => (45º, 45º, 45º)

And below you can see that the algorithm thinks the current result is the best because each number is smaller:

• Difference input-current => (20.8º, 21.3º, 14.5º)
• Difference input-best => (24.2º, 23.7º, 14.5º)

But the angle ends up being bigger:

• Actual angle from input to current => 30.5º
• Actual angle from input to best => 24.5º

Additionally, I need this code to be as performant as possible and avoid branching as I may need to port this to GPU later on.

Spoiler: Original post

I'm trying to snap a direction vector into a few axis-aligned directions, pointing to the vertices, edges and faces of an axis-aligned cube, or in other words, the components of the resulting vector can be either -1, 0 or 1. Hopefully this picture can make it clearer:



I have a rudimentary knowledge of trigonometry but it's never a strength of mine. I've tried a few naive ideas but none of them get correct results.

Can someone point me in the right direction?
Thanks!

 

(Reposting because I couldn't edit my post)

Use the dot product (https://en.wikipedia.org/wiki/Dot_product). When you get the dot product of two unit vectors (vectors whose magnitudes are equal to 1), you get the cosine of the angle between them.

The idea would be would be that you would get the dot product of your direction with each of the 3 main axes (forward, right, up), and go with the axis that gave you the largest dot product (smallest angle between the two vectors).

You don't need to do their negative counterparts if you just use the absolute values of the dot products for your conditions.

Something like this:

Code (CSharp):

1.  
2. Vector3[] directions=new Vector3[]{Vector3.forward, Vector3.right, Vector3.up};
3. int closest axis;
4. float largest_dot;
5. float dot;
6.  
7. for(int axis = 0 ; axis < 3; axis++){
8.   dot=Vector3.dot(direction, directions[axis])
9.   if (abs(dot) > abs(largest_dot){
10.     largest_dot=dot;
11.     closest axis=axis;
12.   }
13. }
14.  
15. Vector3 new_direction=Vector3.zero
16.  
17. new_direction[closest axis]=Mathf.sign(direction[closest axis]);

You could also achieve it without using the dot product by just finding the axis in your direction with the largest value, but dot products are a useful thing to know for Vector math.

 

Hi!

Thank you for answering. This was in fact one of the first things I tried. But you made me realize that actually the vector components already represent the cosine of the angle with each axis, so the dot product is not needed. (the dot operations return the exact same numbers)

Actually an older try was really close but I was using the sine instead of cosine, so with this fix now works perfectly!

Code (CSharp):

1. const float GridAngle = radians(45);
2.  
3.     public static float3 EncodeDirection(float3 dir)
4.     {
5.         float3 angle = acos(dir);
6.         float3 rounded = round(angle / GridAngle) * GridAngle;
7.         float3 vector = normalize(cos(rounded));
8.  
9.         return vector;
10.     }

Now I have to profile and optimize it. There's probably a lot of room for that.

Thank you!

 

Well I spoke too fast. There's some rounding errors with the code above

In the following image you can see a case where the input (white) is rounded to a direction further away (green) than the optimal (red). The optimal is found by looking up a table and calculating the minimum angle with all the possibilities.

Spoiler

I initially thought the error was because I was rounding raw cosines instead of linear values, but I'm comparing the angle (acos).

Any Ideas?

 

How about this:

Code (CSharp):

1.  
2. float MyRound(float f)
3. {
4.    if (Mathf.Abs(f) < Mathf.Cos(Mathf.PI / 8))
5.       return 0;
6.    else
7.       return Mathf.Sign(f);
8. }
9.  
10.  
11. Vector3 n = dir.Normalize();
12. return new Vector3(MyRound(n.x), MyRound(n.y), MyRound(n.z));
13.  

 

Thanks for answering!
There are 2 problems to overcome here:
- First, rounding (middle 0.5f) using the raw cosine has a big error. Alternatively, calculating the arcosine and then the cos is costly. You solution solved this! With some little changes I'm getting the same results as my previous algorithm for far cheaper!
- The second issue (rounding) comes from the fact that I'm rounding each axis independently, and so I get wrong results in the cases where the angle for an axis is longer but the overall angle is shorter. Your answer suffers from the same issue unfortunately.

 

I'm intrigued. What are those changes? I thought my solution was exactly what you were asking for, so where did I misunderstand you?

I think I understand what you mean but do not know, yet, why that would occur. Can you give me a concrete example input and desired output where the approach fails?

 

What I changed was the calculation of the midpoint and the API (I'm using the new mathematics one). Pi / 8 was snapping the angle to 90º so after some tinkering I arrived to this:

cosMidpoint = cos(radians(GridAngle * 1.5f));

And I computed it as a const to avoid the perf cost.

Here's an example of the error (ignore the degree/radian details, I've made everything in degrees for easier following)

Spoiler: Example

What I expect it to output is the option with the smallest angle to the input. The angle is defined as

math.acos(math.dot(input, output))

.

Consider the following case:

• Input => (0.3546994, 0.3626332, 0.861792)
• Current output => (0, 0, 1)
• Best output => (0.5773503, 0.5773503, 0.5773503)

Those vector correspond to the following angles:

• Input angles => ~(69.2º, 68.7º, 30.5º)
• Current output angles => (90º, 90º 45º)
• Best output angles => (45º, 45º, 45º)

And below you can see that the algorithm rounds each component correctly and thinks the current result is the best because the individual differences are smaller:

• Difference input-current => (20.8º, 21.3º, 14.5º)
• Difference input-best => (24.2º, 23.7º, 14.5º)

But the angle ends up being bigger:

• Actual angle from input to current => 30.5º
• Actual angle from input to best => 24.5º

 

Thanks, I now have a theory where I misunderstood you, though I have not verified it yet. You initially wanted to snap to vectors pointing to the surface of a cube. So I assumed you wanted to minimize the distance between the point on the surface given by the (possibly scaled) input vector and the 'snap' point. But you now want to minimize angles which is similar to distances on a sphere. That is a different thing altogether and might explain the error. I'll have to compute the distance on the cube surface for your example to harden the theory... Until I find time for this, can I ask you which of the two things you want, cube or sphere?

 

Yes I realize I've been unclear. My original use case is related to cubes but I should've formulated the problem better. I want to snap angles / points on a sphere.

I've now edited my post with a better explanation of the problem, my current solution and an example of the error I need to solve. I hope it will be clearer for anyone reading this from now on.

Sorry for the confusion.

 

I am not sure if you'll get anywhere with this generalized problem. The requirement of fixing the angle on all of the three axes in general won't give you any solutions. Each angle that you fix on an axis defines a plane. Doing this on the three axes will give you three planes that need to intersect to give you a valid point. In general, three planes only intersect in a single point and for our setup we know that this is the origin. So only the cases where we get valid solutions for this specification is when the line we get by intersecting two of the planes also happens to lie in the third plane. (We then intersect with the unit sphere to get the valid point.) I think this case might only occur with angles of 45º and 90º, but I am not sure. So we might as well solve this only for 45º degrees... or you actually want something else.

 

So, I've finally gotten around to set up a test for this. My initial suggestion for the cube case was just broken. But the error was not really in the value for the angle but in the use of cos where tan would have been the correct choice. And when using tan I also need to actually position the point on the cube. I've created a visualization by generating random points on the unit sphere and coloring all points that end up snapping to the same vector with the same color.

This is the updated code:

Code (CSharp):

1.     private Vector3 findNearest(Vector3 dir)
2.     {
3.         // scale vector to unit cube
4.         float scaleDivisor = Mathf.Abs(dir.x);
5.         if (Mathf.Abs(dir.y) > scaleDivisor)
6.             scaleDivisor = Mathf.Abs(dir.y);
7.         if (Mathf.Abs(dir.z) > scaleDivisor)
8.             scaleDivisor = Mathf.Abs(dir.z);
9.         dir = dir / scaleDivisor;
10.  
11.         return new Vector3(MyRound(dir.x), MyRound(dir.y), MyRound(dir.z));
12.     }
13.  
14.     float MyRound(float f)
15.     {
16.         if (Mathf.Abs(f) < Mathf.Tan(Mathf.PI / 8))
17.             return 0;
18.         else
19.             return Mathf.Sign(f);
20.     }

And this is how it looks:


This is still not the case you are interested in but it's what my original code was supposed to do.

 

Wow! Gonna take a look at that as soon as I can. Also I will upload my own code an test/visualization to avoid any further misunderstanding

 

Well, that "as soon as I can" surely took linger than I wanted. Work is draining all my energy

Anyways, I tried your algorithm and it does what I want, and with lower error too! It's great!

As promised, here's a test script comparing both your and my algorithms, and a few fancy gizmos to visualize the whole thing: https://pastebin.com/vNNWMQyL

My algorithm is plotted green, yours is magenta and the "best possible result" is in green. I get this from searching a lookup table. It's entirely possible that I'm mistaken in my heuristic and that's actually not the best possible result.

In any case, If you're still interested take a look and let me know (btw the algorithm depends on the mathematics package and was made using unity 2020.1 beta, but should work in older versions just fine).

Again thanks a lot for taking the time for answering!