Searches for exclusive Higgs and Z boson decays into J /ψ γ , ψ( 2 S ) γ , and ϒ ( nS ) γ at √ s = 13 TeV with the ATLAS detector

Searches for the exclusive decays of the Higgs and Z bosons into a J /ψ , ψ( 2 S ) , or ϒ ( nS ) ( n = 1 , 2 , 3) meson and a photon are performed with a pp collision data sample corresponding to an integrated luminosity of 36 . 1 fb − 1 collected at √ s = 13 TeV with the ATLAS detector at the CERN Large Hadron Collider. No signiﬁcant excess of events is observed above the expected backgrounds, and 95% conﬁdence- level upper limits on the branching fractions of the Higgs boson decays to J /ψ γ , ψ( 2 S ) γ , and ϒ ( nS ) γ of 3 . 5 × 10 − 4 , 2 . 0 × 10 − 3 , and ( 4 . 9 , 5 . 9 , 5 . 7 ) × 10 − 4 , respectively, are obtained assuming Standard Model production. The corresponding 95% conﬁdence-level upper limits for the branching fractions of the Z boson decays are 2 . 3 × 10 − 6 , 4 . 5 × 10 − 6 and ( 2 . 8 , 1 . 7 , 4 . 8 ) × 10 − 6 , respectively.


Introduction
Following the observation of a Higgs boson H with a mass of approximately 125 GeV by the ATLAS and CMS collaborations [1,2], detailed measurements of its properties show no deviations from the SM.However, its role in fermion mass generation is still to be shown experimentally.In the Standard Model (SM), this mass generation is implemented through Yukawa interactions, and many theories beyond the SM predict substantial modifications of the relevant Higgs boson couplings to fermions.The ATLAS and CMS collaborations have reported measurements of the Higgs coupling to a third-generation fermion with a significance greater than five standard deviations in the H → τ + τ − channel [3][4][5].In addition, progress has been made in the third-generation quark sector with indirect evidence of the coupling of the Higgs boson to the top quark [3].This was recently complemented by direct observation of the associated production of the Higgs boson with a top-quark pair (t tH) [6,7].Evidence of Higgs boson decays into b b has also been found with a significance in excess of three standard deviations by both ATLAS and CMS [8][9][10].No experimental evidence of Higgs boson decays into the first-and second-generation fermions has yet been found, but direct searches were recently performed by the ATLAS Collaboration for H → c c [11] and H → µ + µ − [12,13] and by the CMS Collaboration for H → µ + µ − and H → e + e − decays [14].
Deviations of the cand b-quark Yukawa couplings from the SM expectations can lead to significant increases in the branching fractions for exclusive decays.These deviations can arise in beyond-the-SM theories; for example, the quark masses might not originate entirely from the Higgs mechanism but could also be induced by other subdominant sources of electroweak symmetry breaking [21].Other scenarios include the minimal flavour violation framework [22], the Froggatt-Nielsen mechanism [23], the Higgsdependent Yukawa couplings model [24], the Randall-Sundrum family of models [25], and the possibility of the Higgs boson being a composite pseudo-Goldstone boson [26].An overview of relevant models of physics beyond-the-SM is provided in Ref. [27].
Decays of the Higgs boson into a heavy, vector quarkonium state, Q ≡ J/ψ or Υ(nS), and a photon were searched for by the ATLAS Collaboration with up to 19.2 fb −1 of data collected at √ s = 8 TeV [28], resulting in 95% confidence level (CL) upper limits of 1.5 × 10 −3 for B (H → J/ψ γ) and (1.3, 1.9, 1.3) × 10 −3 for B (H → Υ(nS) γ) (n = 1, 2, 3).The former decay mode was also searched for by the CMS covers the range |η| < 2.4 with resistive plate chambers in the barrel and thin gap chambers in the endcap regions.
A two-level trigger and data acquisition system is used to record events for offline analysis [36].The level-1 trigger is implemented in hardware and uses a subset of detector information to reduce the event rate to at most 100 kHz.It is followed by a software-based high-level trigger which filters events using the full detector information and outputs events for permanent storage at an average rate of 1 kHz.

