Search for Low-Mass Leptoquarks using the $llqq$ final states in Pb-Pb Ultra-Peripheral Collisions

After a review over past experiments and theoretical requirements from low-energy flavour physics, we argue that the possibility of low-mass leptoquarks (LQ) cannot be fully excluded due to the assumptions made in the measurements. Therefore we propose to search for the pair production of low-mass leptoquarks in Pb-Pb ultra-peripheral collisions with the least model dependence (assuming a small coupling constant $\lesssim0.1$). There are a couple of advantages: 1) high photon flux provides high production rate for low mass LQs; 2) the background contamination is much lower than that in $p$-$p$ collisions. The analysis strategy permits a leptoquark to decay to all possible lepton-plus-quark modes. Taking the scalar LQ, $S_3$, with an electric charge $|q|=\frac{4}{3}e$, as example, the mass point of 100 GeV can be excluded at the 95~\% confidence level using a dataset of 4~pb$^{-1}$ Pb-Pb ultra-peripheral collisions at $\sqrt{s}=5.02$ TeV and the performance of the ATLAS detector in Run~2. The proposed method also applies to searching for high-mass LQs in $p$-$p$ collisions as long as the LQ pair production mechanism dominates.

Various LQs have been searched for in the LEP, Tevatron, LHC and other experiments. Ref. [32] is a good review on the LQ physics and experimental searches up to 2016. Table 1 summarizes a list of experimental searches from different colliders. From the table, it seems that LQs are excluded up to 1 TeV or more. However, all searching results are interpreted with assumptions. One common assumption is that a LQ couples only to the same generation of fermions (this assumption is denoted by "SG" in Table 1). For example, a so-called "third-generation" LQ with a charge + 2 3 e only couples to bτ + or tν τ . This assumption is not appropriate because of more and more strong evidences for LUV. Taking the decay B + → K + l + l − as in the LHCb measurement [24] as example, a LQ must be able to couple to µ/e and b/s quark if LQ is the right solution. A much looser assumption is that a LQ couples to a single lepton plus all possible quarks (this assumption is denoted by "SL" in Table 1). This assumption was used in the OPAL experiment [41,42]. Inspired by their strategy, a model-independent method is proposed in Sec. 3. For the measurements using single LQ production processes, we have to interpret the results as a function of the LQ mass and the coupling constant between a LQ and the fermions (denoted by λ or λ lq for the lepton l and the quark q). One common way is to set constraints on the LQ mass assuming a coupling strength which is equal to that in the electromagnetic interaction, namely, λ = √ 4πα em in Table 1, where α em ≈ 1 137 is the fine structure constant. A combined fit to all measurements is desirable to avoid the assumptions, but not available for the moment. In Sec. 2, we will show that 80 GeV is a robust lower limit for the LQ with an electric charge |q| = 4 3 e based on the OPAL measurements [41,42]. Considering a LQ with an electric charge |q| = 4 3 e, it has 9 decay modes (LQ → l + q with l = e, µ, τ and q = d, s, b). Except the OPAL experiment, there is no search for the final state with a tau lepton plus light quarks. Therefore, it seems that we cannot deny the possibility of a low-mass LQ based on previous searches without a global fit. On the other hand, we will show that low-energy flavour physics [67][68][69][70] can impose strong constraint on combinations of λ lq and m LQ in Sec. 2. We can only set upper limits on λ lq for a fixed LQ mass or vice versa. Given these difficulties, model-independent methods to search for LQs are highly necessary.
For the low mass LQs (less than 200 GeV), we find that the nuclear ultra-peripheral collision (UPC) is an ideal place to perform such searches. In UPCs, relativistic nuclei collide at impact parameters (b) that are so large that there are no hadronic interactions. Thus it provides a unique opportunity to study photon-photon interactions. Ref. [72,73] are two good reviews about the physics of UPCs. Given high photon flux at low energies (proportional to the 4-th power of the nucleus charge) and much cleaner collision environment, we propose to search for the LQ pair production in Pb-Pb UPC events. Most of previous searches use the property that LQs carry color charge. This proposal uses the property that LQs have electric charge, and is much less model-dependent compared to the existing searching strategies. Taking the scalar LQ with an electric charge |q| = 4 3 e, denoted as S 3 as in the review [32], as example, a feasibility study using all possible llqq (l = e, µ, τ , q = d, s, b and charge conjugation is always assumed) final states is conducted based on the performance of the ATLAS detector in Run 2 and Pb-Pb UPCs at the center-of-mass (c.m.) energy √ s = 5.02 TeV. In Sec. 2, the constraints on the LQ mass and coupling constant from the experiments and the low-energy flavour physics are reviewed. In Sec. 3, a much less model-dependent method is proposed to search for LQ pair production, which process dominates if the coupling strength of the LQ to the fermions is not big ( 0.1 as shown in Fig. 1). It is worth emphasizing that the proposed method itself is also applicable to search for high-mass LQs Table 1. Overview of the experimental searches for the leptoquarks. "SG" represents the assumption that a LQ couples to the same generation of fermions. "SL" represents the assumption that a LQ couples to single lepton and all possible quarks. It should be noted that the presented exclusion region may be just part of the experiment results.

