Collider production of Electroweak resonances from photon-photon states

We estimate production cross sections for 2-body resonances of the Electroweak Symmetry Breaking sector (in $W_LW_L$ and $Z_LZ_L$ rescattering) from $\gamma\gamma$ scattering. We employ unitarized Higgs Effective Field Theory amplitudes previously computed coupling the two photon channel to the EWSBS. We work in the Effective Photon Approximation and examine both $e^-e^+$ collisions at energies of order 1-2 TeV (as relevant for future lepton machines) and $pp$ collisions at LHC energies. Dynamically generating a spin-0 resonance around 1.5 TeV (by appropriately choosing the parameters of the effective theory) we find that the differential cross section per unit $s$, $p_t^2$ is of order 0.01 fbarn/TeV$^4$ at the LHC. Injecting a spin-2 resonance around 2 TeV we find an additional factor 100 suppression for $p_t$ up to 200 GeV. The very small cross sections put these $\gamma\gamma$ processes, though very clean, out of reach of immediate future searches.


Introduction
Accelerator-based particle physics is making progress in the exploration of the TeV energy range at the LHC. At a minimum, one may make headway in understanding the sector of the Standard Model (SM) responsible for Electroweak Symmetry Breaking (EWSBS), composed of the new Higgs boson h and the longitudinal components of gauge boson pairs W L W L and Z L Z L . These are equivalent to the ω a Goldstone bosons of electroweak symmetry breaking, in the sense of the Equivalence Theorem [1]. Under its hypothesis, that the energy of longitudinal gauge boson scattering is large E 2 ωω = s ωω M 2 W , M 2 Z , m 2 h , the scattering amplitudes involving the W L and Z L (that come to dominate W and Z scattering anyway at high energy) can be exchanged for the scattering amplitudes of the scalar ω a . Employing the latter is advantageous because of the absence of spin complications and because many of their couplings are related, in a transparent manner, by the pattern of symmetry breaking, SU (2) L × SU (2) R → SU (2) c .
Much of the LHC strategy so far has focused on hard collisions, with multiple tracks in the central rapidity region of the detectors, triggering for various high-p t (transverse momentum) scenarios. To reduce noise produced by hadron remainders, and also to directly access quartic gauge couplings, the isolation of γγ initiated events is an interesting additional alley of investigation.
In fact, run-I of the LHC has already found some events corresponding to the reaction γγ → W + W − , initially with low p t below 100 GeV [2], and now up to 200-300 GeV [3]. This later publication presents marginal (3.4σ) evidence with approximately 20 inverse femtobarn of integrated luminosity taken at 7 and 8 TeV in pp collisions. They have a total of 15 reconstructed events in both sets of data (with expected backgrounds summing about 5 events). The data is used to constrain coefficients of the linear realization of the Standard Model Effective Theory (SMEFT), following earlier Tevatron studies [5], but not the nonlinear Higgs EFT (HEFT) that we employ.
Encouraged by this success, CMS and Totem have joined [6] into the CMS-Totem Precision Proton Spectrometer (CTPPS) that will employ the LHC bending magnets to curve the trajectory of slightly deflected protons and detect them off-beam. The ATLAS collaboration is also working in at least two subprojects [7], AFP and ALFA, that allow to identify one or even the two elastically scattered protons a couple hundred meters down the beampipe from the main detector. Tagging of the outgoing protons with these detectors will allow rather exclusive measurements, among others, of γγ initiated reactions, efficiently exploiting the LHC as a photon-photon collider.
Meanwhile, a new generation of e − e + colliders is in very advanced design stages. CLIC [8] and the ILC [9] would naturally run in the 350-500 GeV region (just above the tt threshold, but in a second stage they could reach up to 1-5 to 3 TeV (CLIC) and 1 TeV (ILC) which would allow many interesting new physics studies with W W pairs [10]. The lepton colliders can also easily be adapted to perform γγ physics, and LEP was indeed used this way [11].
Therefore, it is sensible to carry out theoretical studies of the EWSBS in photon-photon collisions since the experimental prospects are reasonably good. Since no clear direction for new physics searches is emerging yet from the LHC [12], there has been a revival of the electroweak chiral Lagrangian -now including an explicit Higgs boson, in what has been called [13] the Higgs Effective Field Theory (HEFT)-and other effective theory formulations.
HEFT is valid to about 4πv 3 TeV (or 4πf in the presence of a new physics scale such as in Composite Higgs Models). Because we use the Equivalence Theorem that requires high energies, we address the 500 GeV-3 TeV region (other groups have examined the lower-energy γγ production of new resonances). In this energy range, m h is negligible, and we thus consistently neglect the Higgspotential self-couplings of order m 2 h . Except for this small assumption, a feature of many BSM (Beyond the Standard Model) approaches, our setup is rather encompassing, as several BSM theories may be cast, at moderate energy, in HEFT. Several groups [14][15][16][17][18][19][20][21][22] have studied in detail this EFT and its derived scattering amplitudes. Since those EFTs violate unitarity (see subsec. 2.2 below for a summary), we [23][24][25] and others [14,[26][27][28][29][30] have pursued methods of unitarization that are sensible in the resonance region.
In a recent contribution [31] we have coupled the EWSBS, well studied in HEFT+unitarity in that body of work, to the γγ channel. The motivation is clear: now we are prepared to address the production cross section of ωω bosons via γγ intermediate states. That is the thrust of the present document.
The electric field of a fast charge is Lorentz contracted in the longitudinal direction and thus practically transverse, appearing as an electromagnetic wave travelling parallel to the particle's momentum, as observed by Fermi [32]; the theory was further developed by Weizsäcker and Williams [33,34] (at a classical level) while Pomeranchuk and Shmushkevitch [35] offered a consistent covariant formulation. The resulting "Equivalent Photon Approximation" whereby the moving charge is accompanied by a quantized radiation field is reviewed and detailed in [36,37], from which we will draw all needed material.
Because we are working under kinematic conditions that make the Equivalence Theorem a good approximation, throughout the article we will use interchangeably the notations W L W L and ωω for the charged, longitudinal gauge bosons and Z L Z L or zz for the neutral ones, computing all amplitudes in terms of the Goldstone bosons.

Partial waves in perturbation theory
The lowest-order γγ partial waves that do not vanish (which we denote by a (0) superindex) are given next in Eq. (1). They are Next to Leading Order (NLO) for J = 0 while Leading Order (LO) suffices for J = 2.
We obtained them in terms of the fine structure constant α = e 2 /4π and the parameters of the EWSBS (that the LHC is constraining) in [31], from earlier work on the effective Lagrangian and the invariant amplitude involving two photons in [38]. They read where the combinations A C and A N refer to the charged basis W + W − and ZZ, which here appear mixed because we employ the custodial isospin basis that characterizes the final state, since the photon coupling is isospin violating and can yield both I = 0 and I = 2. I = 1 is discarded because the ωω state must be Bose symmetric, entailing J = 1, and the γγ state cannot be arranged with one unit of angular momentum as per Landau-Yang's theorem. A C and A N can be written as For completeness, let us quote also the scalar partial wave yielding the scalarisoscalar hh final state, which only couples with positive parity states The scalar partial waves P (0) I0 at this order, and all waves at higher orders, grow polynomially with Mandelstam s according to the chiral counting, if there is BSM physics in the EWSBS, until the new scale of that physics is approached. Therefore, chiral perturbation theory (ChPT) eventually breaks down; the amplitudes can still be represented from first principles (unitarity and causality) by a dispersive analysis, with chiral perturbation theory supplying the low-energy behavior (subtraction constants for the dispersion relations) which gives rise to the well-known unitarized EFT. In the next subsection we quickly recall the application of this unitarization to amplitudes involving two photons.
If no new physics is within reach at the LHC, the corresponding SM expressions are a = 1, c γ = a i = 0, b = a 2 and thus R (0) 0 = 0, as well as A N = A C = 0, so that P

Unitarity and resonances
In this article we do not consider the final hh state, and for simplicity we also assume that it is decoupling from ωω W L W L so we set a 2 = b (as well as the other parameters coupling both channels, d = e = 0).
The scattering amplitude linking ωω and γγ is then a three by three matrix [31] due to custodial isospin. The two-photon state can couple to both I = 0, 2 breaking custodial symmetry, though the presumed BSM interactions in the BSM do not connect the two channels. For each of them, angular momentum can be 0 or 2. This matrix is where the A IJ (s) are the elastic partial waves ωω → ωω from [24,39], and the P IJ (s) photon-photon amplitudes are taken from subsec. 2.1. The two zeroes in the upper left box encode isospin symmetry in the EWSBS; the zero in the lower right corner arises because we work at LO in α, so that γγ F (0) γγ 0. The unitarity condition for this matrix amplitude is not satisfied by the perturbative amplitude because of the derivative couplings growing with s, so unitarization is needed. But since α is a small parameter, it can be taken at leading order. Then, Eq. (6) can be satisfied, in very good approximation, to all orders in s but only to LO in α. Substituting Eq. (5) in Eq. (6) yields In the second equation, the γγ → γγ amplitude has been neglected as it would exceed first order in the α expansion. The elastic ωω → ωω amplitude may be expanded in the HEFT (as recounted in [24]) by This amplitude violates exact elastic unitarity |A| 2 = ImA, satisfying it only in perturbation theory |A (0) | 2 = ImA (1) , which is an important handicap of EFTs and leads to large separations from data at mid-energy (few-hundred MeV above threshold) in hadronic physics. However, if it is employed as the low-energy limit of aÃ satisfying exact unitarity and obtained from dispersion relations, it gives rise to successful methods (such as the IAM, N/D, Improved-K matrix, large-N unitarization, etc.). These methods differ in numerical accuracy but not in substance [23,24], as they all reproduce the same resonances in each elastic IJ channel for similar values of the chiral parameters. The P amplitudes, by Watson's theorem, need to have the same phase asÃ due to strong rescattering. This we guarantee by satisfying Eq. (7). Observing that at low energies, P ≈ P (0) , and enforcing the correct analytical structure in the complex s plane, we proposed [31] the following unitarization method for the γγ → ωω scalar amplitudes, which implements the IAM philosophy; here,Ã(s) = A (0) (s)/(1 − A (1) (s) A (0) (s) ) is the elastic IAM. Now, for J = 2, the IAM cannot be employed, and then we resort to the well-known N/D method (we have also checked that employing the N/D for both J = 0 and J = 2 leads to little material difference). Then, a formula similar to Eq. (9) can be used Here, the N/D elastic amplitude has been employed; this is somewhat more complicated than the IAM, and requires giving further detail on Eq. (8), as the quantities are built from the B and the D factors defined by These are computed in perturbation theory and have been reported earlier in [24]. The amplitudes are µ-independent because B(µ) runs in such a way as to absorb the dependence coming from the logarithms.

Invariant amplitude and differential cross section
The non-vanishing matrix elements can be reconstructed from the (unitarized) partial waves bỹ where I ∈ {0, 2}.T I0 andR 00 are related with the positive parity state (|++ + |−− )/ √ 2 by means of the definitioñ Since we have 4 possible γγ initial states, the differential cross section for γγ → ωω will be And, for γγ → hh, In implementing these two equations, which are a backbone of the computation, we have employed the Inverse Amplitude Method extension in equation (9) for the J = 0 channels, as is it is the one which has been more extensively studied in low-energy chiral perturbation theory and its uncertainties are well understood. For the J = 2 resonances, the Inverse Amplitude Method cannot be used as a parametrization as it would require knowing the NNLO amplitude in the HEFT. As this is not the case, we have compromised and used the N/D method as laid out in Eq. (11).
By using the change of basis from the ωω isospin one, |I, M I , to the charge one, {|ω + ω − , |ω − ω + , |zz }, and taking into account that γγ states do not couple with J = 1 gammagamma states, the unpolarized γγ → {ω + ω − , zz} differential cross section can be written as If we take the SM limit as laid out at the end of subsec. 2.1, we find 3 , respectively. A seeming puzzle with this expression is that the tree-level perturbative expression for γγ → ππ (discussed at length in chiral perturbation theory in [40]), a pure scalar electrodynamics result, is given by which is independent of the polar angle, and does not contain any factor |Y 2 2 | 2 . This difference is an artifact of our partial wave expansion: if we wanted to recover the Born-like result of Eq. (20) we would need to resum the partial wave series. For example, the first few P 0J with even J = 2 . . . 12 are α/(6 √ 2) (given in Eq. (1a)), α/(6 √ 30), α/(6 √ 140), α/(6 √ 420), α/(6 √ 990), α/(6 √ 2002), and the first few P 2J are α/12 (given in Eq. (1b)), α/(12 √ 15), α/(12 √ 70), α/(12 √ 210), α/(12 √ 495), α/(12 √ 1001). Each of these quantities multiplies the corresponding spherical harmonic in reconstructing the perturbative amplitude. The series is well behaved for any fixed angle θ, but in truncating it, we introduce a spurious angle dependence. We have not pursued the issue further since our aim is not to present precise off-resonance cross-sections for production of the EWSBS particles; this can be best computed by standard means (Feynman amplitudes not expanded in J). Both methods can also work together and part of us have recently assessed it, in a separate collaboration [41], to implement in LHC Monte Carlo simulations.
Our goal here is to produce the resonance cross-sections; and near a BSM resonance, the dominance of its corresponding partial wave over all the other, perturbative ones, is warranted in the presence of experimental angular acceptance cuts that avoid any forward Coulomb divergence. Thus, in the figures that follow, one should pay attention to the differential cross-sections near the peak, and not take too seriously the background cross-sections that are affected by factors of order 1. The effect is lesser in directions perpendicular to the beam axis (low rapidity).

Inverse process ωω → γγ
As an aside, and for completeness, we also give expressions for the process ωω → γγ (and for hh → γγ) that may be useful in the study of resonances decaying by the two photon channel. Assuming time reversal invariance i | T | j = j | T | i , and considering that we have (2+1) 2 = 9 possible initial states |I, I z , we obtain Finally, since there is only 1 possible initial state, ω + ω − → γγ and zz → γγ can be written as The aim of this section is to study the differential cross section dσ dsγγ dp 2 T to photoproduce pairs of longitudinal W L electroweak bosons in e − e + → e − e + + γγ → e − e + + W L W L at an energy of 1 TeV, the top of the energy range of the International Linear Collider and above. dσ dsγγ dp 2 T is obtained through the convolution of photon flux functions derived from the Equivalent Photon Approximation and the γγ → ωω cross section described in section 2. Fig. 1 shows the characteristic Feynman diagram to be evaluated.
We will work in the ee center of mass frame so that p e − = (E; 0, 0, E), p e + = (E; 0, 0, −E) and Mandelstam's invariant squared energy will be s ee = 4E 2 . The only exception to this massless-electron kinematics will be in the regulation of x-integrals such as Eq. (26) below; thus, we actually work in the leading me √ s approximation which leads to finite answers. In the Equivalent Photon Approximation [36], the charged leptons can radiate collinear photons [42] (since their boosted Coulomb field is practically transverse to the lepton direction of motion) so that we may take the photon momenta as p γ1 = (ω 1 ; 0, 0, ω 1 ) and p γ2 = (ω 2 ; 0, 0, −ω 2 ); the corresponding invariant is s γγ = s W L W L = 4ω 1 ω 2 . If these two photons do enter into the EWSBS through a resonance, this is not produced in its rest frame as their momenta are not opposite in the laboratory. Instead, each photon carries a different fraction of its parent lepton momentum, Substituting these two into s γγ and eliminating E for s ee we find the constraint The variables x and y are bound above (by the maximum energy available from the electron), Then Eq. (24) also gives a lower bound so that integration over the photon momentum fractions never hits the end points and is regular. Neglecting all masses and photon virtualities we can interpret the Mandelstam γγ variables in the center of mass of the γγ → ωω subsystem as usual and can trade t γγ (the variable in terms of which our Feynman amplitudes are expressed) for the more directly measurable p 2 t as Then, it immediately follows that dp 2 so the Ω angular dependence of the cross-section can be traded for one in p 2 t . Then, Eq. (16) becomes dσ γγ→ωω dp 2 The photon virtualities are also bound [36] and these bounds limit the interval of validity of the photon number density per unit energy and virtuality, Eq. (D.4) from [36], where C and D are two constants that parametrize the internal structure of the charged particle (and, as usual, Q 2 = −q 2 ). In the case of pointlike elementary particles such as e − e + , C = D = 1, the photon flux can be integrated over virtuality to interpret it in a manner analogous to a parton distribution function, The integral over Q 2 can be performed analytically, yielding which is represented in figure 2. With the photon fluxes at hand, we may now mount the cross section for the ωω production process, by means of (35) or, in differential form,

Some numerical examples
We exemplify Eq. (36) with a set of parameters characteristic of the EWSBS in the presence of new physics. For simplicity we will decouple the hh channel setting b a 2 . We keep the LO parameter a = 0.81 fixed and separating from its SM value (that would be a = 1). This particular value is chosen because it is just under the 2σ recently proposed exclusion bound [43]. Those authors report a = 0.98 ± 0.09 (1σ) from current LHC data. We generate elastic ωω resonances by means of the a 4 and a 5 NLO parameters, fixing all others to zero at NLO and higher; in this way, the entire coupling to the γγ sector is provided by the electron squared charge in α em . We have chosen the sets a 4 = 10 −3 , a 5 = 0 and a 4 = 10 −3 , a 5 = 10 −3 , that have increasing BSM strength at NLO and generate resonances at decreasing s γγ . All these parameters are understood to be taken at the renormalization scale µ = 3 TeV (their running to other scales can be found in our earlier work [23,24]). They Both sets induce resonances around 1 TeV; the one at higher mass in the right plot is narrower (f 2 -like) and the one common to both plots is an f 0 -like, broader structure. We show both fixed-angle (in the CM) and fixed-p t scattering.
are basically unconstrained except for the current absence of BSM resonances. The sets we use do provide resonance in the energy region just above 1 TeV.
We have chosen p t = 50, 100, 200 GeV that would pass typical experimental cuts 1 and future e − e + machines will similarly impose p t cuts at trigger time. Both plots show similar features. The fixed-angle cross section is larger and falls slower with the energy than the fixed p t one. Resonances are however clearer at fixed p t , and their line shape is the better resolved the lighter they are (with the parameters chosen, a strong scalar resonance appears around 1 TeV). Fig. 4 shows the convolution of the cross section γγ → ωω in Fig. 3 with the photon flux factors, to yield the e − e + production cross section that can be readily obtained in experiment.
This we present in doubly differential form, respect to E e − e + and respect to p 2 t (of each produced ω). If the resonance of the EWSBS is above 1 TeV, as shown in the left plot of the figure, the resonance shape is not so distinct (especially in the presence of experimental errors), but the line shape exposes a clear knee with a change of power-law slope and is shifted to higher values after the resonance. On the other hand, a resonance below 1 TeV is more clearly visible and can be better reconstructed if p t is modest. For larger p t ∼ 200 GeV, the behavior of the line shape is similar to that of a higher energy resonance. The reader may be intrigued by the growth of the cross-section with p t for small s. One should remember that the underlying Lagrangian is Effective Field Theory-based and thus, coupled derivatively. Therefore, an increase of the transverse momentum provides larger amplitudes at the γγ → W L W L level. For larger s such that unitarity is saturated and for larger p t this effect diminishes and the usual kinematic effects lower the cross-section.

Production in pp collisions
In this section we revisit γγ → W + L W − L at hadron colliders, focusing on the LHC (pp initial state), so the complete reaction is pp → ppW + W − (through γγ).

We estimate the cross section
from that for the inelastic dσ(γγ→ωω) dp 2 t in analogy to the case of lepton colliders; the only difference is that now the photon flux is not directly calculable with arbitrary precision, but partly fit to data. Additionally, one or both of the outgoing protons might be disrupted, giving several situations: elastic (pp in the final state), semielastic (only one of the two protons remains), totally inelastic (both protons break up). Additionally, the inelasticity can be few-body, with momentum transfer in the resonance zone, or deeply inelastic, with large Q 2 and accepting a parton-level description.
Theoretical studies of photoinduced production at the LHC date more than a decade [44], at least for hadron resonances; the cross sections found for those are of course huge in comparison to the electroweak sector. Though the experimental identification of this process is rather difficult in the noisy environment of a hadron collider, taking into account that the photons are not detected and must be inferred, at least 13 events have already been reported by CMS [3]. ATLAS also sees candidate events [4] and while some are expected for the pro-posed Future Circular Collider [45], its energy is too low to be of interest for the resonance region.
Moreover, two projects, CT-PPS (CMS-TOTEM Precision Proton Spectrometer [46]) and AFP (the ATLAS Forward Proton detector [47]) aim at detecting elastically scattered protons near the beampipe. These experiments take data at 200 meters (CT-PPS) and 210 meters (AFP) from their respective interaction points. Both employ the LHC magnets to separate the scattered protons by their momentum slightly different from that of the beam, and detect them downstream.
The kinematics for the process closely parallels the discussion in section 3 through Eq. (24), just substituting e → p as necessary.
Given the photon flux in the proton beam as, once more, f (x), and in the collinear photon approximation, we may write [48] where in the second step we have used the relation analogous to Eq. (24). If the proton is left intact (elastic photon emission) then the flux factor f (x) is calculable from the electromagnetic form factor of the proton. In the deep inelastic regime, we can speak of the photon as a parton of the proton; and in the intermediate region, the proton is left in an excited state (one of several resonances), f (x) then being a nontrivial structure function.
As the photons are collinear with the proton, the angular dependence of the W W emission comes from the reaction γγ → W W , as in e − e + collisions. Again, p 2 t is, unlike Ω, invariant under longitudinal boosts (and is easily measured) so we take it as second variable and write dσ γγ→ω + ω − dp 2 where the first factor stems from the variable change t → p 2 t and the second from Ω → t.
We may then write Eq. (37) in double differential form,

Photon flux in the proton
We need to convolute photon-level cross-sections with the collinear photon flux in the proton, with f γ|p (E γ ) ≡ f (x) (x = E γ /E p ) computed under two kinematic regimes which can be distinguished by experimental triggers.
In the first, we take the absorption cross-section for real photons as not too different from that for virtual photons (small virtuality). On top of this approximation, there is the mild assumption that the cross-section must fall quickly after a certain energy. Then one can find, for elastically scattered protons, an expression in terms of the Sachs electromagnetic form factors G E and G M [36,49] f with lower integration limit Q 2 min = (M p x) 2 /(1 − x). Alternatively, a second kinematic regime is the deep inelastic kinematics; f (x) is directly taken as the parton distribution function. We now show parametrizations of both elastic and deeply inelastic photon fluxes.
First to mention is the very crude parametrization employed by Drees and Zeppenfeld [50] and also recently adopted in [51]. In this high-energy application, the authors neglect the lower limit Q 2 min , and the magnetic form factor G M (Q 2 ), and parametrize the electric one G E (Q 2 ) by a simple dipole form We will plot the resulting photon flux in Fig. 6, and further include (a) the simple improvement of [44] that considers the minimum Q 2 and (b) the parametrization of Kniehl [49] that now includes both G E and G M at the proton-photon vertex. In addition to those classic works, we will also work with more modern parametrizations that reflect progress in hadron physics in the last two decades and help better characterize systematic errors. We will try the low-energy parametrization of Lorenz and Meißner [52] based on a conformal coordinate change from Q 2 to z, with t cut = 4m 2 π the charged pion-pair threshold. This allows a Taylor power series expansion in z of the form factors in terms of free a k constants that were adjusted to low-Q 2 data (up to about 1 GeV 2 ). The expansion converges rather well. This low energy form should provide a very accurate fit only below that scale, but since the form factor is integrated in Eq. (41) to obtain the flux, contributing all the way to Q 2 6 GeV, we need to supplement this parametrization with a high-energy contribution. Thus, for Q 2 > 1 GeV we adopt the simplest Brodsky-Lepage [53] form factor that follows the power-law counting of QCD [54], that yields for large Q 2 where Λ = 200 MeV and C 2 is a constant that we use to match continuously with the low-Q 2 parametrization; this we do at Q 2 match = 850 GeV 2 . As for the electric form factor, in this simple parametrization it is obtained from G E = G M /µ p . Fig. 5 represents the form factor obtained by matching these low-energy data fit and asymptotic parametrization. The use of this form factor does not significantly change the results obtained with those of [49,50], and [44]. Finally, we also include in the figure two more contemporary parametrizations of these form factors. One is the data-oriented fit of Kelly [55] and the other, a theoretical computation by Segovia et al. based on the Dyson-Schwinger equations [56].
We now turn to the deeply inelastic cross sections. We have at our disposal several different photon distribution functions in the proton, published respectively by the collaborations CT14QED (or, for shortness, CTQ14 in what follows) [57], NNPDF3.0QED [58], NNPDF2.3QED [59] and MRST2004QED (or just MRST) [60]. Additionally, the LUX photon [61] "pdf" has also been included in the comparisons. Note that LUX is not only a proper pdf when the photon virtuality is large, but also an effective way of encoding the Weizsäcker-William photon flux of the proton at large energy (even at moderate and low  Two of these sets, MRST and CTQ14, are obtained with a similar analysis and their results are consequently also similar, with the difference between the two sets falling with Q 2 . Therefore we will show our results for the CTQ14 set only, with those obtained from the MRST pdf sets being numerically close. The CTQ14 collaboration has fit isolated photon production in DIS in the interval 10 GeV 2 < Q 2 < 350 GeV 2 , and we expect that the pdf parametrizations will be usable in this momentum squared range. The earlier NNPDF photon distributions were rather different from those of CTQ14 (and MRST) as can be appreciated from Fig. 7, especially so at low x. The difference might have been attributable to NNPDF excluding the direct DIS information on the photon as discussed in [62] which perhaps makes its uncertainties unnecessarily large. Newly for NNPDF3.0, this difference with CTQ14 almost vanishes on the low-x region, as can be seen in the bottom panel of Fig. 7 (where the newest set, the dotted line, is very close to the LuxQED and CTQ14 solid lines) and also in appendices B and C. In any case we employ all these sets so we can explore the systematic uncertainty in the cross section estimates.

Some numerical examples
We now have all ingredients needed to estimate proton-proton cross sections that produce W L W L or Z L Z L by means of intermediate γγ states. In Figs. 8 and 9 we put to use the elastic photon fluxes computed above and shown in Fig. 6 to compute the cross sections for TeV-EWSBS resonance production with intermediate photon states, leaving the protons unharmed.
From the figures, it appears that the cross section is small and since it increases slightly with p t , not much harm is done by imposing experimental cuts thereon that exclude low-lying quarkonia or τ τ production.
It also appears (see Fig. 9) that if a resonance would exist below 1 TeV (which we can achieve by increasing a 4 or adding a contribution from a 5 as done in the figure), the cross section would increase significantly.
As we do not find very strong signals, we need to be comprehensive and increase the kinematic range with the inelastic regime but not DIS (that is, lift any restrictions on the fate of the final state protons which we will perform in the next subsection). A very easy computation can be carried out with the DIS Figure 9: Computation with the elastic photon flux, similar to Fig. 8, but adding the a 5 NLO parameters. pdfs in Fig. 7, where both protons dissociate (there is no difficulty in computing, for example, the instances in which one proton is dispersed elastically and the other dissociated, by combining the different fluxes, all at hand). In Figs. 10 and 11 we show just this computation. The cross section obtained from the NNPDF set is quite larger than that from the CT14 one as NNPDF is the largest of the two for higher x (a 1 − 2 TeV resonance in a 13 TeV collider requires x ∼ 0.1 − 0.2).
The cross section for these DIS events can easily be 5 times larger than the elastic one, but they are very difficult to reconstruct as they can leave charged tracks in the central tracker that would not pass the cuts to reduce background.  fig. 8, that this can be up to one order of magnitude more likely than elastic scattering. The (effective) PDF energy scale is µ 2 = s γγ .
Therefore, an interesting strategy would be to search for inelastic, but not deeply inelastic, events where one or both protons are dissociated in the 1 − 3 GeV resonance region.

Inelastic regime (not necessarily DIS)
The cross sections reported so far in pp collisions, elastic and deeply inelastic are rather small, and there is small hope of measuring the later in pp because it would probably leave activity in the central barrel, so that it would not be easily identifiable over background. See [63] for extensive discussion on how to incorporate various rapidity-gap cuts that assist event identification into theoretical calculations. Actually, we expect most of the cross section not to be in those extreme regimes, but rather correspond to an intermediate, inelastic but not deeply inelastic proton recoil (in the baryon resonance energy region).
A full theory description of that region (1−2.5 GeV) is beyond our ability, as many resonances of various spins populate it and likely contribute. Therefore, we resort once more to a data-driven description, adopting a photon flux f that incorporates information from Jefferson laboratory and other mid-energy facilities. A convenient parametrization of the photon content of the proton useful for pp → γγ + X is provided by the LUX photon [61] "pdf" that is precisely the photon flux that we need 2 , describing low-Q 2 data from A1, CLAS and Hermes GD11-P. In addition to low-energy baryon resonances, the authors of [61] also incorporate the elastic form factors and DIS functions that we have examined above into their photon flux. The high-Q 2 flux 3 is a proper pdf for the photon evolved at NNLO and fit to standard data. Casting form factors and inelasticities in the language of parton distribution functions makes all the pieces fit into the standard Monte Carlo collider machinery.
Note that, in this framework, the energy scale µ 2 at which the pdf is set (and that enters into the well-known LHAPDF library [64]) differs from the virtuality of the actual γ emission process, Q 2 . This can be checked in Eq. (6) of Ref. [61], where µ 2 appears as a cutoff of the integration over Q 2 . Indeed, Fig. 2 of Ref. [61] exposes that, for µ = 100 GeV and x > 0.05, more than half of the effective pdf comes from physics at an energy scale Q 2 < (1 GeV) 2 . Hence, the requirement µ > 10 GeV of the pdf 4 is meant to limit applicability of the photon flux to collider phenomenology at a center of mass energy of s > (10 GeV) 2 (in the spirit of the Weizsäcker-Williams approximataion), but not as a limitation on the virtuality of the emission process that can be soft as in p → γ * p * . The parameter µ should be set at the scale of the large γγ scattering energy.
We have also examined an alternative work [65] that also parametrizes effective PDFs (more properly, photon fluxes to be used with the Weizsäcker-Williams approximation) via elastic (and ∆(1232)-inelastic) form factors where Q t is the photon transverse momentum, t = −(Q 2 t + x 2 m 2 p )/(1 − x), and F 1 the Dirac electromagnetic proton form factor (multiplying γ µ at the photonproton vertex). Note the curious absence of the Pauli (helicity non-conserving) form factor F 2 5 . The contribution of F 2 is included in [66], but that work is limited to the elastic contribution whereas [65] gives an analytical expression accounting for the lowest possible proton excitation, ∆(1232).
In any case, we employ these works for cross-checks and show the outcome produced with the newer and more complete LUX NNLO γ-flux [61].
Proceeding then as in Eq. (40), we obtain the cross section reported in Fig. 12. Because the LUX photon flux requires a µ scale (as it incorporates inelastic structure functions of the proton), we vary this in the graph over a reasonable range. We have taken as parameters a 2 = b = 0.81 2 , c γ = 10 −4 , and a 4 = 10 −3 , consistently with our previous sets. All the other NLO parameters from both the EWSBS and the photon sector are set to zero, namely g = d = e = a 1 = a 2 = a 3 = a 5 = 0. This set yields a typical resonance around 1.5 TeV.   [61] photon flux factor is evaluated), differential respect to the produced s W L W L and the squared transverse momentum. The NLO parameters, visible in the plots, are chosen so that a resonance in the EWSBS amplitudes is present around 1.2-1.5 TeV. After convolution with the photon flux, only a broad shoulder is visible.
The figure shows what happens to this resonance of the EWSBS after convolution with the inelastic photon fluxes: it becomes a broad shoulder, experimentally challenging after accounting for statistical data uncertainties.
In Fig. 13, we scan over c γ (top) and a 1 (bottom left), also with a 2 = b = 0.81 2 and a 4 = a 5 = 0. For completeness, we have also included a case with a 2 = b = 0.95 2 and a 4 = a 5 = 10 −3 (bottom right graph of Fig. 13), consistent with LHC constraints [75]. This set of values introduces a clear resonance at √ s ∼ 1.8 TeV that is narrower (and thus, dominated by NLO parameters).  Finally, in Fig. 14 we compare the signal with a 2 = b = 0.95 2 , a 4 = a 5 = 10 −3 , with the SM background pp → γγp * p * → W + L W − L + X. Note that γγ → ZZ vanishes at LO in the SM. The SM computations have been taken from Refs. [73,74]. Note the big backround coming from the transverse modes. However, such a background can be decreased by looking for events at high p T .
The cross sections that we find are not larger than those in subsec. 4.2, in spite of including further kinematic windows. It may be that the older parametrizations from MRST, CT14QED or NNPDF overestimated the photon flux.

Discussion and outlook
Photon-photon induced production of electroweak resonances, if they exist, is an alternative to their production from gluon-gluon interactions. While the cross sections are naturally small, as we have quantified, they are very clean if the outgoing elastically scattered protons can be tagged (see figures 8 and 9). Such searches are complementary to two-photon reconstruction in the final state of a central collision, a method that is already constraining the EWSBS below about 900 GeV [67].  Our approach to assess the EWSBS from two-photon collisions is based on the EFT formalism supplemented with unitarity. We are not able to distinguish specific models [68] as long as they are cast at low energy in the symmetry mold of the Standard Model: new physics enters through the low-energy EFT couplings only.
We have shown the Standard Model background production of W + W − from expressions available in the literature [70] and incorporated into standard Monte Carlo programs (see figure 14). The Leading Order production is easy to understand in our setup: basically, use Eq. (20) as opposed to our expressions. To take the limit continuously is less trivial as we are relying in a truncated partial wave expansion, which fails at forward/backward angles for which it is not designed, and at low energies where terms of order m h , M W are not negligible (e.g. the SM Higgs potential if it is at work). Still, taking the limit of M W → 0 of the background calculated in figure 14 and of the analytical expressions of [74] would eventually allow to match with our calculations with all parameters taken at the SM values, if the partial wave expansion is put aside and Feynman amplitudes are used for the comparison.
We have computed the elastic-elastic cross section (both protons intact), which is the cleanest experimental channel. The number of events to be found increases with p t , for modest values thereof. If a new resonance was around E γγ = 1 TeV, we have shown in Fig. 8, for example, that the cross section would be rather flat in energy and around 10 −2 fbarn/TeV 4 or somewhat more. This means that an integrated luminosity of 300 fbarn −1 at the LHC run II would prove insufficient to gather enough events at this high invariant boson-boson mass, specially when only certain diboson decay channels are experimentally reconstructed, further reducing the cross section by their branching fraction. Further small reductions are due to absorption effects [69] in the photon debris.
Thus, looking for inelastic processes to increase the cross section seems mandatory. We have shown the deeply inelastic cross sections in which both protons dissociate (figure 8), but those events are difficult to isolate because the non-photon-initiated background is too large, leaving activity in the central silicon trackers. The resonance-mediated inelastic (but not deeply inelastic) events where the proton dissociates but mostly in the forward direction are, therefore, more promising. But precise predictions are here difficult because we find quite some systematic difference due to the chosen pdf set; one can opt for the newest LUXQED set.
The situation is a bit better for resonances below 1 TeV, that may be detectable with this method as the cross sections are an order of magnitude larger. Additionally, for resonances of larger mass there may be hope in collisions involving heavy ions: for example, lead-lead induced γγ collision cross-sections are enhanced by a factor (Z = 82) 2 if collisions are incoherent which is unfortunately diminished by a factor 2000 smaller luminosity than in proton-proton collisions with the current LHC machine, so perhaps p-Pb collisions are the optimal search channel. At small momentum transfer, the entire nucleus can interact collectively and then the P b − P b reaction is enhanced by Z 4 which is more promising; but EFT interactions grow derivatively, so this strategy works only at somewhat large s or the underlying scattering amplitude is in turn decreasing the cross sections.
We have similarly predicted example cross sections for a future electronpositron collider operating in the TeV region. If we take as reference the proposed luminosity of the CLIC collider, that could conceivably accumulate about 650 fbarn −1 per year, our resonance cross sections of order 10 −3 fbarn will only yield a couple of events per year. Thus we find that, while CLIC may be apt for exploring vector resonances that couple in an s-wave to e − e + , it will fall short in luminosity to be a practical tool for scalar or tensor resonances in γγ physics.
The inclusion of all computations reported in Monte Carlo simulations of the LHC (or ILC) detectors by interested collaborations should be possible and is encouraged.
neglected [72]. L 2 and L 4 can be written as and the covariant derivative of the U field is defined as This Lagrangian leads to the Feynman rules computed in Ref. [38]. The amplitude elements A(γγ → W + L W − L , Z L Z L ) have been computed both with the spherical and linear representations of the U field discussed next in appendix A.1, and yielding the same result [38]. The unitarized partial waves can be found on section 2.1 of the present work.

A.1 Spherical (or square-root) parametrization of the coset
Here we remind the reader of two possible choices of the coset parametrization for SU (2) L × SU (2) R /SU (2) L+R . The coordinates on that coset, three Goldstone boson fields, are of course not unique but S matrix elements (on-shell amplitudes) do not depend on their choice. Very often one finds an exponential parametrization withπ = τ a π a (x) and τ a (a = 1, 2, 3) being Pauli matrices. This choice is well suited for SU (3) chiral perturbation theory with three flavors.
However, in the electroweak sector (as well as in two-flavor ChPT in QCD) the coset is just the space SU (2), isomorphic to the S 3 three-dimensional sphere. This suggests the use of simpler "spherical" coordinates: where againω = τ a ω a (x) and ω 2 = a (ω a ) 2 =ω 2 . The resulting Feynman rules and Feynman diagrams are less numerous than for the exponential parametrization and thus, calculations are a bit simpler: yet the final answers are identical to the exponential parametrization, as we showed in [38]. There, we recalled how to change between the two sets of coordinates by rewriting the exponential as where π = √ π 2 with π 2 = a (π a ) 2 and then comparing to the spherical parametrization to recover which implies ω 2 = v 2 sin 2 (π/v). An expansion (formally, in powers of π 2 /v 2 ) yields the series ω a = π a 1 − 1 6 π v 2 + 1 120 The "eaten" Goldstone bosons that provide the longitudinal components of the W ± and Z gauge bosons are then ω ± = (ω 1 ∓ iω 2 )/ √ 2, ω 0 = ω 3 (= z). The Feynman rules involving less than four WBGBs are exactly the same in both parametrizations since they differ in terms at least quadratic in the WBGBs. However the vertices with four WBGBs are indeed different in both parametrizations if the WBGBs are off-shell (but they coincide for on-shell amplitudes). The next section quotes the Feynman rules, that we do not rederive here.

B Uncertainty bands for the NNPDFs
In this brief paragraph we plot the uncertainty bands for the NNPDF sets extracting a photon from the proton; we show that, within that uncertainty, the new sets are compatible with the CT14 and LUXQED pdf sets. Figures 15  and 16 display all the sets. This gives us some confidence in their use to predict photon-initiated cross-sections; the uncertainty bands for these are shown in figure 17.

C Uncertainty bands for the CT14qed
Here, we plot the uncertainty bands for the CT14qed sets extracting a photon from the proton. The uncertainty of LUXQED has been found to be smaller than the size of the line. The CT14qed band is computed according to [57] at 90% CL. That is, the error band includes initial inelastic momentum fraction of the electron up to 0.30%. The PDF line uses 0% initial inelastic momentum fraction. Error bands for LUXqed happen to be too small for representing. Figures 18  and 19 display all the sets. Then, figure 20 propagates these uncertainties to the production cross-section.