Data and simulated data
The search is performed with a sample of pp collision data recorded at a centre-of-mass energy √ s = 13 TeV corresponding to an integrated luminosity of 36.1 fb −1 .The integrated luminosity of the data sample has an uncertainty of 2.1% derived using the method described in Ref. [37].Events are retained for further analysis only if they were collected under stable LHC beam conditions and all relevant detector components were fully functional.
The data samples used in this analysis were collected by a combination of two triggers [36].The first required an isolated photon with transverse momentum p γ T greater than 35 GeV and at least one muon with p µ T greater than 18 GeV.The second required an isolated photon with p γ T > 25 GeV and a muon with p µ T > 24 GeV.Higgs boson production was modelled using the P -B v2 Monte Carlo (MC) event generator [38][39][40][41][42] with CT10 parton distribution functions [43] for the gluon-gluon fusion (ggH) and vector-boson fusion (VBF) processes.Both processes were calculated up to next-to-leading order (NLO) in α S .P -B was interfaced with P 8.186 [44,45] to model the parton shower, hadronisation, and underlying event with the AZNLO set of tuned parameters [46].Additional contributions from the associated production of a Higgs boson and a W or Z boson (denoted W H and Z H, respectively) were modelled by the P 8.186 event generator [44,45] with NNPDF23LO parton distribution functions [47] and the A14 set of tuned parameters for hadronisation and the underlying event [48].The production rates for the SM Higgs boson with m H = 125 GeV, obtained from the compilation in Ref. [27], are assumed throughout this analysis.The ggH production is normalised such that it reproduces the total cross section predicted by a nextto-next-to-next-to-leading-order QCD calculation with NLO electroweak corrections applied [49][50][51][52].The VBF production is normalised to an approximate next-to-next-to-leading-order (NNLO) QCD cross section with NLO electroweak corrections applied [53][54][55].The production of W H and Z H is normalised to cross sections calculated at NNLO in QCD with NLO electroweak corrections [56,57] including the NLO QCD corrections [58] for gg → Z H.The production of a Higgs boson in association with t t (t tH) or b b (b bH) is accounted for by scaling the total cross section used to normalise the ggH signal sample, assuming the signal efficiency of these processes to be equal to that for ggH.The addition of t tH and b bH to the ggH signal changes the acceptance by less than 1%.
The P -B v2 event generator was also used to model inclusive Z boson production.P -8.186 [44,45] with CTEQ6L1 parton distribution functions [59] and the AZNLO set of hadronisation and underlying-event parameters [46] was used to simulate parton showering and hadronisation.The prediction is normalised to the total cross section obtained from the measurement in Ref. [60], which has an uncertainty of 2.9%.This conservatively assumes the luminosity component of the uncertainty to be completely uncorrelated with that of this search.The Higgs and Z boson decays were simulated as a cascade of two-body decays, accounting for effects of the quarkonium helicity on the µ + µ − kinematics.
The branching fractions for the decays Q → µ + µ − are taken from Ref. [61].The simulated events were passed through the detailed G 4 simulation of the ATLAS detector [62,63] and processed with the same software used to reconstruct the data.
Selected Q → µ + µ − candidates are subjected to isolation and vertex-quality requirements.In this case, the primary pp vertex is defined as the reconstructed vertex with the highest Σ i p T 2 i of all associated tracks used to form the vertex.The sum of the p T of the tracks within ∆R = (∆φ) 2 + (∆η) 2 = 10 GeV/p µ T (with a maximum ∆R of 0.3) of the leading muon is required to be less than 6% of the muon p T [65].To mitigate the effects of multiple pp interactions in the same or neighbouring bunch crossings, only ID tracks that originate from the primary vertex are considered.The transverse momentum of the ID track associated with the leading muon is subtracted from the sum, and the subleading muon is also subtracted if it falls within the isolation cone.To reject backgrounds from b-hadron decays, the signed projection of the Q candidate's flight distance between the primary pp vertex and the dimuon vertex onto the direction of its transverse momentum is required to be less than three times its uncertainty.
Photons are reconstructed from clusters of energy in the electromagnetic calorimeter.Clusters without matching ID tracks are classified as unconverted photon candidates while clusters matched to ID tracks consistent with the hypothesis of a photon conversion into e + e − are classified as converted photon candidates [66].Reconstructed photon candidates are required to have transverse momentum p γ T > 35 GeV and pseudorapidity |η γ | < 2.37, excluding the barrel/endcap calorimeter transition region 1.37 < |η γ | < 1.52, and to satisfy the "tight" photon identification criteria [66].Isolation requirements are imposed to further suppress the contamination from jets.The sum of the transverse momenta of all tracks originating from the primary vertex, within ∆R = 0.2 of the photon direction, excluding those associated with reconstructed photon conversions, is required to be less than 5% of p γ T .In addition to this track isolation criterion, a calorimeter isolation criterion is applied where the sum of the transverse momenta of calorimeter energy clusters within ∆R = 0.4 of the photon direction, excluding the transverse energy of the reconstructed photon, is required to be less than 2.45 GeV + 0.022p γ T .The effects of multiple pp interactions in the events are also accounted for in the calorimeter isolation measurements.
Combinations of a Q → µ + µ − candidate and a photon satisfying ∆φ(Q, γ) > π/2 are retained for further analysis.When multiple combinations are possible, a situation that arises in only a few percent of the events, the combination of the highest-p T photon and the Q candidate with an invariant mass closest to the respective quarkonium mass is retained.To improve the sensitivity of the Υ(nS)γ analysis in resolving the individual Υ(nS) states, the events are classified into two exclusive categories based upon the pseudorapidity of the muons.Events where both muons are within the region |η µ | < 1.05 constitute the "barrel" (B) category.Events where at least one of the muons is outside the region |η µ | < 1.05 constitute the "endcap" (EC) category.
To maintain a single search region, while ensuring near-optimal sensitivity for both the Higgs and Z boson analyses, the transverse momentum of the Q candidate p Q T is required to be greater than a value that varies as a function of the invariant mass of the three-body system m Qγ .For the ψ(nS)

