Probing heavy Majorana neutrinos and the Weinberg operator through vector boson fusion processes in proton-proton collisions at $\sqrt{s}$ = 13 TeV

The first search exploiting the vector boson fusion process to probe heavy Majorana neutrinos and the Weinberg operator at the LHC is presented. The search is performed in the same-sign dimuon final state using a proton-proton collision data set recorded at $\sqrt{s}$ = 13 TeV, collected with the CMS detector and corresponding to a total integrated luminosity of 138 fb$^{-1}$. The results are found to agree with the predictions of the standard model. For heavy Majorana neutrinos, constraints on the squared mixing element between the muon and the heavy neutrino are derived in the heavy neutrino mass range 50 GeV-25 TeV; for masses above 650 GeV these are the most stringent constraints from searches at the LHC to date. A first test of the Weinberg operator at colliders provides an observed upper limit at 95% confidence level on the effective $\mu\mu$ Majorana neutrino mass of 10.8 GeV.


1
Despite the enormous success of the standard model (SM) of particle physics over the last fifty years, the observation of neutrino oscillations [1,2] has confirmed that at least two SM neutrinos have small, but nonzero, masses. This provides a compelling hint of physics beyond the SM. Efforts to determine each neutrino mass are being pursued by many experiments within different research fields, including cosmology and tritium beta decay [3][4][5], and various theoretical explanations of these observations have been proposed [6][7][8][9][10][11][12][13].
A natural formalism for generating neutrino masses is through a dimension-five operator, proposed by Weinberg [10], which extends the SM Lagrangian with where , are the flavors of the leptons, which can be electrons, muons, or taus; Λ is the scale at which the particles responsible for neutrino masses become relevant degrees of freedom; C 5 is a flavor-dependent Wilson coefficient; L = (ν , ) is the left-handed lepton doublet; and Φ is the SM Higgs doublet with a vacuum expectation value v = √ 2 Φ ≈ 246 GeV. This operator, known as the Weinberg operator, generates the Majorana neutrino masses as and introduces lepton number violation (LNV), as the Majorana neutrino is its own antiparticle. The LNV processes have been tested extensively in searches for neutrinoless double beta decays of heavy nuclei. These experiments search for events with two same-sign (SS) charge electrons and the absence of neutrinos in the final state, and set stringent upper limits on the associated Majorana neutrino mass (|m ee |) at 90% confidence level (CL) in the ranges 61-165 meV [14] and 79-180 meV [15]. Since muons and tau leptons are much heavier than electrons, SS final states involving these heavier leptons are kinematically forbidden in low-energy nuclear experiments and can only be studied at colliders [16]. As the study of tau leptons is experimentally more challenging, the final state with SS muons is investigated in this search. The observed upper limit on the effective µµ Majorana neutrino mass (|m µµ |), translated from the searches for LNV decays of charged kaons at the NA62 experiment [17], is around 55 GeV at 90% CL [16].
The Weinberg operator can be realized in the context of "seesaw" models [6][7][8][9], assuming the existence of hypothetical heavy states, for example a heavy Majorana neutrino (N) in the type-I seesaw model. The small mass of the SM neutrinos can thus be explained by a suppression due to the high mass of new particles, because m ν ∼ y 2 ν v 2 /m N . In this expression, m ν is the SM neutrino mass, m N is the mass of the heavy Majorana neutrino, and y ν is the Yukawa coupling. This ability to explain the small mass of the SM neutrinos motivates the search for heavy Majorana neutrinos. Heavy Majorana neutrinos can only couple to the SM through mixing with SM neutrinos, which is characterized by the mixing element, V N , between a SM neutrino in its lefthanded interaction state and a heavy Majorana neutrino in its mass eigenstate. The searches that have set limits on V N were based on the production of a heavy Majorana neutrino through either quark-antiquark annihilation or Wγ, with different analysis strategies based on various decay modes, such as many leptons in the final state [18][19][20], and in final states containing LNV lepton pairs [21,22]. For m N at a few GeV, the constraints on |V N | 2 , where = {e, µ}, are around 10 −7 -10 −5 , while at around 1 TeV, the constraints approach unity. This Letter reports the first search probing Majorana neutrinos and the Weinberg operator in vector boson fusion (VBF) processes [16,23] at the CERN LHC, studying events with SS muon pairs. As shown in Fig. 1 by the two representative Feynman diagrams, in the t channel, SS W boson pairs may convert to SS lepton pairs via a TeV-scale Majorana neutrino or through a Weinberg operator process. Because the cross section of t-channel processes is less sensitive to the mass of the intermediate particle compared with s-channel quark-antiquark annihilation processes, the VBF process can complement searches for heavy Majorana neutrinos at the TeV mass scale as its cross section decreases more slowly with increasing N mass compared with the values from s-channel production. Although the various seesaw mechanisms can be considered as realizations of the Weinberg operator, an alternative implementation is considered in this Letter, as shown in Fig. 1 (right). The two VBF processes considered in this work differ primarily in that the Majorana process is assumed to be mediated by a heavy t-channel neutrino whereas the Weinberg operator process is mediated by a lighter t-channel neutrino. The CMS apparatus [25] is a multipurpose, nearly hermetic detector, designed to trigger on [26,27] and identify electrons, muons, photons, and hadrons [28-30]. A global "particle-flow" algorithm [31] aims to reconstruct all individual particles in an event, combining information provided by the all-silicon inner tracker and by the crystal electromagnetic and brass-scintillator hadron calorimeters, operating inside a 3.8 T superconducting solenoid, with data from gasionization muon detectors embedded in the flux-return yoke outside the solenoid. The reconstructed particles are used to build charged leptons, jets, and missing transverse momentum Samples of signal events are simulated with next-to-leading order (NLO) precision using the MADGRAPH5 aMC@NLO 2.6.5 generator [35], based on model implementations from Refs. [16,23,36,37]. For the heavy Majorana neutrino, we assume only one heavy state and scan its mass in a range from 50 GeV to 25 TeV. The process induced by the Weinberg operator is simulated using a novel approach, established in Ref. [16], which approximates the internal neutrino lines and vertex insertion with an effective light Majorana neutrino. The Wilson coefficient C 5 is set to 1, and the energy scale Λ to 200 TeV. The production cross sections scale with V 2 N V 2 N and |C 5 /Λ| 2 for the VBF processes with heavy Majorana neutrinos [23] and the Weinberg operator [16], respectively. Since the kinematic shape of the signal is independent of |V µ N |, we set |V µ N | = 1 for the signal simulation to ensure efficient generation of signal events.
The SM electroweak W ± W ± and W ± Z processes, where both bosons decay leptonically, are simulated using MADGRAPH5 aMC@NLO 2.4.2 [35] at leading order (LO) with six electroweak O(α 6 ) and zero quantum chromodynamics (QCD) vertices. The QCD-induced W ± W ± process is also simulated using MADGRAPH5 aMC@NLO 2.4.2 at LO. Contributions with an initialstate b quark are excluded from the electroweak W ± Z simulation because they are included in the production of a Z boson in association with a single top quark, known as the tZq process, which is simulated at NLO using MADGRAPH5 aMC@NLO 2.3.3. The ZZ → 4 and gg → ZZ processes are simulated with the POWHEG v2 [38] and MCFM 7.0.1 [39] generators, respectively. The production of ttW ± , ttZ, and triple vector boson (VVV) background events is simulated at NLO in QCD using the MADGRAPH5 aMC@NLO 2.2.2 (2.4.2) generator [35, 40,41] for the samples corresponding to the 2016 (2017-2018) data-taking period. The production of W ± W ± events through double-parton scattering is generated at LO using PYTHIA 8.226 (8.230) [42] for 2016 (2017-2018).
The simulated samples for the 2016 (2017-2018) data-taking period use the NNPDF 3.0 (3.1) parton distribution functions [43,44], and are interfaced with PYTHIA to model the fragmentation and hadronization of partons in the initial and final states, along with the underlying event. The CUETP8M1 (CP5) tune [45,46] is used in the simulation for the 2016 (2017-2018) data-taking period. In LO (NLO) simulations performed with MADGRAPH5 aMC@NLO, jets are matched to the parton shower produced by PYTHIA 8.2 following the MLM [47] (FxFx [41]) prescription. The interactions of all final-state particles with the CMS detector are simulated using GEANT4 [48]. Simulated events are mixed with the contribution of particles from additional pp interactions within the same or nearby bunch crossings (pileup); the pileup distribution in the simulated samples matches that observed in data.
Collision events are collected using single-muon triggers [27] that require the presence of an isolated muon with transverse momentum p T > 24 or 27 GeV, depending on the data-taking period. In addition, a set of dimuon triggers [27] with lower p T thresholds is used to ensure a trigger efficiency above 90% for events with muons. Depending on the subsequent off-line selection, the dimuon triggers increase the event acceptance by up to 4%. The primary vertex (PV) is taken to be the vertex corresponding to the hardest scattering in the event, evaluated using tracking information alone, as described in Section 9.4.1 of Ref. [49]. The jets are clustered using the anti-k T algorithm [50,51] with a distance parameter of 0.4, with the tracks assigned to the candidate vertices as inputs. To mitigate the impact from pileup, tracks identified to be originating from pileup vertices are discarded and calorimeter energies are corrected to account for charged and neutral particles originating from those vertices. The p miss T is then defined as the negative vector p T sum of all reconstructed particles in an event. Its magnitude is referred to as p miss T . Jets arising from b quark hadronization and decay (b jets) are identified using a deep neural network algorithm, termed DEEPCSV, whose inputs are tracks displaced from the primary interaction vertex, details of secondary vertices, and the kinematical variables of jets [52]. To classify jets as b tagged, we use a requirement on the score of the DEEPCSV algorithm for which the b jet identification efficiency ranges from 65 to 75% and the rate at which gluon or light-flavor jets are misidentified as b jets is around 1%. The fraction of jets originating from the hadronization of c quarks that are incorrectly identified as b jets is about 10%. Starting from the constituents of reconstructed jets, hadronically decaying τ leptons are identified via the "hadrons-plus-strips" algorithm [32].
The final states analyzed contain two well-identified SS muons with p T > 30 GeV and pseudorapidity |η| < 2.4, and at least two jets with p T > 30 GeV and |η| < 4.7. The two selected jets with the largest p T are considered as the VBF jet candidates, which must satisfy a large pseudorapidity separation |∆η jj | > 2.5 and a large dijet invariant mass m jj > 750 GeV. Events produced via VBF are expected to contain rapidity gaps between the muons and the VBF jets.
To utilize this signature, we require max(Z ) < 0.75, where Z is the Zeppenfeld variable [53] of a lepton : Z = |η − (η j1 + η j2 )/2|/|∆η jj |, where η is the pseudorapidity of the lepton, and η j1 and η j2 are of the two VBF jets. Events with an additional loosely identified muon (electron) with p T > 10 GeV and |η| < 2.5(2.4) or an identified hadronically decaying tau lepton with p T > 20 GeV and |η| < 2.3 are rejected to suppress backgrounds such as W ± Z production. Because the signal processes do not produce final state neutrinos other than those produced in hadron decays, p miss T is a useful discriminant to separate the signal and backgrounds. For the analysis targeting processes with heavy Majorana neutrinos, the azimuthal separation variable ∆φ is better for improving the sensitivity to the signal than p miss T because the two muons in the final state tend to be more back-to-back than in background processes [23]. However, p miss T is more discriminating than azimuthal separation for the Weinberg operator process because of the different event topology resulting from the discrepancy between the two t-channel propagators. To separate the signal and the SM electroweak W ± W ± events, signal regions (SRs) for the heavy Majorana neutrino and Weinberg operator analyses are defined in bins of ∆φ and p miss T , respectively. For the heavy Majorana neutrino analysis, we define a high-∆φ bin (∆φ > 2) and a low-∆φ bin (∆φ < 2) in the SR. For the Weinberg operator analysis, a low-p miss T bin (< 50 GeV) and a high-p miss T bin (> 50 GeV) are defined in the SR.
The background from opposite-sign dilepton events where one of the lepton charges is misidentified is negligible for events with muon pairs [54, 55]. The so-called nonprompt-lepton backgrounds, originating from leptonic decays of heavy quarks or hadrons misidentified as leptons, are suppressed by identification and isolation requirements. Events from a tt pair where only one W boson decays leptonically are the main source of nonprompt-lepton backgrounds. To reduce this contribution, events with at least one b-tagged jet with p T > 20 GeV and |η| < 2.4 are rejected. The remaining contribution of the nonprompt-lepton background is estimated directly from a data sample by applying weights to events containing muon candidates that fail the nominal selection criteria while passing a less stringent isolation requirement. These weights, called fake-lepton factors, are obtained from the probability of a jet being misidentified as a lepton and the efficiency of correctly reconstructing and identifying a lepton. More details about this method are given in Refs. [56,57].
To constrain the reducible background contributions, which include W ± Z production and backgrounds from nonprompt leptons, several control regions (CRs) in data are used. These CRs are defined to select samples of events enriched in those background processes, following the procedure described in Refs. [58][59][60]. A WZ CR, requiring the presence of three muons in an event, is used to estimate W ± Z background contributions. The opposite-sign dimuon combination with invariant mass closest to the Z boson mass (m Z ) is considered as the Z boson decay product, whereas the remaining muon is assumed to originate from the decay of the W boson. The selection criteria defining the SR and the WZ CR are summarized in Table 1. To select event samples enriched in nonprompt leptons, a b-tagged CR is defined by using the same selection as for the SR, however requiring the presence of a b-tagged jet. Similarly, the WZb CR is defined by requiring the same selection as for the WZ CR, but requiring at least one b-tagged jet. The dominant backgrounds in the SR are SM electroweak W ± W ± production and the contribution from nonprompt leptons.
The discriminating variable used to search for a signal for both analyses is H T /p µ 1 T . Here, H T is the scalar sum of jet p T for all jets with p T > 30 GeV and |η| < 4.7, µ 1 indicates the highest p T muon. This discriminating variable measures the amount of hadronic activity relative to leptonic activity in a given event, and is sensitive to the color structure of the hard scattering processes as discussed in Ref. [16]. As hadronic activity is typically suppressed in VBF production processes, H T /p µ 1 T tends to be smaller for signal events than for background events [16,23].
The statistical analysis is performed with a profile likelihood ratio test statistic in which systematic uncertainties are modeled as constrained nuisance parameters. Limits are set using the , and the overall impact on the simulated samples is less than 1%. Theoretical uncertainties arising from the choice of the renormalization and factorization scales are estimated by varying these scales independently up and down by a factor of 2 from their nominal values, excluding the two extreme variations [70]. The largest cross section variations from this procedure are taken as the uncertainty [71], which is about 6%, and which is applied to both the signal and the SM electroweak W ± W ± processes. The parton distribution function uncertainty of the signal and the SM electroweak W ± W ± processes is evaluated according to the procedure described in Ref. [72], and is found to be about 1%. These uncertainties are not applied to the subdominant W ± Z background, where they are not significant compared to other uncertainties considered. For the nonprompt-lepton background prediction, systematic uncertainties arise from the limited size of the sample used to measure the efficiency and from the difference in the flavor composition of the jets that are misidentified as leptons between the measured sample and the signal region, estimated to amount to at most 5 and 7%, respectively. An additional normalization uncertainty of 30% is assigned to the nonprompt-lepton background, based on studies of simulated events [57]. The statistical uncertainty in the signal and background templates is introduced as a systematic uncertainty via the Barlow-Beeston-lite approach [73]. This is the systematic uncertainty that has the largest impact on the results derived from the simultaneous fits described below. Uncertainties related to pileup, theoretical calculations, and the muon trigger and identification are correlated across the three data-taking years. The luminosity and b tagging uncertainties are partially correlated through some of the individual uncertainty sources, and the remaining uncertainties are uncorrelated across the three datasets. The total uncertainties on the results are dominated by the statistical uncertainties.
Two separate fits are performed, based on the SR and the b-tagged, WZ, and WZb CRs: one for the heavy Majorana neutrino analysis using ∆φ bins in the SR, and a second for the Weinberg operator analysis with the p miss T bins in the SR. The normalization factors for the W ± W ± , W ± Z, and tZq background processes, affecting both the SRs and CRs, are included as free parameters in the fit together with the signal strength. The predicted yields are shown in Fig. 2 with their best fit normalizations from the simultaneous fits for the background-only hypothesis, i.e. assuming no contributions from processes involving heavy Majorana neutrinos and the Weinberg operator. The bin boundaries are chosen to optimize the signal sensitivity while leaving the major background contributions in each bin at a reasonable level. The data are found to be consistent with the background expectation. Using the relationship between the cross section and the mixing elements for the heavy Majorana neutrino analysis, upper limits at the 95% CL are derived on |V µ N | 2 as shown in Fig. 3. These results surpass those obtained in previous searches by the ATLAS and CMS Collaborations [18,19,21] for m N 650 GeV, and set the first direct limits for m N > 2 TeV. Previous upper limits on the mixing elements from the CMS Collaboration [18,21] are plotted in Fig. 3 for comparison. These limits on |V µ N | 2 approach unity at around 1 TeV. As the mixing element describes the oscillation probability between the SM neutrinos and the heavy Majorana neutrino, results with |V µ N | 2 > 1 are not physically meaningful.
Equation (2) is used to convert the limit on |C 5 /Λ| 2 extracted from the Weinberg operator analysis into a limit on |m µµ |. The observed (expected) 95% CL upper limit on |m µµ | is found to be 10.8 (12.8) GeV. This is the first upper limit from a collider experiment, and it improves upon a previous limit set in Ref. [16] using the results from the NA62 experiment [17]. In summary, this Letter presents the first search for Majorana neutrinos at several TeV and a first probe of the Weinberg operator at the LHC. These achievements were made possible by considering for the first time vector boson fusion processes resulting in a same-sign dimuon final state. The results are consistent with the predictions from the standard model. For heavy Majorana neutrinos, upper limits on the mixing element |V µ N | 2 are set for the mass range 50 GeV < m N < 25 TeV and the best sensitivity is reached for m N 650 GeV. The phase space explored exceeds the center-of-mass energy of the LHC, improving previous limits for direct production of Majorana neutrinos. The highest mass for which |V µ N | 2 = 1 is excluded is around 23 TeV. The observed (expected) 95% confidence level upper limit on the effective µµ Majorana mass associated with the Weinberg operator is 10.8 (12.8) GeV, exceeding the current best limit from high intensity Kaon experiments.

Acknowledgments
We congratulate our colleagues in the CERN accelerator departments for the excellent performance of the LHC and thank the technical and administrative staffs at CERN and at other CMS institutes for their contributions to the success of the CMS effort. In addition, we gratefully acknowledge the computing centers and personnel of the Worldwide LHC Computing Grid and other centers for delivering so effectively the computing infrastructure essential to our analyses. Finally, we acknowledge the enduring support for the construction and operation of the LHC, the CMS detector, and the supporting computing infrastructure provided by the following funding agencies:  [27] CMS Collaboration, "The CMS trigger system", JINST 12 (2017) P01020, doi:10.1088/1748-0221/12/01/P01020, arXiv:1609.02366.