Search for heavy resonances decaying into two Higgs bosons or into a Higgs boson and a W or Z boson in proton-proton collisions at 13 TeV

A search is presented for massive narrow resonances decaying either into two Higgs bosons, or into a Higgs boson and a W or Z boson. The decay channels considered are HH→bb¯τ+τ−\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \mathrm{H}\mathrm{H}\to \mathrm{b}\overline{\mathrm{b}}{\tau}^{+}{\tau}^{-} $$\end{document} and VH→qq¯τ+τ−\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \mathrm{V}\mathrm{H}\to \mathrm{q}\overline{\mathrm{q}}{\tau}^{+}{\tau}^{-} $$\end{document}, where H denotes the Higgs boson, and V denotes the W or Z boson. This analysis is based on a data sample of proton-proton collisions collected at a center-of-mass energy of 13 TeV by the CMS Collaboration, corresponding to an integrated luminosity of 35.9 fb−1. For the TeV-scale mass resonances considered, substructure techniques provide ways to differentiate among the hadronization products from vector boson decays to quarks, Higgs boson decays to bottom quarks, and quark- or gluon-induced jets. Reconstruction techniques are used that have been specifically optimized to select events in which the tau lepton pair is highly boosted. The observed data are consistent with standard model expectations and upper limits are set at 95% confidence level on the product of cross section and branching fraction for resonance masses between 0.9 and 4.0 TeV. Exclusion limits are set in the context of bulk radion and graviton models:spin-0 radion resonances are excluded below a mass of 2.7 TeV at 95% confidence level. In the spin-1 heavy vector triplet framework, mass-degenerate W′ and Z′ resonances with dominant couplings to the standard model gauge bosons are excluded below a mass of 2.8 TeV at 95% confidence level. These are the first limits for massive resonances at the TeV scale with these decay channels at s=13\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \sqrt{s}=13 $$\end{document} TeV.


Introduction
Heavy resonances that decay to HH, VV, or VH, where H denotes the Higgs boson, and V denotes a W or Z boson, are motivated by theories beyond the standard model (SM) that address the large difference between the electroweak and gravitational scales.These heavy particles arise as Kaluza-Klein (KK) excitations of spin-0 radions [1][2][3], and as spin-2 gravitons predicted in models based on Randall-Sundrum warped extra dimensions [4,5], with the gravitons propagating in the entire five-dimensional bulk [6][7][8].Heavy spin-1 W and Z particles that decay to VV and VH are also postulated in composite Higgs models [9][10][11][12], little Higgs models [13,14], and in the sequential SM (SSM) [15].The models containing new spin-1 states are generalized in the heavy vector triplet (HVT) framework [16].All of the new hypothetical particles with spins of 0, 1, or 2 can be produced at the CERN LHC, via the processes depicted in the Feynman diagrams of Fig. 1.
Feynman diagrams for the production of a spin-0 radion or a spin-2 graviton X that decays to two Higgs bosons (left), and the production of a heavy vector boson V (W or Z ) that decays to a vector boson and a Higgs boson (right).
The bulk graviton model has two free parameters: the mass of the first KK excitation of the spin-2 boson, denoted as the KK bulk graviton, and the ratio k ≡ k/M Pl , where k is the un- known curvature scale of the extra dimension and M Pl ≡ M Pl / √ 8π is the reduced Planck mass.Searches for radions in this model can be described in terms of the radion mass and the ultraviolet cutoff of the theory Λ R [17].The HVT model is formulated in terms of four parameters: the mass of the new vector bosons, their coupling coefficient to fermions c F , their coupling coefficient to the Higgs boson and longitudinally-polarized SM vector bosons c H , and the strength of the new vector boson interaction g V .In the HVT framework two scenarios are considered and henceforth referred to as model A and model B, depending on the couplings of the new physics resonance to the SM particles.In model A (g V = 1, c H = −0.556,c F = −1.316), the coupling strengths to the SM bosons and fermions are comparable and the new particles decay primarily to fermions.In model B (g V = 3, c H = −0.976,c F = 1.024), the coupling to the SM fermions are small, and the branching fraction of the new resonance to SM bosons is nearly 100% [16].Searches for diboson resonances have previously been performed in several final states, placing lower limits on the masses of such resonances above the TeV scale [18][19][20][21][22][23][24][25][26][27][28][29][30][31][32].This paper presents a search for resonances with masses above 900 GeV decaying into HH or VH.The final states analyzed are HH → bbτ + τ − and VH → qqτ + τ − .Events are classified as "semileptonic", denoted as τ h , if one τ lepton decays into a lighter lepton ( ), denoting either a muon or an electron, and two neutrinos, and the other decays hadronically (τ h ) into hadrons and a neutrino.Events are categorized as "fully-hadronic", denoted as τ h τ h , if both τ leptons decay hadronically.The analysis aims to reconstruct the diboson decay products in order to search for a narrow local enhancement in the diboson invariant mass spectrum with a width insignificant with respect to the experimental resolution.This search complements and significantly extends the reach of the CMS search for ZH → qqτ + τ − resonances with data collected at √ s = 8 TeV [27], by using a larger data sample collected at √ s = 13 TeV, a more complex categorization of the final states, and an increased number of signal models.
Since the resonances under study have masses of the order of a few TeV, the bosons resulting from their decays have transverse momenta (p T ) of at least several hundred GeV.As a consequence, the final decay products are collimated such that the hadronically decaying bosons cannot be resolved using standard jet algorithms.Dedicated techniques, called V tagging and H tagging, are applied to exploit the substructure of these large-cone jets to resolve the hadronically decaying vector and Higgs bosons.
For the Higgs boson decaying to a pair of τ leptons, the decay products are also close in angular separation.The τ lepton reconstruction and identification techniques described in Ref. [27] are therefore adopted to achieve optimal signal significance for this particular event topology.

