Monte Carlo study of a single SST-1M prototype for the Cherenkov Telescope Array

TL;DR

This paper validates a Monte Carlo model of the SST-1M prototype telescope for the Cherenkov Array, focusing on gamma/hadron separation and energy/direction reconstruction using machine learning.

Contribution

It presents the first validation of the Monte Carlo model of the SST-1M prototype and assesses its expected performance in real observing conditions.

Findings

Monte Carlo model accurately predicts telescope performance.

Machine learning improves gamma/hadron separation.

Expected energy and direction reconstruction performance is quantified.

Abstract

The SST-1M telescope was developed as a prototype of a Small-Size-Telescope for the Cherenkov Telescope Array observatory and it has been extensively tested in Krakow since 2017. In this contribution we present validation of the Monte Carlo model of the prototype and expected performance in Krakow conditions. We focus on gamma/hadron separation and mono reconstruction of energy and gamma photon arrival direction using Machine learning methods.

Tables

Table 1. Table 1: Expected event rates for the full telescope FoV and energy thresholds for Crab observation in Krakow.

	Energy	Crab	Proton
	threshold	event rate	event rate
	[TeV]	[mHz]	[Hz]
All	2.692	5.722	8.957
triggered
After	3.641	2.447	2.678
cleaning
Quality	5.009	1.237	1.558
cuts
g/h	4.966	1.174	0.605
separation

Full text

1 Introduction

SST-1M was developed as a prototype of a Small-Sized Telescope for the Cherenkov Telescope Array [1], designed for observations of the gamma-ray induced atmospheric showers for energies above 3 TeV. The SST-1M design is based on Davies-Cotton concept with a 4-m multi-segment mirror dish composed of 18 hexagonal facets [2]. The telescope is equipped with an innovative camera which features a fully digital readout and trigger system, called DigiCam [3] and adopts silicon photomultipliers (SiPMs) as light sensing technology. Photo-detection plane of the camera consists of 1296 SiPM pixels and the whole system has a large 9-degree diameter field of view (FoV).

The first fully operational SST-1M prototype is located in Krakow and has been extensively tested since 2017. In order to understand the data from the prototype, to estimate its performance in Krakow conditions and to perform a high level analysis as gamma/hadron separation or image reconstruction, a precise Monte Carlo (MC) simulation of the telescope is neccessary.

In the following sections a brief description of MC model validation, which is also mandatory for the MC validation process of CTA, is presented (the full description can be found in [4]) together with estimation of the prototype performance in mono regime and Krakow atmospheric conditions.

2 Monte-Carlo model of the SST-1M prototype

The MC model of the prototype was created based on the data taken at IFJ in Krakow and the laboratory measurements of individual components like mirrors, window or camera electronics. The gain and electronic noise of individual SiPMs were also measured in the lab, but to check validity of the MC model, raw ADC distributions were compared for a dark run with camera lid closed with sim$\_$telarray [5] simulations of pedestal events for a wide parameter space of gain, noise and dark count rate (DCR). The best matching distribution for a single camera pixel is shown in the left Fig. 1. It turned out that the parameters which need to be set in sim$\_$telarray to reproduce the data are slightly different from the laboratory measurements. While the difference in DCR can be explained by temperature differences between the lab and the site, the difference in amplitude and electronic noise ($15\%$ and $8\%$ respectively) is still a matter of investigation [4].

As a high level test of MC model validity, the measured rate scans were compared with simulations, which is shown in the left Fig. 2. The shape of simulated rate scan corresponds well with the data in Night Sky Background (NSB) dominated part of the dependency. In the left Fig. 2 there is also shown the proton trigger rate for three zenith angles multiplied by 1.5 as a correction for cosmic rays (CR) composition. Detailed analysis showed that the real trigger rate in that part of the rate scan tends to be about $2\times$ higher than in simulations. The right plot of Fig. 2, however, shows that the slope of a linear fit of $\log_{10}{(\mathrm{trigger\;\,rate})}=f(\mathrm{threshold})$ dependency for CR trigger dominated part of the rate scans is consistent with the slope for simulated protons.

The MC model of the prototype was used to build a small library of CORSIKA [6] and sim$\_$telarray simulations of the prototype for Krakow atmospheric conditions for 20 deg zenith angle. The trigger threshold was set on 260 ADC (45.6 p.e.) for the simulations to get about 300 Hz trigger rate for typical NSB level in Krakow (about 300 MHz222Compare with a typical NSB level at CTA-S site of about 40 MHz.), which corresponds with the adopted strategy of threshold settings during data taking.

For calibration, cleaning and Hillas parameters [7] extraction, the SST-1M pipeline digicampipe333https://github.com/cta-sst-1m/digicampipe, which is based on ctapipe444https://github.com/cta-observatory/ctapipe, was used. As an image cleaning method, standard tailcut cleaning is adopted in digicampipe, using two cuts optimized by minimizing the miss Hillas parameter on point source gamma simulations.

After processing data from scientific runs and proton simulations, another important high level test of MC validity can be done by comparing the Hillas parameter distributions. One night of ON source observation555OFF data wasn’t taken, but ON data are still dominated by CR as only quality cuts were applied in this case. was chosen for this test and the distribution of the Hillas parameter size is shown in the right Fig. 1. The protons were simulated with -2.0 spectral index and therefore the distribution of parameters had to be re-weighted with respect to the expected trigger rates for the real CR spectrum (see Sec. 3.4). After re-weighting, both distributions show very good match.