2
The bounds from the experimental searches and low-energy flavour physics

A robust lower mass bound from the OPAL experiments
For a LQ, both its mass and coupling constant are unknown. The partial width of LQ → l+q is where C = 16 for a scalar LQ and 24 for a vector LQ. The LQ would decay outside a detector for a mass of 100 GeV and a coupling constant of less than 10 −8 . Since a LQ carries both electric and color charge, a long-lived one would interact with the detector matter and may appear like a b-tagged jet. There is no dedicated search for such objects. If the lifetime is not too long, some direct searches are summarized in Table 1. Among these results, two of the OPAL measurements [41,42] seem to use the weakest assumption, which forbids a LQ to couple to different generations of leptons, but allow to couple to all possible quarks. In addition, it used the pair production in e + e − collisions. The pair production dominates if λ 10 −2 and its cross section does not depend upon λ. Collecting the branching fractions, we can naively interpret their results in a different way (taking the LQ with the electric charge |q| = 4 3 e and the measurement [42] as example).
where σ(m) is the production cross section of a pair of LQs with a mass m in the process e + e − → LQLQ; B lx is the branching fraction of LQ → lx with x = d, s, b for a fixed lepton l. Using B ex + B µx + B τ x = 1 and the Cauchy-Schwarz inequality, we have Using the cross section formula in Ref. [42] and the highest collision energy, a conservative bound is obtained. m > 80 GeV (2.8) This interpretation neither considers the shape information of the invariant mass reconstructed from a lepton candidate and a jet, nor considers the fact that there are multiple collision energy points with different luminosities. It does not take into account the contribution from LQ → τ + q with a subsequent leptonical τ decay in the eeqq or µµqq signal region, either. A full consideration shall give a better bound (but would not exceed ∼100 GeV anyway). Thus we claim this is a robust bound with little model dependence.
A recent ATLAS measurement looked at the final states eq and µq with q being all possible quarks except the top quark and searched for the LQ mass above 400 GeV [57]. Exclusion regions as a function of the LQ mass and the branching fraction is reported (more details see [57]). Unfortunately, there seems no search for the final state with a tau lepton and light quarks except the OPAL experiment. It is difficult to exclude LQs with higher mass without assumptions.
Other experiments as summarized in Table 1 either focused on only specific final state or made very strong assumptions. No further bounds can be obtained without a combined fit.