The CMS detector
The central feature of the CMS apparatus is a superconducting solenoid of 6 m internal diameter, providing a magnetic field of 3.8 T. Within the solenoid volume are a silicon pixel and strip tracker, a lead tungstate crystal electromagnetic calorimeter (ECAL), and a brass and scintillator hadron calorimeter (HCAL), each composed of a barrel and two endcap sections.Forward calorimeters extend the pseudorapidity (η) coverage provided by the barrel and endcap detectors.Muons are measured in gas-ionization detectors embedded in the steel flux-return yoke outside the solenoid.The CMS two-level trigger system [33] reduces the event rate from the bunch crossing rate of 40 MHz down to around 1 kHz for data storage.
A more detailed description of the CMS detector, together with a definition of the coordinate system used and the kinematic variables, can be found in Ref. [34].

Data sample and simulation
The data sample analyzed in this search corresponds to an integrated luminosity of 35.9 fb −1 , collected in 13 TeV proton-proton collisions with the CMS detector during the 2016 data taking period.The signal processes pp → X → VH → qqτ + τ − and pp → X → HH → bbτ + τ − are simulated at leading order (LO) using the MADGRAPH5 aMC@NLO v2.2.2 [35] Monte Carlo (MC) event generator, for resonance masses between 900 and 4000 GeV, where the Higgs boson is forced to decay to τ pairs and the other boson to a pair of quarks.The signal processes where pp → X → HH → bbVV and pp → X → VH → qqVV are also considered, in which VV → 2 2ν or 2τ2ν, as they can yield final states similar to those of the primary signal process.The natural width of the resonance is assumed to be smaller than the experimental resolution of its reconstructed mass, as consistent with the benchmark radion, graviton and HVT models.
The SM background processes are generated using MC simulation.The Z/γ * +jets events and the W+jets events are simulated at LO with the MADGRAPH5 aMC@NLO generator.The POWHEG v2 generator is used to simulate tt and single top quark production at next-to-leading order [36][37][38][39].The LO PYTHIA 8.205 [40] generator is used for SM diboson (WW, WZ, or ZZ) and multijet events.For all signal and background samples, showering and hadronization are modeled using PYTHIA, τ lepton decays are described using TAUOLA 1.1.5[41], and the response of the detector is simulated using GEANT4 [42].Additional collisions in the same or adjacent bunch crossings (pileup) are superimposed onto the hard scattering processes, with the pileup vertex multiplicity distribution adjusted to match that of data.

