{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "
\n", "\n", "**This is a fixed-text formatted version of a Jupyter notebook**\n", "\n", "- Try online [![Binder](https://mybinder.org/badge.svg)](https://mybinder.org/v2/gh/gammapy/gammapy-webpage/v0.10?urlpath=lab/tree/light_curve.ipynb)\n", "- You can contribute with your own notebooks in this\n", "[GitHub repository](https://github.com/gammapy/gammapy/tree/master/tutorials).\n", "- **Source files:**\n", "[light_curve.ipynb](../_static/notebooks/light_curve.ipynb) |\n", "[light_curve.py](../_static/notebooks/light_curve.py)\n", "
\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Light curves\n", "\n", "## Introduction\n", "\n", "This tutorial explain how to compute a light curve with Gammapy.\n", "\n", "We will use the four Crab nebula observations from the [H.E.S.S. first public test data release](https://www.mpi-hd.mpg.de/hfm/HESS/pages/dl3-dr1/) and compute per-observation fluxes. The Crab nebula is not known to be variable at TeV energies, so we expect constant brightness within statistical and systematic errors.\n", "\n", "The main classes we will use are:\n", "\n", "* [gammapy.time.LightCurve](..\/api/gammapy.time.LightCurve.rst)\n", "* [gammapy.time.LightCurveEstimator](..\/api/gammapy.time.LightCurveEstimator.rst)\n", "\n", "## Setup\n", "\n", "As usual, we'll start with some setup..." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "import matplotlib.pyplot as plt\n", "import astropy.units as u\n", "from astropy.coordinates import SkyCoord, Angle\n", "from regions import CircleSkyRegion\n", "from gammapy.utils.energy import EnergyBounds\n", "from gammapy.data import DataStore\n", "from gammapy.spectrum import SpectrumExtraction\n", "from gammapy.spectrum.models import PowerLaw\n", "from gammapy.background import ReflectedRegionsBackgroundEstimator\n", "from gammapy.time import LightCurve, LightCurveEstimator" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Spectrum\n", "\n", "The `LightCurveEstimator` is based on a 1d spectral analysis within each time bin.\n", "So before we can make the light curve, we have to extract 1d spectra." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "data_store = DataStore.from_dir(\"$GAMMAPY_DATA/hess-dl3-dr1/\")" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Observations\n", "Number of observations: 4\n", "Info for OBS_ID = 23523\n", "- Start time: 53343.92\n", "- Pointing pos: RA 83.63 deg / Dec 21.51 deg\n", "- Observation duration: 1687.0 s\n", "- Dead-time fraction: 6.240 %\n", "Info for OBS_ID = 23526\n", "- Start time: 53343.95\n", "- Pointing pos: RA 83.63 deg / Dec 22.51 deg\n", "- Observation duration: 1683.0 s\n", "- Dead-time fraction: 6.555 %\n", "Info for OBS_ID = 23559\n", "- Start time: 53345.96\n", "- Pointing pos: RA 85.25 deg / Dec 22.01 deg\n", "- Observation duration: 1686.0 s\n", "- Dead-time fraction: 6.398 %\n", "Info for OBS_ID = 23592\n", "- Start time: 53347.91\n", "- Pointing pos: RA 82.01 deg / Dec 22.01 deg\n", "- Observation duration: 1686.0 s\n", "- Dead-time fraction: 6.212 %\n", "\n" ] } ], "source": [ "mask = data_store.obs_table[\"TARGET_NAME\"] == \"Crab\"\n", "obs_ids = data_store.obs_table[\"OBS_ID\"][mask].data\n", "observations = data_store.get_observations(obs_ids)\n", "print(observations)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "# Target definition\n", "target_position = SkyCoord(ra=83.63308, dec=22.01450, unit=\"deg\")\n", "on_region_radius = Angle(\"0.2 deg\")\n", "on_region = CircleSkyRegion(center=target_position, radius=on_region_radius)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "CPU times: user 978 ms, sys: 18.7 ms, total: 996 ms\n", "Wall time: 1.01 s\n" ] } ], "source": [ "%%time\n", "bkg_estimator = ReflectedRegionsBackgroundEstimator(\n", " on_region=on_region, observations=observations\n", ")\n", "bkg_estimator.run()" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "CPU times: user 878 ms, sys: 21.3 ms, total: 899 ms\n", "Wall time: 903 ms\n" ] } ], "source": [ "%%time\n", "ebounds = EnergyBounds.equal_log_spacing(0.7, 100, 50, unit=\"TeV\")\n", "extraction = SpectrumExtraction(\n", " observations=observations,\n", " bkg_estimate=bkg_estimator.result,\n", " containment_correction=False,\n", " e_reco=ebounds,\n", " e_true=ebounds,\n", ")\n", "extraction.run()\n", "spectrum_observations = extraction.spectrum_observations" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Light curve estimation\n", "\n", "OK, so now that we have prepared 1D spectra (not spectral models, just the 1D counts and exposure and 2D energy dispersion matrix), we can compute a lightcurve.\n", "\n", "To compute the light curve, a spectral model shape has to be assumed, and an energy band chosen.\n", "The method is then to adjust the amplitude parameter of the spectral model in each time bin to the data, resulting in a flux measurement in each time bin." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "# Creat list of time bin intervals\n", "# Here we do one time bin per observation\n", "def time_intervals_per_obs(observations):\n", " for obs in observations:\n", " yield obs.tstart, obs.tstop\n", "\n", "\n", "time_intervals = list(time_intervals_per_obs(observations))" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "# Assumed spectral model\n", "spectral_model = PowerLaw(\n", " index=2, amplitude=2.0e-11 * u.Unit(\"1 / (cm2 s TeV)\"), reference=1 * u.TeV\n", ")" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "energy_range = [1, 100] * u.TeV" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "CPU times: user 862 ms, sys: 20.2 ms, total: 883 ms\n", "Wall time: 885 ms\n" ] } ], "source": [ "%%time\n", "lc_estimator = LightCurveEstimator(extraction)\n", "lc = lc_estimator.light_curve(\n", " time_intervals=time_intervals,\n", " spectral_model=spectral_model,\n", " energy_range=energy_range,\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Results\n", "\n", "The light curve measurement result is stored in a table. Let's have a look at the results:" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "['time_min', 'time_max', 'flux', 'flux_err', 'flux_ul', 'is_ul', 'livetime', 'alpha', 'n_on', 'n_off', 'measured_excess', 'expected_excess']\n" ] } ], "source": [ "print(lc.table.colnames)" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "text/html": [ "Table length=4\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
time_mintime_maxfluxflux_err
1 / (cm2 s)1 / (cm2 s)
float64float64float64float64
53343.9223400925953343.941865555561.839834552434914e-111.8993229068363982e-12
53343.9542150925953343.973694259261.9975853809663366e-112.0202859435674313e-12
53345.9619812962953345.981495185182.193293401955398e-112.5526997556899865e-12
53347.91319657407453347.932710462962.3179879073655363e-112.5670181770343347e-12
" ], "text/plain": [ "\n", " time_min time_max ... flux_err \n", " ... 1 / (cm2 s) \n", " float64 float64 ... float64 \n", "------------------ ----------------- ... ----------------------\n", " 53343.92234009259 53343.94186555556 ... 1.8993229068363982e-12\n", " 53343.95421509259 53343.97369425926 ... 2.0202859435674313e-12\n", " 53345.96198129629 53345.98149518518 ... 2.5526997556899865e-12\n", "53347.913196574074 53347.93271046296 ... 2.5670181770343347e-12" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "lc.table[\"time_min\", \"time_max\", \"flux\", \"flux_err\"]" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAERCAYAAACU1LsdAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAGNlJREFUeJzt3XuYZHV95/H3R0EBQUGZeOGyo4TxEpWLo0jYJCAiCAZ1dXcwLgpB0cQLsm5WwvNsNoY8G7wSTRZhBMVsjKJAWBREWQMqEdEBR25jEIXgwLiMclXAOPLdP87pQ1N0d1UPXVXd0+/X8/QzVef8qs63z3SdT/3O5XdSVUiSBPCocRcgSZo/DAVJUsdQkCR1DAVJUsdQkCR1DAVJUmdBhkKSTyS5Lck1c/R+Fya5M8kXe6a/PckNSSrJ9nOxLEmazxZkKABnAAfN4ft9ADh8iun/DLwU+Nc5XJYkzVsLMhSq6uvA7ZOnJdml/cZ/RZJvJHnWLN7vq8A9U0z/blXd9IgLlqQFYrNxFzCHVgJvraofJNkLOBl4yZhrkqQFZZMIhSRbA78NfD7JxOTHtvP+A/AXU7zslqo6cDQVStLCsEmEAs1usDuravfeGVV1DnDO6EuSpIVnQR5T6FVVdwM3JvmPAGnsNuayJGnBWZChkOQzwGXAM5OsTXIU8HrgqCTfA64FXjmL9/sG8Hlg//b9DmynvzPJWmBH4Kokp8317yJJ80kcOluSNGFB9hQkScOx4A40b7/99rV06dJxlyFJC8oVV1zx06pa0q/dgguFpUuXsmrVqnGXIUkLSpKBRmYY2u6jJDsluTjJmiTXJjlmmnb7JlndtvnasOqRJPU3zJ7CBuDdVXVlkm2AK5JcVFXXTTRIsi3NlccHVdXNSX5jiPVIkvoYWk+hqtZV1ZXt43uANcAOPc3+ADinqm5u2902rHokSf2N5OyjJEuBPYDLe2YtA7ZLckk7kN0bpnn90UlWJVm1fv364RYrSYvY0EOhHZfobOBd7ZXHk20GvAA4BDgQ+O9JlvW+R1WtrKrlVbV8yZK+B88lSRtpqGcfJdmcJhA+3Y5B1Gst8NOq+gXwiyRfB3YDrh9mXZKkqQ3z7KMApwNrqurD0zT7P8DvJNksyVbAXjTHHiRJYzDMnsI+NHczuzrJ6nba8cDOAFV1SlWtSXIhcBXwAHBaVc3JLTYlSbM3tFCoqkuBDNDuAzS3w5QkTWHFqZcBcOZb9h76shz7SJLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQRqCFadexopTLxt3GdKsGQqSpI6hIEnqGAqSpI6hIEnqGAqSpI6hIEnqDC0UkuyU5OIka5Jcm+SYGdq+MMmvk7x2WPVIkvrbbIjvvQF4d1VdmWQb4IokF1XVdZMbJXk08D7gy0OsRZI0gKH1FKpqXVVd2T6+B1gD7DBF03cAZwO3DasWSdJgRnJMIclSYA/g8p7pOwCvBk4ZRR2SpJkNPRSSbE3TE3hXVd3dM/uvgfdU1a/7vMfRSVYlWbV+/fphlSpJi94wjymQZHOaQPh0VZ0zRZPlwGeTAGwPHJxkQ1WdO7lRVa0EVgIsX768hlmzJC1mQwuFNFv604E1VfXhqdpU1dMntT8D+GJvIEiSRmeYPYV9gMOBq5OsbqcdD+wMUFUeR5CkeWZooVBVlwKZRfsjhlWLJGkwXtEsSeoYCpKkjqEgSeoYCpKkjqEgSeoYCpKkjqEgSeoYCpKkjqEgSeoYCpK0AKy9496RLMdQkKQF4JY77x/JcgwFSZrnrlvXeyua4Rnq/RQkSRvvpIuu5yNf/UH3fOlx5wNwzP67cuwBy4ayTENBkuapYw9YxrEHLGPFqZdx+Y23c9OJhwx9me4+kiR1DAVJWgB22HaLkSzHUJCkBWDH7bYayXIMBUlSx1CQJHUMBUlSx1CQJHUGuk4hyW8A+wBPA+4DrgFWVdUDQ6xNkjRiM4ZCkv2A44AnAt8FbgO2AF4F7JLkLOBDVTW6a7AlSUPTr6dwMPDmqrq5d0aSzYBXAAcAZw+hNknSiM0YClX1JzPM2wCcO+cVSZLGZqMPNCc5ci4LkSSN3yM5++i9c1aFJGle6Heg+arpZgFPnvtyJEnj1O9A85OBA4E7eqYH+OZQKpIkjU2/UPgisHVVre6dkeSSoVQkSRqbfmcfHTXDvD+Y+3IkSeM06wPNSY4eRiGSpPHbmLOP3jrnVUiS5oWNCYXMeRWSpHlhY0Lh9wdplGSnJBcnWZPk2iTHTNHm9Umuan++mWS3jahHkjRHBgqFJMckeXySAO9NcmWSl/V52Qbg3VX1bODFwNuSPKenzY3A71XV84ETgJWzrF+SNIcG7Sn8YTsS6suAJcCRwIkzvaCq1lXVle3je4A1wA49bb5ZVRPXQHwL2HEWtUvz2to77h13CdKsDRoKE8cRDgY+WVXfYxbHFpIsBfYALp+h2VHAl6Z5/dFJViVZtX79+kEXO6MVp17GilMvm5P3kqZyy533j7sEadYGDYUrknyFJhS+nGQbYKAb7CTZmmZo7XdNd9+F9r4NRwHvmWp+Va2squVVtXzJkiUDliyNz3XrvMWIFqaB7rxGs8HeHfhRVd2b5Ek0u5BmlGRzmkD4dFWdM02b5wOnAS+vqp8NWI80L5100fV85Ks/6J4vPe58AI7Zf1eOPWDZuMqSBjZQKLS33bxy0vOfATNuwNuD0qcDa6rqw9O02Rk4Bzi8qq4ftGhpvjr2gGUce8AyVpx6GZffeDs3nXjIuEuSZmXQnsLG2Ac4HLg6ycTYSccDOwNU1SnAnwFPAk5uMoQNVbV8iDVJkmYwtFCoqkvpczC6qt4EvGlYNUjjtMO2W4y7BGnWHslNdiTNYMftthp3CdKszRgKSZ6X5FtJfpxkZZLtJs379vDLkySNUr+ewseAPweeB1wPXJpkl3be5kOsS5I0Bv2OKWxdVRe2jz+Y5ArgwiSHAzXc0iRJo9YvFJLkCVV1F0BVXZzkNTTXHjxx6NVJkkaq3+6j9wHPnjyhqq4C9qe5vkCStAnpdzvOf+idluQpVXUz8OahVSVJGouNOSX1gjmvQpI0L3jnNUlSZ2NC4eNzXoUkaV4YeJiL9sK1nYBvJdkTYOImOpKkTcNAoZDkBOAI4Ic8eH1CAS8ZTlmSpHEYtKfwn4BdqurfhlmMJGm8Bj2mcA2w7TALkSSN36A9hb8CvpvkGuCXExOr6tChVCVJGotBQ+FTNFc3X82A92aWJC08g4bCT6vqo0OtRJI0doOGwhVJ/go4j4fuPvKUVEnahAwaCnu0/7540jRPSZWkTcxAoVBV+w27EEnS+A10SmqS/5lk20nPt0vyl8MrS5I0DoNep/Dyqrpz4klV3QEcPJySJEnjMmgoPDrJYyeeJNkSeOwM7SVJC9CgB5r/Hvhqkk/SHGD+Q5prFyRJm5BBDzS/P8lVwEtp7qdwQlV9eaiVSZIAOPMte49sWTOGQpJUVQFU1YXAhTO1kSQtbP2OKVyc5B1Jdp48McljkrwkyaeANw6vPEnSKPXbfXQQzfGDzyR5OnAnsCVNmHwFOKmqVg+3REnSqMwYClV1P3AycHKSzYHtgfsmn54qSdp0DHw7zqr6FbBuiLVIksZs0OsUJEmLwKIKhRWnXsaKUy8bdxmSNG8NOvbRc6aYtu+cVyNJGqtBewqfS/KeNLZM8jc0t+iUJG1CBg2FvYCdgG8C3wFuBfaZ6QVJdkpycZI1Sa5NcswUbZLko0luSHJVkj1n+wtIkubOoGcf/Qq4j+YahS2AG6uq372aNwDvrqork2xDc/e2i6rqukltXg7s2v7sBXys/Xck1t5x76gWJUkLwqA9he/QhMILgX8PvC7JWTO9oKrWTdyus6ruAdYAO/Q0eyXwd9X4FrBtkqfO5hd4JG658/5RLUqSFoRBewpHVdWq9vFPgFcmOXzQhSRZSnNLz8t7Zu0A/HjS87XttIdcD5HkaOBogJ13fsiIGxvtunV3z8n7SNKmZNBQuK13/CPga4O8MMnWwNnAu6qqd0ucKV7ysMH1qmolsBJg+fLlj2jwvbV33MvS487vnk88Pmb/XTn2gGWP5K0lacEbNBTOp9lYh+aYwtOBfwF+a6YXtUNjnA18uqrOmaLJWpoD2BN2pDmIPTQ7brcV/3zc/qw49TIuv/F2bjrxkGEuTpIWlIGOKVTV86rq+e2/uwIvAi6d6TVJApwOrKmqD0/T7DzgDe1ZSC8G7qoqh9KQpDEZeOyjydozil7Yp9k+wOHA1UkmRlI9Hti5fY9TgAto7vV8A3AvcOTG1LOxdth2i1EuTpLmvYFCIcl/mfT0UcCewPqZXlNVlzL1MYPJbQp42yA1DMOO2201rkVL0rw0aE9hm0mPN9AcYzh77suRJI3ToPdofu+wC5EkjV+/ezR/gSlOEZ1QVYfOeUWSpLHp11P44EiqkCTNC/1C4caqunkklUiSxq7fdQrnTjxI4oFlSdrE9QuFyaeUPmOYhUiSxq/f7qOa5rGkGZz5lr3HXYK0UfqFwm5J7qbpMWzZPqZ9XlX1+KFWJ0kaqRlDoaoePapCJEnjN+hNdiRJi4ChIEnqGAqSpI6hIEnqGAqSpI6hIEnqGAqSpI6hIEnqGAqSpI6hIEnqGAqSpI6hIEnqGAqSpI6hIEnqGAqSpM6iC4W1d9w77hIkad5adKFwy533j7sESZq3FlUoXLfu7v6NJGkR63eP5k3CSRddz0e++oPu+dLjzgdgh223GFdJkjQvLYqewrEHLOOmEw9hr6c/EaB7vON2W425MkmaXxZFKEiSBrPoQsFdRpI0vUUXCu4ykqTpLbpQkCRNb2ihkOQTSW5Lcs0085+Q5AtJvpfk2iRHDqsWSdJghtlTOAM4aIb5bwOuq6rdgH2BDyV5zBDrkST1MbRQqKqvA7fP1ATYJkmArdu2G4ZVjySpv3EeU/hb4NnArcDVwDFV9cBUDZMcnWRVklXr168fZY2StKiM84rmA4HVwEuAXYCLknyjqh42FkVVrQRWAixfvrzmYuFnvmXvuXgbSdqkjLOncCRwTjVuAG4EnjXGeiRp0RtnKNwM7A+Q5MnAM4EfjbEeSVr0hrb7KMlnaM4q2j7JWuB/AJsDVNUpwAnAGUmuBgK8p6p+Oqx6JEn9DS0Uqup1febfCrxsWMuXJM2eVzRLkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSps9m4CxilM9+y97hLkKR5zZ6CJKljKEiSOoaCJKljKEiSOoaCJKljKEiSOoaCJKljKEiSOoaCJKmTqhp3DbOSZD3wr3P8ttsDP53j95wL1jU71jU71jU7C72uf1dVS/o1WnChMAxJVlXV8nHX0cu6Zse6Zse6Zmex1OXuI0lSx1CQJHUMhcbKcRcwDeuaHeuaHeuanUVRl8cUJEkdewqSpI6hIEnqLLpQSPLoJN9N8sUp5j02yZlJbkhyeZKl86SuI5KsT7K6/XnTiGq6KcnV7TJXTTE/ST7arq+rkuw5T+raN8ldk9bXn42orm2TnJXk+0nWJNm7Z/641le/uka+vpI8c9LyVie5O8m7etqMfH0NWNe4/r6OTXJtkmuSfCbJFj3z52T7tahux9k6BlgDPH6KeUcBd1TVbyY5DHgfsGIe1AVwZlW9fUS1TLZfVU13YczLgV3bn72Aj7X/jrsugG9U1StGVMuEjwAXVtVrkzwG2Kpn/rjWV7+6YMTrq6r+Bdgdmi9EwC3AP/Y0G/n6GrAuGPH6SrID8E7gOVV1X5LPAYcBZ0xqNifbr0XVU0iyI3AIcNo0TV4JfKp9fBawf5LMg7rmq1cCf1eNbwHbJnnquIsahySPB34XOB2gqv6tqu7saTby9TVgXeO2P/DDquodqWDcf1/T1TUumwFbJtmMJthv7Zk/J9uvRRUKwF8D/w14YJr5OwA/BqiqDcBdwJPmQV0Ar2m70Gcl2WkENQEU8JUkVyQ5eor53fpqrW2njbsugL2TfC/Jl5L81ghqegawHvhkuxvwtCSP62kzjvU1SF0w+vU12WHAZ6aYPq6/rwnT1QUjXl9VdQvwQeBmYB1wV1V9pafZnGy/Fk0oJHkFcFtVXTFTsymmDfWc3QHr+gKwtKqeD/xfHvw2MGz7VNWeNN34tyX53Z75I19frX51XUkzzstuwN8A546gps2APYGPVdUewC+A43rajGN9DVLXONYXAO3urEOBz081e4ppIzmHvk9dI19fSbaj6Qk8HXga8Lgk/7m32RQvnfX6WjShAOwDHJrkJuCzwEuS/H1Pm7XATgBtF+0JwO3jrquqflZVv2yffhx4wZBrmljure2/t9HsV31RT5NufbV25OFd2pHXVVV3V9XP28cXAJsn2X7IZa0F1lbV5e3zs2g2xr1tRr2++tY1pvU14eXAlVX1/6aYN5a/r9a0dY1pfb0UuLGq1lfVr4BzgN/uaTMn269FEwpV9adVtWNVLaXpFv5TVfUm7XnAG9vHr23bDPWbySB19exHPZTmgPRQJXlckm0mHgMvA67paXYe8Ib2LJEX03Rp1427riRPmdiXmuRFNH/nPxtmXVX1E+DHSZ7ZTtofuK6n2cjX1yB1jWN9TfI6pt9FM/L1NUhdY1pfNwMvTrJVu+z9efh2YE62X4vx7KOHSPIXwKqqOo/mYNz/TnIDTcIeNk/qemeSQ4ENbV1HjKCEJwP/2P7tbwb8Q1VdmOStAFV1CnABcDBwA3AvcOQ8qeu1wB8l2QDcBxw27HBvvQP4dLvr4UfAkfNgfQ1S11jWV5KtgAOAt0yaNvb1NUBdI19fVXV5krNodl1tAL4LrBzG9sthLiRJnUWz+0iS1J+hIEnqGAqSpI6hIEnqGAqS9Agl+fMkt+TBQfIOnqLNFkm+3V4JfW2S906ad3o7fWLUgq17XvvaJJVkoHsxJ3l8W8/fzvZ3MRS0SUvypEkf1J/0fHC/OaRl7pHktPbxEe2Hef9J81/dTntt+/ySiQ97HhwB9uok1yX5yySPbectSXLhMGrW4NKMknrGFLNOqqrd258Lppj/S+Al7ZXQuwMHtddfABxbVbu1oxbcDHSDX7bX5bwTuLz3DWdwAvC1WbTvGArapLVXg+9eVbsDp/DQD27vFaFz5Xia4Q8mXE1zMdSEw4DvzfD6/arqeTRXaj+D9naLVbUeWJdkn7ktV6PQDuz38/bp5u1PtfPuhma4cGBLHjo8xQnA+4H7JyakGWr/A0m+0/YuJl9T8QKa63l6x0YaiKGgRSvJz9t/903ytSSfS3J9khOTvL7t6l+dZJe23ZIkZ7cfxO9MtXFuv9U9v6omb/S/AbwoyebtboHfBFb3q6/dgLwVeFWSJ7aTzwVe/4h+cQ3L29sN9CfSjFX0MO3GfDVwG3DRpOFHSPJJ4CfAs2i/VCTZA9ipqnrvs3IUzRXeLwReCLw5ydOTPAr4EPAnG/tLGApSYzeae1o8DzgcWFZVL6IZzvwdbZuP0PQ0Xgi8hqmHOl/Ow4cDKZqBDA+kGdTsvEGLar9B3khzTwGAVcDvDPp6zZ00N65ZTfP/fuik3ZAH0tzrYRea3ULraDbMD1NVv257rTvSfFF47qR5R9IMdrcGWNFu4E8C3j3FW72MZgiQ1TS7lZ5E8zfyx8AFVfXjKV4zkEU/zIXU+s7EuDpJfsiDXe+rgf3axy8FnpMHh6h/fJJtquqeSe/zVJqhqnt9lma/8BNoPuTHz6K2yaNf3kaz4dCIVdVe0PQsgSOq6oip2iX5OPCwOyj2vNedSS4BDmLSl4iq+nWSM2m+6Z8DPBe4pP2bewpwXjvkTYB3VNWXe5b9RuB3kvwxsDXwmCQ/r6rekXGnZShIjV9OevzApOcP8ODn5FHA3lV13wzvcx+wRe/Eqvp2+63wvqq6PgPe+6TdHbUUuL6dtEW7DM0jSZ46abC+V/Pw3iJJlgC/agNhS5ovGe9rjyPsUlU3tI9/H/h+Vd0FbD/p9ZcA/7WqViX5Ms34S/9UVb9Ksgy4papeP6n9EcDy2QQCGArSbHyF5qyQDwAk2b2qeo8NrGHq7j7AnzLpYGE/7fGHk4Fzq+qOdvIyptjgaOzen2R3ml2FN9EOppfkacBpVXUwTS/yU2lu8/ko4HNV9cV2N9Gn0twlLzQnIfxRn+WdRvNl4co2SNYDr5qLX8RQkAb3TuB/JbmK5rPzdZoDwZ2q+n6SJ0yxW4mq+tI077sZD+2pXNx+0B9Fc7+IEybN2w84/5H9GnokquoS4JKeaYdP0/ZWmpFeqaqrgD2maPMAzX1V+i13357XHM8MuyGr6gweeg/ngThKqjTHkhwL3FNVfe+53V6DcAPw3HZ3Qb/2XwdeOannIM0pzz6S5t7HeOg3/ym1F6ytBk4eMBCWAB82EDRM9hQkSR17CpKkjqEgSeoYCpKkjqEgSeoYCpKkzv8HUnb63axJJSkAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "lc.plot();" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$$2.0740149 \\times 10^{-11} \\; \\mathrm{\\frac{1}{s\\,cm^{2}}}$$" ], "text/plain": [ "" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Let's compare to the expected flux of this source\n", "from gammapy.spectrum import CrabSpectrum\n", "\n", "crab_spec = CrabSpectrum().model\n", "crab_flux = crab_spec.integral(*energy_range).to(\"cm-2 s-1\")\n", "crab_flux" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAERCAYAAACU1LsdAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAGX5JREFUeJzt3X20XXV95/H3R4gaC4oKVQk4sVRRK0I0PrSZdhCsUbCIDzPVOiiUlqG1Ch2HMbhmOrXMGlFa29qOuiIqdmp9gjSlYkGngo8jNUAkYIpDhSIBh2iNYJtaHr7zx953c3O4uefccPc59+a+X2vdlX32/p2zv3fn3PM5v/3w26kqJEkCeMikC5AkLRyGgiSpYyhIkjqGgiSpYyhIkjqGgiSpsyhDIckHk9yR5Lp5er1Lk+xI8qmB+b+R5MYkleTA+ViXJC1kizIUgAuAF8/j650HnDTD/C8DLwT+fh7XJUkL1qIMhar6AvAP0+clOaz9xn9Vki8meeocXu+vgbtmmH9NVd38oAuWpEVi30kXMI/WA6dX1f9N8jzgPcAxE65JkhaVvSIUkuwH/AzwySRTsx/WLnsF8DszPG1bVa0dT4WStDjsFaFAsxtsR1UdNbigqjYAG8ZfkiQtPovymMKgqroTuCnJvwVI48gJlyVJi86iDIUkHwX+D3B4kluTnAq8Fjg1ydeB64GXzeH1vgh8Eji2fb217fw3JbkVOAS4Nsn58/27SNJCEofOliRNWZQ9BUlSPxbdgeYDDzywVq5cOekyJGlRueqqq75bVQcNa7foQmHlypVs2rRp0mVI0qKSZKSRGXrbfZTk0CSXJ9ma5PokZ+ym3dFJNrdtPt9XPZKk4frsKdwDvLmqrk6yP3BVks9W1TemGiQ5gObK4xdX1S1JfrzHeiRJQ/TWU6iq26vq6nb6LmArsGKg2S8BG6rqlrbdHX3VI0kabixnHyVZCawCrhxY9BTg0UmuaAeye91unn9akk1JNm3fvr3fYiVpCes9FNpxiS4CzmyvPJ5uX+DZwPHAWuC/JnnK4GtU1fqqWl1Vqw86aOjBc0nSHur17KMky2gC4SPtGESDbgW+W1X/CPxjki8ARwLf7LMuSdLM+jz7KMAHgK1V9a7dNPsL4GeT7JvkEcDzaI49SJImoM+ewhqau5ltSbK5nfdW4IkAVfW+qtqa5FLgWuA+4PyqmpdbbEqS5q63UKiqLwEZod15NLfDlCTNYOW6SwC4+dzje1+XYx9JkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIkjqGgiSpYyhIPVi57hJWrrtk0mVIc2YoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqdNbKCQ5NMnlSbYmuT7JGbO0fU6Se5O8qq96JEnD7dvja98DvLmqrk6yP3BVks9W1TemN0qyD/AO4LIea5EkjaC3nkJV3V5VV7fTdwFbgRUzNH0jcBFwR1+1SJJGM5ZjCklWAquAKwfmrwBeDrxvHHVIkmbXeygk2Y+mJ3BmVd05sPgPgLdU1b1DXuO0JJuSbNq+fXtfpUrSktfnMQWSLKMJhI9U1YYZmqwGPpYE4EDguCT3VNXG6Y2qaj2wHmD16tXVZ82StJT1FgppPuk/AGytqnfN1KaqnjSt/QXApwYDQZI0Pn32FNYAJwFbkmxu570VeCJAVXkcQZIWmN5Coaq+BGQO7U/uqxZJ0mi8olmS1DEUJEkdQ0GS1DEUJEkdQ0GS1DEUJEkdQ0GS1DEUJEkdQ0GS1DEUJGkB23jNtm56zbmf2+VxHwwFSVqgNl6zjbM3bOkeb9uxk7M3bOk1GAwFSVqgzrvsBnbevevtZnbefS/nXXZDb+s0FCRpgbptx845zZ8PhoIkLVAHH7B8TvPng6EgSQvUWWsPZ/myfXaZt3zZPpy19vDe1tnr7TglSXvuxFUrADjz4819ylYcsJyz1h7eze+DPQVJWsCmB8CX1x3TayCAoSBJmsZQkCR1DAVJUsdQkCR1Rjr7KMmPA2uAg4GdwHXApqq6r8faJEljNmsoJHkBsA54DHANcAfwcOBE4LAkFwK/V1V39l2oJKl/w3oKxwG/WlW3DC5Isi/wUuDngYt6qE2SNGazhkJVnTXLsnuAjfNekSRpYvb4QHOSU+azEEnS5D2Ys4/eNm9VSJIWhGEHmq/d3SLgcfNfjiRpkoYdaH4csBb4/sD8AF/ppSJJ0sQMC4VPAftV1ebBBUmu6KUiSdLEDDv76NRZlv3S/JcjSZqkOR9oTnJaH4VIkiZvT84+On3eq5AkLQh7EgqZ9yokSQvCnoTCL4zSKMmhSS5PsjXJ9UnOmKHNa5Nc2/58JcmRe1CPJGmejBQKSc5I8sgkAd6W5OokLxrytHuAN1fV04DnA29I8vSBNjcB/6aqngmcA6yfY/2SpHk0ak/hl9uRUF8EHAScApw72xOq6vaqurqdvgvYCqwYaPOVqpq6BuKrwCFzqF1akDZes62bXnPu53Z5LC10o4bC1HGE44APVdXXmcOxhSQrgVXAlbM0OxX4q908/7Qkm5Js2r59+6irndXKdZewct0l8/Ja0pSN12zj7A1busfbduzk7A1bDAYtGqOGwlVJPkMTCpcl2R8Y6QY7SfajGVr7zN3dd6G9b8OpwFtmWl5V66tqdVWtPuigg0YsWRq/8y67gZ1337vLvJ1338t5l90woYqkuRnpzms0H9hHAd+qqn9K8liaXUizSrKMJhA+UlUbdtPmmcD5wEuq6nsj1iMtSLft2Dmn+dJCM1JPoaruq6qrq2pH+/h7VbW7wfIAaA9KfwDYWlXv2k2bJwIbgJOq6ptzK11aeA4+YPmc5ksLzYMZOnuYNcBJwDFJNrc/xyU5PcnUBXC/BTwWeE+7fFOP9Ui9O2vt4Sxfts8u85Yv24ez1h4+oYqkuRl199GcVdWXGHIwuqp+BfiVvmqQxu3EVc0Jdmd+vBlDcsUByzlr7eHdfGmh67OnIC1J0wPgy+uOMRC0qMwaCkmOSPLVJN9Osj7Jo6ct+5v+y5MkjdOwnsJ7gd8GjgC+CXwpyWHtsmU91iVJmoBhxxT2q6pL2+nfTXIVcGmSk4DqtzRJ0rgNC4UkeVRV/QCgqi5P8kqaaw8e03t1kqSxGrb76B3A06bPaK9POJbm+gJJ0l5k2O04/2xwXpLHV9UtwK/2VpUkaSL25JTUT897FZKkBcE7r0mSOnsSCu+f9yokSQvCyMNctBeuHQp8NcmzAKZuoiNJ2juMFApJzgFOBv6O+69PKOCYfsqSJE3CqD2FfwccVlX/0mcxkqTJGvWYwnXAAX0WIkmavFF7Cm8HrklyHfCjqZlVdUIvVUmSJmLUUPgwzdXNWxjx3sySpMVn1FD4blW9u9dKJEkTN2ooXJXk7cDF7Lr7yFNSJWkvMmoorGr/ff60eZ6SKkl7mZFCoape0HchkqTJG/Xitf8BvLOqdrSPHw28uar+S5/Fzbejjz66m/7Ot77XzPvqeROqRnsz31+aT1PvJ849vvd1jXqdwkumAgGgqr4PHNdPSZKkSRn1mMI+SR5WVT8CSLIceFh/ZfXjiiuu6KZXrrukmTeG5NXS4/tL82nq/TQOo4bCnwJ/neRDNAeYf5nm2gVJ0l5k1APN70xyLfBCmvspnFNVl/VamSQJgJvH2OOcNRSSpKoKoKouBS6drY0kaXEbdqD58iRvTPLE6TOTPDTJMUk+DLy+v/IkSeM0bPfRi2mOH3w0yZOAHcBymjD5DPD7VbW53xIlSeMyayhU1T8D7wHek2QZcCCwc/rpqZKkvcfIt+OsqruB23usRZI0YaNevCZJWgKWVCisXHfJWC8CkaTFZqRQSPL0GeYdPe/VSJImatSewieSvCWN5Un+iOYWnZKkvcioofA84FDgK8DXgNuANbM9IcmhSS5PsjXJ9UnOmKFNkrw7yY1Jrk3yrLn+ApKk+TNqKNwN7KS5RuHhwE1VNexezffQDK/9NJqb87xhht1QLwGe3P6cBrx31MIfjI3XbOum15z7uV0eS9JSNmoofI0mFJ4D/GvgNUkunO0JVXX71O06q+ouYCuwYqDZy4A/qcZXgQOSPGEuv8CeOHvDlm56246dnL1hi8EgSYweCqdW1W9V1d1V9Z2qehnwF6OuJMlKmlt6XjmwaAXw7WmPb+WBwUGS05JsSrJp+/bto652t3befe8DHp932Q0P+nUlabEb9eK1OwbHPwI+P8oTk+wHXAScWVV3Di6e4SkPGFyvqtYD6wFWr17dy+B7t+3Y2cfLStKiMmooXELzYR2aYwpPAm4Afmq2J7VDY1wEfKSqNszQ5FaaA9hTDqE5iD12Bx+wfBKrlaQFZaTdR1V1RFU9s/33ycBzgS/N9pwkAT4AbK2qd+2m2cXA69qzkJ4P/KCqeh9KY/myfR7w+Ky1h/e9Wkla8Pboiub2APJzhjRbA5wEHJNkc/tzXJLTk5zetvk08C3gRuD9wK/vST1z9fZXHNFNrzhgOW9/xRGcuOoBhzIkackZafdRkv847eFDgGcBsx7xraovMfMxg+ltCnjDKDXMpxNXreDMjzcjfn953THjXr0kLVijHlPYf9r0PTTHGC6a/3IkSZM06j2a39Z3IZKkyRt2j+a/ZIZTRKdU1QnzXpEkaWKG9RR+dyxVSJIWhGGhcFNV3TKWSiRJEzfslNSNUxNJPLAsSXu5YaEw/ZTSn+izEEnS5A3bfVS7mZY0i5vPPX7SJUh7ZFgoHJnkTpoew/J2mvZxVdUje61OkjRWs4ZCVe0z23JJ0t5lj8Y+kiTtnQwFSVLHUJAkdQwFSVLHUJAkdQwFSVLHUJAkdQwFSVLHUJAkdQwFSVLHUJAkdQwFSVLHUJAkdQwFSVLHUJAkdZZMKGy8Zls3vebcz02wEklauJZEKGy8Zhtnb9jSPd62Y+cEq5GkhWtJhMJ5l93AzrvvnXQZkrTgLYlQuM2egSSNZEmEwsEHLJ90CZK0KCyJUDhr7eEsX7bPpMuQpAVvSYTCiatW8PZXHNE9XmHPQZJmtCRCAZpgmPLldcdMsBJJWriWTChIkobrLRSSfDDJHUmu283yRyX5yyRfT3J9klP6qkWSNJo+ewoXAC+eZfkbgG9U1ZHA0cDvJXloj/VIkoboLRSq6gvAP8zWBNg/SYD92rb39FWPJGm4SR5T+GPgacBtwBbgjKq6b6aGSU5LsinJpu3bt4+zRklaUvad4LrXApuBY4DDgM8m+WJV3TnYsKrWA+sBVq9eXfOx8pvPPX4+XkaS9iqT7CmcAmyoxo3ATcBTJ1iPJC15kwyFW4BjAZI8Djgc+NYE65GkJa+33UdJPkpzVtGBSW4F/huwDKCq3gecA1yQZAsQ4C1V9d2+6pEkDddbKFTVa4Ysvw14UV/rlyTNnVc0S5I6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqWMoSJI6hoIkqbPvpAsYp5vPPX7SJUjSgmZPQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUMRQkSR1DQZLUSVVNuoY5SbId+Pt5ftkDge/O82vOB+uaG+uaG+uam8Ve17+qqoOGNVp0odCHJJuqavWk6xhkXXNjXXNjXXOzVOpy95EkqWMoSJI6hkJj/aQL2A3rmhvrmhvrmpslUZfHFCRJHXsKkqSOoSBJ6iy5UEiyT5JrknxqhmUPS/LxJDcmuTLJygVS18lJtifZ3P78yphqujnJlnadm2ZYniTvbrfXtUmetUDqOjrJD6Ztr98aU10HJLkwyd8m2ZrkpweWT2p7Datr7NsryeHT1rc5yZ1JzhxoM/btNWJdk3p//WaS65Ncl+SjSR4+sHxePr+W1O04W2cAW4FHzrDsVOD7VfWTSV4NvAP4xQVQF8DHq+o3xlTLdC+oqt1dGPMS4Mntz/OA97b/TrougC9W1UvHVMuUPwQurapXJXko8IiB5ZPaXsPqgjFvr6q6ATgKmi9EwDbgzweajX17jVgXjHl7JVkBvAl4elXtTPIJ4NXABdOazcvn15LqKSQ5BDgeOH83TV4GfLidvhA4NkkWQF0L1cuAP6nGV4EDkjxh0kVNQpJHAj8HfACgqv6lqnYMNBv79hqxrkk7Fvi7qhocqWDS76/d1TUp+wLLk+xLE+y3DSyfl8+vJRUKwB8A/xm4bzfLVwDfBqiqe4AfAI9dAHUBvLLtQl+Y5NAx1ARQwGeSXJXktBmWd9urdWs7b9J1Afx0kq8n+askPzWGmn4C2A58qN0NeH6SHxtoM4ntNUpdMP7tNd2rgY/OMH9S768pu6sLxry9qmob8LvALcDtwA+q6jMDzebl82vJhEKSlwJ3VNVVszWbYV6v5+yOWNdfAiur6pnA/+b+bwN9W1NVz6Lpxr8hyc8NLB/79moNq+tqmnFejgT+CNg4hpr2BZ4FvLeqVgH/CKwbaDOJ7TVKXZPYXgC0u7NOAD450+IZ5o3lHPohdY19eyV5NE1P4EnAwcCPJfn3g81meOqct9eSCQVgDXBCkpuBjwHHJPnTgTa3AocCtF20RwH/MOm6qup7VfWj9uH7gWf3XNPUem9r/72DZr/qcweadNurdQgP7NKOva6qurOqfthOfxpYluTAnsu6Fbi1qq5sH19I82E82Gbc22toXRPaXlNeAlxdVf9vhmUTeX+1dlvXhLbXC4Gbqmp7Vd0NbAB+ZqDNvHx+LZlQqKqzq+qQqlpJ0y38XFUNJu3FwOvb6Ve1bXr9ZjJKXQP7UU+gOSDdqyQ/lmT/qWngRcB1A80uBl7XniXyfJou7e2TrivJ46f2pSZ5Ls37/Ht91lVV3wG+neTwdtaxwDcGmo19e41S1yS21zSvYfe7aMa+vUapa0Lb6xbg+Uke0a77WB74OTAvn19L8eyjXST5HWBTVV1MczDufyW5kSZhX71A6npTkhOAe9q6Th5DCY8D/rx97+8L/FlVXZrkdICqeh/waeA44Ebgn4BTFkhdrwJ+Lck9wE7g1X2He+uNwEfaXQ/fAk5ZANtrlLomsr2SPAL4eeA/TJs38e01Ql1j315VdWWSC2l2Xd0DXAOs7+Pzy2EuJEmdJbP7SJI0nKEgSeoYCpKkjqEgSeoYCpL0ICX57STbcv8gecfN0ObhSf6mvRL6+iRvm7bsA+38qVEL9ht47quSVJKR7sWc5JFtPX8819/FUNBeLcljp/2hfmfgD/crPa1zVZLz2+mT2z/mY6ctf3k771Xt4yum/thz/wiwW5J8I8l/T/KwdtlBSS7to2aNLs0oqRfMsOj3q+qo9ufTMyz/EXBMeyX0UcCL2+svAH6zqo5sRy24BegGv2yvy3kTcOXgC87iHODzc2jfMRS0V2uvBj+qqo4C3seuf7iDV4TOl7fSDH8wZQvNxVBTXg18fZbnv6CqjqC5UvsnaG+3WFXbgduTrJnfcjUO7cB+P2wfLmt/ql12JzTDhQPL2XV4inOAdwL/PDUjzVD75yX5Wtu7mH5NxbNprucZHBtpJIaClqwkP2z/PTrJ55N8Isk3k5yb5LVtV39LksPadgcluaj9Q/zaTB/O7be6Z1bV9A/9LwLPTbKs3S3wk8DmYfW1HyCnAycmeUw7eyPw2gf1i6svv9F+QH8wzVhFD9B+mG8G7gA+O234EZJ8CPgO8FTaLxVJVgGHVtXgfVZOpbnC+znAc4BfTfKkJA8Bfg84a09/CUNBahxJc0+LI4CTgKdU1XNphjN/Y9vmD2l6Gs8BXsnMQ52v5oHDgRTNQIZraQY1u3jUotpvkDfR3FMAYBPws6M+X/MnzY1rNtP8v58wbTfkWpp7PRxGs1vodpoP5geoqnvbXushNF8UnjFt2Sk0g91tBX6x/YD/feDNM7zUi2iGANlMs1vpsTTvkV8HPl1V357hOSNZ8sNcSK2vTY2rk+TvuL/rvQV4QTv9QuDpuX+I+kcm2b+q7pr2Ok+gGap60Mdo9gs/iuaP/K1zqG366Jd30HxwaMyq6nnQ9CyBk6vq5JnaJXk/8IA7KA681o4kVwAvZtqXiKq6N8nHab7pbwCeAVzRvuceD1zcDnkT4I1VddnAul8P/GySXwf2Ax6a5IdVNTgy7m4ZClLjR9Om75v2+D7u/zt5CPDTVbVzltfZCTx8cGZV/U37rXBnVX0zI977pN0dtRL4Zjvr4e06tIAkecK0wfpezgN7iyQ5CLi7DYTlNF8y3tEeRzisqm5sp38B+Nuq+gFw4LTnXwH8p6ralOQymvGXPldVdyd5CrCtql47rf3JwOq5BAIYCtJcfIbmrJDzAJIcVVWDxwa2MnN3H+Bsph0sHKY9/vAeYGNVfb+d/RRm+MDRxL0zyVE0uwpvph1ML8nBwPlVdRxNL/LDaW7z+RDgE1X1qXY30YfT3CUvNCch/NqQ9Z1P82Xh6jZItgMnzscvYihIo3sT8D+TXEvzt/MFmgPBnar62ySPmmG3ElX1V7t53X3ZtadyefuH/hCa+0WcM23ZC4BLHtyvoQejqq4ArhiYd9Ju2t5GM9IrVXUtsGqGNvfR3Fdl2HqPHnjOW5llN2RVXcCu93AeiaOkSvMsyW8Cd1XV0Htut9cg3Ag8o91dMKz9F4CXTes5SPPKs4+k+fdedv3mP6P2grXNwHtGDISDgHcZCOqTPQVJUseegiSpYyhIkjqGgiSpYyhIkjqGgiSp8/8Bg4xaWw7XMzoAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "ax = lc.plot(marker=\"o\", lw=2)\n", "ax.hlines(\n", " crab_flux.value,\n", " xmin=lc.table[\"time_min\"].min(),\n", " xmax=lc.table[\"time_max\"].max(),\n", ");" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Exercises\n", "\n", "* Change the assumed spectral model shape (e.g. to a steeper power-law), and see how the integral flux estimate for the lightcurve changes.\n", "* Try a time binning where you split the observation time for every run into two time bins.\n", "* Try to analyse the PKS 2155 flare data from the H.E.S.S. first public test data release.\n", " Start with per-observation fluxes, and then try fluxes within 5 minute time bins for one or two of the observations where the source was very bright." ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [], "source": [] } ], "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.6.0" }, "nbsphinx": { "orphan": true } }, "nbformat": 4, "nbformat_minor": 2 }