The bounds from low-energy flavour physics
Low-energy flavour physics like rare meson decays and lepton decays can provide stringent bounds on combinations of the coupling constant and mass [67][68][69][70]. In most cases, a LQ propagator would inevitably couple to a lepton and a quark on both ends. The extra amplitude induced by a LQ would be proportional to a factor of λ lq λ l q /m 2 LQ if the involved energy scale is much smaller than the LQ mass. This factor is invariant under the scaling transformation λ → kλ and m LQ → km LQ . Here we quote a few bounds from the literature. For K L → e ± µ ∓ , Ref. [68] obtained For the µ − → e − conversion on the Titanium [71], Ref. [69] obtained Interestingly, the neutral meson (K 0 , D 0 , B 0 ) mixing is found to be able to constrain λ 2 /m LQ . A recent calculation [70] gives roughly (the exact bound is model-dependent) Due to the scaling invariance above, one cannot provide definite bounds on a LQ's mass due to the mixing with the coupling constant. The existing anomalies in B-meson decays need significant contribution from new physics. The bounds based on those anomalies are highly model-dependent. We shall not expand the discussion here. In order to effectively bound a LQ's mass without touching its coupling to the fermions, we shall consider the physics processes where only LQs and gauge bosons are involved. For example, we can measure the differential cross section of γγ → γγ, γγ → gg or gg → γγ. Both the SM particles and LQs contribute at loop level. The first process is the cleanest, but the first and the second ones shall only be sensitive to low-mass LQs in heavy-ion collisions. The third one might be sensitive to high-mass LQs in p-p collisions, though it suffers from the irreducible background qq → γγ. The resonance effect near the threshold of a pair of LQs shall enhance the sensitivity. Anyway, the method presented in next section is more available and straightforward though more experimental efforts are needed.

A model-independent searching method using the LQ pair production
The LQ, S 3 , considered in this feasibility study has a charge |q| = 4 3 e and is allowed to decay to all possible lq modes (l = e, µ, τ and q = d, s, b). The pair production in p-p collisions and Pb-Pb UPCs dominate if the coupling constant is not too big (e.g. λ 0.1 as shown in Fig. 1). This is consistent with the bounds from the low-energy physics in last section. The cross section of pair production does not depend on the coupling constant, but on the branching fractions of the decay modes. Thus we consider all possible final states, namely, eejj, eµjj, eτ had jj, µµjj, µτ had jj, and τ had τ had jj, where τ had denotes a hadronic tau decay candidate and j denotes a jet candidate. We do not attempt to identify the jet flavour for simplicity. For each final state, all LQ's contributions are considered. For example, all LQ→ eq/µq/τ q decays contribute in the eµjj final state.
Three floating parameters, κ, B ex and B µx , are introduced in the statistical interpretation. The branching fraction B τ x is constrained to be 1 − B ex − B µx . κ is the ratio of the measured cross section to the model prediction. We further impose the conditions 0 ≤ B ex , B µx ≤ 1 and B ex +B µx ≤ 1. For a given LQ mass, the measurement can determine the upper limit at 95 % confidence level (C.L.) of κ. If the upper limit is less than 1, then a LQ with the mass is excluded.
The method is model-independent in the sense that it allows a LQ to couple to all possible fermions.

Cross sections in UPC
For a nucleus A, the cross section for the process where σ γγ→X (k 1 , k 2 ) is the production cross section of the final state X in two-photon interactions with energy k 1 and k 2 ; d 2 L γγ /dk 1 dk 2 is the corresponding two-photon luminosity.
where R A is the nucleus radius; and d 3 N γ /dkd 2 b is the photon flux from a nucleus at a distance b. For a charge Z nucleus with velocity β and Lorentz boost factor γ, the photon flux is where x = kb βγ and K 0 (x) and K 1 (x) are modified Bessel functions. The photon flux decreases exponentially above the cut off k max ≈ γ R A . In Eq. 4.2, we have assumed that the probability of not having a hadronic interaction, denoted by P 0had (b), is 0 if b < R A and 1 if b > R A . In general, we can replace the integra- Using a more complicated description of P 0had (b) leads to a 10 − 15% difference in the flux [72].

