Expected performance of air-shower measurements with the radio-interferometric technique
PPrepared for submission to JINST
Expected performance of air-shower measurements withthe radio-interferometric technique
F. Schlüter 𝑎,𝑏, T. Huege 𝑎,𝑐 𝑎 Karlsruher Institut für Technologie, Institut für Astroteilchenphysik, Karlsruhe, Germany 𝑏 Universidad Nacional de San Martín, Instituto de Tecnologías en Detección y Astropartículas,Buenos Aires, Argentina 𝑐 Vrije Universiteit Brussel, Astrophysical Institute, Brussels, Belgium
E-mail: [email protected]
Abstract: Interferometric measurements with arrays of radio antennas are a powerful and widelyused technique in astronomy. Recently, this technique has been revisited for the reconstructionof extensive air showers [1]. This “radio-interferometric technique” exploits the coherence inthe radio emission emitted by billions of secondary shower particles to reconstruct the showerparameters, in particular the shower axis and depth of the shower maximum 𝑋 max . The accuracypreviously demonstrated on simulations with an idealized detector is very promising. The prospectof an accurate 𝑋 max measurement for inclined air showers combined with measurements of theelectromagnetic energy (also with radio antennas) and the muonic shower content (via groundparticle detectors) is very intriguing as it would provide a high sensitivity to the mass of cosmicrays, key information to study their origin. In this article we evaluate the potential of interferometric 𝑋 max measurements using (simulated) inclined air showers with sparse antenna arrays under realisticconditions. To determine prerequisites for the application of the radio-interferometric techniquewith various antenna arrays, the influence of inaccuracies in the time synchronisation betweenantennas and its inter-dependency with the antenna density is investigated in detail. We find a strongcorrelation between the antenna multiplicity (per event) and the maximum acceptable time jitter, i.e.,inaccuracy in the time synchronisation. For data recorded with a time synchronisation accurate towithin 1 ns in the commonly used frequency band of 30 MHz to 80 MHz, an antenna multiplicityof >
50 is needed to achieve an 𝑋 max resolution of 𝜎 𝑋 max (cid:46)
20 g cm − . For data recorded with 2 nsaccuracy, already (cid:38)
200 antennas are needed to achieve this 𝑋 max resolution. Furthermore, we findno advantage reconstructing 𝑋 max from data simulated at higher observation frequencies, i.e., upto several hundred MHz. Finally, we provide a generalisation of our results from very inclined airshowers to vertical geometries.Keywords: Radio Detection, Interferometry, Extensive Air Showers, Ultra-high-energy cosmic rays Corresponding author. a r X i v : . [ a s t r o - ph . I M ] F e b ontents Radio signals from extensive air showers exhibit wave phenomena. An example for such phenomenais the nowadays well-established Cherenkov cone which results from the interference of the signalsemitted by billions of shower particles in an atmosphere with a refractive index gradient.Interferometric techniques expose this coherence in the radio emission. Thereby, both thesignal’s amplitude and phase information is used, while traditional reconstruction methods ofextensive air showers rely on the amplitude information only. Interferometric techniques are standardin radio astronomy, where sources are at infinity and hence all antennas receive the same signal witha planar wavefront. Application to radio emission from extensive air showers is more challenging,as the source is typically nearby, is extended, and the emission from different parts of the showerpropagates through different refractive index gradients.Interferometric techniques have previously been used successfully for cosmic-ray radio detectionin the LOPES experiment to identify coherent air-shower radio pulses amongst strong and time-correlated radio-frequency interference, to estimate the energy of the primary particle, and to provide– 1 –n image of the intensity distribution on the sky from which the arrival direction can be determined[2]. They have also been employed to determine the depth of shower maximum from LOPES data[3, 4] with an experimental accuracy worse than 100 g/cm but potential in pure simulations to reachan accuracy as good as 30 g/cm . Attempts to apply interferometric techniques to ground-basedradio arrays with a larger extension than the small-scale LOPES experiment, for example within theAuger Engineering Radio Array, had not been successful [5], presumably because the then-madeassumption that antennas see identical signals no longer holds for larger arrays. Another experimentroutinely using interferometric techniques to identify and reconstruct air-shower radio emissionis ANITA [6]. Finally, real-time interferometric triggering is also being investigated for particleshowers in ice [7] and air [8].In [1], the so-called radio-interferometric technique (RIT) is developed and successfully appliedto air shower simulations with an idealized detector (zenith-angle dependent dense antenna array,perfect time synchronisation between antennas and perfectly known antenna locations) to reconstructthe shower axis and depth of the shower maximum 𝑋 max with high accuracy. A resolution of betterthan 0 . ° (< 0 . ° ) in the arrival direction and 3 g cm − (10 g cm − ) in 𝑋 max for inclined (vertical)showers is demonstrated.The limited size of the footprint illuminated by the strongly forward-beamed radio emissionin vertical air showers demands the use of comparatively dense and small antenna arrays and thusrestricts the observation of cosmic rays to energies around and below 1 EeV [9, 10]. In inclined airshowers the radio-emission footprint is spread over large areas thus enabling the observation of airshowers with sparse antenna arrays [11, 12]. This allows one to instrument large areas (> 1000 km )to detect ultra-high energy cosmic rays (UHECRs) with energies up to ∼
100 EeV, soon to be realizedwith the AugerPrime Radio Detector [13]. An accurate reconstruction of 𝑋 max for inclined airshowers using RIT in addition to the measurement of the energy content of the electromagneticcascade by the same radio antennas and the mounic content by ground-particle detectors wouldprovide excellent sensitivity to the mass composition of cosmic rays [14, 15] and could thus providekey information in the quest for the origin of UHECRs .Here, we investigate whether the promising results achieved in [1] for simulations with anidealized detector (zenith-angle dependent dense antenna array, perfect time synchronisation betweenantennas) can be confirmed for air showers measured with realistically dimensioned air showerdetector arrays, i.e., coarse discretely spaced antenna arrays, as needed to instrument the requiredlarge fiducial areas, and an imperfect time synchronisation between antennas not connected bycables. The primary objective of this study is to investigate the application of RIT for inclined airshowers, where the potential is largest in terms of achievable 𝑋 max resolution and complementarityto ground-based measurements, and to formulate prerequisites for the application of RIT with sparseantenna arrays which can cover the required large fiducial areas.For interferometry, the signal arrival times and the antenna positions have to be knownvery accurately to preserve the coherence within the measured signals. In [1], the authors quotethat the timing accuracy has to be better than a quarter of the signals’ oscillation period, e.g., With “antenna” we refer to an antenna(-station) consisting of at least two orthogonally aligned antennas allowing todetermine the full 3-dimensional electric field of the incoming radio emission. In the context of air-shower simulations“antenna” refers to a location at which the radio pulse is sampled. – 2 – 𝑡 = √︃ 𝜎 𝑡 signal + ( 𝜎 (cid:174) 𝑥 antenna / 𝑐 ) < ( · 𝜈 ) − ∼ 𝜈 =
50 MHz. Furthermore theyreport that a maximum inaccuracy of 𝜎 𝑡 = 𝜎 𝑡 < .
667 ns at 50 MHz). Air shower experiments whichaim to instrument large areas rely on self-powering detector stations with wireless communication.Thus the time synchronisation between those stations, achieved with GPS clocks, is typically of theorder of a few nanoseconds ( 𝜎 𝑡 ∼ 𝜎 𝑡 (cid:46) ∼
10 cm with differential GPS surveys.Thus for measurements of the radio emission below (cid:46)
100 MHz the contribution of 𝜎 (cid:174) 𝑥 antenna to 𝜎 𝑡 canbe ignored. However, for frequencies of several hundred MHz the 𝜎 (cid:174) 𝑥 antenna can become significant.Thus verifying which coherence criterion is sufficient is crucial for the design and planning of anexperiment which aims to employ interferometric reconstructions.The investigation presented here mainly refers to the frequency band of the radio emissionfrom 30 MHz to 80 MHz. This frequency band, also used in [1], is used by most current-generationlarge scale radio detector arrays [19–21] as well as the upcoming AugerPrime Radio Detector [13].Additionally, we investigate the performance achievable with higher frequency bands, in particular50 MHz to 200 MHz as proposed for the GRAND experiment [22] and 150 MHz to 350 MHz foreven higher frequencies such as those accessible by the upcoming SKA-Low [23, 24] array or theIceCube Radio Surface Array [25]. Furthermore, we investigate how an inaccurate knowledge of theatmospheric refractivity profile affects the reconstruction. We do not consider ambient noise, i.e.,radio-frequency-interference, in our study. However, we briefly discuss this matter in Sec. 6.This article is structured as follows. First, we elaborate on the shower simulations used inthis work. In section 3 we describe the reconstruction of the shower axis and 𝑋 max with RIT.Furthermore in Sec. 3.3 the effect of inaccuracies in the knowledge of the atmospheric refractivityon the reconstruction is shown. In section 4.1 we evaluate RIT for inclined air showers with differentzenith angles measured with a 1 . We evaluate the potential of RIT using CoREAS [26] simulations. The simulations used in this workcan be divided into two sets. The first set contains a total of 1902 showers, half of which is inducedby proton and the other half by iron primaries. The showers are simulated with antennas situated ona hexagonal grid with 1 . ≡ atmosphericprofile, as listed in [27, p. 162], with a refractive index at sea level of 𝑛 = + . · − ), amagnetic field with an inclination of ∼ − ° and a strength of 0 . 𝜇 G, and an altitude of the detector– 3 – igure 1 . Top : Radio-emission footprint of a 77 . ° zenith-angle air shower coming from east measured witha dense 250 m array. The energy fluence 𝑓 , i.e., energy deposit per square meter, is color coded. The footprintexhibits the typical Cherenkov cone. Middle : Same shower measured with a sparse 1500 m (sub-)array withthe same central antenna.
Bottom : Same shower measured on a 1500 m (sub-)array with a different centralantenna. of 1400 m a.s.l.. The chosen refractive index at sea level reflects the yearly average for the locationof the Pierre Auger Observatory. The yearly fluctuations in refractivity are of the order of 7 % [28,p. 51]. We use QGSJetII-04 [29] and UrQMD [30] as high- and low-energy hadronic interactionmodels and set a thinning level of 1 × − with optimized weight limitation [31].The showers simulated with the 1 . ( 𝐸 / eV ) = . ( 𝐸 / eV ) = . ( 𝐸 / eV ) . The arrivaldirections, i.e., the azimuth 𝜙 and zenith 𝜃 angles, are uniformly randomized in 𝜙 from 0 ◦ to360 ◦ and in sin 𝜃 from 65 ◦ to 85 ◦ . The shower impact point at ground (in the following called“core”) is randomly distributed within a finite 3000 km detector array. For each shower all antennasare simulated within a maximum distance to the shower axis, beyond which the signals becomenegligible and are typically dominated by ambient radio-frequency background. The average numberof simulated antennas per shower and the maximum antenna-axis distance binned in zenith angle arelisted in Tab. 1. For (actual) measurements the antenna multiplicity is, in addition to the detector– 4 – able 1 . Average number of antennas simulated and maximum antenna-axis distance (measured perpendicularto the shower axis, i.e., in the shower plane) for the 1.5 km hexagonal grid as a function of the zenith angle in2 . ° -bins. (cid:104) 𝜃 (cid:105)/ ◦ (cid:104) 𝑛 ant (cid:105) ± 𝜎 ant ± ± ± ± ± ±
11 87 ±
21 173 ± 𝑟 maxant / m 1500 1500 1508 1822 2230 2785 3563 4707 Table 2 . The number of reconstructions 𝑛 rec performed on the dense simulations, i.e., the amount of all uniquesub-arrays for all 50 showers, and average number of antennas on each sub-array (cid:104) 𝑛 ant (cid:105) for the different arrayspacings. spacing / m 250 500 750 1000 1250 1500 𝑛 rec
50 200 450 800 1250 1800 (cid:104) 𝑛 ant (cid:105) dense grid, which can be divided in several sub-arrays with larger antennaspacings, are suitable. Since the computational cost for each shower scales almost linearly withthe number of simulated pulses we need to limit our phase space of densely sampled, simulatedshowers. Thus we simulate 50 proton showers with only one energy log ( 𝐸 / eV ) = .
4, one zenithangle 𝜃 = . ◦ and two azimuth angles 𝜙 = ◦ (arriving from geomagnetic east) and 𝜙 = ◦ (arriving from north of east), for each of which we simulate 25 showers. For a hexagonal arraywhich is invariant for rotations of 60 ° , showers from 0 ° and 30 ° cover the two extreme cases ofa shower falling into the array exactly parallel to a line of antennas and with the largest possibleangle between two lines of antennas. Pulses are simulated on a grid with 250 m spacing and amaximum axis distance of 2235 . (cid:46) . To study the reconstructionperformance for different array spacings we define various sub-arrays. The following arrays areinvestigated: 250 m, 500 m, 750 m, 1000 m, 1250 m and 1500 m. For each spacing (except 250 m)several unique sub-arrays can be defined, each of them corresponding to a different (relative) coreposition for a given simulation. Thus, for example, one single shower can be reconstructed on 36unique sub-arrays with a spacing of 1500 m. In Fig. 1 an example shower measured with the full250 m grid (top) and two different sub-arrays with a spacing of 1500 m (middle, bottom) is shown.In Tab. 2 the number of all unique sub-arrays for all 50 showers and the average number of antennason these sub-arrays for each spacing are summarized. Unlike for the simulation set with the 1 . – 5 – igure 2 . Left : Cross section of the coherent energy fluence profile of a 1 EeV, 77 . ° proton shower sampledwith 1335 antennas on a 250 m grid. The cross section shows the longitudinal profile along the y-axis ing/cm and the lateral profile in (cid:174) 𝑣 × (cid:174) 𝐵 -direction along the x-axis. The coherent energy fluence is color-coded.The vertical black line indicates the shower axis. The grey lines illustrate 2-dimensional lateral cross-sections(dashed lines) along an initial guessed axis (solid grey line) used for the shower axis reconstruction (theshown line does not reflect the typical deviation of a guessed axis with a Gaussian resolution in 𝜙 and 𝜃 of0 . ° from the true shower axis). Right : Longitudinal profile of the coherent energy fluence 𝑓 𝐵 𝑗 for the sameshower. First, to narrow down the position of the maximum, the longitudinal profile is sampled with a coarse100 g cm − sampling (large blue circles). Around the found maximum, in a 200 g cm − window, the profile isthen sampled more precisely with 10 g cm − steps (small orange circles). In this window, a Gaussian curve isfitted to the profile (black curve). The vertical blue line shows the found maximum as determined from thefitted Gaussian parameters. In this section, we describe the reconstruction of the shower axis and the depth of the shower maximum 𝑋 max with RIT. The algorithms, developed in [1], make use of 3-dimensional interferometric mapsproviding information about the longitudinal development of air showers. From these maps thecosmic-ray properties, in particular the arrival direction and depth of the shower maximum 𝑋 max ,can be inferred. The algorithms described below are adapted from [1], however, their actualimplementation is independent and has, in parts, changed.RIT exploits the coherence in the radio emission from air showers and one searches for animaginary point source for which the coherent signal becomes maximal. The time-dependentcoherent (beam-formed) signal 𝐵 𝑗 ( 𝑡 ) originating at an arbitrary location in the atmosphere (cid:174) 𝑗 iscalculated by the sum over all time-shifted antenna signals 𝑆 𝑖 ( 𝑡 − Δ 𝑖, 𝑗 ) at positions (cid:174) 𝑖𝐵 𝑗 ( 𝑡 ) = 𝑛 ant ∑︁ 𝑖 𝑆 𝑖 ( 𝑡 − Δ 𝑖, 𝑗 ) . (3.1)The time shift between an antenna location (cid:174) 𝑖 and the source location (cid:174) 𝑗 is Δ 𝑖, 𝑗 = 𝑑 𝑖, 𝑗 · 𝑛 𝑖, 𝑗 𝑐 (3.2)– 6 –ith the geometrical distance 𝑑 𝑖, 𝑗 and effective (averaged) refractive index 𝑛 𝑖, 𝑗 between the positions (cid:174) 𝑖 and (cid:174) 𝑗 , and the vacuum speed of light 𝑐 . That means, Δ 𝑖, 𝑗 corresponds to the light propagation timebetween positions (cid:174) 𝑖 and (cid:174) 𝑗 . To calculate the effective refractive indices between source locations (cid:174) 𝑗 and antenna positions (cid:174) 𝑖 a model for the refractivity 𝑁 ( ℎ ) = 𝑛 ( ℎ ) − 𝑁 ( ℎ ) = 𝑁 ( ) · 𝜌 ( ℎ )/ 𝜌 ( ) . (3.3)In [32] it is shown that this approximation is adequate for the frequency band of 30 MHz to 80 MHz,for higher frequencies the Global Data Assimilation System (GDAS) can be used to refine therefractivity model by then also including the influence of humidity. The practical calculation of theeffective refractivity between two positions, 𝑁 𝑖, 𝑗 , which cannot be calculated analytically in a curvedatmosphere, is explained in appendix 8.1. The calculation of the light propagation time (using theeffective refractivity) along straight lines corresponds to the algorithm adopted in CoREAS. Innature, the emission between sources and observers propagates on slightly bent trajectories dueto refraction in the atmospheric refractive index gradient. In [33] we found that the calculationon straight lines reproduces the relative propagation times between two different sources in theatmosphere better than within 0 . 𝐵 𝑗 ( 𝑡 ) , the electric field values of the time-shifted signals 𝑆 𝑖 ( 𝑡 − Δ 𝑖, 𝑗 ) are linearlyinterpolated to fit the finite time binning Δ 𝑡 of 𝐵 𝑗 ( 𝑡 ) . For each trace 𝐵 𝑗 ( 𝑡 ) we determine atime-independent signal, namely the sum over the squared amplitudes in a 100 ns signal windowaround the peak amplitude 𝑓 𝐵 𝑗 = 𝜖 𝑐 Δ 𝑡 𝑡 peak + ∑︁ 𝑡 peak − 𝐵 𝑗 ( 𝑡 ) (3.4)where 𝜖 is the vacuum permittivity and 𝑐 the speed of light in vacuum. The peak amplitude andthe peak time 𝑡 peak in 𝐵 𝑗 ( 𝑡 ) are determined from the maximum of the absolute Hilbert envelope of 𝐵 𝑗 ( 𝑡 ) . The quantity 𝑓 𝐵 𝑗 can be understood as the coherent energy fluence received by the array ofobservers 𝑖 from a given location (cid:174) 𝑗 .Eqs. (3.1), (3.2), and (3.4) allow us now to calculate the coherent energy fluence received fromany position in the atmosphere. In Fig. 2 ( Left ) a cross section of the coherent energy fluence froman example shower sampled at 1335 antenna locations is shown. The longitudinal profile along theshower axis (vertical black line) is expressed in g cm − (y-axis) while the lateral profile is shownperpendicular to the shower axis along the (cid:174) 𝑣 × (cid:174) 𝐵 -direction (x-axis) in meters. It is apparent that theprofile of the coherent energy fluence correlates with the particle cascade of the air shower, i.e., 𝑓 𝐵 is strongest around the shower axis and exhibits a maximum. It has been shown that this maximum,defined as 𝑋 RIT , correlates linearly with the shower maximum of the particle cascade [1]. Thus RITallows to reconstruct the shower properties, e.g., the depth of the shower maximum and shower axis. A linear interpolation is not strictly physically correct. Application of a phase gradient to the Fourier spectrum oradequate up-sampling of 𝑆 𝑖 would be more physically motivated. However, linear interpolation is computationally moreefficient and we validated that the reconstruction accuracy is independent of this procedure for Δ 𝑡 = Δ 𝑡 = .
33 ns for frequencies up to 350 MHz. – 7 – igure 3 . Histograms of the opening angle distribution between true and reconstructed arrival directionfor the different antenna spacings (colors) for showers with 𝜙 = ◦ resp. 𝜙 = ◦ . The legend shows theresolution of the arrival direction reconstruction in terms of the 68 %-quantile for all shower and all showerwith 𝜙 = ◦ resp. 𝜙 = ◦ . As in [1], only the signal in the (cid:174) 𝑣 × (cid:174) 𝐵 polarisation ( (cid:174) 𝑣 : direction of the primary particle trajectory,i.e., shower axis, (cid:174) 𝐵 : direction of the Earth’s magnetic field), which is obtained by rotating the Electricfield vector simulated in the North-South, West-East, Vertical polarisations using the true arrivaldirection, is used for reconstruction ( ≡ 𝑆 𝑖 ( 𝑡 ) ). It seems natural to separate the radio emission basedon its emission mechanisms, i.e., separate between geomagnetic and charge-excess emission, as anyphase-shift in the signals between both mechanisms would reduce the signals coherence. Such phaseshifts correspond to a small degree of circular polarization observed both in simulations and data,see [9]. In inclined air showers the geomagnetic emission, which constitutes most of the signal in the (cid:174) 𝑣 × (cid:174) 𝐵 polarisation, is dominant while the signal in the (cid:174) 𝑣 × (cid:174) 𝑣 × (cid:174) 𝐵 polarisation is completely comprisedby the sub-dominant charge-excess emission. In fact, determining the longitudinal profiles 𝑓 𝐵 ( 𝑋 ) with signals in the (cid:174) 𝑣 × (cid:174) 𝑣 × (cid:174) 𝐵 polarisation yields no well-defined maxima which can be correlated tothe depth of the showers. The shower axis, i.e., the extrapolated trajectory of the primary particle, is reconstructed with RITby searching for an axis along which the longitudinal profile of the coherent energy fluence ismaximal. For this propose, the lateral profile of the coherent radio emission, i.e., the cross-section of 𝑓 𝐵 𝑗 ( 𝑋 = const ) , is sampled at several depths along the shower’s development. For each cross-sectionthe location of its maximum is determined and interpreted as its intersection with the shower– 8 –xis. Given these intersections, a straight line is fitted minimizing the distance between line andintersections, weighted by the signal strength of each maximum.Each cross-section is sampled in a plane perpendicular to an initial (guessed) axis which isdetermined given the true arrival direction, but smeared in zenith and azimuth angle with a Gaussianresolution of 0 . ° each, and an intersection point at ground given by the intersection of the trueMonte-Carlo (MC) shower axis smeared in a perpendicular plane with a Gaussian resolution of100 m. This accommodates for the imperfect knowledge of the shower axis from a traditionalreconstruction as starting point for a RIT reconstruction under practical circumstances.The following procedure is applied to find the maximum in each lateral cross-section at depthsof 500, 600, 700, 800, 900, 1000 and 1100 g cm − : In a first iteration, the maximum is searchedon a quadratic grid which is characterized by its overall size and grid spacing. Here we chose agrid spacing of 60 m. The grid covers an area of 1 km and is set such that the location of the MCshower axis is within the search grid. This is due to performance reasons, not all MC shower axeswould be contained in a 1 km -grid around the initial guessed axis under the starting conditionsmentioned above. The 68 % quantile of the distance between MC and guessed shower axis foran MC zenith angle of 77 . − is 675 m. However the area is sufficientlydimensioned to realistically model interferometric maps containing grating lobes, i.e., local maxima.For experimental measurements one has to ensure to make the search region sufficiently scaled, ofcourse at the expense of computational effort. In a second iteration, the cross-section is sampled on arefined quadratic grid around the previously found maximum, i.e., zoomed-in around the previouslyfound maximum. This process is repeated until the grid spacing becomes smaller than 0 . ° .In Fig. 3 the opening angle distribution between true and reconstructed arrival direction forthe dense simulations reconstructed on arrays with different antenna spacings with perfect timesynchronisation is shown. The histogram is separated between showers arriving with azimuth anglesof 𝜙 = ◦ ( Left ) and 𝜙 = ◦ (Right). The resolution in terms of the 68 %-quantile is shown in therespective legend. The overall accuracy, especially for antenna spacings ≤ . ◦ for almost all configurations. For larger antenna spacings a bigger differencebetween showers from the two different incoming directions is evident. While the worsening of theresolution as function of the antenna spacing for showers from 𝜙 = ◦ is moderate and just becomessignificant for the 1500 m spacing, the degeneration for showers from 𝜙 = ◦ is much more dramatic.For those showers a footprint where all high-signal antennas are aligned on a straight line parallel tothe shower axis projected on the ground is likely as the antenna grid gets too coarse to sample theCherenkov cone along the whole plane (cf. middle panel of Fig. 1). Inferring the correct arrivaldirection is more difficult for such geometries. To reconstruct the depth of the shower maximum, we determine the maximum 𝑋 RIT of the longitudinalprofile of the coherent signal 𝑓 𝐵 𝑗 ( 𝑋 ) along the (Monte-Carlo or reconstructed) shower axis. Aprofile of 𝑓 𝐵 𝑗 ( 𝑋 ) along the MC shower axis as function of the slant depth 𝑋 is shown in Fig. 2( Right ). To find 𝑋 RIT we employ the following algorithm: We sample the longitudinal profile 𝑓 𝐵 𝑗 ( 𝑋 ) along the shower axis in steps of 100 g cm − between 500 and 1000 g cm − . If the maximum isfound at an edge, the sampling range is dynamically extended. Once the maximum is well-confined,a 200 g cm − window around the found maximum is sampled with a refined step size of 10 g cm − .– 9 – igure 4 . Top : Reconstructed 𝑋 RIT as a function of 𝑋 MCmax for the dense simulations. The black line indicatesthe calibration curve according to Eq. (3.5). Reconstruction along the MC shower axis with perfect timesynchronisation between the antennas. The different colors refer to reconstructions with (sub-)arrays ofdifferent spacings (transparency increases with number of reconstructions). The legend illustrates the bias andresolution for different array spacing.
Bottom : Residuals between reconstructed and true depth of the showermaximum 𝑋 max . 𝑋 RIT is then determined by the maximum of a Gaussian curve fitted to this 200 g cm − window (cf.Fig. 2, Right ).In Fig. 4 (top) the reconstruction of 𝑋 RIT for all dense simulations with a zenith angle of 𝜃 = . ◦ and the different aforementioned array spacings (color coded) as a function of the trueshower depth 𝑋 MCmax is shown. The reconstruction is performed with a perfect time synchronisationbetween the different antennas, i.e., the signal arrival times are exactly known, and along the MCshower axis. A good, linear correlation is found between the reconstructed 𝑋 RIT and 𝑋 MCmax . Thuswith a linear equation the shower maximum can be reconstructed as a function of 𝑋 RIT : 𝑋 max ( 𝑋 RIT ) = . · 𝑋 RIT + .
15 g cm − . (3.5)The resulting residuals for the different spacings are shown on the bottom panel of the same figure(bias and resolutions of this residual in the legend of the top panel). It can be seen that regardless ofthe array spacing and thus the number of pulses used in the reconstruction (antenna multiplicity) anaccurate reconstruction is achieved.In Fig. 5 the reconstruction of 𝑋 max for the simulations on the 1 . 𝑋 max and 𝑋 RIT for thehere considered zenith angles range. A fit to showers with 𝜃 ≥ ◦ yields the following calibration– 10 – igure 5 . Left : Reconstructed 𝑋 max as a function of 𝑋 MCmax for the simulations on the 1.5 km grid along theMC shower axis with perfect time synchronisation between the antennas. The reconstruction of 𝑋 max is afunction of 𝑋 RIT and the zenith angle 𝜃 according to Eq. (3.6). The black dashed line indicates the identity,the color code shows the zenith angle. Right : Residuals between reconstructed and true depth of the showermaximum 𝑋 max as function of the zenith angle. Bottom panel shows the profile, i.e., mean 𝜇 and standarddeviation 𝜎 of the above residual. function: 𝑋 max ( 𝑋 RIT , 𝜃 ) = . · 𝑋 RIT + (cid:18) . − 𝜃 − . ◦ . ◦ (cid:19) g cm − . (3.6)In Fig. 5 ( Left ) the reconstructed 𝑋 max as function of the true 𝑋 MCmax is shown. The comparisonexhibits a significant scatter, only for showers with higher zenith angles (color coded) is a goodcorrelation achieved. The residual of the reconstructed 𝑋 max as function of the zenith angle and itsprofile (mean and standard deviation binned in 2 . ° zenith-angle bins) is also shown ( Right ). It isapparent that the reconstruction accuracy strongly depends on the zenith angle. The dominant effecthere is the insufficient antenna multiplicity for lower zenith angles (cf. Tab. 1). The dependence ofthe reconstruction accuracy on the antenna multiplicity is investigated in more detail in section 4.Comparing Eqs. (3.5) and (3.6), evaluated for 𝜃 = . ◦ , reveals no significant deviationbetween each other. Furthermore no significant bias between the 𝑋 max reconstructions with differentantenna spacings is evident (cf. legend in Fig. 4). Hence it seems that the calibration between 𝑋 max and 𝑋 RIT is independent of the antenna spacing and the different detector layouts covered in thiswork. The calibration found here is also in good agreement with [1] which found an average depth 𝑋 RIT ∼
615 g cm − for a true depth of the shower maximum of 𝑋 max ∼
700 g cm − for showers with 𝜃 = ◦ . In this section, we evaluate the impact of an inaccurately known atmospheric refractive index profileon the reconstruction. To that end, we reconstruct the showers with different atmospheric refractivity– 11 – igure 6 . Deviation in 𝑋 RIT reconstructed with different refractive index profiles along the MC shower axis.A deviation of (cid:46) − which cor response to (cid:46)
150 m is small compared to the absolute distance ∼
55 km.Here, all 50 dense showers are reconstructed once for each array spacing. profiles than used in the CoREAS simulations. We use the refractivity profiles ( ≡ refractivity at sealevel & atmospheric density profile) as adequate for the site of the Pierre Auger Observatory for themonths of February and June for reconstruction, while the October profile was used for the CoREASsimulations. These two months represent the extrema in the yearly fluctuation of the refractivity atground at the location of the Pierre Auger Observatory (for which the simulated October atmosphereresembles a good yearly average) [28, Fig. 3.21]. The yearly fluctuation is on the order of 7 % andthus larger than for the locations of other radio air-shower experiments such as LOFAR or Tunka-Rexwith a yearly fluctuation of 4 % and 3 % [28, p. 51], respectively. Thus using the refractivity profilesfor February and June implies minimal knowledge of the true refractivity in the atmosphere.Using a mismatching atmospheric density profile for reconstruction will yield a wrong atmo-spheric depth even if the point of origin of the maximally coherent emission is correctly determined.Thus, the atmospheric depth of a maximum, reconstructed with an inaccurate refractive index profile,is determined using the correct atmospheric density profile. The deviation in 𝑋 RIT between thereconstruction with different refractive index profiles shown in Fig. 6 is (cid:46) − .The uncertainty due to an inaccurate knowledge of the atmospheric density profile is identicalto the corresponding uncertainty of 𝑋 max measurements with the fluorescence technique which is onthe order of 2 g cm − to 4 g cm − [34] for higher energies at the Pierre Auger Observatory.For the following investigation the (correct) simulated atmosphere density profile and refractiveindex at sea level is used for the reconstruction.– 12 – Interferometric reconstruction of the depth of the shower maximum under realisticconditions
Having evaluated the achievable performance under idealized conditions, largely confirming theresults reported in [1], in the following sections we will evaluate the 𝑋 max reconstruction with RITfor a more practical scenario, i.e., with imperfect time synchronisation and along the reconstructedshower axis. Before we examine the interferometric reconstruction for simulations with varying detector layoutswe evaluate the technique on simulations with the finite 1 . 𝑋 max .In Fig. 7 the resolution of the 𝑋 max reconstruction via Eq. (3.6) binned as function of the antennamultiplicity is shown. The reconstructions along the MC and reconstructed shower axes are shownwith solid and dashed lines, respectively. Different Gaussian time jitters are shown in different colorsand markers. The horizontal error bars indicate the bin size, the vertical bars correspond to thestatistical fluctuation of the resolution determined via a bootstrapping procedure. The reconstructionquality depends on both the antenna multiplicity and the accuracy of the time synchronisation. Evenwith a perfect time synchronisation, a minimum number of antennas (cid:38)
12 is required to keep theresolution below 40 g cm − . An imperfect time synchronisation limits the achievable accuracy. Veryaccurate results ( 𝜎 𝑋 max <
20 g cm − ) are only achieved for a time jitter of 1 ns or less and (cid:38) 𝜎 𝑡 = (cid:104) 𝑛 ant (cid:105) ∼
100 and 3 in the last bin).In section 6 the results obtained here are discussed and compared to other measurements of thedepth of the shower maximum.
Now we study the effect of imperfect time synchronisation between antennas on the reconstructionof showers measured with different array spacings / antenna multiplicities. To that end, we repeatthe reconstruction of the 50 simulated showers several times on various different sub-arrays afterintroducing random Gaussian time jitters mimicking an inaccurate time synchronisation betweenthe antennas. Figure 8 shows the resolution in 𝑋 max as a function of the antenna spacing and fordifferent time jitters. The average number of antennas (cid:104) 𝑛 ant (cid:105) per event and spacing is shown onthe top x-axis. The figure demonstrates that again the resolution worsens in the presence of a timejitter. This deterioration is amplified for showers reconstructed with a low antenna multiplicity.With a time jitter of 3 ns the reconstruction on a very dense array with > igure 7 . Resolution of the 𝑋 max reconstruction for the simulations on the 1.5 km grid along the true andreconstructed shower axes (solid and dashed lines, respectively), and for different Gaussian time jitters. Theresolution is binned as a function of the antenna multiplicity, the horizontal error bars indicate the bin size,the vertical bars correspond to the statistical fluctuation of the resolution. The 𝑋 max reconstruction along thereconstructed axes can contain a few outliers where the axis reconstruction was inaccurate, causing the visiblefluctuations in the 𝑋 max resolution. very accurate with a resolution of (cid:46)
10 g cm − . However, when the antenna multiplicity is (cid:46) (cid:38)
40 g cm − . Reconstructing 𝑋 max along an imperfectly reconstructed axis seems to have nosignificant implications for data taken with a Gaussian time jitter of 0 ns to 2 ns and only a little effectfor the 3 ns time jitter. In contrast to the simulations with the 1 . array and thus their footprints are eventually notevenly sampled, the showers with the dense 250 m-spaced array are always evenly sampled.Incorrect time synchronisation between antennas also affects the arrival direction reconstruction.In case of a 3 ns Gaussian time jitter the resolution in the direction reconstruction worsens by afactor of ∼ 𝜙 = ◦ and30 ° (cf. Fig. 3) decreases.It is important to stress that the primary factor governing the 𝑋 max resolution, in case of imperfecttime synchronisation between antennas, is the antenna multiplicity. The average number of antennasfor each array spacing as listed in table 2 refers to an instrumented area of 72 . (cf. Fig. 1).Smaller but denser arrays still need to accommodate a sufficiently high number of antennas and/orvery good time synchronisation to allow accurate reconstructions. In addition, our tests have shownthat a complete and symmetric sampling of the radio-emission footprint, i.e., inside, on top, andoutside the Cherenkov cone, is needed to ensure accurate reconstruction.In section 6 we use an analytic description of the radio-emission induced area as function of thezenith angle to generalize the results acquire here with showers with 𝜃 = . ◦ and an instrumented– 14 – igure 8 . Reconstruction resolution in 𝑋 max of the 50 densely sampled showers with a zenith angle of 77 . ° .Resolution is shown for different time jitter scenarios (different colors & markers) and along the MC showeraxis (solid line) or reconstructed axis (dashed line) as a function of the antenna spacing. The average numberof antennas per spacing is given on the top x-axis as reference. Figure 9 . Reconstruction resolution in 𝑋 max for the 30 MHz to 80 MHz, 50 MHz to 200 MHz, and 150 MHzto 350 MHz frequency bands (different colors) along the MC or reconstructed shower axis (solid or dashedline) for perfect time synchronisation as a function of the antenna spacing. area of 72 . to lower zenith angles and smaller and denser arrays.– 15 – able 3 . Parameters of Eq. (3.5), i.e., 𝑋 max = 𝑎 · 𝑋 RIT + 𝑏 for the different frequency bands. a b30 to 80 MHz 1.029 76.15 g cm −
50 to 200 MHz 1.027 76.97 g cm −
150 to 350 MHz 1.024 92.91 g cm − Many next-generation radio-detection experiments aim to observe extensive air showers with broaderfrequency bands and at higher frequencies. Here, we test the interferometric reconstruction of 𝑋 max for two additional frequency bands: 50 MHz to 200 MHz and 150 MHz to 350 MHz. Applying theinterferometric algorithm to data recorded at higher frequencies requires more stringent coherencecriteria and thus a more accurate time synchronisation between antennas. The considered time jitterscenarios do not reflect equivalent phase-accuracy across the frequency bands but were chosen froma practical point of view, i.e., what time synchronisation accuracy an experiment has to achieve toemploy RIT for higher frequencies. Furthermore, we found that the 3-dimensional profile of thecoherent signal 𝑓 𝐵 𝑗 around the shower axis is increasingly narrow for higher frequencies. Hence theresolution with which the lateral cross-sections are sampled to infer the shower axis needs to berefined. For the axis reconstruction of showers recorded in the frequency band from 150 MHz to350 MHz an overall size of 0 .
16 km and a search grid spacing of 20 m was used.In Fig. 9 the 𝑋 max resolution with perfect time synchronisation and along the MC andreconstructed shower axes for the different frequency bands (color coded) are compared. For thehigher frequencies the reconstruction accuracy (along the MC shower axis) slightly decreases forsparser antenna arrays (150 MHz to 350 MHz). This is even more prominent for the reconstructionalong the reconstructed axis (50 MHz to 200 MHz and 150 MHz to 350 MHz). One reason for thiscould be an insufficient sampling of the Cherenkov cone with these sparse arrays. As the Cherenkovcone itself is more dominating but also more narrow for showers measured with higher frequencies,the reconstruction accuracy depends more strongly on the antenna spacing. Furthermore, an offsetin 𝑋 RIT between the different frequency bands is found, e.g., the 𝑋 RIT reconstructed for 150 MHzto 350 MHz is ∼ . − smaller as for 30 MHz to 80 MHz. This offset has been taken intoaccount when determining 𝑋 max . For this purpose we repeated the parametrisation of Eq. (3.5) forthe different frequency band, see Table 3.Figure 10 shows the achieved 𝑋 max resolution for the two higher frequency bands ( Left : 50 MHzto 200 MHz,
Right : 150 MHz to 350 MHz) and different Gaussian time jitter scenarios. As expected,for the higher frequency bands already more modest time jitters, e.g., 2 ns for 50 MHz to 200 MHzand 1 ns for 150 MHz to 350 MHz, worsen the results significantly. The antenna multiplicity has tohigher than ∼
70 (100) and the time synchronisation better than 𝜎 𝑡 = 𝜎 𝑋 max <
40 g cm − .As mentioned before, the axis reconstruction for 150 MHz to 350 MHz requires a finer samplingof the lateral cross-sections. With the refined search-grid spacing of 20 m a similar resolution of 𝑋 max , reconstructed along the reconstructed axis, compared to the lower frequency bands with asearch-grid spacing of 60 m is observed. However, the overall grid size used to reconstruct the– 16 – igure 10 . Reconstruction resolution in 𝑋 max for the 50 MHz to 200 MHz ( Left ) and 150 MHz to 350 MHz(
Right ) frequency bands and different time jitters (different colors and markers) along the MC or reconstructedshower axis (solid or dashed line) as a function of the antenna spacing. shower maximum for data recorded with 150 MHz to 350 MHz does not reflect the reconstructionunder practical circumstances as the searched area is too small to reliably contain the true maximumfor an axis (starting point) which is known with a accuracy of 0 . ° in zenith and azimuth each.A, more sophisticated, gradient-descent based algorithm could reduce the computing timesignificantly and would allow the shower axis reconstruction under practical circumstances also forhigher frequency bands. Such an algorithm has to be robust against grating lobes, i.e., local maximain the interferometric maps.In [1], in addition to an algorithm similar to the one described above, the shower axis isrefined by maximizing the integrated longitudinal profile along the axis. Also in LOPES the arrivaldirection is inferred by a two-folded approach, first applying a raster-search algorithm and upon thisa gradient-descent algorithm. Here, we discuss our results and also mention a few aspects which could not be studied in detailwithin the scope of this work.The study presented here for different detector layouts, especially those with an antenna spacingof < 1 . 𝜃 = . ◦ . In [1], a modest improvementof the reconstruction accuracy with increasing zenith angle, i.e., with increasing distance betweenobserver and source region, is found for simulations with a constant antenna multiplicity. The drasticimprovement in resolution with increasing zenith angle found in this work (cf. Fig. 5) cannot beexclusively explained by this. In fact, the antenna multiplicity is identified as the crucial factor foran accurate reconstruction. Showers with a lower zenith angle illuminate smaller areas at groundand thus the antenna multiplicity decreases for showers measured with a constantly-spaced array.Assuming that the size 𝐴 of the radio-emission footprint at ground scales with the radius of theCherenkov cone 𝑟 che (defined in the shower plane) yields the following relation: 𝐴 ∼ 𝜋𝑟 cos 𝜃 . (6.1)– 17 –f an air shower is approximated as a point source moving with the speed of light, we can approximatethe Cherenkov radius with the radius of a cone with its apex at the shower maximum and anopening angle equal the Cherenkov angle at the location defined by the refractive index 𝑛 ( ℎ ) ,cos ( 𝛼 che ) = / 𝑛 ( ℎ ) : 𝑟 che = tan (cid:20) cos − (cid:18) 𝑛 ( ℎ max ( 𝜃, 𝑋 max )) (cid:19) (cid:21) · 𝑑 max ( 𝜃, 𝑋 max , ℎ obs ) (6.2)with the height of the shower maximum above sea level ℎ max , the distance along the shower axisbetween ground and shower maximum 𝑑 max and the altitude of the observation plane ℎ obs . Theantenna multiplicity is proportional to the footprint area 𝐴 . To instrument a given area with a certainnumber of antennas the antenna spacing Δ ant has to satisfy the relation: √ Δ ant ∼ 𝐴 . Figure 11 showsthe antenna spacing as a function of the zenith angle necessary to satisfy the antenna multiplicity (cid:104) 𝑛 ant (cid:105) =
84, i.e., the antenna multiplicity of the reference point at 𝜃 = . ◦ (black marker) forshowers measured with a 1000 m hexagonal grid (cf. Tab. 2). Ignoring any additional zenith-anglerelated effects, this curve indicates the necessary antenna spacing to achieve a reconstruction asaccurate as for the reference, i.e., 𝜎 𝑋 max ( 𝜎 𝑡 = ) =
16 g cm − or 𝜎 𝑋 max ( 𝜎 𝑡 = ) =
29 g cm − (cf.Figs. 7, 8). The different colors refer to different observation heights and the shaded areas correspondto 𝑋 max values ranging from 550 g cm − to 950 g cm − with a nominal value of 750 g cm − (solidlines). These values resemble the mean and range of 𝑋 max values for a cosmic ray composition ofhalf proton and half iron primaries with energies around 10 EeV. It is apparent that an accuratereconstruction for more vertical showers with 𝜃 < ◦ is only achievable with antenna spacings below100 m. For showers with zenith angles below 25 ° the antenna spacing cannot be larger than ∼ tensof meters. Those showers, when measured at an altitude of 1000 m a.s.l., can reach ground beforedeveloping the full maximum (this happens with a depth of the maximum of >
850 g cm − ) andthus no lower limit can be calculated. This emphasizes that for vertical air showers the observationaltitude matters. This is further underlined by the finding in [1] that the accuracy in 𝑋 max deteriorateswhen the distance between observer and source becomes smaller. The area associated to a givenantenna spacing and (cid:104) 𝑛 ant (cid:105) =
84 is shown on the second (right) y-axis.Furthermore, in this study the issue of triggering the readout of the radio signals, be it based onradio signals or measured particles, has not been considered. The interferometric reconstructionprofits also from low signals and thus a readout of all antennas for a given event is optimal (assimulated in this study), even if no measurable radio pulse or particles are present. However mostcosmic-ray experiments employ a trigger based on the signal strength per detector station to reducethe amount of data recorded. This might lead to a reduction of recorded radio pulses and thus limitthe accuracy of interferometric measurements. Experiments with an accompanying particle detectorcan profit from a lower trigger threshold and thus recorded more radio pulses for vertical showers. Forair showers with zenith angles beyond 75 ° to 80 ° , however, the size of the radio-emission footprinteventually exceeds the size of the particle footprint [12], and a trigger relying on information ofparticle detectors alone will limit the amount of radio pulses recorded. For instance, the AugerPrimeRadio Detector will only record radio pulses of antennas for which the water-Cherenkov detector(WCD) beneath has triggered. Thus the number of pulses recorded is governed by the particlefootprint, i.e., the footprint for which the water-Cherenkov detectors will trigger, and as such is afunction of the primary energy and the zenith angle. To trigger more than 12 WCD, a minimum– 18 – igure 11 . Antenna spacing required to achieve an antenna multiplicity of ∼
84. This is the meanantenna multiplicity for showers with 𝜃 = 77 . ° measured with a 1000 m hexagonal grid (black dot) overan instrumented area of 72 . . The second y-axis shows the area associated to a given antenna spacingand the aforementioned antenna multiplicity. For these showers an 𝑋 max resolution below 20 g cm − whenmeasured with a 1 ns Gaussian time jitter in the 30 MHz to 80 MHz band is achieved (cf. Fig. 8). The differentcolors refer to different observation heights and the shaded areas correspond to 𝑋 max values ranging from550 g cm − to 950 g cm − with a nominal value of 750 g cm − (solid lines). The curves are calculated usingEqs. (6.1) and (6.2). For the calculation of the radius of the Cherenkov cone the US standard model of theatmospheric density profile and a refractivity at ground of 𝑁 = · − are used. For showers with a zenithangles below 25 ° measured at an observation altitude of 1000 m the shower maximum can lie underground,thus no lower limit can be calculated. zenith angle of (cid:38) ° or energy of (cid:38)
10 EeV has to be reached [35]. More than 50 stations are almostnever triggered. Given these limitations in the triggering for the AugerPrime Radio Detector, theapplication of a RIT reconstruction unfortunately does not seem very promising even if the timesynchronisation can be improved to 𝜎 𝑡 ∼ (cid:174) 𝑣 × (cid:174) 𝐵 polarisation, is accessible from the experimentalmeasurements.Besides detector effects, ambient noise is a crucial aspect for the detection of air showers withradio antennas. However, the ambient noise conditions can change dramatically between differentlocations around the Earth, their significance depends on the observed range of cosmic-ray energiesand the frequency-band of choice, and their modeling taking into account different contributions,e.g., narrow- and broadband radio-frequency-interference, is not straightforward. Furthermore it isanticipated that the effect of ambient noise is attenuated for interferometric measurements scalingwith the square root of the number of antennas [1]. Dedicated studies for specific experiments areneeded to determine the impact of noise at their specific location. Such a study can also determine if– 19 –he interferometric detection threshold can be lowered when measuring in a higher frequency band[36].Judging the required 𝑋 max resolution to be achieved with a large sparse antenna arrays is acomplex question as it depends on several factors such as the scientific objective, e.g., measuringthe average mass composition or aiming for a light-heavy particle discrimination, the availablestatistics, and the astrophysical scenario, i.e., the actual mass composition of cosmic rays. To simplify,we compare the achievable 𝑋 max resolution with RIT to different experimental results. Recently,the Pierre Auger Collaboration has demonstrated that an accurate 𝑋 max reconstruction with the1 . − andimproves to 25 g cm − for energies above 20 EeV. The resolution achieved by the Auger FluorescenceDetector is 25 g cm − (15 g cm − ) for energies above 1 EeV (10 EeV) [34]. LOFAR, a radio airshower experiment, measures vertical showers with hundreds of antennas in the energy range from10 to 10 eV with a typical accuracy of 17 g cm − [38]. Tunka-Rex, another radio air showerexperiment, measures 𝑋 max with a low, typical multiplicity of 7 antennas with an accuracy of25 g cm − [39]. However, both results achieved with these radio experiments rely on the extensiveuse of very time-consuming and computing-intensive Monte-Carlo simulations and are not applicablefor larger antenna arrays (with higher event statistics) or more inclined air showers.The results obtained in this work show that the application of RIT with large, sparse antennaarrays relying on wireless communication is very challenging. Even with specialized hardwarewhich improves the time synchronisation to ∼ (cid:38)
50 has to be achieved to obtain competitive results. These requirements will likelynot be met by existing or currently planned experiments such as the Pierre Auger Observatory orGRAND.More suitable for the application of RIT seem smaller ultra-dense antenna arrays with cabledcommunication such as the Square Kilometer Array SKA-Low [23], which in fact is designed as aninterferometer and thus will meet the required timing accuracy for interferometric analyses.
This study explores the potential for interferometric measurements of the depth of maximum ofextensive air showers 𝑋 max with large arrays of radio antennas under realistic conditions. It hasbeen shown that in addition to a very good time synchronisation of 1 ns for the frequency band of30 MHz to 80 MHz, also a sufficiently large number of antennas per shower ( (cid:38) 𝑋 max . Given the size of the radio-emission footprint as a function of thezenith angle, this constrains the maximal suitable array spacing to <
100 m for vertical showers withzenith angles 𝜃 < ◦ and to a few hundred meters for showers with zenith angles 𝜃 (cid:46) ◦ . Only forhigher zenith angles, arrays with an antenna spacing of 1000 m or larger accommodate a sufficientantenna multiplicity. However, any kind of trigger based on the signal strength per detector stationwill reduce the amount of recorded radio pulses. Thus it seems very challenging to accommodatesuch requirements for (existing) air shower arrays with spacings (cid:38) 𝑋 max reconstructionis achieved when the geometry and signals’ arrival times are exactly known. Thus no advantage isfound when applying the interferometric reconstruction to data recorded with higher frequencies.Experiments which facilitate a large number of antennas combined with a very accuratetime synchronisation such as the Square Kilometer Array SKA have great potential to exploitinterferometric measurements of 𝑋 max . If combined with a muon detector, this approach couldyield very valuable information to study the physics of extensive air showers and their hadronicinteractions, as well as the mass composition of cosmic rays, with unprecedented detail. Acknowledgements
We are grateful to F.G. Schröder for his valuable comments on our manuscript. Felix Schlüter issupported by the Helmholtz International Research School for Astroparticle Physics and EnablingTechnologies (HIRSAP) (grant number HIRS-0009). Simulations for this work were performed onthe supercomputer BwUniCluster 2.0 and ForHLR II at KIT funded by the Ministry of Science,Research and the Arts Baden-Württemberg and the Federal Ministry of Education and Research.The authors acknowledge support by the state of Baden-Württemberg through bwHPC.
From equation (3.3) follows that the effective refractivity between two positions (cid:174) 𝑖 and (cid:174) 𝑗 is calculatedvia the integral along the line of sight with length 𝑙 𝑖, 𝑗 : 𝑁 𝑖, 𝑗 = 𝑁 ( ) 𝜌 ( ) ∫ 𝑗𝑖 𝜌 ( ℎ ( 𝑙 )) d 𝑙𝑙 𝑖, 𝑗 . (8.1)For sources with a zenith angle 𝜃 (cid:46) ◦ the atmosphere can be approximated to be flat and theintegral over d 𝑙 can be substituted with d 𝑙 = d ℎ / cos ( 𝜃 ) . This simplifies the equation to 𝑁 𝑖, 𝑗 = 𝑁 ( ) 𝜌 ( ) 𝑇 𝑖 − 𝑇 𝑗 Δ ℎ 𝑖, 𝑗 (8.2)with the analytically described mass-overburden 𝑇 𝑥 = 𝑇 ( ℎ 𝑥 ) . For more inclined geometries thecurvature of the Earth has to be taken into account and the density along 𝑙 𝑖, 𝑗 becomes a function ofthe zenith angle and distance from ground. In that case the integral in Eq. (8.1) cannot be solvedanalytically anymore and thus has to be solved numerically. This is computationally very demandingand thus we use a pre-calculated table of the integrated refractivity . The table comprises theintegrated refractivity as a function of the zenith angle a given line of sight makes with the Earth’ssurface (which is not identical to the zenith angle measured at higher altitudes along the line of This is equivalent to the treatment included in CoREAS since version v7.7000. – 21 – igure 12 . Top : Total light propagation time for (9) different source positions and numerous differentobservers. The x-axis denotes the zenith angle under which a source is seen by an observer. The propagationtime is calculated along straight lines with the effective refractivity being calculated by very fine-grainedpiece-wise numerical integration or with the pre-calculated tables.
Middle : Absolute difference of the lightpropagation time between numerical and tabulated calculation.
Bottom : Relative agreement of the effectiverefractivity between numerical and tabulated calculation. sight) and the distance 𝑑 between the Earth’s surface and a point along the line of sight. For any twopoints in the atmosphere for which the line through both points also intersects with the sphericalEarth the integrated refractivity between those points can then be determined directly from thepre-tabulated values. The grid points of the table are spaced in tan 𝜃 and equidistant in distance 𝑑 .The integrated refractivity for any arbitrary point in the atmosphere is determined by a bi-linearinterpolation within this table. A python implementation of this table and the interpolation has beenmade publicly available at [40].In Fig. 12 we compare the effective refractivity determined with this model (pre-calculatedtable) to the exact numerical solution. The total light propagation time can be accurately calculatedusing the tabulated effective refractivity (cf. Fig. 12 top panel). A minor zenith-angle-dependentbias is visible in the absolute residual (cf. middle panel), however, the deviation is negligible since itis below any coherence criteria for signals in the MHz regime. In addition to the bias, small wiggles ,which originate from the unpyhsical, linear interpolation between zenith angle bins, are visible butalso negligible. – 22 – eferences [1] H. Schoorlemmer and W. R. Carvalho Jr. Radio interferometry applied to the observation of cosmic-rayinduced extensive air showers, 2020.[2] H. Falcke et al. Detection and imaging of atmospheric radio flashes from cosmic ray air showers. Nature , 435:313–316, 2005. doi: 10.1038/nature03614.[3] W. D. Apel et al. The wavefront of the radio signal emitted by cosmic ray air showers.
JCAP , 09:025,2014. doi: 10.1088/1475-7516/2014/09/025.[4] W. D. Apel et al. Final results of the LOPES radio interferometer for cosmic-ray air showers.
TheEuropean Physical Journal C , 81(2):176, Feb 2021. ISSN 1434-6052. doi:10.1140/epjc/s10052-021-08912-4. URL https://doi.org/10.1140/epjc/s10052-021-08912-4 .[5] I. Jandt. Beamforming with AERA. Master’s thesis, Bergischen Universität Wuppertal, 1 2012. URL https://astro.uni-wuppertal.de/fileadmin/physik/astro/mainpage/publications/theses/Diplom/Jandt-Diplom.pdf .[6] A. Romero-Wolf et al. An interferometric analysis method for radio impulses from ultra-high energyparticle showers.
Astropart. Phys. , 60:72–85, 2015. doi: 10.1016/j.astropartphys.2014.06.006.[7] A. G. Vieregg et al. A ground-based interferometric phased array trigger for ultra-high energyneutrinos.
PoS , ICRC2017:1013, 2018. doi: 10.22323/1.301.1013.[8] K. Hughes et al. Towards interferometric triggering on air showers induced by tau neutrino interactions.
PoS , ICRC2019:917, 2020. doi: 10.22323/1.358.0917.[9] T. Huege. Radio detection of cosmic ray air showers in the digital era.
Phys. Rept. , 620:1–52, 2016. doi:10.1016/j.physrep.2016.02.001.[10] F. G. Schröder. Radio detection of cosmic-ray air showers and high-energy neutrinos.
Prog. Part. Nucl.Phys. , 93:1–68, 2017. doi: 10.1016/j.ppnp.2016.12.002.[11] T. Huege and A. Haungs.
Radio Detection of Cosmic Rays: Present and Future , volume 9 of
JPSConference Proceedings . Journal of the Physical Society of Japan, Apr 2016. doi:10.7566/JPSCP.9.010018. URL https://doi.org/10.7566/JPSCP.9.010018 .[12] A. Aab et al. Observation of inclined EeV air showers with the radio detector of the Pierre AugerObservatory.
JCAP , 10:026, 2018. doi: 10.1088/1475-7516/2018/10/026.[13] B. Pont on behalf of the Pierre Auger Collaboration. A large radio detector at the Pierre AugerObservatory - measuring the properties of cosmic rays up to the highest energies.
PoS(ICRC2019)395 ,Jul 2019.[14] E. M. Holt, F. G. Schröder, and A. Haungs. Enhancing the cosmic-ray mass sensitivity of air-showerarrays by combining radio and muon detectors.
Eur. Phys. J. , C79(5):371, 2019. doi:10.1140/epjc/s10052-019-6859-4.[15] F. Schroeder et al. Radio Detection of Cosmic Rays.
SnowMass 2021 , CF7, 2020. URL .[16] F. G. Schröder et al. New method for the time calibration of an interferometric radio antenna array.
Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers,Detectors and Associated Equipment , 615(3):277 – 284, 2010. ISSN 0168-9002. doi: – 23 – ttps://doi.org/10.1016/j.nima.2010.01.072. URL .[17] P. Allison et al. Timing calibration and synchronization of Surface and Fluorescence Detectors of thePierre Auger Observatory. In ,volume 8 of
International Cosmic Ray Conference , page 307, jan 2005.[18] A. Aab et al. Nanosecond-level time synchronization of autonomous radio detector stations forextensive air showers.
Journal of Instrumentation , 11(01):P01018–P01018, jan 2016. doi:10.1088/1748-0221/11/01/p01018. URL https://doi.org/10.1088/1748-0221/11/01/p01018 .[19] E. M. Holt on behalf of the Pierre Auger Collaboration. Recent results of the Auger Engineering RadioArray (AERA).
PoS(ICRC2017)492 , Jul 2017.[20] P. Schellart et al. Detecting cosmic rays with the LOFAR radio telescope.
Astron. Astrophys. , 560:A98,2013. doi: 10.1051/0004-6361/201322683.[21] P.A. Bezyazeekov et al. Measurement of cosmic-ray air showers with the Tunka Radio Extension(Tunka-Rex).
Nucl. Instrum. Meth. A , 802:89–96, 2015. doi: 10.1016/j.nima.2015.08.061.[22] J. Álvarez Muñiz et al. The Giant Radio Array for Neutrino Detection (GRAND): Science and Design.
Sci. China Phys. Mech. Astron. , 63(1):219501, 2020. doi: 10.1007/s11433-018-9385-7.[23] T. Huege et al. High-precision measurements of extensive air showers with the SKA.
PoS , ICRC2015:309, 2016. doi: 10.22323/1.236.0309.[24] E. de Lera Acedo et al. SKALA, a log-periodic array antenna for the SKA-low instrument: design,simulations, tests and system considerations.
Experimental Astronomy , 39(3):567–594, Oct 2015. ISSN1572-9508. doi: 10.1007/s10686-015-9439-0. URL https://doi.org/10.1007/s10686-015-9439-0 .[25] F. G Schröder. Science case of a scintillator and radio surface array at IceCube. arXiv preprintarXiv:1908.11469 , 2019.[26] T. Huege, M. Ludwig, and C.W. James. Simulating radio emission from air showers with CoREAS.
AIP Conf. Proc. , 1535(1):128, 2013. doi: 10.1063/1.4807534.[27] D. Heck and T. Pierog. Extensive air shower simulation with CORSIKA: A user’s guide (version 7.7400from may 27, 2020). URL https://web.ikp.kit.edu/corsika/usersguide/usersguide.pdf .[28] C. Glaser.
Absolute energy calibration of the Pierre Auger Observatory using radio emission ofextensive air showers . PhD thesis, RWTH Aachen, 2017.[29] S. Ostapchenko. Monte carlo treatment of hadronic interactions in enhanced pomeron scheme:QGSJET-II model.
Phys. Rev. D , 83:014018, 1 2011. doi: 10.1103/PhysRevD.83.014018.[30] M. Bleicher et al. Relativistic hadron-hadron collisions in the ultra-relativistic quantum moleculardynamics model.
J. Phys. G , 25(9):1859–1896, 9 1999. doi: 10.1088/0954-3899/25/9/308.[31] M. Kobal. A thinning method using weight limitation for air-shower simulations.
Astropart. Phys. , 15:259–273, 6 2001. doi: 10.1016/S0927-6505(00)00158-4.[32] P. Mitra et al. Reconstructing air shower parameters with LOFAR using event specific GDASatmosphere.
Astroparticle Physics , 123:102470, Dec 2020. ISSN 0927-6505. doi:10.1016/j.astropartphys.2020.102470. URL http://dx.doi.org/10.1016/j.astropartphys.2020.102470 .[33] F. Schlüter et al. Refractive displacement of the radio-emission footprint of inclined air showerssimulated with CoREAS.
The European Physical Journal C , 80(7), 6 2020. doi: – 24 – https://doi.org/10.1140%2Fepjc%2Fs10052-020-8216-z .[34] A. Aab et al. Depth of maximum of air-shower profiles at the Pierre Auger Observatory. I.Measurements at energies above 10 . eV. Phys. Rev. D , 90:122005, Dec 2014. doi:10.1103/PhysRevD.90.122005. URL https://link.aps.org/doi/10.1103/PhysRevD.90.122005 .[35] A. Aab et al. Reconstruction of inclined air showers detected with the Pierre Auger Observatory.
Journal of Cosmology and Astroparticle Physics , 2014(08):019–019, aug 2014. doi:10.1088/1475-7516/2014/08/019. URL https://doi.org/10.1088/1475-7516/2014/08/019 .[36] A. Balagopal V., A. Haungs, T. Huege, and F. G. Schroeder. Search for PeVatrons at the Galactic Centerusing a radio air-shower array at the South Pole.
Eur. Phys. J. C , 78(2):111, 2018. doi:10.1140/epjc/s10052-018-5537-2. [Erratum: Eur.Phys.J.C 78, 1017 (2018)].[37] A. Aab et al. Deep-Learning based reconstruction of the shower maximum 𝑋 max using thewater-Cherenkov detectors of the Pierre Auger Observatory, 2021. URL https://arxiv.org/abs/2101.02946 .[38] S. Buitink et al. Method for high precision reconstruction of air shower 𝑋 max using two-dimensionalradio intensity profiles. Phys. Rev. D , 90:082003, Oct 2014. doi: 10.1103/PhysRevD.90.082003. URL https://link.aps.org/doi/10.1103/PhysRevD.90.082003 .[39] P. A. Bezyazeekov et al. Reconstruction of cosmic ray air showers with Tunka-Rex data using templatefitting of radio pulses.
Phys. Rev. D , 97:122004, Jun 2018. doi: 10.1103/PhysRevD.97.122004. URL https://link.aps.org/doi/10.1103/PhysRevD.97.122004 .[40] A tool package for cosmic-ray and neutrino radio detectors. URL https://github.com/nu-radio/radiotools . Revision from 11.08.2020.. Revision from 11.08.2020.