The algorithm samples a random distance based on the highest extinction (or scattering or absorption) coefficient along the ray (or line segment in my case), or any coefficient greater than or equal to that (such as the highest in the entire medium, which is often easier to find). Then it probabilistically takes another leap based on the actual coefficient at the sampled location.
|Algorithm from Unbiased Global Illumination with Participating Media.|
Using the lowest max coefficient allows larger steps to be taken, which decreases the number of steps needed and increases performance. Since Rayleigh and ozone scattering and absorption coefficients vary predictably based on wavelength, the lowest max coefficient can be found for given wavelength, and with Rayleigh scattering, the lowest coefficient along the line segment can be found by using the coefficient at the point on the line segment that is closest to the center of the Earth (assuming density is monotonically decreasing with altitude).
After implementing this algorithm and the related optimizations, not only are my renders unbiased, but they are also 5 to 10 times faster—a pretty huge speed increase.