1 Introduction

Many extensions to the Standard Model (SM) predict the existence of heavy resonances that decay into pairs of vector bosons (WW, WZ, and ZZ, collectively referred to as VV with \(V=W,Z\)). These theoretically well-motivated extensions include the two-Higgs-doublet model [1], composite Higgs models [2, 3], technicolour [4,5,6] models, and warped extra dimensions [7, 8]. The Large Hadron Collider (LHC), as the world’s highest-energy proton–proton (pp) collider, is a unique facility for the search for these heavy resonances. Indeed, both the ATLAS and CMS collaborations have reported searches for diboson resonances in various production modes and in a variety of decay final states of the vector bosons [9,10,11,12,13,14,15,16].

Depending on the assumed model, the predicted diboson resonances can be produced through gluon–gluon fusion (ggF), Drell–Yan (DY), or vector-boson fusion (VBF) processes. Representative Feynman diagrams of these processes are shown in Fig. 1.

Fig. 1
figure 1

Representative Feynman diagrams for the production of heavy resonances X with their decays into a pair of vector bosons. The hashed circles represent direct or effective couplings

This paper reports on a search for heavy resonances X in the mass range 300 GeV to 5 TeV in the \(X\rightarrow VV\) diboson decay in pp collisions at \(\sqrt{s}=13\) TeV. Three types of diboson resonances are considered in the search. The first is a neutral scalar resonance, the radion (R) [17, 18] which appears in some Randall–Sundrum (RS) models and which can decay into WW or ZZ. The second is the heavier versions of the SM W and Z bosons, \(W^\prime \) and \(Z^\prime \) bosons, as parameterised in the Heavy Vector Triplet (HVT) framework [19], which can decay through \(W^\prime \rightarrow WZ\) and \(Z^\prime \rightarrow WW\). The third diboson resonance is a spin-2 graviton (\(G_{\mathrm {KK}}\)) of the first Kaluza–Klein (KK) excitation in a bulk RS model [7, 20, 21] and decays into WW or ZZ.

Semileptonic VV final states in which one vector boson decays leptonically (\(V_\ell \): \(W\rightarrow \ell \nu \), \(Z\rightarrow \ell \ell \) or \(Z\rightarrow \nu \nu \)) while the other decays hadronically (\(V_h\): \(V \rightarrow qq\)) are considered, leading to three distinct channels: \(ZV\rightarrow \nu \nu qq\) (0-lepton ), \(WV\rightarrow \ell \nu qq\) (1-lepton ), and \(ZV\rightarrow \ell \ell qq\) (2-lepton ). Here \(\ell \) denotes either an electron (e) or a muon (\(\mu \)). The hadronic \(V \rightarrow qq\) decays are reconstructed either as two separate small-radius jets (small-R jet, or j) or as one large-radius jet (large-R jet, or J) depending on the transverse momentum (\(p_{\text {T}}\)) of the boson. The reconstructed transverse mass (\(m_{\mathrm {T}}\)) of the VV system for the 0-lepton channel and VV invariant mass (\(m_{VV}\)) for the 1-lepton and 2-lepton channels are used for signal–background discrimination via maximum-likelihood fits to their observed distributions.

Compared with previous searches in these final states [22, 23], the current search is performed with a data set approximately four times larger. Several improvements are made which include utilising a multivariate technique to identify and distinguish production processes, using tracking information in the large-R jet reconstruction, and introducing b-quark jet tagging for large-R jets.

2 Detector and data sample

The ATLAS experiment [24, 25] at the LHC is a multipurpose particle detector with a forward–backward symmetric cylindrical geometry and a near \(4\pi \) coverage in solid angle.Footnote 1 It consists of an inner tracking detector surrounded by a thin superconducting solenoid providing a 2 T axial magnetic field, electromagnetic and hadron calorimeters, and a muon spectrometer.

The inner tracking detector (ID) covers the pseudorapidity range \(|\eta | < 2.5\). It consists of silicon pixel, silicon microstrip, and transition radiation tracking detectors. Lead/liquid-argon (LAr) electromagnetic calorimeters (ECAL) provide electromagnetic (EM) energy measurements with high granularity. A steel/scintillator-tile hadron calorimeter (HCAL) covers the central pseudorapidity range (\(|\eta | < 1.7\)). The endcap and forward regions are instrumented with LAr calorimeters for EM and hadronic energy measurements up to \(|\eta | = 4.9\). The muon spectrometer (MS) surrounds the calorimeters and is based on three large air-core toroidal superconducting magnets with eight coils each. The field integral of the toroids ranges between 2.0 and 6.0 Tm across most of the detector. The muon spectrometer includes a system of precision tracking chambers and fast detectors for triggering. A two-level trigger system [26] selects events to be recorded at a reduced rate. The first level is a hardware implementation aiming to reduce the rate to around 100 kHz, while the software-based high-level trigger provides the remaining rate reduction to approximately 1 kHz.

This search uses the pp collision data at \(\sqrt{s}=13\) \(\text {TeV}\) recorded by the ATLAS detector during the data-taking between 2015 and 2018 with a total integrated luminosity of \(139.0\pm 2.4\) \(\mathrm {fb^{-1}} \)[27].

A combination of multiple single-lepton and missing transverse momentum (\(E_{\text {T}}^{\text {miss}}\)) triggers with varying thresholds, as well as lepton quality and isolation requirements is used [26, 28]. During data-taking, as the instantaneous luminosity increased, the thresholds for unprescaled single-lepton triggers with tight isolations were increased in stages: the electron transverse energy (\(E_{\text {T}} \)) threshold was increased from 24 to 26 GeV, and the muon transverse momentum (\(p_{\text {T}} \)) threshold was increased from 20 to 26 GeV. Similarly, the threshold of the \(E_{\text {T}}^{\text {miss}}\) triggers increased from 70 to 110 GeV. Lepton triggers with tight isolations were complemented by those with looser isolations but higher \(E_{\text {T}} \) or \(p_{\text {T}} \) thresholds. The search uses the \(E_{\text {T}}^{\text {miss}}\) triggers in the 0-lepton channel and single-lepton triggers in the 2-lepton channel. The trigger efficiencies are greater than 90% for signal events targeted by these two channels, independent of the resonance mass. For the 1-lepton channel, the single-electron triggers were used in the electron case, but the single-muon triggers were used only for \(p_{\text {T}} (\mu \nu )<150\) GeV in the muon case. For \(p_{\text {T}} (\mu \nu )>150\) GeV, since the calculation of \(E_{\text {T}}^{\text {miss}}\) at the trigger level does not account for the presence of minimum ionising muons, the \(E_{\text {T}}^{\text {miss}}\) triggers were used instead. Using the unprescaled \(E_{\text {T}}^{\text {miss}}\) triggers miminises the impact of the efficiency loss from the limited geometric coverage of the muon triggers. The trigger efficiency for the 1-lepton channel increases from approximately 80% at 300 GeV to be above 90% at a resonance mass of 500 GeV.

Events were retained for analysis if they were recorded with all detector systems operating normally and pass data-quality requirements [29]. Collision vertices are formed from tracks with \(p_{\text {T}} >500\) MeV. The vertex candidate with the highest \(\sum p_{\text {T}} ^2\) of its associated tracks is selected as the primary vertex. All events are required to contain a primary vertex with at least two associated tracks.

3 Simulation of signal and background processes

Monte Carlo (MC) simulations were used for background modellings and estimations, evaluations of signal efficiencies, optimisations of event selections, and estimations of systematic uncertainties. Generated signal and background events were processed through the full ATLAS detector simulation program [30] based on \(\textsc {Geant} 4\) [31]. Multiple overlaid pp collisions (pile-up) were simulated with the soft QCD processes of Pythia 8.186 [32] using the A3 set of tuned parameters [33] and the NNPDF23lo parton distribution function (PDF) set [34]. All simulated events are processed with the same trigger and reconstruction algorithm as the data. Scale factors were used to correct differences between the data and simulations.

3.1 Signal models and simulation

Three types of resonances corresponding to different spins are considered in the search. The first one is a scalar neutral radion, introduced in some bulk RS models to stabilise the radius of the compactified extra dimension \(r_\mathrm {c}\) [17, 18]. The coupling of the RS radion field to SM fields is inversely proportional to \(\Lambda _R = \mathrm {e}^{-k\pi r_\mathrm {c}} \sqrt{6M^3_5/k}\) [35,36,37], where \(M_5\) is the five-dimensional Planck mass, and k is the curvature factor. The RS radion events were simulated with \(k\pi r_\mathrm {c}=35\) and \(\Lambda _R=3\) \(\text {TeV}\)  [36]. The RS radion couples to SM fermions with a strength proportional to the fermion mass and to SM vector bosons with a strength proportional to the square of the boson mass, similarly to a heavy Higgs boson. However, the RS radion has a much narrower width due to its overall weaker couplings to SM particles. For example, the intrinsic width of a 3 \(\text {TeV}\) RS radion is approximately 3% of its mass, assuming \(\Lambda _R=3\) \(\text {TeV}\). RS radions can be produced through both the ggF and VBF processes at the LHC as shown in Fig. 1.

The second type considered comprises two heavy vector bosons described in the HVT framework [19]: an electrically charged \(W^\prime \) boson and an electrically neutral \(Z^\prime \) boson produced through the DY and VBF processes. The new heavy vectors couple to the Higgs and the SM gauge bosons via a combination of parameters \(g_V c_H\) and to the fermions via the combination \(g^2/g_V~c_F\). The parameter \(g_V\) represents the typical strength of the vector boson interaction, while the parameters \(c_H\) and \(c_F\) are expected to be of the order of unity in most models. Benchmark Model A [19] (\(g_V =1\)) is representative of a model of weakly coupled vector resonances in an extension of the SM gauge group where the HVT bosons have comparable decay branching ratios into SM fermions and vector bosons. Model B [19] with \(g_V = 3\), is representative of a composite model scenario where the HVT boson couplings to fermions are suppressed. In Model C, \(g_V=c_H=1\) and the HVT boson coupling to fermions was set to zero, so that only VBF production is possible. The \(W^\prime \rightarrow WZ\) and \(Z^\prime \rightarrow WW\) decays were considered in this search.