Event reconstruction
The particle flow (PF) event algorithm [43] reconstructs and identifies each individual particle through an optimized combination of information from the various elements of the CMS detector.The energy of each electron is determined from the electron momentum as determined by the tracker, the energy of the corresponding ECAL cluster, and the energy sum of all bremsstrahlung photons spatially compatible with originating from the electron track.The energies of muons are obtained from the curvature of the corresponding tracks.The energies of charged hadrons are determined from a combination of their momentum measured in the tracker and the matching ECAL and HCAL energy deposits, corrected for zero-suppression effects and for the response function of the calorimeters to hadronic showers.Finally, the energies of neutral hadrons are obtained from the corresponding corrected ECAL and HCAL energies.
The identified particles are clustered into jets using the anti-k T algorithm [44], implemented in the FASTJET package [45].Two distance parameters are used in the analysis, 0.4 and 0.8, yielding jet collections referred to as AK4 and AK8 jets, respectively.The AK4 jets are used primarily to reject or select events with top quarks, while the larger AK8 jets are used to identify and contain hadronically decaying W, Z, and Higgs boson candidates.The charged-hadron subtraction (CHS) algorithm for mitigating pileup [43] discards charged particles not originating from the primary vertex (PV).The PV is defined as the vertex with the largest p 2 T sum of the physics objects.Here, the physics objects are the AK4 jets, and the associated missing transverse momentum, which is taken as the negative vector sum of the p T of the jets.The residual contamination from neutral pileup particles is estimated to be proportional to the event energy density and the jet area, and is removed from the jet energy calculation.The jet momentum is determined as the vectorial sum of all particle momenta in the jet, and is found from simulation to be within 5 to 10% of the true momentum over the entire p T spectrum and detector acceptance.Jet energy corrections are obtained from simulation, and are confirmed with in situ measurements of the energy balance in dijet, multijet, γ+jets, and leptonically decaying Z+jets events [46].Additional selection criteria are applied to each event to remove spurious jet-like features originating from isolated noise patterns in certain HCAL regions [47].The AK4 and AK8 jets must have p T > 20 and p T > 200 GeV, respectively, and |η| < 2.4 to be considered in the subsequent steps of the analysis.
To determine the jet mass and the substructure variables used in the identification of the hadronic decays of bosons, the so-called pileup per particle identification (PUPPI) algorithm [48] is applied to AK8 jets, instead of the CHS algorithm, to retain better stability of the substructure variables in events with a large amount of pileup.The PUPPI algorithm uses the local distribution of particles, event-pileup properties, and tracking information to compute a weight describing the likelihood for each particle to originate from a pileup interaction.The weights are used to rescale the particle four-momenta, superseding the need for further jet-based corrections.The PUPPI constituents are subsequently clustered with the same algorithm used for CHS jets, and then matched to the AK8 jets clustered with the CHS constituents that are then used for the kinematic selections.
Subsequently, the soft-drop algorithm [49,50], designed to remove contributions from soft radiation, is applied to the AK8 PUPPI jets.The soft-drop jet mass m j is defined as the invariant mass associated with the four-momentum of the soft-drop jet.Dedicated mass corrections, derived from data in a region enriched with tt events containing merged hadronic Wdecays, are applied to m j to remove any dependence on jet p T , and to match the jet mass and resolution observed in data.After the application of these corrections, the W(qq) jet mass resolution is measured to be 10% in the jet p T range considered.
The two-prong hadronic decays of Wand Z boson candidates are used to discriminate against jets initiated from single quarks and gluons.The constituents of the jet are clustered again using the k T algorithm [44], and the procedure is stopped when N subjets are obtained.Subsequently, the N-subjettiness as defined in Ref. [51] is calculated on the PUPPI-corrected jet for the one and two-subjet hypotheses as The normalization factor is given by the factor d 0 = ∑ k p T,k R 0 , where R 0 is the radius of the original jet, the index k increments over the jet constituents, and ∆R J,k = √ (∆η J,k ) 2 + (∆φ J,k ) 2 are the angular distances in terms of the differences in η and azimuth (φ) calculated between the kth jet constituent and the axes of the Jth subjet.Small values of the ratio of 2-subjettiness to 1-subjettiness, τ 21 = τ 2 /τ 1 , correspond to a high compatibility with the hypothesis that the jet is produced by two partons from the decay of a massive object, rather than arising from a single parton.The efficiency of the τ 21 selection is measured from data in a tt-enriched sample [52].
Jets originating from the dominant bb decays of Higgs bosons are likely to have two displaced vertices because of the long lifetime and large mass of the b quarks.Following the procedure above, jet clustering using the anti-k T algorithm is halted when two subjets are identified.The inclusive, combined secondary-vertex b tagging algorithm [53] is applied to the two subjets, which are considered as b-tagged if they pass a working point that provides a misidentification rate of ≈10% while maintaining an 85% efficiency.To remove backgrounds containing top quark decays, events with AK4 jets that do not overlap with the AK8 jet are subjected to a veto based on the same b tagging algorithm, but with a working point corresponding to an efficiency of ≈70% for identifying B hadrons and a ≈1% misidentification rate.The ratio of the b tagging efficiency determined from data and simulation is used as a scale factor to correct the simulated events.
A dedicated algorithm is used to reconstruct Higgs bosons decaying to one or two τ h candidates [54].The procedure begins by using the Cambridge-Aachen algorithm [55] with a distance parameter of 0.8 to identify jets with a large cone size, called CA8 jets.For each CA8 jet with p T > 100 GeV, the last step of the clustering is retracted, obtaining two subjets.If these subjets are found to have p T > 10 GeV and satisfy the mass drop-condition, which requires max(m subjet1 , m subjet2 )/m CA8jet < 2/3, the two subjets are used as seeds in the standard τ lepton reconstruction and the "hadron-plus-strips" algorithm [56] is applied to them to identify τ h candidates.If the p T and mass drop conditions are not met, the unclustering and identification procedures are repeated (iteratively) for the most energetic subjet.The τ h candidates selected through the hadron-plus-strips algorithm are then required to have |η| < 2.3 and p T > 20 GeV, and to satisfy a multivariate-discriminator threshold, obtained using a boosted decision tree technique [57].This algorithm is trained to discriminate between genuine τ h and generic jets, using variables related to energy deposits and track impact parameters that are correlated with the τ lepton lifetime.If no τ h candidates are identified with this method, then the procedure is repeated using AK4 jets as seeds, with similar selection requirements.The τ h candidate of highest p T is required to satisfy an isolation requirement that corresponds to a 50-60% efficiency in the considered topology.If two τ h candidates are identified, as in the τ h τ h channel, the isolation requirement on the τ h of second highest p T is relaxed to achieve a 70-80% efficiency.The probability to misidentify a CA8 jet as an H → τ + τ − decay is below 0.1%, after these selection criteria.
Electrons are reconstructed in the region |η| < 2.5 by matching energy deposits in the ECAL with tracks reconstructed in the tracker [58].Electron identification is based on the distribu-tion of energy deposited along the electron trajectory and the direction and momentum of the track in the inner tracker.Additional requirements are applied to remove electrons produced through photon conversions.Electrons are also required to be isolated from other particles in the detector, by imposing an upper threshold on the isolation parameter.The electron isolation parameter is defined as the magnitude of the p T sum of all the PF candidates (excluding the electron) within ∆R < 0.3 around the electron direction, after the contributions from pileup and particles associated with reconstructed τ h candidates within the isolation cone are removed.
Muons are reconstructed within the acceptance of the CMS muon system, |η| < 2.4, using information from both the muon spectrometer and the silicon tracker [59].Muon candidates are identified based on the compatibility of tracks reconstructed in the silicon tracker with tracks reconstructed from a combination of hits in both the tracker and the muon detector.In addition, the trajectory is required to be compatible with originating from the primary vertex, and to have a sufficient number of hits in the tracker and muon systems.Muons are required to be isolated by imposing a limit on the magnitude of the p T sum of all the PF candidates (excluding the muon) within ∆R < 0.4 around the muon direction, after the contributions particles associated with reconstructed τ h candidates within the isolation cone are removed.
The missing transverse momentum vector p miss T is defined as the negative vectorial sum of the momenta of all PF candidates associated with the primary vertex projected onto the plane perpendicular to the beam direction.The missing transverse momentum p miss T is defined as the magnitude of p miss T .The observable H miss T is defined as the magnitude of the vectorial p T sum of all AK4 jets with p T > 30 GeV and |η| < 3.0.