Signal modelling
For the ψ(nS) γ → µ + µ − γ final state, the total signal efficiencies (kinematic acceptance, trigger, reconstruction, identification, and isolation efficiencies) are 19% and 11% for the Higgs and Z boson decays, respectively.The corresponding efficiencies for the Υ(nS) γ → µ + µ − γ final states are 22-23% and 15-16%, respectively.The difference in efficiency between the Higgs and Z boson decays arises primarily from the softer p γ T and p µµ T distributions associated with Z → Q γ production, as seen by comparing Figures 1(a) and 1(b) for the J/ψ γ case and Figures 1(c) and 1(d) for the Υ(nS) γ case.
Accounting for quarkonium helicity effects in the Q → µ + µ − decays leads to a 3-4% decrease of the expected Higgs boson efficiency and a corresponding 6-9% increase of the expected Z boson efficiency.
The m µ + µ − γ resolution is 1.6-1.8%for both the Higgs and Z boson decays.For each of the final states, a two-dimensional (m µ + µ − γ and m µ + µ − ) probability density function (pdf) is used to model the signal.The Higgs boson signals are modelled with two-dimensional multivariate Gaussian distributions, which retain the correlation between m µ + µ − γ and m µ + µ − of the final states.For the Z boson decays, the m µ + µ − γ distributions of the signal are modelled with Voigtian pdfs (a convolution of Breit-Wigner and Gaussian pdfs) corrected with mass-dependent efficiency factors, and the m µ + µ − distributions are modelled as Gaussian pdfs.
The m µ + µ − distribution for selected ψ(nS) γ candidates, which pass the background generation region criteria described in Section 6, is shown in Figure 2(a) and exhibits clear peaks at the J/ψ and ψ(2S) masses.In Figures 2(b) and 2(c), the corresponding distributions for the selected Υ(nS) γ candidates are shown, where the Υ(1S, 2S, 3S) peaks can be observed.The ψ(nS) and Υ(nS) peaks are fitted with Gaussian pdfs and are used to cross-check the parameters obtained from the fit to simulated signal event samples, while the background is modelled with a Chebychev polynomial function.The experimental resolution in m µ + µ − is approximately 54 MeV for the J/ψ γ candidates (43 MeV for events in the barrel category and 64 MeV for events in the endcap category).

