3dcoordinate systemsrotationstrigonometry
I'm looking for the formula to calculate the cartesian coordinates of a point A(Xa, Ya, Za) after a rotation of an angle δ (in radians) around an axis B(Xb, Yb, Zb)C(Xc, Yc, Zc).
I found how to do it using rotation matrices or quaternions.
But how can I do that using the cartesian coordinates of the A, B and C points?
You can work with the rotation operators in $\mathbb{C}^2$. This method is well-known in quantum mechanics (rotation on the Bloch sphere).
I assume the following spherical coordinates system ($\theta$: colatitude, $\phi$: longitude):
Given the coordinate angles $(\theta, \phi)$, define a vector $Q(\theta, \phi)$ in $\mathbb{C}^2$ by $$ Q(\theta, \phi) = \begin{pmatrix} \cos \frac{\theta}{2} \\ e^{i\phi}\sin\frac{\theta}{2}\end{pmatrix}. $$ Conversely, given a vector $\begin{pmatrix} z_1 \\ z_2 \end{pmatrix} \in \mathbb{C}^2$, define the coordinate angles $(\theta, \phi)$ by $$ Q^{-1}\begin{pmatrix} z_1 \\ z_2 \end{pmatrix} = \left( 2\, \textrm{atan} \frac{\textrm{Mod}(z_2)}{\textrm{Mod}(z_1)}, \textrm{Arg}(z_2)-\textrm{Arg}(z_1)\right). $$
Then the rotation to the angle $\alpha$ around a Cartesian axis maps the coordinate angles $(\theta,\phi)$ to the coordinate angles $$ (\theta', \phi') = (Q^{-1} \circ R_{\textrm{axis}}(\alpha) \circ Q)(\theta,\phi) $$ where the rotation operator $R_{\textrm{axis}}(\alpha)$ is, depending on the rotation axis: \begin{gather} R_x(\alpha) = \begin{pmatrix} \cos\frac{\alpha}{2} & -i\sin\frac{\alpha}{2} \\ -i\sin\frac{\alpha}{2} & \cos\frac{\alpha}{2} \end{pmatrix}, \\ R_y(\alpha) = \begin{pmatrix} \cos\frac{\alpha}{2} & -\sin\frac{\alpha}{2} \\ \sin\frac{\alpha}{2} & \cos\frac{\alpha}{2} \end{pmatrix}, \\ R_z(\alpha) = \begin{pmatrix} e^{-i\frac{\alpha}{2}} & 0 \\ 0 & e^{i\frac{\alpha}{2}} \end{pmatrix}. \end{gather}
My interpretation of the problem statement is as follows (please correct me if I've misunderstood):
First, note that the vector from $A$ to $B$ can be expressed as $$\vec{AB} = B-A = (x'-x, \ y'-y, \ z'-z)$$
There are actually infinitely many different orientations that "point along" this vector $\vec{AB}$. To see why, imagine an aircraft doing an aileron roll while flying from point $A$ to point $B$. The aircraft's direction of motion is always pointed towards the destination point $B$, but its orientation is constantly changing throughout the maneuver.
So, we only need to find one possible orientation meeting the above criteria.
We can think of an orientation as a 3D rotation applied to the default coordinate axes. There are many distinct ways to represent rotations, as outlined in http://en.wikipedia.org/wiki/Rotation_formalisms_in_three_dimensions. (Conveniently, this article also defines conversions to and from the various formalisms.)
Let's work with the rotation matrix representation. First, let's say we want the unit vector $\hat{x} = (1,0,0)$ to map onto the unit vector $$\hat{u} = \frac{\vec{AB}}{|\vec{AB}|},$$
which means that $\hat{u}$ gives us the first column of the rotation matrix.
Since the transformed axes ($\hat{u}$, $\hat{v}$, and $\hat{w}$) must be mutually orthogonal, we know that the second column $\hat{v}$ of the rotation matrix must satisfy $\hat{u}\cdot\hat{v}=0.$ That is, $$v_1 u_1 + v_2 u_2 + v_3 u_3 = 0.$$ And, since $\hat{v}$ must be a unit vector, $$|\hat{v}|^2 = v_1^2+v_2^2+v_3^2 = 1.$$
With two equations governing three unknowns, we have one degree of freedom remaining. (As with the aileron-roll example, we've fixed the direction in which one transformed axis will point, but the coordinate system is still free to rotate around the transformed axis.)
Let's arbitrarily pick a value of $\hat{v}$ with zero $z$ component and a positive $y$ component; that is, let's set $v_3=0$, solve the (quadratic) system of equations, and pick the solution that yields $v_2 > 0$. (In rare cases where $\hat{u}$ lies exactly along the $y$ axis, forcing $v_2 = 0$, we can choose a slightly different, but still essentially arbitrary, constraint.)
This gives us two columns of the rotation matrix. We can solve for the final column, $\hat{w}$, by again using the unit-vector and orthogonality properties. The quadratic equations should have two solutions, but we'll pick the one unique solution that satisfies the same "right-hand rule" that our initial coordinate system satisfied (typically $\hat{u}\times\hat{v}=\hat{w}$).
Now, with a full rotation matrix, we can translate this into a roll, pitch, and yaw using the conversions listed in the "Rotation formalisms" article cited above.
Note: in order to get a unique answer, we had to introduce an arbitrary constraint along the way, when solving for $\hat{v}$. Other constraints might be equally valid, depending on the use case. For instance, one approach would be to specify that the roll $R$ should remain constant ($R = R'$).
Best Answer
$\newcommand{\Vec}[1]{\mathbf{#1}}$Here is an algorithm (not a formula per se, but can be turned into one):
Geometrically, $(\Vec{e}_{1}', \Vec{e}_{2}', \Vec{e}_{3})$ is the result of rotating the "original" orthonormal basis $(\Vec{e}_{1}, \Vec{e}_{2}, \Vec{e}_{3})$ through angle $\delta$ about $\Vec{e}_{3}$. The formula in item 4. decomposes $A - B$ (the displacement of point $A$ from the "origin" $B$) into components with respect to the original basis, and uses those as components with respect to the rotated basis.
Note: If $\Vec{e}_{3} = (0, 0, \pm 1)$ we can pick $\Vec{e}_{1} = (1, 0, 0)$ and $\Vec{e}_{2} = (0, 1, 0)$. Otherwise, we have $\Vec{e}_{3} = (X, Y, Z)$ with $X^{2} + Y^{2} \neq 0$, and may pick $$ \Vec{e}_{1} = \frac{(-Y, X, 0)}{\sqrt{X^{2} + Y^{2}}}. $$ This is not numerically robust if the segment $\overline{BC}$ is nearly parallel to the third Cartesian axis. In practice, find the two largest (in absolute value) components of $\Vec{e}_{3}$ and pick $\Vec{e}_{1}$ using those instead of the first two.