Explore Calibrated Data

import numpy as np
from astropy import units as u
from matplotlib import pyplot as plt

import ctapipe
from ctapipe.calib import CameraCalibrator
from ctapipe.image import hillas_parameters, tailcuts_clean
from ctapipe.io import EventSource
from ctapipe.utils.datasets import get_dataset_path
from ctapipe.visualization import ArrayDisplay, CameraDisplay

# %matplotlib inline
plt.style.use("ggplot")

print(ctapipe.__version__)
print(ctapipe.__file__)

0.27.1.dev15+g1af6977ab
/home/docs/checkouts/readthedocs.org/user_builds/ctapipe/envs/2868/lib/python3.10/site-packages/ctapipe/__init__.py

Let’s first open a raw event file and get an event out of it:

filename = get_dataset_path("gamma_prod5.simtel.zst")
source = EventSource(filename, max_events=2)

for event in source:
    print(event.index.event_id)

Downloading gamma_prod5.simtel.zst:   0%|          | 0.00/17.5M [00:00<?, ?B/s]
Downloading gamma_prod5.simtel.zst:   0%|          | 10.2k/17.5M [00:00<03:27, 84.2kB/s]
Downloading gamma_prod5.simtel.zst:   0%|          | 41.0k/17.5M [00:00<01:34, 184kB/s]
Downloading gamma_prod5.simtel.zst:   1%|          | 92.2k/17.5M [00:00<00:59, 291kB/s]
Downloading gamma_prod5.simtel.zst:   1%|          | 215k/17.5M [00:00<00:30, 572kB/s]
Downloading gamma_prod5.simtel.zst:   3%|▎         | 440k/17.5M [00:00<00:16, 1.03MB/s]
Downloading gamma_prod5.simtel.zst:   5%|▌         | 891k/17.5M [00:00<00:08, 1.93MB/s]
Downloading gamma_prod5.simtel.zst:  10%|█         | 1.79M/17.5M [00:00<00:04, 3.69MB/s]
Downloading gamma_prod5.simtel.zst:  21%|██        | 3.60M/17.5M [00:00<00:01, 7.21MB/s]
Downloading gamma_prod5.simtel.zst:  38%|███▊      | 6.60M/17.5M [00:01<00:00, 12.6MB/s]
Downloading gamma_prod5.simtel.zst:  52%|█████▏    | 9.09M/17.5M [00:01<00:00, 14.9MB/s]
Downloading gamma_prod5.simtel.zst:  72%|███████▏  | 12.5M/17.5M [00:01<00:00, 19.0MB/s]
Downloading gamma_prod5.simtel.zst:  89%|████████▉ | 15.6M/17.5M [00:01<00:00, 20.8MB/s]
Downloading gamma_prod5.simtel.zst: 100%|██████████| 17.5M/17.5M [00:01<00:00, 11.8MB/s]
4009
5101

filename

PosixPath('/home/docs/.cache/ctapipe/minio-cta.zeuthen.desy.de/dpps-testdata-public/data/ctapipe-test-data/v1.1.0/gamma_prod5.simtel.zst')

source

SimTelEventSource

Read events from a SimTelArray data file (in EventIO format).

ctapipe makes use of the sim_telarray metadata system to fill some information not directly accessible from the data inside the files itself. Make sure you set this parameters in the simulation configuration to fully make use of ctapipe. In future, ctapipe might require these metadata parameters.

This includes:

• Reference Point of the telescope coordinates. Make sure to include the

  user-defined parameters LONGITUDE and LATITUDE with the geodetic coordinates of the array reference point. Also make sure the ALTITUDE config parameter is included in the global metadata.

• Names of the optical structures and the cameras are read from

  OPTICS_CONFIG_NAME and CAMERA_CONFIG_NAME, make sure to include these in the telescope meta.

• The MIRROR_CLASS should also be included in the telescope meta

  to correctly setup coordinate transforms.

If these parameters are not included in the input data, ctapipe will fallback guesses these based on available data and the list of known telescopes for ctapipe.instrument.guess_telescope.

