Extending background docs

docs/src/background/spectral_line_formation.ipynb → docs/src/background/forward_problem.ipynb

@@ -4,15 +4,39 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "# Spectral line formation"
+    "# Forward problem\n",
+    "---\n",
+    "\n",
+    "The forward problem describes how the physical properties of an (astrophysical) medium give rise to a spectral line observation.\n",
+    "We collectively represent the physical properties of an object by a model vector, $\\boldsymbol{m}$, and we represent an observation by a vector, $\\boldsymbol{o}$.\n",
+    "Hence, the forward problem can be thought of as a function, $f$, that maps a model vector, $\\boldsymbol{m}$, to the corresponding observation vector,\n",
+    "\\begin{equation*}\n",
+    "    \\boldsymbol{o} \\ = \\ f\\left( \\boldsymbol{m} \\right) .\n",
+    "\\end{equation*}\n",
+    "We distinguish two stages in the forward function, such that, $f\\left( \\boldsymbol{m} \\right) \\ = \\ f_{2}\\left( f_{1}\\left( \\boldsymbol{m} \\right) \\right)$.\n",
+    "The first part, $f_{1}$, describes the spectral line formation, i.e.\\ the radiative processes that determine the amount of electromagnetic radiation that emanates from the medium and travels towards our detectors.\n",
+    "This is described [here](##Spectral_line_formation).\n",
+    "The second part, $f_{2}$, describes the observation and includes the instrumentation effects on the observed signal. This is described in [here](##Observation)."
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {
+    "vscode": {
+     "languageId": "plaintext"
+    }
+   },
+   "source": [
+    "## Spectral line formation\n",
+    "---"
    ]
   },
   {
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "## Physical model\n",
-    "\n",
+    "### Physical model\n",
+    "---\n",
     "In pomme, we use Cartesian coordinates, $\\boldsymbol{x}=(x,y,z)$, and assume the line of sight to be along the positive $z$-axis, such that the plane of the sky corresponds to the $(x,y)$-plane.\n",
     "Furthermore, we will assume that the observer is located at $z=0$, and that the model box has a depth, $L$, along the $z$-axis.\n",
     "Then, in the absence of radiation scattering, considering a single ray through the model (resulting in a single pixel in the observed image), the observed intensity at $\\boldsymbol{x}=(x,y,0)$, at frequency, $\\nu$, reads,\n",
@@ -30,7 +54,14 @@
     "\\tau_{\\text{obs}}(\\nu; x, y, z)\n",
     "\\ \\equiv \\\n",
     "\\int_{0}^{z} \\text{d}z' \\ \\chi\\big( \\nu_{\\text{com}}(\\nu; \\, x, y, z'); \\, x, y, z' \\big) .\n",
-    "\\end{equation*}\n",
+    "\\end{equation*}"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "\n",
     "Since we look along the $z$-axis, the non-relativistic Doppler-shifted frequency in the co-moving frame is given by,\n",
     "\\begin{equation*}\n",
     "    \\nu_{\\text{com}}(\\nu, \\boldsymbol{x}) \\ = \\ \\left( 1 + \\frac{v_{z}(\\boldsymbol{x})}{c} \\right) \\nu ,\n",
@@ -92,14 +123,16 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "## Numerical implementation"
+    "### Numerical implementation\n",
+    "---"
    ]
   },
   {
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "### Line optical depth\n",
+    "#### Line optical depth\n",
+    "---\n",
     "Velocity gradients play a key role in line radiative transfer. Since spectral lines are narrowly peaked in frequency space, they are very sensitive to Doppler shifts, and thus motion (gradients), along the line-of-sight.\n",
     "Therefore, when numerically solving a line transfer problem, it is key to properly trace the velocity (gradient) along the line-of-sight. \n",
     "Since we assume to know the line profile function analytically, we can take care of this sharp frequency dependence by resolving its dependence analytically."
@@ -158,7 +191,8 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "### Line radiative transfer"
+    "#### Line radiative transfer\n",
+    "---"
    ]
   },
   {
@@ -204,6 +238,67 @@
     "The implementation of this intensity increment can be found in [src/pomme/forward.py](https://github.com/Magritte-code/pomme/blob/main/src/pomme/forward.py). For the same reaseon as with the line optical depth, we compute both cases for all values, to minimise the use of masks, and we add a small number ($10^{-30}$) to the denominator to avoid NaNs in gradients."
    ]
   },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## Observation\n",