3 Mono reconstruction methods and the prototype performance

3.1 Arrival direction reconstruction

For the arrival direction reconstruction of the showers, several implementations of DISP method have been tested. The DISP parameter, which is the angular distance between predicted source position and the cleaned image center of gravity, depends mostly on the length and width Hillas parameters. First, we tested several formulas for DISP [8, 9, 10] and minimized parameters in the equations with respect to squared distance of reconstructed source position from the center of the FOV ($\theta^{2}$) using simulations of on-axis point source gammas. It turned out, that the resolution ($R_{68}$) we can reach using this approach varies from 0.45 deg to 0.60 deg, depending on the formula used. Another approach we tested consist in filling multidimensional look-up tables DISP = f[length,width, size, time$\_$gradient666A slope of photon arrival times in each pixel along main axis of Hillas elipse.], which leads to $R_{68}$ = 0.16 deg. Finally, we tested K-Neighbors (KN) and Random Forest (RF) regressors from scikit-learn Python library [11] for DISP determination and RF classifier to determine the correct sign of DISP. The simulation data set was split into train, validation and test dataset and optimal parameters for each method were grid-searched and cross validated. Both methods tested perform similarly according to coefficient of determination (KN: $R^{2}_{\mathrm{test}}=0.945$, RF: $R^{2}_{\mathrm{test}}=0.956$) and the image resolution after Random Forrest regression, which was finally adopted as a primary method of arrival direction reconstruction of the prototype data, is $R_{68}$ = 0.15 deg. The angular resolution as a function of reconstructed energy is shown in Fig. 3.

3.2 Energy reconstruction

The number of emitted Cherenkov photons is proportional to energy of primary gamma ray. Therefore, the energy of primary photon can be reconstructed from the number of photons (or photoelectrons) detected, which is given by the Hillas parameter size, and distance of shower core from the telescope, which is expressed by the impact parameter. In case of multi-telescope detection, impact parameter can be calculated from trivial geometry, which is not possible for mono reconstruction. In that case, however, the impact parameter for distant showers can be derived from the time$\_$gradient of a shower image.

For the energy reconstruction, we tested KN and RF regressors from scikit-learn library, trained on simulated point source gammas. Optimal parameters for each method were found by grid-searching and cross validated. The RF regressor performs slightly better ($R^{2}_{\mathrm{test}}=0.943$) than KN regressor ($R^{2}_{\mathrm{test}}=0.925$) and it was therefore used as a primary method for energy reconstruction in this study. The correlation between simulated and reconstructed energy is shown in right Fig. 4

The energy resolution and the bias were calculated from a gaussian fit of the relative error distribution in each true energy bin (see left Fig. 4). One can see that energy resolution which can be reached in our mono reconstruction is approximately $20\,\%$ above 10 TeV with $-3\%$ energy bias between 10 and 100 TeV.

3.3 Gamma-hadron separation

In IACT observations, the trigger rate from diffuse CR background is about $1000\times$ higher than the trigger rate from gamma ray photons and therefore a strong background suppression is necessary. In this analysis, we used RF classifier from scikit-learn, trained on diffuse protons and point source gammas after quality cuts777width / length ¿ $10^{-3}$, n$\_$islands = 1, n$\_$pixels ¿ 2, size ¿ 50 p.e., leakage2 ¡ 0.3. Optimal parameters were grid-searched using area under ROC curve (AUC) as a measure of classification performance. The optimal cut on hadronness was found by maximizing the F1-score [13]. In the left Fig. 5 one can see that the most important features for classification are Reduced Scaled Width and $\log(\mathsf{var}(E))$, which is the inter-tree variance of reconstructed energy, higher in general for protons than for gammas as a consequence of the fact that RF for energy reconstruction was trained on gamma events only [12]. The optimal performance of the separation gives AUC = 0.905, precision = 0.762 and recall = 0.908. Having false positive rate of $28\%$ the separation is weaker than we expected, but our preliminary tests with the use of more powerful boosted decision trees methods lead to similar performance, which suggests that such performance is close to mono reconstruction limits of the telescope.

3.4 Expected trigger rates and differential sensitivity

In order to estimate the sensitivity of the prototype for Krakow conditions with high NSB level, the effective areas and the expected event rates for a gamma point source with Crab spectrum and background diffuse protons with CR spectrum were calculated. In Fig. 6, effective areas and differential event rates are plotted for each analysis step: all triggered events, events that survived cleaning, events after cleaning and quality cuts and finally, events after gamma/hadron separation. The energy thresholds for gamma ray showers detection and the integrated event rates for the full telescope FoV are listed in Tab. 1. The sensitivity, calculated with respect to three criteria on the flux ($5\sigma$ detection $\&$ $N_{\gamma}>5\%$ of the residual background $\&$ $N_{\mathrm{excess}}>10$) with an additional cut on distance from the source of 0.15 deg, shows that even in Krakow a solid $5\sigma$ detection of Crab should be possible for 50h observation.

4 Acknowledgments

We gratefully acknowledge financial support from the agencies and organizations listed here: http://www.cta-observatory.org/consortium$\_$acknowledgments. The work is supported by the projects of Ministry of Education, Youth and Sports: MEYS LM2015046, LTT17006 and EU/MEYS CZ.02.1.01/0.0/0.0/16_013/0001403, Czech Republic, and by the grant Nr. DIR/WK/2017/12 from the Polish Ministry of Science and Higher Education. We greatly acknowledge financial support form the Swiss State Secretariat for Education Research and Innovation SERI.