"Euler Angles" you can think of as a function $(S^1)^3 \to SO_3$ or $\mathbb R^3 \to SO_3$. The derivative of this function does not always have rank 3, so you have degenerate submanifolds where the function is many-to-one. In this special case that's called "gimbal lock".
One formalism that avoids this is quaternions. You can of course use other formalisms, and many other formalisms are naturally related to the quaternion version, so people tend to gravitate to the quaternion version. One version that's closely related to quaternions would be to use the exponential map for the unit quaternion group. But this also has "gymbal lock" but of a different kind. But it does have the rather appealing interpretation as rotations about arbitrary axis -- this is perhaps more useful if you're only interested in rotations that differ from the identity matrix (or some given matrix) by a small amount, they're very natural coordinates on "small scales" in $SO_3$.
Are there any special properties you'd like for coordinates on $SO_3$? That might give a sense for where you want to go with this.
Edit in response to 1st comment: There's various conventions for Euler angles. Let's use what Wikipedia calls "Proper Euler Angles" http://en.wikipedia.org/wiki/Euler_angles
Given a unit vector $\vec v \in \mathbb R^3$ and a number $\theta \in \mathbb R$, let $R_{\vec v, \theta}$ be the rotation matrix which fixes $\vec v$ and which rotates counter-clockwise by an angle $\theta$ in the plane orthogonal to $\vec v$. One version of "Proper Euler angles" would be:
$$R_{(1,0,0), \theta_1} R_{(0,1,0), \theta_2} R_{(1,0,0), \theta_3}$$
The "Tait-Bryan" variant would use the three distinct vectors, notice above we only used the standard vectors for the xy-plane. The Tait-Bryan variant is the one appropriate to pitch/roll/yaw of an aeroplane. Wikipedia has an excellent picture:
In this case you'd have
$$R_{(1,0,0),\theta_1} R_{(0,1,0), \theta_2} R_{(0,0,1),\theta_3} $$
$\theta_3$ would be the roll, $\theta_2$ the pitch and $\theta_1$ the yaw. To relate my coordinates to the picture, $(1,0,0)$ is the direction the plane is pointing. $(0,1,0)$ is the direction of the left wing. $(0,0,1)$ is the direction of the yellow axis sticking out of the top of the plane.
You can write out the matrices quite explicitly:
$$ R_{(1,0,0), \theta} =
\pmatrix{ 1 & 0 & 0 \cr
0 & \cos \theta & -\sin \theta \cr
0 & \sin \theta & \cos \theta }$$
So in this case, one occurance of "gimbal lock" is $\theta_2 = 0$. In the Tait-Bryan variant it would be when the aeroplane is either pointing straight up or down, which is $\theta_2 = \pm \pi/2$.
After a lot of searching, I have finally found a reference that has what I needed. Although it doesn't appear in a 'Quaternion to Euler' or 'Quaternion to Tait-Bryan' google search, when I gave up and started looking for 'Quaternion to Axis-Angle' with the intention of going through that representation as an intermediate step, I came across the following wonderful document:
Technical Concepts
Orientation, Rotation, Velocity and Acceleration,
and the SRM
Version 2.0, 20 June 2008
Author: Paul Berner, PhD
Contributors: Ralph Toms, PhD, Kevin Trott, Farid Mamaghani, David Shen,
Craig Rollins, Edward Powell, PhD
http://www.sedris.org/wg8home/Documents/WG80485.pdf
It covers a lot of the formalisms, but most importantly, shows derivations and solutions for 3-1-3 and 3-2-1 Euler angle representation. It also seems to cover inter-conversion between pretty much every other rotation representation I'm aware of, and so I would also recommend it as a good general reference.
Oh, and the actual solution for a 3-2-1 ($z-y-x$) Tait-Bryan rotation convention from that reference:
$$
\phi = \operatorname{arctan2}\left(q_2 q_3 + q_0 q_1,\frac{1}{2}-(q_1^2 + q_2^2)\right) \\
\theta = \arcsin(-2(q_1 q_3 - q_0 q_2)) \\
\psi = \operatorname{arctan2}\left(q_1 q_2 + q_0 q_3,\frac{1}{2}-(q_2^2 + q_3^2)\right)
$$
Note that the gimbal-lock situation occurs when $2(q_1 q_3 + q_0 q_2) = \pm1$ (which gives a $\theta$ of $\pm \frac{\pi}{2} $), so it can be clearly identified before you attempt to evaluate $\phi$ and $\psi$.
(Convention for arctan2 is $\operatorname{arctan2}(y, x)$, as hinted on page 3.)
Best Answer
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'$).