.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "tutorials/analysis-3d/event_sampling.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. or to run this example in your browser via Binder .. rst-class:: sphx-glr-example-title .. _sphx_glr_tutorials_analysis-3d_event_sampling.py: Event sampling ============== Learn to sampling events from a given sky model and IRFs. Prerequisites ------------- To understand how to generate a model and a `~gammapy.datasets.MapDataset` and how to fit the data, please refer to the `~gammapy.modeling.models.SkyModel` and :doc:`/tutorials/analysis-3d/simulate_3d` tutorial. Context ------- This tutorial describes how to sample events from an observation of a one (or more) gamma-ray source(s). The main aim of the tutorial will be to set the minimal configuration needed to deal with the Gammapy event-sampler and how to obtain an output photon event list. The core of the event sampling lies into the Gammapy `~gammapy.datasets.MapDatasetEventSampler` class, which is based on the inverse cumulative distribution function `(Inverse CDF) `__. The `~gammapy.datasets.MapDatasetEventSampler` takes in input a `~gammapy.datasets.Dataset` object containing the spectral, spatial and temporal properties of the source(s) of interest. The `~gammapy.datasets.MapDatasetEventSampler` class evaluates the map of predicted counts (`npred`) per bin of the given Sky model, and the `npred` map is then used to sample the events. In particular, the output of the event-sampler will be a set of events having information about their true coordinates, true energies and times of arrival. To these events, IRF corrections (i.e. PSF and energy dispersion) can also further be applied in order to obtain reconstructed coordinates and energies of the sampled events. At the end of this process, you will obtain an event-list in FITS format. Objective --------- Describe the process of sampling events from a given Sky model and obtain an output event-list. Proposed approach ----------------- In this section, we will show how to define an observation and create a Dataset object. These are both necessary for the event sampling. Then, we will define the Sky model from which we sample events. In this tutorial, we propose examples for sampling events of: - `a point-like source <#sampling-the-source-and-background-events>`__ - `a time variable point-like source <#time-variable-source-using-a-lightcurve>`__ - `an extended source using a template map <#extended-source-using-a-template>`__ - `a set of observations <#simulate-multiple-event-lists>`__ We will work with the following functions and classes: - `~gammapy.data.Observations` - `~gammapy.datasets.Dataset` - `~gammapy.modeling.models.SkyModel` - `~gammapy.datasets.MapDatasetEventSampler` - `~gammapy.data.EventList` .. GENERATED FROM PYTHON SOURCE LINES 79-84 Setup ----- As usual, let’s start with some general imports… .. GENERATED FROM PYTHON SOURCE LINES 84-115 .. code-block:: Python from pathlib import Path import numpy as np import astropy.units as u from astropy.coordinates import Angle, SkyCoord from astropy.io import fits from astropy.time import Time from regions import CircleSkyRegion import matplotlib.pyplot as plt from IPython.display import display from gammapy.data import ( DataStore, FixedPointingInfo, Observation, observatory_locations, ) from gammapy.datasets import MapDataset, MapDatasetEventSampler from gammapy.irf import load_irf_dict_from_file from gammapy.makers import MapDatasetMaker from gammapy.maps import MapAxis, WcsGeom from gammapy.modeling.models import ( ExpDecayTemporalModel, FoVBackgroundModel, Models, PointSpatialModel, PowerLawNormSpectralModel, PowerLawSpectralModel, SkyModel, TemplateSpatialModel, ) .. GENERATED FROM PYTHON SOURCE LINES 116-118 Check setup ----------- .. GENERATED FROM PYTHON SOURCE LINES 118-122 .. code-block:: Python from gammapy.utils.check import check_tutorials_setup check_tutorials_setup() .. rst-class:: sphx-glr-script-out .. code-block:: none System: python_executable : /home/runner/work/gammapy-docs/gammapy-docs/gammapy/.tox/build_docs/bin/python python_version : 3.9.19 machine : x86_64 system : Linux Gammapy package: version : 1.3.dev283+g0ce322691 path : /home/runner/work/gammapy-docs/gammapy-docs/gammapy/.tox/build_docs/lib/python3.9/site-packages/gammapy Other packages: numpy : 1.26.4 scipy : 1.13.0 astropy : 5.2.2 regions : 0.8 click : 8.1.7 yaml : 6.0.1 IPython : 8.18.1 jupyterlab : not installed matplotlib : 3.8.4 pandas : not installed healpy : 1.16.6 iminuit : 2.25.2 sherpa : 4.16.0 naima : 0.10.0 emcee : 3.1.6 corner : 2.2.2 ray : 2.20.0 Gammapy environment variables: GAMMAPY_DATA : /home/runner/work/gammapy-docs/gammapy-docs/gammapy-datasets/dev .. GENERATED FROM PYTHON SOURCE LINES 123-135 Define an Observation --------------------- You can firstly create a `~gammapy.data.Observations` object that contains the pointing position, the GTIs and the IRF you want to consider. Hereafter, we chose the IRF of the South configuration used for the CTA DC1 and we set the pointing position of the simulated field at the Galactic Center. We also fix the exposure time to 1 hr. Let’s start with some initial settings: .. GENERATED FROM PYTHON SOURCE LINES 135-149 .. code-block:: Python path = Path("$GAMMAPY_DATA/cta-caldb") irf_filename = "Prod5-South-20deg-AverageAz-14MSTs37SSTs.180000s-v0.1.fits.gz" # telescope is pointing at a fixed position in ICRS for the observation pointing = FixedPointingInfo( fixed_icrs=SkyCoord(0.0, 0.0, frame="galactic", unit="deg").icrs, ) livetime = 1 * u.hr location = observatory_locations["cta_south"] irfs = load_irf_dict_from_file(path / irf_filename) .. GENERATED FROM PYTHON SOURCE LINES 150-152 Now you can create the observation: .. GENERATED FROM PYTHON SOURCE LINES 152-162 .. code-block:: Python observation = Observation.create( obs_id=1001, pointing=pointing, livetime=livetime, irfs=irfs, location=location, ) print(observation) .. rst-class:: sphx-glr-script-out .. code-block:: none Observation obs id : 1001 tstart : 51544.00 tstop : 51544.04 duration : 3600.00 s pointing (icrs) : 266.4 deg, -28.9 deg deadtime fraction : 0.0% .. GENERATED FROM PYTHON SOURCE LINES 163-185 Define the MapDataset --------------------- Let’s generate the `~gammapy.datasets.Dataset` object (for more info on `~gammapy.datasets.Dataset` objects, please checkout :doc:`/tutorials/api/datasets` tutorial): we define the energy axes (true and reconstructed), the migration axis and the geometry of the observation. *This is a crucial point for the correct configuration of the event sampler. Indeed, the spatial and energetic binning should be treated carefully and… the finer the better. For this reason, we suggest to define the energy axes (true and reconstructed) by setting a minimum binning of least 10-20 bins per decade for all the sources of interest. The spatial binning may instead be different from source to source and, at first order, it should be adopted a binning significantly smaller than the expected source size.* For the examples that will be shown hereafter, we set the geometry of the dataset to a field of view of 2degx2deg and we bin the spatial map with pixels of 0.02 deg. .. GENERATED FROM PYTHON SOURCE LINES 185-200 .. code-block:: Python energy_axis = MapAxis.from_energy_bounds("0.1 TeV", "100 TeV", nbin=10, per_decade=True) energy_axis_true = MapAxis.from_energy_bounds( "0.03 TeV", "300 TeV", nbin=20, per_decade=True, name="energy_true" ) migra_axis = MapAxis.from_bounds(0.5, 2, nbin=150, node_type="edges", name="migra") geom = WcsGeom.create( skydir=pointing.fixed_icrs, width=(2, 2), binsz=0.02, frame="galactic", axes=[energy_axis], ) .. GENERATED FROM PYTHON SOURCE LINES 201-206 In the following, the dataset is created by selecting the effective area, background model, the PSF and the Edisp from the IRF. The dataset thus produced can be saved into a FITS file just using the `write()` function. We put it into the `evt_sampling` sub-folder: .. GENERATED FROM PYTHON SOURCE LINES 208-220 .. code-block:: Python empty = MapDataset.create( geom, energy_axis_true=energy_axis_true, migra_axis=migra_axis, name="my-dataset", ) maker = MapDatasetMaker(selection=["exposure", "background", "psf", "edisp"]) dataset = maker.run(empty, observation) Path("event_sampling").mkdir(exist_ok=True) dataset.write("./event_sampling/dataset.fits", overwrite=True) .. GENERATED FROM PYTHON SOURCE LINES 221-228 Define the Sky model: a point-like source ----------------------------------------- Now let’s define a sky model for a point-like source centered 0.5 deg far from the Galactic Center and with a power-law spectrum. We then save the model into a yaml file. .. GENERATED FROM PYTHON SOURCE LINES 228-250 .. code-block:: Python spectral_model_pwl = PowerLawSpectralModel( index=2, amplitude="1e-12 TeV-1 cm-2 s-1", reference="1 TeV" ) spatial_model_point = PointSpatialModel( lon_0="0 deg", lat_0="0.5 deg", frame="galactic" ) sky_model_pntpwl = SkyModel( spectral_model=spectral_model_pwl, spatial_model=spatial_model_point, name="point-pwl", ) bkg_model = FoVBackgroundModel(dataset_name="my-dataset") models = Models([sky_model_pntpwl, bkg_model]) file_model = "./event_sampling/point-pwl.yaml" models.write(file_model, overwrite=True) .. GENERATED FROM PYTHON SOURCE LINES 251-257 Sampling the source and background events ----------------------------------------- Now, we can finally add the `~gammapy.modeling.models.SkyModel` we want to event-sample to the `~gammapy.datasets.Dataset` container: .. GENERATED FROM PYTHON SOURCE LINES 257-261 .. code-block:: Python dataset.models = models print(dataset.models) .. rst-class:: sphx-glr-script-out .. code-block:: none DatasetModels Component 0: SkyModel Name : point-pwl Datasets names : None Spectral model type : PowerLawSpectralModel Spatial model type : PointSpatialModel Temporal model type : Parameters: index : 2.000 +/- 0.00 amplitude : 1.00e-12 +/- 0.0e+00 1 / (cm2 s TeV) reference (frozen): 1.000 TeV lon_0 : 0.000 +/- 0.00 deg lat_0 : 0.500 +/- 0.00 deg Component 1: FoVBackgroundModel Name : my-dataset-bkg Datasets names : ['my-dataset'] Spectral model type : PowerLawNormSpectralModel Parameters: norm : 1.000 +/- 0.00 tilt (frozen): 0.000 reference (frozen): 1.000 TeV .. GENERATED FROM PYTHON SOURCE LINES 262-270 The next step shows how to sample the events with the `~gammapy.datasets.MapDatasetEventSampler` class. The class requests a random number seed generator (that we set with `random_state=0`), the `~gammapy.datasets.Dataset` and the `~gammapy.data.Observations` object. From the latter, the `~gammapy.datasets.MapDatasetEventSampler` class takes all the meta data information. .. GENERATED FROM PYTHON SOURCE LINES 272-275 .. code-block:: Python sampler = MapDatasetEventSampler(random_state=0) events = sampler.run(dataset, observation) .. GENERATED FROM PYTHON SOURCE LINES 276-283 The output of the event-sampler is an event list with coordinates, energies (true and reconstructed) and time of arrivals of the source and background events. `events` is a `~gammapy.data.EventList` object (for details see e.g. :doc:`/tutorials/data/cta` tutorial.). Source and background events are flagged by the MC_ID identifier (where 0 is the default identifier for the background). .. GENERATED FROM PYTHON SOURCE LINES 283-287 .. code-block:: Python print(f"Source events: {(events.table['MC_ID'] == 1).sum()}") print(f"Background events: {(events.table['MC_ID'] == 0).sum()}") .. rst-class:: sphx-glr-script-out .. code-block:: none Source events: 138 Background events: 15319 .. GENERATED FROM PYTHON SOURCE LINES 288-290 We can inspect the properties of the simulated events as follows: .. GENERATED FROM PYTHON SOURCE LINES 290-294 .. code-block:: Python events.select_offset([0, 1] * u.deg).peek() plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_event_sampling_001.png :alt: event sampling :srcset: /tutorials/analysis-3d/images/sphx_glr_event_sampling_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 295-299 By default, the `~gammapy.datasets.MapDatasetEventSampler` fills the metadata keyword `OBJECT` in the event list using the first model of the SkyModel object. You can change it with the following commands: .. GENERATED FROM PYTHON SOURCE LINES 299-302 .. code-block:: Python events.table.meta["OBJECT"] = dataset.models[0].name .. GENERATED FROM PYTHON SOURCE LINES 303-306 Let’s write the event list and its GTI extension to a FITS file. We make use of `fits` library in `astropy`: .. GENERATED FROM PYTHON SOURCE LINES 306-314 .. code-block:: Python primary_hdu = fits.PrimaryHDU() hdu_evt = fits.BinTableHDU(events.table) hdu_gti = fits.BinTableHDU(dataset.gti.table, name="GTI") hdu_all = fits.HDUList([primary_hdu, hdu_evt, hdu_gti]) hdu_all.writeto("./event_sampling/events_0001.fits", overwrite=True) .. rst-class:: sphx-glr-script-out .. code-block:: none /home/runner/work/gammapy-docs/gammapy-docs/gammapy/.tox/build_docs/lib/python3.9/site-packages/astropy/io/fits/hdu/table.py:374: AstropyDeprecationWarning: The following keywords are now recognized as special column-related attributes and should be set via the Column objects: TCTYPn, TCUNIn. In future, these values will be dropped from manually specified headers automatically and replaced with values generated based on the Column objects. warnings.warn( .. GENERATED FROM PYTHON SOURCE LINES 315-325 Time variable source using a lightcurve --------------------------------------- The event sampler can also handle temporal variability of the simulated sources. In this example, we show how to sample a source characterized by an exponential decay, with decay time of 2800 seconds, during the observation. First of all, let’s create a lightcurve: .. GENERATED FROM PYTHON SOURCE LINES 325-332 .. code-block:: Python t0 = 2800 * u.s t_ref = Time("2000-01-01T00:01:04.184") times = t_ref + livetime * np.linspace(0, 1, 100) expdecay_model = ExpDecayTemporalModel(t_ref=t_ref.mjd * u.d, t0=t0) .. GENERATED FROM PYTHON SOURCE LINES 333-338 where we defined the time axis starting from the reference time `t_ref` up to the requested exposure (`livetime`). The bin size of the time-axis is quite arbitrary but, as above for spatial and energy binning, the finer the better. .. GENERATED FROM PYTHON SOURCE LINES 340-343 Then, we can create the sky model. Just for the sake of the example, let’s boost the flux of the simulated source of an order of magnitude: .. GENERATED FROM PYTHON SOURCE LINES 343-361 .. code-block:: Python spectral_model_pwl.amplitude.value = 2e-11 sky_model_pntpwl = SkyModel( spectral_model=spectral_model_pwl, spatial_model=spatial_model_point, temporal_model=expdecay_model, name="point-pwl", ) bkg_model = FoVBackgroundModel(dataset_name="my-dataset") models = Models([sky_model_pntpwl, bkg_model]) file_model = "./event_sampling/point-pwl_decay.yaml" models.write(file_model, overwrite=True) .. GENERATED FROM PYTHON SOURCE LINES 362-365 For simplicity, we use the same dataset defined for the previous example: .. GENERATED FROM PYTHON SOURCE LINES 365-369 .. code-block:: Python dataset.models = models print(dataset.models) .. rst-class:: sphx-glr-script-out .. code-block:: none DatasetModels Component 0: SkyModel Name : point-pwl Datasets names : None Spectral model type : PowerLawSpectralModel Spatial model type : PointSpatialModel Temporal model type : ExpDecayTemporalModel Parameters: index : 2.000 +/- 0.00 amplitude : 2.00e-11 +/- 0.0e+00 1 / (cm2 s TeV) reference (frozen): 1.000 TeV lon_0 : 0.000 +/- 0.00 deg lat_0 : 0.500 +/- 0.00 deg t0 : 2800.000 +/- 0.00 s t_ref (frozen): 51544.001 d Component 1: FoVBackgroundModel Name : my-dataset-bkg Datasets names : ['my-dataset'] Spectral model type : PowerLawNormSpectralModel Parameters: norm : 1.000 +/- 0.00 tilt (frozen): 0.000 reference (frozen): 1.000 TeV .. GENERATED FROM PYTHON SOURCE LINES 370-372 And now, let’s simulate the variable source: .. GENERATED FROM PYTHON SOURCE LINES 374-380 .. code-block:: Python sampler = MapDatasetEventSampler(random_state=0) events = sampler.run(dataset, observation) print(f"Source events: {(events.table['MC_ID'] == 1).sum()}") print(f"Background events: {(events.table['MC_ID'] == 0).sum()}") .. rst-class:: sphx-glr-script-out .. code-block:: none Source events: 1523 Background events: 15246 .. GENERATED FROM PYTHON SOURCE LINES 381-385 We can now inspect the properties of the simulated source. To do that, we adopt the `select_region` function that extracts only the events into a given `SkyRegion` of a `~gammapy.data.EventList` object: .. GENERATED FROM PYTHON SOURCE LINES 385-393 .. code-block:: Python src_position = SkyCoord(0.0, 0.5, frame="galactic", unit="deg") on_region_radius = Angle("0.15 deg") on_region = CircleSkyRegion(center=src_position, radius=on_region_radius) src_events = events.select_region(on_region) .. GENERATED FROM PYTHON SOURCE LINES 394-396 Then we can have a quick look to the data with the `peek` function: .. GENERATED FROM PYTHON SOURCE LINES 396-400 .. code-block:: Python src_events.peek() plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_event_sampling_002.png :alt: event sampling :srcset: /tutorials/analysis-3d/images/sphx_glr_event_sampling_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 401-404 In the right figure of the bottom panel, it is shown the source lightcurve that follows a decay trend as expected. .. GENERATED FROM PYTHON SOURCE LINES 406-416 Extended source using a template -------------------------------- The event sampler can also work with a template model. Here we use the interstellar emission model map of the Fermi 3FHL, which can be found in the `$GAMMAPY_DATA` repository. We proceed following the same steps showed above and we finally have a look at the event’s properties: .. GENERATED FROM PYTHON SOURCE LINES 416-438 .. code-block:: Python template_model = TemplateSpatialModel.read( "$GAMMAPY_DATA/fermi-3fhl-gc/gll_iem_v06_gc.fits.gz", normalize=False ) # we make the model brighter artificially so that it becomes visible over the background diffuse = SkyModel( spectral_model=PowerLawNormSpectralModel(norm=5), spatial_model=template_model, name="template-model", ) bkg_model = FoVBackgroundModel(dataset_name="my-dataset") models_diffuse = Models([diffuse, bkg_model]) file_model = "./event_sampling/diffuse.yaml" models_diffuse.write(file_model, overwrite=True) dataset.models = models_diffuse print(dataset.models) .. rst-class:: sphx-glr-script-out .. code-block:: none DatasetModels Component 0: SkyModel Name : template-model Datasets names : None Spectral model type : PowerLawNormSpectralModel Spatial model type : TemplateSpatialModel Temporal model type : Parameters: norm : 5.000 +/- 0.00 tilt (frozen): 0.000 reference (frozen): 1.000 TeV lon_0 (frozen): 0.000 deg lat_0 (frozen): -0.062 deg Component 1: FoVBackgroundModel Name : my-dataset-bkg Datasets names : ['my-dataset'] Spectral model type : PowerLawNormSpectralModel Parameters: norm : 1.000 +/- 0.00 tilt (frozen): 0.000 reference (frozen): 1.000 TeV .. GENERATED FROM PYTHON SOURCE LINES 439-445 .. code-block:: Python sampler = MapDatasetEventSampler(random_state=0) events = sampler.run(dataset, observation) events.select_offset([0, 1] * u.deg).peek() plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_event_sampling_003.png :alt: event sampling :srcset: /tutorials/analysis-3d/images/sphx_glr_event_sampling_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 446-457 Simulate multiple event lists ----------------------------- In some user case, you may want to sample events from a number of observations. In this section, we show how to simulate a set of event lists. For simplicity, we consider only one point-like source, observed three times for 1 hr and assuming the same pointing position. Let’s firstly define the time start and the livetime of each observation: .. GENERATED FROM PYTHON SOURCE LINES 457-461 .. code-block:: Python tstarts = Time("2020-01-01 00:00:00") + [1, 5, 7] * u.hr livetimes = [1, 1, 1] * u.hr .. GENERATED FROM PYTHON SOURCE LINES 462-486 .. code-block:: Python n_obs = len(tstarts) irf_paths = [path / irf_filename] * n_obs events_paths = [] for idx, tstart in enumerate(tstarts): irfs = load_irf_dict_from_file(irf_paths[idx]) observation = Observation.create( obs_id=idx, pointing=pointing, tstart=tstart, livetime=livetimes[idx], irfs=irfs, location=location, ) dataset = maker.run(empty, observation) dataset.models = models sampler = MapDatasetEventSampler(random_state=idx) events = sampler.run(dataset, observation) path = Path(f"./event_sampling/events_{idx:04d}.fits") events_paths.append(path) events.table.write(path, overwrite=True) .. rst-class:: sphx-glr-script-out .. code-block:: none /home/runner/work/gammapy-docs/gammapy-docs/gammapy/.tox/build_docs/lib/python3.9/site-packages/gammapy/utils/random/inverse_cdf.py:36: RuntimeWarning: invalid value encountered in divide pdf = pdf.ravel() / pdf.sum() /home/runner/work/gammapy-docs/gammapy-docs/gammapy/.tox/build_docs/lib/python3.9/site-packages/astropy/units/quantity.py:673: RuntimeWarning: invalid value encountered in divide result = super().__array_ufunc__(function, method, *arrays, **kwargs) /home/runner/work/gammapy-docs/gammapy-docs/gammapy/.tox/build_docs/lib/python3.9/site-packages/gammapy/utils/random/inverse_cdf.py:36: RuntimeWarning: invalid value encountered in divide pdf = pdf.ravel() / pdf.sum() /home/runner/work/gammapy-docs/gammapy-docs/gammapy/.tox/build_docs/lib/python3.9/site-packages/astropy/units/quantity.py:673: RuntimeWarning: invalid value encountered in divide result = super().__array_ufunc__(function, method, *arrays, **kwargs) /home/runner/work/gammapy-docs/gammapy-docs/gammapy/.tox/build_docs/lib/python3.9/site-packages/gammapy/utils/random/inverse_cdf.py:36: RuntimeWarning: invalid value encountered in divide pdf = pdf.ravel() / pdf.sum() /home/runner/work/gammapy-docs/gammapy-docs/gammapy/.tox/build_docs/lib/python3.9/site-packages/astropy/units/quantity.py:673: RuntimeWarning: invalid value encountered in divide result = super().__array_ufunc__(function, method, *arrays, **kwargs) .. GENERATED FROM PYTHON SOURCE LINES 487-490 You can now load the event list and the corresponding IRFs with `DataStore.from_events_files`: .. GENERATED FROM PYTHON SOURCE LINES 490-497 .. code-block:: Python path = Path("./event_sampling/") events_paths = list(path.rglob("events*.fits")) data_store = DataStore.from_events_files(events_paths, irf_paths) display(data_store.obs_table) .. rst-class:: sphx-glr-script-out .. code-block:: none OBS_ID TSTART TSTOP ONTIME ... TIMEUNIT TIMESYS TIMEREF s s s ... ------ ----------- ----------- ----------------- ... -------- ------- ------- 1 631170005.0 631173605.0 3600.000000000002 ... s utc LOCAL 0 631155605.0 631159205.0 3600.000000000006 ... s utc LOCAL 2 631177205.0 631180805.0 3599.999999999997 ... s utc LOCAL .. GENERATED FROM PYTHON SOURCE LINES 498-502 Then you can create the observations from the data store and make your own analysis following the instructions in the :doc:`/tutorials/starting/analysis_2` tutorial. .. GENERATED FROM PYTHON SOURCE LINES 502-507 .. code-block:: Python observations = data_store.get_observations() observations[0].peek() plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_event_sampling_004.png :alt: Effective area, Background rate, Energy dispersion, Point spread function, Events :srcset: /tutorials/analysis-3d/images/sphx_glr_event_sampling_004.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 508-516 Exercises --------- - Try to sample events for an extended source (e.g. a radial gaussian morphology); - Change the spatial model and the spectrum of the simulated Sky model; - Include a temporal model in the simulation .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 22.916 seconds) .. _sphx_glr_download_tutorials_analysis-3d_event_sampling.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: binder-badge .. image:: images/binder_badge_logo.svg :target: https://mybinder.org/v2/gh/gammapy/gammapy-webpage/main?urlpath=lab/tree/notebooks/dev/tutorials/analysis-3d/event_sampling.ipynb :alt: Launch binder :width: 150 px .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: event_sampling.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: event_sampling.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_