Explore Calibrated Data

[1]:
import ctapipe
from ctapipe.utils.datasets import get_dataset_path
from ctapipe.io import EventSource, EventSeeker
from ctapipe.visualization import CameraDisplay
from ctapipe.instrument import CameraGeometry
from matplotlib import pyplot as plt
from astropy import units as u
import numpy as np

%matplotlib inline
plt.style.use("ggplot")
[2]:
print(ctapipe.__version__)
print(ctapipe.__file__)
0.16.1.dev2+g02147768
/home/runner/work/ctapipe/ctapipe/ctapipe/__init__.py

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

[3]:
filename = get_dataset_path("gamma_prod5.simtel.zst")
source = EventSource(filename, max_events=2)

for event in source:
    print(event.index.event_id)
4009
5101
[4]:
filename
[4]:
PosixPath('/home/runner/.cache/ctapipe/cccta-dataserver.in2p3.fr/data/ctapipe-extra/v0.3.4/gamma_prod5.simtel.zst')
[5]:
source
[5]:
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 avaible 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)
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 `CameraFrame` attached to the `CameraGeometry` instances in the `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 GainSelector to use. (default: ThresholdGainSelector)
input_url /home/runner/.cache/ctapipe/cccta-dataserver.in2p3.fr/data/ctapipe-extra/v0.3.4/gamma_prod5.simtel.zst Path to the input file containing events. (default: traitlets.Undefined)
max_events 2 Maximum number of events that will be read from the file (default: None)
skip_calibration_events True Skip calibration events (default: True)
[6]:
event
[6]:
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
                    pointing.*: Array and telescope pointing positions with
                                default None
                 calibration.*: Container for calibration coefficients for the
                                current event with default None
                         mon.*: container for event-wise monitoring data (MON)
                                with default None
[7]:
print(event.r1)
{'tel': {26: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         27: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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 ]], dtype=float32)},
         29: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         47: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         53: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         59: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         69: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         73: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         121: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         122: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         124: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         148: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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 ]], dtype=float32)},
         162: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], dtype=float32)},
         166: {'selected_gain_channel': array([0, 0, 0, ..., 0, 0, 0], 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]], 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:

R0R1DL0DL1

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).

[8]:
from ctapipe.calib import CameraCalibrator

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)

[9]:
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  : (1764, 25)
  - dl1 image shape : (1764,)
TEL027: MST_MST_FlashCam
  - r0  wave shape  : (1, 1764, 25)
  - r1  wave shape  : (1764, 25)
  - dl1 image shape : (1764,)
TEL029: MST_MST_FlashCam
  - r0  wave shape  : (1, 1764, 25)
  - r1  wave shape  : (1764, 25)
  - dl1 image shape : (1764,)
TEL047: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (2048, 128)
  - dl1 image shape : (2048,)
TEL053: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (2048, 128)
  - dl1 image shape : (2048,)
TEL059: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (2048, 128)
  - dl1 image shape : (2048,)
TEL069: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (2048, 128)
  - dl1 image shape : (2048,)
TEL073: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (2048, 128)
  - dl1 image shape : (2048,)
TEL121: MST_MST_NectarCam
  - r0  wave shape  : (2, 1855, 60)
  - r1  wave shape  : (1855, 60)
  - dl1 image shape : (1855,)
TEL122: MST_MST_NectarCam
  - r0  wave shape  : (2, 1855, 60)
  - r1  wave shape  : (1855, 60)
  - dl1 image shape : (1855,)
TEL124: MST_MST_NectarCam
  - r0  wave shape  : (2, 1855, 60)
  - r1  wave shape  : (1855, 60)
  - dl1 image shape : (1855,)
TEL148: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (2048, 128)
  - dl1 image shape : (2048,)
TEL162: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (2048, 128)
  - dl1 image shape : (2048,)
TEL166: SST_ASTRI_CHEC
  - r0  wave shape  : (1, 2048, 128)
  - r1  wave shape  : (2048, 128)
  - dl1 image shape : (2048,)

Some image processing:

Let’s look at the image

[10]:
from ctapipe.visualization import CameraDisplay

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
[11]:
disp = CameraDisplay(geometry, image=image)
../_images/tutorials_calibrated_data_exploration_15_0.png
[12]:
from ctapipe.image import tailcuts_clean, hillas_parameters
[13]:
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)
../_images/tutorials_calibrated_data_exploration_17_0.png
[14]:
params = hillas_parameters(geometry, cleaned)
print(params)
params
{'intensity': 201.13648319244385,
 'kurtosis': 1.7664608543531797,
 'length': <Quantity 0.0862792 m>,
 'length_uncertainty': <Quantity 0.00266303 m>,
 'phi': <Angle -0.42975131 rad>,
 'psi': <Angle -0.30256932 rad>,
 'r': <Quantity 0.78989421 m>,
 'skewness': 0.012113377631296642,
 'width': <Quantity 0.02976784 m>,
 'width_uncertainty': <Quantity 0.00125979 m>,
 'x': <Quantity 0.71806864 m>,
 'y': <Quantity -0.32910527 m>}
[14]:
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]
[15]:
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)
[15]:
(-0.8291052663127615, 0.1708947336872385)
../_images/tutorials_calibrated_data_exploration_19_1.png
[16]:
source.metadata
[16]:
{'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:

[17]:
tels_in_event = set(
    event.dl1.tel.keys()
)  # use a set here, so we can intersect it later
tels_in_event
[17]:
{26, 27, 29, 47, 53, 59, 69, 73, 121, 122, 124, 148, 162, 166}
[18]:
cam_ids = set(sub.get_tel_ids_for_type("MST_MST_NectarCam"))
cam_ids
[18]:
{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,
 128,
 129,
 130}
[19]:
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:

[20]:
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

[21]:
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))
[21]:
Text(0.5, 1.0, 'Sum of 3x MST_MST_NectarCam')
../_images/tutorials_calibrated_data_exploration_28_1.png

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)

[22]:
from ctapipe.visualization import ArrayDisplay
[23]:
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()
<Figure size 432x288 with 0 Axes>
../_images/tutorials_calibrated_data_exploration_31_1.png
[ ]: