DL2_to_DL3.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "eb8526e8",
   "metadata": {},
   "source": [
    "## DL2 to DL3"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f52f6f1f",
   "metadata": {},
   "source": [
    "In this notebook, we will combine the DL2 *.h5* files with the instrument response function (IRF) to achieve the DL3 with magic-cta-pipe (MCP). We will use a small data sample due to time constraints.\n",
    "\n",
    "Let's start importing some basic libraries:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "e6187f24",
   "metadata": {},
   "outputs": [],
   "source": [
    "import glob\n",
    "import logging\n",
    "import os\n",
    "import subprocess\n",
    "import sys\n",
    "log = logging.getLogger()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b0146a19",
   "metadata": {},
   "source": [
    "### Processing DL2 to DL3\n",
    "\n",
    "Here we run a script that:\n",
    "1) produces DL3 data (i.e.: DL2 data after gammaness and theta cuts). To get the cut values for every run, it interpolates the cut tables produced with the MC test sample, using the run average pointing coordinates.\n",
    "2) interpolates IRFs (actually in this tutorial we don't interpolate, but we show you how to do it). It also interpolates the effective area values calculated with the MC test samples for every run using its average pointing coordinates.\n",
    "\n",
    "The DL3 and IRFs are packed together in the output files (*.fits.gz*). \n",
    "\n",
    "In the configuration file, we set:\n",
    "\n",
    "1) `interpolation_method: \"nearest\"`, where we use the nearest MC point to estimate the IRF of our target. Ideally, we should use the other options, i.e. \"linear\", where we do a linear interpolation of 3 MC points around our target, or \"cubic\", where we do a cubic interpolation instead.\n",
    "2) Maximum angular distance for interpolation.\n",
    "3) Target name or coordinates.\n",
    "\n",
    "Let's take a look at the configuration file section suited for this analysis:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "a5137b3c",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "dl2_to_dl3:\n",
      "    interpolation_method: \"nearest\"  # select \"nearest\", \"linear\" or \"cubic\"\n",
      "    interpolation_scheme: \"cosZd\" # select \"cosZdAz\" or \"cosZd\"\n",
      "    max_distance: \"45. deg\"  # angle type Quantity, or comment out to remove the cut\n",
      "    source_name: \"Crab\"\n",
      "    source_ra: null  # used when the source name cannot be resolved\n",
      "    source_dec: null  # used when the source name cannot be resolved"
     ]
    }
   ],
   "source": [
    "%%capture\n",
    "os.system(\"sed -n '260,266p' /fefs/aswg/workspace/analysis-school-2024/MCP/dl3/conf/config.yaml\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e5220f27",
   "metadata": {},
   "source": [
    "Let's setup the configuration and data paths:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "d9223a4d",
   "metadata": {},
   "outputs": [],
   "source": [
    "# The output directory where the DL3 data will be saved:\n",
    "#dir_dl3 = '/fefs/aswg/workspace/analysis-school-2024/MCP/dl3/output'\n",
    "dir_dl3='/fefs/aswg/workspace/YOUR_WORKSPACE/DL3/output'\n",
    "f=open(f'{dir_dl3}/DL3.log','w')\n",
    "\n",
    "# The input directory with the IRF:\n",
    "dir_irf='/fefs/aswg/workspace/analysis-school-2024/MCP/irf/output'\n",
    "\n",
    "# The input directory with the DL2 files:\n",
    "input_data=('/fefs/aswg/workspace/analysis-school-2024/MCP/dl3/input/*')\n",
    "input_file_data = glob.glob(input_data)\n",
    "input_file_data.sort()\n",
    "\n",
    "# The configuration file directory:\n",
    "config=('/fefs/aswg/workspace/analysis-school-2024/MCP/dl3/conf/config.yaml')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d3ba9faa",
   "metadata": {},
   "source": [
    "Here we use subprocess to run the function that converts DL2 to DL3 and get a log file:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "a4561159",
   "metadata": {},
   "outputs": [],
   "source": [
    "for input_file in input_file_data:  \n",
    "    y = subprocess.run(['lst1_magic_dl2_to_dl3', f'-d{input_file}', f'-i{dir_irf}',\\\n",
    "         f'-o{dir_dl3}', f'-c{config}'], stdout=f, stderr=f)        "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bdd07e5c",
   "metadata": {},
   "source": [
    "### Creating Index files\n",
    "\n",
    "Then we run the command below to create **Header Data Unit (HDU)** and **Observation** index files, needed by Gammapy to process DL3 data. Also these files are saved in the DL3 output directory."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "2f17d3a1",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "subprocess.run(['create_dl3_index_files', f'-i{dir_dl3}'],stdout=f, stderr=f);    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "eafa9362",
   "metadata": {},
   "outputs": [],
   "source": [
    "more $dir_dl3/DL3.log"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c65a7788",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "a02bc283",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "DL2_to_DL3.ipynb  SED_and_LC.ipynb\r\n"
     ]
    }
   ],
   "source": [
    "ls"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a65828d5",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.11.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}