allowed_tels	None	list of allowed tel_ids, others will be ignored. If None, all telescopes in the input stream will be included (default: None)
atmosphere_profile_choice	AtmosphereProfileKind.AUTO	Which type of atmosphere density profile to load from the file, in case more than one exists. If set to AUTO, TABLE will be attempted first and if missing, FIVELAYER will be loaded. (default: AtmosphereProfileKind.AUTO)
back_seekable	False	Require the event source to be backwards seekable. This will reduce in slower read speed for gzipped files and is not possible for zstd compressed files (default: False)
calib_scale	1.0	Factor to transform ADC counts into number of photoelectrons. Corrects the DC_to_PHE factor. (default: 1.0)
calib_shift	0.0	Factor to shift the R1 photoelectron samples. Can be used to simulate mis-calibration. (default: 0.0)
focal_length_choice	FocalLengthKind.EFFECTIVE	If both nominal and effective focal lengths are available in the SimTelArray file, which one to use for the `~ctapipe.coordinates.CameraFrame` attached to the `~ctapipe.instrument.CameraGeometry` instances in the `~ctapipe.instrument.SubarrayDescription`, which will be used in CameraFrame to TelescopeFrame coordinate transforms. The 'nominal' focal length is the one used during the simulation, the 'effective' focal length is computed using specialized ray-tracing from a point light source (default: FocalLengthKind.EFFECTIVE)
gain_selector_type	ThresholdGainSelector	The name of a GainSelector subclass. Possible values: ['ManualGainSelector', 'ThresholdGainSelector'] (default: ThresholdGainSelector)
input_url	/home/docs/.cache/ctapipe/minio-cta.zeuthen.desy.de/dpps-testdata-public/data/ctapipe-test-data/v1.1.0/gamma_prod5.simtel.zst	Path to the input file containing events. (default: None)
max_events	2	Maximum number of events that will be read from the file (default: None)
override_obs_id	None	Use the given obs_id instead of the run number from sim_telarray (default: None)
select_gain	None	Whether to perform gain selection. The default (None) means only select gain of physics events, not of calibration events. (default: None)
skip_calibration_events	True	Skip calibration events (default: True)
skip_r1_calibration	False	Skip the simtel category-A calibration. The R1 waveforms (in ADC counts) will be the same as the R0 waveforms. (default: False)

event

ctapipe.containers.ArrayEventContainer:
                       index.*: event indexing information with default None
                          r0.*: Raw Data with default None
                          r1.*: R1 Calibrated Data with default None
                         dl0.*: DL0 Data Volume Reduced Data with default None
                         dl1.*: DL1 Calibrated image with default None
                         dl2.*: DL2 reconstruction info with default None
                  simulation.*: Simulated Event Information with default None
                                with type <class
                                'ctapipe.containers.SimulatedEventContainer'>
                     trigger.*: central trigger information with default None
                         count: number of events processed with default 0
                  monitoring.*: container for event-wise monitoring data (MON)
                                with default None
                        muon.*: Container for muon analysis results with default
                                None

print(event.r1)