The third benchmark resonance searched for is a spin-2 bulk RS graviton \(G_{\mathrm {KK}}\) which appears as the first KK excitation of the gravitational field in a bulk RS graviton model [7, 20, 21]. The \(G_{\mathrm {KK}}\) couplings to light fermions are suppressed and therefore decays into final states involving heavy quarks, Higgs or vector bosons are favoured. The strength of the coupling depends on \(k/\overline{M}_{\mathrm {Pl}}\), where k corresponds to the curvature of the warped extra dimension and \(\overline{M}_{\mathrm {Pl}}\) is the effective four-dimensional Planck scale. The value of \(k/\overline{M}_{\mathrm {Pl}}\) is typically of \({{\mathcal {O}}}(1)\), and this and the \(G_{\mathrm {KK}}\) mass are the only two free parameters. The \(G_{\mathrm {KK}}\) has a mass-dependent width, which is 3.7% of its mass at 500 GeV and 6.4% at 5 \(\text {TeV}\) for \(k/\overline{M}_{\mathrm {Pl}}=1\). It can be produced through the ggF and VBF processes and decays into WW and ZZ with sizeable branching ratios. The \(G_{\mathrm {KK}}\) samples were generated with \(k/\overline{M}_{\mathrm {Pl}}=1\).

Signal events for the HVT and bulk RS graviton (radion) models were generated with MadGraph5_aMC@NLO v2.2.2 (v2.6.1) [38] at leading order (LO) using the NNPDF23lo PDF set. For the production of resonances in the HVT model, both the DY and VBF mechanisms were simulated, and the RS radion and \(G_{\mathrm {KK}}\) resonances were produced via both the ggF and VBF mechanisms. For all signal models and production mechanisms, the generated events were interfaced to Pythia 8.186 (8.230 for the RS radion model) [39] for parton showering, hadronisation, and the underlying event. This interface relied on the A14 set of tuned parameters [40] for events generated with MadGraph5_aMC@NLO at LO.

As examples, Table 1 shows the theoretical cross-sections, the diboson decay branching ratios, and the total widths of the resonances for two different mass values.

Table 1 List of benchmark signal models. Predictions of cross-section \(\sigma \), branching ratio \(\mathcal {B}\) into WW, WZ, or ZZ, and intrinsic width divided by the resonance mass \(\Gamma /m\), for the given hypothetical new particle at \(m=800\,\hbox {GeV}\) and \(3\,\text {TeV}\) are summarised

3.2 Background process simulation

Background processes include W and Z boson production in association with jets (\(W\!+\!\mathrm {jets}\) and \(Z\!+\!\mathrm {jets}\), collectively denoted by \(V\!+\!\mathrm {jets}\)), top-quark production (both top-quark pair, \(t\bar{t}\), and single-top-quark), non-resonant diboson production (WWWZ and ZZ), and multijet production. MC samples were produced to model these background processes with the exception of multijet production, for which data were used to estimate its contribution.

The production of \(V+\)jets was simulated with the Sherpa v2.2.1 [41] generator using the matrix elements (ME) with next-to-leading order (NLO) accuracy for up to two jets, and with leading-order (LO) accuracy for up to four jets, calculated with the Comix [42] and OpenLoops [43, 44] libraries. They were matched with the Sherpa parton shower [45] using the MEPS@NLO prescription [46,47,48,49] using the set of tuned parameters developed by the Sherpa authors. The NNPDF30nnlo set of PDFs [50] was used and the samples were normalised to a next-to-next-to-leading-order (NNLO) prediction [51] with a flat K-factor. Simulated V+jets events from MadGraph5_aMC@nlo v2.2.2 [38] using LO-accurate ME with up to four final-state partons were used to estimate the possible mismodelling of the Sherpa sample. The ME calculation employed the NNPDF30nlo set of PDFs [50] or the NNPDF23lo set of PDFs. Events were interfaced to Pythia 8.186 for the modelling of the parton shower, hadronisation, and underlying event. The A14 tune [40] of Pythia was used with the NNPDF23lo PDF set. The decays of bottom and charm hadrons were performed by EvtGen v1.2.0.

Samples of \(t\bar{t}\) and single-top-quark events were generated with Powheg-Box [52,53,54,55] v2 at NLO with the NNPDF30nlo PDF set. The parameter \(h_{\mathrm {damp}}\), which regulates the high-\(p_{\text {T}}\) radiation in the Powheg, was set to 1.5 \(m_{t}\) to obtain good data–MC agreement at high \(p_{\text {T}}\) [56], where \(m_{t} = 172.5~\hbox {GeV}\) was the top-quark mass used in the simulation. The parton shower, fragmentation, and underlying event were simulated using Pythia 8.230 [39] with the NNPDF23lo PDF set and the A14 tune. The decays of bottom and charm hadrons were performed by EvtGen v1.6.0.

Diboson processes were simulated with Sherpa v2.2.1 using the ME at NLO accuracy in QCD for up to one additional parton and at LO accuracy for up to three additional parton emissions, including off-shell effects and Higgs boson contributions. The NNPDF30nnlo PDF set was used. The electroweak VVjj samples were generated by MadGraph5_aMC@NLO 2.4.3 [38] and were used, together with the Sherpa diboson sample, for the VBF analysis. The NNPDF30lo PDF set was used. The parton showers and hadronisation were modelled with Pythia 8.186 using the A14 tune.

Theoretical cross-sections were used to normalise background contributions. The cross-sections of single-top-quark t- and s-channel production were calculated with the Hathor v2.1 program [57, 58], while the Wt-channel followed the prescriptions from Refs. [59, 60]. Cross-sections for diboson production were calculated at NLO [55, 61]. The normalisations of \(V\!+\!\mathrm {jets}\) and \(t\bar{t}\) contributions were estimated from data using the control regions as described in Sect. 6.1.

4 Object reconstruction and identification

Leptons, jets, and \(E_{\text {T}}^{\text {miss}}\) are basic building blocks for this search. Their identification requirements are summarised briefly in this section.

4.1 Leptons

Electrons are reconstructed from energy clusters that are consistent with EM showers in the ECAL and are matched to tracks in the ID [62]. They are required to have transverse energy \(E_{\text {T}} >7\) GeV and pseudorapidity \(|\eta |<2.47\), excluding the ECAL barrel-endcap transition region: \(1.37<|\eta |<1.52\). To reduce backgrounds from misidentification and non-prompt sources, electrons must meet a likelihood-based criterion [62]. The likelihood is used to classify electrons as having either Loose, Medium, or Tight quality.

Muons are identified by matching MS tracks with those in the ID and are required to have transverse momentum \(p_{\text {T}} >7\) GeV and pseudorapidity \(|\eta |<2.5\). An identification requirement based on information from the ID and MS systems is applied to reduce backgrounds from misreconstruction and muons originating from hadron decays in flight. Similarly to electrons, muons are classified as having either Loose, Medium, or Tight quality [63].

Leptons are required to have associated tracks satisfying \(|d_0/\sigma _{d_0}|<5~(3)\) and \(|z_0\times \sin \theta |<0.5\) mm for electrons (muons), where \(d_0\) is the transverse impact parameter relative to the beam line, \(\sigma _{d_0}\) is its uncertainty, and \(z_0\) is the distance between the longitudinal position of the track along the beam line at the point where \(d_0\) is measured and the longitudinal position of the primary vertex.

Leptons from W and Z boson decays are required to have \(p_{\text {T}} >30\) GeV. They are expected to be isolated from other energy deposits in the detector. Thus, isolation criteria based on the sum of track \(p_{\text {T}} \), the sum of calorimeter \(E_{\text {T}} \), or both, in small cones around the lepton direction are used to further reduce backgrounds from non-isolated sources. Leptons of Loose quality with \(p_{\text {T}} <100\) GeV are required to pass a FixedCutLoose isolation requirement and no isolation requirement is applied for \(p_{\text {T}} >100\) GeV so that the leptons from high-\(p_{\text {T}}\) \(Z \rightarrow \ell \ell \) decays are not removed in the presence of nearby leptons. Details can be found in Refs. [62, 63].

4.2 Jets

Small-R jets are reconstructed from calorimeter energy clusters using the anti-\(k_t\) algorithm [64, 65] with a radius parameter of \(R=0.4\). Energy- and \(\eta \)-dependent correction factors derived from MC simulations are applied in order to correct jets back to the particle level [66]. Jets are required to have \(p_{\text {T}} >30\) GeV and \(|\eta |<4.5\). To suppress jets from pile-up interactions, a jet vertex tagger [67] is applied to jets with \(p_{\text {T}} <120\) GeV and \(|\eta |<2.5\), based on information about tracks associated with the primary vertex and pile-up vertices. For forward jets, the uncertainty on pileup modelling is taken into account.

A multivariate algorithm for the identification of small-R jets containing b-hadrons (b-tagging) [68] is used. The algorithm is based on information such as track impact-parameter significances and positions of reconstructed secondary decay vertices. The identified jets, called b-jets, are restricted up to \(|\eta |<2.5\) due to the ID coverage. The b-tagging algorithm has an efficiency of 85% for b-hadrons in simulated \(t\bar{t} \) events, a light-flavour jet rejection factor of 33 and a c-jet rejection of about 3 [68].

Large-R jets are reconstructed from track-calo clusters [69] with the anti-\(k_t\) algorithm, but with the radius parameter increased to \(R=1.0\). The track-calo clusters are formed by combining information from the calorimeter and the ID, utilising the excellent angular resolution of the ID and the improving energy resolution of the calorimeter at high energies. A trimming algorithm [70] is applied to reduce the impact of pile-up and soft radiation overlapping with the jet. The constituents of each jet are reclustered with the \(k_t\) algorithm [71] into smaller \(R = 0.2\) subjets and those subjets are removed if \(p_{\mathrm {T}}^{\mathrm {subjet}} / p_{\text {T}}^{J}< 0.05\), where \(p_{\mathrm {T}}^{\mathrm {subjet}}\) and \(p_{\text {T}}^{J}\) are the transverse momenta of the subjet and the large-R jet, respectively. The large-R jets are required to have \(p_{\text {T}} >200\) GeV, \(|\eta |<2.0\), and a jet mass (\(m_J\)) greater than 50 GeV.

Variable-radius (VR) jets are used to identify b-jets from boosted hadronic \(V \rightarrow qq\) decays that are reconstructed as large-R jets. They are reconstructed from ID tracks associated with large-R jets by using the anti-\(k_t\) algorithm with a \(p_{\text {T}} \)-dependent radius R parameter between 0.02 and 0.4 and a \(\rho \)-parameter of 30 GeV  [72]. They are required to have \(p_{\text {T}} >10\) GeV and \(|\eta |<2.5\). The same b-tagging algorithm which is used for small-R jets is applied to identify variable-radius jets from b-hadrons.

4.3 Overlap removal

