import numpy as np
import hipecta.core
import astropy.units as u
from ctapipe.core import Component
from traitlets import Int
'''
TO:
-1 Set a config by telescope_type
-Fill ctapipe event container r1, and dl0
'''
class TelescopeInfo:
def __init__(self, reco_temporary):
"""
Attributes
----------
reco_temporary: hipecta.core.PRecoTemporary
Parameters
----------
reco_temporary: hipecta.core.PRecoTemporary containing:
focalLength, matNeighbourQuadSum, nbSliceBeforPeak, matCalibratedSignal,
matNeighbourSlice, nbSliceWindow, matKeepSignalQuad, matSignal,
newMatCalibratedSignal, matNeighbourPixelSum, matSignalQuad,
newMatKeepSignalQuad
"""
self.reco_temporary = reco_temporary
[docs]class MCZfitsReco(Component):
"""
Process a raw telescope waveform (R1) to DL1 (Hillas parameters)
"""
wavelet_threshold = Int(default_value=3, help="wavelet cleaning threshold").tag(config=True)
hillas_threshold_signal_tel = Int(default_value=500, help="hillas threshold").tag(config=True)
def __init__(self, config=None, tool=None, **kwargs):
"""
handle the r0 to dl1 reconstruction. Calibration/Integration/Cleaning/Hillas.
Attributes
----------
telescope_info : dictionary
Python dict. key is telescope id. Value is an instance of Telscope_info
Parameters
----------
config : traitlets.loader.Config
Configuration specified by config file or cmdline arguments.
Used to set traitlet values.
Set to None if no configuration to pass.
tool : ctapipe.core.Tool or None
Tool executable that is calling this component.
Passes the correct logger to the component.
Set to None if no Tool to pass.
waveletThreshold : int
Wavelet cleaning threshold.
hillasThresholdSignalTel: int
The hillas threshold
kwargs
"""
super().__init__(config=config, parent=tool, **kwargs)
self.telescope_info = dict()
self.cut_config = hipecta.core.PConfigCut()
self.cut_config.waveletThreshold = self.wavelet_threshold
self.cut_config.hillasThresholdSignalTel = self.hillas_threshold_signal_tel
def _add_mc_telescope(self, event, number_samples, pixels_position, reference_shape):
#telescope_description: ctapipe.instrument.TelescopeDescription,
#camera: ctapipe.io.containers.MCCameraEventContainer):
"""
Add a telescope information and create gain, pedestal, ref_shape and pixel pos arrays
Parameters
----------
telescope_id: int
number_samples: int
pixels_position :
reference_shape:
#telescope_description: `ctapipe.instrument.telescope.TelescopeDescription`
#camera: ctapipe.io.containers.MCCameraEventContainer
Returns
-------
"""
gain_high = telescope.tabGainHi
gain_low = telescope.tabGainLo
pedestal_high = telescope.tabPedHi
pedestal_low = telescope.tabPedLo
reco_temporary = hipecta.core.createRecoTemporary(number_samples,
0,
pixels_position,
gain_high,
gain_low,
pedestal_high,
pedestal_low,
)
hipecta.core.updateRecoTemporaryWithRefShape(reco_temporary, reference_shape, 0.1)
self.telescope_info[telescope.telescopeId] = TelescopeInfo(reco_temporary)
def _process(self, telescope_id, waveform):
"""
Process a waveform and return its Hillas parameters
Parameters
----------
waveform: `numpy.ndarray` of shape (number_slice,
Returns
-------
Hillas parameters
"""
reco_temporary = self.telescope_info[telescope_id].reco_temporary
hillas_parameters, event_reco = hipecta.core.fullAnalysis(waveform, self.cut_config, reco_temporary)
return hillas_parameters, event_reco
[docs] def process(self, telescope, event, reference_shape):
"""
Process a ctapipe r0 event and return its Hillas parameters
Parameters
----------
telescope_id: int
event: ctapipe.io.containers.DataContainer
fill_container: bool
Returns
-------
Hillas parameters
"""
telescopeId = event.telescopeID
if telescopeId not in self.telescope_info:
pixels_position = telescope.tabPos.tabPixelPosXY
number_samples = len(event.matSliceHi)
self._add_mc_telescope(telescope, number_samples, pixels_position, reference_shape)
matslice = event.matSliceHi
return self._process(telescopeId, matslice)
def get_hillas_parameters_container(hillas_param):
"""
convert tuple of hillas parameter (from HiPeCTA to ctapipe HillasParametersContainer
Parameters
----------
hillas_param: tuples containaing hillas paramters
Returns
-------
ctapipe HillasParametersContainer
"""
return HillasParametersContainer(intensity=hillas_param[hipecta.core.getHillasImageAmplitude()],
x=hillas_param[hipecta.core.getHillasGx()] * u.m,
y=hillas_param[hipecta.core.getHillasGy()] * u.m,
width=np.sqrt(2) * hillas_param[hipecta.core.getHillasLength()] * u.m,
length=np.sqrt(2) * hillas_param[hipecta.core.getHillasWidth()] * u.m,
r=np.sqrt(hillas_param[hipecta.core.getHillasGx()] ** 2
+ hillas_param[hipecta.core.getHillasGy()] ** 2) * u.m,
phi=Angle(hillas_param[hipecta.core.getHillasPhi()] * u.rad),
psi=Angle(
(hillas_param[hipecta.core.getHillasDirection()] + np.pi / 2) * u.rad),
skewness=hillas_param[hipecta.core.getHillassSewness()],
kurtosis=hillas_param[hipecta.core.getHillasKurtosis()]
)