- Email: [email protected]
A calibration of WFCTA prototype telescopes using N2 laser
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Contents lists available at ScienceDirect
Nuclear Inst. and Methods in Physics Research, A journal homepage: www.elsevier.com/locate/nima
A calibration of WFCTA prototype telescopes using 𝑁2 laser J.L. Liu a, *, F.R. Zhu b , H.Y. Jia b , Y.X. Bai c , Z. Cao c , M.J. Chen c , S.Z. Chen c , S.H. Feng c , B. Gao c , M.H. Gu c , H.H. He c , C. Hou c , J. Liu c , L.L. Ma c , G. Xiao c , M. Zha c , B.K. Zhang d , S.S. Zhang c , Y. Zhang c, *, J. Zhao c , B. Zhou c , X. Zuo c , for the LHAASO collaboration a b c d
Physics Department, Kunming University, Kunming, 650214, PR China School of Physical Science and Technology, Southwest Jiaotong University, Chengdu, 610031, PR China Key Laboratory of Particle Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences, Beijing 100049, PR China Physics Department, Normal College of Fuyang, Fuyang, 236037, PR China
a r t i c l e
i n f o
Keywords: Detector calibration Effective area Light flux Uncertainty
a b s t r a c t An end-to-end calibration method, which is different from the traditional way, has been employed in WFCTA prototype experiment based on 𝑁2 laser. As we known, the calculation of the expected and the observed light fluxes at the front of the aperture of the telescope are crucial in the such method. To calculate predicted light flux accurately, the propagation and scattering of photons in the atmosphere have been considered in detail. In order to get the detector effective areas, which will be used in the transformation of ADC signals to the observed light flux, a ray trace procedure, considering all details of detector response, has been performed for all incident direction of photons in field of view of the detector. By matching the observed and expected light flux at the entrance of the telescope, the overall calibration gains of two telescopes can be obtained. During the test run of the calibration, a overall calibration gain around 𝐺 = 0.38 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 and 𝐺 = 0.36 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 for WFCT01 and WFCT02 are given, respectively. Furthermore, in view of all uncertainty sources, systematic uncertainties for the calibration have also been estimated. They are around the level of (−7.3%, +7.8%). The calibration results achieved by end-to-end method and the traditional way are in consistent at the level of 10%, which is within the range of uncertainties. Since the end-to-end and traditional calibration methods have their own advantages, the performance of the combination of these two calibration methods has also been discussed in the end. © 2017 Elsevier B.V. All rights reserved.
1. Introduction The energy spectrum of cosmic ray (CR) contains important clues concerning the production and acceleration mechanisms in galactic and extragalactic particle sources. Over the past few decades, many ground-based experiments have been constructed focusing on its measurement in different energy regions above several TeV. However, because the dynamic range of individual experiments is limited and each experiment has its own energy scale and uncertainty, the observed results are not consistent among them. The precise measurement of the spectrum at a same observatory with sufficient overlap between different energy regions, which allows to obtain a complete and selfconsistent observation over the whole energy range from several TeV to several EeV, is necessary. The project of the Wide Field of View (FOV) Cherenkov/Fluorescence Telescope Array (WFCTA), composed
of 24 telescopes, which is a main component of Large High Attitude Air Shower Observatory (LHAASO) [1–3], has been designed to reach this goal. Using the atmosphere as the absorber and the scintillator of a large calorimeter, WFCTA prototypes [4] observe CR by measuring Cherenkov lights induced in air showers. The key issue in CR energy spectrum measurement in such an experiment is the determination of the total energy of reconstructed air showers which depend on an evaluation of the convertor from ADC counts recorded by electronics to a light flux at the aperture of the telescope. That is, it requires an absolute detector calibration. Various techniques to calibrate the response of such telescopes have been performed or proposed. One is using CRs to give an overall calibrated response of the telescope, including the effect of atmosphere transmission [5]. One can adjust a global calibration factor
* Corresponding authors.
E-mail addresses: [email protected] (J.L. Liu), [email protected] (Y. Zhang). https://doi.org/10.1016/j.nima.2017.09.041 Received 11 February 2017; Received in revised form 29 August 2017; Accepted 17 September 2017 Available online 6 October 2017 0168-9002/© 2017 Elsevier B.V. All rights reserved.
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
by comparing the MC-predicted CR counting rate and the measured counting rate. One main drawback of this technique is that the predicted rates depend on the CR flux and especially the composition. This would be involved large uncertainties. Another technique is comparing Cherenkov ring images in the focal plane produced by predicted and measured CR muons [6,7]. However, complications may arise from the difference in the wavelength spectrum of Cherenkov light from muons and from air showers. In addition, the predicted feature of Cherenkov rings by muons depending on the muon spectrum. Therefore one only can measure a spectrum weighted average response to the detector. The using of a fast pulsed calibrated light source would make the absolute calibration of these telescopes more straightforward and reliable. The traditional way is to calibrate all components of the detector piece by piece by using a calibrated light source, such as the mirror reflectivity, the transmission of the filter at the entrance, the quantum efficiency, photon–electron collecting efficiency, amplifier gain of photon multiplier tube (PMT), and the electronic gain. Details of all components can be found in [4]. The product of these factors provides the overall calibration factor. This is the way in current calibration and data analysis for WFCTA prototype experiment [8,9]. This method is effective on investigating the detector response in somewhat details, for instance the wavelength dependence of the quantum efficiency and the filter transmission. Difficulties for this method are (i) one cannot make sure whether or not all effects are considered, and (ii) every single piece introduces a certain of uncertainty and there is no well understanding way to estimate the combination effect of these errors to the overall calibration factor, due to the interaction between components of the system. Another approach that is so called end-toend method, which gives a overall convertor factor from ADC count to the light flux at the telescope entrance, has been applied here. One can estimate the uncertainty well in this method. The shortage is that it cannot provide a full understanding of the wavelength dependence of the detector gain. Auger experiment has applied this method to calibrate the telescope. A drum illuminator consists of a pulsed UV LED, emitting around 375 nm has been positioned at the entrance of the telescope under calibration. Absolute calibration constants were obtained from the known pulsed flux of photons emitted by the drum, and the corresponding ADC pulse integrals of the camera pixels [10]. In this paper, another end-to-end calibration method which uses a 𝑁2 laser mounted at several kilometers away from the telescope has been performed. In 2011, a 𝑁2 laser has been set up in the FOV of WFCTA prototype telescopes. Under the same condition as telescopes were operated in CR observation and some well chosen clear nights, the laser pulse were shot with several fixed directions respectively. After propagating in the atmosphere, telescopes recorded signals from the laser pulse along the track. Based on the calculation of the light flux at the front of the detector, one can directly calibrate the detector by analyzing the measured ADC counts. Uncertainties of this method are mainly from the laser energy probing, the light calculation related with aerosol and Rayleigh scattering models, as well as ADC counts. The remainder of this paper is arranged as follows. The light scattering and propagation model in the atmosphere is described in Section 2. Section 3 is the briefly description of the detector and the laser system. The status of data sets, the calculation of light flux, the detector response to the photons, as well as the estimation of the effective detector area, together with the calibration are presented in Section 4. Section 5 is the discuss about the calibration uncertainties and the performance of the combination of two calibration methods. Section 6 is the summary, as well as the prospect of detector calibration in LHAASO-WFCTA experiment in the future.
Fig. 1. Ozone level distribution as a function of altitude.
photons are emitted in a direction close to the incidence direction of the primary particle. Attenuation of the atmosphere for these lights is due to both absorptive process and scattering process on their way to the detector. Absorptive process is dominant at ultraviolet wavelengths where ozone component of the atmosphere begins to absorb the lights, and above 800 nm where vapor and CO2 absorptions begin [11,12]. In the wavelength range of 280 nm to 700 nm where the ground telescope is sensitive, scattering process initiated by molecules and aerosols in the atmosphere are dominant. That is, the contribution of absorption to the light attenuation can be negligible except for ozone absorption. Ozone absorption effect can be expressed by its transmission coefficient 𝑇𝑂 [13]: 𝑇𝑂 = exp(−𝐴𝑂3 𝐶𝑂3 (𝑙2 − 𝑙1 ))
(1)
where 𝐴𝑂3 is the ozone absorptive coefficient depending on the wavelength, and 𝐴𝑂3 = 0.117∕(cm atm) at 𝜆 = 337 nm. 𝐶𝑂3 represents the ozone level varying with the altitude, which is shown in Fig. 1 [13]. The molecule and aerosol scattering processes that contribute to the overall attenuation can be treated separately expressed by transmission coefficients 𝑇𝑅 and 𝑇𝐴 , respectively [14]. The molecule scattering is the most straightforward in two processes and can be described by Rayleigh Scattering [15]. According to this theory, the number of photons scattered out of the beam per unit length is 𝑁𝜌 400 4 𝑑𝑁 =− ( ) 𝑑𝑙 𝑋𝑅 𝜆
(2)
where 𝜆 is the wavelength of the photon and in an unit of nm. 𝑋𝑅 is the mean free path for scattering and is 2970 g/cm2 at 𝜆 = 400 nm [12]. And 𝜌 = 𝜌0 exp(−ℎ∕𝐻)
(3) 0.00129 g/cm3 , which is the
is the density profile of the atmosphere. 𝜌0 = atmospheric density at 0 ◦ C at sea level for an isothermal atmosphere. 𝐻 is the scale height of the atmosphere and is usually assumed to be 7.5 km. This parameter varies with the pressure and the temperature of the atmosphere. Thus it will be considered as an uncertainty source in light propagation. By solving function (2), one can find the transmission coefficient for lights passing from a point at slant depth 𝑋1 to a point at slant depth 𝑋2 can be expressed as 𝑇𝑅 = exp[−
2. Light propagation in the atmosphere
𝑋2 − 𝑋1 400 4 ⋅( ) ]. 𝑋𝑅 𝜆
(4)
Aerosol scattering is more complicated and can be described by Mie scattering theory [16,17]. However it relies on the simplifying assumption of spherical scatters, a condition which is not always fulfilled in the field. Moreover, it also depends on the particle composition, which changes rapidly relying on the wind and the whether condition. Based
When CRs enter the atmosphere, they induce extensive air showers (EAS) of secondary particles. While charged particles in an EAS cross the atmosphere with a speed higher than the light’s in the air, Cherenkov 279
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 2. Phase functions of aerosol scattering and Rayleigh scattering.
on Mie scattering theory, when photons traveling from a point at height ℎ1 and length 𝑙1 to a point at height ℎ2 and length 𝑙2 , the average aerosol transmission coefficient can be described by 𝑙2 − 𝑙1 ) (5) Λ Λ is the aerosol scattering length. While the height is lower than the mixing layer height (about 0.1 km), Λ is almost a constant and is assumed to be 25 km. Above the mixing layer, Λ increases with the height exponentially and can be read as
𝑇𝐴 = exp(−
Λ=
𝐿𝐴 ℎ2 − ℎ1 𝐻𝐴 𝑒ℎ2 ∕𝐻𝐴 − 𝑒ℎ1 ∕𝐻𝐴
(6)
Fig. 3. Typical pulse signals recorded by the PMT electronics of WFCTA prototype in a laser event. The 𝑋 axis is the FADC channel and each channel is 20 ns. The 𝑌 axis is the FADC count in each channel.
where 𝐿𝐴 is the aerosol attenuation length describing the light attenuation due to the aerosol at the ground level. 𝐻𝐴 is the aerosol scale height (about 1.2 km) accounting for its dependence on the height [17,18]. When lights arrive at a position in the FOV of the telescope, the number of photons scattered into an solid angle from the laser beam both by molecules and aerosols depend on the angular distribution of scattered lights, which is referred as phase function. The phase function for Rayleigh scattering is well understood and can be estimated analytically as function (7). With its key feature of the symmetry in the forward and backward directions, it is proportional to (1+cos2 𝜃). A small factor 𝛾 = 0.0155 is used to correct the effect of anisotropy of the 𝑂2 and 𝑁2 molecules [19]: 3 𝑃𝑅 (𝜃) = [(1 + 3𝛾) + (1 − 𝛾) cos2 𝜃]. 16𝜋(1 + 2𝛾)
Table 1 Summary of reflector parameters in WFCTA prototypes. Mirror facet length Mirror facet area Curvature radius Reflective material Reflectivity (>240 nm)
298.4 mm 0.23 𝑚2 4740 mm 𝐴𝑙 + 𝑀𝑔𝐹2 >83%
Table 2 Summary of PMT parameters in WFCTA prototypes. Manufacturer Model Photocathode Voltage Typical voltage Wavelength Gain Nonlinearity Rise time FWHM
(7)
For the phase function of aerosol scattering 𝑃𝐴 (𝜃), it is quite complex and model dependent. Two models referred as Elbert’s model and Longtin’s wind-dependent-desert model [20] are available to describe its gross feature on the light scattering probability distribution, which is sufficient for the purpose of Cherenkov detection. Elbert’s model is used as default in the study. In general, the angular distribution of lights scattered by the aerosol is very strong peaked in the forward direction and reaches a minimum near 90 ◦ . All these different phase functions are displayed in Fig. 2 and all of them are normalized over 4𝜋 solid angle.
Photonis XP3062/FL 40 mm 900 V ∼ 1300 V 1100 V 280 nm ∼ 700 nm 2 ×105 < 2% 3 ns 4.5 ns
an area of 5.0 m2 . For the detailed parameter information, see Table 1. The camera is constituted by 256 hexagonal PMTs. Detailed information of PMTs are shown in Table 2. Each PMT has a FOV of 1◦ × 1◦ and all PMTs are arranged as a 16 × 16 array, forming a total FOV of 16◦ × 14◦ . Details of the telescope can be found in reference [4,8]. When a laser pulse passes the FOV of the detector, scattering lights will trigger tubes watching the direction of the laser track and a FADC electronics system will record pulse signals every 20 ns. Fig. 3 displays such pulse signals recorded by PMTs in a laser event. A 𝑁2 laser has been set up at the west by north of ARGO-YBJ experiment, 2520 m away from the center of the experimental hall. Its parameters are summarized in Table 3. The laser is mounted on a rotating platform, enabling the laser beam has an angular resolution
3. WFCTA prototypes and the laser set Two telescopes have been set up in WFCTA prototype experiment. They are located at the west by north (denoted as WFCT01) and the south by east (denoted as WFCT02) of the ARGO-YBJ experiment [21], about 99 m and 79 m away from the center of the ARGO hall respectively. Their pointing directions are set as 32◦ and 31◦ in zenith direction, 170◦ and 171◦ in azimuth direction respectively. Each telescope mainly consists of a UV-light collecting mirror and an image camera. The reflected mirror is constituted by 20 hexagonal spherical facets, forming 280
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 6. The geometry of lights propagation and scattering in the calibration. The shaded area represents the FOV of the telescope.
Fig. 4. The distribution of laser energies during a run.
the calculation of the angular width of the bin along the track and the correlation of the effective area of the detector. In this section, firstly, data sets in the analysis will be introduced. secondly, the calculation of the light flux at the entrance of the telescope will be described. Thirdly, the binning and the estimation of detector effective area will be given. Finally, based on above studies, the calibration procedure and the result will be presented in detail.
Table 3 𝑁2 laser parameters. Manufacturer Model Wavelength Spectrum width Energy Pulse duration Stability Beam size Divergence
Stanford research systems NL100 337.1 nm 0.1 nm 170 μJ @ 10 Hz 3.5 ns 3% 3 × 7 mm 3 × 8 mrad
4.1. Data sets To avoid the influence of the aerosol to the calibration, the data that observed at a very clear night is chosen. Furthermore, as indicated in Fig. 2, values of phase functions of aerosol scattering both in two models are relatively smaller while the scattering angle is around 90◦ . Therefore the zenith of the laser pulse was set to around 60◦ , which makes the scattering angle being around 90◦ . Thirdly, to avoid the influence of marginal effect to the calculation of the effective area, the data with image on the center area of FOV are chosen. After data selection, three data sets, that is, 410 events with (𝜃, 𝜑) = (60◦ , −8◦ ), 732 events with (𝜃, 𝜑) = (50◦ , −8◦ ) and 684 events with (𝜃, 𝜑) = (60◦ , −10◦ ) (where 𝜃 is the zenith and 𝜑 is the azimuth), are available in the calibration. The former two are for WFCT01 and the latter is for WFCT02.
of 0.1◦ and a spatial resolution of 0.1 mm. To minimize the influence of the aerosol, the calibration has been carried out at very clear nights. The laser shots 337.1 nm photons with a pulse duration of 3.5 ns, elevations of 30◦ ∼ 60◦ and azimuths of −22◦ ∼ 4◦ . To match the laser pulse energy with the CR energy and to maintain a precise beam direction, a collimator has been mounted at the entrance of the laser. The absolute energy of the laser beam is measured by an energy probe, which has a calibration uncertainty of less than 5%. Energies of laser beams during a run shown in Fig. 4 indicate that the average energy of laser pulses is 90.8 μJ with a small fluctuation of 1.2%.
4.2. Light flux calculation 4. Absolute calibration of WFCTA prototypes The calculation of the light flux at the aperture of the telescope is crucial for end-to-end calibration. The calculation of the predicted light flux will be described in detail in this section. For given wavelength 𝜆 and energy 𝐸0 , the total number of photons in a laser pulse is 𝑁0 = 𝐸0 ∕(ℎ𝑐∕𝜆). ℎ is Planck constant and 𝑐 is light velocity here. For this study with 𝐸0 = 90.8 μJ and 𝜆 = 337.1 nm, photons emitted into the air in a pulse is 𝑁0 = 1.54×1014 . As shown in Fig. 6, these
As mentioned in Section 1, the matching of predicted light flux and observed light flux at the aperture of the detector is the key to the end-toend calibration method. As displayed in Fig. 5, the image of an observed laser event is the assemble of trigged PMTs. To transform recorded ADC counts for each PMT in the observation to the light flux in front of the detector, one not only needs the calibration factor, but also requires
Fig. 5. The typical laser event observed by WFCTA prototypes. The left panel is for WFCT01 and the right panel is for WFCT02. Green hexagons are PMTs on the image plane. Red solid circles represent signals collected by triggered PMTs. The diameter of circles is proportional to the logarithm of the amplitude of signals.
281
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Table 4 The flux of photons in front of the telescope, which are scattered from a 90.8 μJ 𝑁2 laser pulse at 1240 m above the ground. The unit of the flux is 𝑝ℎ𝑜𝑡𝑜𝑛𝑠∕𝑚2 ∕𝑑𝑒𝑔𝑟𝑒𝑒. The aerosol scale height 𝐻𝐴 is assumed to be 1.2 km for all combinations of aerosol scattering parameters here. Aerosol
𝑃𝐴 (𝜃)
Elbert Longtin Rayleigh
In order to include enough tubes in a bin so as to avoid a large fluctuation for the effective area due to the marginal tubes, a angular size of 1.5◦ for a bin is chosen in the binning. Once all bins are defined, all none-noised tubes are sorted into bins. If a tube’s center in SDP lies within 0.25◦ of a bin center, it is considered wholly within the bin (i.e. the fraction 𝑓 of the contribution for the tube to this bin is 1). If a tube’s center lies greater than 0.25◦ and less than 1.25◦ from a bin center, it is considered partly in the bin and a fraction 𝑓 is calculated according to formula (11) to determine how much of the tube is contained into the bin.
𝐿𝐴 12 km (standard)
25 km (average)
1000 km (Rayleigh)
345 226 1040
182 119 1145
5 3 1249
1 𝑓 =1− [𝑎 cos(1.5◦ −2 △ 𝜃 ′ )−(1.5◦ −2 △ 𝜃 ′ ) sin(𝑎 cos(1.5◦ −2 △ 𝜃 ′ ))] (11) 𝜋 lights will undergo a path from the laser (𝐿) to a scatter (𝑆) in the FOV of the detector at first. After scattering to the direction of the detector (𝐷), they propagate to the detector through the atmosphere. Propagating along the path 𝐿𝑆 and 𝑆𝐷, lights suffer the same scattering as at position 𝑆. While arrived at position 𝑆, the number of photons decrease to 𝑁1 for attenuations of Rayleigh scattering 𝑇𝑅1 , aerosol scattering 𝑇𝐴1 and ozone absorption 𝑇𝑂1 .
where △𝜃 ′ is the angular distance between the tube center and the bin center. The deduction of this function is presented in Appendix. If a tube’s center lies greater than 1.25◦ from the bin center, it is not included in the bin. Once all tubes are sorted into bins, the ADC count for each bin is calculated following 𝜙0 =
where the calculation of 𝑇𝑅1 , 𝑇𝐴1 and 𝑇𝑂1 follows Eqs. (4), (5) and (1), respectively. At location 𝑆 in the FOV of the detector, the number of photons scattered to a solid angle corresponding to the mirror size from the laser beam depends on the angular distribution of scattered lights, which is referred as phase function. The number of photons scattered from an unit angular track length and reached to an unit area of the mirror can be read as 𝑁2 = 𝑁1 [(1 − 𝑇𝑅 (𝑆))𝑃𝑅 (𝜃) + (1 − 𝑇𝐴 (𝑆))𝑃𝐴 (𝜃)]
1 1 2 Δ𝜃 𝑙𝑆𝐷
(12)
4.4. Effective area of the detector As mentioned at the beginning of this section, in order to transform the ADC signal to the light flux at the front of the detector, one needs to know the effective light collecting area of the detector. The mirror area of 5.0 m2 mentioned in Section 3 is only the geometric area. The true effective area not only depends on the geometry effect associated with the direction of incident rays, but also depends on the location and the size of the spot on the camera. That is because PMTs have different sensitivity varying with the location of the photocathode and there has insensitive area between hexagonal PMTs. Therefore the ratio of sensitive area to insensitive area changes as the spot moves around on the image plane. The effective area has to be determined using simulation for every position and direction of the light source as the source moves over the FOV of the detector. As described in Section 4.3, when a bin of the track has been defined, the vector of the bin center has been calculated simultaneously. Five thousands trial incident rays (𝑁𝑡𝑟𝑖𝑎𝑙 ) hit on the mirror from directions stepping down the track in a range of 3.5◦ centered at the bin center vector. A simulation procedure named as ray tracing has been employed for each incident ray to calculate the effective area for each bin. In the ray tracing procedure, it has considered the shape and the volume of the system container and the PMT cluster, and the detailed structure of the reflector and the camera, as well as shadows of the support structure of the filter in front of the telescope entrance and the cluster in the ray tracing procedure. By this ray tracing, one can calculate the effective area 𝑁 following 𝐴𝑒𝑓 𝑓 = 𝐴0 𝑁𝑠𝑢𝑐𝑐 of the bin, where 𝐴0 is the mirror geometric 𝑡𝑟𝑖𝑎𝑙 area on which the trial rays uniformly hit and 𝑁𝑠𝑢𝑐𝑐 is the number of rays which successfully reach to the effective area of PMTs. When a tube collects a trial ray, we check whether the tube is triggered during the measurement. If the tube is not triggered, the ray will not be counted into 𝑁𝑠𝑢𝑐𝑐 . As described in Section 4.1, events in three data sets have different geometry. Thus the geometry of photons received by the telescope are also different. Taking the angle (denoted as 𝑑𝑡ℎ𝑒) between the incident direction of photons and the primary optical axis of the telescope as an example, average values for three data sets are 2.5◦ , 2.8◦ and 4.5◦ , respectively. 𝐴𝑒𝑓 𝑓 distributions in Fig. 7 show that average values of 𝐴𝑒𝑓 𝑓 for three data sets are 2.29 m2 , 2.26 m2 and 2.11 m2 respectively. It is in agreement with the result in bottom-right panel of Fig. 7. That is 𝐴𝑒𝑓 𝑓 decreases with the increase of 𝑑𝑡ℎ𝑒. To study the uncertainty of the calculation of 𝐴𝑒𝑓 𝑓 , measured errors of all parameters associated with the detector in ray trace simulation have been estimated. They are curvature radius of mirror, 𝑅𝑚𝑖𝑟 , diameter
(9)
where Δ𝜃 is the angular bin size around the position 𝑆, 𝑙𝑆𝐷 is the path length from 𝑆 to 𝐷. Photons scattered to the direction of the detector undergo a same attenuation as propagating along the path 𝐿𝑆. Thus the photons (i.e. the light flux) decrease to 𝑁3 while at location 𝐷. 𝑁3 = 𝜙exp = 𝑁2 𝑇𝑅2 𝑇𝐴2 𝑇𝑂2
𝑓 𝑁𝐴𝐷𝐶 (𝑖)
where 𝑁𝐴𝐷𝐶 (𝑖) is the total ADC count recorded by the 𝑖th tube.
(8)
𝑁1 = 𝑁0 𝑇𝑅1 𝑇𝐴1 𝑇𝑂1
∑
(10)
where the calculation of 𝑇𝑅2 , 𝑇𝐴2 and 𝑇𝑂2 are same as 𝑇𝑅1 , 𝑇𝐴1 and 𝑇𝑂1 except for replacing the integral path 𝐿𝑆 with 𝑆𝐷. To understand the uncertainty associated with different atmospheric parameters, such as phase function models, scattering length, attenuation length, aerosol distribution etc., the flux at the front of the detector has been calculated under different atmospheric conditions. As an example, Table 4 gives the results of light flux scattered from a given position 𝑆 on the laser beam which is 1240 m above the ground. From the results listed in Table 4, the contribution of aerosol scattering can be as small as 0.24%–0.40% according to different aerosol scattering phase function as long as the atmosphere is clear. At other hand, even under the ‘‘average’’ atmosphere condition, that is 𝐿𝐴 = 25 km and 𝐻𝐴 = 1.2 km [12,17,18], the uncertainty of light flux associated with phase function is about 4.7%, although the aerosol scattering contribution can be as large as 13.7%. 4.3. Angular binning for the track There are several major steps in binning of the laser track. (i) Remove noise tubes according to triggered PMTs’ location on the image plane and their trigger times. (ii) Determine the Laser-Detector-Plane (LDP) based on the direction of the laser pulse and locations of the laser and the detector. (iii) Project all none-noised tubes into the LDP. (iv) Find the topmost tube and mark the first bin center. (v) Step down along the laser track in LDP with 1.5𝑜 interval and calculate the vector which is from detector pointing to each bin center. (vi) Sort tubes into bins, possibly splitting a tube among different bins by assigning a weighting factor for the tube. 282
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 7. The distribution of effective area for WFCT01 (top) and WFCT02 (bottom-left), and as a function of 𝑑𝑡ℎ𝑒 for both two telescopes (bottom-right). Table 5 Nominal values and estimated uncertainties for all parameters used in ray trace simulation. The unit for all parameters is mm. Parameter
𝑅𝑚𝑖𝑟
𝐷𝑚𝑖𝑟
𝐷𝑚𝑐
𝑊𝑐
𝐻𝑐
𝑇𝑐
𝑅𝑠𝑝𝑜𝑡
Nominal value Uncertainty
4740 ±10
596.8 ±3
2305 ±3
711 ±5
635 ±5
355 ±5
15 ±3
337.1 nm. 𝑄𝑒 , 𝜀 and 𝐺𝑃 𝑀𝑇 represent the quantum efficiency, photon electron collecting efficiency and gain of PMT respectively, and 𝐺𝑒 is the electronics gain. Measuring these four effects item by item is difficult. In the calibration, a UV-LED at 375 nm mounted at the center of the mirror is applied to obtain the overall response for each PMT [4]. For the calibration is done at 375 nm, a further correction is applied according to the wavelength dependence of the quantum efficiency for the laser calibration at 337.1 nm. The other way of the calibration is to measure the overall factors using a laser beam, which is so called end-to-end method. In such method, for the given energy and geometries of the laser and the detector, one can calculate the expected light flux 𝜙exp at the front of the detector following the method described in Section 4.2. Then one can determine the factor 𝐺 by matching the predicted light flux in Eq. (10) and the observed light flux in Eq. (13). That is
of a mirror segment, 𝐷𝑚𝑖𝑟 , distance from mirror center to cluster center, 𝐷𝑚𝑐 , dimensions of cluster box, 𝑊𝑐 (width), 𝐻𝑐 (height), 𝑇𝑐 (thickness), average spot size at the cluster, 𝑅𝑠𝑝𝑜𝑡 . Details of nominal value and estimated error for these parameters are listed in Table 5. The overall uncertainty to 𝐴𝑒𝑓 𝑓 caused by these parameters is ±2%. 4.5. Detector calibration Once the total ADC signal 𝜙0 and the effective area 𝐴𝑒𝑓 𝑓 for each bin along the track are calculated, the observed light flux at the front of the detector can be expressed as 𝜙𝑜𝑏𝑠 =
𝜙0 𝐴𝑒𝑓 𝑓 △ 𝜃𝐺
.
𝐺=
.
(15)
In a laser event, a gain factor can be obtained from each bin along the track. To avoid the influence of the fluctuation, only a gain factor is given in a event by averaging calibration factors obtained from all bins. All calibration results are presented in Fig. 8. The analysis of both two data sets of WFCT01 (top panels) give a consistent calibration constant by end to end calibration method, that is 𝐺 = 0.38 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 in average. For WFCT02 (the bottom panel), the calibration result is 𝐺 = 0.36 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 in average. Calibration results obtained from piece by piece method are also shown to give a comparison. Average gain factors from piece by piece method are about 10% larger than ones from end-to-end method, which is within the level of the calibration uncertainties.
(13)
The goal of the calibration is to determine the value of 𝐺. As described in Section 1, one way is of the calibration to directly measure all factors in Eq. (14) piece by piece, which is so called traditional method. 𝐺𝑝 = 𝑅𝑚𝑖𝑟𝑟𝑜𝑟 𝑇𝑓 𝑖𝑙𝑡𝑒𝑟 𝑄𝑒 𝜀𝐺𝑃 𝑀𝑇 𝐺𝑒
𝜙0 𝜙exp 𝐴𝑒𝑓 𝑓 △ 𝜃
(14)
where 𝑅𝑚𝑖𝑟𝑟𝑜𝑟 = 0.83 is the mirror reflectivity, 𝑇𝑓 𝑖𝑙𝑡𝑒𝑟 = 0.865 is the transmission of the filter at the entrance, both two values are at 𝜆 = 283
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 8. The distribution of calibration gain by end-to-end (solid curves) and piece by piece (dashed curves) methods. Top panels are for WFCT01 and the bottom panel is for WFCT02.
5. Discuss
5.2. Uncertainties on end-to-end calibration
5.1. Calibration method
As stated at the end of Section 1, sources of systematic uncertainty in end-to-end calibration are (i) laser energy probing, ±5%, (ii) aerosol scattering model, see Table 4, (iii) the calculation of Rayleigh scattering due to the uncertainty of the temperature and pressure of the atmosphere, ±2%, (iv) the calculation of 𝐴𝑒𝑓 𝑓 , ±2%, see Table 5, and (v) ADC count error due to the overall effect of the fluctuation and the systematic error of electronics initiated by the variation of the temperature of the external environment, ±2%. (i) and (ii) are the most significant effects. To study the influence of aerosol scattering model, we have changed the horizontal attenuation and phase function displayed in Table 4, respectively. The error contribution of this source is (−4.1%, +4.9%). The overall uncertainty to predicted light flux introduced both by (i), (ii) and (iii) is (−6.8%, +7.3%), see Fig. 9. (iv) and (v), together with the uncertainty of 𝐺𝑝 (±7%) [4], contribute to the uncertainty of observed light flux. The overall contribution is (−7.5%, +7.5%), also see Fig. 9. For the calibration gain factor 𝐺 under end-to-end method, the overall uncertainty caused both by above sources is (−7.3%, +7.8%). And for 𝐺𝑃 , the overall uncertainty is ±7%. From Fig. 8, one can see that the difference of average gain factors obtained by two calibration method is around 10%. So calibration gains achieved by two methods are consistent within the level of uncertainties.
As mentioned in Section 1, in the calibration, end-to-end method and the traditional way have their own advantages and disadvantages. The combination of two methods provides the best way to do the calibration. The gain 𝐺𝑝 obtained by the traditional method is applied in Eq. (13) to get the observed light flux at the front of the detector. The expected light flux in Eq. (10) is used to find an adjustment 𝑓𝑎𝑑𝑗 for the current gain 𝐺𝑝 of the detector. The predicted and observed light fluxes along the track at the front of the detector obtained by above method are presented in Fig. 9. The topmost and the bottommost bins are affected seriously by the marginal effect, thus these two bins are removed. They indicate that the observed fluxes are about 10% lower than predicted. To avoid the influence of the ∑ fluctuation for the observed light flux, an integral flux 𝜙 is taking for each laser event. Results of integral fluxes for three data sets are shown in Fig. 10. The adjustment of the gain in current data analysis at 337.1 nm can be obtained by ∑ 𝜙 𝑓𝑎𝑑𝑗 = ∑ 𝑜𝑏𝑠 . (16) 𝜙exp
6. Summary
In this study, only one night data is available. Thus we just present the calibration method and procedure, and cannot give the conclusion of 𝑓𝑎𝑑𝑗 measurement in this paper.
Through this study, first, a complete analysis procedure and an end-to-end calibration method have been established for the absolute 284
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 9. The distribution of predicted (open circles) and observed (solid circles) light fluxes as a function of the light source height for WFCT01 (top) and WFCT02 (bottom) telescopes. The shades around the points represent the uncertainties.
Appendix. Derivation of function (11)
calibration of WFCTA telescopes, including the detailed calculation of light scattering and propagation in the atmosphere, the shower track binning, and the ray trace which is used to calculate the effective area of detector upon different incident directions. Second, by using this method, in the test run, a overall calibration gain around 𝐺 = 0.38 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 and 𝐺 = 0.36 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 for WFCT01 and WFCT02 are given. Third, both the advantage and disadvantage of traditional and end-to-end calibration methods are discussed. Based on the discussion, the performance of the combination of two methods in the calibration is investigated. Finally, the uncertainty of the calibration considering all error contributions is estimated. They are in the level of (−7.3%, +7.8%). The calibration results obtained by end-to-end method and the traditional way are in consistent at the order of 10%, which is within the level of uncertainties. In the future, we will monitor the aerosol condition hourly, as well as monitor the temperature, the pressure and humidity of the atmosphere in the calibration of LHAASOWFCTA telescopes. This will help to decrease the calibration uncertainty caused by atmosphere model significantly.
Eq. (11) in the article is used to calculate the fraction of a PMT contributing to the bin according their central vectors. Two cases of the binning are shown in Fig. 11. That is, the center of a PMT is inner the FOV of the bin (left) and is outer the FOV of the bin (right). Since the whole derivations of Eq. (11) in such two cases are the same, for concise, only the left case is provided as follows: For convenience, we denote by 𝑅 ∶ 𝑡ℎ𝑒 𝑟𝑎𝑑𝑖𝑢𝑠 𝑜𝑓 𝑡ℎ𝑒 𝑐𝑖𝑟𝑐𝑙𝑒 𝐴𝐵𝐶𝐷. 𝑠 ∶ 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑡ℎ𝑒 𝑐𝑖𝑟𝑐𝑙𝑒 𝐴𝐵𝐶𝐷. 𝑠0 ∶ 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑡ℎ𝑒 𝑠𝑒𝑔𝑚𝑒𝑛𝑡 𝐴𝐷𝐶. 𝑠1 ∶ 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑡ℎ𝑒 𝑠𝑒𝑔𝑚𝑒𝑛𝑡 𝐴𝐵𝐶. 𝑠2 ∶ 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑡ℎ𝑒 𝑠𝑒𝑐𝑡𝑜𝑟 𝐴𝑀𝐶𝐵. 𝑠3 ∶ 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑡ℎ𝑒 𝑡𝑟𝑖𝑎𝑛𝑔𝑙𝑒 𝐴𝑀𝐶. 𝑓 ∶ 𝑡ℎ𝑒 𝑝𝑟𝑜𝑝𝑜𝑟𝑡𝑖𝑜𝑛 𝑜𝑓 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑠𝑒𝑔𝑚𝑒𝑛𝑡 𝐴𝐷𝐶 𝑖𝑛 𝑡ℎ𝑒 𝑐𝑖𝑟𝑐𝑙𝑒 𝐴𝐵𝐶𝐷. 𝑓1 ∶ 𝑡ℎ𝑒 𝑝𝑟𝑜𝑝𝑜𝑟𝑡𝑖𝑜𝑛 𝑜𝑓 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑠𝑒𝑔𝑚𝑒𝑛𝑡 𝐴𝐵𝐶 𝑖𝑛 𝑡ℎ𝑒 𝑐𝑖𝑟𝑐𝑙𝑒 𝐴𝐵𝐶𝐷. △𝜃 ′ ∶ 𝑡ℎ𝑒 𝑎𝑛𝑔𝑙𝑒 𝑏𝑒𝑡𝑤𝑒𝑒𝑛 𝑡ℎ𝑒 𝑡𝑢𝑏𝑒 𝑐𝑒𝑛𝑡𝑒𝑟 𝑎𝑛𝑑 𝑡ℎ𝑒 𝑏𝑖𝑛 𝑐𝑒𝑛𝑡𝑒𝑟. Noting that 𝑠1 , 𝑠 in which 𝑠1 = 𝑠2 − 𝑠3 . By the area formula of a sector, we obtain that
Acknowledgments This work is supported by NSFC, Yunnan provincial Science and Technology Department and Kunming University through grants of 11563004, 2014FB154 and YJL17004, respectively. It is also funded by Program for Innovative Research Team (in Science and Technology) in University of Yunnan Province. We acknowledge the essential support of W.Y. Chen, G. Yang, X.F. Yuan, C.Y. Zhao. We are grateful to Prof. Xinhua Bai for his valuable comments and professional revision suggestions (SDSM&T). We are also thankful to Dr. Chang-Bo Yang for the discussions of mathematical problems with him.
𝑓 = 1 − 𝑓1 = 1 −
(17)
1 2 𝜃𝑅 , 2
(18)
𝑠2 =
where 𝜃 = 2𝑎 cos( 𝑀𝐸 ), 𝑀𝐸 = 0.75 − △𝜃 ′ and 𝑅 = 0.5. 𝑅 By simple computation, we have
285
𝜃 = 2𝑎 cos(1.5 − 2 △ 𝜃 ′ ),
(19)
𝑠2 = 𝑎 cos(1.5 − 2 △ 𝜃 ′ )𝑅2 ,
(20)
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 10. The distribution of integral predicted (dashed curves) and observed (solid curves) light fluxes for WFCT01 (top) and WFCT02 (bottom) telescopes.
References [1] Z. Cao, A future project at Tibet: The large high altitude air shower observatory (LHAASO), Chin. Phys. C 34 (02) (2010) 249–252. [2] M. Zha, the ARGO-YBJ collabaration, the LHAASO collabaration, Status of the large high altitude air shower observatory project, Nucl. Instrum. and Meth. in Phys. Res. A 692 (2012) 77–82. [3] Z. Cao, the LHAASO Collaboration, Status of LHAASO updates from ARGO-YBJ, Nucl. Instrum. and Meth. in Phys. Res. A 742 (2014) 95–98. [4] S.S. Zhang, Y.X. Bai, Z. Cao, et al., Properties and performance of two wide field of view Cherenkov/fluorescence telescope array prototypes, Nucl. Instrum. and Meth. in Phys. Res. A 629 (2011) 57–65. [5] A. Konopelkoa, et al., the HEGRA Collaboration, Detection of gamma rays above 1 TeV from the Crab Nebula by the second HEGRA imaging atmospheric Cherenkov telescope at La Palma, Astropart. Phys. 4 (3) (1996) 199–215. [6] T.B. Humenskya, for the VERITAS Collaboration, Calibration of VERITAS Telescope 1 via Muons, 2005. arXiv:astro-ph/0507449v1. [7] G. Vacanti, P. Fleury, Y. Jiang, et al., Muon ring images with an atmospheric Cherenkov telescope, Astropart. Phys. 2 (1) (1994) 1–11. [8] S.S. Zhang, et al., the ARGO Collaboration, LHAASO Collaboration, Energy spectrum of cosmic protons and helium nuclei by a hybrid measurement at 4300 m a.s.l., Chin. Phys. C 38 (4) (2014) 045001. [9] S.S. Zhang, et al., the ARGO Collaboration, LHAASO Collaboration, The knee of the Cosmic Hydrogen and Helium Spectrum below 1 PeV measured by ARGO-YBJ and a Cherenkov telescope of LHAASO, Phys. Rev. D 92 (092005) (2015). [10] A.A. Watson, the Auger collabaration, Properties and performance of the prototype instrument for the Pierre Auger Observatory, Nucl. Instrum. and Meth. in Phys. Res. a 523 (2004) 50–95. [11] A.S. Jursa, et al., Handbook of Geophysics and the Space Environment, Air Force Geophysics Laboratory, 1985. [12] P. Sokosky, Introduction to Ultralhigh Energy Cosmic Ray Physics, Addison Welsley, New York, 2004, p. 201. [13] T. Watanabe, T. Ogawa, Precision measurements of stratospheric ozone profiles by rocket-borne optical ozonesondes, Adv. Space Res. 7 (9) (1987) 123–126. [14] J.J. DeLuisi, B.G. Schuster, R.K. Sato, Separation of dust and molecular scattering contributions to the lidar observation: A method, Appl. Opt. 14 (8) (1975) 1917– 1923.
Fig. 11. Binning diagram. Two cases of binning are presented. Left: the center of a PMT is inner the FOV of the bin. Right: the center of a PMT is outer the FOV of the bin. Both in two cases, opened circle is the FOV of the PMT. Two long dashed lines represent the FOV of the bin.
and 𝑀𝐸 )]𝑀𝐸 𝑅 = 𝑅 (1.5 − 2 △ 𝜃 ′ ) sin[𝑎 cos(1.5 − 2 △ 𝜃 ′ )].
𝑠3 = 𝑅 sin[𝑎 cos( 2
(21)
Therefore, we have 𝑎 cos(1.5−2 △ 𝜃 ′ )𝑅2 −𝑅2 (1.5−2 △ 𝜃 ′ ) sin[ 𝑐𝑎 cos(1.5−2 △ 𝜃 ′ ) 𝜋𝑅2 1 ′ ◦ = 1− [𝑎 cos(1.5−2 △ 𝜃 )−(1.5 −2 △ 𝜃 ′ ) sin(𝑎 cos(1.5−2 △ 𝜃 ′ ))]. (22) 𝜋
𝑓 = 1−
286
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287 [19] A. Bucholtz, Rayleigh-scattering calculations for the terrestrial atmosphere, Appl. Opt. 34 (15) (1995) 2765–2773. [20] D.R. Longtin, E.P. Shettle, J.R. Hummel, J.D. Pryce, A Wind Dependent Desert Aerosol Model: Radiative Properties, (No. OMI-221), Optimetrics Inc, Burlington MA, 1988. [21] G. Aiellia, et al., the ARGO-YBJ collaboration, The status of the ARGO experiment at YBJ, Nucl. Phys. B (Proc. Suppl.) 166 (2007) 96–102.
[15] G. Placzek, The Rayleigh and Raman Scattering, Lawrence Radiation Laboratory, 1959. [16] G. Mie, Contributions to the optics of turbid media, particularly of colloidal metal solutions, Ann. Phys 25 (3) (1908) 377–445. [17] L. Elterman, UV, Visible, and IR Attenuation for Altitudes to 50 km, 1968, Air Force Cambridge Researsh Labs., 1968. [18] J.D. Spinhirne, J.A. Reagan, B.M. Herman, Vertical distribution of aerosol extinction cross section and inference of aerosol imaginary index in the troposphere by lidar technique, J. Appl. Meteorol. 19 (4) (1980) 426–438.
287
Contents lists available at ScienceDirect
Nuclear Inst. and Methods in Physics Research, A journal homepage: www.elsevier.com/locate/nima
A calibration of WFCTA prototype telescopes using 𝑁2 laser J.L. Liu a, *, F.R. Zhu b , H.Y. Jia b , Y.X. Bai c , Z. Cao c , M.J. Chen c , S.Z. Chen c , S.H. Feng c , B. Gao c , M.H. Gu c , H.H. He c , C. Hou c , J. Liu c , L.L. Ma c , G. Xiao c , M. Zha c , B.K. Zhang d , S.S. Zhang c , Y. Zhang c, *, J. Zhao c , B. Zhou c , X. Zuo c , for the LHAASO collaboration a b c d
Physics Department, Kunming University, Kunming, 650214, PR China School of Physical Science and Technology, Southwest Jiaotong University, Chengdu, 610031, PR China Key Laboratory of Particle Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences, Beijing 100049, PR China Physics Department, Normal College of Fuyang, Fuyang, 236037, PR China
a r t i c l e
i n f o
Keywords: Detector calibration Effective area Light flux Uncertainty
a b s t r a c t An end-to-end calibration method, which is different from the traditional way, has been employed in WFCTA prototype experiment based on 𝑁2 laser. As we known, the calculation of the expected and the observed light fluxes at the front of the aperture of the telescope are crucial in the such method. To calculate predicted light flux accurately, the propagation and scattering of photons in the atmosphere have been considered in detail. In order to get the detector effective areas, which will be used in the transformation of ADC signals to the observed light flux, a ray trace procedure, considering all details of detector response, has been performed for all incident direction of photons in field of view of the detector. By matching the observed and expected light flux at the entrance of the telescope, the overall calibration gains of two telescopes can be obtained. During the test run of the calibration, a overall calibration gain around 𝐺 = 0.38 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 and 𝐺 = 0.36 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 for WFCT01 and WFCT02 are given, respectively. Furthermore, in view of all uncertainty sources, systematic uncertainties for the calibration have also been estimated. They are around the level of (−7.3%, +7.8%). The calibration results achieved by end-to-end method and the traditional way are in consistent at the level of 10%, which is within the range of uncertainties. Since the end-to-end and traditional calibration methods have their own advantages, the performance of the combination of these two calibration methods has also been discussed in the end. © 2017 Elsevier B.V. All rights reserved.
1. Introduction The energy spectrum of cosmic ray (CR) contains important clues concerning the production and acceleration mechanisms in galactic and extragalactic particle sources. Over the past few decades, many ground-based experiments have been constructed focusing on its measurement in different energy regions above several TeV. However, because the dynamic range of individual experiments is limited and each experiment has its own energy scale and uncertainty, the observed results are not consistent among them. The precise measurement of the spectrum at a same observatory with sufficient overlap between different energy regions, which allows to obtain a complete and selfconsistent observation over the whole energy range from several TeV to several EeV, is necessary. The project of the Wide Field of View (FOV) Cherenkov/Fluorescence Telescope Array (WFCTA), composed
of 24 telescopes, which is a main component of Large High Attitude Air Shower Observatory (LHAASO) [1–3], has been designed to reach this goal. Using the atmosphere as the absorber and the scintillator of a large calorimeter, WFCTA prototypes [4] observe CR by measuring Cherenkov lights induced in air showers. The key issue in CR energy spectrum measurement in such an experiment is the determination of the total energy of reconstructed air showers which depend on an evaluation of the convertor from ADC counts recorded by electronics to a light flux at the aperture of the telescope. That is, it requires an absolute detector calibration. Various techniques to calibrate the response of such telescopes have been performed or proposed. One is using CRs to give an overall calibrated response of the telescope, including the effect of atmosphere transmission [5]. One can adjust a global calibration factor
* Corresponding authors.
E-mail addresses: [email protected] (J.L. Liu), [email protected] (Y. Zhang). https://doi.org/10.1016/j.nima.2017.09.041 Received 11 February 2017; Received in revised form 29 August 2017; Accepted 17 September 2017 Available online 6 October 2017 0168-9002/© 2017 Elsevier B.V. All rights reserved.
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
by comparing the MC-predicted CR counting rate and the measured counting rate. One main drawback of this technique is that the predicted rates depend on the CR flux and especially the composition. This would be involved large uncertainties. Another technique is comparing Cherenkov ring images in the focal plane produced by predicted and measured CR muons [6,7]. However, complications may arise from the difference in the wavelength spectrum of Cherenkov light from muons and from air showers. In addition, the predicted feature of Cherenkov rings by muons depending on the muon spectrum. Therefore one only can measure a spectrum weighted average response to the detector. The using of a fast pulsed calibrated light source would make the absolute calibration of these telescopes more straightforward and reliable. The traditional way is to calibrate all components of the detector piece by piece by using a calibrated light source, such as the mirror reflectivity, the transmission of the filter at the entrance, the quantum efficiency, photon–electron collecting efficiency, amplifier gain of photon multiplier tube (PMT), and the electronic gain. Details of all components can be found in [4]. The product of these factors provides the overall calibration factor. This is the way in current calibration and data analysis for WFCTA prototype experiment [8,9]. This method is effective on investigating the detector response in somewhat details, for instance the wavelength dependence of the quantum efficiency and the filter transmission. Difficulties for this method are (i) one cannot make sure whether or not all effects are considered, and (ii) every single piece introduces a certain of uncertainty and there is no well understanding way to estimate the combination effect of these errors to the overall calibration factor, due to the interaction between components of the system. Another approach that is so called end-toend method, which gives a overall convertor factor from ADC count to the light flux at the telescope entrance, has been applied here. One can estimate the uncertainty well in this method. The shortage is that it cannot provide a full understanding of the wavelength dependence of the detector gain. Auger experiment has applied this method to calibrate the telescope. A drum illuminator consists of a pulsed UV LED, emitting around 375 nm has been positioned at the entrance of the telescope under calibration. Absolute calibration constants were obtained from the known pulsed flux of photons emitted by the drum, and the corresponding ADC pulse integrals of the camera pixels [10]. In this paper, another end-to-end calibration method which uses a 𝑁2 laser mounted at several kilometers away from the telescope has been performed. In 2011, a 𝑁2 laser has been set up in the FOV of WFCTA prototype telescopes. Under the same condition as telescopes were operated in CR observation and some well chosen clear nights, the laser pulse were shot with several fixed directions respectively. After propagating in the atmosphere, telescopes recorded signals from the laser pulse along the track. Based on the calculation of the light flux at the front of the detector, one can directly calibrate the detector by analyzing the measured ADC counts. Uncertainties of this method are mainly from the laser energy probing, the light calculation related with aerosol and Rayleigh scattering models, as well as ADC counts. The remainder of this paper is arranged as follows. The light scattering and propagation model in the atmosphere is described in Section 2. Section 3 is the briefly description of the detector and the laser system. The status of data sets, the calculation of light flux, the detector response to the photons, as well as the estimation of the effective detector area, together with the calibration are presented in Section 4. Section 5 is the discuss about the calibration uncertainties and the performance of the combination of two calibration methods. Section 6 is the summary, as well as the prospect of detector calibration in LHAASO-WFCTA experiment in the future.
Fig. 1. Ozone level distribution as a function of altitude.
photons are emitted in a direction close to the incidence direction of the primary particle. Attenuation of the atmosphere for these lights is due to both absorptive process and scattering process on their way to the detector. Absorptive process is dominant at ultraviolet wavelengths where ozone component of the atmosphere begins to absorb the lights, and above 800 nm where vapor and CO2 absorptions begin [11,12]. In the wavelength range of 280 nm to 700 nm where the ground telescope is sensitive, scattering process initiated by molecules and aerosols in the atmosphere are dominant. That is, the contribution of absorption to the light attenuation can be negligible except for ozone absorption. Ozone absorption effect can be expressed by its transmission coefficient 𝑇𝑂 [13]: 𝑇𝑂 = exp(−𝐴𝑂3 𝐶𝑂3 (𝑙2 − 𝑙1 ))
(1)
where 𝐴𝑂3 is the ozone absorptive coefficient depending on the wavelength, and 𝐴𝑂3 = 0.117∕(cm atm) at 𝜆 = 337 nm. 𝐶𝑂3 represents the ozone level varying with the altitude, which is shown in Fig. 1 [13]. The molecule and aerosol scattering processes that contribute to the overall attenuation can be treated separately expressed by transmission coefficients 𝑇𝑅 and 𝑇𝐴 , respectively [14]. The molecule scattering is the most straightforward in two processes and can be described by Rayleigh Scattering [15]. According to this theory, the number of photons scattered out of the beam per unit length is 𝑁𝜌 400 4 𝑑𝑁 =− ( ) 𝑑𝑙 𝑋𝑅 𝜆
(2)
where 𝜆 is the wavelength of the photon and in an unit of nm. 𝑋𝑅 is the mean free path for scattering and is 2970 g/cm2 at 𝜆 = 400 nm [12]. And 𝜌 = 𝜌0 exp(−ℎ∕𝐻)
(3) 0.00129 g/cm3 , which is the
is the density profile of the atmosphere. 𝜌0 = atmospheric density at 0 ◦ C at sea level for an isothermal atmosphere. 𝐻 is the scale height of the atmosphere and is usually assumed to be 7.5 km. This parameter varies with the pressure and the temperature of the atmosphere. Thus it will be considered as an uncertainty source in light propagation. By solving function (2), one can find the transmission coefficient for lights passing from a point at slant depth 𝑋1 to a point at slant depth 𝑋2 can be expressed as 𝑇𝑅 = exp[−
2. Light propagation in the atmosphere
𝑋2 − 𝑋1 400 4 ⋅( ) ]. 𝑋𝑅 𝜆
(4)
Aerosol scattering is more complicated and can be described by Mie scattering theory [16,17]. However it relies on the simplifying assumption of spherical scatters, a condition which is not always fulfilled in the field. Moreover, it also depends on the particle composition, which changes rapidly relying on the wind and the whether condition. Based
When CRs enter the atmosphere, they induce extensive air showers (EAS) of secondary particles. While charged particles in an EAS cross the atmosphere with a speed higher than the light’s in the air, Cherenkov 279
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 2. Phase functions of aerosol scattering and Rayleigh scattering.
on Mie scattering theory, when photons traveling from a point at height ℎ1 and length 𝑙1 to a point at height ℎ2 and length 𝑙2 , the average aerosol transmission coefficient can be described by 𝑙2 − 𝑙1 ) (5) Λ Λ is the aerosol scattering length. While the height is lower than the mixing layer height (about 0.1 km), Λ is almost a constant and is assumed to be 25 km. Above the mixing layer, Λ increases with the height exponentially and can be read as
𝑇𝐴 = exp(−
Λ=
𝐿𝐴 ℎ2 − ℎ1 𝐻𝐴 𝑒ℎ2 ∕𝐻𝐴 − 𝑒ℎ1 ∕𝐻𝐴
(6)
Fig. 3. Typical pulse signals recorded by the PMT electronics of WFCTA prototype in a laser event. The 𝑋 axis is the FADC channel and each channel is 20 ns. The 𝑌 axis is the FADC count in each channel.
where 𝐿𝐴 is the aerosol attenuation length describing the light attenuation due to the aerosol at the ground level. 𝐻𝐴 is the aerosol scale height (about 1.2 km) accounting for its dependence on the height [17,18]. When lights arrive at a position in the FOV of the telescope, the number of photons scattered into an solid angle from the laser beam both by molecules and aerosols depend on the angular distribution of scattered lights, which is referred as phase function. The phase function for Rayleigh scattering is well understood and can be estimated analytically as function (7). With its key feature of the symmetry in the forward and backward directions, it is proportional to (1+cos2 𝜃). A small factor 𝛾 = 0.0155 is used to correct the effect of anisotropy of the 𝑂2 and 𝑁2 molecules [19]: 3 𝑃𝑅 (𝜃) = [(1 + 3𝛾) + (1 − 𝛾) cos2 𝜃]. 16𝜋(1 + 2𝛾)
Table 1 Summary of reflector parameters in WFCTA prototypes. Mirror facet length Mirror facet area Curvature radius Reflective material Reflectivity (>240 nm)
298.4 mm 0.23 𝑚2 4740 mm 𝐴𝑙 + 𝑀𝑔𝐹2 >83%
Table 2 Summary of PMT parameters in WFCTA prototypes. Manufacturer Model Photocathode Voltage Typical voltage Wavelength Gain Nonlinearity Rise time FWHM
(7)
For the phase function of aerosol scattering 𝑃𝐴 (𝜃), it is quite complex and model dependent. Two models referred as Elbert’s model and Longtin’s wind-dependent-desert model [20] are available to describe its gross feature on the light scattering probability distribution, which is sufficient for the purpose of Cherenkov detection. Elbert’s model is used as default in the study. In general, the angular distribution of lights scattered by the aerosol is very strong peaked in the forward direction and reaches a minimum near 90 ◦ . All these different phase functions are displayed in Fig. 2 and all of them are normalized over 4𝜋 solid angle.
Photonis XP3062/FL 40 mm 900 V ∼ 1300 V 1100 V 280 nm ∼ 700 nm 2 ×105 < 2% 3 ns 4.5 ns
an area of 5.0 m2 . For the detailed parameter information, see Table 1. The camera is constituted by 256 hexagonal PMTs. Detailed information of PMTs are shown in Table 2. Each PMT has a FOV of 1◦ × 1◦ and all PMTs are arranged as a 16 × 16 array, forming a total FOV of 16◦ × 14◦ . Details of the telescope can be found in reference [4,8]. When a laser pulse passes the FOV of the detector, scattering lights will trigger tubes watching the direction of the laser track and a FADC electronics system will record pulse signals every 20 ns. Fig. 3 displays such pulse signals recorded by PMTs in a laser event. A 𝑁2 laser has been set up at the west by north of ARGO-YBJ experiment, 2520 m away from the center of the experimental hall. Its parameters are summarized in Table 3. The laser is mounted on a rotating platform, enabling the laser beam has an angular resolution
3. WFCTA prototypes and the laser set Two telescopes have been set up in WFCTA prototype experiment. They are located at the west by north (denoted as WFCT01) and the south by east (denoted as WFCT02) of the ARGO-YBJ experiment [21], about 99 m and 79 m away from the center of the ARGO hall respectively. Their pointing directions are set as 32◦ and 31◦ in zenith direction, 170◦ and 171◦ in azimuth direction respectively. Each telescope mainly consists of a UV-light collecting mirror and an image camera. The reflected mirror is constituted by 20 hexagonal spherical facets, forming 280
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 6. The geometry of lights propagation and scattering in the calibration. The shaded area represents the FOV of the telescope.
Fig. 4. The distribution of laser energies during a run.
the calculation of the angular width of the bin along the track and the correlation of the effective area of the detector. In this section, firstly, data sets in the analysis will be introduced. secondly, the calculation of the light flux at the entrance of the telescope will be described. Thirdly, the binning and the estimation of detector effective area will be given. Finally, based on above studies, the calibration procedure and the result will be presented in detail.
Table 3 𝑁2 laser parameters. Manufacturer Model Wavelength Spectrum width Energy Pulse duration Stability Beam size Divergence
Stanford research systems NL100 337.1 nm 0.1 nm 170 μJ @ 10 Hz 3.5 ns 3% 3 × 7 mm 3 × 8 mrad
4.1. Data sets To avoid the influence of the aerosol to the calibration, the data that observed at a very clear night is chosen. Furthermore, as indicated in Fig. 2, values of phase functions of aerosol scattering both in two models are relatively smaller while the scattering angle is around 90◦ . Therefore the zenith of the laser pulse was set to around 60◦ , which makes the scattering angle being around 90◦ . Thirdly, to avoid the influence of marginal effect to the calculation of the effective area, the data with image on the center area of FOV are chosen. After data selection, three data sets, that is, 410 events with (𝜃, 𝜑) = (60◦ , −8◦ ), 732 events with (𝜃, 𝜑) = (50◦ , −8◦ ) and 684 events with (𝜃, 𝜑) = (60◦ , −10◦ ) (where 𝜃 is the zenith and 𝜑 is the azimuth), are available in the calibration. The former two are for WFCT01 and the latter is for WFCT02.
of 0.1◦ and a spatial resolution of 0.1 mm. To minimize the influence of the aerosol, the calibration has been carried out at very clear nights. The laser shots 337.1 nm photons with a pulse duration of 3.5 ns, elevations of 30◦ ∼ 60◦ and azimuths of −22◦ ∼ 4◦ . To match the laser pulse energy with the CR energy and to maintain a precise beam direction, a collimator has been mounted at the entrance of the laser. The absolute energy of the laser beam is measured by an energy probe, which has a calibration uncertainty of less than 5%. Energies of laser beams during a run shown in Fig. 4 indicate that the average energy of laser pulses is 90.8 μJ with a small fluctuation of 1.2%.
4.2. Light flux calculation 4. Absolute calibration of WFCTA prototypes The calculation of the light flux at the aperture of the telescope is crucial for end-to-end calibration. The calculation of the predicted light flux will be described in detail in this section. For given wavelength 𝜆 and energy 𝐸0 , the total number of photons in a laser pulse is 𝑁0 = 𝐸0 ∕(ℎ𝑐∕𝜆). ℎ is Planck constant and 𝑐 is light velocity here. For this study with 𝐸0 = 90.8 μJ and 𝜆 = 337.1 nm, photons emitted into the air in a pulse is 𝑁0 = 1.54×1014 . As shown in Fig. 6, these
As mentioned in Section 1, the matching of predicted light flux and observed light flux at the aperture of the detector is the key to the end-toend calibration method. As displayed in Fig. 5, the image of an observed laser event is the assemble of trigged PMTs. To transform recorded ADC counts for each PMT in the observation to the light flux in front of the detector, one not only needs the calibration factor, but also requires
Fig. 5. The typical laser event observed by WFCTA prototypes. The left panel is for WFCT01 and the right panel is for WFCT02. Green hexagons are PMTs on the image plane. Red solid circles represent signals collected by triggered PMTs. The diameter of circles is proportional to the logarithm of the amplitude of signals.
281
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Table 4 The flux of photons in front of the telescope, which are scattered from a 90.8 μJ 𝑁2 laser pulse at 1240 m above the ground. The unit of the flux is 𝑝ℎ𝑜𝑡𝑜𝑛𝑠∕𝑚2 ∕𝑑𝑒𝑔𝑟𝑒𝑒. The aerosol scale height 𝐻𝐴 is assumed to be 1.2 km for all combinations of aerosol scattering parameters here. Aerosol
𝑃𝐴 (𝜃)
Elbert Longtin Rayleigh
In order to include enough tubes in a bin so as to avoid a large fluctuation for the effective area due to the marginal tubes, a angular size of 1.5◦ for a bin is chosen in the binning. Once all bins are defined, all none-noised tubes are sorted into bins. If a tube’s center in SDP lies within 0.25◦ of a bin center, it is considered wholly within the bin (i.e. the fraction 𝑓 of the contribution for the tube to this bin is 1). If a tube’s center lies greater than 0.25◦ and less than 1.25◦ from a bin center, it is considered partly in the bin and a fraction 𝑓 is calculated according to formula (11) to determine how much of the tube is contained into the bin.
𝐿𝐴 12 km (standard)
25 km (average)
1000 km (Rayleigh)
345 226 1040
182 119 1145
5 3 1249
1 𝑓 =1− [𝑎 cos(1.5◦ −2 △ 𝜃 ′ )−(1.5◦ −2 △ 𝜃 ′ ) sin(𝑎 cos(1.5◦ −2 △ 𝜃 ′ ))] (11) 𝜋 lights will undergo a path from the laser (𝐿) to a scatter (𝑆) in the FOV of the detector at first. After scattering to the direction of the detector (𝐷), they propagate to the detector through the atmosphere. Propagating along the path 𝐿𝑆 and 𝑆𝐷, lights suffer the same scattering as at position 𝑆. While arrived at position 𝑆, the number of photons decrease to 𝑁1 for attenuations of Rayleigh scattering 𝑇𝑅1 , aerosol scattering 𝑇𝐴1 and ozone absorption 𝑇𝑂1 .
where △𝜃 ′ is the angular distance between the tube center and the bin center. The deduction of this function is presented in Appendix. If a tube’s center lies greater than 1.25◦ from the bin center, it is not included in the bin. Once all tubes are sorted into bins, the ADC count for each bin is calculated following 𝜙0 =
where the calculation of 𝑇𝑅1 , 𝑇𝐴1 and 𝑇𝑂1 follows Eqs. (4), (5) and (1), respectively. At location 𝑆 in the FOV of the detector, the number of photons scattered to a solid angle corresponding to the mirror size from the laser beam depends on the angular distribution of scattered lights, which is referred as phase function. The number of photons scattered from an unit angular track length and reached to an unit area of the mirror can be read as 𝑁2 = 𝑁1 [(1 − 𝑇𝑅 (𝑆))𝑃𝑅 (𝜃) + (1 − 𝑇𝐴 (𝑆))𝑃𝐴 (𝜃)]
1 1 2 Δ𝜃 𝑙𝑆𝐷
(12)
4.4. Effective area of the detector As mentioned at the beginning of this section, in order to transform the ADC signal to the light flux at the front of the detector, one needs to know the effective light collecting area of the detector. The mirror area of 5.0 m2 mentioned in Section 3 is only the geometric area. The true effective area not only depends on the geometry effect associated with the direction of incident rays, but also depends on the location and the size of the spot on the camera. That is because PMTs have different sensitivity varying with the location of the photocathode and there has insensitive area between hexagonal PMTs. Therefore the ratio of sensitive area to insensitive area changes as the spot moves around on the image plane. The effective area has to be determined using simulation for every position and direction of the light source as the source moves over the FOV of the detector. As described in Section 4.3, when a bin of the track has been defined, the vector of the bin center has been calculated simultaneously. Five thousands trial incident rays (𝑁𝑡𝑟𝑖𝑎𝑙 ) hit on the mirror from directions stepping down the track in a range of 3.5◦ centered at the bin center vector. A simulation procedure named as ray tracing has been employed for each incident ray to calculate the effective area for each bin. In the ray tracing procedure, it has considered the shape and the volume of the system container and the PMT cluster, and the detailed structure of the reflector and the camera, as well as shadows of the support structure of the filter in front of the telescope entrance and the cluster in the ray tracing procedure. By this ray tracing, one can calculate the effective area 𝑁 following 𝐴𝑒𝑓 𝑓 = 𝐴0 𝑁𝑠𝑢𝑐𝑐 of the bin, where 𝐴0 is the mirror geometric 𝑡𝑟𝑖𝑎𝑙 area on which the trial rays uniformly hit and 𝑁𝑠𝑢𝑐𝑐 is the number of rays which successfully reach to the effective area of PMTs. When a tube collects a trial ray, we check whether the tube is triggered during the measurement. If the tube is not triggered, the ray will not be counted into 𝑁𝑠𝑢𝑐𝑐 . As described in Section 4.1, events in three data sets have different geometry. Thus the geometry of photons received by the telescope are also different. Taking the angle (denoted as 𝑑𝑡ℎ𝑒) between the incident direction of photons and the primary optical axis of the telescope as an example, average values for three data sets are 2.5◦ , 2.8◦ and 4.5◦ , respectively. 𝐴𝑒𝑓 𝑓 distributions in Fig. 7 show that average values of 𝐴𝑒𝑓 𝑓 for three data sets are 2.29 m2 , 2.26 m2 and 2.11 m2 respectively. It is in agreement with the result in bottom-right panel of Fig. 7. That is 𝐴𝑒𝑓 𝑓 decreases with the increase of 𝑑𝑡ℎ𝑒. To study the uncertainty of the calculation of 𝐴𝑒𝑓 𝑓 , measured errors of all parameters associated with the detector in ray trace simulation have been estimated. They are curvature radius of mirror, 𝑅𝑚𝑖𝑟 , diameter
(9)
where Δ𝜃 is the angular bin size around the position 𝑆, 𝑙𝑆𝐷 is the path length from 𝑆 to 𝐷. Photons scattered to the direction of the detector undergo a same attenuation as propagating along the path 𝐿𝑆. Thus the photons (i.e. the light flux) decrease to 𝑁3 while at location 𝐷. 𝑁3 = 𝜙exp = 𝑁2 𝑇𝑅2 𝑇𝐴2 𝑇𝑂2
𝑓 𝑁𝐴𝐷𝐶 (𝑖)
where 𝑁𝐴𝐷𝐶 (𝑖) is the total ADC count recorded by the 𝑖th tube.
(8)
𝑁1 = 𝑁0 𝑇𝑅1 𝑇𝐴1 𝑇𝑂1
∑
(10)
where the calculation of 𝑇𝑅2 , 𝑇𝐴2 and 𝑇𝑂2 are same as 𝑇𝑅1 , 𝑇𝐴1 and 𝑇𝑂1 except for replacing the integral path 𝐿𝑆 with 𝑆𝐷. To understand the uncertainty associated with different atmospheric parameters, such as phase function models, scattering length, attenuation length, aerosol distribution etc., the flux at the front of the detector has been calculated under different atmospheric conditions. As an example, Table 4 gives the results of light flux scattered from a given position 𝑆 on the laser beam which is 1240 m above the ground. From the results listed in Table 4, the contribution of aerosol scattering can be as small as 0.24%–0.40% according to different aerosol scattering phase function as long as the atmosphere is clear. At other hand, even under the ‘‘average’’ atmosphere condition, that is 𝐿𝐴 = 25 km and 𝐻𝐴 = 1.2 km [12,17,18], the uncertainty of light flux associated with phase function is about 4.7%, although the aerosol scattering contribution can be as large as 13.7%. 4.3. Angular binning for the track There are several major steps in binning of the laser track. (i) Remove noise tubes according to triggered PMTs’ location on the image plane and their trigger times. (ii) Determine the Laser-Detector-Plane (LDP) based on the direction of the laser pulse and locations of the laser and the detector. (iii) Project all none-noised tubes into the LDP. (iv) Find the topmost tube and mark the first bin center. (v) Step down along the laser track in LDP with 1.5𝑜 interval and calculate the vector which is from detector pointing to each bin center. (vi) Sort tubes into bins, possibly splitting a tube among different bins by assigning a weighting factor for the tube. 282
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 7. The distribution of effective area for WFCT01 (top) and WFCT02 (bottom-left), and as a function of 𝑑𝑡ℎ𝑒 for both two telescopes (bottom-right). Table 5 Nominal values and estimated uncertainties for all parameters used in ray trace simulation. The unit for all parameters is mm. Parameter
𝑅𝑚𝑖𝑟
𝐷𝑚𝑖𝑟
𝐷𝑚𝑐
𝑊𝑐
𝐻𝑐
𝑇𝑐
𝑅𝑠𝑝𝑜𝑡
Nominal value Uncertainty
4740 ±10
596.8 ±3
2305 ±3
711 ±5
635 ±5
355 ±5
15 ±3
337.1 nm. 𝑄𝑒 , 𝜀 and 𝐺𝑃 𝑀𝑇 represent the quantum efficiency, photon electron collecting efficiency and gain of PMT respectively, and 𝐺𝑒 is the electronics gain. Measuring these four effects item by item is difficult. In the calibration, a UV-LED at 375 nm mounted at the center of the mirror is applied to obtain the overall response for each PMT [4]. For the calibration is done at 375 nm, a further correction is applied according to the wavelength dependence of the quantum efficiency for the laser calibration at 337.1 nm. The other way of the calibration is to measure the overall factors using a laser beam, which is so called end-to-end method. In such method, for the given energy and geometries of the laser and the detector, one can calculate the expected light flux 𝜙exp at the front of the detector following the method described in Section 4.2. Then one can determine the factor 𝐺 by matching the predicted light flux in Eq. (10) and the observed light flux in Eq. (13). That is
of a mirror segment, 𝐷𝑚𝑖𝑟 , distance from mirror center to cluster center, 𝐷𝑚𝑐 , dimensions of cluster box, 𝑊𝑐 (width), 𝐻𝑐 (height), 𝑇𝑐 (thickness), average spot size at the cluster, 𝑅𝑠𝑝𝑜𝑡 . Details of nominal value and estimated error for these parameters are listed in Table 5. The overall uncertainty to 𝐴𝑒𝑓 𝑓 caused by these parameters is ±2%. 4.5. Detector calibration Once the total ADC signal 𝜙0 and the effective area 𝐴𝑒𝑓 𝑓 for each bin along the track are calculated, the observed light flux at the front of the detector can be expressed as 𝜙𝑜𝑏𝑠 =
𝜙0 𝐴𝑒𝑓 𝑓 △ 𝜃𝐺
.
𝐺=
.
(15)
In a laser event, a gain factor can be obtained from each bin along the track. To avoid the influence of the fluctuation, only a gain factor is given in a event by averaging calibration factors obtained from all bins. All calibration results are presented in Fig. 8. The analysis of both two data sets of WFCT01 (top panels) give a consistent calibration constant by end to end calibration method, that is 𝐺 = 0.38 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 in average. For WFCT02 (the bottom panel), the calibration result is 𝐺 = 0.36 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 in average. Calibration results obtained from piece by piece method are also shown to give a comparison. Average gain factors from piece by piece method are about 10% larger than ones from end-to-end method, which is within the level of the calibration uncertainties.
(13)
The goal of the calibration is to determine the value of 𝐺. As described in Section 1, one way is of the calibration to directly measure all factors in Eq. (14) piece by piece, which is so called traditional method. 𝐺𝑝 = 𝑅𝑚𝑖𝑟𝑟𝑜𝑟 𝑇𝑓 𝑖𝑙𝑡𝑒𝑟 𝑄𝑒 𝜀𝐺𝑃 𝑀𝑇 𝐺𝑒
𝜙0 𝜙exp 𝐴𝑒𝑓 𝑓 △ 𝜃
(14)
where 𝑅𝑚𝑖𝑟𝑟𝑜𝑟 = 0.83 is the mirror reflectivity, 𝑇𝑓 𝑖𝑙𝑡𝑒𝑟 = 0.865 is the transmission of the filter at the entrance, both two values are at 𝜆 = 283
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 8. The distribution of calibration gain by end-to-end (solid curves) and piece by piece (dashed curves) methods. Top panels are for WFCT01 and the bottom panel is for WFCT02.
5. Discuss
5.2. Uncertainties on end-to-end calibration
5.1. Calibration method
As stated at the end of Section 1, sources of systematic uncertainty in end-to-end calibration are (i) laser energy probing, ±5%, (ii) aerosol scattering model, see Table 4, (iii) the calculation of Rayleigh scattering due to the uncertainty of the temperature and pressure of the atmosphere, ±2%, (iv) the calculation of 𝐴𝑒𝑓 𝑓 , ±2%, see Table 5, and (v) ADC count error due to the overall effect of the fluctuation and the systematic error of electronics initiated by the variation of the temperature of the external environment, ±2%. (i) and (ii) are the most significant effects. To study the influence of aerosol scattering model, we have changed the horizontal attenuation and phase function displayed in Table 4, respectively. The error contribution of this source is (−4.1%, +4.9%). The overall uncertainty to predicted light flux introduced both by (i), (ii) and (iii) is (−6.8%, +7.3%), see Fig. 9. (iv) and (v), together with the uncertainty of 𝐺𝑝 (±7%) [4], contribute to the uncertainty of observed light flux. The overall contribution is (−7.5%, +7.5%), also see Fig. 9. For the calibration gain factor 𝐺 under end-to-end method, the overall uncertainty caused both by above sources is (−7.3%, +7.8%). And for 𝐺𝑃 , the overall uncertainty is ±7%. From Fig. 8, one can see that the difference of average gain factors obtained by two calibration method is around 10%. So calibration gains achieved by two methods are consistent within the level of uncertainties.
As mentioned in Section 1, in the calibration, end-to-end method and the traditional way have their own advantages and disadvantages. The combination of two methods provides the best way to do the calibration. The gain 𝐺𝑝 obtained by the traditional method is applied in Eq. (13) to get the observed light flux at the front of the detector. The expected light flux in Eq. (10) is used to find an adjustment 𝑓𝑎𝑑𝑗 for the current gain 𝐺𝑝 of the detector. The predicted and observed light fluxes along the track at the front of the detector obtained by above method are presented in Fig. 9. The topmost and the bottommost bins are affected seriously by the marginal effect, thus these two bins are removed. They indicate that the observed fluxes are about 10% lower than predicted. To avoid the influence of the ∑ fluctuation for the observed light flux, an integral flux 𝜙 is taking for each laser event. Results of integral fluxes for three data sets are shown in Fig. 10. The adjustment of the gain in current data analysis at 337.1 nm can be obtained by ∑ 𝜙 𝑓𝑎𝑑𝑗 = ∑ 𝑜𝑏𝑠 . (16) 𝜙exp
6. Summary
In this study, only one night data is available. Thus we just present the calibration method and procedure, and cannot give the conclusion of 𝑓𝑎𝑑𝑗 measurement in this paper.
Through this study, first, a complete analysis procedure and an end-to-end calibration method have been established for the absolute 284
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 9. The distribution of predicted (open circles) and observed (solid circles) light fluxes as a function of the light source height for WFCT01 (top) and WFCT02 (bottom) telescopes. The shades around the points represent the uncertainties.
Appendix. Derivation of function (11)
calibration of WFCTA telescopes, including the detailed calculation of light scattering and propagation in the atmosphere, the shower track binning, and the ray trace which is used to calculate the effective area of detector upon different incident directions. Second, by using this method, in the test run, a overall calibration gain around 𝐺 = 0.38 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 and 𝐺 = 0.36 𝐴𝐷𝐶 𝑐𝑜𝑢𝑛𝑡𝑠∕𝑝ℎ𝑜𝑡𝑜𝑛 for WFCT01 and WFCT02 are given. Third, both the advantage and disadvantage of traditional and end-to-end calibration methods are discussed. Based on the discussion, the performance of the combination of two methods in the calibration is investigated. Finally, the uncertainty of the calibration considering all error contributions is estimated. They are in the level of (−7.3%, +7.8%). The calibration results obtained by end-to-end method and the traditional way are in consistent at the order of 10%, which is within the level of uncertainties. In the future, we will monitor the aerosol condition hourly, as well as monitor the temperature, the pressure and humidity of the atmosphere in the calibration of LHAASOWFCTA telescopes. This will help to decrease the calibration uncertainty caused by atmosphere model significantly.
Eq. (11) in the article is used to calculate the fraction of a PMT contributing to the bin according their central vectors. Two cases of the binning are shown in Fig. 11. That is, the center of a PMT is inner the FOV of the bin (left) and is outer the FOV of the bin (right). Since the whole derivations of Eq. (11) in such two cases are the same, for concise, only the left case is provided as follows: For convenience, we denote by 𝑅 ∶ 𝑡ℎ𝑒 𝑟𝑎𝑑𝑖𝑢𝑠 𝑜𝑓 𝑡ℎ𝑒 𝑐𝑖𝑟𝑐𝑙𝑒 𝐴𝐵𝐶𝐷. 𝑠 ∶ 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑡ℎ𝑒 𝑐𝑖𝑟𝑐𝑙𝑒 𝐴𝐵𝐶𝐷. 𝑠0 ∶ 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑡ℎ𝑒 𝑠𝑒𝑔𝑚𝑒𝑛𝑡 𝐴𝐷𝐶. 𝑠1 ∶ 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑡ℎ𝑒 𝑠𝑒𝑔𝑚𝑒𝑛𝑡 𝐴𝐵𝐶. 𝑠2 ∶ 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑡ℎ𝑒 𝑠𝑒𝑐𝑡𝑜𝑟 𝐴𝑀𝐶𝐵. 𝑠3 ∶ 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑡ℎ𝑒 𝑡𝑟𝑖𝑎𝑛𝑔𝑙𝑒 𝐴𝑀𝐶. 𝑓 ∶ 𝑡ℎ𝑒 𝑝𝑟𝑜𝑝𝑜𝑟𝑡𝑖𝑜𝑛 𝑜𝑓 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑠𝑒𝑔𝑚𝑒𝑛𝑡 𝐴𝐷𝐶 𝑖𝑛 𝑡ℎ𝑒 𝑐𝑖𝑟𝑐𝑙𝑒 𝐴𝐵𝐶𝐷. 𝑓1 ∶ 𝑡ℎ𝑒 𝑝𝑟𝑜𝑝𝑜𝑟𝑡𝑖𝑜𝑛 𝑜𝑓 𝑡ℎ𝑒 𝑎𝑟𝑒𝑎 𝑜𝑓 𝑠𝑒𝑔𝑚𝑒𝑛𝑡 𝐴𝐵𝐶 𝑖𝑛 𝑡ℎ𝑒 𝑐𝑖𝑟𝑐𝑙𝑒 𝐴𝐵𝐶𝐷. △𝜃 ′ ∶ 𝑡ℎ𝑒 𝑎𝑛𝑔𝑙𝑒 𝑏𝑒𝑡𝑤𝑒𝑒𝑛 𝑡ℎ𝑒 𝑡𝑢𝑏𝑒 𝑐𝑒𝑛𝑡𝑒𝑟 𝑎𝑛𝑑 𝑡ℎ𝑒 𝑏𝑖𝑛 𝑐𝑒𝑛𝑡𝑒𝑟. Noting that 𝑠1 , 𝑠 in which 𝑠1 = 𝑠2 − 𝑠3 . By the area formula of a sector, we obtain that
Acknowledgments This work is supported by NSFC, Yunnan provincial Science and Technology Department and Kunming University through grants of 11563004, 2014FB154 and YJL17004, respectively. It is also funded by Program for Innovative Research Team (in Science and Technology) in University of Yunnan Province. We acknowledge the essential support of W.Y. Chen, G. Yang, X.F. Yuan, C.Y. Zhao. We are grateful to Prof. Xinhua Bai for his valuable comments and professional revision suggestions (SDSM&T). We are also thankful to Dr. Chang-Bo Yang for the discussions of mathematical problems with him.
𝑓 = 1 − 𝑓1 = 1 −
(17)
1 2 𝜃𝑅 , 2
(18)
𝑠2 =
where 𝜃 = 2𝑎 cos( 𝑀𝐸 ), 𝑀𝐸 = 0.75 − △𝜃 ′ and 𝑅 = 0.5. 𝑅 By simple computation, we have
285
𝜃 = 2𝑎 cos(1.5 − 2 △ 𝜃 ′ ),
(19)
𝑠2 = 𝑎 cos(1.5 − 2 △ 𝜃 ′ )𝑅2 ,
(20)
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287
Fig. 10. The distribution of integral predicted (dashed curves) and observed (solid curves) light fluxes for WFCT01 (top) and WFCT02 (bottom) telescopes.
References [1] Z. Cao, A future project at Tibet: The large high altitude air shower observatory (LHAASO), Chin. Phys. C 34 (02) (2010) 249–252. [2] M. Zha, the ARGO-YBJ collabaration, the LHAASO collabaration, Status of the large high altitude air shower observatory project, Nucl. Instrum. and Meth. in Phys. Res. A 692 (2012) 77–82. [3] Z. Cao, the LHAASO Collaboration, Status of LHAASO updates from ARGO-YBJ, Nucl. Instrum. and Meth. in Phys. Res. A 742 (2014) 95–98. [4] S.S. Zhang, Y.X. Bai, Z. Cao, et al., Properties and performance of two wide field of view Cherenkov/fluorescence telescope array prototypes, Nucl. Instrum. and Meth. in Phys. Res. A 629 (2011) 57–65. [5] A. Konopelkoa, et al., the HEGRA Collaboration, Detection of gamma rays above 1 TeV from the Crab Nebula by the second HEGRA imaging atmospheric Cherenkov telescope at La Palma, Astropart. Phys. 4 (3) (1996) 199–215. [6] T.B. Humenskya, for the VERITAS Collaboration, Calibration of VERITAS Telescope 1 via Muons, 2005. arXiv:astro-ph/0507449v1. [7] G. Vacanti, P. Fleury, Y. Jiang, et al., Muon ring images with an atmospheric Cherenkov telescope, Astropart. Phys. 2 (1) (1994) 1–11. [8] S.S. Zhang, et al., the ARGO Collaboration, LHAASO Collaboration, Energy spectrum of cosmic protons and helium nuclei by a hybrid measurement at 4300 m a.s.l., Chin. Phys. C 38 (4) (2014) 045001. [9] S.S. Zhang, et al., the ARGO Collaboration, LHAASO Collaboration, The knee of the Cosmic Hydrogen and Helium Spectrum below 1 PeV measured by ARGO-YBJ and a Cherenkov telescope of LHAASO, Phys. Rev. D 92 (092005) (2015). [10] A.A. Watson, the Auger collabaration, Properties and performance of the prototype instrument for the Pierre Auger Observatory, Nucl. Instrum. and Meth. in Phys. Res. a 523 (2004) 50–95. [11] A.S. Jursa, et al., Handbook of Geophysics and the Space Environment, Air Force Geophysics Laboratory, 1985. [12] P. Sokosky, Introduction to Ultralhigh Energy Cosmic Ray Physics, Addison Welsley, New York, 2004, p. 201. [13] T. Watanabe, T. Ogawa, Precision measurements of stratospheric ozone profiles by rocket-borne optical ozonesondes, Adv. Space Res. 7 (9) (1987) 123–126. [14] J.J. DeLuisi, B.G. Schuster, R.K. Sato, Separation of dust and molecular scattering contributions to the lidar observation: A method, Appl. Opt. 14 (8) (1975) 1917– 1923.
Fig. 11. Binning diagram. Two cases of binning are presented. Left: the center of a PMT is inner the FOV of the bin. Right: the center of a PMT is outer the FOV of the bin. Both in two cases, opened circle is the FOV of the PMT. Two long dashed lines represent the FOV of the bin.
and 𝑀𝐸 )]𝑀𝐸 𝑅 = 𝑅 (1.5 − 2 △ 𝜃 ′ ) sin[𝑎 cos(1.5 − 2 △ 𝜃 ′ )].
𝑠3 = 𝑅 sin[𝑎 cos( 2
(21)
Therefore, we have 𝑎 cos(1.5−2 △ 𝜃 ′ )𝑅2 −𝑅2 (1.5−2 △ 𝜃 ′ ) sin[ 𝑐𝑎 cos(1.5−2 △ 𝜃 ′ ) 𝜋𝑅2 1 ′ ◦ = 1− [𝑎 cos(1.5−2 △ 𝜃 )−(1.5 −2 △ 𝜃 ′ ) sin(𝑎 cos(1.5−2 △ 𝜃 ′ ))]. (22) 𝜋
𝑓 = 1−
286
J.L. Liu et al.
Nuclear Inst. and Methods in Physics Research, A 877 (2018) 278–287 [19] A. Bucholtz, Rayleigh-scattering calculations for the terrestrial atmosphere, Appl. Opt. 34 (15) (1995) 2765–2773. [20] D.R. Longtin, E.P. Shettle, J.R. Hummel, J.D. Pryce, A Wind Dependent Desert Aerosol Model: Radiative Properties, (No. OMI-221), Optimetrics Inc, Burlington MA, 1988. [21] G. Aiellia, et al., the ARGO-YBJ collaboration, The status of the ARGO experiment at YBJ, Nucl. Phys. B (Proc. Suppl.) 166 (2007) 96–102.
[15] G. Placzek, The Rayleigh and Raman Scattering, Lawrence Radiation Laboratory, 1959. [16] G. Mie, Contributions to the optics of turbid media, particularly of colloidal metal solutions, Ann. Phys 25 (3) (1908) 377–445. [17] L. Elterman, UV, Visible, and IR Attenuation for Altitudes to 50 km, 1968, Air Force Cambridge Researsh Labs., 1968. [18] J.D. Spinhirne, J.A. Reagan, B.M. Herman, Vertical distribution of aerosol extinction cross section and inference of aerosol imaginary index in the troposphere by lidar technique, J. Appl. Meteorol. 19 (4) (1980) 426–438.
287