A practical guide to writing a radiative transfer code

Using our decades-long experience in radiative transfer (RT) code development for Earth science, we endeavor to reduce the knowledge gap of bringing RT from theory to code quickly. Despite numerous classic and recent literature, it is still hard to develop an RT code from scratch within a few weeks. It is equally hard to understand, not to mention modify, an existing “monster” RT code, for which the developer is either located remotely or has retired. Following the format of “Numerical Recipes” by Press et al., we collocate in this paper small pieces of necessary theory with corresponding small pieces of RT code. These are arranged in an order that is natural for code development, which is often opposite of the natural order for laying out the theoretical basis. We focus on the transfer of unpolarized monochromatic solar radiation in a plane-parallel atmosphere over a reflecting surface. Both the surface and the atmosphere are homogeneous (uniform) at all directions. The multiple scattering is numerically solved using the deterministic method of Gauss-Seidel iterations. Except for the presented Python-Numba open-source RT code gsit, the paper does not report any new scientific results, but rather serves as an academic demonstration. If development time is an issue or the reader is familiar with basic concepts of RT theory, we recommend proceeding directly to Sec. 3 “RT code development”. Program title: gsit (pronounced “jeezit”) CPC Library link to program files:https://doi.org/10.17632/d3zt5zhx49.1 Developer's repository link:https://github.com/korkins/gsit Licensing provisions: MIT Programming language: Python 3 Nature of problem: We present a tutorial code in Python for deterministic (non-stochastic) numerical simulation of multiple scattering of monochromatic solar light in a plane-parallel Earth atmosphere bounded from below by a reflecting surface. The problem is solved in a simplified form (i.e., uniform atmosphere, no polarization, uniform surface reflectance, etc.) to better explain numerical features, rather than physics, of propagation of light in the atmosphere. Solution method: The method of Gauss-Seidel iterations. It relies on the Fourier decomposition of the Radiative Transfer Equation over azimuth, Gauss quadrature for numerical integration over the zenith and iterative process for integration over height (optical depth) with analytical (hence known) single scattering approximation being the starting point. The method is relatively simple