example/Tutorial_example-GRS.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# GaussPy+ tutorial for a test field of the Galactic Ring Survey\n",
    "\n",
    "This notebook is intended to guide users through the typical GaussPy+ procedures for the decomposition of a position-position-velocity (PPV) dataset. \n",
    "\n",
    "For more information exceeding this tutorial we recommend taking a look at the following papers and documents:\n",
    "\n",
    "- For a description about the GaussPy+ decomposition package see:\n",
    "> - [Riener et al. 2019](https://arxiv.org/abs/1906.10506)\n",
    "\n",
    "- For a description about the GaussPy algorithm see: \n",
    "> - [GaussPy documentation](https://gausspy.readthedocs.io/en/latest/)\n",
    "> - [Lindner et al. 2015](https://arxiv.org/abs/1409.2840)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# import necessary libraries and plotting functions\n",
    "\n",
    "import os\n",
    "\n",
    "from astropy.io import fits\n",
    "\n",
    "%matplotlib inline\n",
    "import matplotlib\n",
    "import matplotlib.pyplot as plt\n",
    "from pylab import cm\n",
    "\n",
    "from astropy.wcs import WCS\n",
    "\n",
    "from gausspyplus.plotting import get_points_for_colormap, shiftedColorMap\n",
    "\n",
    "\n",
    "def get_cmap_rchi2(vmin, vmax):\n",
    "    orig_cmap = matplotlib.cm.RdBu_r\n",
    "    start, stop = get_points_for_colormap(vmin, vmax, central_val=1.)\n",
    "    midpoint = (1 - vmin) / (vmax - vmin)\n",
    "    return shiftedColorMap(orig_cmap, start=0., midpoint=midpoint, stop=stop)\n",
    "\n",
    "\n",
    "def add_style(ax):\n",
    "    ax.set_xlabel('Galactic Longitude')\n",
    "    ax.set_ylabel('Galactic Latitude')\n",
    "\n",
    "\n",
    "if not os.path.exists('decomposition_grs'):\n",
    "    !mkdir decomposition_grs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Overview\n",
    "\n",
    "This tutorial consists of the following steps:\n",
    "\n",
    "|        |  |\n",
    "| :-      |:-                                                |\n",
    "| 0.       | [Create a configuration file (optional)](#step0) |\n",
    "| 1.       | [Create a training set from the data](#step1) |\n",
    "| 2.       | [Find the optimal values for the smoothing parameters $\\alpha_{1}$ and $\\alpha_{2}$](#step2) |\n",
    "| 3.       | [Prepare the data for the decomposition](#step3) |\n",
    "| 4.       | [Decomposition of the data](#step4) |\n",
    "| 5.       | [Spatially coherent refitting - phase 1](#step5) |\n",
    "| 6.       | [Spatially coherent refitting - phase 2](#step6) |\n",
    "| 7.       | [Finalize the decomposition results](#step7) |\n",
    "|          | [Example: Create FITS files of fit parameters](#example) |\n",
    "\n",
    "The directory `gausspyplus/data` contains the FITS cube `grs-test_field.fits`, which is a subset of the Galactic Ring Survey ([Jackson et al. 2006](https://ui.adsabs.harvard.edu/abs/2006ApJS..163..145J/abstract)) that we will use for this tutorial. This is the same dataset that was used as a test field in Riener+ 2019."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"step0\"></a>\n",
    "## 0. Create a configuration file (optional)\n",
    "\n",
    "We start by generating a default configuration file. The next code creates a configuration file called 'gausspy+.ini' in the current working directory."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import gausspyplus.config_file as cf\n",
    "cf.make()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This configuration file only contains the most essential parameters, which in many cases should already be sufficient to create first good results. \n",
    "\n",
    "To take a look or access the full range of keywords that can be changed by the user change the above command to `cf.make(all_keywords=True)`.\n",
    "\n",
    "It is not necessary to create a configuration file for running `GaussPy+`. If no configuration file is supplied `GaussPy+` will resort to the default value for the parameters. \n",
    "\n",
    "Parameters can also be supplied or changed later on as will be done in the scripts we call further below. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"step1\"></a>\n",
    "## 1. Create a training set from the data\n",
    "\n",
    "Now we create a training set from the FITS cube `grs-test_field.fits`. For this we execute the `step_1-training_set--grs.py` script contained in the `example` directory. See the `step_1-training_set--grs.py` script for more comments.\n",
    "\n",
    "<div class=\"alert alert-block alert-warning\">\n",
    "NOTE: Running this script will use 75% of all CPUs on the machine you are running it unless the `use_ncpus` parameter is specified.\n",
    "</div>\n",
    "\n",
    "Depending on the number of available CPUs the execution of this script might take a couple of minutes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#  run the script\n",
    "!python step_1-training_set--grs.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Results\n",
    "After this script was successfully executed, the `decomposition_grs` directory should contain a folder named `gpy_training` with two files:\n",
    "\n",
    "- `grs-test_field-training_set_100_spectra.pickle`: a pickled dictionary of the decomposition results of the training set\n",
    "\n",
    "- `grs-test_field-training_set_100_spectra_plots.pdf`: plots of the decomposition results. We recommend to always plot a subsample of the training set to check whether the fitting of the training set worked out well.\n",
    "\n",
    "Since this example serves only to illustrate how training sets can be created we kept the number of spectra of the training set deliberately low. We recommend to include at least 200 spectra in the training set to get good results for the smoothing parameters."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"step2\"></a>\n",
    "## 2. Find the optimal values for the smoothing parameters $\\alpha_{1}$ and $\\alpha_{2}$\n",
    "\n",
    "After we checked that the decomposition of the training set gave good results, we can supply it to the machine learning procedure that GaussPy employs to find the best smoothing parameters values $\\alpha_{1}$ and $\\alpha_{2}$. For this we execute the `step_2-train--grs.py` script contained in the `example` directory. See the `step_2-train--grs.py` script for more comments.\n",
    "\n",
    "<div class=\"alert alert-block alert-warning\">\n",
    "NOTE: Running this script will use 75% of all CPUs on the machine you are running it unless the `use_ncpus` parameter is specified.\n",
    "</div>\n",
    "\n",
    "Depending on the number of available CPUs the execution of this script might take a couple of minutes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#  run the script\n",
    "!python step_2-train--grs.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Results\n",
    "\n",
    "In our case, we needed 47 iterations until a stable convergence was achieved. This yielded values for the smoothing parameters of $\\alpha_{1} = 2.58$ and $\\alpha_{2} = 5.14$.\n",
    "\n",
    "Note that the `decomposition_grs` directory now contains a new folder named `gpy_log`, which contains a log of the training results. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"step3\"></a>\n",
    "## 3. Prepare the data for the decomposition\n",
    "\n",
    "Next, we have to prepare our data cube for the decomposition. Executing the script `step_3-prepare--grs.py` will automatically:\n",
    "\n",
    "- estimate the root-mean-square noise $\\sigma_{\\mathrm{rms}}$ of the spectra\n",
    "\n",
    "- determine regions in the spectra that are likely to contain signal peaks\n",
    "\n",
    "- mask out negative noise spikes (by default all data peaks with a minimum $< -5 \\times \\sigma_{\\mathrm{rms}}$)\n",
    "\n",
    "- produce a pickled dictionary containing all necessary information for the decomposition\n",
    "\n",
    "<div class=\"alert alert-block alert-warning\">\n",
    "NOTE: Running this script will use 75% of all CPUs on the machine you are running it unless the `use_ncpus` parameter is specified.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#  run the script\n",
    "!python step_3-prepare--grs.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Results\n",
    "After this script was successfully executed, the `decomposition_grs` directory should contain a folder named `gpy_prepared` that contains `grs-test_field.pickle`, which is a pickled dictionary of the prepared data cube.\n",
    "\n",
    "The `decomposition_grs` directory should also contain a new folder named `gpy_maps` that contains the file `grs-test_field_noise_map.fits`, which is a map of the estimated rms noise values of our GRS test field.\n",
    "\n",
    "Lets take a look at the noise map:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "filepath = os.path.join('decomposition_grs', 'gpy_maps', 'grs-test_field_noise_map.fits')\n",
    "noise = fits.getdata(filepath)\n",
    "wcs = WCS(fits.getheader(filepath))\n",
    "\n",
    "fig, ax = plt.subplots(figsize=(6, 4), subplot_kw=dict(projection=wcs))\n",
    "\n",
    "img_noise = ax.imshow(noise, cmap='plasma_r', vmin=0.075, vmax=0.375)\n",
    "fig.colorbar(img_noise, ax=ax, extend='max')\n",
    "ax.set_title('Noise map')\n",
    "add_style(ax)\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `step_3-prepare--grs.py` script also plotted all spectra within pixel ranges $30 \\leq \\text{X} \\leq 34$ and $25 \\leq \\text{Y} \\leq 29$ and saved it as `grs-test_field_plots.pdf` in the `gpy_plots` directory.\n",
    "\n",
    "The spectra are plotted in the correct spatial order to aid in comparing differences between neighboring spectra. The shaded red areas indicate the regions of the spectrum that were identified to contain signal (goodness of fit calculations are restricted to these areas). The dotted horizontal red lines indicate the estimated noise values of $\\pm \\sigma_{\\mathrm{rms}}$ and the horizontal dashed line marks a S/N ratio of 3. The title of each subplot contains information about the location in terms of pixels and the index of the spectrum in the `grs-test_field.pickle` dictionary."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"step4\"></a>\n",
    "## 4. Decomposition of the data\n",
    "\n",
    "After the successful preparation of the data, we can proceed to the decomposition of the data. \n",
    "\n",
    "The following script will run an improved fitting routine on top of the original GaussPy decomposition routine; for more details see Sect. 3.2. in [Riener et al. 2019](https://arxiv.org/abs/1906.10506).\n",
    "\n",
    "<div class=\"alert alert-block alert-warning\">\n",
    "NOTE: Running this script will use 75% of all CPUs on the machine you are running it unless the `use_ncpus` parameter is specified.\n",
    "</div>\n",
    "\n",
    "Depending on the number of available CPUs the execution of this script might take a couple of minutes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#  run the script\n",
    "!python step_4-decompose--grs.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Results\n",
    "After this script was successfully executed, the `decomposition_grs` directory should contain a folder named `gpy_decomposed` that contains `grs-test_field_g+_fit_fin.pickle`, which is a pickled dictionary of the decomposition results.\n",
    "\n",
    "In addition, we also produced a map showing the number of fitted components and a map showing the $\\chi_{\\mathrm{red}}^{2}$ values of the fit. \n",
    "\n",
    "Lets take a look at both of these maps:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "filepath = os.path.join('decomposition_grs', 'gpy_maps', 'grs-test_field_g+_rchi2_map.fits')\n",
    "rchi2 = fits.getdata(filepath)\n",
    "wcs = WCS(fits.getheader(filepath))\n",
    "\n",
    "fig, axes = plt.subplots(ncols=2, figsize=(12, 4), subplot_kw=dict(projection=wcs))\n",
    "\n",
    "ax = axes.flatten()[0]\n",
    "\n",
    "vmin = min(rchi2.flatten())\n",
    "vmax = 2.5\n",
    "new_cmap = get_cmap_rchi2(vmin, vmax)\n",
    "\n",
    "img_rchi2 = ax.imshow(rchi2, cmap=new_cmap, vmin=vmin, vmax=vmax)\n",
    "fig.colorbar(img_rchi2, ax=ax, extend='max')\n",
    "ax.set_title('$\\chi_{\\mathrm{red}}^{2}$ map')\n",
    "add_style(ax)\n",
    "\n",
    "ax = axes.flatten()[1]\n",
    "\n",
    "ncomps = fits.getdata(os.path.join(\n",
    "    'decomposition_grs', 'gpy_maps', 'grs-test_field_g+_component_map.fits'))\n",
    "\n",
    "vmax = 6\n",
    "new_cmap = cm.get_cmap('Spectral_r', vmax + 1)\n",
    "\n",
    "img_ncomps = ax.imshow(ncomps, cmap=new_cmap, vmin=0, vmax=vmax)\n",
    "fig.colorbar(img_ncomps, ax=ax)\n",
    "ax.set_title('Number of fitted components')\n",
    "add_style(ax)\n",
    "\n",
    "fig.suptitle('After first decomposition with improved fitting routine')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"step5\"></a>\n",
    "## 5. Spatially coherent refitting - phase 1\n",
    "\n",
    "Next we will try to improve upon the decomposition results obtained in the last step. For this, we try to refit spectra that were flagged as not satisfying our chosen quality criteria. For more details see Sect.$\\,$3.3.1. in [Riener et al. 2019](https://arxiv.org/abs/1906.10506) for more details.\n",
    "\n",
    "For this we execute the `step_5-spatial_refitting-p1--grs.py` script contained in the `example` directory. See the `step_5-spatial_refitting-p1--grs.py` script for more comments.\n",
    "\n",
    "<div class=\"alert alert-block alert-warning\">\n",
    "NOTE: Running this script will use 75% of all CPUs on the machine you are running it unless the `use_ncpus` parameter is specified.\n",
    "</div>\n",
    "\n",
    "Depending on the number of available CPUs the execution of this script might take a couple of minutes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#  run the script\n",
    "!python step_5-spatial_refitting-p1--grs.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Results\n",
    "After this script was successfully executed, the `gpy_decomposed` folder contains the new file `grs-test_field_g+_fit_fin_sf-p1.pickle`, which is a pickled dictionary of the new decomposition results.\n",
    "\n",
    "In addition, we also produced new corresponding mapa of the number of fitted components and $\\chi_{\\mathrm{red}}^{2}$ values of the fit. \n",
    "\n",
    "Lets take a look again at both of these maps:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "filepath = os.path.join(\n",
    "    'decomposition_grs', 'gpy_maps', 'grs-test_field_g+_fit_fin_sf-p1_rchi2_map.fits')\n",
    "rchi2 = fits.getdata(filepath)\n",
    "header = fits.getheader(filepath)\n",
    "wcs = WCS(header)\n",
    "\n",
    "fig, axes = plt.subplots(ncols=2, figsize=(12, 4), subplot_kw=dict(projection=wcs))\n",
    "\n",
    "ax = axes.flatten()[0]\n",
    "\n",
    "vmin = min(rchi2.flatten())\n",
    "vmax = 2.5\n",
    "new_cmap = get_cmap_rchi2(vmin, vmax)\n",
    "\n",
    "img_rchi2 = ax.imshow(rchi2, cmap=new_cmap, vmin=vmin, vmax=vmax)\n",
    "fig.colorbar(img_rchi2, ax=ax, extend='max')\n",
    "ax.set_title('$\\chi_{\\mathrm{red}}^{2}$ map')\n",
    "add_style(ax)\n",
    "\n",
    "ax = axes.flatten()[1]\n",
    "\n",
    "ncomps = fits.getdata(os.path.join(\n",
    "    'decomposition_grs', 'gpy_maps', 'grs-test_field_g+_fit_fin_sf-p1_component_map.fits'))\n",
    "\n",
    "vmax = 6\n",
    "new_cmap = cm.get_cmap('Spectral_r', vmax + 1)\n",
    "\n",
    "img_ncomps = ax.imshow(ncomps, cmap=new_cmap, vmin=0, vmax=vmax)\n",
    "fig.colorbar(img_ncomps, ax=ax)\n",
    "ax.set_title('Number of fitted components')\n",
    "add_style(ax)\n",
    "\n",
    "fig.suptitle('After spatially coherent refitting - phase 1')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see that the $\\chi_{\\mathrm{red}}^{2}$ values of the fit did not change significantly; however the map of the number of fitted components shows already more spatial coherence. We can also check the new fit results for our 25 neighboring spectra that are plotted in `grs-test_field_g+_fit_fin_sf-p1_plots.pdf`. These plots also confirm that we gained a significant improvement in terms of spatial coherence of the fit results."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"step6\"></a>\n",
    "## 6. Spatially coherent refitting - phase 2\n",
    "\n",
    "In the last step, we try to further improve upon the decomposition results obtained in the last step by checking the coherence of the centroid positions of the fitted Gaussian components. See Sect.$\\,$3.3.2. in [Riener et al. 2019](https://arxiv.org/abs/1906.10506) for more details.\n",
    "\n",
    "For this we execute the `step_6-spatial_refitting-p2--grs.py` script contained in the `example` directory. See the `step_6-spatial_refitting-p2--grs.py` script for more comments.\n",
    "\n",
    "<div class=\"alert alert-block alert-warning\">\n",
    "NOTE: Running this script will use 75% of all CPUs on the machine you are running it unless the `use_ncpus` parameter is specified.\n",
    "</div>\n",
    "\n",
    "Depending on the number of available CPUs the execution of this script might take a couple of minutes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#  run the script\n",
    "!python step_6-spatial_refitting-p2--grs.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Results\n",
    "After this script was successfully executed, the `gpy_decomposed` folder contains the new file `grs-test_field_g+_fit_fin_sf-p2.pickle`, which is a pickled dictionary of the new decomposition results.\n",
    "\n",
    "In addition, we also produced new corresponding maps of the number of fitted components and $\\chi_{\\mathrm{red}}^{2}$ values of the fit. \n",
    "\n",
    "Lets take a look again at both of these maps:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "filepath = os.path.join(\n",
    "    'decomposition_grs', 'gpy_maps', 'grs-test_field_g+_fit_fin_sf-p2_rchi2_map.fits')\n",
    "rchi2 = fits.getdata(filepath)\n",
    "wcs = WCS(fits.getheader(filepath))\n",
    "\n",
    "fig, axes = plt.subplots(ncols=2, figsize=(12, 4), subplot_kw=dict(projection=wcs))\n",
    "\n",
    "ax = axes.flatten()[0]\n",
    "\n",
    "vmin = min(rchi2.flatten())\n",
    "vmax = 2.5\n",
    "new_cmap = get_cmap_rchi2(vmin, vmax)\n",
    "\n",
    "img_rchi2 = ax.imshow(rchi2, cmap=new_cmap, vmin=vmin, vmax=vmax)\n",
    "fig.colorbar(img_rchi2, ax=ax, extend='max')\n",
    "ax.set_title('$\\chi_{\\mathrm{red}}^{2}$ map')\n",
    "add_style(ax)\n",
    "\n",
    "ax = axes.flatten()[1]\n",
    "\n",
    "ncomps = fits.getdata(os.path.join(\n",
    "    'decomposition_grs', 'gpy_maps', 'grs-test_field_g+_fit_fin_sf-p2_component_map.fits'))\n",
    "\n",
    "vmax = 6\n",
    "new_cmap = cm.get_cmap('Spectral_r', vmax + 1)\n",
    "\n",
    "img_ncomps = ax.imshow(ncomps, cmap=new_cmap, vmin=0, vmax=vmax)\n",
    "fig.colorbar(img_ncomps, ax=ax)\n",
    "ax.set_title('Number of fitted components')\n",
    "add_style(ax)\n",
    "\n",
    "fig.suptitle('After spatially coherent refitting - phase 2')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Judging by the map showing the number of fitted components we seem to have gained more spatial coherence. The new and final fit results for our 25 neighboring spectra are plotted in `grs-test_field_g+_fit_fin_sf-p2_plots.pdf`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"step7\"></a>\n",
    "## 7. Finalize the decomposition results\n",
    "\n",
    "Finally, we can produce a table containing the final fit results from the second phase of the spatially coherent refitting stage.\n",
    "\n",
    "<div class=\"alert alert-block alert-warning\">\n",
    "NOTE: Running this script will use 75% of all CPUs on the machine you are running it unless the `use_ncpus` parameter is specified.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#  run the script\n",
    "!python step_7-finalize--grs.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Results\n",
    "After this script was successfully executed, the `gpy_decomposed` folder contains the new file `grs-test_field_g+_fit_fin_sf-p2_finalized.dat`, which is a table containing the entire decomposition results.\n",
    "\n",
    "By default the table is similar to Table 1 presented in [Riener et al. 2020](https://arxiv.org/abs/1911.06860) (see Sect. 3.1 in [Riener et al. 2020](https://arxiv.org/abs/1911.06860) for a description of the table columns)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<a id=\"example\"></a>\n",
    "## Example: Create FITS files of fit parameters\n",
    "\n",
    "Finally, we show an example of how it is possible to extract information on individual fit components. GaussPy+ saves all the relevant data as pickled dictionaries in the ``gpy_prepared`` and ``gpy_decomposed`` directories. By executing the following script, we create FITS files of the fit parameters (amplitude, mean position, FWHM) of the component with the lowest fitted mean position per spectrum."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#  run the script\n",
    "!python example_fits_first_component--grs.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Results\n",
    "After this script was successfully executed, the `gpy_maps` directory contains the following new files:\n",
    "\n",
    "- `grs-test_field_g+_fit_fin_sf-p2-map_amp_fc.fits`: a map of the amplitude values of the first fit components of all spectra\n",
    "- `grs-test_field_g+_fit_fin_sf-p2-map_vel_fc.fits`: a map of the mean velocity values of the first fit components of all spectra\n",
    "- `grs-test_field_g+_fit_fin_sf-p2-map_fwhm_fc.fits`: a map of the FWHM values of the first fit components of all spectra"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Note\n",
    "We could have combined the preparation and decomposition steps in a single script, but for illustrative purposes we have executed them seperately in this tutorial."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [conda env:py3p8]",
   "language": "python",
   "name": "conda-env-py3p8-py"
  },
  "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.8.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}