MC sample preparation
The MadGraph 5 [74] is used to compute cross sections by convolving the Weizäcker-Williams equivalent photon flux (Eq. 4.3) of the nucleus with the elementary two-photon cross section, and to produce the signal and background samples. The parton showers are simulated with Pythia 8 [75] and the detector response is simulated with Delphes 3 [76] based on the Run 2 performance of the ATLAS detector.
The signal model is obtained from [77] and documented in Ref. [78]. The LQ of interest in this work is the scalar boson S 3 , which is allowed to couple to all charged leptons and down-type quarks. Six signal samples corresponding to the processes γγ → S 3S3 → eeqq/eµqq/eτ qq/µµqq/µτ qq/τ τ qq are produced with a mass 100 GeV and all coupling constant fixed at 10 −2 . The possible interference effect with the SM process is not considered. The fiducial region at generator level is defined as where p T , η and ∆R represent the transverse momentum, rapidity and the angular distance 1 , respectively. Figure 1 shows the fiducial cross section of Pb+Pb→Pb+Pb+e + e − dd as a function of the coupling constant λ ed . Little dependence is seen for λ ed 0.1. The cross section used in the normalization does not consider the branching fractions, which will be treated as the parameters of interest in the statistical interpretation, and is scaled by a factor of 10 −2 in the plots in Sec. 5 just for better illustration.  The dominant background is from the SM process γγ → l + l − qq. It is irreducible. The background events with multiple leptons or multiple jets are also considered and their contribution is minor. In addition, the background from Pb-Pb hadronic interactions is inevitable. It is very difficult to simulate this background in reality. We usually design tight cuts to suppress it and estimate it in a data-driven way. We assume this background is negligible after the event selection in next section.

Event selection and Observable reconstruction
Full reconstruction of all 4 objects is not desirable because of low efficiency. We require two lepton candidates and at least one jet candidates to be reconstructed. A charged lepton candidate (e/µ) is required to have transverse momentum p T > 15 GeV and rapidity |η| < 2.5. A hadronic tau candidate τ had is required to have p T > 20 GeV and |η| < 2.5. The associated number of charged tracks is required to be 1 or 3 with unit sum of charges. A jet candidate denoted by j is required to have p T > 20 GeV and |η| < 5.0.
An overlapping removal is performed according to the priority order µ → e → τ had → j. We first select muon candidate and then select electron candidate. The electron candidate is required not to overlap with the muon candidate within the angular distance ∆R = 0.4. τ had candidate is required not to overlap with electron and muon candidates within the same angular distance. Similarly the j candidate is required not to overlap with µ, e and τ had candidates. 12 signal regions (SR) are designed below (x = 1, 2).
• "eexj" signal region: 2 electron candidates with opposite charge sign and x jets.
• "emuxj" signal region: 1 electron candidate and 1 muon candidate with opposite charge sign and x jets.
• "etaxj" signal region: 1 electron candidate and 1 τ had candidates with opposite charge sign and x jets.
• "mumuxj" signal region: 2 muon candidates with opposite charge sign and x jets.
• "mutaxj" signal region: 1 muon candidate and 1 τ had candidate with opposite charge sign and x jets.
• "tataxj" signal region: 2 τ had candidate and x jets. For the signal region "tata2j", the missing transverse energy 2 E miss T is required to be greater than 10 GeV to suppress the multi-jet background.
If we have real data samples, we can apply further requirements to suppress the background due to the hadronic Pb-Pb interactions: 1) no additional charged tracks in the inner detector and no additional energy deposit from neutral particles in the calorimeters as used in Ref. [79,80]; 2) no significant energy deposit in the Zero Degree Calorimeters [81] for the ATLAS detector ( because the Pb nuclei are intact in UPCs) as used in Ref.
[82]; 3) adoption of multi-variate-analysis (MVA) methods [86,87] to separate signal events from background events. In Appendix A, it is demonstrated that the signal purity in Pb-Pb collisions is much better than that in p-p collisions.
The signal is identified as a peak structure in the invariant mass spectrum of a lepton and a jet candidates. For the events with only one jet, the jet is used twice to reconstruct the LQ mass. For the events with two jets, the best alignment is determined by minimizing their relative mass difference, namely, The subscript associated with the lepton candidate represents the charge sign, namely, l 1,2 = l ± , in the signal regions with the same flavour leptons (eexj, mumuxj and tataxj).
Otherwise the subscript represents the lepton flavour. For the etaxj and mutaxj signal regions, we assume that the missing transverse energy is due to the neutrinos from the tau decays and thus added to the τ had candidate, namely, where φ miss is the azimuthal angle of the missing transverse energy; and the collinear approximation, which assumes that the flight direction of the invisible neutrino(s) in a tau decay is the same as that of the visible tau decay products [83,84], is used in Eq. 5.4. Assuming a luminosity of 1 fb −1 , Fig. 2-7 show the distribution of the invariant mass of a lepton and a jet candidates in all signal regions. We can see that the S 3 resonance is well reconstructed, which also confirms that the resonant production of a S 3 pair dominates. Other features in these distributions are also expected: 1) the mass resolution of m(ej) and m(µj) is better than that of m(τ j); 2) The resolution of m(τ had j) is even worse because it is only a partial reconstruction of the S 3 resonance; 3) the background level in the emu2j signal region is quite low because lepton flavour is conserved in the SM.
In order to improve the sensitivity to the mode LQ → τ q, two mass reconstruction methods are tried in the tata2j signal region. One method is using the collinear approximation (technical details can be found in Ref. [83,84]), the other method is from Ref. [84] which reconstructs the mass by sampling the neutrino(s)' momentum distribution in the center-of-mass frame of the tau lepton. The two mass variables, m coll and m xlg , are shown in Fig. 8. The latter mass m xlg has a better resolution and accurately reconstructs the S 3 mass. Thus it is used as the final observable for the tata2j signal region.

