{ "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://static.mybinder.org/badge.svg)](https://mybinder.org/v2/gh/gammapy/gammapy-webpage/v0.18.2?urlpath=lab/tree/light_curve_flare.ipynb)\n", "- You can contribute with your own notebooks in this\n", "[GitHub repository](https://github.com/gammapy/gammapy/tree/master/docs/tutorials).\n", "- **Source files:**\n", "[light_curve_flare.ipynb](../_static/notebooks/light_curve_flare.ipynb) |\n", "[light_curve_flare.py](../_static/notebooks/light_curve_flare.py)\n", "
\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Light curve - Flare\n", "\n", "## Prerequisites:\n", "\n", "- Understanding of how the light curve estimator works, please refer to the [light curve notebook](light_curve.ipynb).\n", "\n", "## Context\n", "\n", "Frequently, especially when studying flares of bright sources, it is necessary to explore the time behaviour of a source on short time scales, in particular on time scales shorter than observing runs.\n", "\n", "A typical example is given by the flare of PKS 2155-304 during the night from July 29 to 30 2006. See the [following article](https://ui.adsabs.harvard.edu/abs/2009A%26A...502..749A/abstract).\n", "\n", "**Objective: Compute the light curve of a PKS 2155-304 flare on 5 minutes time intervals, i.e. smaller than the duration of individual observations.**\n", "\n", "## Proposed approach:\n", "\n", "We have seen in the general presentation of the light curve estimator, see [light curve notebook](light_curve.ipynb), Gammapy produces datasets in a given time interval, by default that of the parent observation. To be able to produce datasets on smaller time steps, it is necessary to split the observations into the required time intervals. \n", "\n", "This is easily performed with the `select_time` method of `~gammapy.data.Observations`. If you pass it a list of time intervals it will produce a list of time filtered observations in a new `~gammapy.data.Observations` object. Data reduction can then be performed and will result in datasets defined on the required time intervals and light curve estimation can proceed directly.\n", "\n", "In summary, we have to:\n", "\n", "- Select relevant `~gammapy.data.Observations` from the `~gammapy.data.DataStore`\n", "- Apply the time selection in our predefined time intervals to obtain a new `~gammapy.data.Observations`\n", "- Perform the data reduction (in 1D or 3D)\n", "- Define the source model\n", "- Extract the light curve from the reduced dataset\n", "\n", "Here, we will use the PKS 2155-304 observations from the H.E.S.S. first public test data release. We will use time intervals of 5 minutes duration. The tutorial is implemented with the intermediate level API.\n", "\n", "## Setup \n", "\n", "As usual, we'll start with some general imports..." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "import astropy.units as u\n", "import numpy as np\n", "from astropy.coordinates import SkyCoord\n", "from astropy.time import Time\n", "from regions import CircleSkyRegion\n", "from astropy.coordinates import Angle\n", "\n", "import logging\n", "\n", "log = logging.getLogger(__name__)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now let's import gammapy specific classes and functions" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "from gammapy.data import DataStore\n", "from gammapy.datasets import SpectrumDataset\n", "from gammapy.modeling.models import PowerLawSpectralModel, SkyModel\n", "from gammapy.maps import MapAxis\n", "from gammapy.estimators import LightCurveEstimator\n", "from gammapy.makers import (\n", " SpectrumDatasetMaker,\n", " ReflectedRegionsBackgroundMaker,\n", " SafeMaskMaker,\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Select the data\n", "\n", "We first set the datastore." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "data_store = DataStore.from_dir(\"$GAMMAPY_DATA/hess-dl3-dr1/\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we select observations within 2 degrees of PKS 2155-304. " ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Number of selected observations : 21\n" ] } ], "source": [ "target_position = SkyCoord(\n", " 329.71693826 * u.deg, -30.2255890 * u.deg, frame=\"icrs\"\n", ")\n", "selection = dict(\n", " type=\"sky_circle\",\n", " frame=\"icrs\",\n", " lon=target_position.ra,\n", " lat=target_position.dec,\n", " radius=2 * u.deg,\n", ")\n", "obs_ids = data_store.obs_table.select_observations(selection)[\"OBS_ID\"]\n", "observations = data_store.get_observations(obs_ids)\n", "print(f\"Number of selected observations : {len(observations)}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Define time intervals\n", "We create the list of time intervals. Each time interval is an `astropy.time.Time` object, containing a start and stop time." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[53945.85416667 53945.86111111]\n" ] } ], "source": [ "t0 = Time(\"2006-07-29T20:30\")\n", "duration = 10 * u.min\n", "n_time_bins = 35\n", "times = t0 + np.arange(n_time_bins) * duration\n", "time_intervals = [\n", " Time([tstart, tstop]) for tstart, tstop in zip(times[:-1], times[1:])\n", "]\n", "print(time_intervals[0].mjd)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Filter the observations list in time intervals\n", "\n", "Here we apply the list of time intervals to the observations with `~gammapy.data.Observations.select_time()`.\n", "\n", "This will return a new list of Observations filtered by time_intervals. For each time interval, a new observation is created that convers the intersection of the GTIs and time interval. " ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Number of observations after time filtering: 44\n", "\n", "GTI info:\n", "- Number of GTIs: 1\n", "- Duration: 600.0 s\n", "- Start: 53945.861865555555 MET\n", "- Start: 2006-07-29T20:41:05.184 (time standard: TT)\n", "- Stop: 53945.86881 MET\n", "- Stop: 2006-07-29T20:51:05.184 (time standard: TT)\n", "\n" ] } ], "source": [ "short_observations = observations.select_time(time_intervals)\n", "# check that observations have been filtered\n", "print(\n", " f\"Number of observations after time filtering: {len(short_observations)}\\n\"\n", ")\n", "print(short_observations[1].gti)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As we can see, we have now observations of duration equal to the chosen time step.\n", "\n", "Now data reduction and light curve extraction can proceed exactly as before." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Building 1D datasets from the new observations\n", "\n", "Here we will perform the data reduction in 1D with reflected regions.\n", "\n", "*Beware, with small time intervals the background normalization with OFF regions might become problematic.*" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Defining the geometry\n", "\n", "We define the energy axes. As usual, the true energy axis has to cover a wider range to ensure a good coverage of the measured energy range chosen. \n", "\n", "We need to define the ON extraction region. Its size follows typical spectral extraction regions for HESS analyses." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "# Target definition\n", "e_reco = MapAxis.from_energy_bounds(0.4, 20, 10, \"TeV\")\n", "e_true = MapAxis.from_energy_bounds(0.1, 40, 20, \"TeV\", name=\"energy_true\")\n", "\n", "on_region_radius = Angle(\"0.11 deg\")\n", "on_region = CircleSkyRegion(center=target_position, radius=on_region_radius)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Creation of the data reduction makers\n", "\n", "We now create the dataset and background makers for the selected geometry." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "dataset_maker = SpectrumDatasetMaker(\n", " containment_correction=True, selection=[\"counts\", \"exposure\", \"edisp\"]\n", ")\n", "bkg_maker = ReflectedRegionsBackgroundMaker()\n", "safe_mask_masker = SafeMaskMaker(methods=[\"aeff-max\"], aeff_percent=10)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Creation of the datasets\n", "\n", "Now we perform the actual data reduction in the time_intervals." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "CPU times: user 18.3 s, sys: 576 ms, total: 18.9 s\n", "Wall time: 19.1 s\n" ] } ], "source": [ "%%time\n", "datasets = []\n", "\n", "dataset_empty = SpectrumDataset.create(\n", " e_reco=e_reco, e_true=e_true, region=on_region\n", ")\n", "\n", "for obs in short_observations:\n", " dataset = dataset_maker.run(dataset_empty.copy(), obs)\n", "\n", " dataset_on_off = bkg_maker.run(dataset, obs)\n", " dataset_on_off = safe_mask_masker.run(dataset_on_off, obs)\n", " datasets.append(dataset_on_off)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Define the Model\n", "\n", "The actual flux will depend on the spectral shape assumed. For simplicity, we use the power law spectral model of index 3.4 used in the [reference paper](https://ui.adsabs.harvard.edu/abs/2009A%26A...502..749A/abstract).\n", "\n", "Here we use only a spectral model in the `~gammapy.modeling.models.SkyModel` object." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "spectral_model = PowerLawSpectralModel(\n", " index=3.4, amplitude=2e-11 * u.Unit(\"1 / (cm2 s TeV)\"), reference=1 * u.TeV\n", ")\n", "spectral_model.parameters[\"index\"].frozen = False\n", "\n", "sky_model = SkyModel(\n", " spatial_model=None, spectral_model=spectral_model, name=\"pks2155\"\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Assign to model to all datasets\n", "\n", "We assign each dataset its spectral model" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "for dataset in datasets:\n", " dataset.models = sky_model" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Extract the light curve\n", "\n", "We first create the `~gammapy.time.LightCurveEstimator` for the list of datasets we just produced. We give the estimator the name of the source component to be fitted." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "lc_maker_1d = LightCurveEstimator(\n", " energy_edges=[0.7, 20] * u.TeV,\n", " source=\"pks2155\",\n", " time_intervals=time_intervals,\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can now perform the light curve extraction itself. To compare with the [reference paper](https://ui.adsabs.harvard.edu/abs/2009A%26A...502..749A/abstract), we select the 0.7-20 TeV range." ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "CPU times: user 26 s, sys: 358 ms, total: 26.4 s\n", "Wall time: 30.8 s\n" ] } ], "source": [ "%%time\n", "lc_1d = lc_maker_1d.run(datasets)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Finally we plot the result for the 1D lightcurve:" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "No handles with labels found to put in legend.\n" ] }, { "data": { "text/plain": [ "" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYkAAAERCAYAAACO6FuTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3dfbRcdX3v8ffHcDCnBT1IouAJIYg0PrHg0BTIzWoXRmsAuZJL6QW7iopUwKu92mpSsK4qdnmJ1tYroEBUBNoK2BgjRUBZBipSnvIACYiFCDWcBJsYjZhrICH53j/2HjKZzJzZc2bvM7NnPq+1zsrMfprfnpns7/yevlsRgZmZWT0v6XQBzMysezlImJlZQw4SZmbWkIOEmZk15CBhZmYNOUiYmVlDpQwSkq6RtEnSIzkd73ZJWyXdUrP8CEn3S3pC0k2S9s/j9czMyqKUQQK4Fjg5x+P9HXBOneWfAT4fEUcBvwTOy/E1zcy6XimDRET8APhF9TJJR6Y1gpWS7pb0uhaO933g1zXHEzAXWJIuug6Y317JzczKZb9OFyBHi4ELI+IJSScAXyK5yI/XwcDWiHghfT4KDLdZRjOzUumJICHpAOC/Af+SVAAAeGm67gzgU3V22xAR88Y6bJ1lzmFiZn2lJ4IESbPZ1og4tnZFRCwFlo7jmD8HhiTtl9YmpgEb2yummVm5lLJPolZEPAs8JemPIelPkHRMm8cM4E7gzHTRu4Fvt1VQM7OSURmzwEq6ATgJmAL8F/AJYDlwJXAoMADcGBH1mpnqHe9u4HXAAcAW4LyI+K6k1wA3Aq8AVgN/GhHP53s2Zmbdq5RBwszMJkZPNDeZmVkxStdxPWXKlJgxY0ani2FmViorV678eURMbXW/0gWJGTNmsGLFik4Xw8ysVCT9dDz7ubnJzMwacpAwM7OGHCTMzKyh0vVJmJn1u507dzI6Ospzzz23z7rJkyczbdo0BgYGcnktBwkzs5IZHR3lwAMPZMaMGVTlqyMi2LJlC6OjoxxxxBG5vJabm8zMSua5557j4IMP3itAAEji4IMPrlvDGC8HCTOzEqoNEM2Wj5eDhJmZNeQg0YPOuvpezrr63k4Xw8x6gIOEmVkJNUrOmnfSVgcJM7OSmTx5Mlu2bNknIFRGN02ePDm31/IQWDOzkpk2bRqjo6Ns3rx5n3WVeRJ5cZAwMyuZgYGB3OZBNOPmJjMza8hBwszMGiosSEiaLOkBSQ9LelTSJXW2OUnSryQ9lP79TVHlsb15mKyZZVFkn8TzwNyI2CZpAPihpNsi4r6a7e6OiNMKLIeZmY1TYUEikrFZ29KnA+lfvgN4zcysUIX2SUiaJOkhYBNwR0TcX2ez2WmT1G2S3tjgOOdLWiFpRb0hX7bHstUbWL1+K/c/9QvmLFrOstUbOl0kMyuxQoNEROyKiGOBacDxkt5Us8kq4PCIOAa4HFjW4DiLI2JWRMyaOrXl+3j3jWWrN3Dx0rXs2LUbgA1bt3Px0rUOFGY2bhMyTyIitkq6CzgZeKRq+bNVj2+V9CVJUyLi5xNRrl5R6YBevX7riwGiYvvOXSxcsob5I8OdKJqZlVyRo5umShpKHw8CbwV+XLPNIUrz2ko6Pi3PlqLK1OtqA0Sz5WZmzRRZkzgUuE7SJJKL/zci4hZJFwJExFXAmcD7Jb0AbAfOjryzU/WBmy6YDcCcRcvZsHX7PuuHhwYnukhm1iOKHN20Bhips/yqqsdXAFcUVYZ+s2DeTC5eupbtO3e9uGxwYBIL5s3sYKnMrMw847qHzB8Z5tIzjmb/ScnHOjw0yKVnHL1Pf4RHQJlZVk7w12PmjwxzwwPrgT3NUNUajYCq7GtmVs1Bok94BJSZjYebm/qMR0CZWStck+gTHgFlZuPhmkSfWTBvJoMDk/Za5hFQZtaIaxJ9ptLvsHDJGnbs2s3w0CAL5s10f4SZ1eUg0YeajYAyM6twc5OZmTXkmkQPcu3AzPLimoSZmTXkIGFmZg05SJiZWUMOEmZm1pCDhJmZNeTRTX3KI6DMLAvXJMzMrCEHiS5y1tX3vpjS28ysGzhIWEMOWmbmIGFmZg05SJiZWUMOEmZm1lBhQULSZEkPSHpY0qOSLqmzjSRdJmmdpDWSjiuqPGZm1roi50k8D8yNiG2SBoAfSrotIu6r2uYU4Kj07wTgyvRfMzPrAoXVJCKxLX06kP5FzWanA9en294HDEk6tKgymZlZawrtk5A0SdJDwCbgjoi4v2aTYeDpquej6bLa45wvaYWkFZs3by6uwGZmtpdCg0RE7IqIY4FpwPGS3lSziertVuc4iyNiVkTMmjp1ahFF7bhlqzewev1W7n/qF8xZtJxlqzd0ukhmZhMzuikitgJ3ASfXrBoFDqt6Pg3YOBFl6ibLVm/g4qVr2bFrNwAbtm7n4qVruz5QeLKdWe8rcnTTVElD6eNB4K3Aj2s2uxl4VzrK6UTgVxHxTFFl6jaVi+zCJWvYvnPXXuu279zFwiVrOlQy12zMLFHk6KZDgeskTSIJRt+IiFskXQgQEVcBtwKnAuuA3wDnFlierlWpQWRdXrRGNRuA+SP7dBmZWQ8rLEhExBpgpM7yq6oeB/CBosrQ7SrpuucsWs6Grdv3WT88NDih5ak0Ha1ev3WfAFWp2ThImPUXz7juAgvmzWRwYNJeywYHJrFg3syOlKfbajZm1jm+6VAXqPw6X7hkDTt27WZ4aJAF82ZO+K/2bqvZmFnnuSbRJeaPDDMyfYgTjngF91w0t6PNOt1WszGzznFNwvbRLTUbM+s8Bwmra/7IMDc8sB6ofz/syhDZHbt2M2fRcgcRsx7l5iZrWVkn/5lZ6xwkrCXNJv95BrZZb3GQsJZ5iKxZ/3CQsJbcdMHshkNhh4cG6/ZfmFl5OUh0kZsumF2Ki6yHyJr1D49uspZ5iKxZ/3CQsHFpNkTWzHqDm5vMzKwh1ySsIdcQzMw1CTMza8hBwszMGnKQMDOzhjL1SUh6JTAHeDWwHXgEWBERnmJrZtbDxgwSkt4MXAS8AlgNbAImA/OBIyUtAf4+Ip4tuqDWfdyxbdb7mtUkTgXeFxHra1dI2g84DfhD4JsFlM3MzDpszCAREQvGWPcCsCz3EpmZWdcYd8e1pHPzLIiZmXWfdkY3XTLWSkmHSbpT0mOSHpX0oTrbnCTpV5IeSv/+po3ymJlZzpp1XK9ptAp4VZNjvwB8JCJWSToQWCnpjoj4Uc12d0fEadmKa2ZmE6lZx/WrgHnAL2uWC/j3sXaMiGeAZ9LHv5b0GDAM1AYJa0Hlzm8eWWRmE6FZkLgFOCAiHqpdIemurC8iaQYwAtxfZ/VsSQ8DG4GPRsSjdfY/HzgfYPr06Vlf1szM2tRsdNN5Y6z7kywvIOkAkiGyH64zn2IVcHhEbJN0KsloqaPqvNZiYDHArFmzIsvrmplZ+1ruuE5/1WfddoAkQPxzRCytXR8Rz0bEtvTxrcCApCmtlsnMzIoxntFNF2bZSJKArwKPRcQ/NNjmkHQ7JB2flmfLOMpkZmYFGM/9JJRxuznAOcBaSZU+jY8B0wEi4irgTOD9kl4gyQl1dkS4OcnMrEuMJ0j89ywbRcQPaRJQIuIK4IpxlMHMzCZApuYmSR+S9LK0aegSSaskva3gspmZWYdl7ZN4bzoy6W3AVOBcYFFhpbK6lq3ewOr1W7n/qV8wZ9Fylq3e0OkimVmPyxokKs1GpwJfi4iHyd43YTlYtnoDFy9dy45dyS08NmzdzsVL1zpQmFmhsvZJrJT0PeAI4OI0zYZvODQBKjOsV6/f+mKAqNi+cxcLl6xh/shwJ4pmZn0ga5A4DzgWeDIifiPpYJImJ5sgtQGi2XIzszxkChLpbUpXVT3fguczvChLPqXx5lyqbD9n0XI2bN2+z/rhocGWjmdm1op2UoXbBFowbyaDA5P2WjY4MIkF82Z2qERm1g/GM0/COqDS77BwyRp27NrN8NAgC+bNdH+EmRXKQaJE5o8Mc8MDye3GnSrczCbCmM1Nko6WdJ+kpyUtlnRQ1boHii+emZl1UrM+iSuBTwJHA48DP5R0ZLpuoMBymZlZF2jW3HRARNyePv6cpJXA7ZLOAZyIz8ysxzWrSUjSyytPIuJO4I+AfwQOL7JgvcTpNMysrJoFic8Ar69eEBFrgLcA+9xEyPbldBpmVmbNbl/69dplkg6JiPXA+worVYlUagk7du1mzqLlLw5LdTqNxHgnEZpZdxjPENhbgePyLkgZNaolVMs7nYYvtmY2kYq8M13PylJLePzTpwBOp2Fm5TaetBxfzr0UJZWlluB0GmZWZplrEulEusOA+yQdBxARq8beqze1knTP6TTMrMwyBQlJfwu8B/gJe+ZHBDC3mGKVw4J5M7l46Vq279z14rJ6tQSn0zCzsspak/ifwJERsaPIwpSNawlm1uuy9kk8Agy1cmBJh0m6U9Jjkh6V9KE620jSZZLWSVpTacYqk/kjw4xMH+KEI17BPRfNdYCo4kmEZuWXtSZxKbBa0iPA85WFEfGOMfZ5AfhIRKxKb3e6UtIdEfGjqm1OAY5K/04gyRV1QisnYN1prOHBDqRm5ZE1SFxHMvt6LRnvbR0RzwDPpI9/LekxYBioDhKnA9dHRJB0iA9JOjTd10rqrKvvHXN48A0PrHffjFlJZA0SP4+Iy8b7IpJmACPA/TWrhoGnq56PpsscJEquzPfk9ixxsz2yBomVki4Fbmbv5qamQ2AlHQB8E/hwRDxbu7rOLvtkl5V0PnA+wPTp0zMWubv00wXnpgtmjzk8uJ/eC7OyyxokRtJ/T6xa1nQIrKQBkgDxzxFRLyHgKMnci4ppwMbajSJiMbAYYNasWU5RXgJZhwfnyTUAs/xlChIR8eZWDyxJwFeBxyLiHxpsdjPwQUk3knRY/8r9Eb3Bw4PNekPWyXT/B/hsRGxNnx9EMnLp42PsNgc4B1gr6aF02ceA6QARcRVJssBTgXXAb4Bzx3MS1p08idCs/LI2N50SER+rPImIX0o6FWgYJCLihzRJBpiOavpAxjJ0LV8AzaxXZZ1MN0nSSytPJA0CLx1jezMz6wFZg8Q/Ad+XdJ6k9wJ3kMydMJsQZ11974sd00XyLHGzvWXtuP6spDXAW0makP42Ir5baMnMWtDoDoGtHsOzxM32NmaQkKS034CIuB24faxtzGpNRH9NHhd3zxI3q69Zc9Odkv5c0l4z2CTtL2mupOuAdxdXPLOxnXX1vSxcsmav+Riw5+LeShNVmWeJmxWlWZA4GdgF3CBpo6QfSXoKeAJ4J/D5iLi24DJaj2u3vyGPi/tNF8xueEtZzxK3fjZmkIiI5yLiSxExBzgceAswEhGHR8T7IuKhsfY3K1qeF3ffatZsX5nvcR0ROyPimcqEOrOJ0mzEUV4X9/kjw1x6xtHsPyn5bzE8NMilZxztTmvra5nvcW3WCVk6pfNMAeJZ4mZ7c5CwrtXKiCNf3M2Kkam5SdIb6iw7KffSmNXwiCOzzsraJ/ENSX+V3pN6UNLlJLc0NSuMRxyZdV7WIHECyX0f/h14kOSeD3OKKpT1j4nqlDaz8cnaJ7ET2A4MApOBpyLC9X1ry0R3SmflGorZHsqSUUPSw8C3gb8FDgauBnZGxJnFFm9fs2bNihUrVkz0y1rOGnVKA+w/6SWMTB/a62Ltu861z+9hf5O0MiJmtbpf1prEeRFRuTL/DDhd0jmtvphZNXdKm3W/rH0SmyRNr/4D/q3Igllvc6e0WTlkDRLfAW5J//0+8CRwW1GFsv7QiU7pibovhVmvyHo/iaOrn0s6DrigkBJZ3+hEp7SZtWZcM64jYpWk38u7MN3InX3FyjpT2u+/WWdkChKS/rLq6UuA44DNhZTIzMy6RtaaxIFVj18g6Zv4Zv7FMTOzbpK1T+KSVg8s6RrgNGBTRLypzvqTSOZePJUuWhoRn2r1dczMrDjN7nH9r0DD2XYR8Y4xdr8WuAK4foxt7o6I08Yqg1leKilAduzazZxFy/uqk7yfz93a06wm8bnxHjgifiBpxnj3N8tTlhQgvaqfz93a1yxIPBUR6wt8/dlpyo+NwEcj4tF6G0k6HzgfYPr06QUWxzqh6JFLrdyXopdURuaNde4OEtZMs8l0yyoPJOXdUb0KODwijgEur36tWhGxOCJmRcSsqVOn5lwM6wf9nAKkn8/d2tcsSKjq8WvyfOGIeDYitqWPbwUGJE3J8zXMoH9TgNx0weym527WTLMgEQ0et03SIZKUPj4+LcuWPF+jXc3udWDl0c/3pejnc7f2NeuTOEbSsyQ1isH0MenziIiXNdpR0g3AScAUSaPAJ4ABkh2vAs4E3i/pBZJ7VZwdWfKWTxB39vWWfk4B0s/nbu3LdD+JblL0/STG6uyD5F4Hj3/6lMJe34rVz2lW+vncbfz3k8iaBbbvuLPPzGycCf56Rb1fVpXHcxYtZ8PW7fvs484+M+snrkk04M4+M7M+r0mMxZ19Nha371u/cJAYQ9Z7HVh5+HM0a42bm8zMrCHXJMz6hGtRNh6uSZi1yDPxrZ84SJi1oNFMfAcK61V9GyT8a9BaddbV97JwyRq279y11/JK2u3KiCezXtKXQcK/Bm28PBPf+k1f5W5yXiZr11gz8e+5aG4HSmSWjXM3tcC/Bm28PBPf+k1fDYF1XiZrl2fiW7/py5qEfw1aO+aPDDMyfYgTjngF91w01wGiw866+l4PGihQX9UkKvxr0Mwsm74MEuC8TGZmWfRlc5OZWav6tVmrb2sSZu1w7bO3OPV7Y65JmFlpOXNC8RwkzKyUnDlhYri5ycxKZazMCZU8Wh6pmJ/CahKSrpG0SdIjDdZL0mWS1klaI+m4ospiZr0nr8wJWZqs+rlZq8iaxLXAFcD1DdafAhyV/p0AXJn+O2HcSWVWPnlmTmjUZAV75lNl2aaXFRYkIuIHkmaMscnpwPWRZBi8T9KQpEMj4pmiymRmvWPBvJlcvHTtXqnbW8mccNbV947ZZFWZR9Vsm17/sdnJjuth4Omq56Ppsn1IOl/SCkkrNm/ePCGFM7PuNn9kmEvPOJr9JyWXseGhQS494+iWft1nabLq94SgnQwSqrOsbt7yiFgcEbMiYtbUqVMLLpZZ+eQ10atsE8bayaN10wWzGzZNDQ8NctMFszNt0+s6GSRGgcOqnk8DNnaoLGbWwxoFvyzJPvs9IWgng8TNwLvSUU4nAr9yf4RZbyhLjSRLk1UezVplVljHtaQbgJOAKZJGgU8AAwARcRVwK3AqsA74DXBuUWUxs97VbpNPlmSf/ZwQtMjRTe9ssj6ADxT1+mbWmspcgB27djNn0XKnzzfAaTnMClOWJhfo7RQX/TwRLg9Oy2FWcu3UALo9xUW72VmzToTLcvx+a2aqcJAwK7G8ZgPnORegG5qtuj34NdNNqcsdJMxKKK+LYJ4pLqD7Ulj0+0S4PLhPwqzE8roItjsXoNL/snDJmr3SZMCeoNWqdvoSskyEs2wcJMwKUHRnad4XwbzmAuSZmTWPjvR+nwiXBzc3meUs7yaXsdqn201yV62duQB5NVvl3ZdQ2XbhkjXs2LWb4aHBrh/a2w19OtVckzDLUbMml7yHxHbbbOC8frnn2ZfQTn6nidaNQ5FdkzDL2UR3luY5GzjL/mPVbNr95Z53R/pEa2dUUpbU5Z0Y7eSahFmO8s4aWsaJYHn8cs+7L6HSh9PtunE0loOEWc7yusB1Y9PDROm2ZrSJ0K1pyR0kzHKWxwUu7+GkeZnImk2Z+hIgn/emG0djuU/CrAB59BN0W9NDK6O2ytC0k6e8RrR142gsJclYy2PWrFmxYsWKThfDrKksnZhjbTNW5+09F83NqZTNjTUsFWD/SS/h8U+fMmHl6TaNOpwheW9Gpg+1HDSLSMshaWVEzGp1Pzc3mRWk3c7Sbmt66LaaTTfp5ffGQcKsQ5q1YXdL561TXIytWzuc8+IgYdYBWUcudVPnbbfVbLpJLw/ZdZ+E2QQrog17oixbvaGrOlW7Sdb3plNpwMfbJ+HRTWYdUNY27Dxnd3fTPRPy0Kv3wXaQMJtgN10we8yRS91+gen28nVSs/em25L3ZeE+CbMO6Pf2/TKmG2lXWWfQuyZh1gHdOGlqonTb3esmQrcm78ui0CAh6WTgC8Ak4CsRsahm/UnAt4Gn0kVLI+JTRZbJrFv0aht2I2W/73S7ytoPVViQkDQJ+CLwh8Ao8KCkmyPiRzWb3h0RpxVVDrNu1g/BoVZZL5btKHM/VJF9EscD6yLiyYjYAdwInF7g65lZF+v3SXll7YcqMkgMA09XPR9Nl9WaLelhSbdJemO9A0k6X9IKSSs2b95cRFnNbIKU9WLZrm6ZQd+qIvskVGdZ7cy9VcDhEbFN0qnAMuCofXaKWAwshmQyXd4FNbOJ08+d9mXshyoySIwCh1U9nwZsrN4gIp6tenyrpC9JmhIRPy+wXGbWYWW8WOalbOdbZJB4EDhK0hHABuBs4E+qN5B0CPBfERGSjidp/tpSYJnMrEuU7WLZrwoLEhHxgqQPAt8lGQJ7TUQ8KunCdP1VwJnA+yW9AGwHzo6yJZMyM+thTvBnZtYHfNMhMzPLnYOEmZk15CBhZmYNOUiYmVlDDhJmZtaQg4SZmTXkIGFmZg2Vbp6EpM3AT5tsNgXo99Qe/f4e+Pz7+/zB70Ht+R8eEVNbPUjpgkQWklaMZ9JIL+n398Dn39/nD34P8jp/NzeZmVlDDhJmZtZQrwaJxZ0uQBfo9/fA52/9/h7kcv492SdhZmb56NWahJmZ5cBBwszMGipdkJD0n5LWSnpI0j43lpD0ckn/KulhSY9KOjfrvmWQ4fwPkvQtSWskPSDpTVXrTpb0H5LWSbpoYkuenzbfg174DgxJWiLpx5IekzS7Zr0kXZZ+zmskHVe1rvTfgTbPv/SfP2R6D14n6V5Jz0v6aM261r4DEVGqP+A/gSljrP8Y8Jn08VTgF8D+WfYtw1+G8/874BPp49cB308fTwJ+ArwG2B94GHhDp89nIt+DHvoOXAf8Wfp4f2CoZv2pwG2AgBOB+3vpOzDe8++Vzz/je/BK4PeATwMfrVre8negdDWJDAI4UJKAA0iCxAudLdKEegPwfYCI+DEwQ9KrgOOBdRHxZETsAG4ETu9cMQvV6D0oPUkvA/4A+CpAROyIiK01m50OXB+J+4AhSYfSA9+BNs+/J2R5DyJiU0Q8COys2b3l70AZg0QA35O0UtL5ddZfAbwe2AisBT4UEbsz7lsGzc7hYeAMAEnHA4cD04Bh4Omq7UbTZWU03vcgy77d7jXAZuBrklZL+oqk367ZptFn3QvfgXbOH8r/+UO296CRlr8DZQwScyLiOOAU4AOS/qBm/TzgIeDVwLHAFWnkzbJvGTQ7h0XAQZIeAv4cWE1Sk1KdY5V1/PN434Ms+3a7/YDjgCsjYgT4f0Btu3Kjz7oXvgPtnD+U//OHbO9BIy1/B0oXJCJiY/rvJuBbJNWnaucCS9Oq5jrgKZJ26Sz7dr1m5xARz0bEuRFxLPAukn6Zp0h+MRxWtek0ktpW6bTxHvTCd2AUGI2I+9PnS0guGLXb1Puse+E70M7598LnD9neg7H2bek7UKogIem3JR1YeQy8DXikZrP1wFvSbV4FzASezLhvV8tyDumoh/3Tp38G/CAingUeBI6SdES6/mzg5okrfT7aeQ964TsQET8DnpY0M130FuBHNZvdDLwrHeVzIvCriHiGHvgOtHP+vfD5Q+b3oJHWvwOd7qVvsUf/NSTtzQ8DjwJ/nS6/ELgwffxq4Hsk/RGPAH861r5l+st4/rOBJ4AfA0uBg6r2PxV4nGR0Q+nOv933oBe+A+l5HAusANYAy4CDas5fwBfTz3ktMKvHvgPjOv9e+fwzvgeHkNQangW2po9fNp7vgNNymJlZQ6VqbjIzs4nlIGFmZg05SJiZWUMOEmZm1pCDhJlZGyR9UtKGNGngQ5JOrbPNZCXJJiuJRy+pWndMmoxvrZLkpC+r2Xe6pG21ifqalOlySdvaO7OEg4T1NEkHV/3n/VnVf+Ztkr5U0Gt+WNK70sfXSvpNZXx+uuwLkkLSlPT5tvTfGZK2p6kWHksvKu+u2u+06ouLTTxJJ0m6ts6qz0fEsenfrXXWPw/MjYhjSIavnpzO4QD4CnBRRBxNMsFvQe2xSRIWZi3jLGAo6/bNOEhYT4uILZX/vMBV7PnPfEBE/K+8X0/SfsB7ga9XLV5HmkRN0kuANwMbGhziJxExEhGvJ5no9Bfak+7+O8A7JP1W3uW2YkWi8st+IP2rzD+YCfwgfXwH8EeV/STNB54kmddB1fK3pbWPVZL+RdIB6fJJJFmQF+ZVdgcJ60vpL8Jb0seflHSdpO8pud/AGZI+m1b/b5c0kG73u5L+TUlyuO+qfmbRucCqiKjOPHwDcFb6+CTgHjJkJo6IJ4G/BP53+jyAu4DTxnXSVqQPKrl3xTWSDqq3gaRJSvKJbQLuiD1pNR4B3pE+/mPStBnprPC/Ai6pOc4U4OPAWyPJQ7WC5HsC8EHg5khm2OfCQcIscSTwdpJf/P8E3JlW/7cDb08DxeXAmRHxu8A1JLn6a80BVtYsewKYml483kmSnjmrVaS5x1IrgN9vYX/LgaT70wv8V0hqc5UmzHnAlSTfn2OBZ4C/r3eMiNiV1minAcdrz82w3kuSbHAlcCCwI11+CUnNt7Zv4USSdPj3pGV6N3C4pFeTBJnL8znrxH55HsysxG6LiJ2S1pLcmOX2dPlaYAZJk8CbgDskkW5T79faocBjdZYvJWk+OgG4oIVy1Wbt3ESSesYmUEScAEkNFHhPRLyn3naSvgzc0uRYWyXdBZwMPBLJPU/elu7/OyQ/ViD5rpwp6bMkfQy7JT0H/JSkJvLOmtd+O/BaYF36Hf0tSesi4rUtn3AVBwmzxPMAEbFb0s7Yk69mN8n/EwGPRsTsRgdIbQcm11l+I0mt4Lr0NbKWa4S9g87k9DWsS0g6tKp5539QJ2mgpKnAzjRADAJvBT6TrntlRGxK+6s+TtJ3RkT8ftX+nwS2RcQV6bG+KOm1EbEu7aOaFhHfIcnZVNlnW7sBAtzcZJbVf5A0Gb2R6SIAAAEySURBVM0GkDQg6Y11tnuM5NfcXiJiPfDXQOYRVZJmAJ9j7+aD36GEmUt7XKX/ag3JoIS/AJD0akmVkU6HAnem2zxIUhOo1DjeKelxkoSUG4GvjfViEbEZeA9wQ3q8+9i7STJXrkmYZRAROySdCVwm6eUk/3f+LzWjTkiGKv5jg2NcXbssHQ31fNWiIyWtJqkx/Bq4PCKqLxpvBi4e94lYWyLiLpLBA9XLzmmw7UaSjKtExBqSWmG97b4AfKHJ636y5vlykntYj7XPAWOtz8pZYM1yJulbwMKIeCLDtscAX46Ipje/UXJ/lK9HxFtyKKZZJm5uMsvfRSTNC2OSdCHJ8NiPZzzudOAjbZTLrGWuSZiZWUOuSZiZWUMOEmZm1pCDhJmZNeQgYWZmDTlImJlZQ/8fhz85TKOeZs8AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "lc_1d.plot(marker=\"o\")" ] }, { "cell_type": "code", "execution_count": null, "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.7.0" }, "nbsphinx": { "orphan": true } }, "nbformat": 4, "nbformat_minor": 4 }