{'tel': {26: {'calibration_monitoring_id': None,
              'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
              'event_type': <EventType.SUBARRAY: 32>,
              'first_cell_id': None,
              'module_hires_local_clock_counter': None,
              'pedestal_intensity': None,
              'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(1764,), dtype=uint8),
              'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(1764,), dtype=int8),
              'waveform': array([[[ 0.5449172 ,  0.7336623 ,  0.85949236, ..., -0.08423293,
          0.10451212,  0.37714386],
        [ 0.26182497,  0.15617064,  0.26182497, ..., -0.66793317,
         -0.71019495, -0.71019495],
        [-0.5535298 , -0.15927683,  0.27647644, ...,  0.62922907,
          0.64997923,  0.64997923],
        ...,
        [ 0.00745648, -0.11931664, -0.28834748, ...,  0.00745648,
          0.21874502,  0.366647  ],
        [ 0.06618547, -0.12305467, -0.24921475, ...,  0.12926552,
          0.08721215,  0.12926552],
        [-0.23363174, -0.19188231, -0.15013288, ..., -0.08750875,
         -0.04575931, -0.12925817]]], shape=(1, 1764, 25), dtype=float32)},
         27: {'calibration_monitoring_id': None,
              'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
              'event_type': <EventType.SUBARRAY: 32>,
              'first_cell_id': None,
              'module_hires_local_clock_counter': None,
              'pedestal_intensity': None,
              'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(1764,), dtype=uint8),
              'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(1764,), dtype=int8),
              'waveform': array([[[ 0.8074397 ,  0.88952863,  0.78691745, ...,  0.31490615,
          0.04811714,  0.10968383],
        [-0.07342623, -0.07342623,  0.0933422 , ...,  0.00995799,
         -0.03173412, -0.19850256],
        [ 0.09717706,  0.05602828, -0.06741804, ..., -0.02626926,
         -0.21143875, -0.19086435],
        ...,
        [-0.09182719, -0.1328717 , -0.11234944, ...,  0.35966238,
          0.48279592,  0.5238404 ],
        [ 0.25864878,  0.31945735,  0.23837928, ...,  0.7045781 ,
          0.82619524,  0.7451171 ],
        [-0.16240191, -0.225794  , -0.18353261, ..., -0.07787914,
         -0.16240191, -0.2046633 ]]], shape=(1, 1764, 25), dtype=float32)},
         29: {'calibration_monitoring_id': None,
              'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
              'event_type': <EventType.SUBARRAY: 32>,
              'first_cell_id': None,
              'module_hires_local_clock_counter': None,
              'pedestal_intensity': None,
              'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(1764,), dtype=uint8),
              'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(1764,), dtype=int8),
              'waveform': array([[[-0.04938218,  0.03202369,  0.3372957 , ..., -0.13078806,
         -0.02903072, -0.06973365],
        [ 0.59601194,  0.32254168,  0.02803528, ..., -0.24543495,
         -0.24543495, -0.37165198],
        [-0.24841212, -0.1858834 , -0.01914014, ..., -0.35262665,
         -0.39431247, -0.4568412 ],
        ...,
        [ 0.05627692,  0.09927555,  0.12077487, ..., -0.5027053 ,
         -0.58870256, -0.61020184],
        [ 0.3565684 ,  0.47816813,  0.5592346 , ..., -0.17036383,
         -0.23116371, -0.3122302 ],
        [ 0.16858962,  0.08948522,  0.10926133, ...,  0.46523112,
          0.40590283,  0.36635062]]], shape=(1, 1764, 25), dtype=float32)},
         47: {'calibration_monitoring_id': None,
              'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
              'event_type': <EventType.SUBARRAY: 32>,
              'first_cell_id': None,
              'module_hires_local_clock_counter': None,
              'pedestal_intensity': None,
              'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(2048,), dtype=uint8),
              'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(2048,), dtype=int8),
              'waveform': array([[[-0.04045208,  0.04424386, -0.12514801, ..., -0.04045208,
          0.00189589,  0.00189589],
        [-0.07548809, -0.02983804, -0.02983804, ..., -0.02983804,
         -0.02983804, -0.12113815],
        [-0.00396833, -0.04704034, -0.17625633, ..., -0.17625633,
         -0.04704034, -0.13318433],
        ...,
        [-0.04439949, -0.04439949,  0.04724088, ..., -0.18186006,
         -0.22768025, -0.09021968],
        [-0.04938671,  0.04159823, -0.04938671, ..., -0.18586414,
         -0.00389424, -0.00389424],
        [-0.13006851, -0.17539828, -0.039409  , ..., -0.039409  ,
         -0.08473875, -0.08473875]]], shape=(1, 2048, 128), dtype=float32)},
         53: {'calibration_monitoring_id': None,
              'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
              'event_type': <EventType.SUBARRAY: 32>,
              'first_cell_id': None,
              'module_hires_local_clock_counter': None,
              'pedestal_intensity': None,
              'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(2048,), dtype=uint8),
              'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(2048,), dtype=int8),
              'waveform': array([[[-0.0082417 , -0.0082417 ,  0.03595864, ...,  0.03595864,
         -0.0082417 ,  0.08015899],
        [ 0.00310245,  0.00310245,  0.08867562, ...,  0.00310245,
          0.1314622 ,  0.04588903],
        [ 0.00668868,  0.09877883,  0.09877883, ...,  0.05273375,
         -0.13144656, -0.0393564 ],
        ...,
        [-0.0738365 , -0.16169758, -0.0738365 , ...,  0.10188565,
          0.10188565,  0.01402458],
        [ 0.19253932,  0.14814104,  0.14814104, ...,  0.14814104,
         -0.07385035,  0.10374276],
        [-0.16682361,  0.0083629 , -0.07923036, ..., -0.03543373,
         -0.03543373, -0.16682361]]], shape=(1, 2048, 128), dtype=float32)},
         59: {'calibration_monitoring_id': None,
              'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
              'event_type': <EventType.SUBARRAY: 32>,
              'first_cell_id': None,
              'module_hires_local_clock_counter': None,
              'pedestal_intensity': None,
              'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(2048,), dtype=uint8),
              'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(2048,), dtype=int8),
              'waveform': array([[[-0.1164972 , -0.03435972, -0.07542846, ...,  0.00670901,
          0.04777775,  0.17098396],
        [-0.19899125, -0.02344115, -0.11121619, ..., -0.11121619,
         -0.06732867, -0.06732867],
        [-0.0580824 ,  0.03068624,  0.11945489, ...,  0.07507057,
          0.03068624,  0.11945489],
        ...,
        [-0.02629879, -0.0717934 ,  0.06469042, ..., -0.117288  ,
         -0.02629879, -0.02629879],
        [-0.15135607, -0.05747958, -0.05747958, ...,  0.03639691,
          0.1302734 , -0.05747958],
        [-0.00660102, -0.14327039, -0.09771393, ...,  0.08451191,
         -0.00660102, -0.09771393]]], shape=(1, 2048, 128), dtype=float32)},
         69: {'calibration_monitoring_id': None,
              'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
              'event_type': <EventType.SUBARRAY: 32>,
              'first_cell_id': None,
              'module_hires_local_clock_counter': None,
              'pedestal_intensity': None,
              'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(2048,), dtype=uint8),
              'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(2048,), dtype=int8),
              'waveform': array([[[ 0.00841863,  0.04984551,  0.00841863, ...,  0.00841863,
          0.04984551,  0.09127238],
        [ 0.01125244, -0.06966043, -0.06966043, ..., -0.15057331,
          0.05170888, -0.029204  ],
        [ 0.2630286 ,  0.21625374,  0.21625374, ...,  0.12270404,
          0.07592919,  0.21625374],
        ...,
        [-0.03710477,  0.05580388,  0.00934956, ...,  0.24162118,
          0.14871253,  0.19516686],
        [-0.01875117, -0.06278115, -0.06278115, ..., -0.06278115,
         -0.01875117, -0.06278115],
        [-0.11524381, -0.11524381, -0.20668873, ...,  0.25053585,
          0.02192356,  0.15909094]]], shape=(1, 2048, 128), dtype=float32)},
         73: {'calibration_monitoring_id': None,
              'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
              'event_type': <EventType.SUBARRAY: 32>,
              'first_cell_id': None,
              'module_hires_local_clock_counter': None,
              'pedestal_intensity': None,
              'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(2048,), dtype=uint8),
              'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(2048,), dtype=int8),
              'waveform': array([[[-0.03057504,  0.01178788, -0.07293797, ...,  0.01178788,
          0.09651372,  0.01178788],
        [-0.08434547, -0.08434547, -0.1278464 , ...,  0.00265641,
         -0.08434547,  0.00265641],
        [ 0.09269284,  0.13629547,  0.00548759, ...,  0.00548759,
          0.04909021,  0.1798981 ],
        ...,
        [-0.05265228, -0.18097262, -0.09542572, ..., -0.09542572,
         -0.13819917, -0.05265228],
        [ 0.07313758,  0.07313758,  0.03256791, ..., -0.0891411 ,
         -0.00800176, -0.0891411 ],
        [-0.1561605 ,  0.06265412, -0.1561605 , ...,  0.06265412,
         -0.02487173, -0.06863465]]], shape=(1, 2048, 128), dtype=float32)},
         121: {'calibration_monitoring_id': None,
               'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
               'event_type': <EventType.SUBARRAY: 32>,
               'first_cell_id': None,
               'module_hires_local_clock_counter': None,
               'pedestal_intensity': None,
               'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(1855,), dtype=uint8),
               'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(1855,), dtype=int8),
               'waveform': array([[[-0.21347916, -0.25365788, -0.2938366 , ...,  0.52982706,
         -0.03267494, -0.11303237],
        [ 0.01376745,  0.0734081 ,  0.25233003, ..., -0.34407645,
         -0.18503472, -0.12539406],
        [-0.35374358, -0.25606352, -0.17791945, ..., -0.02163133,
          0.2909449 , -0.00209532],
        ...,
        [ 0.1588768 , -0.01432226,  0.35132024, ..., -0.22601001,
         -0.01432226,  0.4860306 ],
        [-0.18021354, -0.14094928,  0.19279698, ..., -0.29800636,
         -0.37653488, -0.396167  ],
        [-0.12698945, -0.0883472 , -0.06902608, ..., -0.32020068,
         -0.10766833,  0.06622178]]], shape=(1, 1855, 60), dtype=float32)},
         122: {'calibration_monitoring_id': None,
               'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
               'event_type': <EventType.SUBARRAY: 32>,
               'first_cell_id': None,
               'module_hires_local_clock_counter': None,
               'pedestal_intensity': None,
               'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(1855,), dtype=uint8),
               'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(1855,), dtype=int8),
               'waveform': array([[[ 0.0766639 ,  0.32452905,  0.3054625 , ..., -0.19026779,
         -0.11400159,  0.15293011],
        [-0.01299111,  0.04359313,  0.08131595, ..., -0.23932804,
         -0.1827438 , -0.14502099],
        [ 0.01621819, -0.1590356 ,  0.07463612, ..., -0.27587143,
         -0.23692615, -0.21745351],
        ...,
        [ 0.25392857,  0.34937772, -0.03241886, ..., -0.2614968 ,
         -0.1851375 , -0.29967648],
        [-0.19018087, -0.0539299 ,  0.00446337, ..., -0.19018087,
          0.04339222,  0.10178549],
        [ 0.13077052,  0.8790263 ,  1.5121658 , ..., -0.15702018,
         -0.31050855, -0.34888065]]], shape=(1, 1855, 60), dtype=float32)},
         124: {'calibration_monitoring_id': None,
               'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
               'event_type': <EventType.SUBARRAY: 32>,
               'first_cell_id': None,
               'module_hires_local_clock_counter': None,
               'pedestal_intensity': None,
               'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(1855,), dtype=uint8),
               'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(1855,), dtype=int8),
               'waveform': array([[[-0.11797108, -0.11797108, -0.09868825, ..., -0.0601226 ,
          0.11342283,  0.11342283],
        [-0.10063321,  0.10081094, -0.06400701, ..., -0.0273808 ,
          0.02755852, -0.10063321],
        [-0.19428314, -0.2135887 , -0.0591442 , ..., -0.23289426,
         -0.02053308, -0.09775533],
        ...,
        [-0.04796196, -0.04796196,  0.30063868, ...,  0.01013814,
          0.4362056 ,  0.33937207],
        [-0.10645148, -0.12572832, -0.12572832, ...,  0.31763902,
          0.39474636,  0.10559376],
        [ 0.50643456,  0.26373467,  0.35708076, ..., -0.128319  ,
         -0.14698821, -0.07231133]]], shape=(1, 1855, 60), dtype=float32)},
         148: {'calibration_monitoring_id': None,
               'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
               'event_type': <EventType.SUBARRAY: 32>,
               'first_cell_id': None,
               'module_hires_local_clock_counter': None,
               'pedestal_intensity': None,
               'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(2048,), dtype=uint8),
               'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(2048,), dtype=int8),
               'waveform': array([[[ 0.00710569, -0.03808828, -0.08328226, ...,  0.00710569,
          0.09749365,  0.00710569],
        [ 0.09487705,  0.09487705,  0.01112803, ..., -0.07262099,
         -0.03074648,  0.05300254],
        [-0.11470292, -0.02782934,  0.10248102, ..., -0.02782934,
         -0.02782934, -0.1581397 ],
        ...,
        [ 0.07964301,  0.0340572 , -0.05711443, ...,  0.07964301,
         -0.01152861,  0.07964301],
        [ 0.00197998,  0.00197998,  0.13376029, ...,  0.08983352,
          0.04590675, -0.08587357],
        [ 0.03801603, -0.00807014,  0.03801603, ..., -0.00807014,
         -0.1002425 , -0.1002425 ]]], shape=(1, 2048, 128), dtype=float32)},
         162: {'calibration_monitoring_id': None,
               'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
               'event_type': <EventType.SUBARRAY: 32>,
               'first_cell_id': None,
               'module_hires_local_clock_counter': None,
               'pedestal_intensity': None,
               'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(2048,), dtype=uint8),
               'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(2048,), dtype=int8),
               'waveform': array([[[-5.6908205e-03,  8.2228467e-02,  8.2228467e-02, ...,
         -4.9650464e-02, -5.6908205e-03,  3.8268823e-02],
        [-3.4884274e-02,  8.1552938e-03, -7.7923842e-02, ...,
         -1.2096341e-01, -7.7923842e-02,  8.1552938e-03],
        [-3.2151666e-02, -7.7254020e-02,  1.0315538e-01, ...,
          5.8053035e-02, -3.2151666e-02,  1.2950684e-02],
        ...,
        [-4.1048851e-02, -8.1529431e-02, -8.1529431e-02, ...,
         -5.6826987e-04, -5.6826987e-04, -1.6249059e-01],
        [ 4.2881422e-02, -2.7393538e-04,  4.2881422e-02, ...,
          8.6036779e-02,  2.1550286e-01,  2.5865820e-01],
        [-3.9726321e-02, -8.2472593e-02, -3.9726321e-02, ...,
          3.0224383e-01,  1.7400502e-01,  1.3125876e-01]]],
      shape=(1, 2048, 128), dtype=float32)},
         166: {'calibration_monitoring_id': None,
               'event_time': <Time object: scale='tai' format='mjd' value=0.0>,
               'event_type': <EventType.SUBARRAY: 32>,
               'first_cell_id': None,
               'module_hires_local_clock_counter': None,
               'pedestal_intensity': None,
               'pixel_status': array([4, 4, 4, ..., 4, 4, 4], shape=(2048,), dtype=uint8),
               'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], shape=(2048,), dtype=int8),
               'waveform': array([[[-0.03838527, -0.08015627, -0.03838527, ...,  0.00338574,
          0.00338574, -0.03838527],
        [-0.08626862, -0.08626862, -0.04286342, ..., -0.08626862,
         -0.08626862, -0.12967381],
        [ 0.05237084, -0.21012151, -0.07887534, ..., -0.07887534,
         -0.03512662, -0.07887534],
        ...,
        [-0.11608443, -0.19859755, -0.07482787, ...,  0.00768524,
         -0.03357131, -0.15734099],
        [-0.00771839,  0.03798387, -0.14482519, ..., -0.09912293,
         -0.09912293, -0.05342066],
        [ 0.17404978,  0.1298405 ,  0.1298405 , ..., -0.04699665,
         -0.00278737, -0.04699665]]], shape=(1, 2048, 128), dtype=float32)}}}