Systematical uncertainties
The systematical uncertainties can be categorized into experimental uncertainties and theoretical uncertainties. Here is a summary of experimental uncertainties based on the ATLAS detector in Run 2 and main theoretical uncertainties.   Figure 2.
The distribution of the invariant mass of a lepton candidate and a jet candidate in the signal region "ee1j" (Top) and "ee2j" (Bottom). The left one is m(e + j) and the right one is m(e − j). The black points in the upper pads show the asimov data which is just the sum of all background events with uncertainties completely suppressed for better illustration. The ratio of the asimov data and the total background is shown in the lower pads, where the error bars on the black points represent the expected data uncertainty while the green bands represent the total MC statistical uncertainty.
• Object triggering, reconstruction and identification efficiencies: the performance of final object triggering, reconstruction and identification can found in Ref.    Figure 3.
The distribution of the invariant mass of a lepton candidate and a jet candidate in the signal region "emu1j" (Top) and "emu2j" (Bottom). The left one is m(ej) and the right one is m(µj). The black points in the upper pads show the asimov data which is just the sum of all background events with uncertainties completely suppressed for better illustration. The ratio of the asimov data and the total background is shown in the lower pads, where the error bars on the black points represent the expected data uncertainty while the green bands represent the total MC statistical uncertainty. uncertainty is about the probability of not having a hadronic interaction, P 0had (b). As emphasized in the end of Sec. 4.1, the different choices of P 0had (b) show a photon flux uncertainty of 10-15% [72]. More discussions about these theoretic factors can be found in Ref. [72,73]. Photon flux uncertainty is directly related to the uncertainty of the signal production cross section. Here we present a data-driven estimation of this uncertainty. Because we have found a robust mass bound m LQ > 80 GeV in Sec. 2.1 and the irreducible background with m(lj) < 80 GeV share the same final objects as  Figure 4. The distribution of the invariant mass of a lepton candidate and a jet candidate in the signal region "eta1j"(Top) and "eta2j"(Bottom). The left one is m(ej) and the right one is m(τ j). The black points in the upper pads show the asimov data which is just the sum of all background events with uncertainties completely suppressed for better illustration. The ratio of the asimov data and the total background is shown in the lower pads, where the error bars on the black points represent the expected data uncertainty while the green bands represent the total MC statistical uncertainty.
the signal, thus the region with m(lj) < 80 GeV can be used as a control region to constrain the photon flux uncertainty.
• Uncertainty of the background due to Pb-Pb hadronic interactions: the MC simulation of this background is unavailable and we have assumed it is negligible in the signal regions. One of the key features of photon-photon interactions in UPCs is that the net transverse momentum of the final state is very small. This feature can be used to estimate the Pb-Pb hadronic background in the signal regions "ee2j", "emu2j", and "mumu2j" where no neutrinos are in the signal events. In other signal regions, we may use empirical functions to model it based on the fact that the background shape is expected to be smooth compared to the signal resonance structure.
• Jet flavour: no attempt is made to identify the jet flavour. The jet reconstruction efficiency may differ for different flavours. This can be estimated by MC and the systematical uncertainty can be controlled by studying the process of γγ → qq. A more detailed analysis by categorizing the events into signal regions with light jets  Events / 10 GeV  Figure 5. The distribution of the invariant mass of a lepton candidate and a jet candidate in the signal region "mumu1j"(Top) and "mumu2j"(Bottom). The left one is m(µ + j) and the right one is m(µ − j). The black points in the upper pads show the asimov data which is just the sum of all background events with uncertainties completely suppressed for better illustration. The ratio of the asimov data and the total background is shown in the lower pads, where the error bars on the black points represent the expected data uncertainty while the green bands represent the total MC statistical uncertainty.
(from u, d, s quarks), c-tagged jets and b-tagged jets will surely improve the sensitivity. But the methodology is the same.

