Quaternion code needs clarification

[CoChrists]

The C++ (and also the CUDA) code that generates quaternions from three floats needs comments. Could anybody comment on the role of the parameters r1, r2, r3 in https://github.com/tjlane/thor/blob/master/src/scatter/_cpuscatter.cpp#L22 please?

What is the meaning of these three numbers in the context of rotations? Are these Euler angles? Or how can I control the rotation that the resulting quaternion represents?

Thanks a lot!

[tjlane]

This is a good point, punctuated by the fact that I can't remember if I wrote that code or stole it from somewhere... I remember basically what it does, and it looks like you should be able to modify it to do precise (in addition to random) rotations if you'd like.

I suggest we work in the new branch ("scatterdoc") I referenced in the PR above until you're happy and then merge it into master.

With regards to "controlling" the rotation, I think your answer is here: http://en.wikipedia.org/wiki/Rotation_formalisms_in_three_dimensions#Quaternions

If you want & feel capable, feel free to modify the code to suit your needs without breaking anything (should be possible here). Don't hesitate to hit me up if you need help or want to talk through a plan. Not 100% sure what you're trying to do, so if you need more help just elaborate a bit on how I can assist!

[CoChrists]

Thanks for writing some comments in the code. However, the comments and the references you gave do not explain the role of the three random float parameters, The wiki reference is not even related to your implementation, e.g. the scaler element is saved in q1 if I can trust your internal variable naming, i.e. q=(w,xi,yj,zk).

Remaining questions: 1) What's the role of r1,r2,r3? 2) Would I have to rewrite the implementation if I wanted to do "precise" rotation? Or could I keep the function and just feed it some precise r1,r2,r3 in the way that the answer to question 1) says?

[tjlane]

OK, so first, I'm pretty sure I didn't write this code, so I was hedging a little (ps @dermen did you write this function? If not I probably stole it from somewhere).

That said by cross-referencing the wiki article, it is possible to figure out what is going on, and understand how r1/r2/r3 are used to generate a random rotation.

The wiki article shows how a quaterion embeds an Euler axis and rotation theta around that axis.

The function defines two random angles (theta1, theta2) and also two random directions (sig1, sig2). It uses these to compute a random Euler axis and rotation about that axis. The code shows how r1/r2/r3 get mapped onto the axis and angle.

Depending on what you are doing, it either makes sense to modify this function or write your own version. B/c the code is so small, the second option is not necessarily a bad one. That's why I asked exactly what you're trying to do so I could possibly help guide things forward.

[CoChrists]

In your code comments you write: "that quaternion represents a random rotation draw uniformly from SO(3)". It would be nice to have a reference that proofs this.

Further more you say: "[The function] uses these to compute a random Euler axis and rotation about that axis." Could you point me to the lines that compute the Euler axis and the rotation angle, please? I cannot see them... Maybe you meant that implicitly these are defined by the code? Genuinely, I want to understand how this works.

[CoChrists]

To give some background information: What I want to do eventually is to simulate scattering from a fixed orientation. Of course, I would like to modify the exiting code/interfaces as little as necessary only. But in order to do so I need to understand this function. Otherwise, i could write a wrapper, that would either call this function or myGenerateFixedQuaternion()...

[tjlane]

So if you force all the random numbers (r1, r2, r3) to be the same for every simulation, you'll get the same rotation every time, eliminating the randomness. If you pass (0,0,0), you should get no rotation.

I acknowledge that you want to understand how this algorithm works, but I think a description of quaternion algebra is beyond the scope of the documentation here.

The code is pretty explicit IMO. This line: https://github.com/tjlane/thor/blob/master/src/scatter/_cpuscatter.cpp#L34

Corresponds directly to the formulae for q{1,2,3,4} in the wiki article: http://en.wikipedia.org/wiki/Rotation_formalisms_in_three_dimensions#Quaternions

I can add that as a comment? I'm a little unclear of what would satisfy you here, as things are clear to me. Feel free to modify the code/add documentation as you see fit.

[rmcgibbo]

Not to get too involved in something I'm not a part of, but I think a comment -- i.e. what you would put in a python docstring for the same function -- would go a long way here. It's not to explain how quaterion algebra works, its just to explain what the semantics of the input and output arguments are.

[CoChrists]

@tjlane The wiki article does not have a theta2 nor does it have a sqrt(1-s). The formula in the wiki is: w = cos(theta / 2), which appears nowhere in the code...

@rmcgibbo I totally agree with you.

But it is unclear which properties r1,r2,r3 must fulfill in order for this function to result in what it is supposed to, namely a random rotation on SO(3). Ideally I would like to see a reference here, which explains the algorithm. There one could look up the details. That's all. The wiki article and the 3dprogramming just explain basic Quaternion definitions. That's not needed as I can look up general definitions myself. The reference we need should explain how to create a random quaternion rotation in SO(3) given three independent random numbers uniformly distributed over [0,1).

[rmcgibbo]

Here is some python code for random quaterions that I ported from yank, I think: https://github.com/rmcgibbo/mdtraj/blob/master/MDTraj/utils/rotation.py#L90

From comparing the two, I am pretty sure r1, r2 and r3 are supposed to be uniformly distributed in [0,1]. This page also has some info: http://planning.cs.uiuc.edu/node198.html

[CoChrists]

"This page also has some info: http://planning.cs.uiuc.edu/node198.html"

This is great! Exactly what I was looking for and it even has reference to a book! That you very much RMCGIBBO.

[CoChrists]

@tjlane Could you include the link from rmcgibbo in the code comments, please?

[tjlane]

@rmcgibbo thanks for the link. fwiw I did include a doc string in the PR linked above.

@CoChrists I'll commit the link in a second, but part of my badgering you a bit is to try and share the burden of maintaining this code base, so that it doesn't always fall on me. So in the future feel free to just go ahead and contribute when and where the need is clear.

[CoChrists]

@tjlane Thank you! but really no need for the two other irrelevant citations.