An overlap-removal procedure is applied to the selected leptons and jets. If two electrons share the same track, or the separation between their two energy clusters satisfies \(|\Delta \eta |<0.075\) and \(|\Delta \phi |<0.125\), then the lower-\(p_{\text {T}}\) electron is discarded. Electrons that fall within \(\Delta R=0.02\) of a selected muon are also discarded. For nearby electrons and small-R jets, the jet is removed if the separation between the electron and jet satisfies \(\Delta R<\) 0.2; the electron is removed if the separation satisfies 0.2 \(<\Delta R<\) 0.4. For nearby muons and small-R jets, the jet is removed if the separation between the muon and jet satisfies \(\Delta R<\) 0.2 and if the jet has less than three tracks or the energy and momentum differences between the muon and the jet are small; otherwise the muon is removed if the separation satisfies \(\Delta R<\) 0.4. To prevent double-counting of energy from an electron inside the large-R jet, the large-R jet is removed if the separation between the electron and the large-R jet satisfies \(\Delta R < 1.0\).

4.4 Missing transverse quantities

The missing transverse momentum (\(\mathbf {E}_{\mathrm{T}}^{{\mathrm{miss}}}\)) is calculated as the negative vectorial sum of the transverse momenta of calibrated electrons, muons, small-R jets, and unassociated tracks. Large-R jets are not included in the \(\mathbf {E}_{\mathrm{T}}^{{\mathrm{miss}}}\) calculation to avoid double-counting of energy between the small-R jets and large-R jets. Energy depositions due to the underlying event and other soft radiation are taken into account by constructing a ‘soft term’ from ID tracks associated with the primary vertex but not with any reconstructed object [73]. Similarly, the track-based missing transverse momentum, \(\mathbf {p}_{{\mathrm{T}}}^{{\mathrm{miss}}}\), is the negative vectorial sum of the transverse momenta of all good-quality inner-detector tracks that are associated with the primary vertex.

5 Event classification and selections

The search begins with the selection of the leptonically decaying boson \({V_\ell }\). Candidate events are first selected according to the number of Loose leptons and assigned to 0-lepton (\({V_\ell }=Z,\,Z\rightarrow \nu \nu \)), 1-lepton (\({V_\ell }=W,\,W\rightarrow \ell \nu \)) and 2-lepton (\({V_\ell }=Z,\,Z\rightarrow \ell \ell \)) channels. Other lepton multiplicities are excluded from the analysis. Although specific selections differ, the three channels follow the same analysis flow as illustrated in Fig. 2. For each channel, events are further classified into two exclusive VBF and ggF/DY categories as described in Sect. 5.1, targeting their corresponding production processes for heavy resonances.

The selection proceeds to identify the hadronically decaying boson \({V_h}\). Depending on the \(V_h\)-boson momentum, the energy deposits of the two jets from the hadronically decaying V bosons can be well separated or can largely overlap in the detector. Thus the \(V \rightarrow qq\) decay, including \(Z\rightarrow bb\), can be either reconstructed from two resolved small-R jets (\(V \rightarrow jj\)) for low-energy bosons or identified as one merged large-R jet (\(V \rightarrow J\)) for energetic bosons. The \({V_h}\) candidates are identified first through the merged \(V \rightarrow J\) identification and then, if it fails, through the resolved \(V \rightarrow jj\) reconstruction.

Selections specific to each channel are presented in Sect. 5.3. Multiple signal regions (SRs) are defined in order to enhance search sensitivities, as described in Sect. 5.4. Section 5.5 discusses the mass variables used as the final discriminants. The analysis flow is run twice, once for \({V_h}=W\) and once for \({V_h}=Z\), which involves selecting different ranges of \(m_{jj}\) or \(m_J\).

Fig. 2
figure 2

Illustration of the selection flow and signal regions of the \(X\rightarrow VV\rightarrow {V_\ell }{V_h}\) search. The VBF category targets VBF production whereas the ggF/DY category is for the rest. Three signal regions (high purity, low purity and resolved) are selected for each category, based on the \(V \rightarrow qq\) reconstruction. For the 0-lepton channel, no resolved selection is considered. For final states with hadronically decaying Z bosons, the three signal regions in the ggF/DY category are each further split into tagged and untagged according to the b-tagging information about jets from \(Z \rightarrow qq\) decays

5.1 Categorisation of production processes

For the three production processes shown in Fig. 1, the ggF and DY processes have the same final states while the VBF process possesses two additional jets, called VBF-tag jets. The kinematics of these jets differ from those from the V-boson hadronic decays. They are typically well separated in pseudorapidity and usually have large dijet invariant mass. These characteristics were used in previous searches [22, 23] to separate VBF production from ggF/DY production. In this search, a recurrent neural network (RNN) [74, 75] is used to classify the VBF and ggF/DY event topologies. It is built with the Keras [76] library using the Theano python library [77] as a back end for mathematical computations. The RNN has 2 hidden layers with 25 recurrent cells to exploit the hidden correlation of the input sequence.

The RNN uses the four-momenta of small-R jets as input. It is well suited for a variable-length input sequence such as the jet information. The RNN allows to recover events with only one VBF-tag jet reconstructed (\(\sim 30\%\) of signal events), and those events were not selected in previous searches where two VBF-tag jets were required [22, 23]. Although the RNN permits to deal with a large number of input jets, a maximum of two input jets is chosen to minimise the impact of systematic uncertainties associated with additional jets. Only a small increase (2–3%) in the tagging efficiency of VBF events is observed if more than two jets are used as inputs.

For events with large-R jets, small-R jets with angular separations of \(\Delta R<1\) from the leading large-R jets are removed. If there is no large-R jet in an event, the pair of small-R jets with dijet invariant mass closest to the V-boson mass is removed. Up to two remaining small-R jets with the highest \(p_{\text {T}} \) are chosen as the input to the RNN. Events with no small-R jets left are automatically classified as ggF/DY events.

The RNN score distributions depend on the assumed model of a heavy resonance, its mass and decay mode. The RNN trained with the 1 \(\text {TeV}\) scalar resonance in the \(X\rightarrow ZZ\rightarrow \ell \ell qq\) decay is applied for the three leptonic channels, the three resonance models, and all resonance masses.

Figure 3a compares the RNN score of simulated events from VBF and ggF/DY production of a 1 \(\text {TeV}\) resonance in the signal models considered in this search. The RNN discrimination power increases with the resonance mass. An event is classified as a VBF event if its RNN score is above 0.8 and otherwise as a ggF/DY event. The threshold is chosen to maximise the sensitivity to VBF signals. Figure 3b shows the fractions of simulated signal events passing the RNN requirement as functions of the resonance mass for different signal models. The RNN correctly classifies VBF events more than 40% of the time for a diboson resonance heavier than 1 \(\text {TeV}\) with a ggF/DY contamination of about 2–5%. It yields a relative increase in the VBF event selection efficiency from 10% (5%) at 0.5 TeV to 60% (50%) at 3 TeV for a scalar resonance (spin-1 or spin-2 resonance) compared to the previous cut-based selection [22, 23], with similar background rejections.

Fig. 3
figure 3

a RNN score distributions for the production of a 1 \(\text {TeV}\) resonance in the signal models considered for this search; b the fractions of signal events passing the VBF requirement on the RNN score as functions of the resonance mass for both VBF and ggF production

5.2 Reconstruction and identification of the \(V \rightarrow qq\) decay

The \(V \rightarrow J\) candidates are identified from the highest-\(p_{\text {T}}\) large-R jet in an event by requiring its mass \(m_J\) to be in a \(p_{\text {T}}\)-dependent window centred around the expected value of the V-boson mass from simulations [78, 79], as shown in Fig. 4a. The mass window depends on the jet mass resolution [69] and is approximately 30 GeV wide at \(p_{\text {T}} =500\) GeV and increases to about 60 GeV at \(p_{\text {T}} =2.5\) \(\text {TeV}\). A jet substructure variable \(D_{2}^{(\beta =1)}\) is used to assess the quality of the \(V \rightarrow J\) candidates. The variable \(D_{2}^{(\beta =1)}\) is defined as the ratio of three-point to two-point energy correlation functions [80, 81] based on the energies and pairwise angular distances of particles within a large-R jet. The variable is optimised to distinguish between jets originating from a single parton and those originating from \(V \rightarrow qq\) decays. A \(p_{\text {T}}\)-dependent upper (lower) requirement on \(D_{2}^{(\beta =1)}\), shown in Fig. 4b, is employed to select high (low)-purity signal regions as described in Sect. 5.4. Efficiencies for the \(m_J\) requirement alone and for the combined \(m_J\) and \(D_{2}^{(\beta =1)}\) requirements as functions of the large-R jet \(p_{\text {T}} \) are shown in Fig. 5. The efficiency for tagging \(V \rightarrow qq\) decay varies from approximately 40% at low \(p_{\text {T}} \) to 70% at high \(p_{\text {T}} \). The background rejection factor of the W (Z) tagger is estimated using the simulated \( W \rightarrow \ell \nu \left( Z \rightarrow \ell \ell \right) \)+jets events, and is approximately 5 (6) at \(p_{\text {T}} =200~\)GeV and 35 (30) at \(p_{\text {T}} >700\) GeV.

Fig. 4
figure 4

a The upper and lower bounds of \(m_J\) and b the upper (lower) requirements on \(D_{2}^{(\beta =1)}\) selecting the high (low)-purity signal regions as functions of the large-R jet \(p_{\text {T}} \) for the \(V \rightarrow J\) tagging for the W boson and the Z boson

Fig. 5
figure 5

Efficiencies of the \(m_J\) and \(D_{2}^{(\beta =1)}\) requirements as functions of the large-R jet \(p_{\text {T}} \) for the \(V \rightarrow J\) tagging for a the W boson and b the Z boson

The \(V \rightarrow jj\) candidates are reconstructed from two small-R jets within \(|\eta |<2.5\). The leading jet is required to have \(p_{\text {T}} >60\) GeV and the subleading jet is required to have \(p_{\text {T}} >30\ (45)\) GeV in the 2-lepton (1-lepton ) channel. No resolved \(V \rightarrow jj\) reconstruction is considered for the 0-lepton channel due to the large multijet background. The two highest-\(p_{\text {T}}\) small-R jets in \(|\eta |<2.5\) are chosen to form the \(V \rightarrow jj\) candidate except for the \(Z\rightarrow bb\) reconstruction, for which events are required to have exactly two b-tagged jets, and in which case they are used. The invariant mass of the two jets, \(m_{jj}\), must be consistent with that of the V boson by satisfying \(62<m_{jj}<97\) GeV for \(W \rightarrow jj\) and \(70<m_{jj}<105\) GeV for \(Z \rightarrow jj\). Fixed mass windows are applied because the dijet mass resolution is largely independent of the dijet \(p_{\text {T}}\) for the resonance masses to which the resolved analysis is sensitive.

5.3 Event selections for individual leptonic channels

Event selections for all three leptonic channels consist of the selections for the leptonically and hadronically decaying V bosons and an event-level selection designed to reduce backgrounds specific to each channel. The selection of hadronically decaying V bosons is common to all three channels. It requires a \(V \rightarrow qq\) candidate identified by either the merged or resolved technique as discussed above. The other selections are specific to individual leptonic channels and are described below. An overview of the selections is shown in Table 2.