Perform basic calibration:

Here we will use a CameraCalibrator which is just a simple wrapper that runs the three calibraraton and trace-integration phases of the pipeline, taking the data from levels:

R0 → R1 → DL0 → DL1

You could of course do these each separately, by using the classes R1Calibrator, DL0Reducer, and DL1Calibrator. Note that we have not specified any configuration to the CameraCalibrator, so it will be using the default algorithms and thresholds, other than specifying that the product is a “HESSIOR1Calibrator” (hopefully in the near future that will be automatic).

calib = CameraCalibrator(subarray=source.subarray)
calib(event)

Now the r1, dl0 and dl1 containers are filled in the event

• r1.tel[x]: contains the “r1-calibrated” waveforms, after gain-selection, pedestal subtraciton, and gain-correction

• dl0.tel[x]: is the same but with optional data volume reduction (some pixels not filled), in this case this is not performed by default, so it is the same as r1

• dl1.tel[x]: contains the (possibly re-calibrated) waveforms as dl0, but also the time-integrated image that has been calculated using a ImageExtractor (a NeighborPeakWindowSum by default)

for tel_id in event.dl1.tel:
    print("TEL{:03}: {}".format(tel_id, source.subarray.tel[tel_id]))
    print("  - r0  wave shape  : {}".format(event.r0.tel[tel_id].waveform.shape))
    print("  - r1  wave shape  : {}".format(event.r1.tel[tel_id].waveform.shape))
    print("  - dl1 image shape : {}".format(event.dl1.tel[tel_id].image.shape))