Background modelling
The dominant source of background exhibits a non-resonant m µ + µ − γ distribution and is composed of two distinct contributions: genuine Q → µ + µ − decays and non-resonant dimuon production.The decay Z → µ + µ − γ with final-state radiation (Z FSR) constitutes a further, smaller background contribution exhibiting a characteristic resonant structure in the m µ + µ − γ distribution.The m µ + µ − γ continuum background is modelled with a non-parametric data-driven approach using templates to describe the kinematic distributions.The normalisation of the background is extracted directly from a fit to the data.The shape of the background model in the final discriminant variable is also profiled in the fit.A similar procedure was used in the earlier search for Higgs and Z boson decays into J/ψ γ and Υ(nS) γ [28] and the searches for Higgs and Z boson decays into φ γ and ρ γ [30,31].
The background model uses a sample of 5500 ψ(nS) γ and 2300 Υ(nS) γ candidate events.These events pass all of the kinematic selection requirements described previously, except that the photon and T requirement of 30 GeV.These events constitute the background "generation region" defined in Section 6.
Q candidates are not required to satisfy the nominal isolation requirements and a looser minimum p Q T requirement of 30 GeV is imposed.These events define the background-dominated "generation region" (GR).From these events, pdfs are constructed to describe the distributions of the relevant kinematic and isolation variables and their most important correlations.The data control samples are corrected for contamination from Z → µ + µ − γ decays.
The pdfs of these kinematic and isolation variables are sampled to generate an ensemble of pseudocandidates, each with complete Q and γ four-vectors and the associated Q and photon isolation values.The important correlations among the kinematic and isolation variables of the background events, in particular between p Q T and p γ T , are retained in the generation of the pseudocandidates through the following sampling scheme: • Initially, a value for p Q T is sampled from the corresponding pdf.
• The distribution of p γ T is parameterised in bins of p Q T , and a value is sampled from the corresponding bin given the previously sampled value of p Q T .The muon isolation variables are parameterised in bins of p Q T and p γ T and sampled according to the previously selected values.• The distributions of the pseudorapidity difference between Q and γ, ∆η(Q, γ), the photon calorimeter isolation variable, and their correlations are parameterised in a binned two-dimensional distribution in the same bins of p Q T used to describe the p γ T and muon isolation variables.
• Given the selected values of relative photon calorimeter isolation and p Q T , a value for the relative photon track isolation is sampled.
• Given the selected values of ∆η(Q, γ) and p Q T , a value is sampled for the azimuthal angular separation between Q and γ, ∆φ(Q, γ).
• Values for η Q and φ Q are sampled from a binned histogram of the corresponding distributions in the data control sample.These are combined with ∆η(Q, γ) and ∆φ(Q, γ) to give the values of η γ and φ γ .
• A value for m Q is sampled from within the required region of m µ + µ − .Separate pdfs are used to describe the m µ + µ − distributions of resonant ψ(nS) and Υ(nS) production and non-resonant dimuons, which is referred to as "combinatoric" in the following.
The use of this procedure ensures a good description of the background and avoids any reliance on the accuracy or limited sample size of simulated background events.
The nominal selection requirements are imposed and the surviving pseudocandidates are used to construct templates for the m Q γ distributions which are then smoothed using a Gaussian kernel density estimation [67].Potential contamination of the GR sample from signal events is expected to be negligible, and it was verified, through signal injection tests, that such a potential signal contamination would not affect the shape of the background model.
The normalisation of the exclusive background from Z FSR is determined directly from the fit to the data.The shape of this background in the m µ + µ − γ distribution is modelled with a Voigtian pdf, while the m µ + µ − distribution is modelled with a first order polynomial.The parameters of the former pdf are derived from the simulated Z → Q γ signal samples.The parameters of the latter pdf are determined directly in the fit to data.
To validate this background model with data, the m Q γ distributions in several regions defined by kinematic and isolation requirements looser than the nominal signal requirements are used to compare the prediction of the background model with the data.Three validation regions are defined using the GR selection as the basis and adding the p Q T requirement (VR1), the muon isolation requirements (VR2), or the photon isolation requirement (VR3).The m Q γ distributions in these validation regions are shown in Figure 3 for the ψ(nS) γ and the Υ(nS) γ final states.The background model is found to describe the data well in all validation regions.