Statistical interpretation and sensitivity results
A simultaneous binned likelihood fit is performed to the distributions shown in Figs. 2-8. The likelihood function is where the index i represents the i-th signal region; the index j represent the j-th bin; the index k represent the k-th signal mode; N i bins represents the number of bin in the i-th signal region; N i sig represent the number of signal modes in the i-th signal region; and s i j,k (B ex , B µx ) represents the yield contribution of the k-th signal mode to the j-th bin in the i-th signal region. The yield of each signal mode should be associated to the   Figure 6. The distribution of the invariant mass of a lepton candidate and a jet candidate in the signal region "muta1j"(Top) and "muta2j"(Bottom). The left one is m(µj) and the right one is m(τ j). The black points in the upper pads show the asimov data which is just the sum of all background events with uncertainties completely suppressed for better illustration. The ratio of the asimov data and the total background is shown in the lower pads, where the error bars on the black points represent the expected data uncertainty while the green bands represent the total MC statistical uncertainty. Table 3. Expected sensitivity to search for S 3 with a mass of 100 GeV using 1 fb −1 of Pb-Pb UPCs at √ s = 5.02 TeV. Allowing B ex and B µx to float, the left plot of Fig. 9 shows the −2∆ ln L(κ) ≡ −2 ln L(κ) Lmax as a function of κ assuming the absence of new physics. We can see that κ > 1.1 × 10 −2 can be excluded at 95 % confidence level for a dataset with the luminosity L = 1 fb −1 . To exclude κ > 1, a dataset of L ≈ 4.0 pb −1 is needed as shown in right plot of Fig. 9. The sensitivity results are summarized in Table 3.   Figure 7. The distribution of the invariant mass of a lepton candidate and a jet candidate in the signal region "tata1j"(Top) and "tata2j"(Bottom). The left one is m(τ + had j) and the right one is m(τ − had j). The black points in the upper pads show the asimov data which is just the sum of all background events with uncertainties completely suppressed for better illustration. The ratio of the asimov data and the total background is shown in the lower pads, where the error bars on the black points represent the expected data uncertainty while the green bands represent the total MC statistical uncertainty.