TEL026: MST_MST_FlashCam
  - r0  wave shape  : (1, 1764, 25)
  - r1  wave shape  : (1, 1764, 25)
  - dl1 image shape : (1764,)
TEL027: MST_MST_FlashCam
  - r0  wave shape  : (1, 1764, 25)
  - r1  wave shape  : (1, 1764, 25)
  - dl1 image shape : (1764,)
TEL029: MST_MST_FlashCam
  - r0  wave shape  : (1, 1764, 25)
  - r1  wave shape  : (1, 1764, 25)
  - dl1 image shape : (1764,)
TEL047: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (1, 2048, 128)
  - dl1 image shape : (2048,)
TEL053: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (1, 2048, 128)
  - dl1 image shape : (2048,)
TEL059: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (1, 2048, 128)
  - dl1 image shape : (2048,)
TEL069: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (1, 2048, 128)
  - dl1 image shape : (2048,)
TEL073: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (1, 2048, 128)
  - dl1 image shape : (2048,)
TEL121: MST_MST_NectarCam
  - r0  wave shape  : (2, 1855, 60)
  - r1  wave shape  : (1, 1855, 60)
  - dl1 image shape : (1855,)
TEL122: MST_MST_NectarCam
  - r0  wave shape  : (2, 1855, 60)
  - r1  wave shape  : (1, 1855, 60)
  - dl1 image shape : (1855,)
