Let $\mathcal{I}^3\subset\mathbb{R}^4$
be the standard hyperboloid model for hyperbolic $3$-space and consider the usual $\mathrm{SO}(3,1)$
action of $\mathrm{PSL}(2,\mathbb{C})$
on $\mathcal{I}^3$.
Given a matrix
$\gamma=\begin{pmatrix}
a & b\\
c & d
\end{pmatrix}\in\mathrm{PSL}(2,\mathbb{C})$,
there is a unique geodesic $g\in\mathcal{I}^3$
preserved by the action of $\gamma$.
Particularly,
$\gamma$
translates along and rotates around $g$
(and these quantities can be written in terms of the trace $a+d$).

I want to know how to determine the axis $g$
in terms of the entries of $\gamma$.
One way to go about this would be to use the fact that $g$
is the intersection of $P\cap\mathcal{I}^3$
where $P$
is a uniquely determiend Euclidean plane passing through the origin in $\mathbb{R}^4$.
Perhaps we can write down two vectors in terms of $a,b,c,d$
that span this plane,
using the fact that the extended action of $\gamma$
to $\mathbb{R}^4$
preserves the plane.

By comparison,
in the upper half-space model $\mathcal{H}^3$,
the fixed godesic will be the Euclidean half-circle
or half-line
orthogonal to $\partial\mathcal{H}^3$
that connects the points $\dfrac{a-d\pm\sqrt{(a-d)^2+4bc}}{2c}$,
where we interpret these as $0$
and $\infty$
when $c=0$.
But in the hyperboloid model,
points on the boundary are not so accessible.
So how do we describe the fixed geodesic there?