Can we exclude a 100 GeV LQ unambiguously?
In this section, we allow a LQ to decay to all possible fermions and investigate the possibility of a 100 GeV LQ with an electric charge 4 3 e considering all available experimental constraints. The results from the H1 experiment [43][44][45], CDF experiment [46][47][48], and ATLAS experiment [56] are used. We do not include other experiments which did not look at this mass region although they may potentially provide stronger constraints. Table I summarizes the available constraints with the orignal interpretation with the corresponding assumption. We need to reinterpret them to remove the assumption. If an experiment does not explicitly identify the flavour of a jet, we will always include either all light downtype quarks (denoted by q = d, s) or all down-type quarks (denoted by Q = d, s, b) in the reinterpretation.
In the H1 experiment, a LQ is produced via the single production mode e+p → LQ+X. For LQ → lq, the final yield is proportional to λ 2 ed B lq (we ignore the possible interference  Figure 8. The distribution of the reconstructed mass of a lepton candidate and a jet candidate in the signal region "tata2j" based on the collinear approximation method (Top) and the method in Ref. [84] (Bottom). The left one is m(τ + j) and the right one is m(τ − j). The black points in the upper pads show the asimov data which is just the sum of all background events with uncertainties completely suppressed for better illustration. The ratio of the asimov data and the total background is shown in the lower pads, where the error bars on the black points represent the expected data uncertainty while the green bands represent the total MC statistical uncertainty.    Fig. 3 in [45] λ ed 0.006 (B ed = 1) λ 2 ed B eq (0.006) 2 CDF [46] Table I in [46] σ < 10 pb (B τ Q = 1) σ th B 2 τ Q < 10 pb, σ th = 9.8 pb CDF [47] Table I in [47] σ < 1.1 pb (B eq = 1) σ th B 2 eq < 1.1 pb, σ th = 11.13 pb CDF [48] Table II in [48] σ < 1.35 pb (B µq = 1) σ th B 2 µq < 1.35 pb, σ th = 11.13 pb ATLAS [56] Fig. 4 in [56] Exclude 200-800 GeV (SG) see text with the SM process). The constraints are then re-expressed as which are from Ref. [43][44][45], respectively. In the CDF experiment, a pair of LQs is produced in the reaction of pp → LQLQ + X. The final yield is proportional to the product of the pair-production cross section and the branching fractions. Using the cross section from the theory, we can obtain the constraints on the branching fractions, namely, which are from Ref. [46][47][48], respectively. Based on the constraints above, it is impossible to exclude a 100 GeV LQ definitely because all couplings are not covered. The ATLAS and CMS experiments have also performed comprehensive searches for various LQs via either single-or pair-production mechanism. However, we find only Ref. [56] from ATLAS has a public dataset which can be used to probe the existence of a 100 GeV LQ with an electric charge 4 3 e. The original motivation of Ref. [56] is to search for a LQ decaying to tτ /bν or bτ /tν. We will use their data with the final state τ lep τ had plus one or more b-tagged jets, namely, Fig. 4 (middle left) from Ref. [56] and Fig. 6b, Fig. 6c and Fig. 7b in its auxiliary material [96]. The data with the τ had τ had final state is not used due to the complicated trigger requirements which depend upon the real-time data-taking status. Figure 10 shows the distribution of the pair mass of a tau-lepton candidate and its matched jet candidate in the signal region with one b-tagged jets (denoted by τ lep τ had 1j b ) and signal region with two b-tagged jets (denoted by τ lep τ had 2j b ). The matching criterion is to minimize the mass difference of the LQ pair candidate (see Ref. [56]).
We allow all possible decay modes, and do not differentiate e and µ or d and s. We thus have four decay modes, LQ → bτ, bl, qτ, ql where l = e, µ and q = d, s and 10 individual signal topologies from a LQ pair. They are also shown in Fig. 10. In preparing signal samples, we have tried our best to apply the same selection criteria as described in Ref. [56]. Given that no significant deviation from the SM prediction is seen and following the methodology presented in this work, the upper limit at 95 % C.L. of the pair-production cross section modifier κ is determined to be 0.83 as shown in Fig. 11, for which the CDF constraints in Eq. 8.5-8.6 have been considered. It means that a LQ with a mass 100 GeV is excluded at 95 % C.L. under the general condition that a LQ can decay to all possible fermions.
It is worthwhile to make two comments: 1) the sensitivity is dominated by the modes qτ qτ and qτ ql (marked as squares and open circles in Fig. 10) for which we should design a signal region without any b-tagged jet in reality; 2) Ref. [56] claims to exclude a LQ up to 1030 GeV assuming B bτ = 1. It is likely that the exclusion region would shrink much dropping the assumption. Given the observation in this section, this work is better seen as a proof of concept in order to search for LQs in a least model-dependent way.     Figure 10. The distribution of the pair mass of a tau-lepton candidate and its mass-matched jet candidate in the signal region with one b-tagged jet (Top) and two b-tagged jets (Bottom). The left one is m(lj) and the right one is m(τ j). The black points represent the data from Ref. [56,96]. The open circles, black squares and histograms represent different signal topologies for a 100 GeV LQ.