Table 2 Overview of the main \(X\rightarrow VV\rightarrow {V_\ell }{V_h}\) selection criteria; the text gives more details. \(\mathcal{R}_{p_{\text {T}}/m}\) stands for \(\min (p_{\text {T}} ^{{V_\ell }},p_{\text {T}} ^{{V_h}})/m_{VV}\)

For the merged selection, since the leading large-R jet is considered as the \(V \rightarrow J\) candidate, any small-R jet within an angular radius \(R=1\) around it is removed. For the resolved selection, large-R jets are ignored and no small-R jets are removed.

An event veto based on b-tagging information is also applied. For signal events, b-jets can arise from the \(Z\rightarrow bb\) decays. For \({V_h}= Z\), classification based on number of b-tagged jets in \(Z \rightarrow qq\) candidates is applied, as described in Sect. 5.4. For \({V_h}= W\), events are required to have at most one b-tagged jets in \(W\rightarrow qq\) candidates, considering mis-identification rate for charm hadrons from \(W \rightarrow sc\) decay. In the merged selection, events are vetoed if there are more than two b-tagged variable-radius jets associated with the leading large-R jets.

Unless specifically noted, the same selections are applied for the VBF and ggF/DY categories. The merged selection is applied first and the resolved selection is applied only to events failing the merged selection.

5.3.1 0-lepton : \(ZV\rightarrow \nu \nu qq\)

