Yes. You must have made an error somewhere. To simulate identifying a rotated ellipsoid, you must first create points on the unrotated ellipsoid, then rotate and shift them under the rotation matrix $R$ and a shift vector $r_0$, thus generating a new set of points.
Now, the original (basic) ellipsoid had the equation
$$ r^T Q_0 r = 1 $$
The image of coordinate vector $r$ is $ r' = R r + r_0 $. It follows that the equation of the rotated/shifted ellipsoid is
$$ (r' - r_0)^T R Q_0 R^T (r' - r_0) = 1 $$
Let set $Q = R Q_0 R^T $ then the equation is
$$ (r - r_0)^T Q (r - r_0 ) = 1 $$
Note the $3 \times 3$ matrix $Q$ is of the form
$$ Q = \begin{bmatrix} a && d && e \\ d && b && f \\ e && f && c \end{bmatrix} $$
Thus, by multiplying out the above quadratic form, we obtain
$ a (x - x_0)^2 + b (y - y_0)^2 + c (z - z_0)^2 \\+ 2 d (x - x_0)(y - y_0) + 2 e (x - x_0)(z - z_0) + 2 f (y - y_0) (z - z_0) = 1 $
Thus the equation of the ellipsoid is of the form
$$ A x^2 + B y^2 + C z^2 + D x y + E x z + F y z + G x + H y + I z + J = 0 $$
All the points generated by the simulator will satisfy this equation, thus we can write
$$ A x_i^2 + B y_i^2 + C z_i^2 + D x_i y_i + E x_i z_i + F y_i z_i + G x_i + H y_i + I z_i + J = 0 $$
The unknowns here are the constants $A$ through $J$. So we'll define a vector $\theta$ as follows
$ \theta = [ A, B, C, D, E, F, G, H, I, J ]^T $
and we'll define the matrix $\Phi$ to be the $N \times 10$ matrix whose $i$-th row is
$ \Phi_i = \begin{bmatrix}x_i^2 && y_i^2 && z_i^2 && x_i y_i && x_i z_i && y_i z_i && x_i && y_i && z_i && 1 \end{bmatrix} $
Now our equations can be written as
$ \Phi \theta = 0 $
This is an $N$ dimensional vector. To be able to solve for $\theta$, we'll pre-multiply both sides of the equation by $ \Phi^T$, to get
$ \Phi^T \Phi \theta = 0 $
Which is a set of $10$ equation in $10$ unknowns. To get a non-zero solution, we'll truncate the system of equations, and take the first $9$ equations, thus disregarding the last equation. Now, using Gauss-Jordan elimination, we can solve the $ 9 \times 10$ system of equations for the the direction of vector $\theta$.
The solution will be of the form $ \theta = t \varphi $
where $\varphi$ is a fixed vector and $ t $ is an arbitrary scalar.
From here, we now have matrix $Q$ correct up to a scaling factor, as follows
$ Q = t \begin{bmatrix} \varphi_1 && \dfrac{1}{2} \varphi_4 && \dfrac{1}{2} \varphi_5 \\ \dfrac{1}{2} \varphi_4 && \varphi_2 && \dfrac{1}{2} \varphi_6 \\
\dfrac{1}{2} \varphi_5 && \dfrac{1}{2} \varphi_6 && \varphi_3 \end{bmatrix} $
But first, from our equation of the rotated/shifted ellipsoid, we have
$$ (r - r_0)^T Q (r - r_0) = 1 $$
which when expanded becomes
$$ r^T Q r + r^T a + b = 0 $$
where $ a = - 2 Q r_0 $ and $ b = r_0^T Q r_0 - 1 $
It follows that $r_0 = -\dfrac{1}{2} Q^{-1} a $
Now from our identification we have
$ a = t \begin{bmatrix} \varphi_7 && \varphi_8 && \varphi \end{bmatrix} $
Therefore, we can determine the shift vector $r_0$ because the $t$'s cancel out.
Once this is done, we can determine the value of $t$ from
$ b = r_0^T Q r_0 - 1 $
which is linear in $t$. Then our matrix $Q$ is fully specified.
But we're not done yet, to find the semi-axes of the ellipsoid we have to diagonalize $Q$ into $ Q = R D R^T $ and this can done using, for example, MATLAB, or Mathematica. The diagonal matrix $D$ has on its diagonal the reciprocal of the semi-axes lengths squared, i.e.
$ D = \text{diag}(\dfrac{1}{a^2}, \dfrac{1}{b^2}, \dfrac{1}{c^2} ) $
Hence, using this reverse-engineering process (i.e. identification) we were able to retrieve all the parameters of the original ellipsoid, its rotation matrix, and it shift vector.