Event selection
Events are selected using a set of triggers that require p miss T or H miss T larger than 90 GeV, in combination with additional trigger requirements, such as the presence of a jet with p T > 80 GeV.The efficiency of the trigger for events subsequently satisfying the offline event selection, measured in an independent sample of events selected with muon triggers, is verified to be > 95%, with an uncertainty of 2%.
All events are required to contain one Higgs boson candidate decaying to τ h or τ h τ h .The other boson candidate is reconstructed as a jet, using the same kinematic criteria in all categories.Its soft-drop jet mass must be in the interval of 65-135 GeV.If the mass is in the range 65-85 GeV, the candidate is classified as a Wboson, if it is in the range 85-105 GeV it is classified as a Z boson, and if it is the range 105-135 GeV it is considered to be a Higgs boson.To discriminate against backgrounds, the jets are required to have small values of τ 21 , and events are divided into categories of high-purity (HP) if τ 21 < 0.4, and low-purity (LP) if 0.4 < τ 21 < 0.75.A jet is V tagged if it fulfills the soft-drop jet mass and τ 21 requirements.The normalization scale factors 0.99 ± 0.06 for the HP category and 0.96 ± 0.11 for the LP category [52] are applied to simulated events with genuine hadronic boson decays.Higgs boson jet candidates are classified according to the number of subjets (1 or 2) that pass the b tagging selection.Subjet b tagging is not used for jets compatible with Wor Z candidates and no N-subjettiness requirement is applied to the Higgs boson candidate jet.If both the N-subjettiness and the b tagging requirements are satisfied, the event is considered further.
Events are divided into categories depending on the number of identified τ h candidates (1 or 2), and on the classification of the large jet cone, i.e. either HP or LP in τ 21 , or either 1 or 2 b-tagged subjets.
Since the undetected neutrinos carry a significant fraction of the momentum in the ττ system, signal events are expected to have a large p miss T , thereby justifying the use of triggers that require large p miss T or H miss T .A stringent offline requirement of greater than 200 GeV is applied to the reconstructed p miss T , to ensure a stable trigger efficiency and to suppress the background contribution from multijet events.Events with top quark pairs or single top quarks are suppressed by removing events in which any AK4 jet not overlapping with the AK8 jet is b-tagged.
Several selection requirements are applied to remove SM backgrounds, such as meson an baryon resonances, Z+jets, W+jets and tt and single top quark production.The angular distance ∆R ττ should be smaller than 1.5, in order to reject W+jets events in which a jet misidentified as a τ lepton is typically spatially well-separated from the genuine lepton.To further increase the signal purity, the di-tau mass, as estimated from the SVFIT procedure [60][61][62], should be between 50 and 150 GeV.The SVFIT algorithm, based on a likelihood approach, estimates the di-τ system mass using the measured momenta of the visible decay products of both τ leptons, the reconstructed p miss T , and the p miss T resolution.
Finally, the resonance candidate mass m X , defined as the invariant mass of the H → ττ candidate and the hadronically decaying boson jet, is required to be larger than 750 GeV in order to ensure full trigger and reconstruction efficiencies.
After these selection requirements, the selection efficiency of a radion signal in the 1 and 2 b tag categories is 1-6% in the τ h τ h channel, and 3-10% in the τ h channel, for low and high resonance masses, respectively.The efficiencies for a W signal passing the V-tagging selection are 2-9% and 8-19% in the τ h τ h and τ h channels, respectively.Events in which H → VV are found contribute up to an additional 20-30% of the total signal yield in the τ h channels, and less than 10% in the τ h τ h channels.

