{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# How precise must be the input velocity parameter ?\n",
    "\n",
    "The input velocity (set through the parameter ``pos_cov``) must be given with a good precision. The error should not generally be larger than a spectrum channel (since the FWHM is around 1.5 channels)\n",
    "\n",
    "Here's some ways you can determine the size of a channel."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from orcs.process import SpectralCube\n",
    "# load spectral cube\n",
    "cube = SpectralCube('/home/thomas/M31_SN3.merged.cm1.1.0.hdf5')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compute the channel width (in cm-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "channel width (cm-1) 2.0983247457061225\n"
     ]
    }
   ],
   "source": [
    "# you can simply get it from the axis of the cube (which is in cm-1)\n",
    "print('channel width (cm-1)', cube.params.base_axis[1] - cube.params.base_axis[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "line width (cm-1) 2.6229059321332664\n",
      "number of steps at the right of ZPD (contributing to the resolution) 672\n",
      "ratio of the step number / the number of steps at the right of the ZPD 1.25\n",
      "channel width (cm-1) 2.098324745706613\n"
     ]
    }
   ],
   "source": [
    "# you can relate the channel size to the line width as the FWHM of a sinc line is equal to 1.20671 x width \n",
    "# and the width of a sinc line is equal to 1.25 the channel size by construction (because there are 25% more \n",
    "# steps made than necessary on the left size of the ZPD to be able to compute the phase)\n",
    "# the number of steps really contributing to the resolution are the steps made at the right of the ZPD which can be\n",
    "# esaily calculated\n",
    "print('line width (cm-1)', cube.params.line_fwhm / 1.20671)\n",
    "print('number of steps at the right of ZPD (contributing to the resolution)', cube.params.step_nb - cube.params.zpd_index)\n",
    "print('ratio of the step number / the number of steps at the right of the ZPD', cube.params.step_nb / float(cube.params.step_nb - cube.params.zpd_index))\n",
    "print('channel width (cm-1)', cube.params.line_fwhm / 1.20671 / 1.25)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Spectral resolution:  4733.531668752227\n",
      "Mean wavenumber 14980.989521972366\n",
      "Line width (cm-1) 3.1648651726294244\n",
      "channel width (cm-1) 2.0981778041978103\n"
     ]
    }
   ],
   "source": [
    "# Finally, you can relate the line width with the resolution at the mean wavenumber of the cube \n",
    "# and then compute the channel size\n",
    "print('Spectral resolution: ', cube.params.resolution)\n",
    "mean_cm1 = (cube.params.base_axis[-1] + cube.params.base_axis[0]) / 2\n",
    "print('Mean wavenumber', mean_cm1)\n",
    "line_fwhm = mean_cm1 / cube.params.resolution \n",
    "print('Line width (cm-1)',  line_fwhm)\n",
    "print('channel width (cm-1)', line_fwhm / 1.20671 / 1.25)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## In terms of velocity\n",
    "\n",
    "The channel width can be related to an uncertainty on the velocity of the lines at a given wavelength"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "channel width in km/s at the mean wavenumber of the axis -41.98773900973879\n"
     ]
    }
   ],
   "source": [
    "import orb.utils.spectrum\n",
    "channel_width = cube.params.base_axis[1] - cube.params.base_axis[0]\n",
    "mean_cm1 = (cube.params.base_axis[-1] + cube.params.base_axis[0]) / 2\n",
    "print('channel width in km/s at the mean wavenumber of the axis', orb.utils.spectrum.compute_radial_velocity(mean_cm1 + channel_width, mean_cm1, wavenumber=True))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "channel width in km/s at the mean wavenumber of the axis 42.019749282152794\n"
     ]
    }
   ],
   "source": [
    "# it could be directly computed with the following formula\n",
    "c = 3e5 # lightspeed in km/s\n",
    "print('channel width in km/s at the mean wavenumber of the axis', c * channel_width / mean_cm1)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.7"
  },
  "toc": {
   "colors": {
    "hover_highlight": "#DAA520",
    "navigate_num": "#000000",
    "navigate_text": "#333333",
    "running_highlight": "#FF0000",
    "selected_highlight": "#FFD700",
    "sidebar_border": "#EEEEEE",
    "wrapper_background": "#FFFFFF"
   },
   "moveMenuLeft": true,
   "nav_menu": {
    "height": "80px",
    "width": "252px"
   },
   "navigate_menu": true,
   "number_sections": true,
   "sideBar": true,
   "threshold": 4,
   "toc_cell": false,
   "toc_section_display": "block",
   "toc_window_display": false,
   "widenNotebook": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}