Because the refractive index can be expressed as a function of $\theta$ (the angle between the wave vector and the ambient magnetic field) in the case of a cold plasma, it is possible to solve for the ray refractive index ($n_r$) analytically. This would, of course, reduce the computational overhead of both calculating emission coefficients and ray tracing. This is not important in the present iteration of the code since $n_r$ is only calculated once per emission coefficient. The general form of the equation is:
$$
n_r^2=\left|n^2\sin{\theta}\frac{y}{\frac{\partial}{\partial\theta}\left(\frac{\cos{\theta}+x\sin{\theta}}{y}\right)}\right|,
$$
where $x \equiv \left( \frac{1}{n}\frac{\partial n}{\partial\theta}\right)_\omega$ and $y \equiv \sqrt{1+x^2}$ are my own definitions for simplicity of computation.1 We can then plug in an analytical expression for $n$, that being the real part of the refractive index defined by the Appleton-Hartree equation for a uniformly magnetized cold plasma.
Now, we can implicitly differentiate the dispersion relation to determine $\frac{\partial n}{\partial \theta}$. This yields
$$A'n^4 + 4An^3n' - (B'n^2 + 2Bnn') = 0,$$
where $X' \equiv \frac{\partial X}{\partial \theta}$, $A = S\sin^2\theta + P\cos^2\theta$, $B = RL\sin^2\theta + PS(1+\cos^2\theta)$, and $R$, $L$, $S$ and $P$ are all constant with respect to $\theta$.2 Differentiating again,
$$A' = (S-P)\sin(2\theta),$$
$$B' = (RL-PS)\sin(2\theta).$$
Solving the differentiated dispersion relation for $n'$, we find
$$\frac{\partial n}{\partial \theta} = \frac{n}{2}\left(\frac{B' - A'n^2}{2An^2 - B}\right).$$
Now, let’s solve the remaining expressions.
$$\frac{\partial}{\partial\theta}\left(\frac{\cos\theta + x\sin\theta}{y}\right) = \frac{(x'\sin\theta + x\cos\theta - \sin\theta)y - y'(x\sin\theta + \cos\theta)}{y^2},$$
where
$$x' = \frac{n'' - x}{n},$$
$$y' = \frac{xx'}{y}.$$
This necessitates one further unwieldy derivative calculation for $n''=q/d$, where
$$q \equiv n'n^2(3A'B - 2AB') - 2AA'n'n^4 + n^3(A''B - 3A'B' + 2AB'') + 2n^5(A'^2 - AA'') - BB'n' + n(B'^2 - BB''),$$
$$d \equiv 2(B - 2An^2)^2,$$
$$A'' = 2(S-P)\cos(2\theta),$$
$$B'' = 2(RL-PS)\cos(2\theta).$$
The final expression for $n_r$ can then be readily assembled via analytical computation.
This calculation has been implemented here, however, the results do not display symmetry about the magnetic field axis and therefore either the derivation above or the implementation are incorrect.