+    "---\n",
+    "In the previous section, we presented a theoretical model for spectral line formation.\n",
+    "In practice, however, spectral line observations are affected by instrumentation effects, such as binning and noise, the particular form of which highly depends on the way the object is observed.\n",
+    "\n",
+    "High spatial and spectral resolution observations are typically obtained using interferometry, for which the situation is even more complicated, since there are spatial scale-dependent effects which cause the smallest and largest structures to be unresolved in the resulting images.\n",
+    "These effects must carefully be taken into account, as they (in part) determine which information about the spatial distribution of the physical properties is encoded in the observations and which is not.\n",
+    "At this point, we should clarify what we mean with interferometric data and distinguish between two types.\n",
+    "The first type are the visibilities, i.e. the correlations between the coherent electromagnetic waves, measured in several frequency bins for the different pairs of antennas in the array.\n",
+    "These visibilities are often further processed into images that map the actual brightness in the plane of the sky for each frequency bin.\n",
+    "This is what we distinguish as the second type of interferometric data.\n",
+    "The conversion of the sky brightness maps from the visibilities can be achieved, for instance, with the *clean* algorithm [(Hogbom 1974)](https://ui.adsabs.harvard.edu/abs/1974A%26AS...15..417H/abstract), as implemented, for instance, in [CASA](https://casa.nrao.edu).\n",
+    "This conversion is in itself already an inverse problem.\n",
+    "Therefore, when reconstructing a model from interferometric observations, it makes sense to start directly from the visibilities.\n",
+    "However, since visibilities are more difficult to interpret than sky brightness maps, we consider both options.\n",
+    "\n",
+    "\n",
+    "The conversion of the observed intensity in the plane of the sky, $I_{\\text{obs}}(\\nu; x, y)$, to visibilities, $V(\\nu;u,v)$, can be expressed as,\n",
+    "\\begin{equation*}\n",
+    "    V(\\nu; u, v) \n",
+    "    \\ = \\\n",
+    "    \\frac{1}{d^{2}}\n",
+    "    \\int_{-\\infty}^{+\\infty} \\text{d}x\n",
+    "    \\int_{-\\infty}^{+\\infty} \\text{d}y \\\n",
+    "    \\hat{A}(x,y) \\, I_{\\text{obs}}(\\nu; x, y) \\exp\\big( -2\\pi i \\left( u x + v y \\right) \\big)\n",
+    "\\end{equation*}\n",
+    "in which $\\hat{A}(x,y) \\equiv A(x,y) / A_{0}$, is the antenna response function, normalised by, $A_{0}$, the response at the centre of the beam, $d$, is the distance between the object and the observer, and $u$ and $v$ are the components of the projected distance vectors, often referred to as baselines, between pairs of antennas in the array.\n",
+    "Since there are only a finite number of antennas and thus a finite number of antenna pairs, the visibilities are only known at a finite number of baseline samples, $\\{(u_{d},v_{d})\\}_{d=1}^{N_{d}}$.\n",
+    "This sampling, which depends on the particular antenna configuration, causes a loss of information about the spatial distribution of the source in an intricate way.\n",
+    "We can see that the transformation of the intensity in the plane of the sky into visibilities is essentially a 2D Fourier transform.\n",
+    "Although the Fourier transform is invertible, the fact that we can only sample the visibilities for a finite number of baselines, makes that the combined transformation from a sky brightness map to a finite set of visibilities is not invertible.\n",
+    "\n",
+    "In principle, we could use the exact pipeline to model the instrumentation effects in our forward model, for instance, using [CASA](https://casa.nrao.edu).\n",
+    "However, later, in our reconstruction algorithm,\n",
+    "we want the implementation of the forward function, $f$, to be automatically differentiable using the *autograd* functionality in [PyTorch](https://pytorch.org).\n",
+    "Therefore, we implemented this map ourselves in [PyTorch](https://pytorch.org) using a fast Fourier transform (FFT), following the [Galario](https://github.com/mtazzari/galario) package for visibility modelling [(Tazzari et al. 2018)](https://ui.adsabs.harvard.edu/abs/2018MNRAS.476.4527T/abstract).\n",
+    "\n",
+    "\n",
+    "Alternatively, one may want to reconstruct 3D physical models from previously obtained sky brightness maps.\n",
+    "In that case, it is important to consider the beam, i.e.\\ a kernel that models the spatial spread of the intensity in the plane of the sky, and the antenna response function.\n",
+    "\\begin{equation*}\n",
+    "    \\tilde{I}_{\\text{obs}}(\\nu; x, y)\n",
+    "    \\ = \\\n",
+    "    \\frac{1}{d^{2}}\n",
+    "    \\int_{-\\infty}^{+\\infty} \\text{d}x'\n",
+    "    \\int_{-\\infty}^{+\\infty} \\text{d}y' \\\n",
+    "    \\hat{A}(x,y) \\,\n",
+    "    I_{\\text{obs}}(\\nu; x', y') \\,\n",
+    "    B(x-x',y-y') ,\n",
+    "\\end{equation*}\n",
+    "in which $B(x, y)$ is the beam kernel, which we assume to be a 2D Gaussian, centred around $(0,0)$, and, $d$, is again the distance between the object and the observer.\n",
+    "\n",
+    "It should be emphasised that in both processes of line formation and observation information is lost.\n",
+    "Mathematically, this implies that the function, $f$, is not invertible.\n",
+    "In the next section, we describe how, with a probabilistic approach, we can circumvent this problem and obtain a probabilistic inverse."
+   ]
+  },
   {
    "cell_type": "markdown",
    "metadata": {},

docs/src/background/index.rst

@@ -2,10 +2,12 @@ Background
 ##########
 
 Some background on the physics of spectral line observations, the Bayesian reconstruction algorithm, and their numerical implementation.
+This is all based `this paper<https://ui.adsabs.harvard.edu/abs/2024arXiv240218525D/abstract>`_.
+Please cite `the paper<https://ui.adsabs.harvard.edu/abs/2024arXiv240218525D/abstract>`_ if you use any of this information.
 
 .. toctree::
-   :maxdepth: 2
+   :maxdepth: 4
    :caption: Contents:
 
-   spectral_line_formation
-   bayesian_reconstruction
+   forward_model
+   inverse_model