Background estimation
The main sources of background originate from top quark pair production and from the production of a vector boson in association with jets (Z+jets and W+jets), while minor contributions arise from single top quark, diboson, and multijet production.These background contributions are split into either tt and single top quark production (tt, t), or into V+jets production.The latter includes Z+jets and W+jets, multijet, and SM diboson production.
The shape of the distribution of the top quark pair and single top quark background is determined from simulation, while the normalization is determined from data through dedicated control regions that are enriched in top quark events.Control regions having a purity larger than 80% for top quarks are selected by inverting the b tag veto on the AK4 jets and tightening the b tagging criteria.Events are separated according to the requirements of large-cone jet identification.Data are found to be well-described by simulations in terms of the jet and dijet resonance-mass distributions.Multiplicative scale factors are used to correct for the difference in normalization of data and simulation in the control regions, after subtracting the other background contributions.Scale factors obtained in control regions of τ h events are applied also to the τ h τ h channels, where there are fewer events.The normalization of top quark production processes in each region is also corrected using the scale factors reported in Table 1.
The estimated contribution from V+jets backgrounds is based on data, in regions defined by applying the complete signal selection apart from the jet-mass requirements.Data are divided into the τ h τ h and the τ h channels.Two jet-mass sidebands (SB) are defined with jet masses in the range of 30-65 GeV for the low sideband (LSB), or above 135 GeV for the high sideband (HSB), and used to predict the background contribution in the signal region (SR).Analytic functions are fitted to the simulated distributions of the jet mass, separately for V+jets and for top quark processes.The former is modeled by a third-order polynomial or an exponential convoluted with an error function, and the latter with a sum of one or two Gaussian functions to model the W and top quark distributions.The background model, composed of the sum of the functions modeling the two backgrounds, is fitted to the data in the jet mass sidebands.In this unbinned likelihood fit, only the V+jets distribution shape and normalization parameters are left free to float.The number of expected events in the SR is then obtained by integrating the background components in the jet mass window between 105 and 135 GeV.An indicative fit is shown in Fig. 2 (left) for the HP τ 21 category of τ h events.The procedure is repeated with an alternative function used for modeling the V+jets jet mass, with the observed difference in the normalization taken to be the associated systematic uncertainty.The expected number of background events in each signal region is reported in Table 2.
The distribution of the V+jets background resonance mass (m X ) in the SR is estimated from the SBs through a transfer function α(m X ), which accounts for the small kinematical differences and the correlations involved in the interpolation from the SBs to the SR, and does not depend on the systematic uncertainties that affect the simulated V+jets spectra, since they cancel out in the ratio.The α function is determined from simulations as where N MC,V+jets SR (m X ) and N MC,V+jets SB (m X ) are analytic functions fitted to the simulated m X distributions in the SR and SB regions, respectively.Depending on the category, the function can either be an exponential using one or two parameters, or a power law with one parameter.The distribution of the V+jets background in the SR is then estimated by fitting an analytic function to data in the SBs, after subtracting the top quark background estimated from simulation, and multiplying by the α(m X ) transfer function.The normalization of the V+jets is determined from the fit to the jet mass, as reported in Table 2.The resonance mass distribution is shown in Fig. 2 (right) for τ h events in the HP category.
The overall background in the SR is then expected to be: where N tt,t SB and N tt,t SR are the distributions for the top quark process in the SB and SR, respectively.The shape and the normalization of the distribution are fixed from simulation, with the latter corrected using the appropriate scale factors in Table 1.The data in the SR and the background predictions before and after the fit in the SR are shown in Figs. 3 and 4.