Systematic uncertainties
The systematic uncertainties in the expected signal yields are summarised in Table 1.Uncertainties in the Higgs production cross sections are obtained from Ref. [27].The Z boson production cross section uncertainty is taken from the measurement in Ref. [60].
The uncertainty in the acceptance of the Higgs boson signal due to missing higher orders in the generator, parton distribution functions, underlying-event set of tuned parameters, and parton showering is estimated by studying the effect of the different variations on the acceptance at generator level.The total uncertainty in the acceptance due to these effects is estimated to be 1.8%.For the Z boson, the respective uncertainty is determined to be 6% by comparing the generator-level acceptance predictions of the nominal simulated sample with those of M 5_ MC@NLO v2.2.2 [68] and S 2.2.1 [69].
Trigger efficiencies for photons are determined from samples enriched with Z → e + e − events in data [70].The systematic uncertainty in the expected signal yield associated with the trigger efficiency is estimated to be 2.0%.Photon identification efficiencies are determined using the enriched Z → e + e − events, as well as inclusive photon events and Z → γ events [66,71].The impact on the yields of the photon identification efficiency uncertainties, for both the converted and unconverted photons, are 1.4% for the Higgs and Z boson signals.The effect of the muon reconstruction and identification efficiency uncertainty is 2.8% [65].
The photon energy scale uncertainty, determined from Z → e + e − events and validated using Z → γ events [72,73], is propagated through the simulated signal samples as a function of η γ and p γ T .The uncertainty associated with the determination of the photon energy scale and resolution in the simulation is 0.3% in the yield of the Higgs and Z boson signal.Similarly, the systematic uncertainty associated with the scale of the muon momentum measurement is 0.2% [65].The effect of the energy/momentum scale and resolution uncertainties on the signal shape is negligible for the Higgs and Z boson signals.For the Z FSR background shape, the effect is 0.2% on the value of the mean.Including the corresponding nuisance parameter results in a 0.5% change in the expected limit of the Z → Υ(nS)γ channels with negligible impact in the other channels.
The shape of the background model is allowed to vary around the nominal shape, and the parameters controlling these systematic variations are treated as nuisance parameters in the maximum-likelihood fit described in Section 8. Three such shape variations are implemented through scale variations on the p γ T distribution in the model, linear distortions of the shape of the ∆φ(Q, γ) distribution, and a global tilt of the three-body mass distribution around a pivot point.The first two variations are straightforward alterations to the underlying kinematics of the pseudocandidates, which cause corresponding changes in the three-body mass.The third variation is applied directly to the final three-body mass template.These systematic uncertainties allow the background template to adjust to the observed data.

