{ "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.11?urlpath=lab/tree/spectrum_pipe.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", "[spectrum_pipe.ipynb](../_static/notebooks/spectrum_pipe.ipynb) |\n", "[spectrum_pipe.py](../_static/notebooks/spectrum_pipe.py)\n", "
\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Spectrum analysis with Gammapy (run pipeline)\n", "\n", "In this tutorial we will learn how to perform a 1d spectral analysis.\n", "\n", "We will use a \"pipeline\" or \"workflow\" class to run a standard analysis. If you're interested in implementation detail of the analysis in order to create a custom analysis class, you should read ([spectrum_analysis.ipynb](spectrum_analysis.ipynb)) that executes the analysis using lower-level classes and methods in Gammapy. \n", "\n", "In this tutorial we will use the folling Gammapy classes:\n", "\n", "- [gammapy.data.DataStore](..\/api/gammapy.data.DataStore.rst) to load the data to \n", "- [gammapy.scripts.SpectrumAnalysisIACT](..\/api/gammapy.scripts.SpectrumAnalysisIACT.rst) to run the analysis\n", "\n", "We use 4 Crab observations from H.E.S.S. for testing." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Setup\n", "\n", "As usual, we'll start with some setup for the notebook, and import the functionality we need." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", "import astropy.units as u\n", "from astropy.coordinates import SkyCoord\n", "from regions import CircleSkyRegion\n", "\n", "from gammapy.utils.energy import EnergyBounds\n", "from gammapy.data import DataStore\n", "from gammapy.scripts import SpectrumAnalysisIACT\n", "from gammapy.catalog import SourceCatalogGammaCat\n", "from gammapy.maps import Map\n", "from gammapy.spectrum.models import LogParabola\n", "from gammapy.spectrum import CrabSpectrum" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Select data\n", "\n", "First, we select and load some H.E.S.S. data (simulated events for now). In real life you would do something fancy here, or just use the list of observations someone send you (and hope they have done something fancy before). We'll just use the standard gammapy 4 crab runs." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[23523 23526 23559 23592]\n" ] } ], "source": [ "data_store = DataStore.from_dir(\"$GAMMAPY_DATA/hess-dl3-dr1/\")\n", "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(obs_ids)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Configure the analysis\n", "\n", "Now we'll define the input for the spectrum analysis. It will be done the python way, i.e. by creating a config dict containing python objects. We plan to add also the convenience to configure the analysis using a plain text config file." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "crab_pos = SkyCoord.from_name(\"crab\")\n", "on_region = CircleSkyRegion(crab_pos, 0.15 * u.deg)\n", "\n", "model = LogParabola(\n", " alpha=2.3,\n", " beta=0.01,\n", " amplitude=1e-11 * u.Unit(\"cm-2 s-1 TeV-1\"),\n", " reference=1 * u.TeV,\n", ")\n", "\n", "flux_point_binning = EnergyBounds.equal_log_spacing(0.7, 30, 5, u.TeV)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQ8AAAEHCAYAAACwfMNTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAE7hJREFUeJzt3X+QXWV9x/H3h034kWyALAGkIXUvtAwEqkSDUmOJoDJU8FdVBhSEarG1iERjtUVHZaQdBUVtp3XsyM8SRTQRkVEhpikYMWkDhIQYlMpGi2QEgiYBY3DJt3+c55qbzWb37LN7du/Z/bxm7txzn3Puud8nyf3k/LjnOYoIzMyGap+xLsDM6snhMQIkaaxrGI661w/170Md63d4jIynx7qAYap7/VD/PtSufoeHmWVxeJhZFtX5bMuMGTOiu7t7rMugp6eHRqMx1mVkq3v9UP8+tHP9PT09sXnz5j02NCaNRTEjpbu7m9WrV491GWbj2ty5c/s9mOvdFjPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPL4vAwsyyVhoekWZKWS9ogab2kS1vmXSLpx6n9ypb2qyStljS/ytrMbHiqHs+jF1gYEfdJmgbcK2kpcDjweuAFEbFD0mEAko5N7zsFuB64q+L6zCxTpeEREZuATWl6m6QNwEzgIuCTEbEjzXs8vaUD2AkEULvRpM0mklE75iGpG5gDrAKOAf5M0ipJd0k6CSAi1gNTgBXAF/aynsmjUrCZ7aH1+zcqwxBK6gQWAwsiYqukScB04GTgJOAWSUdF4ZIB1iPgV817XHR1dY1C9WYTW09PD5KeSS9D0rSIiMq3PFJSLQYWRcSS1PwosCSFxX9T7KrMGGxdafnOiJgaEVPbdcBYs/Gk0WjQ/M6l719A9WdbBFwDbIiIq1tm3QqclpY5BtgXeLLKWsxsZFW92zIPOB9YJ2lNarsMuBa4VtKDwLPABc00M7N6qPpsywr2ftbkvCo/28yq5V+YmlkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4mFmWSsND0ixJyyVtkLRe0qV95n9AUkiakV7vI+lGSfdIOr7K2sxseKre8ugFFkbEccDJwMWSZkMRLMCrgZ+3LH86sAp4I7Cw4trMbBgqDY+I2BQR96XpbcAGYGaa/Vngg0C0vKUD2JkeqrI2MxueUTvmIakbmAOskvQ64BcR8UCfxe4A5gO3AVfvZT0dFZZpZgOQNLk5PWmUPrATWAwsoNiV+TDFLspuIqIXOGeA9QjYkp7p6uqqpF4z26WnpwdJz6SXIWlaRETlWx4pqRYDiyJiCXA00AAekLQROBK4T9LzBltXFDojYmpETG00GlWWbmZAo9Gg+Z1L37+Airc80hbCNcCGiLgaICLWAYe1LLMRmBsRT1ZZi5mNrKq3POYB5wOnSVqTHq+p+DPNbBRUuuURESsY5KxJRHRXWYOZVcO/MDWzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLKXG85A0D/g48Pz0HlGMCnhUdaWZWTsrOxjQNcD7gHuB56orx8zqomx4bImI71RaiZnVStnwWC7pKmAJsKPZ2Lyhk5lNPGXD46XpeW5LWwCnjWw5ZlYXpcIjIk6tuhAzq5dSp2olHSTpakmr0+Mzkg6qujgza19lf+dxLbANODs9tgLXVVWUmbW/ssc8jo6IN7W8vlzSmioKMrN6KLvlsV3Sy5sv0o/GtldTkpnVQdktj3cDN6TjHAKeAi6sqigza39lz7asAV4o6cD0emulVZlZ2xswPCSdFxE3SXp/n3YAIt353swmnsG2PKam52n9zIsRrsXMamTA8IiIL6bJ70XED1rnpYOmZjZBlT3b8i8l28xsghjsmMefAi8DDu1z3ONAoKPKwsysvQ12zGNfoDMt13rcYyvw5qqKMrP2N9gxj7uAuyRdHxE/G6WazKwGyh7z+I2kqyR9W9J/Nh+DvUnSLEnLJW2QtF7Span9KkkPSVor6RuSDm55z1Xp4rv5mX0ys1FQNjwWAQ8BDeByYCPwPyXe1wssjIjjgJOBiyXNBpYCJ0TEC4CfAP8AIOnY9L5TgItL1mZmY6BseBwSEdcAv4uIuyLiHRRhMKCI2NQcbSwitgEbgJkRcWdE9KbFVgJHpukOYCfFb0g0hH6Y2SgrGx6/S8+bJJ0paQ67vvClSOoG5gCr+sx6B/AdgIhYD0wBVgBf2Mt6fJanhiTt9rB6av3+lb0w7op0UdxCit93HEgxmnrZD+wEFgMLWq+LkfRhil2bRc22iLhkgPUI2JKe6erqKluCjaH+wkISEf6Rch309PQg6Zn0MiRNi4goe2Hc7WlyCzCkIQklTaYIjkURsaSl/QLgLOCVUfJfUVqus/l67ty5/tdnVrFGo8HmzZun9m0ve9OnQ4GLgO7W96RjHwO9TxT3fNnQehGdpDOADwHzI+I3ZWows/ZSdrflm8D3ge8xtJs+zQPOB9a1jDx2GfDPwH7A0rRJuzIi/mYI6zWzMVY2PKZExIeGuvKIWEH/Z02+PdR1WX1FxB7HPXy8o/7Knm25XdJrKq3ExrWI2O1h9Vc2PC6lCJDtkrZK2ibJo4mZTWBlz7b0NxiQmU1gg12Sf2xEPCTpRf3N971qzSauwbY8FlKcov1MP/N8r1qzCWywS/IvSs++V62Z7Waw3Za/GGh+6y9GzWxiGWy35bUDzAvA4WE2QQ222/KXo1WImdVLqd95SPqnPqN9TZd0RXVlmVm7K/sjsT+PiF83X0TErwD/4tRsAisbHh2S9mu+kHQAxYVtZlZC38GQxsOASGUvjLsJWCbpOooDpe8AbqisKrMJoO4DIpX9efqVktYCr6K4SvYTEXFHpZWZjRPjYSujP2W3PKAYvLg3Ir4naUoaimxbVYWZjRf9DUkwHpQ923IR8HWgeePrmcCtVRVlZu2v7AHTiylGBdsKEBEPA4dVVZTZeNPfsY06H++A8rstOyLi2eaml6RJFAdOzaykuodFX2W3PO6SdBlwgKRXA18DvlVdWWbW7sqGx98DTwDrgL+mGIP0I1UVZWbtr+yp2p2SbgVujYgnKq7JzGpgwC0PFT4u6UmKG13/WNITkj46OuWZWbsabLdlAcVZlpMi4pCI6AJeCsyTVPp2k2Y2/gwWHm8Hzo2InmZDRDwCnJfmmdkENVh4TI6IJ/s2puMek6spyczqYLDweDZznpmNc4OdbXnhXm7uJGD/Cuoxs5oYbBjCjtEqxMzqpeyPxMzMduPwMLMsDg8zy+LwMLMsDg8zy+LwMLMsDg8zy1JpeEiaJWm5pA2S1ku6NLV3SVoq6eH0PD217yPpRkn3SDq+ytrMbHiq3vLoBRZGxHHAycDFkmZTDC60LCL+GFiWXgOcDqwC3ggsrLg2MxuGSsMjIjZFxH1pehvF7RtmAq9n102jbgDekKY7gJ3pMf7GqjcbR0btmIekbmAOxZbF4RGxCYqAYddI7HcA84HbgKv3sh7/ZN5sjLR+/4Zy06fhfGAnsBhYEBFb93YDnIjoBc4ZYD0CtqRnurq6KqjWzFr19PQg6Zn0MtIN36LyLQ9JkymCY1FELEnNv5R0RJp/BPB4mXVFoTMipkbE1EajUU3RZvZ7jUaD5ncuff8Cqj/bIuAaYENEtO6G3AZckKYvAL5ZZR1mNvKq3m2ZB5wPrJO0JrVdBnwSuEXSO4GfA2+puA4zG2GVhkdErGDvZ01eWeVnm1m1/AtTM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLA4PM8vi8DCzLJWGh6RrJT0u6cGWthMlrZS0RtJqSS9J7ftIulHSPZKOr7IuMxu+qrc8rgfO6NN2JXB5RJwIfDS9BjgdWAW8EVhYcV1mNkyVhkdE3A081bcZODBNHwQ8lqY7gJ3poSrrMrPhmzQGn7kAuEPSpynC62Wp/Q7gJuDtwLvGoC4zG4KxCI93A++LiMWSzgauAV4VEb3AOYO9WVJHRDxXdZFmtqfW758iouoP6wZuj4gT0ustwMEREZIEbImIAwdYReu6BDzdfN3V1TVl8+bNI1+0mf3eIYccwlNPPfWb9DKAaRERY3Gq9jFgfpo+DXi47BujMLX5aDQalRRoZrs0Gg1avnedkbY4Kt1tkfQV4BXADEmPAh8DLgI+L2kS8Ft8fMOslioNj4g4dy+zXlzl55pZ9fwLUzPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPL4vAwsywODzPLUvmtF6o0Y8aM6O7uHrXP6+npYTyN2O7+tL926FNPT09s3rx5jw2NWofHaJP0TERMHes6Ror70/7auU/ebTGzLA4PG098g/RR5N2WIZCkGEd/YO5P+2vnPjk8zCyLd1sSSRdKOmus6zCri0pvN9muJG0EtgHPAb0RMTfNOlvSGcAvI+ITAy2blvs80AF8KSI+mdovBE4FtgObgMnACcDZEfHsaPSvlaQOYDXwi4g4a2+1t1vdfUnaH7gb2I/i3+3XI+JjaV4d+zMLuBF4HrAT+PeI+HyaV4/+RMSEewAbgRl92i4E3pamvzrIsh3AT4GjgH2BB4DZLet5a5pelp4vA+aMUV/fD3wZuH2g2tut7n76IaAzTU8GVgEn17g/RwAvStPTgJ+kumvTH++27G5Leh7sQNBLgP+NiEeiSP2bgde3zN+anp9Iz89S/I85qiQdCZwJfKmleaDa26Lu/kTh6fRycnoE9e3Ppoi4L01vAzYAM6lRfyZqeARwp6R7Jb0rY9mZwP+1LPNoams3nwM+SLFZ3FSX2vcgqUPSGuBxYGlErKLG/WmS1A3Modiaqk1/JuQxD2BeRDwm6TBgqaSHIuL65syIOGegZen/9wSR3rvHeiLi0xX0YUDp4O/jEXGvpFe0zupn8WiXugcSEc8BJ0o6GPiGpBOocX8AJHUCi4EFEbFVUm36MyG3PCLisfT8OPANik3FoSz7KDCrZbEjgceqqjfTPOB16YDvzcBpkm6iHrUPKCJ+DfwXcAY17o+kyRTBsSgilqTm+vRnrA8cjfYDmApMa5m+BzhjKMtSbLE9AjTYdVDr+LHu2wB9fgW7DpjWqvaWPhwKHJymDwC+D5xV4/6I4mzL5/q016Y/E3G35XCKTV4o/qK+HBHfHeqykt4D3EFxdPzaiFhfdeEjISJ6a1r7EcAN6dTzPsAtEXE71PbvYh5wPrAuHccBuCwivl2X/vgXpmaWZUIe8zCz4XN4mFkWh4eZZXF4mFkWh4eZZXF4mFkWh4eZZXF4jEOSnpO0RtKDkr6VrgVB0h9I+nqJ9z+9l/Y3SJo9yHsfkPSVvMpHRtl+2vA4PMan7RFxYkScADwFXAzFdToR8eZhrPcNFGNL9EvScRT/pk6RNGa3CxiBfloJDo/x74ekS7oldUt6ME1PkXSLpLWSvipplaTmiGpI+se0FbFS0uGSXga8DrgqbdUc3c9nvRX4D+DOtGxzXe+V9KP0WTentk5J10lal9rflNpPl/RDSfdJ+lq66hRJGyVdntrXSTo2tc9P9ayRdL+kaX36uX/L59wv6dTUfqGkJZK+K+lhSVeO8J/7+DfWF9f4MfIP4On03AF8jXThH9ANPJimPwB8MU2fAPQCc9PrAF6bpq8EPpKmrwfePMDn/gR4PnA6cFtL+2PAfmm6eXHbp2i5KAyYDsygGGpwamr7EPDRNL0RuCRN/y3F8HwA36IYNgGgk+IapNZ+LgSuS9PHAj8H9qcYmesR4KD0+mfArLH+u6vTw1se49MB6WKrzUAXsLSfZV5Ocak+EfEgsLZl3rPA7Wn6Xoov44AknQQ8ERE/A5YBL5I0Pc1eCyySdB5FSAG8CvjX5vsj4lcUwwrOBn6Q6r+AIoyampett9b0A+BqSe+lCKZedvdyiq0hIuIhipA4Js1bFhFbIuK3wI/6fJYNwuExPm2PiBMpvgz7ko559DHQDZJ+F+m/aoqBn8tcfX0ucGwaP+SnwIHAm9K8MymC4sXAvZImpc/ve1WmKEYIOzE9ZkfEO1vm7+hbUxQDT/8VxWX6K5u7MyX7uaNlumw/LXF4jGMRsQV4L/CBNPBMqxXA2QDpDMqflFjlNorBencjaR/gLcALIqI7Iropxt08N82bFRHLKYZEPJhi9+JO4D0t65gOrATmSfqj1DZF0jEMQNLREbEuIj5FMUp83/C4G3hbWvYY4A+BH5foqw3C4THORcT9FAPKnNNn1r8Bh0paS3FsYS27BoDem5uBv0sHHlsPmJ5CcWuHX7S03U2xCzITuEnSOuB+4LNRjAR2BTA9nU5+ADg1Ip6gOBbxlVTXSvYMg74WtKxjO/CdfvrZkT7/q8CFEbGj70ps6DyexwSVBtWZHBG/TUGwDDgm2uCeJlYP3sebuKYAy9PujIB3OzhsKLzlYWZZfMzDzLI4PMwsi8PDzLI4PMwsi8PDzLL8P0hjCDuCzcLHAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "exclusion_mask = Map.create(skydir=crab_pos, width=(10, 10), binsz=0.02)\n", "\n", "gammacat = SourceCatalogGammaCat(\"$GAMMAPY_DATA/gamma-cat/gammacat.fits.gz\")\n", "\n", "regions = []\n", "for source in gammacat:\n", " if not exclusion_mask.geom.contains(source.position):\n", " continue\n", " region = CircleSkyRegion(source.position, 0.15 * u.deg)\n", " regions.append(region)\n", "\n", "exclusion_mask.data = exclusion_mask.geom.region_mask(regions, inside=False)\n", "exclusion_mask.plot();" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "config = dict(\n", " outdir=\".\",\n", " background=dict(\n", " on_region=on_region,\n", " exclusion_mask=exclusion_mask,\n", " min_distance=0.1 * u.rad,\n", " ),\n", " extraction=dict(containment_correction=False),\n", " fit=dict(\n", " model=model,\n", " stat=\"wstat\",\n", " forward_folded=True,\n", " fit_range=flux_point_binning[[0, -1]],\n", " ),\n", " fp_binning=flux_point_binning,\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Run the analysis\n", "\n", "TODO: Clean up the log (partly done, get rid of remaining useless warnings)" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
FCN = 108.84841996225438TOTAL NCALL = 118NCALLS = 118
EDM = 2.5716206464251214e-06GOAL EDM = 1e-05\n", " UP = 1.0
\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
ValidValid ParamAccurate CovarPosDefMade PosDef
TrueTrueTrueTrueFalse
Hesse FailHasCovAbove EDMReach calllim
FalseTrueFalseFalse
" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
+NameValueHesse ErrorMinos Error-Minos Error+Limit-Limit+Fixed?
0par_000_amplitude3.329060.221997No
1par_001_alpha2.323550.19305No
2par_002_beta18.64439.95249No
\n", "
\n",
       "\n",
       "
" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/html": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "analysis = SpectrumAnalysisIACT(observations=observations, config=config)\n", "analysis.run(optimize_opts={\"print_level\": 1})" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Results\n", "\n", "Let's look at the results, and also compare with a previously published Crab nebula spectrum for reference." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Fit result info \n", "--------------- \n", "Model: LogParabola\n", "\n", "Parameters: \n", "\n", "\t name value error unit min max frozen\n", "\t--------- --------- --------- -------------- --- --- ------\n", "\tamplitude 3.329e-11 2.220e-12 cm-2 s-1 TeV-1 nan nan False\n", "\treference 1.000e+00 0.000e+00 TeV nan nan True\n", "\t alpha 2.324e+00 1.930e-01 nan nan False\n", "\t beta 1.864e-01 9.952e-02 nan nan False\n", "\n", "Covariance: \n", "\n", "\t name amplitude reference alpha beta \n", "\t--------- ---------- --------- ---------- ----------\n", "\tamplitude 4.928e-24 0.000e+00 2.247e-13 -6.318e-14\n", "\treference 0.000e+00 0.000e+00 0.000e+00 0.000e+00\n", "\t alpha 2.247e-13 0.000e+00 3.727e-02 -1.744e-02\n", "\t beta -6.318e-14 0.000e+00 -1.744e-02 9.905e-03 \n", "\n", "Statistic: 39.256 (wstat)\n", "Fit Range: [ 0.87992254 27.82559402] TeV\n", "\n" ] } ], "source": [ "print(analysis.fit.result[0])" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAADhCAYAAADiQzMhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3X2cZFV95/HPt6q6ZwiPI4wm8jSwQZAEIqEFjEQFow4RHYMig2gU0RF9oe5mzQaiiWbZBPOgG1EUJjCAiiCyKCiDKApMTGBDj7ICInFCmNAiOwPyzEx3Vd1f/ri3empququrq+pWdXV/37zqVfeee+45p3qK/vU5995zFBGYmZm1q9DvBpiZ2WBzIDEzs444kJiZWUccSMzMrCMOJGZm1hEHEjMz64gDiZmZdcSBxMzMOuJAYmZmHSn1uwG9sNdee8WyZcv63Qwzs4Gyfv36RyNi6Uz5FkQgWbZsGaOjo/1uhpnZQJG0sZV8HtoyM7OOOJDMZbec1+8WmJnNyIFkLrvtk/1uQe85eJoNnAVxjaRTT15yIscWNwIFFAKEKEDUv6fbk8ejsF3aZI7JdFGIyRw77BdCrAXe8Pk3UwQUoiBRjOw4UEQUs7zFbL9Aul1KW0sJKFE7nqYXJ9PTY0MKigHDytKO/TDDxRIllRgqlhguFikVSgwXhxgqFhkuDjFcLDJULFGUKBaVvhdEqSAK2XuxsC1NUms/7Ns+Cced091/QDPL1ZwPJJIOBD4K7B4Rb5kuLU8FxK8+dSQJCQEECaFgu/8UO6YrIQgSJTBFPlQlCMjyoYQzn/oF73/6/0/W/c1NNwPw+d2W8PklzwMSUIKU4zoyd72zpWwRgigAhfQ9isTkdhrqYnK7kAbXKGwXhLcF1wKFEKPA0Z9bSWEyYBa2C56lWhqiFNuC4g7BMdsfEgxFMCQYBoZqL8EwCYsEwxEsesUfsag0xHCpxOLSEIuKwywqlVhcGma4VKRUEMOlAkPFAqWiGC4WJoNly0FyOrect/CC50L8zPOYplvYStKPWzh/c0S8etrCpTXAicCmiPjNuvTlwGdI/0i+OCJmHMORdE1j0JgqbSojIyMxkHdtfWJ3+MSTUx6KCKpRpRpVJiplKlGlklQZr5SpJFXKSSVNT6qUkyrlajl7r1CuVqgkVSaSMuVKZTLPRHWCSlKl+qMvUY6ggqgAFdLtKkEZqE6mQQWoElQIqtl2lSARWVqQ1NKzIFvLEwoSgnc99TDvqQueNRfuticX7vG8LCDXgm+S2497OpPBcLvgWJzcJ2p9wbowpzS0FVSkmL0XVKSoEsXsvVQoUVSJGzas5s2HfGiyBzhUKFEqDGXbaXAbLgwxXCyxqDiU9ghLQ2nAKw2xqDjE4tIwi0vp+05DQ/zK8CJ2Ghpmp2x/qDhEqZD2MouFYs9/hjto8t22uUPS+ogYmSlfsx5JEfj9ZnUA189Q/mXA54Av1jWsCFwAvAYYA+6UdH1WX+MA+bsjYtMMdSxIkiipRIkSi4qLulv4Ubl38qZX9wvmzOxVrxZAK0ll8r1xu5JUqEQaMMcrZcrVClurZcYr2ataYbwywUS1wkSWXrnrK5QRE8AEUAkoAxOTAVOUs+1yFjRrwbQWJKskVKlS1QRJBAlJ+q6gHAnJZG80qXuvAvCvW6+DPHuZ2/0QAYooC46FrLdYoEAhCtnQamFy6DTt/WXvtZ4e2ta7I+3tDWe9vsUKhoFFBDspWKRgMcFiJVmvUJwA3HLpcVnPMS2r1FD+tvcd04q02Qs8/YbOfnY2pWaB5H0R0fQeYkkfaHY8ItZJWtaQfBSwISIeyMq4ClgREeeR9l6s5pVn97sFc85kAC10eVT2J1eQ/YadxlTH1PA+C09shCcfmty9+8H0f7Xq7vuydfd9mUA8FwWeC7E1xNYosCUKbAkxjibfxyPdnoB0HzERMBHKAmEaELd/xWRPsxYQ00CYvipKB3GT7D2UkGRDqqgKqiK2baMqUjXrLVZRoTrlR37/40/wgSeemtw/buMPAfj8HrvxhSV7zO7nF7VrhIX0GmBou0AzVYArIYa+90FKhbSnV+ulDRWmf69tN+7PdH5jvslyVOp8KHQOmvb/xoj4wUwnt5JnCnsDD9XtjwFHT5dZ0p7AXwJHSDonIs6bKm2K81YBqwD222+/Npo5ByzEMeR+Bc9+/qVa1wsrAjtnryX9a9F2kiRIIqgkwZZyhfFywnglYUu5yng5YWulyni5ynglYWs5YctEhfFKhefK4zxbnmC8XGZreYIHdy/z4coE45UyF/70VN52wCVMVCeYqJbZZ0uFickeZZlypO/V2NbLTKggJaDKtqA1GcwSVEtXFVTZLthJCYVCFT3wrxSyPPXnQoVQdXLgNU+tBJwdjtUdn+68WtoRzz+Clzz/Jbl+hkZN/6yT9DLg7cDvAr8GbAHuAW4AvhwR7QxyThWOp/1TMCIeo2GEY6q0Kc5bDayG9BrJ7JtpfbEQg+ccVyhkV4CKsHioS9dXPgFfeef2Q6gRQRKkw4ERREA1C2JJQKWaBqrnJipsLVfToFWu7rA9UamypZwwXq6ypVxlopIGvolKwkS1bruc7ZerjFeztEoaGMuVMuPVCRLqeltZoEoDVhXYPnhp8niVoVLCUCkoFROGSgmlYrpdLCYUI6FAQpEqBRJElUIkKNIgVlb6CsokPDcZ3CpRJols2DaZmAyw5WqZSlQmf44f+K0PzJ1AIulG4GHgOtK//jcBi4EXAccB10n6dETMdJ2k0Riwb93+Plk9ZgvTQhzCnOIzS6Io2r/+0USSBNUsQCVJGqyqEURCXXotz7Ze2NYsGG2ZSNhSrrBlIg1atWC0tVzNglHCeGVb0BqvVLP3hNJ//BNbo8QTz3sJW+vybM2C3NZKQjWZ3d+6AhYPFVg8VExfwwUWl2DxcLDb+P5d//nNpFmP5B0R8WhD2jPAD7PXpyTt1UaddwIHSToA+DmwEnhbG+WYzQ8LsRfW489c61V1Q0RkPaUsINUCUFK/zWTa0mv+jADGVryTau3chpsPy1lPqTZEOF6pDRvWb2//Pl5Jg1r9/ng5oZr0/q68ZtdIHgWQtDOwJSISSS8CDgFujIjyFIFmO5KuBF4F7CVpDPh4RFwi6SzgJtIh4TURcW93Po6ZWb4kUSrOIihlw4EHvWDX7ZKrSVBJksmgU022vWrXpJIkfa82vJo5YOnOs/5MnWrl1pd1wO9KWgJ8DxgFTgFOm+nEiDh1mvS1wNpZtNPMbF5JZ35or/dQC0LVZMeAMzSbINclrQQSRcRzks4APhsRfyPpR3k3zMzMptZJEMpDK5M2Krt76zTSu7VgAKZWMTOz3mglkHwYOAf4ekTcm81zdUu+zTIzs0ExY88iItaRXiep7T8AfCjPRpmZ2eDweiRmZtYRBxIzM+uIA4mZmXWkaSCR9DpJZzTO4Cvp3Xk2yszMBse0gUTSX5GuQngY8D1JH6w7fFbeDTMzs8HQrEfyBuD4iPivwJHACZL+d3Zs/k2ob2ZmbWkWSEoR6dzEEfEEaWDZTdLXSNeJMTMzaxpI/k3SK2s7EVGNiDOA+4EX594yMzMbCM0CycnAvzQmRsTH2H49ETMzW8CaTSO/pbYt6XBgWUP+a/NrlpmZDYoZp0iRtAY4HLgXqC3HEjiQmJkZrc3ie0xEHJp7S6aRTRL5UWD3iHhLllYAzgV2A0Yj4vJ+tc/MbKFr5cn22yW1FUgkrZG0SdI9DenLJd0vaYOkpgtWR8QD2UX+eiuAvYEy6RrwZmbWJ630SC4nDSaPAOOkz5BERBzewrmXAZ8DvlhLkFQELgBeQxoE7pR0Pemyu+c1nP/uiNg0RbkHA7dHxEWSriFdudHMzPqglUCyBngHcDfbrpG0JCLWNU6vAhwFbMimo0fSVcCKiDgPOLHFoseAiWy7OlUGSauAVQD77bffbJptZmaz0MrQ1n9ExPUR8e8RsbH26qDOvYGH6vbHsrQpSdpT0oXAEZLOyZKvBV4n6bPUrZVSLyJWR8RIRIwsXbq0g+aamVkzrfRIfirpK8A3SYe2AIiIdu/ammp6lZguc0Q8BpzZkPYc0HjdxMzM+qCVQLITaQB5bV1aJ7f/jrH9A437AA+3WZaZmfVZK0vtnt7lOu8EDpJ0APBzYCXwti7XYWZmPTLjNRJJl0vao25/SfaQ4owkXQncDhwsaUzSGdlEkGcBNwH3AVdHxL3tNd/MzPqtlaGtw7PZfwGIiMclHdFK4RFx6jTpa4G1rTXRzMzmslbu2ipIWlLbkfQ8WgtAZma2ALQSED4F/HP24F8AbwX+MtdWmZnZwGjlYvsXJY0Cx5PeuntSRPwk95aZmdlAmDaQSNolIp4ByALHDsGjPo+ZmU3jiU6e4Z77ml0juU7SpyS9QtLOtURJB0o6Q9JNwPL8m2hmNuCefGjmPAOs2cJWr5b0+8D7gJdnF9nLpEvt3gC8MyIe6U0zzcy64NLX977OR37cv7oBTr8h9yqaXiPxbbpmZm16YuP2PZGNP0jfd98X9ti/P23KiW/jNbOFowd/ne/g0tenQeQTT/a+7h5p5TkSMzOzaTmQmJnlbfd9Z84zwKYNJJJGJX0mWxZ3cS8bZWY2r8yzayKNmvVIjgG+DrwKuE3SWkkflvSinrTMzMwGQrPbfyvArdkLSb8GnAD8L0m/DtwRER/oQRvNzGwOa/murYj4Ben67WskFYCX5dYqMzMbGG1dbI+IJCL+qduNmYqkN0n6B0nXSXptlnagpEuyiSTNzKyPcr1rS9IaSZsk3dOQvlzS/ZI2SDq7WRkR8Y2IeC/wLuCULO2BiPCa7WZmc0DeDyReBnwO+GItQVIRuAB4Den67XdKuh4oAuc1nP/uiNiUbX8sO8/MzOaQGQOJpPOnSH4SGI2I65qdGxHrJC1rSD4K2BARD2TlXwWsiIjzgBOnqF/AJ4EbI+KHM7W37rxVwCqA/fbbr9XTzMxslloZ2loMvAT4WfY6HHgecIakv2+jzr2B+qkwx7K06XwQ+D3gLZLOBJC0p6QLgSMknTPVSRGxOiJGImJk6dKlbTTTzMxa0crQ1q8Dx2e3AyPpC8B3SIem7m6jTk2RFtNljojzgfMb0h4DzmyjbjMz67JWeiR7AzvX7e8MvDAiqsB4G3WOAfXzBewDPNxGOWZmNge00iP5G+AuSbeS9iZeAfxVttjVzW3UeSdwkKQDgJ8DK4G3tVGOmZnNAU0DSXah+zuka5IcRRpI/jQiaj2IP57h/CtJp1jZS9IY8PGIuETSWcBNpHdqrYmIezv6FGZm1jczLWwVkr4REUcCTe/Qmub8U6dJ94JZZmbzRCvXSO6Q9NLcW2JmZgOplWskxwFnSnoQeJZ0eCsi4vA8G2ZmZoOhlUByQu6tMDOzgTXj0FZEbCS9Xff4bPu5Vs4zM7OFYcaAIOnjwJ8AtSfIh4Av59koMzMbHK30LP4AeCPp9RGyW393zbNRZmY2OFoJJBMREWTTmGQPIpqZmQGtBZKrJV0E7CHpvaRPs/9Dvs0yM7NBMeNdWxHxd5JeAzwFHAz8eUR8N/eWmZnZQGhpYasscDh4mJnZDqYd2pL0rZlObiWPmZnNb816JMdmS+BOR8ChXW6PmZkNmGaBZEUL5090qyFmZjaYpg0kEXFbLxtiZmaDqaWL7f0k6U3A64HnAxdExHckFYBzgd2A0Yi4vJ9tNDNbyHKdM0vSGkmbJN3TkL5c0v2SNkg6u1kZEfGNiHgv8C7glCx5BekSwGXSpXvNzKxPWplr6/lTpB3cYvmXAcsbzi0CF5DOKnwocKqkQyUdJulbDa/6uj+WnQfp8yy3R8QfAe9vsS1mZpaDVnok/yjprbUdSf8d+HorhUfEOuCXDclHARsi4oGImACuAlZExN0RcWLDa5NSfw3cGBE/zMoYAx7PtqtT1S1plaRRSaObN29upblmZtaGVgLJq4B3SPqapHXAi0iDQbv2Bh6q2x/L0qbzQeD3gLdIOjNLuxZ4naTPAuumOikiVkfESESMLF26tIPmmplZM61MkfILSd8mnUY+Ac6JiGc6qFNTVdOk/vOB8xvSngPO6KANZmbWJTMGEknfBX4B/CawD7BG0rqI+EibdY6RLpRVsw/wcJtlmZlZn7UytHVBRPxhRDwREfcAvwM82UGddwIHSTpA0jCwEmj2BL2Zmc1hrSy1+42G/UpEnNtK4ZKuBG4HDpY0JumMiKgAZwE3AfcBV0fEvbNvupmZzQWtDG09zbZrGMOkS+0+ExG7z3RuRJw6TfpaYO0s2mlmZnNUKxfbt1tWN3vSvJO7tszMbB6Z9ZPt2VDX8Tm0xczMBlArQ1sn1e0WgBGa3K5rZmYLSyuTNr6hbrsCPEhrU8ybmdkC0Mo1ktN70RAzMxtM0waSbPqRZk+cfyiXFpmZ2UBp1iMZ7VkrzMxsYDULJFdkDw+amZlNq9ntv/9S28iGuczMzHbQLJDUz9L78rwbYmZmg6lZIPGzImZmNqNm10gOkfRj0p7Jf8m2yfYjIg7PvXVmZjbnNQskL+5ZK8zMbGBNG0giYmMvG2JmZoNp1pM29pqkF0u6UNI1kt6fpR0o6RJJ1/S7fWZmC12ugUTSGkmbJN3TkL5c0v2SNkg6u1kZEXFfRJwJvJV0wkgi4oGI8JrtZmZzQNNAIunw7P2wNsu/DFjeUGYRuAA4ATgUOFXSoZIOk/Sthtfzs3PeCPwA+F6b7TAzs5zM1CN5t6SDgLb++o+IdcAvG5KPAjZkvYoJ4CpgRUTcHREnNrw2ZeVcHxG/A5zWTjvMzCw/0wYSSR/Pjt8BFCT9eZfq3Bt4qG5/LEubrh2vknS+pIvIlueVtKekC4EjJJ0zzXmrJI1KGt28eXOXmm5mZo2a3bX1F9mQUgm4OSKu71KdmiKt2SzDtwK3NqQ9BpzZrJKIWA2sBhgZGfHDlWZmOZlpaOvoiPgA8NIu1jkG7Fu3vw/wcBfLNzOzHmoaSCLio9n7n3WxzjuBgyQdIGkYWAl0q7djZmY9lvftv1cCtwMHSxqTdEY2Nf1ZwE3AfcDVEXFvnu0wM7P8tLJme9si4tRp0teSXTg3M7PBNuefbDczs7mt2e2/RUnvk3SupJc3HPtY/k0zM7NB0KxHchHwSuAx4HxJn647dlKurTIzs4HRLJAcFRFvi4i/B44GdpF0raRFTP0siJmZLUDNAslwbSMiKhGxCrgL+D6wS94NMzOzwdAskIxK2m7CxYj4n8ClwLI8G2VmZoNj2kASEW+PiG9PkX5xRAzl2ywzMxsUze7a+h912yc3HPurPBtlZmaDo9nQ1sq67cYZdpdjZmZG80Ciaban2jczswWq2RQpMc32VPtmZjaV02/odwty1yyQ/Jakp0h7Hztl22T7i3NvmZmZDYRmC1sVe9kQMzMbTJ600czMOuJAYmZmHXEgMTOzjihi/t+AJWkzsLEheXfgyVkW1eo53cq3F/BoC+UMonZ+/oNSf7fKbrec2Z43m/yt5G0lj7/bg1H//hGxdMZcEbEgX8DqvM7pVj5gtN8/p7n08x+U+rtVdrvlzPa82eRvJW+Lefzdnkf1L+ShrW/meE63881H/f7sedbfrbLbLWe2580mfyt5+/1v22/9/vw9rz/3oa1sdcW7IuJZSW8Hfhv4TEQ0DjVZA0mjETHS73aYdZu/2/NLLwLJj4HfAg4HvgRcApwUEa/MteI6e+21VyxbtqxX1ZmZzQvr169/NFq4RtLsyfZuqURESFpB2hO5RNI7e1DvpGXLljE6OtrLKs3MBp6klkaOehFInpZ0DvB24BWSioDXMzFboE656HYAvvq+l/W5JdYtvbjYfgowDpwREY8AewN/24N6zcysB3LvkWTB49N1+/8BfDHves3MrDdyCySSnmbq6eYFRETsllfdZmbWO7kFkojYNa+yzWxwPb21zFNbK6zf+DhH7r+k382xLujZA4mSni9pv9qrV/Wa2dyxfuPj/PSRpxl7fAunXXwH6zc+3u8mWRfkHkgkvVHSz4B/B24DHgRuzLteM5t77njgMZJswLtcSbjjgcf62yDril70SM4FjgH+NSIOAF4N/FMP6jWzOeaYA/ekoHR7qFTgmAP37G+DrCt6EUjKEfEYUJBUiIhbgJf0oF4zm2OO3H8Jh/zqruyzZCeueM8xvkYyT/TigcQnJO0CrAOukLQJqPSgXjObg3ZdPMSui4ccROaRXvRIVgBbgP8GfBv4N+ANPajXBszvf2Ydx/71930B1mzA5B5IIuLZiKhGRCUiLo+I87Ohro5IWi7pfkkbJJ3djbZa//huHrPB1Yu7tp6W9FT22iqpKumpDsssAhcAJwCHAqdKOrQb7Z0rFtpf576bx2xw9XypXUlvAo6KiD/toIyXAZ+IiNdl++cARMR5U+Xfdddd48gjj2y3up7bussL+cUhJ0OhSCESXnDf1Sx+5uF+NytXC/Ezm811t9122/pW1o3pxcX27UTEN7owFLU38FDd/hhwdH0GSauAVQCLFi3qsLre2rrbvlAookKRiHR/vv9SXfzMw/zaT7+WftanHpr3n9dsPsk9kEg6qW63AIww9Rxcsyp2irTtyoyI1cBqgJGRkbj11ls7rLJ31m98nNMuvoNyJWFoaIhLzjvbd7jMU55S3eYyaapftTvqRY+k/g6tCumT7Ss6LHMM2Ldufx9g3vwJe+T+S7jiPcdwxwOPccyBezqIzGOed8rmg15MI396DsXeCRwk6QDg58BK4G051NM3R+6/xL9Y5rnanWpJwGkX3+EH9Gxg5TmN/GdpMoQVER9qt+yIqEg6C7gJKAJrIuLedssz64ep7lRzILFBlGePpLZI+stJb9H9arZ/MrC+08IjYi2wttNyzPqlNu9UEp53ygZbnuuRXA4g6V3AcRFRzvYvBL6TV71mg6I279RTWyt8ZuUR7o3YwOrFxfYXArsCv8z2d8nSzBa8tR9+Rb+bYNaxXgSSTwI/knRLtv9K4BM9qNfMzHqgF3dtXSrpRrY9MHh2RDySd71mZtYbuc21JemQ7P23SYeyHspeL8zSzMxsHsizR/JHpFOUfGqKYwEcn2PdZmbWI3netbUqez8urzrMzKz/ejGN/MmSds22PybpWklH5F2vmZn1Ri9WSPyziHha0rHA64DLgQt7UK+ZmfVALwJJNXt/PfCFiLgOGO5BvWZm1gO9CCQ/l3QR8FZgraRFParXzMx6oBe/0N9KOrni8oh4Ange8Mc9qNfMzHog90ASEc8Bm4Bjs6QK8LO86zUzs97oxV1bHwf+BDgnSxoCvtxBeSdLuldSImnGtYTNzCxfvRja+gPgjcCzABHxMOkkju26BzgJWNd508zM8rV+4+NccMsG1m98vN9NyU0vJm2ciIiQFACSdu6ksIi4LyunG20zM8vN+o2Pc9rFdzBRSRguFebtKpi96JFcnd21tYek9wI3AxfnXamkVZJGJY1u3rw57+rMzHZwxwOPMVFJSGLbKpjzUS9m//07Sa8BngIOBv48Ir7b7BxJNwO/OsWhj2bPobRS72pgNcDIyMi0S/6ameXlmAP3ZLhUoFxJ5vUqmIro7e9YSUVgZURc0WE5twIfiYjRFvJuBjY2JO8OPDnLals9p1v59gIebaGcQdTOz39Q6u9W2e2WM9vzZpO/lbyt5Fkw320N77RzYfhXdk0mnns6JrY82+v6O7R/RCydMVdE5PICdiO9U+tzwGsBAWeR/kK/rgvl3wqMdHD+6rzO6VY+YDSvf59+v9r5+Q9K/d0qu91yZnvebPK3krfFPP5uz6P687xG8iXSoay7gfeQrtN+MrAiIla0W6ikP5A0BrwMuEHSTW0W9c0cz+l2vvmo3589z/q7VXa75cz2vNnkbyVvv/9t+63fn7/n9ec2tCXp7og4LNsuknZj94uIp3OpcB6SNBoRflbG5h1/t+eXPHsk5dpGRFSBf3cQmbXV/W6AWU783Z5H8uyRVMkeQiS9PrIT8Fy2HRGxWy4Vm5lZT/X8ri0zM5tfPJ27mZl1xIHEzMw64kBiZmYdcSAZIJJ2lnS5pH+QdFq/22PWDZIOlHSJpGv63RZrjwNJn0laI2mTpHsa0pdLul/SBklnZ8knAddExHtJp+Y3m5Nm872OiAci4oz+tNS6wYGk/y4DltcnZA9wXgCcABwKnCrpUGAf4KEsW7WHbTSbrcto/XttA86BpM8iYh3wy4bko4AN2V9qE8BVwApgjDSYgP/tbA6b5ffaBpx/Gc1Ne7Ot5wFpANkbuBZ4s6Qv0P/5fMxma8rvtaQ9JV0IHCHpnKlPtbmsFysk2uxNtfxjRMSzwOm9boxZl0z3vX4MOLPXjbHucY9kbhoD9q3b3wd4uE9tMesWf6/nKQeSuelO4CBJB0gaBlYC1/e5TWad8vd6nnIg6TNJVwK3AwdLGpN0RkRUSBcBuwm4D7g6Iu7tZzvNZsPf64XFkzaamVlH3CMxM7OOOJCYmVlHHEjMzKwjDiRmZtYRBxIzM+uIA4mZmXXEgcQWPElVSXfVvc6e+az8SXpQ0t2SRiR9PWvbBklP1rX1d6Y59z2SvtSQ9oJsavchSV+V9EtJb+rNp7H5zM+R2IIn6ZmI2KXLZZayB/A6KeNBYCQiHq1LexXwkYg4cYZzlwA/A/aJiK1Z2lnAYRHxvmz/y6Tr23yjk3aauUdiNo2sR/AXkn6Y9QwOydJ3zhZuulPSjyStyNLfJelrkr4JfEdSQdLnJd0r6VuS1kp6i6RXS/p6XT2vkXRtB+18qaTbJK2XdKOkF0TE48A/A6+vy7oSuLLdesym40BiBjs1DG2dUnfs0Yj4beALwEeytI8C34+IlwLHAX8raefs2MuAd0bE8aQrWi4DDgPekx0D+D7wYklLs/3TgUvbabikRcBngDdHxJHAl4Fzs8NXkgYPJO2btWVdO/WYNeNp5M1gS0S8ZJpjtZ7CetLAAPBa4I2SaoFlMbBftv3diKgt6HQs8LWISIBHJN0C6bzp2fWLt0u6lDTA/GGbbX8x8BvAzZIAiqSz7EI6IeL5knYBTiGd2yppsx6zaTmQmDU3nr1X2fb/i0h7APfXZ5R0NPBsfVKTci8lXZxsK2mwafd6ioAfR8TvNh6IiGfzfX4mAAABFklEQVQl3Uy6CuFK4P1t1mHWlIe2zGbvJuCDyroAko6YJt8PSFe0LEh6AfCq2oGIeJh0LY6Pka5v3q6fkK4yeFTWlmFJv1F3/Ergj4E9IuLODuoxm5YDidmO10g+OUP+c4Eh4MeS7mHbNYlG/4d0mOke4CLg/wJP1h2/AngoIn7SbsMjYhx4C/BpSf8P+BFwdF2Wb5MOu13Vbh1mM/Htv2Y5krRLRDwjaU/gX4CXR8Qj2bHPAT+KiEumOfdBGm7/7XLbfPuvdYV7JGb5+paku4B/BM6tCyLrgcNJ77Kazmbge5JGut0oSV8FXk56jcasI+6RmJlZR9wjMTOzjjiQmJlZRxxIzMysIw4kZmbWEQcSMzPryH8CY3a4x6OEWDQAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "opts = {\n", " \"energy_range\": analysis.fit.fit_range,\n", " \"energy_power\": 2,\n", " \"flux_unit\": \"erg-1 cm-2 s-1\",\n", "}\n", "axes = analysis.spectrum_result.plot(**opts)\n", "CrabSpectrum().model.plot(ax=axes[0], **opts)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Exercises\n", "\n", "Rerun the analysis, changing some aspects of the analysis as you like:\n", "\n", "* only use one or two observations\n", "* a different spectral model\n", "* different config options for the spectral analysis\n", "* different energy binning for the spectral point computation\n", "\n", "Observe how the measured spectrum changes." ] }, { "cell_type": "code", "execution_count": 9, "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 }