Two-sided constraints on Lorentz Invariance violation from Tibet-AS$\gamma$ and LHAASO Very-High-Energy photon observations

We present new two-sided constraints on Lorentz Invariance violation energy scale for photons with quartic dispersion relation from recent gamma ray observations by Tibet-AS$\gamma$ and LHAASO experiments. The constraints are based on the consideration of the processes of photon triple splitting (superluminal scenario) and the suppression of shower formation (subluminal). The constraints in subluminal scenario are better than the pair production constraints and are the strongest in the literature.

Introduction Since 2019, a new window in Very-High-Energy (VHE) astronomy focusing on gamma-rays with energy from 100 TeV to several PeV has been opened. The photon spectra of galactic sources which continue beyond 100 TeV have been measured by the collaborations Tibet-ASγ [1], HAWC [2,3,4,5], LHAASO [6,7] and Carpet-3 [8]. In addition, the Tibet-ASγ collaboration has reported the detection of a diffuse photon flux from the Galactic disk [9]. Two of these collaborations, Tibet-ASγ [9] and LHAASO [7], have reported single photon events with energy near PeV -957 TeV and 1.42 PeV, respectively. It is expected that these recent experimental data will shed a lot of light on the mechanism of cosmic rays acceleration in galactic sources. In addition, the physics of these gamma rays may be sensitive to modification of some fundamental laws of physics, such as violation of Lorentz Invariance (LI).
The motivation for considering broken LI was born in several different attempts to construct quantum gravity theory (see [10,11] and references therein for a review). Lorentz Invariance Violation (LIV) in the matter sector is usually phenomenologically considered in the framework of Standard Model Extension [12] (SME). The main feature of LIV is the modification of dispersion relation for particles [13]. This results in kinematical effects such as the modification of thresholds for some processes in particle physics [13,14]. However, the dynamical effects of LIV including the modification of cross-sections for non-threshold processes are also important. The processes important for LIV photon propagation and detection in the atmosphere include photon decay [13,14,15], photon triple splitting [16,17], pair production on soft photon background [18,19] and in the Coulomb field of nuclei [20,15]. The last two processes hold in the standard case and modified in a LIV scenario. The modern constraints on LIV in the photon sector are gathered in the Data Tables [21], see also the recent review [22].
On this date, several constraints on LIV parameters in the photon sector based on LHAASO data have been already obtained in [23] and [24] 1 . The article [23] is devoted to constraining cubic corrections to the photon dispersion relation while the article [24] is devoted to quartic correction, focusing on the processes of photon decay and photon splitting. In this Letter we complement these works with Tibet-ASγ data as well as with the constraints from air shower formation.
The Model. In this Letter we focus on quartic corrections to the photon dispersion relation, Here M LV denotes the energy scale of LIV quantum gravity 2 . The sign plus (minus) in the dispersion relation (1) refers to a superluminal (subluminal) LIV scenario. The dispersion relation (1) can be obtained from the Lagrangian, Note that Feynman rules for photons in the model (2) are modified, see [15] for details.
Photon decay. The simplest process characteristic to a superluminal LIV scenario, is the photon decay to an electron-positron pair γ → e + e − . The process has a certain energy threshold: below the threshold it is forbidden by the energy-momentum conservation; for energies larger than the threshold one the decay rate is very fast 3 [15]. The detection of even a single photon of a given energy means that the energy is larger than the threshold one and so establishes the constraint, The significance of the constraint (3) coincides with the significance of the corresponding photon event.
1 After the first version of the current paper the LHAASO collaboration article on constraining LIV appeared [25].
2 In the literature the notation E QG, (2) is also used. 3 The mean free path for TeV photons is about one hundred meters.
Photon splitting. In the superluminal scenario γ → e + e − is not the only channel of the photon decay. The process of photon splitting to three photons γ → 3γ occurs without threshold if the photon dispersion relation is superluminal, so if the photon decay is kinematically forbidden the splitting may be the most effective process. The splitting width reads [17], The splitting process starts making a visible effect when the photon mean free path associated with splitting is of the same order as the distance to the source, (Γ γ→3γ ) −1 ∼ L source . The probability for a given photon to pass the distance L source not being splitted is given by the exponential distribution, Considering the initial photon flux instead of a single photon, P splitting determines the suppression of the photon flux.
Shower formation. In the subluminal LIV scenario any types of spontaneous photon decay are forbidden. Moreover, some photon decay processes in an external field, allowed in the standard LI case, may become suppressed. The most important of these processes is the Bethe-Heitler process -photon decay to an electron-positron pair in the Coulomb field of a nuclei. This process is crucial for the formation of atmospheric showers initiated by high-energy photons. The suppression of the Bethe-Heitler process in a subluminal LIV scenario was firstly estimated by [20] and calculated in the model in [15]. As a result, photon-induced air showers in the subluminal scenario would start in the atmospheric deeper than in the standard case. Very deep showers cannot be registered by the experiment. As the simplest criterium for the registration we compare the depth X 0 of the photon first interaction in the atmosphere with the total atmosphere depth X atm . If X 0 > X atm , the mean shower starts near the ground surface, so the event will not be detected. Following this criterium, the probability for pair production of a photon in the atmosphere reads, where the mean depth of the first interaction for the LIV case X 0 LV is expressed via the LI mean depth X 0 LI = 57 g cm −2 and the ratio of the LI and the LIV Bethe-Heitler cross-sections. Taking into account the expression for LV Bethe-Heitler cross-section [15], one obtains here the limit Suppression of the observed photon flux. One can note that the different effects of photon splitting and suppression of shower formation appeared in different LIV scenarios (super-or subluminal) lead to the same result: the suppression of the observed photon flux. Thus, the prediction for observed flux in superluminal/subluminal LIV models reads, Here the suppression factor P (E, M LV ) is given by eq. (5) for photon splitting and by eq. (7) for shower formation. The hypothesis of a certain P (E, M LV ) can be tested with experimental data. Note that these tests require a certain model of photon flux dΦ dE LI from the source. In order to set the most robust constraints on M LV , in the case of more than one source models we choose the model with the maximal predicted flux.
Pair Production on CMB Another process important for LIV probing is pair production by a high-energy photon on soft photon backgrounds γγ b → e + e − . This process occurs above a certain threshold, which for the background photon reads th = m 2 e Eγ in the LI scenario. In the presence of LIV the threshold gets modified [26], Here the sign "-" refers to superluminal propagation while the sign"+" to subluminal. The cross-section of the process above the threshold does not significantly depend on LIV [15]. The important outcome of this process is the attenuation of high-energy (TeV energy range) flux of extragalactic photons due the interaction with extragalactic background lights (EBL) [27]. In the subluminal LIV scenario the pair production threshold is shifted to lagrer energies of background photons (related to smaller EBL density), making the extragalactic photon spectrum less attenuated [18,19]. The absence of this effect constrains the subluminal LIV scenario [22]. Considering galactic gamma rays in the 1 − 100 TeV energy range, one can neglect the pair production effect due to small distances inside the galaxy compared to the mean free path of pair production on the infrared background.
For the galactic photons of sub-PeV energies observed by Tibet-ASγ and LHAASO the threshold of pair production decreases to the energy range of cosmic microwave background (CMB), which density exceeds the EBL density by more than two orders of magnitude. The mean free path for a ∼ 1 PeV photon related to the pair production on CMB is ∼ 10 kpc [28,29], which is comparable to galactic scales. The galactic gamma-ray flux is attenuated by the factor P p.p.CM B = e −Lsource/L mf p (Eγ ) . In the subluminal LIV scenario the threshold increases above the energy of the CMB peak, and the attenuation disappears. Thus, in the subluminal LIV scenario the modification of pair production increases the observed gammaray flux while the shower suppression effect leads to the suppression. However, for the spectra observed in [7,9] the factor P p.p.CM B cannot be less than 0.5, so the effect of the modified pair production on CMB seems to be subleading compared to the shower suppression effect.
In the next two paragraphs we will consider the constraints from the observational data of Tibet-ASγ [9] and LHAASO [7].
Tibet-ASγ. Tibet-ASγ is a hybrid experiment devoted to the observations of VHE photons and cosmic rays, which is located at the altitude 4300 meters above sea level in Tibet, China. The whole setup consists of an array of surface air shower (AS) detectors, combined with an underground array of muon detectors. The presence of a muon detector array yields a very high degree of gamma/hadron discrimination. In the recent work [9] the Tibet-ASγ collaboration presented gamma ray events from the dataset of 719 days of observation from the whole Tibet-ASγ field of view with energies larger than 100 TeV, grouped into three energy bins.
In the last energy bin centered at 534 TeV, 16 gamma-ray-like events from the galactic disk (angle |b| < 5 deg from the galactic plane) against 1.35 background events (statistical significance 6.0 σ) have been observed. The maximal energy of the observed gamma-ray event is 957 +166 −141 TeV. All of these events are not connected with any known TeV source, so they are assumed to be diffuse gamma rays. According to modern theoretical models [31,32,33], it is assumed that these diffuse gamma rays are one of the products of hadronic interaction of charged cosmic rays in the interstellar medium. The measured diffuse flux in the last energy bin exceeds the model prediction [31] by the factor of two. However, the authors of a recent work [34] suggested the model of diffuse gamma rays from the outer Galactic disk which is in an agreement with Tibet-ASγ data as well as with Fermi/LAT data on GeV-range diffuse gamma rays. In both superluminal and subluminal LIV scenarious the observed gamma-ray flux is suppressed compared to the model [34]. We test the hypothesis of this suppression with the Tibet-ASγ observational data. For the calculation we have to specify a number of parameters. The maximal atmospheric depth at the Tibet location is X atm = 780 g cm −2 since the maximal zenith angle for Tibet-ASγ events is 40 degrees. Although the exact distance to the source (outer Galactic disk) for each diffuse gamma-ray event is unknown, one can estimate the distance with the maximal/minimum level in order to make the most conservative LIV constraint for a given process. Thus, considering the splitting process we take the minimal distance L source ≈ 1 kpc while for the pair production on CMB we take L source ≈ 5 kpc. The absorption coefficient related to the pair production on CMB for the last energy bin reads P p.p.CM B = 0.73 [29].
Diffuse gamma ray spectra of the outer Galactic disk measured by the Tibet-ASγ collaboration, the theoretical model of the spectrum by Koldobskiy et.al. [34] as well as the predictions of both LIV hypotheses are presented at Fig. 1. First, one can notice that the cut-off of the spectrum in the splitting scenario is very sharp. On the other side, the spectrum steepening associated with the formation of deeper air showers is more smooth. However, even for the mean energy for the pre-ultimate bin the difference between LIV and LI scenarios is negligible. So, we test the hypothesis of a certain LIV mass scale only with the last bin with the χ 2 test. The results are the following. The constraint associated with the absence of the splitting process (superluminal LIV scenario) reads We set the photon decay bound with the lowest energy estimation for the highest energy photon event, E = 705 TeV: at 95% CL the energy of the event is larger than E. Substituting the energy to (3), one obtains (superluminal) M LV > 4.9 × 10 14 GeV, 95% CL.
Note that the photon decay bound (12) is one order of magnitude weaker than the splitting one (10). A hadronic origin of the diffuse gamma rays requires that the direct process of diffuse photon production is neutral pion decay, π 0 → γγ. In the scenario of superluminal LIV photons the pion decay process becomes kinematically forbidden if M LV > 2E 2 γ mπ , where m π is the neutral pion mass (LIV in pions is assumed to be negligible). However, this bound is weaker than the photon decay bound (3).
LHAASO The Large High Altitude Air Shower Observatory (LHAASO) is a new-generation experiment designed for observation of photons and cosmic rays of very high energies located at the altitude 4410 meters in China. The whole installation consists of three components: a Kilometer Square Array (KM2A) of AS detectors (which include surface and underground muon detectors) connected with the Water Cherenkov Detector Array (WCDA) and Wide Field-of-view Cherenkov Telescope Array (WFCTA). The observational data from less than a year of KM2A operation has been presented in a recent paper [7]. The LHAASO collaboration declared [7] the observation of 12 galactic gamma-ray sources above 100 TeV energy. Particularly, the spectra of three gamma-ray sources, J2226+6057, J1908+0621 and J1825-1326, which continue up to 500 TeV, have been provided in [7]. In this paragraph we analyze this data and set the two-sided constraints on the mass scale M LV . The photon decay and photon splitting constraints (both related to the superluminal scenario) from LHAASO data have been already obtained in [24] by a slightly different method; the results will be compared in this paragraph. We complement these constraints with shower formation constraints. The maximal zenith angle of LHAASO events is 50 deg [7], so the maximal atmosphere depth for inclined showers at the LHAASO location is X atm = 890 g cm −2 .
It is shown in [7] that three observed source spectra have a slight steepening above 100 TeV compared to a power-law fit for the spectra, so the best fit for each of these spectra is a log-parabola curve. However, a power-law dependence is still not excluded. In order to estimate the maximal predicted flux at high energies we take the power-law fit. Based on this prediction we test both LIV hypotheses with experimental data (8).
At the Fig. 2 we present differential energy spectra for each LHAASO source: experimental data points [35], the best power-law fit, and the prediction for the flux under the superluminal and subluminal LIV hypotheses. Note that the energy bins of LHAASO data are smaller than Tibet-ASγ ones. The bin selection criterion for the LIV hypotheses test is as follows. On one side, we do not consider energy bins for which the difference between LI and LIV models are negligible. On the other side, we omit bins which statistical significance is less than two sigma since our aim is a 95% CL exclusion limit -the last bin for the source J1908+0621 and the two last bins for J2226+6057 and J1825-1326. It turns out that in order to test the splitting-induced cut-off a single bin for each source is relevant while for the shower suppression effect 4 one has to analyze two bins for each source. Following these prescriptions we perform a χ 2 test and constrain the mass scale M LV at 95% CL for both superluminal and subluminal LIV scenarios. The results are presented in Table 1.
The difference in the superluminal scenario between our constraints and the ones obtained by Chen et al. [24] is less than a factor of two. The possible explanation is that we take the most conservative power-law fit for the spectra while Chen et al. [24] used a log-parabola fit.
Let us also provide here the photon decay bound based on the gamma ray bound with the maximal energy obtained by [24], Here the 95% CL lower energy bound E γ = 1.21 PeV for the maximal energy event has been taken.

