Determining the torque needed to rotate a spacecraft to a given rotation/quaternion Which is the correct algorithm to determine the torque needed to rotate a spacecraft to a given quaternion ?
I have a set of quaternions of a spacecraft and the time difference (delta T) between the quaternions.
I've tried estimating the angular velocity and angular acceleration as: 
1. Find quaternion q so q = q1/q0


*Convert quaternion q to axis and angle:


len=sqrt(q.xq.x+q.yq.y+q.z*q.z);
angle=2*atan2(len, q.w);
axis = q.xyz()/len; 


*angular velocity w = axis * angle / dt. 

*angular acceleration a = w/dt. 


Using this algorithm, the computed angular acceleration (a) is not correct. 
What I'm doing wrong? Is this the correct way to obtain the angular acceleration and torque from a sequence of two quaternions?
 A: From your timeseries of the rotation quaternion $\gamma$, estimate the following derivatives:
$$\Omega = \gamma^{-1} \,\dot{\gamma}$$
$$\dot{\Omega} = \gamma^{-1}\, \ddot{\gamma}-\Omega^2$$
with all the wonted provisos and warnings about calculating numerical derivatives from timeseries. The above two quaternions will be purely imaginary quaternions because $\gamma$ is always a unit quaternion. Then, to calculate the torque, you need to know the inertia tensor of your spacecraft and you use my formula (6) below, at the end of the following analysis.

To work with quaternions ($SU(2)$ representation of the rotation group, with $SU(2)\stackrel{\mathrm{Ad}}{\rightarrow} SO(3)$), you do the following.
As you likely know, a point  in 3-space are represented by the purely imaginary quaterion $X=x\,\mathbf{i}+y\,\mathbf{j}+z\,\mathbf{k}$, and the rotation represented by quaternion $\gamma$ acts on $X$ through 
$$X\mapsto \gamma\,X\,\gamma^{-1}\tag{1}$$
Differentiate this and you'll find that:
$$\dot{X} = \gamma\,[\Omega,X]\,\gamma^{-1}
\tag{2}$$
where $\Omega = \gamma^{-1}\,\dot{\gamma}$ is also a pure quaternion and represents the angular velocity. Indeed, if you write $\Omega = \omega_x\,\mathbf{i}+\omega_y\,\mathbf{j}+\omega_z\,\mathbf{k}$, then the operation represented by the Lie bracket $[\Omega,X] = \Omega\,X-X\,\Omega$ (as matrices) is indeed the cross product $\Omega\times X$ when $\Omega$ and $X$ are thought of as 3-vectors. The moment of the velocity in (2) is:
$$L = \gamma\,[X,\,[\Omega,\,X]]\,\gamma^{-1}\tag{3}$$
At this point, it becomes a little hard to stay in quaternion notation completely because we need to integrate (3) over the whole body to get the inertia tensor. We need to calculate $\int_{body} \rho(X)\,[X,\,[X,\,\Omega]]\mathrm{d} V$; this is simply a $3\times3$ matrix operator that acts on the three components of $\Omega$ and given by:
$$I = \int_{body} \rho(X)\,\mathrm{ad}(X)^2 \mathrm{d} V\tag{4}$$
where $\mathrm{ad}(X)$ is the $3\times 3$ matrix of the linear operation $Y\mapsto [X,\,Y]$ for a pure quaternion $Y$. 
Let's write this linear mapping on $\Omega$ by the inertia tensor $I(\Omega)$.
To calculate the torque, we need to calculate the time derivative of (3); repeating the trick we used to get (2), we find:
$$\tau = \gamma\,\left([\Omega,\,I(\Omega)]+\mathrm{d}_t I(\Omega)\right)\,\gamma^{-1}\tag{5}$$
which, when we rotate with the frame attached to the rigid body, becomes (on transforming $\tau\mapsto \gamma^{-1}\,\tau\,\gamma$):
$$\tau = [\Omega,\,I(\Omega)]+ I(\dot{\Omega})\tag{6}$$
which is the equivalent of obtaining the Euler equations through the formula $\mathrm{d}_t = D_t + \Omega\times$.
A: To find the average angular acceleration $\boldsymbol{\alpha}$ over the time period $\Delta t$, use the formula $$\boldsymbol{\alpha}=\frac{\boldsymbol{\omega}_f-\boldsymbol{\omega}_i}{\Delta t},$$
where $\boldsymbol{\omega}_i$ is the angular velocity at the beginning of time period and $\boldsymbol{\omega}_f$ is the angular velocity at the end. 
Now you say you already know how to compute the torque, but for completeness, I will give the formula for this is well. The average torque $\boldsymbol{\tau}$ is given by
$$ \boldsymbol{\tau}=\frac{\mathbf{I}_f\boldsymbol{\omega}_f-\mathbf{I}_i\boldsymbol{\omega}_i}{\Delta t},$$
where $\mathbf{I}_i$ is the moment of inertia tensor at the beginning and $\mathbf{I}_f$ is the moment of inertia tensor at the end.