TEL124: MST_MST_NectarCam
  - r0  wave shape  : (2, 1855, 60)
  - r1  wave shape  : (1, 1855, 60)
  - dl1 image shape : (1855,)
TEL148: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (1, 2048, 128)
  - dl1 image shape : (2048,)
TEL162: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (1, 2048, 128)
  - dl1 image shape : (2048,)
TEL166: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (1, 2048, 128)
  - dl1 image shape : (2048,)

Some image processing:

Let’s look at the image

tel_id = sorted(event.r1.tel.keys())[1]
sub = source.subarray
geometry = sub.tel[tel_id].camera.geometry
image = event.dl1.tel[tel_id].image

disp = CameraDisplay(geometry, image=image)

mask = tailcuts_clean(
    geometry,
    image,
    picture_thresh=10,
    boundary_thresh=5,
    min_number_picture_neighbors=2,
)
cleaned = image.copy()
cleaned[~mask] = 0
disp = CameraDisplay(geometry, image=cleaned)

params = hillas_parameters(geometry, cleaned)
print(params)
params

{'intensity': np.float64(201.13648319244385),
 'kurtosis': np.float64(1.7664608543531792),
 'length': <Quantity 0.0862792 m>,
 'length_uncertainty': <Quantity 0.00266303 m>,
 'phi': <Angle -0.42975131 rad>,
 'psi': <Angle -0.30256932 rad>,
 'psi_uncertainty': <Angle 0.05077488 rad>,
 'r': <Quantity 0.78989421 m>,
 'skewness': np.float64(0.01211337763129691),
 'transverse_cog_uncertainty': <Quantity 0.00321526 m>,
 'width': <Quantity 0.02976784 m>,
 'width_uncertainty': <Quantity 0.00125979 m>,
 'x': <Quantity 0.71806864 m>,
 'y': <Quantity -0.32910527 m>}

