Spectrograph simulation
=======================

.. autoclass:: qmostetc.Spectrograph
    :members:
    :undoc-members:


Default 4MOST spectrograph properties
-------------------------------------

The spectrograph properties together with the telescope properties are
generated with TOAD in a `separate project
<https://www.4most.eu/readthedocs/is/etc_data>`_. 
All data are stored in the file
``qmostetc/data/TOAD_4FS_ETC_interface_data.fits`` which is a FITS
file with several extensions. The extensions relevant for the
spectrographs are:

**Material**

  The transmission of the system, including telescope mirrors and
  spectrograph. Relevant columns:

  * **spectrograph**: Spectrograph (``lrs`` or ``hrs``)
  * **spectrograph_arm**: Spectrograph arm (one of ``blue``,
    ``green``, ``red``)
  * **design_wave_interval**: Wavelength interval the spectrograph
    arm is designed for [nm].
  * **wavelength_vector**: Wavelength information in [nm].
  * **transmission_vector**: Transmission curve [counts/photon]
  * **extraction_vector**: Efficiency to extract the spectral bins from
    the image.

  The spectrograph takes the product of transmission and extraction to
  convert the photon flux into counted electrons.

  .. exec::

     import sys
     from astropy.table import QTable
     material = QTable.read('qmostetc/data/toad/TOAD_4FS_ETC_interface_data.fits',
                            hdu='Material')
     print("Design wavelength intervals:\n")
     for row in material:
         print("* **{spectrograph} / {spectrograph_arm}**: "
	       "{design_wave_interval[0]:.1f} … {design_wave_interval[1]:.1f}"
	       .format(**row))

  .. plot::

     import matplotlib.pyplot as plt
     from astropy.table import QTable

     plt.title('Transmission × Extraction')
     material = QTable.read('qmostetc/data/toad/TOAD_4FS_ETC_interface_data.fits',
                            hdu='Material')
     lstyle = {'hrs': 'solid', 'lrs': 'dashed'}
     for row in material:
         plt.plot(row['wavelength_vector'],
                  row['transmission_vector']*row['extraction_vector'],
                  label=f"{row['spectrograph']} / {row['spectrograph_arm']}",
                  color=row['spectrograph_arm'],
                  linestyle=lstyle[row['spectrograph']])
     plt.xlabel(f"Wavelength [{row['wavelength_vector'].unit}]")
     plt.ylabel("Transmission [counts/photon]")
     plt.legend()
     plt.show()


**Resolution**

  Spectrograph resolution curve. Relevant columns:

  * **spectrograph**: Spectrograph (``lrs`` or ``hrs``)
  * **spectrograph_arm**: Spectrograph arm (one of ``blue``,
    ``green``, ``red``)
  * **wavelength_vector**: Wavelength vector [nm]
  * **fwhm_vector**: Full Width Half Maximum vector [nm]

  .. plot::

     import matplotlib.pyplot as plt
     from astropy.table import QTable

     plt.title('FWHM')
     resolution = QTable.read('qmostetc/data/toad/TOAD_4FS_ETC_interface_data.fits',
                              hdu='Resolution')
     lstyle = {'hrs': 'solid', 'lrs': 'dashed'}
     for row in resolution:
         plt.plot(row['wavelength_vector'], row['fwhm_vector'],
                  label=f"{row['spectrograph']} / {row['spectrograph_arm']}",
                  color=row['spectrograph_arm'],
                  linestyle=lstyle[row['spectrograph']])
     plt.xlabel(f"Wavelength [{row['wavelength_vector'].unit}]")
     plt.ylabel(f"FWHM [{row['fwhm_vector'].unit}]")
     plt.legend()
     plt.show()

     plt.title('Resolution (wavelength/FWHM)')
     for row in resolution:
         plt.plot(row['wavelength_vector'],
	          row['wavelength_vector']/row['fwhm_vector'],
                  label=f"{row['spectrograph']} / {row['spectrograph_arm']}",
                  color=row['spectrograph_arm'],
                  linestyle=lstyle[row['spectrograph']])
     plt.xlabel(f"Wavelength [{row['wavelength_vector'].unit}]")
     plt.ylabel("Resolution")
     plt.legend()
     plt.show()


**Dispersion**

  This binary table contains the pixel position on the detector. The
  columns which are relevant are:

  * **spectrograph**: Spectrograph (``lrs`` or ``hrs``).
  * **spectrograph_arm**: Spectrograph arm (one of ``blue``, ``green``, ``red``).
  * **wavelength_vector**: Wavelength information in [nm].
  * **pixel_width_vector**: Spectral width of the pixels [nm].

  .. plot::

     import matplotlib.pyplot as plt
     from astropy.table import QTable

     plt.title('Pixel width')
     dispersion = QTable.read('qmostetc/data/toad/TOAD_4FS_ETC_interface_data.fits',
                              hdu='Dispersion')
     lstyle = {'hrs': 'solid', 'lrs': 'dashed'}
     for row in dispersion:
         plt.plot(row['wavelength_vector'], row['pixel_width_vector'],
                  label=f"{row['spectrograph']} / {row['spectrograph_arm']}",
                  color=row['spectrograph_arm'],
                  linestyle=lstyle[row['spectrograph']])
     plt.xlabel(f"Wavelength [{row['wavelength_vector'].unit}]")
     plt.ylabel(f"Pixel width [{row['pixel_width_vector'].unit}]")
     plt.legend()
     plt.show()


**Detector**

  CCD properties. The columns relevant for the ETC are:

  * **spectrograph**: Spectrograph (``lrs`` or ``hrs``)
  * **spectrograph_arm**: Spectrograph arm (one of ``blue``,
    ``green``, ``red``)
  * **dark_current**: The dark current in electrons per pixel per
    second at the operation temperature of the instrument [count/s].
  * **read_noise**: The read noise per pixel, including the amplifier
    electronics [count].
  * **gain**: ADU gain [count/adu]
  * **saturation**: Saturation value [adu]

  .. exec::

     import sys
     from astropy.table import Table
     detector = Table.read('qmostetc/data/toad/TOAD_4FS_ETC_interface_data.fits',
                           hdu='Detector')
     detector['gain'].format = '%.4f'
     detector['dark_current'].format = '%.6f'
     detector['spectrograph',
              'spectrograph_arm',
	      'dark_current',
	      'read_noise',
	      'gain',
	      'saturation'].write(sys.stdout, format="ascii.rst")