Systematic uncertainties
Systematic uncertainties, resulting from experimental and theory sources, may affect both the normalization and shape of the signal and background distributions.Figure 2: Soft-drop jet mass distribution in data in the HP τ h category, together with the background prediction (fitted to the data as explained in the text)(left).Spectrum of the resonance mass in data events in the SBs (right) used for the estimation of the V+jets distribution in the SR.The lower panels depict the pulls in each bin, (N data − N bkg )/σ, where σ is the statistical uncertainty in data, as given by the Garwood interval [63].
Table 2: Predicted number of background events and the observed number in the signal region, for all event categories.The regions denoted by W, Z and H are intervals in the jet soft-drop mass distribution that range from 65 to 85 GeV, from 85 to 105 GeV, and from 105 to 135 GeV, respectively.Separate sources of uncertainty in the expected number are reported as the statistical uncertainty in the V+jets contribution from the fitting procedure (fit), the difference between the nominal and alternative function form chosen for the fit (alt), and the uncertainty in the background from top quarks from the fit to the simulated jet mass spectrum.The principal background is V+jets, and its modeling represents the largest uncertainty in the analysis.The systematic uncertainty in the V+jets background is dominated by the statistical uncertainty associated with the number of events in the jet mass distribution SBs in data and simulation.An additional uncertainty is related to the choice of model used for the jet mass in the V+jets background.It is evaluated from the differences in the expected background yields obtained when using the alternative fitting functions.For the top quark processes, uncertainties from normalization and shape in the parametrization are propagated to the final background estimation.The single top quark and top quark pair production normalization uncertainty arises predominantly from the limited number of events in the CRs.