Discussion
We have obtained two-sided constraints on the LIV mass scale from Tibet-ASγ and LHAASO observations of gamma rays in the sub-PeV energy range. The constraint on M LV for the subluminal LIV scenario, based on Tibet-ASγ data, M LV > 1.7 · 10 13 GeV, (eq. (11)) is the strongest in the literature. The previous constraint of this type [30], based on Tibet-ASγ observation of Crab Nebula spectrum [1], is improved by a factor of 20. The constraint on the same parameter, based on the process of pair production on extragalactic background lights, is M LV > 2.4 · 10 12 GeV [36] -one order of magnitude weaker than (11). The shower formation constraint is based on the model of the given spectrum of the source. Although there is a variety of models for the diffuse gamma ray flux from the Galactic disk, the significantly larger flux than proposed in [34] at sub-PeV energies would contradict to the measurements of the diffuse gamma rays in the GeV energy range as well as to observations of sub-PeV astrophysical neutrinos. It is worth noting that the method for constraining the suppression of shower formation based on a decreasing flux in a given energy bin leads to conservative bounds. Better constraints can be obtained by analysing the shower parameters (basically, the maximal depth X max ) for each photon-induced shower. This method has been proposed in [37] for constraining LIV on photons in 10 19 − 10 20 eV energy range, which have not been detected yet. Inverting this argument, one can propose a unique signal for searching for subluminal LIV in atmospheric showers: the LIV photon shower would start deeper (X max would be larger) than the LI one, but the rest of photon shower parameters such as shower width and number of muons would not change.