Results
The data are compared with background and signal predictions using a two-dimensional (2D) simultaneous unbinned maximum-likelihood fit to the m µ + µ − γ and m µ + µ − distributions.A simultaneous 2D fit is required to distinguish the Z FSR background from the Z → Qγ signal and the non-resonant background.The parameters of interest are the Higgs and Z boson signal normalisations.Systematic uncertainties are modelled using additional nuisance parameters in the fit; in particular, the background normalisations are free parameters.The fit uses the selected events with m Qγ < 300 GeV.
In total, 1033 events were observed in the ψ(nS) γ and 906 in the Υ(nS) γ signal regions.The expected and observed numbers of background events within the m Q γ ranges relevant to the Higgs and Z boson signals are shown in Table 2.The results of the background-only fits for the ψ(nS) γ and Υ(nS) γ analyses are shown in Figure 4.
Table 2: The number of observed events and the mean expected background, with its total uncertainty, for the m Q γ ranges of interest.The expected Z and Higgs boson contributions are shown for branching fraction values of 10 −6 and 10 −3 , respectively.
for for 81-101 120-130 The systematic uncertainties described in Section 7 result in a 1.0% increase of the expected 95% CL upper limit on the branching fraction of the H → ψ(nS) γ decays.For the Z → ψ(nS) γ decays, the effect is larger, 2.6%, mostly due to the systematic uncertainty in the background shape.Similar behaviour is observed in the Υ(nS) γ analysis with systematic uncertainties resulting in a 2.5-2.7%deterioration in the sensitivity to the H → Υ(nS) γ decays and a 2.8-2.9%deterioration in the sensitivity to the Z → Υ(nS) γ decays, also mostly due to the systematic uncertainty in the background shape.
On the basis of the fit to the observed data, the largest excess observed is 2.2σ in the search for Z → J/ψ γ.
GeV) and 54.4 GeV (52.7 GeV) are imposed for the regions m Qγ ≤ 91 GeV and m Qγ ≥ 140 GeV, respectively.The thresholds are varied between their minimum and maximum values as a linear function of m Qγ in the region 91 GeV < m Qγ < 140 GeV.

Figure 1 :
Figure1: Generator-level transverse momentum (p T ) distributions of the photon and of the muons, ordered in p T , for (a) H → J/ψ γ, (b) Z → J/ψ γ, (c) H → Υ(nS) γ and (d) Z → Υ(nS) γ simulated signal events, respectively.The leading muon candidate is denoted by p µ1 T and the subleading candidate by p µ2 T .The hatched histograms denote the full event selection while the dashed histograms show the events at generator level that fall within the analysis geometric acceptance (both muons are required to have |η µ | < 2.5 while the photon is required to have |η γ | < 2.37, excluding the region 1.37 < |η γ | < 1.52).The dashed histograms are normalised to unity, and the relative difference between the two sets of distributions corresponds to the effects of reconstruction, trigger, and event selection efficiencies.

Figure 2 :
Figure 2: Distribution of µ + µ − invariant mass for (a) ψ(nS) γ and Υ(nS) γ ((b) barrel and (c) endcap categories) candidates.The candidates satisfy the event selection but without the nominal isolation requirements and with a looser minimum p QT requirement of 30 GeV.These events constitute the background "generation region" defined in Section 6.

Figure 3 :
Figure 3: The distribution of m µ + µ − γ in data compared to the prediction of the background model for ((a), (b) and (c)) H(Z) → ψ(nS) γ and ((d), (e) and (f)) H(Z) → Υ(nS) γ in the VR1, VR2 and VR3 validation regions.Z FSR refers to the Z → µ + µ − γ background contribution.The background model is normalised to the observed number of events within the region shown.The uncertainty band corresponds to the uncertainty envelope derived from variations in the background modelling procedure.

Figure 4 :Table 3 : 4 −
Figure4: The m µ + µ − γ and m µ + µ − distributions for the selected (a) ψ(nS) γ and Υ(nS) γ ((b) barrel and (c) endcap categories) candidates along with the results of the maximum-likelihood fits with background-only models.Z FSR refers to the Z → µ + µ − γ background contribution.The ratios of the data to the background-only fits are also shown.

