For this answer I will not rely on Cartesian coordinates, but I will not work directly in terms of latitude and longitude either. Instead, it will be convenient to represent a point with spherical coordinates $(\phi,\theta)$ as $z=e^{i\phi}\tan\dfrac\theta2$, corresponding to a stereographic projection of the sphere into the complex plane. We can now apply Mobius transformations to move the complex numbers $\alpha$, $\beta$ (corresponding to $A$, $B$) to more helpful locations.
Specifically, we can use rotations to move $\alpha$ to zero, then map $\beta$ to the real line, and finally rotate $\{\alpha,\beta\}$ to $\{r,r^{-1}\}$ on the positive real line. This amounts to rotating our coordinates so that $A$ and $B$ have zero longitude and opposite latitudes. To make this simpler, we shall assume that $\alpha$ is real without loss of generality (i.e. we take the longitude of $B$ relative to $A$).
We now state these transforms explicitly. First, note that the subgroup of the Mobius group corresponding to rotations consists of mappings of the form $z\mapsto \dfrac{u z-\overline{v}}{v z+\overline{u}}$ for complex $u,v$. A mapping of this form which takes $\alpha$ to the origin is $z\mapsto \dfrac{z-\alpha}{\alpha z+1}$; this also maps $\beta\mapsto \dfrac{\beta-\alpha}{1+\beta\alpha}= R e^{i\psi}$ for some real $R,\psi$. A rotation $z\mapsto e^{-i\phi} z$ then brings $\beta$ to $R$.
For the last transformation, consider the mapping $z\mapsto \dfrac{1+z r}{r-z}$ for real $r>0$ which takes the real line to itself. Then $\{0,R\}\mapsto \left\{\dfrac{1}{r},\dfrac{1+r R}{r-R}\right\}$; since we want the images of $0$ and $R$ to be on equal and opposite sides of the 'equator', we require their product to be unity. From this we deduce $r=R+\sqrt{1+R^2}$. Composing the three transformation gives the net rotation
$$\mathcal{T}:\;z
\mapsto \dfrac{z-\alpha}{1+\beta\alpha}
\mapsto e^{i\phi}\left(\frac{z-\alpha}{1+\beta\alpha}\right)
\mapsto \dfrac{1+e^{i\phi}\left(\frac{z-\alpha}{1+\beta\alpha}\right)r}{r-e^{i\phi}\left(\frac{z-\alpha}{1+\beta\alpha}\right)}
=\dfrac{e^{-i\phi}(1+\beta\alpha)+(z-\alpha)r}{r e^{-i\phi}(1+\beta\alpha)-(z-\alpha)}$$
We now can give a criterion for the original question: a point $z$ lies between $\alpha$ and $\beta$ if $r^{-1}\leq | \mathcal{T}(z)|r$. In principle one can do inverse transformations to map this annulus back to the original coordinates, but I won't do so explicitly. Instead I'll just summarize the result:
Let $\alpha,\beta$ be the stereographic coordinates of two points $A,B$ on the unit sphere, and define $w=R e^{i\phi}:=\dfrac{\beta-\alpha}{1+\beta\alpha}$. Then a point $C$ on the unit sphere lies between $A$ and $B$ iff its sterographic coordinate $z$ satisfies
$$\frac{1}{r}< \left|\dfrac{e^{-i\phi}(1+\beta\alpha)+(z-\alpha)r}{r e^{-i\phi}(1+\beta\alpha)-(z-\alpha)}\right|<r$$ where $r=R+\sqrt{1+R^2}$.
There is a standard formula for this problem
on the Aviation Formulary web page.
There are also many other useful formulas on that page.
Now to consider various special cases.
One special case is when the two points are the same. Then the bearing (also known as a course angle) from one point to the other is not well defined.
Since the bearing is used in the general algorithm for the intersection,
we can't use that algorithm. But the obvious intersection point is
the same as the starting point(s).
If either point is at the north or south pole, every path from it has the same bearing from North, so the bearing doesn't actually tell you which path
is taken from that point. This special case has to be treated as an error,
I think, except perhaps if one point if the north pole and the other is the
south pole (in which case you might choose either one as "the" intersection).
The remaining special cases are for two points not the same,
neither one at a pole.
If the starting points are antipodes of each other, that is,
diametrically opposite each other on the Earth,
then each of the paths reaches the other point after traveling half way
around the Earth.
For example, the point $20$N $30$E is the antipodes of $20$S $150$W.
If the initial bearings of the two paths add up to any whole multiple
of $\pi$ radians ($180$ degrees), then the two paths travel along the
same great circle and there are infinitely many points of intersection.
If you must pick one, it could be one of the two starting points
or one of the two points midway between them; you could decide by
considering the initial bearings at each point.
If the points are not antipodes, it is still possible that the paths from the two points travel along exactly the same great circle.
For example, if the two starting points are $20$N $10$E and $40$N $10$E,
and both start with bearing $0$ degrees from North, then both paths travel the same great circle in the same direction. Using the same two starting points, if one path starts out at bearing $0$ degrees and the other
starts with bearing $180$ degrees, they will travel the same great
circle in opposite directions. In either case, the path from each point
goes through the other point and eventually circles around the Earth to
arrive at its starting point.
One way to tell if the paths travel the same great circle is to look at the bearing or course angle from each point to the other point.
If the direction from point $A$ to point $B$ is the same as the given
bearing of the path starting at $A$, then that path will pass through
both points. If the direction from $A$ to $B$ is exactly opposite
the given bearing (a difference of $\pm\pi$ radians or $\pm180$ degrees),
the path will still go through both points; you just have to go more than
half way around the Earth to reach the other point.
You can write a function to decide whether the given path from $A$ passes
through $B$. You can then use the same function to decide whether the
path from $B$ passes through $A$. If both paths pass through the other
point then they are both on the same great circle.
Technically, any point is equally "the" intersection point, so you just
need to pick one that makes sense to you: the first point, the second point, a point midway between them, or even the antipodes of the midway point.
The choice could depend on whether the initial bearing at each point is
directly toward or away from the other point.
In any other case, the paths will intersect in exactly two points.
The two intersections are antipodes of each other. If the initial bearings
at both points are toward the same point of intersection, I would choose that point as "the" intersection. (The intersection might be up to halfway around the Earth from either starting point.)
But it is possible that each path's initial bearing is toward a different
one of the two points of intersection. Then you have to come up with some criteria for picking one point or the other, depending possibly on which point is closer to which intersection. In the application for which I needed this algorithm, this case should not have come up normally, so I arbitrarily decided that "the" intersection point would be the intersection that was first reached by the first path.
Best Answer
If you project your vectors onto a plane, you can derive the angle from the dot product. $$a \cdot b = ||a||\,||b|| \cos \theta$$
Or you may want the spherical angle, in which case you will also need the radius of your sphere.