Category
The uncertainties in the shape of the V+jets distribution are estimated from the covariance matrix of the fit to m X in the sidebands and the uncertainties in the α(m X ) ratio, which depend on the number of events in data and simulation, respectively.
The uncertainties in the trigger efficiency, and in the electron and muon reconstruction, identification, and isolation efficiencies are obtained by varying the corresponding scale factors by their uncertainty, and each is found to be 1-2% [58,59].For the τ h reconstruction and identification, the uncertainties vary between 6 and 8% and between 10 and 13%, depending on the resonance mass, in the τ h and the τ h τ h channels, respectively [56].A separate uncertainty due to the extrapolation of the reconstruction and identification of τ h leptons at large p T has an impact on the signal normalization of 18% in the τ h , and 30% in the τ h τ h channels, for a 4 TeV signal hypothesis.This uncertainty is responsible for an increase of 1% in the width of the signal distribution.
Jet energy scale and resolution uncertainties affect both the selection efficiencies and the shape of distributions.The corrections to the jet mass scale and resolution are also taken into account, and result in a variation of 1-8% in the expected number of signal events.The jet energy scale uncertainty accounts for a variation in signal efficiency of 1-3%, while the variation in the jet energy resolution has an impact of 1-2%.The effect on the mass distribution is at the level of 1-2% for the mean and the width of the signal distribution.Event migrations between the mass windows due to the effect of jet mass scale and resolution variations are estimated to be between 2 and 15%, depending on the signal and the vector boson mass region.
Scale factors for V tagging and b tagging represent the largest source of normalization uncertainty for the signal.Uncertainties in normalization correspond to 6 and 11% in the HP and LP categories, respectively.An additional uncertainty from the extrapolation of the W tagging from the tt scale to larger values of jet p T is estimated using an alternative HERWIG [64] shower model, and varies from 2 to 18% for 0.9-4 TeV mass hypotheses and the two V tag categories.
In addition, the contribution to the signal normalization uncertainty from the b tagging uncertainty varies between 3 (4)% to 7 (5)% for the 2 (1) b-tagged subjet categories.
The τ lepton energy-scale uncertainties affect both the selection efficiencies and their distribution shapes.The effect on the signal yield amounts to 1% in the τ h channel.In the τ h τ h channel, the effect decreases from 5% for a resonance with mass of 0.9 TeV to 3% for a mass of 4.0 TeV.
Normalization uncertainties from the choice of the parton distribution function (PDF) grow larger with higher resonance mass, and are larger for gluon-initiated processes than for quarkinitiated processes.For W and Z production, which are sensitive to quark PDFs, effects range from 6 to 37%, while radion and graviton production depend on gluon PDFs, and result in a variation of 10 to 64% in the number of expected signal events.Uncertainties of similar magnitude arise from factor of two variations in the factorization and renormalization scales, resulting in 3 to 13% variations for W and Z , and 10 to 19% for radion and graviton production.While normalization uncertainties are not considered in setting limits on production, effects on signal acceptance are propagated to the final fit, amounting to 0.5-2% for PDF uncertainties, depending on resonance mass.
Other systematic uncertainties affecting the normalization of signal and minor backgrounds considered in the analysis include pileup contributions (0.5%) and integrated luminosity (2.5%) [65].A list of the main systematic uncertainties is given in Table 3.
Table 3: Summary of systematic uncertainties for the background and signal events.Uncertainties marked with "shape" are propagated also to the shape of the distributions, and those marked with † are not included in the limit bands, but instead reported in the theory band.The dash symbol is reported where the uncertainty is not applicable to a certain signal or background.The symbols qq' and gg refer to quark-initiated and gluon-initiated processes, respectively.Events / 100 GeV

Results
Results are obtained from a combined fit of the signal and background to the resonance mass distribution in data, based on a profile likelihood, where the systematic uncertainties are considered as "nuisance" parameters [66,67].The background-only hypothesis is tested against the signal hypothesis simultaneously in the different categories.No evidence of significant deviations from the background expectation are found.Assuming that the signals have widths that are negligible relative to the resonance-mass resolution of approximately 7%, the 95% confidence level (CL) upper limits are determined for the signals using the asymptotic frequentist method [66,68,69].Limits are obtained on the product of the cross section and branching fraction for a heavy resonance (X) that decays to HH, WH, or ZH, as reported in Figs.5-6.
Resonance spins of 0 and 2 are considered for the HH final state, while the resonance spin is assumed to be 1 for the WH and ZH final states.For the WH and ZH final states, the W and Z boson mass regions are combined because there are contributions from both signals to the two mass regions.The exclusion limit ranges from 80 to 5 fb for resonances of spin 0 and 2, and from 180 to 5 fb for spin-1 resonances.The τ h and τ h τ h final states, and the one and two b-tagged sub-jet categories are combined, to obtain the limits.The solid red lines and the red dashed areas correspond to the cross sections predicted by the bulk radion and graviton and their corresponding uncertainties, as reported in Table 3.
The predictions from bulk radion and graviton models are superimposed on the exclusion limits in Fig. 5, assuming Λ R = 1 TeV and k = 0.5.With this assumption for the theory parameters, a radion resonance with mass below 2.7 TeV is excluded at 95% CL.For a spin-1 signal, the results are interpreted in the context of the simplified HVT benchmark models A and B. As shown in Fig. 6, a W (Z ) resonance of mass lower than 2.6 (1.8) TeV is excluded at 95% CL in the HVT benchmark model B. The HVT benchmark model A is also reported for completeness.In the mass-degenerate spin-1 triplet hypothesis, the expected and observed limits on the V resonance are shown in Fig. 7 (left).
The exclusion limit shown in Fig. 7 (left) can be interpreted as a limit in the space of the HVT model parameters [g V c H , g 2 c F /g V ].Combining all channels, the excluded region in such a parameter space for narrow resonances is shown in Fig. 7 (right).The region of parameter space where the natural resonance width is larger than the typical experimental resolution of 7%, for which the narrow width assumption is not valid, is shaded.Expected limits are shown with ±1 and ±2 standard deviation uncertainty bands.The τ h and τ h τ h final states, for the HP and LP τ 21 categories, and the W and Z boson mass signal regions are combined, to obtain the limits.The solid lines and the relative dashed areas in magenta and red correspond to the cross sections predicted by the HVT models A and B, respectively, and their corresponding uncertainties, as reported in Table 3. Expected and observed 95% CL upper limit on σB(V → VH) with ±1 and ±2 standard deviation uncertainty bands (left) in the τ h and τ h τ h , τ 21 HP and LP categories, with the W and Z boson mass signal regions combined.Observed exclusion limit (right) in the space of the HVT model parameters [g V c H , g 2 c F /g V ], described in the text, for three different mass hypotheses of 1.5, 2.0, and 3.0 TeV.The region of parameter space where the natural resonance width is larger than the typical experimental resolution of 7%, for which the narrow width assumption is not valid, is shaded in grey.