ctapipe.containers.CameraHillasParametersContainer:
                     intensity: total intensity (size) with default nan
                      skewness: measure of the asymmetry with default nan
                      kurtosis: measure of the tailedness with default nan
                             x: centroid x coordinate with default nan m [m]
                             y: centroid x coordinate with default nan m [m]
                             r: radial coordinate of centroid with default nan m
                                [m]
                           phi: polar coordinate of centroid with default nan
                                deg [deg]
                        length: standard deviation along the major-axis with
                                default nan m [m]
            length_uncertainty: uncertainty of length with default nan m [m]
                         width: standard spread along the minor-axis with
                                default nan m [m]
             width_uncertainty: uncertainty of width with default nan m [m]
                           psi: rotation angle of ellipse with default nan deg
                                [deg]
               psi_uncertainty: uncertainty of psi with default nan deg [deg]
    transverse_cog_uncertainty: uncertainty on the center of gravity along the
                                transverse axis of the image with default nan m
                                [m]

params = hillas_parameters(geometry, cleaned)

plt.figure(figsize=(10, 10))
disp = CameraDisplay(geometry, image=image)
disp.add_colorbar()
disp.overlay_moments(params, color="red", lw=3)
disp.highlight_pixels(mask, color="white", alpha=0.3, linewidth=2)