Table 1 :
Summary of the systematic uncertainties in the expected signal yields.
[77] −1 of √ s = 8 TeV pp collision data with the addition of the first upper limits on the H/Z → ψ(2S) γ decays.MNE/IFA, Romania; MES of Russia and NRC KI, Russian Federation; JINR; MESTD, Serbia; MSSR, Slovakia; ARRS and MIZŠ, Slovenia; DST/NRF, South Africa; MINECO, Spain; SRC and Wallenberg Foundation, Sweden; SERI, SNSF and Cantons of Bern and Geneva, Switzerland; MOST, Taiwan; TAEK, Turkey; STFC, United Kingdom; DOE and NSF, United States of America.In addition, individual groups and members have received support from BCKDF, the Canada Council, CANARIE, CRC, Compute Canada, FQRNT, and the Ontario Innovation Trust, Canada; EPLANET, ERC, ERDF, FP7, Horizon 2020 and Marie Skłodowska-Curie Actions, European Union; Investissements d'Avenir Labex and Idex, ANR, Région Auvergne and Fondation Partager le Savoir, France; DFG and AvH Foundation, Germany; Herakleitos, Thales and Aristeia programmes co-financed by EU-ESF and the Greek NSRF; BSF, GIF and Minerva, Israel; BRF, Norway; CERCA Programme Generalitat de Catalunya, Generalitat Valenciana, Spain; the Royal Society and Leverhulme Trust, United Kingdom.The crucial computing support from all WLCG partners is acknowledged gratefully, in particular from CERN, the ATLAS Tier-1 facilities at TRIUMF (Canada), NDGF (Denmark, Norway, Sweden), CC-IN2P3 (France), KIT/GridKA (Germany), INFN-CNAF (Italy), NL-T1 (Netherlands), PIC(Spain), ASGC(Taiwan), RAL (UK) and BNL (USA), the Tier-2 facilities worldwide and large non-WLCG resource providers.Major contributors of computing resources are listed in Ref.[77].Also at Centre for High Performance Computing, CSIR Campus, Rosebank, Cape Town; South Africa.Aso at CPPM, Aix-Marseille Université, CNRS/IN2P3, Marseille; France.e Aso at Département de Physique Nucléaire et Corpusculaire, Université de Genève, Genève; Switzerland.fAlsoat Departament de Fisica de la Universitat Autonoma de Barcelona, Barcelona; Spain.gAlso at Departamento de Física Teorica y del Cosmos, Universidad de Granada, Granada(Spain); Also at Department of Applied Physics and Astronomy, University of Sharjah, Sharjah; United Arab Emirates. iAlso at Department of Financial and Management Engineering, University of the Aegean, Chios; Also at Department of Physics and Astronomy, University of Louisville, Louisville, KY; United States of America.k Also at Department of Physics and Astronomy, University of Sheffield, Sheffield; United Kingdom.l Also at Department of Physics, California State University, Fresno CA; United States of America.m Also at Department of Physics, California State University, Sacramento CA; United States of America.n Also at Department of Physics, King's College London, London; United Kingdom.o Also at Department of Physics, St. Petersburg State Polytechnical University, St. Petersburg; Russia.
b c Also at CERN, Geneva; Switzerland.d h j Also at Institute for Mathematics, Astrophysics and Particle Physics, Radboud University Nijmegen/Nikhef, Nijmegen; Netherlands.z Also at Institute for Particle and Nuclear Physics, Wigner Research Centre for Physics, Budapest; Hungary.aa Also at Institute of Particle Physics (IPP); Canada.ab Also at Institute of Physics, Academia Sinica, Taipei; Taiwan.ac Also at Institute of Physics, Azerbaijan Academy of Sciences, Baku; Azerbaijan.ad Also at Institute of Theoretical Physics, Ilia State University, Tbilisi; Georgia.ae Also at LAL, Université Paris-Sud, CNRS/IN2P3, Université Paris-Saclay, Orsay; France.a f Also at Louisiana Tech University, Ruston LA; United States of America.ag Also at Manhattan College, New York NY; United States of America.ah Also at Moscow Institute of Physics and Technology State University, Dolgoprudny; Russia.ai Also at National Research Nuclear University MEPhI, Moscow; Russia.a j Also at Near East University, Nicosia, North Cyprus, Mersin; Turkey.ak Also at Physikalisches Institut, Albert-Ludwigs-Universität Freiburg, Freiburg; Germany.al Also at School of Physics, Sun Yat-sen University, Guangzhou; China.am Also at The City College of New York, New York NY; United States of America.an Also at The Collaborative Innovation Center of Quantum Matter (CICQM), Beijing; China.ao Also at Tomsk State University, Tomsk, and Moscow Institute of Physics and Technology State University, Dolgoprudny; Russia.ap Also at TRIUMF, Vancouver BC; Canada.aq Also at Universita di Napoli Parthenope, Napoli; Italy.
y * Deceased