Summary
A search has been conducted for heavy resonances that decay to two bosons, one of which is a W, Z, or Higgs boson that decays to a pair of quarks, and the other is a Higgs boson that decays to a pair of τ leptons.The analyzed data are collected by the CMS experiment at √ s = 13 TeV, corresponding to an integrated luminosity of 35.9 fb −1 .Reconstruction techniques have been developed to select events in which the τ lepton pair is highly boosted.The data are consistent with the standard model expectations and upper limits at 95% confidence level are set on the product of cross section and branching fraction for resonance masses between 0.9 and 4.0 TeV.This search yields the first results at √ s = 13 TeV in the considered mass range and final states.Assuming the ultraviolet cutoff of the theory Λ R = 1, Kaluza-Klein excitations of spin-0 radions with mass smaller than 2.7 TeV are excluded at 95% confidence level.In the heavy vector triplet model B context, a mass-degenerate vector triplet V resonance with mass below 2.8 TeV is excluded at 95% confidence level.

Figure 3 :
Figure 3: Data and expected backgrounds in the τ h channel.The W mass window is shown in the HP (upper left) and LP (upper right) categories, the Z mass window for the HP (middle left) and LP (middle right) categories, and the H mass window for the two b-tagged subjet (lower left) and one b-tagged subjet (lower right) categories.The lower panels depict the pulls in each bin, (N data − N bkg )/σ, where σ is the statistical uncertainty in data, as given by the Garwood interval [63], and provide estimates of the goodness of fit.Signal contributions are shown, assuming benchmark HVT model B for the V and Λ R = 1 for the radion.

Figure 4 :
Figure 4: Data and expected backgrounds in the τ h τ h channel.The W mass window is shown in the HP (upper left) and LP (upper right) categories, the Z mass window for the HP (middle left) and LP (middle right) categories, and the H mass window for the two b-tagged subjet (lower left) and one b-tagged subjet (lower right) categories.The lower panels depict the pulls in each bin, (N data − N bkg )/σ, where σ is the statistical uncertainty in data, as given by the Garwood interval [63], and provide estimates of the goodness of fit.Signal contributions are shown, assuming benchmark HVT model B for the V and Λ R = 1 for the radion.

Figure 5 :
Figure5: Observed 95% CL upper limits on σB(X(spin-0) → HH) (left) and σB(X(spin-2) → HH) (right).Expected limits are shown with ±1 and ±2 standard deviation uncertainty bands.The τ h and τ h τ h final states, and the one and two b-tagged sub-jet categories are combined, to obtain the limits.The solid red lines and the red dashed areas correspond to the cross sections predicted by the bulk radion and graviton and their corresponding uncertainties, as reported in Table3.

Figure 6 :
Figure 6: Observed 95% CL upper limits on σB(W → WH) (left) and σB(Z → ZH) (right).Expected limits are shown with ±1 and ±2 standard deviation uncertainty bands.The τ h and τ h τ h final states, for the HP and LP τ 21 categories, and the W and Z boson mass signal regions are combined, to obtain the limits.The solid lines and the relative dashed areas in magenta and red correspond to the cross sections predicted by the HVT models A and B, respectively, and their corresponding uncertainties, as reported in Table3.

Figure 7 :
Figure 7:  Expected and observed 95% CL upper limit on σB(V → VH) with ±1 and ±2 standard deviation uncertainty bands (left) in the τ h and τ h τ h , τ 21 HP and LP categories, with the W and Z boson mass signal regions combined.Observed exclusion limit (right) in the space of the HVT model parameters [g V c H , g 2 c F /g V ], described in the text, for three different mass hypotheses of 1.5, 2.0, and 3.0 TeV.The region of parameter space where the natural resonance width is larger than the typical experimental resolution of 7%, for which the narrow width assumption is not valid, is shaded in grey.

Table 1 :
Normalization scale factors for top quark production for different event categories, depending on the V tagging and H tagging requirement applied.Uncertainties are due to the limited number of events in the control regions and the uncertainty in the b tagging efficiency.