plt.xlim(params.x.to_value(u.m) - 0.5, params.x.to_value(u.m) + 0.5)
plt.ylim(params.y.to_value(u.m) - 0.5, params.y.to_value(u.m) + 0.5)

(-0.8291052663127615, 0.1708947336872385)

source.metadata

{'is_simulation': False}

More complex image processing:

Let’s now explore how stereo reconstruction works.

first, look at a summed image from multiple telescopes

For this, we want to use a CameraDisplay again, but since we can’t sum and display images with different cameras, we’ll just sub-select images from a particular camera type

These are the telescopes that are in this event:

tels_in_event = set(
    event.dl1.tel.keys()
)  # use a set here, so we can intersect it later
tels_in_event

{162, 69, 122, 166, 59, 73, 47, 148, 53, 121, 26, 27, 124, 29}

cam_ids = set(sub.get_tel_ids_for_type("MST_MST_NectarCam"))
cam_ids

{128, 129, 130, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124}

cams_in_event = tels_in_event.intersection(cam_ids)
first_tel_id = list(cams_in_event)[0]
tel = sub.tel[first_tel_id]
print("{}s in event: {}".format(tel, cams_in_event))

MST_MST_NectarCams in event: {121, 122, 124}

Now, let’s sum those images:

image_sum = np.zeros_like(
    tel.camera.geometry.pix_x.value
)  # just make an array of 0's in the same shape as the camera

for tel_id in cams_in_event:
    image_sum += event.dl1.tel[tel_id].image

And finally display the sum of those images

plt.figure(figsize=(8, 8))

disp = CameraDisplay(tel.camera.geometry, image=image_sum)
disp.overlay_moments(params, with_label=False)
plt.title("Sum of {}x {}".format(len(cams_in_event), tel))

Text(0.5, 1.0, 'Sum of 3x MST_MST_NectarCam')

let’s also show which telescopes those were. Note that currently ArrayDisplay’s value field is a vector by tel_index, not tel_id, so we have to convert to a tel_index. (this may change in a future version to be more user-friendly)

nectarcam_subarray = sub.select_subarray(cam_ids, name="NectarCam")

hit_pattern = np.zeros(shape=nectarcam_subarray.n_tels)
hit_pattern[[nectarcam_subarray.tel_indices[x] for x in cams_in_event]] = 100

plt.set_cmap(plt.cm.Accent)
plt.figure(figsize=(8, 8))

ad = ArrayDisplay(nectarcam_subarray)
ad.values = hit_pattern
ad.add_labels()

• 
•