The 0-lepton channel targets the \(ZV\rightarrow \nu \nu qq\) final state from \(R\rightarrow ZZ\), \(W'\rightarrow ZW\) and \(G_{\mathrm {KK}}\rightarrow ZZ\) decays. Events in this final state have a large \(E_{\text {T}}^{\text {miss}}\) and a \(V \rightarrow qq\) candidate. Due to high \(E_{\text {T}}^{\text {miss}}\) trigger thresholds and the expected large multijet background from mismeasurement at low \(E_{\text {T}}^{\text {miss}}\), events are required to have \(E_{\text {T}}^{\text {miss}} >250\) GeV and no Loose leptons, to suppress background from multijet events and single W bosons respectively.

Additional requirements are applied to further reduce the multijet background. These include \(p_{{\text {T}}}^{\text {miss}}>50\) GeV, and an azimuthal opening angle between \(\mathbf {E}_{\mathrm{T}}^{{\mathrm{miss}}}\) and \(\mathbf {p}_{{\mathrm{T}}}^{{\mathrm{miss}}}\) satisfying \(\Delta \phi (\mathbf {E}_{\mathrm{T}}^{{\mathrm{miss}}},\mathbf {p}_{{\mathrm{T}}}^{{\mathrm{miss}}})<1\). Furthermore, the azimuthal angle between \(\mathbf {E}_{\mathrm{T}}^{{\mathrm{miss}}}\) and the nearest small-R jet must satisfy \(\min \Delta \phi (\mathbf {E}_{\mathrm{T}}^{{\mathrm{miss}}}, j)>0.4\). With these angular requirements along with the \(E_{\text {T}}^{\text {miss}}\) and \(p_{{\text {T}}}^{\text {miss}}\) requirements, the multijet background becomes negligible.

The high \(E_{\text {T}}^{\text {miss}}\) requirement is efficient only for signal events with very heavy resonances. Therefore, only the merged selection is considered for this channel.

5.3.2 1-lepton : \(WV\rightarrow \ell \nu qq\)

The 1-lepton channel is designed for the \(WV\rightarrow \ell \nu qq\) final state from \(R\rightarrow WW\), \(W'\rightarrow WZ\), \(Z'\rightarrow WW\), and \(G_{\mathrm {KK}}\rightarrow WW\). Events in this channel must have exactly one Tight electron or Medium muon, with \(p_{\text {T}} > 30\) GeV, and no other leptons satisfying the Loose quality; \(E_{\text {T}}^{\text {miss}}\) greater than 60 GeV; and a transverse momentum of the lepton-\(\mathbf {E}_{\mathrm{T}}^{{\mathrm{miss}}}\) system (i.e. the reconstructed \({V_\ell }\)), \(p_{\text {T}}^{V_\ell }\), greater than 75 GeV.

For the merged selection, the \(E_{\text {T}}^{\text {miss}}\) and \(p_{\text {T}}^{V_\ell }\) thresholds are raised to 100 GeV and 200 GeV, respectively. Considering that boson \(p_{\text {T}}\) is expected to be approximately \(0.5 m_{VV}\), events are further required to have \(\mathcal{R}_{p_{\text {T}}/m}\), defined as \(\min (p_{\text {T}}^{V_\ell },p_{\text {T}}^{V_h})/m_{VV}\), greater than \(0.35\,(0.25)\) for the ggF/DY (VBF) category. Here \(p_{\text {T}}^{V_h}\) is the \(p_{\text {T}}\) of the \({V_h}\) candidate, i.e. the leading large-R jet, and \(m_{VV}\) is the invariant mass of the VV system reconstructed from the \(\ell \nu \) and the leading large-R jet.Footnote 2 This requirement suppresses background significantly at large \(m_{VV}\) while maintaining high efficiencies for signal events.

For the resolved selection, a set of azimuthal angular requirements are designed and applied to reduce large multijet backgrounds expected at low \(E_{\text {T}}^{\text {miss}}\) and \(p_{\text {T}}^{V_\ell }\). They are \(\Delta \phi (\ell ,\mathbf {E}_{\mathrm{T}}^{{\mathrm{miss}}})<1.5\), \(\Delta \phi (j_1, j_2)<1.5\), \(\Delta \phi (\ell ,j_{1/2})>1.0\) and \(\Delta \phi (\mathbf {E}_{\mathrm{T}}^{{\mathrm{miss}}},j_{1/2})>1.0\). Here \(j_{1/2}\) refers to both \(j_1\) and \(j_2\), which form the \(V \rightarrow jj\) candidate. Similarly to the merged selection, a kinematic criterion of \(\mathcal{R}_{p_{\text {T}}/m}>0.35\,(0.25)\) is imposed for the ggF/DY (VBF) category, where \(p_{\text {T}}^{V_h}\) is the \(p_{\text {T}}\) of the \(V \rightarrow jj\) candidate and \(m_{VV}\) is reconstructed from the \(\ell \nu \) and dijet system.

Additional b-jets can originate from background \(t\bar{t}\) and single-top-quark events. To reduce this background, events are vetoed if there are one or more small-R b-jets beyond those selected as the \(V\rightarrow qq\) candidate.

5.3.3 2-lepton : \(ZV\rightarrow \ell \ell qq\)

The 2-lepton channel is intended for the \(ZV\rightarrow \ell \ell qq\) final state from \(R\rightarrow ZZ\), \(W'\rightarrow ZW\) and \(G_{\mathrm {KK}}\rightarrow ZZ\). The event selection begins with the identification of the \(Z\rightarrow \ell \ell \) decay. The \(Z\rightarrow \ell \ell \) candidates are formed from two same-flavour leptons with \(p_{\text {T}} > 30\) GeV and satisfying the Loose criteria defined in Sect. 4. Muon pairs are required to have opposite charges. Because electrons are more susceptible to charge misidentification, no charge requirement is applied. The dilepton invariant mass, \(m_{\ell \ell }\), must be consistent with the Z boson mass. A fixed \(m_{\ell \ell }\) window of [83, 99] GeV is applied to electron pairs, while the dilepton \(p_{\text {T}}\) is used to define a \(p_{\text {T}}^{\ell \ell }\)-dependent window of \([85.6\;\hbox {GeV}- 0.0117\times p_{\text {T}}^{\ell \ell },\; 94.0\;\hbox {GeV}+ 0.0185\times p_{\text {T}}^{\ell \ell }]\) that is required for the muon pairs because of the deteriorating muon momentum resolution at high \(p_{\text {T}}\). The mass windows are chosen to maintain approximately 95% selection efficiency for \(Z\rightarrow \ell \ell \). Events with additional Loose leptons are removed.

The selected \(Z\rightarrow \ell \ell \) events are required to have \(\mathcal{R}_{p_{\text {T}}/m}>0.35\;(0.25)\) for the ggF/DY (VBF) category for the merged selection and \(\mathcal{R}_{p_{\text {T}}/m}>0.35\) for both the ggF/DY and VBF categories for the resolved selection. Again, \(p_{\text {T}}^{V_\ell }\) is the \(p_{\text {T}}\) of the leptonic V candidate (\(p_{\text {T}}^{V_\ell }=p_{\text {T}}^{\ell \ell }\) in this case), \(p_{\text {T}}^{V_h}\) is the \(p_{\text {T}}\) of the \({V_h}\) candidate, and \(m_{VV}\) is the invariant mass of the \({V_\ell }\) and \({V_h}\) system. This requirement exploits the kinematic feature of highly boosted boson decays expected from heavy resonances to reduce backgrounds.

5.4 Signal region definitions

Multiple signal regions are defined using the properties of the hadronically decaying V boson as illustrated in Fig. 2. Events passing the merged selection are assigned to either high-purity (HP) or low-purity (LP) signal regions according to the quality of their \(V \rightarrow J\) candidates. Those with \(V \rightarrow J\) candidates passing the \(D_{2}^{(\beta =1)}\) requirement of the boson tagger [78] are selected for the HP SR, otherwise for the LP SR. The combined \(m_J\) and \(D_{2}^{(\beta =1)}\) efficiencies for the HR SR are shown in Fig. 5 as functions of V \(p_{\text {T}} \). Events passing the resolved selection form the resolved SR.

For \(V_h=Z\), about 21% of \(Z \rightarrow qq\) decays are \(Z\rightarrow bb\), whereas jets from the dominant background source, \(V\!+\!\mathrm {jets}\), have a smaller heavy-quark content. To exploit this difference, the HP, LP, and resolved SRs are each further split into tagged and untagged SRs in the ggF/DY category if the hadronically decaying boson is a Z boson, i.e. \(V=Z\). The b-tagging is not applied in the VBF category due to the limited number of events. Classification based on the b-tagging is not applied for \(V_h~=~W\). For the merged selection, the splitting is made by applying b-tagging to variable-radius jets associated to the leading large-R jet. Events are tagged only if the two leading variable-radius jets are both b-tagged. For the resolved selection, events are tagged if the \(Z\rightarrow jj\) is formed from two b-jets and untagged otherwise.

Classifications in terms of ggF/DY and VBF categories include: merged and resolved reconstruction of the \(V_h\rightarrow qq\) decay, high and low purity for the merged reconstruction, tagged and untagged identification of the \(Z\rightarrow qq\) decay, and different mass windows for the \(W \rightarrow qq\) and \(Z \rightarrow qq\) decays. This results in 10 SRs for the 0-lepton channel and 15 SRs each for the 1-lepton and 2-lepton channels, for a total of 40 SRs. Because of overlapping mass windows used to select the hadronic decays of the W and Z bosons, these SRs are not orthogonal. In terms of the diboson final states of resonance decays, there are 6 SRs for \(X\rightarrow WW\), 15 for \(X\rightarrow ZZ\), and 19 for \(X\rightarrow WZ\). SRs in each diboson final state are orthogonal. The \(X\rightarrow WW\) and \(X\rightarrow ZZ\) SRs are orthogonal by design, while the \(X\rightarrow WZ\) SRs can overlap with either \(X\rightarrow WW\) or \(X\rightarrow ZZ\) SRs.

5.5 Reconstruction of invariant and transverse resonance mass

Either the invariant mass \(m_{VV}\) or the transverse mass \(m_{\mathrm {T}}\) of the selected VV final states is used as the final discriminant to extract the signal. Heavy resonances would manifest themselves as resonant structures above the SM background in the invariant mass distributions or as broad enhancements in the transverse mass distributions.

For the 0-lepton channel of \(X\rightarrow ZV\rightarrow \nu \nu qq\), no resonance mass reconstruction is possible because of the two undetected neutrinos. Instead a transverse mass defined as:

$$\begin{aligned} m_{\mathrm {T}}= \sqrt{(p_{\text {T}}^{J}+E_{\text {T}}^{\text {miss}})^2-(\mathbf {p}_\text {T}^J+\mathbf {E}_{\mathrm{T}}^{{\mathrm{miss}}})^2} \end{aligned}$$

is used as the discriminant for the merged selection. For the 1-lepton and 2-lepton channels, the VV mass is calculated for both the merged and the resolved reconstruction of the \(V \rightarrow qq\) decay.

Muon momentum resolution deteriorates at high \(p_{\text {T}}\), significantly impacting the \(Z\rightarrow \mu \mu \) mass resolution and consequently the resonance mass resolution in the 2-lepton channel. This deterioration is particularly severe for very heavy resonances, especially in the merged selection. To mitigate the impact, a scale of \(m_Z/m_{\mu \mu }\) is applied to the four-momentum of the dimuon system, effectively fixing the dimuon mass to the Z boson mass [82]. This scaling improves the \(m_{\mu \mu J}\) resolution by about 7% in the merged analysis for a scalar resonance. The scaling is not applied for \(W\rightarrow \mu \nu \) because of the undetected neutrino.

For the resolved selection, the \(m_{VV}\) resolution is improved by \(2\%\) through the scaling of the dijet four-momentum by a factor of \(m_V/m_{jj}\), with \(m_V\) being either the Z or W boson mass [82]. No \(m_V/m_J\) scaling is applied to the merged selection as the improvement in the \(m_{VV}\) resolution is found to be negligible.

5.6 Signal efficiencies and mass resolutions

Signal selection efficiencies depend on the signal model, the production process, and the mass of heavy resonances. Figures 6, 7 and 8 show the acceptance times efficiency (\(A\times \epsilon \)) of the signal events from MC simulations as a function of the resonance mass for (a) ggF/DY and (b) VBF production, combining all SRs of both the ggF/DY and VBF categories of both the resolved and merged analyses. The \(A\times \epsilon \) curves are largely determined by the merged analyses. The resolved analyses contribute only in the low mass region, up to approximately 1 \(\text {TeV}\).

Large differences in \(A\times \epsilon \) shown in the figures for different resonances are due to the different spins of these resonances. The spin-0 RS radions are produced with isotropic angular distributions for both ggF and VBF production. On the other hand, the spin-1 HVT resonances and spin-2 RS gravitons are produced more centrally (more forward) for ggF/DY (VBF) production. These different angular distributions lead to very different efficiencies of the \(\mathcal{R}_{p_{\text {T}}/m}\) requirement. Moreover, the angular requirements between jets and \(E_{\text {T}}^{\text {miss}}\) in the 0-lepton channel are more efficient for DY production of HVT resonances than for ggF production of RS radions and RS gravitons due to the different colour factors for initial-state quarks and gluons.

Signal contributions from \(W\rightarrow \tau \nu \rightarrow \ell \nu \nu \nu \) decays are included in the 1-lepton channel, but not in the 2-lepton channels. Approximately 10–12% of the signal events in SRs are from \(W\rightarrow \tau \nu \rightarrow \ell \nu \nu \nu \) decays in the 1-lepton channel. These events have mass distributions similar to those from \(W\rightarrow \ell \nu \) decays. In the 2-lepton channel, signal contributions from \(Z\rightarrow \tau \tau \rightarrow 2\ell +4\nu \) decays are suppressed by the small \(\tau \tau \rightarrow 2\ell +4\nu \) branching ratio and the Z boson mass requirement. They are found to be negligible. The 0-lepton channel targeting the \(X\rightarrow ZV\rightarrow \nu \nu qq\) signal should also be sensitive to the \(X\rightarrow WV\rightarrow \ell \nu qq,\,\tau \nu qq\) signal due to either the inefficiency of the lepton veto or the lack of a \(\tau \)-veto. This additional ‘cross-channel’ signal contribution is neglected in this search.

Fig. 6
figure 6

Selection acceptance times efficiency for the \(X\rightarrow ZV\rightarrow \nu \nu qq\) signal events from MC simulations as a function of the resonance mass for a ggF/DY and b VBF production, combining HP and LP signal regions. The light shaded band represents the total statistical and systematic uncertainties for the RS radion model, and the total uncertainties are similar for the other signal models

Fig. 7
figure 7

Selection acceptance times efficiency for the \(X\rightarrow WV\rightarrow (e\nu /\mu \nu /\tau \nu )~qq\) signal events from MC simulations as a function of the resonance mass for a ggF/DY and b VBF production, combining all SRs of both the resolved and merged analyses. Signal contributions from \(W\rightarrow \tau \nu \) decays are included in the acceptance calculation. The light shaded band represents the total statistical and systematic uncertainties for the RS radion model, and the total uncertainties are similar for the other signal models. The ‘bump’ structure around 800 GeV is due to the decreasing contribution from the resolved analysis at higher masses

Fig. 8
figure 8

Selection acceptance times efficiency for the \(X\rightarrow ZV\rightarrow \ell \ell qq\) signal events from MC simulations as a function of the resonance mass for a ggF/DY and b VBF production, combining all SRs of both the resolved and merged analyses. The light shaded band represents the total statistical and systematic uncertainties for the RS radion model, and the total uncertainties are similar for the other signal models. The decreases in efficiencies for resonance masses above approximately 2.5 \(\text {TeV}\) are due to the merging of electrons from the highly boosted \(Z\rightarrow ee\) decays. The ‘bump’ structure around 800 GeV is due to the decreasing contribution from the resolved analysis at higher masses

The resonance decays are fully reconstructed in the 1-lepton and 2-lepton channels. In the 1-lepton channel, the \(m_{\ell \nu jj}\) distributions from the resolved \(V \rightarrow jj\) reconstruction have widths of approximately 8% of the resonance mass for narrow resonances, whose intrinsic widths are smaller than the detector mass resolution, with masses of 0.5–1 \(\text {TeV}\).Footnote 3 The width of the \(m_{\ell \nu J}\) distribution from the merged \(V \rightarrow J\) reconstruction varies from 7% at 1 \(\text {TeV}\) to 4% at 5 \(\text {TeV}\). Similarly, in the 2-lepton channel the \(m_{\ell \ell jj}\) resolution is \(\sim 6\%\) for resonance masses of 0.5–1 \(\text {TeV}\) and the \(m_{\ell \ell J}\) resolution varies from approximately 4% at 1 \(\text {TeV}\) to 2% at 5 \(\text {TeV}\).

6 Background estimations

Relevant background sources for the search are \(V\!+\!\mathrm {jets}\), \(t\bar{t}\) and single top, non-resonant diboson, and multijet production. Their relative importance depends on the final state. The largest contributions are from \(Z\!+\!\mathrm {jets}\) and \(W\!+\!\mathrm {jets}\) in the 0-lepton channel and \(W\!+\!\mathrm {jets}\) and \(t\bar{t}\) in the 1-lepton channel. In the 2-lepton channel, the \(Z\!+\!\mathrm {jets}\) background dominates. The multijet background is negligible except in the resolved SRs of the 1-lepton channel. In the tagged SRs, the \(t\bar{t}\) and single-top-quark contributions are enhanced and are in fact dominant in the 1-lepton channel.

MC simulations are used to simulate kinematics of background events except for multijet events. Contributions from diboson and single-top-quark processes are normalised to their theoretical cross-sections, whereas the \(V\!+\!\mathrm {jets}\) and \(t\bar{t}\) contributions are normalised using data through control regions. The multijet background is estimated from data. The definitions of control regions and the method of multijet estimation are described below.

6.1 Control regions for \(W\!+\!\mathrm {jets}\), \(Z\!+\!\mathrm {jets}\), and \(t\bar{t}\)

Control regions (CR) are designed to constrain the normalisations of the \(W\!+\!\mathrm {jets}\), \(Z\!+\!\mathrm {jets}\) and \(t\bar{t} \) background contributions using data, eliminating the reliance on the theoretical cross-sections, which are often less reliable in the phase-space regions covered by this search. Events in CRs are selected from those failing the selections of the SRs, but are otherwise expected to have event topologies similar to those in SRs and small contaminations from potential signals.

CRs for the \(W\!+\!\mathrm {jets}\) background are defined using events in the 1-lepton channel by reversing the \(m_J\) or \(m_{jj}\) requirements of the SR selections, but events are otherwise selected in the same way as those in the corresponding SRs. For the merged selection, \(m_J\) must fall outside of the mass windows of both W and Z boson tagging. For the resolved selection, \(m_{jj}\) is required to be in the range from 50 to 150 GeV, excluding the combined W and Z mass window of 62–105 GeV. The \(W\!+\!\mathrm {jets}\) events are expected to be the dominant contribution in these CRs, except in the b-tagged CRs, where the \(t\bar{t}\) contribution dominates. CRs for the \(Z\!+\!\mathrm {jets}\) background are defined the same way, but using 2-lepton events. The \(Z\!+\!\mathrm {jets}\) events dominate in all \(Z\!+\!\mathrm {jets}\) CRs, even in the tagged CRs.

CRs for the \(t\bar{t}\) background are defined using 1-lepton events, selected in the same way as the 1-lepton SRs except for the requirement of an additional small-R b-jet unassociated with the \(V\rightarrow jj/J\) candidate instead of the b-jet veto. Moreover, \(m_{jj}\) is required to be between 50 and 150 GeV in the case of the resolved selection.

6.2 Multijet background

In the resolved SRs of the 1-lepton channel, multijet production is a non-negligible background source. Multijet events can mimic signal events if there is a lepton from either jet misidentification or heavy-quark decay, along with a large \(E_{\text {T}}^{\text {miss}}\) from energy mismeasurements. The multijet contributions are difficult to model through MC simulations and are therefore estimated from data. A template method is used to estimate the multijet contributions. The method derives the shapes of the \(E_{\text {T}}^{\text {miss}} \) distributions of the multijet contributions from multijet-enriched control regions (MJCR), one for each signal and control region. MJCRs are designed to be orthogonal to both the SRs and CRs as defined above. For the muon channel, MJCRs are defined only for the single-muon trigger, i.e. events with \(p_{\text {T}} (\mu \nu )<150\) GeV, since the multijet contributions to the \(E_{\text {T}}^{\text {miss}}\)-triggered events with \(p_{\text {T}} (\mu \nu )>150\) GeV are found to be negligible.

Events in MJCRs are selected by modifying the lepton requirements used for the SR and CR selections. Electron candidates are required to satisfy the Medium quality criteria and not the Tight quality criteria. Muon candidates must pass a relaxed, but fail the tight, isolation requirement. All other selections remain unchanged. More than 80% of the selected events in MJCRs are estimated to originate from multijet production. These MJCR samples are used to model the kinematics of multijet contributions in their corresponding CRs and SRs, after subtracting contributions from other sources. The multijet scale factors, the ratios of the multijet contributions in the CRs to those in their MJCRs, are extracted through fits to the \(E_{\text {T}}^{\text {miss}} \) distributions in CRs using the multijet \(E_{\text {T}}^{\text {miss}} \) distribution shapes in MJCRs as templates. In the fits, contributions from other sources are constrained to their expectations from MC simulations within their uncertainties. These scale factors are then applied to their corresponding SRs to estimate multijet contributions.

7 Systematic uncertainties

Systematic uncertainties impact the search sensitivity through their effects on background estimations, signal selection efficiencies, and the distributions of the mass discriminants. The sources of these uncertainties can be classified broadly into two groups: (a) those experimental in nature related to the detector and reconstruction performance and (b) those of theoretical origins associated with the MC modelling of both the background and signal processes. The uncertainties and the methods used to evaluate them are discussed below. Unless explicitly stated, the uncertainties quoted are the uncertainties in the quantities themselves, not the impact on the search sensitivity.

7.1 Experimental uncertainties

Experimental uncertainties arise from the luminosity, triggers, and reconstruction and identification of leptons and jets, as well as the calculation of the \(E_{\text {T}}^{\text {miss}}\). They also include uncertainties in the energy and momentum scales and resolutions of leptons and jets.

The uncertainty of the combined 2015–2018 integrated luminosity is 1.7%. It is derived from the calibration of the luminosity scale using xy beam-separation scans, following a methodology similar to that detailed in Ref. [27], and using the LUCID-2 detector for the baseline luminosity measurement [83]. A variation in the pile-up reweighting of MC events is included to cover the uncertainty in the ratio of the predicted and measured inelastic cross-sections [84].

Uncertainties in the efficiencies of lepton triggers are found to be negligible. The modelling of the electron and muon reconstruction, identification and isolation efficiencies is studied with a tag-and-probe method using \(Z\rightarrow \ell \ell \) events in data and simulation [62, 63]. Small corrections are applied to the simulation to better model the performance seen in data. These corrections have associated uncertainties of the order of 1%. Uncertainties in the lepton energy (or momentum) scale and resolution, especially for muon momentum resolution (3%), are also taken into account.

Uncertainties for the energy scale and resolution of the small-R jets are determined using MC simulation and in situ techniques [66]. For central jets, the total relative uncertainty in the jet energy scale varies in the range 1–4% for \(p_{\text {T}} >20\) GeV. For forward jets, additional 2–4% uncertainty depending on \(p_{\text {T}}\) is applied based on \(\eta \)-intercalibration study. The uncertainty in the jet energy resolution ranges from \(20\%\) for jets with a \(p_{\text {T}} \) of 20 GeV to less than \(5\%\) for jets with \(p_{\text {T}} >200\) GeV.

Uncertainties in the scale of the large-R jet \(p_{\text {T}}\) are estimated by comparing the calorimeter- and track-based energy and mass measurements in data and simulation [85]. The precision of the relative jet energy scale is 1–2% for \(200\,\hbox {GeV}<p_T< 2\) \(\text {TeV}\), while that of the mass scale is 2–10%. The jet energy resolution uncertainty is estimated to be approximately 2%. The efficiency of the W/Z boson tagging based on the \(m_J\) and \(D_{2}^{(\beta =1)}\) requirements is estimated using data control samples, following the technique described in Ref. [86]. The efficiency for large-R jets from W/Z boson decays is estimated using \(t{\bar{t}}\) control samples for \(p_{\text {T}} < 600\) GeV. The measurement is extrapolated to the higher \(p_{\text {T}}\) region with additional uncertainties estimated from simulations [87]. The efficiency for background large-R jets from gluons or light quarks is estimated using dijet and \(\gamma \)+jets samples.

Uncertainties in the efficiencies for tagging b-jets and for mis-tagging light-flavour jets are determined from \(t\bar{t} \) control samples [68, 88, 89]. The total uncertainties are 1–10%, 15–50%, and 50–100% for b-jets, c-jets, and light-flavour jets respectively.

Uncertainties in the \(E_{\text {T}}^{\text {miss}}\) trigger efficiencies have negligible impact on the search as the efficiencies for the selected signal events are high. The uncertainty in \(E_{\text {T}}^{\text {miss}}\) is calculated from those in the energy scales and resolutions of leptons and jets as well as those in the energy deposits unassociated with any identified physics objects [73].

Multijet backgrounds are only important for the resolved analysis in the 1-lepton channel and are estimated using data control regions. The dominant uncertainties are from the multijet \(E_{\text {T}}^{\text {miss}}\) and mass templates, obtained from MJCRs after subtracting \(W\!+\!\mathrm {jets}\) and \(t\bar{t}\) contributions. They are estimated by varying the \(W\!+\!\mathrm {jets}\) and \(t\bar{t}\) subtractions and are found to range from a few percent to up to 15%.

7.2 Theoretical uncertainties

Theoretical uncertainties affect the normalisations of diboson and single-top-quark backgrounds, the shapes of mass distributions of background processes, and the signal acceptances. They arise from sources such as the choices of event generators, parton distribution functions (PDFs), parton shower models, and underlying-event tunes. Modelling uncertainties in the shapes of the mass distributions are estimated by varying the renormalisation/factorisation scales, PDF set and \(\alpha _\mathrm {s}\) values used in the nominal MC samples. Alternative generators are used to estimate the uncertainties due to the choices of generators, parton shower models and event tunes.

Background contributions from diboson and single-top-quark processes are estimated from MC simulations and are normalised to their theoretical cross-sections. For the diboson process, the cross-section uncertainty is estimated to be 10% [61, 90]. An additional contribution from electroweak production, simulated with Madgraph5_aMC@NLO+Pythia8, leads to an increase in the normalisation of the diboson background for the VBF process by a factor of 1.60 (1.85) in the resolved (merged) analyses. A uncertainty of 50% is applied to the normalisation of the electroweak diboson contribution. The impact on the ggF/DY analysis is negligible. For the cross-section of single-top-quark processes, an uncertainty of 20% is assumed [91].

Background contributions from \(V\!+\!\mathrm {jets}\) and \(t\bar{t}\) are normalised using data control regions in the 1-lepton and 2-lepton channels. Their overall normalisations are free parameters in the likelihood fit (Sect. 8) and thus only uncertainties in the shapes of discriminant variables are considered. For \(V\!+\!\mathrm {jets}\), the nominal Sherpa samples are compared with samples produced using MadGraph5_aMC@NLO. Moreover, the resummation scale and the CKKW [48, 49] matching scale in the nominal samples are also varied. The shape systematic uncertainty varies the background expectation in each bin and it is typically smaller than 10%, with the Sherpa-MadGraph comparison reaching 25% at the highest mass bin in the merged ggF/DY WZ untagged signal regions for the 1-lepton channel. For \(t\bar{t}\), the default Powheg-Box sample is compared with the alternative MadGraph5_aMC@NLO sample interfaced with Pythia 8.230. The difference is found to be approximately 4% in the merged signal regions, twice the value in the resolved signal regions. The difference between the Pythia 8.230 sample using the A14 tune and the alternative Herwig 7.04 [92, 93] sample using the H7UE set of tuned parameters [93] and the MMHT2014LO PDF set [94] is found to be between 2 and 5% in the various mass bins. The changes resulting from varying the parameter values for the nominal generator are less than 5%. In the 0-lepton channel, there is no pure control region to evaluate the \(V\!+\!\mathrm {jets}\) and \(t\bar{t}\) background, so the normalisation factors for the 0-lepton channel are assumed to be the same as for the 1-lepton channel (\(W\!+\!\mathrm {jets}\) and \(t\bar{t}\)) and 2-lepton channel (\(Z\!+\!\mathrm {jets}\)). Systematic uncertainties in this normalisation are obtained by the data/prediction double ratio between the default and the alternative MC generator and is estimated to be between 10 and 20% for \(V\!+\!\mathrm {jets}\) and up to 30% for \(t\bar{t}\). The \(t\bar{t}\) background is negligible in the 2-lepton channel and therefore its uncertainty is not considered for this channel.

Uncertainties in the signal acceptances are estimated for the choice of PDF set and the modelling of initial- and final-state radiation. The PDF uncertainties are estimated by taking the acceptance difference due to applying internal PDF error sets and the difference due to choosing different PDF sets. The uncertainty due to ISR/FSR modelling is studied by varying parameter values in the tunes used and applied to the HVT, the RS graviton, and the RS radion models. These uncertainties, calculated for several resonant mass points, are retrieved for each model, production process and decay. The PDF uncertainties are evaluated to be under 5% for all models. ISR/FSR uncertainties range from 2% for the merged analysis of ggF HVT production to about 11% for the resolved analysis of VBF HVT production.

Table 3 Dominant relative uncertainties in the best-fit signal-strength parameter \(\mu \) of hypothesised signal production of ggF RS graviton with \(m(G_{\mathrm {KK}})=600\) GeV and \(m(G_{\mathrm {KK}})=2\) \(\text {TeV}\). For this study, the RS graviton production cross-section is assumed to be 100 fb at 600 GeV and 2 fb at 2 \(\text {TeV}\), corresponding to approximately the expected median upper limits at these two mass values. Uncertainties with smaller contributions are not included

7.3 Impact of systematic uncertainties

The effects of systematic uncertainties on the search are studied for hypothesised signals using the signal-strength parameter \(\mu \), the ratio of the extracted cross-section (Sect. 8) to the injected hypothesised signal cross-section. The expected relative uncertainties in the best-fit \(\mu \) value from the leading sources of systematic uncertainties are shown in Table 3 for the ggF production of an RS graviton with \(m(G_{\mathrm {KK}})=600\) GeV and 2 \(\text {TeV}\). Apart from the statistical uncertainties in the data, the uncertainties with the largest impact on the sensitivity of the searches are from the sizes of the MC samples, measurements of small-R and large-R jets, background normalisations and modellings. Uncertainties related to the jet measurements, such as jet energy scale and resolution, affect the search primarily through their impacts on the shapes of the discriminant mass distributions of both signal and background processes. Uncertainties on the normalisations of background contributions estimated using CRs arise from CR statistics as well as MC event generators used to extroplate from CRs to SRs. Background modelling uncertainties include uncertainties on their normalisations, if estimated from MC simulations, as well as on the shapes of the mass distributions. The normalisations are affected by the uncertainties on the theoretical cross sections and on the luminosity. The shapes are affected by, in addition to experimental sources, theoretical sources such as PDF, ISR/FSR, and MC generator etc. For signals with higher mass, the data statistical uncertainty becomes dominant. The effects of systematic uncertainties for the other searches are similar.

8 Results and interpretations

8.1 Statistical procedure

The statistical analysis is based on the framework described in Refs. [95,96,97]. A profile-likelihood-ratio test statistic is used to test the compatibility of the background-only hypothesis and the observed data, and to test the signal-plus-background hypothesis for the production of a heavy resonance X, with its production cross-section in the VV decay mode, \(\sigma (pp\rightarrow X\rightarrow VV)\), as the parameter of interest. Maximum-likelihood fits are made to the observed binned distributions of the final discriminants in SRs, \(m_{\mathrm {T}}\) in 0-lepton , \(m_{\ell \nu J}\) or \(m_{\ell \nu jj}\) in 1-lepton and \(m_{\ell \ell J}\) or \(m_{\ell \ell jj}\) in 2-lepton , and to the numbers of observed events in CRs simultaneously. The mass ranges fitted are 300–3000 GeV for the resolved analysis and 500–6000 GeV for the merged analysis. The normalisations of the \(V\!+\!\mathrm {jets}\) and \(t\bar{t} \) contributions are free parameters in these fits and are constrained by the data in both the CRs and SRs. Systematic uncertainties, described in Sect. 7, and their correlations are incorporated as constraints into the likelihood calculations through nuisance parameters, where each is given a prior distribution based on individual studies or is allowed to float freely, constrained simultaneously by the SRs and CRs.

Fig. 9
figure 9

Comparisons of the observed data and the expected background distributions of \(m_{\mathrm {T}}\) in the 6 ZZ SRs of the 0-lepton channel. The background predictions are obtained through a background-only simultaneous fit to the 6 WW and 15 ZZ SRs and their respective \(V\!+\!\mathrm {jets}\) and \(t\bar{t}\) CRs (see text). The bottom panes show the ratios of the observed data to the background predictions. The blue triangles indicate bins where the ratio is non-zero and outside the vertical range of the plot. The hatched bands represent the uncertainties in the total background predictions, combining statistical and systematic contributions

Fig. 10
figure 10

Comparisons of the observed data and the expected background distributions of \(m_{\ell \nu jj}\) or \(m_{\ell \nu J}\) in the 6 WW SRs of the 1-lepton channel. The background predictions are obtained through a background-only simultaneous fit to the 6 WW and 15 ZZ SRs and their respective \(V\!+\!\mathrm {jets}\) and \(t\bar{t}\) CRs (see text). The bottom panes show the ratios of the observed data to the background predictions. The blue triangles indicate bins where the ratio is non-zero and outside the vertical range of the plot. The hatched bands represent the uncertainties in the total background predictions, combining statistical and systematic contributions

Fig. 11
figure 11

Comparisons of the observed data and the expected background distributions of \(m_{\ell \ell jj}\) or \(m_{\ell \ell J}\) in the 9 ZZ SRs of the 2-lepton channel. The background predictions are obtained through a background-only simultaneous fit to the 6 WW and 15 ZZ SRs and their respective \(V\!+\!\mathrm {jets}\) and \(t\bar{t}\) CRs (see text). The bottom panes show the ratios of the observed data to the background predictions. The blue triangles indicate bins where the ratio is non-zero and outside the vertical range of the plot. The hatched bands represent the uncertainties in the total background predictions, combining statistical and systematic contributions

Two types of fits, referred to as the \(WW+ZZ\) and WZ fits below, are performed. The \(WW+ZZ\) fits include all 21 SRs of the \(X\rightarrow WW\) and \(X\rightarrow ZZ\) searches and the WZ fits includes the 19 SRs of the \(X\rightarrow WZ\) search, along with their respective CRs. Separate fits are performed for the ggF/DY and VBF production modes and for different resonance mass hypotheses, but including SRs and CRs in both the ggF/DY and VBF categories. The \(WW+ZZ\) fits are used to search for the RS radion and RS graviton signals as both the WW and ZZ decay modes are expected from these resonances. The fits are also used to search for HVT \(Z'\rightarrow WW\) production. In this case, the \(X\rightarrow ZZ\) SRs effectively become additional CRs for the search. The WZ fits are used to search for HVT \(W'\rightarrow WZ\) production.

Table 4 The expected background events with contributions from individual sources in 6 WW and 15 ZZ SRs compared with the data. The backgrounds are estimated from a background-only simultaneous fit to all WW and ZZ SRs and their corresponding CRs

8.2 Data and background comparisons

To test the compatibility of the data and the background expectations, the data are first fit to the background-only hypothesis for both the \(WW+ZZ\) and WZ fits. Good agreement is found between the observed mass distributions and the estimated post-fit background contributions in all SRs. As examples, the data are compared with the expected backgrounds from the \(WW+ZZ\) fit in Fig. 9 for the \(m_{\mathrm {T}}\) distributions of the 0-lepton channel, in Fig. 10 for the \(m_{\ell \nu jj}/m_{\ell \nu J}\) distributions of the 1-lepton channel, and in Fig. 11 for the \(m_{\ell \ell jj}/m_{\ell \ell J}\) distributions of the 2-lepton channel. The largest excess is seen at \(m_{\mathrm {T}}\sim 1.5\) \(\text {TeV}\) in Fig. 9a. This excess is estimated to have a local significance of 2.8 standard deviations when modelled using RS radion production. The differences between the pre- and post-fit background estimates are less than 10% for a majority of the bins of the mass distributions. Larger differences, but comparable with the sizes of the statistical uncertainties, are observed for some bins with small numbers of events. Moreover the fits do not significantly alter the shapes of the background distributions for most of the signal regions. The largest change in the shape, up to 15%, is seen in the mass distribution of the VBF merged HP signal region in the 1-lepton channel. The difference is traced primarily to the difference in the simulated \(V\!+\!\mathrm {jets}\) samples produced with the Sherpa and MadGraph5 programs. Table 4 shows the post-fit estimated background event yields from different sources in all WW and ZZ SRs compared with the numbers of observed events in data. A similar level of agreement is obtained for the WZ fit.

Fig. 12
figure 12

Observed (black solid curve) and expected (black dashed curve) 95% CL upper limits on the a ggF and b VBF production cross-section of an RS radion at \(\sqrt{s}=13\) \(\text {TeV}\) in its diboson (WW and ZZ) decay mode as functions of the RS radion mass, combining the \(R\rightarrow WW\) and \(R\rightarrow ZZ\) searches in the three leptonic channels. The theoretical prediction for \(\mathcal {B}(R\rightarrow WW)/\mathcal {B}(R\rightarrow ZZ)\) is assumed for the combination of the WW and ZZ decay modes. The green (inner) and yellow (outer) bands represent \(\pm 1\sigma \) and \(\pm 2\sigma \) uncertainty in the expected limits. Limits expected from individual leptonic channels (dot-dashed curves in blue, magenta, and brown) are also shown for comparison. Limits are calculated in the asymptotic approximation below 3 (1) \(\text {TeV}\) and are obtained from pseudo-experiments above that for ggF (VBF) production. Theoretical predictions (red solid curve) as functions of the RS radion mass are overlaid

8.3 Limits on the production of heavy resonances

Constraints on the production of heavy resonances are derived by repeating the fit to the signal-plus-background hypothesis for different signal models. Upper limits on cross-sections are calculated with a modified frequentist method [98], also known as CL\(_\mathrm {s}\), using the \(\tilde{q}_{\mu }\) test statistic in the asymptotic approximation [99] for resonance masses below 3 (1) \(\text {TeV}\) for ggF/DY (VBF) production. For heavier resonances, the low number of events makes the approximation unreliable and the limits are obtained through pseudo-experiments instead. The observed (expected) limits from pseudo-experiments are approximately 20–30% (10–20%) higher than those from the asymptotic calculations at the highest resonance masses in the search.

8.3.1 Limits on the production of RS radions

Upper limits on the production cross-section of an RS radion in its VV decay modes, \(\sigma (pp\rightarrow R\rightarrow VV)\), are obtained by combining the \(R\rightarrow WW\) and \(R\rightarrow ZZ\) searches in the three leptonic final states. The 1-lepton channel is sensitive to the \(R\rightarrow WW\) decay while the 0-lepton and 2-lepton channels are sensitive to the \(R\rightarrow ZZ\) decay. The limits are derived separately for the ggF and VBF processes through the \(WW+ZZ\) fits for different RS radion mass hypotheses. Figure 12 shows the observed and expected upper limits at 95% confidence level (CL) as functions of its mass for both the ggF and VBF processes. The observed limits for the VBF process are noticeably higher than the expected limits around an RS radion mass of 1.5 \(\text {TeV}\), reflecting the excess seen in the \(m_{\mathrm {T}}\) distribution from the merged HP signal region of the 0-lepton channel. The observed (expected) combined limit on \(\sigma (pp\rightarrow R\rightarrow VV)\) varies from 1.8 (3.3) pb at 300 GeV to 0.38 (0.43) fb at 5 \(\text {TeV}\) for the ggF production process and from 0.60 (1.15) pb at 300 GeV to 0.23 (0.26) fb at 5 \(\text {TeV}\) for the VBF production process. These observed (expected) upper limits exclude the ggF production of an RS radion with a mass below 3.2 (2.9) \(\text {TeV}\) while no mass exclusion can be derived for the VBF production.

Except for masses below approximately 1 \(\text {TeV}\), the 1-lepton channel dominates the combined search sensitivities for both the ggF and VBF processes. This is not surprising as \(\mathcal {B}(R\rightarrow WW)\) is \(\approx 2\times \mathcal {B}(R\rightarrow ZZ)\) in the RS radion model. The 2-lepton channel is the most sensitive at low masses, benefiting from its good mass resolution and small background contributions, and is the least sensitive at high masses because of its small expected signal yields. The 0-lepton channel has the worst sensitivity at low masses and contributes non-negligibly at high masses.

For the RS radion search here as well as the HVT and RS graviton searches presented below, the resolved analysis is important for masses below 600 GeV while the merged analysis dominates the search for higher masses. The LP signal regions improve the cross-section sensitivities of the merged analysis by approximately 5% for the entire mass range of the search. The missing \(X\rightarrow WV\rightarrow \ell \nu qq,\tau \nu qq\) signal contribution in the 0-lepton channel as discussed in Sect. 5.6, if included, would lower the expected cross-section upper limits by up to 4% at 2 \(\text {TeV}\) and up to 10% at 5 \(\text {TeV}\).

Fig. 13
figure 13

Observed (black solid curve) and expected (black dashed curve) 95% CL upper limits on the a DY and b VBF production cross-section of an HVT \(W'\) boson at \(\sqrt{s}=13\) \(\text {TeV}\) in the WZ decay mode as functions of its mass, combining searches in the three leptonic channels. The green (inner) and yellow (outer) bands represent \(\pm 1\sigma \) and \(\pm 2\sigma \) uncertainty in the expected limits. Limits expected from individual leptonic channels (dot-dashed curves in blue, magenta, and brown) are also shown for comparison. Limits are calculated in the asymptotic approximation below 3 (1) \(\text {TeV}\) and are obtained from pseudo-experiments above that for DY (VBF) production. Theoretical predictions as functions of the \(W'\) boson mass are overlaid in a for Model A (red solid curve) and Model B (red dotted curve) and in b for Model C (red solid curve)

Fig. 14
figure 14

Observed (black solid curve) and expected (black dashed curve) 95% CL upper limits on the a DY and b VBF production cross-section of an HVT \(Z'\) boson at \(\sqrt{s}=13\) \(\text {TeV}\) in the WW decay mode as functions of its mass from the search in the 1-lepton channel. The green (inner) and yellow (outer) bands represent \(\pm 1\sigma \) and \(\pm 2\sigma \) uncertainty in the expected limits. Limits are calculated in the asymptotic approximation below 3 (1) \(\text {TeV}\) and are obtained from pseudo-experiments above that for DY (VBF) production. Theoretical predictions as functions of the \(Z'\) boson mass are overlaid in a for Model A (red solid curve) and Model B (red dotted curve) and in b for Model C (red solid curve)

8.3.2 Limits on the production of HVT resonances

Upper limits on the production cross-sections of HVT \(W'\) and \(Z'\) bosons in their WZ and WW decays are obtained through the WZ and \(WW+ZZ\) fits, respectively. All leptonic channels contribute to the \(W'\rightarrow WZ\) search while only the 1-lepton channel contributes to the \(Z'\rightarrow WW\) search. Limits as functions of resonance masses are shown in Figs. 13 and 14 for \(W'\rightarrow WZ\) and \(Z'\rightarrow WW\), respectively, for both DY and VBF processes. The theoretical predictions of the HVT Model A, Model B, and Model C are also shown for comparison. The observed (expected) \(\sigma (pp\rightarrow W'\rightarrow WZ)\) limit ranges from 1.9 (2.5) pb at 300 GeV to 0.16 (0.17) fb at 5 \(\text {TeV}\) for DY production and from 1.3 (1.8) pb at 300 GeV to 0.35 (0.51) fb at 4 \(\text {TeV}\) for VBF production. These observed (expected) limits exclude an HVT \(W'\) boson produced in the DY process lighter than 3.9 (3.8) \(\text {TeV}\) for Model A and 4.3 (4.0) \(\text {TeV}\) for Model B, but fail to exclude any mass region in the VBF process for the benchmark Model C. For both production processes, the 2-lepton channel is the most sensitive for masses up to \(\sim 1.5\) \(\text {TeV}\). At high masses, the sensitivity is dominated by the 0-lepton and 1-lepton channels.

Only the 1-lepton channel contributes to the \(Z'\rightarrow WW\) search. The observed (expected) \(\sigma (pp\rightarrow Z'\rightarrow WW)\) limit ranges from 0.9 (2.7) pb at 300 GeV to 0.18 (0.18) fb at 5 \(\text {TeV}\) for the DY process and from 1.36 (3.15) pb at 300 GeV to 0.25 (0.32) fb at 4 \(\text {TeV}\) for the VBF process. These limits exclude an HVT \(Z'\) boson lighter than 3.5 (3.4) \(\text {TeV}\) for Model A and 3.9 (3.7) \(\text {TeV}\) for Model B in the DY process.

Fig. 15
figure 15

Observed (black solid curve) and expected (black dashed curve) 95% CL upper limits on the a ggF and b VBF production cross-section of an RS graviton at \(\sqrt{s}=13\) \(\text {TeV}\) in its diboson (WW and ZZ) decay mode as functions of its mass, combining the searches for the \(G_{\mathrm {KK}}\rightarrow WW\) and \(G_{\mathrm {KK}}\rightarrow ZZ\) decays in the three leptonic channels. The theoretical prediction for \(\mathcal {B}(G_{\mathrm {KK}}\rightarrow WW)/\mathcal {B}(G_{\mathrm {KK}}\rightarrow ZZ)\) is assumed for the combination of the WW and ZZ decay modes. The green (inner) and yellow (outer) bands represent \(\pm 1\sigma \) and \(\pm 2\sigma \) uncertainty in the expected limits. Limits expected from individual leptonic channels (dashed curves in blue, magenta, and brown) are also shown for comparison. Limits are calculated in the asymptotic approximation below 3 (1) \(\text {TeV}\) and are obtained from pseudo-experiments above that for ggF (VBF) production. Theoretical predictions (red solid curves) for the chosen model as functions of the RS graviton mass are overlaid

Table 5 Observed (expected) 95% CL lower limits on the mass, in \(\text {TeV}\), of different resonances in the benchmark models studied. The symbol “−” means no limit is set

8.3.3 Limits on the production of RS gravitons

Upper limits on the production cross-section of an RS graviton in its VV decays, \(\sigma (pp\rightarrow G_{\mathrm {KK}}\rightarrow VV)\), are obtained following the same procedure used to derive the RS radion limits. The observed and expected upper limits as functions of its mass for both the ggF and VBF processes are shown in Fig. 15. The observed (expected) limit varies from 1.4 (3.6) pb at 300 GeV to 0.26 (0.28) fb at 5 \(\text {TeV}\) for the ggF production process and from 0.40 (0.61) pb at 300 GeV to 0.30 (0.33) fb at 5 \(\text {TeV}\) for the VBF production process. Compared with theory cross-sections, with the benchmark parameters, the observed (expected) upper limits at 95% CL exclude the production of an RS graviton lighter than 2.0 (2.2) \(\text {TeV}\) in the ggF process and lighter than 0.76 (0.77) \(\text {TeV}\) in the VBF process.

Similar to the RS radion case, the 1-lepton channel contributes dominantly to the combined search sensitivities at high masses while the 2-lepton channel is slightly more sensitive than the 1-lepton channel for masses below \(\sim 1\) \(\text {TeV}\). The 0-lepton channel is the least sensitive at low masses, but provides significant contributions at high masses.

8.4 Comparisons of the limits

Table 5 summarises the observed and expected 95% CL lower limits on the masses of the resonances in the benchmark models studied in this paper. These mass limits and the cross-section upper limits presented above are significantly more stringent than those published previously from similar searches. The cross-section upper limits are a factor of three or more lower than those of the search using the same data set, but in the hadronic \(VV\rightarrow JJ\) final state [16]. Compared to the searches with the 36.1 \(\mathrm{fb}^{1}\) data set in the same leptonic final states [22, 23], an improvement of a factor of three or more in the cross-section upper limits is also obtained for most of the searches at the highest masses. The observed lower mass limits of this search for the HVT \(W'\) and \(Z'\) as well as for the RS graviton are also similar to those from the combinations of searches in the fully-leptonic, semi-leptonic, and fully-hadronic final states from smaller datasets of \(\sim 36\) \(\mathrm{fb}^{1}\)  [12, 15].

9 Summary

Searches for the production of heavy diboson resonances are performed using the proton–proton collision data with an integrated luminosity of 139 \(\mathrm {fb^{-1}} \) at \(\sqrt{s}=13\) \(\text {TeV}\). The data were recorded by the ATLAS experiment between 2015 and 2018 at the LHC. The WW, ZZ and WZ decay modes of the heavy resonances in the 0-lepton , 1-lepton and 2-lepton final states are considered. The data are found to be in good agreement with background expectations. Upper limits on the production of heavy resonances in the mass range 300–5000 GeV through gluon–gluon fusion, Drell–Yan or vector-boson fusion processes are derived for Standard Model extensions with an additional neutral scalar, a heavy vector triplet, or warped extra dimensions.

Combining the WW and ZZ decay modes, the observed 95% confidence-level upper limit on \(\sigma (pp\rightarrow X\rightarrow VV)\) for the ggF (VBF) process ranges from 1.8 (0.60) pb at 300 GeV to 0.38 (0.23) fb at 5 \(\text {TeV}\) for an RS radion and from 1.4 (0.40) pb at 300 GeV to 0.26 (0.30) fb at 5 \(\text {TeV}\) for an RS graviton. These observed limits set lower mass limits of 3.2 \(\text {TeV}\) for the ggF production of an RS radion, and 2.0 (0.76) \(\text {TeV}\) for the ggF (VBF) production of an RS graviton.

For the production of \(W'\) and \(Z'\) bosons in the HVT framework, the observed upper limit on \(\sigma (pp\rightarrow W'\rightarrow WZ)\) varies from 1.9 pb at 300 GeV to 0.16 fb at 5 \(\text {TeV}\) for DY production and from 1.3 pb at 300 GeV to 0.35 fb at 4 \(\text {TeV}\) for VBF production. Similarly, the limits on \(\sigma (pp\rightarrow Z'\rightarrow WW)\) are observed to vary from 0.9 pb at 300 GeV to 0.18 fb at 5 \(\text {TeV}\) for DY production and from 1.36 pb at 300 GeV to 0.25 fb at 4 \(\text {TeV}\) for VBF production. In the benchmark Model A (Model B), these cross-section upper limits exclude the ggF production of a \(W'\) boson with \(m(W')<3.9\ (4.3)\) \(\text {TeV}\) and a \(Z'\) boson with \(m(Z')<3.5\ (3.9)\) \(\text {TeV}\).