Summary
The possibility of low-mass leptoquarks cannot be excluded due to the assumptions made in previous experimental measurements. In this work, we propose a least model-dependent method (assuming λ 0.1) to search for them in the Pb-Pb UPCs thanks to the high photon flux and potentially low background contamination. It allows a LQ to couple to all possible lepton-quark pairs. Taking the scalar LQ S 3 as example, the feasibility is investigated using all llqq final states based on the Pb-Pb UPCs at √ s = 5.02 TeV and the performance of the ATLAS detector in Run 2. The LQ with a mass of 100 GeV can be excluded at 95 % confidence level using a dataset of 4.0 pb −1 . This proposal uses the property that LQs carry electric charge and is complementary to the searches in p-p collisions. The method also works for searching for high-mass LQs in p-p collisions as long as the pair production mechanism dominates. Acknowledgments I am greatly indebted to David d'Enterria for support with MadGraph and the theoretical discussions about UPC. I also wish to express my appreciation to Prof. Shan Jin for helping to resume my research career. Especially I would like to thank Fang Dai for her encouragement and financial support.

A Signal purity comparison between in Pb-Pb collisions and in p-p collisions
In this appendix, we compare the signal purity between in Pb-Pb collisions at 5.02 TeV and in p-p collisions at 13 TeV using the bbτ τ final state (an experimental search see Ref. [34]). The final objects are reconstructed in the same way as described in Sec. 5. The b-jet is reconstructed at a working point with an efficiency of 80 %. The events with one charged lepton candidate (e or µ), one hadronic tau decay candidate and one b-jet candidate are selected. The luminosity is 1 fb −1 for Pb-Pb collisions and 36 fb −1 for p-p collisions. Figure 12 shows the distribution of the invariant mass of a hadronic tau candidate and a b-jet candidate. Both the production cross section in Pb-Pb collisions and in p-p collisions are scaled by a factor 10 −2 for better illustration (not including the scaling factor as shown in the legend in Fig. 12). We can see that the signal-to-background ratio in Pb-Pb collisions is around 50 times higher than that in p-p collisions. The reason is that the tt background dominates in p-p collisions while it is negligible in Pb-Pb collisions (the cross section is about 1 % that of bbτ τ production in the SM). This is one of the advantages of searching for low-mass LQs in Pb-Pb UPCs.  Figure 12.
The distribution of the invariant mass of a hadronic tau candidate and a b-jet candidate in the Pb-Pb collisions (L) and the p-p collisions (R). The signal histograms are scaled by a factor of 10 −2 besides the scaling factor in the legends. The black points in the upper pads show the asimov data which is just the sum of all background events with uncertainties completely suppressed for better illustration. The ratio of the asimov data and the total background is shown in the lower pads, where the error bars on the black points represent the expected data uncertainty while the green bands represent the total MC statistical uncertainty.