That's all right and very deep, but I claim that this is
not what Bjorken and Drell used when they solved Dirac equation in spherical coordinates. Instead, they used a dirty but smart trick that bypasses all this deep math. Let me explain it in detail.
Being deep:
In general curved coordinates the Dirac equation has the form
$$[i\gamma^{\mu}(D_{\mu}-ieA_{\mu})+m]\psi=0 \;\;\;\;\;\; (1)$$
where ##A_{\mu}## is the electromagnetic potential and ##D_{\mu}## is the covariant derivative that carries the information about curved coordinates and local tetrads. The equation above can be written as
$$[i\gamma^{\mu}D_{\mu}+e\gamma^{\mu}A_{\mu}+m]\psi=0 \;\;\;\;\;\; (2)$$
The quantity ##D_{\mu}## is quite complicated to compute in practice, which is the price for being deep. Fortunately, there is a trick that in some cases can be used to bypass the computation of ##D_{\mu}##. That trick is used in Bjorken and Drell, which is what I discuss next.
Being smart:
Let us start from the laboratory Lorentz frame with coordinates ##x^{\mu}## in which Eq. (2) takes the form
$$[i\gamma^{\mu}_{\rm fix}\partial_{\mu}+e\gamma^{\mu}_{\rm fix}A_{\mu}+m]\psi_{\rm lab}=0 \;\;\;\;\;\; (3)$$
where ##\gamma^{\mu}_{\rm fix}## are the standard fixed Dirac matrices. This equation is not general covariant, but is valid in one smartly chosen system of coordinates. Now let us introduce some new curvilinear coordinates ##x'^{\nu}##. Instead of properly transforming everything to new coordinates (which would be quite complicated in practice), the smart dirty trick is just to use the identity
$$\partial_{\mu}=\frac{\partial x'^{\nu}}{\partial x^{\mu}} \partial'_{\nu} \;\;\;\;\;\; (4)$$
which implies that (3) can be written as
$$[i\gamma^{\mu}_{\rm fix}\frac{\partial x'^{\nu}}{\partial x^{\mu}} \partial'_{\nu}+e\gamma^{\mu}_{\rm fix}A_{\mu}+m]\psi_{\rm lab}=0 \;\;\;\;\;\; (5)$$
In this way we have written the Dirac equation in terms of curvilinear coordinates without calculating the complicated quantity ##D_{\mu}##. The price paid is that (5) is very very non-covariant. Yet it is correct in the given circumstances. And it is much simpler to deal with than would be the fully covariant equation. Basically, Bjorken and Drell solved Eq. (5) in the special case ##x'^{\nu}=(t,r,\theta,\varphi)##, ##A^{\mu}=(A^0(r),0,0,0)##.
Now let us see how my formalism works in this case. My formalism (just as the formalism used by Bjorken and Drell) is covariant only with respect to Lorentz transformations, not with respect to general coordinate transformations. Nevertheless, it can be used to solve Dirac equation in curvilinear coordinates by using the same smart dirty trick above. In my formalism, the Dirac equation in the laboratory Lorentz frame is
$$[i\Gamma^{\mu}_{\rm lab}\partial_{\mu}+e\Gamma^{\mu}_{\rm lab}A_{\mu}+m]\Psi=0 \;\;\;\;\;\; (6)$$
Having in mind that my formalism is related to the standard Bjorken and Drell formalism via
$$\Psi=\psi_{\rm lab}, \;\;\; \Gamma^{\mu}_{\rm lab}=\gamma^{\mu}_{\rm fix} \;\;\;\;\;\; (7)$$
we see that my Eq. (6) is equivalent to the standard Eq. (3). Therefore one can proceed with (4) and (5) as above.
The moral: Being deep is for mathematical physicists. Practical physicists must be smart.