Exotic diboson $Z'$ decays in the U$\mu\nu$SSM

Searches for new leptophobic resonances at high energy colliders usually target their decay modes into pairs of light quarks, top quarks, or standard model bosons. Additional decay modes may also be present, producing signatures to which current searches are not sensitive. We investigate the performance of generic searches that look for resonances decaying into two large-radius jets. As benchmark for our analysis we use a supersymmetric $\text{U}(1)'$ extension of the Standard Model, the so-called U$\mu\nu$SSM, where all the SM decay modes of the $Z'$ boson take place, plus additional (cascade) decays into new scalars. The generic searches use a generic multi-pronged jet tagger and take advantage of the presence of $b$ quarks in the large-radius jets, and are sensitive to all these $Z'$ decay modes (except into light quarks) at once. For couplings that are well below current experimental constraints, these generic searches are sensitive at the $3\sigma-4\sigma$ level with Run 2 LHC data.


Introduction
Realistic extensions of the Standard Model (SM) involve several new particles in addition to the already discovered ones, which likely appear at different mass scales. Their collider signals, arguably, are often not captured by the simplified scenario framework in which most searches at the Large Hadron Collider (LHC) are performed. The variety of possible signals makes dedicated searches unmanageable, therefore, a generic approach is compulsory.
The simplest scenario where generic searches become very useful is the decay of a very heavy resonance R (at the multi-TeV scale) into two particles at the electroweak scale, which can in turn decay into quarks, or into lighter new particles which ultimately decay into quarks. The signature of such cascade decay is a pair of massive jets, each of them corresponding to one of the decay products of R. Among the possible signatures with this topology, one can have spin, and decay products that can be either scalars or vector bosons, because the generic tagger is so designed.
As framework to study these signals we consider the recently proposed UµνSSM [7], 1 which is a U (1) extension of the 'µ from ν' supersymmetric standard model (µνSSM) [10] (for a recent review, see Ref. [11]). In such a framework the three decay modes (diboson, and exotic di/multiboson) are present. This approach is the opposite to the one adopted in Ref. [6]. There, a minimal model was written down in which the usual diboson decays Z → W W are absent, so that generic searches sensitive to exotic dibosons and multibosons are compulsory. Here, we focus on a model in which all the three types of bosonic decays are present, to assess the sensitivity of a generic search in that situation and see, with a concrete example, how the different Z bosonic decay modes contribute to the combined sensitivity. Other models with a Z boson may have one or more of the decays described above.
In UµνSSM scenarios, the presence of R-parity violating couplings involving righthanded (RH) neutrino superfields,ν c , solves simultaneously the µ problem [12] (for a recent review, see Ref. [13]) of the minimal supersymmetric standard model (MSSM) (for reviews, see e.g. Refs. [14][15][16]) and the ν problem being able to reproduce the neutrino data [7,8,10,[17][18][19][20][21]. In the superpotential of these scenarios, in addition to Yukawa couplings for neutrinos Y νĤ uLν c , the couplings λν cĤ dĤu are allowed generating an effective µ-term when the RH sneutrinos develop electroweak-scale vacuum expectation values (VEVs), The UµνSSM models built in Ref. [7] are constrained by anomaly cancellation conditions and have the attractive properties of forbidding baryon-number-violating operators as well as explicit Majorana masses and µ terms, and avoiding potential domain wall problems. Besides, the extra U (1) gauge symmetry provides the RH neutrinos with a non-vanishing charge, avoiding the uneasy situation from the theoretical viewpoint of being the only fields with no quantum numbers under the gauge group. The presence of several right sneutrinos as well as of the two Higgs doublets H d and H u , provides a large spectrum of scalar particles which could give rise to the aforementioned multiboson signals.
In particular, assuming three families of RH neutrino superfieldsν c i , i = 1, 2, 3, there are eight (seven) neutral scalar (pseudoscalar) states from the mixing between Higgses and sneutrinos. However, the three left sneutrinos are almost decoupled, and upon diagonalisation we are left in the scalar sector with the SM-like Higgs, the heavy doublet-like neutral Higgs and three CP-even singlet-like states, and in the pseudoscalar sector with the heavy doublet-like pseudoscalar and three CP-odd singlet-like states. Therefore, decay channels into dijets from a pair production at the LHC of CP-even and CP-odd right sneutrinos via the Z , are viable in the UµνSSM. In this work, we will carry out a detailed analysis of this signal.
The paper is organised as follows. Sec. 2 will be devoted to the discussion of the UµνSSM benchmarks we use to illustrate the production of exotic multiboson signals. In Sec. 3 we examine the sensitivity to such signals with generic searches. Finally, our conclusions and a discussion of possible generalisations are left for Sec. 4.

UµνSSM benchmarks
In this work we consider scenario 1 of Ref. [7] in which the Z boson is leptophobic. This is the case of interest for resonance searches in hadronic final states: otherwise, the Z decay into lepton pairs gives much cleaner and distinctive signals. There are other possibilities for leptophobic Z bosons in UµνSSM models, and the results are quantitatively quite similar. The U(1) charges of the relevant fields are collected in Table 1.
The superpotential of the UµνSSM and interactions of the Z boson are summarised in Ref. [7], and we refer the reader to that work for details. We restrict ourselves here to the aspects of the model most directly related to our present analysis.
The total Z production cross section for the scenario considered is presented in Fig. 1 as a function of its mass. The Z partial width into fermion pairs is with g Z the U (1) gauge coupling, N c the number of colours and m f the fermion mass. Notice that the Z boson predominantly decays into up-type quarks, given the U (1) charges in Table 1, where z(f ) (−z(f c )) are the charges for the left (right) chiral fermions. The Z partial width into SM boson pairs is where tan β ≡ v u /v d is the usual ratio of neutral scalar VEVs, e.g. in the MSSM, and we have assumed the alignment limit in the latter channel. We have defined the usual kinematical function Notice that for the leptophobic scenario considered z(H d ) = 0. These two partial widths are nearly equal for m Z m W,Z,h .
In addition, we have decays into scalar pairs. The interaction of the new gauge boson with neutral scalar fields, generically denoted as φ, is with z(φ) the corresponding U (1) charge. If one writes φ = (φ R +iφ I )/ √ 2, the interaction reads connecting a CP-even weak eigenstate φ R and a CP-odd one φ I . Assuming CP conservation in the scalar mixing, the mass eigenstates are either CP-even (denoted generically by H) or CP-odd (denoted by A) and the Z decays into pairs H i A j . Considering that additional singlet scalars under the SM gauge group can be present in UµνSSM models [7], we will take the simplifying assumption that they are heavier than the Z and do not have a significant mixing with H, A. Therefore, the partial width into pairs of scalar eigenstates s 1 , s 2 with the same U(1) charges is The presence of the new decay modes in (6) slightly decreases the branching ratio into SM final states, and relaxes the current limits with respect to those in Ref. [7]. In our numerical benchmarks we will assume these two additional channels: In the UµνSSM these dibosons are two CP-even and two CP-odd right sneutrinos, The current limits from the Z boson mass and couplings from dijet [22,23], tt [24,25], diboson [3] and Zh [4] resonance searches are presented in Fig. 2. We also show the cross section predictions for two g Z values. As in Ref. [7], we have assumed tan β = 2. In the following, we will use two benchmarks, • m Z = 2 TeV and g Z = 0.3, with m H 1 ,A 1 = 80 GeV, m H 2 ,A 2 = 200 GeV.
The Z mass and coupling in both cases are well below the current 95% confidence level (CL) upper limits. As mentioned, we have assumed CP conservation, so that the neutral scalar mass eigenstates have a definite parity. For simplicity, we also assume that the (H i , A i ) pairs, i = 1, 2, have the same mass, and that the rest of scalars present in the model are heavier. The Z is narrow in both scenarios, with Γ Z = 1.35 GeV and Γ Z = 5.6 GeV, respectively, and the intrinsic Z width is not relevant for the calculations.
The new scalars H 1 , A 1 are expected to decay predominantly into bb via a small mixing with the neutral scalar and pseudoscalar Higgses (denoted by H d and A d respectively), as depicted in Fig. 3 (left). The heavier CP-even scalar H 2 can decay into W W , ZZ and bb via small mixings with the neutral scalar Higgses H u and H d , but in relevant regions of the parameter space is expected to predominantly decay into H 1 H 1 and A 1 A 1 , as depicted  on the right panel. The partial widths for the different modes are Here g is the SU (2) [7]. The values of |Z H H 2 H d | 2 and |Z H H 2 Hu | 2 are expected to be small, in the range 10 −2 − 10 −4 (see e.g. Ref. [26]) and we can take them equal for simplicity in our discussion. As we can see in Fig. 4 in the upper end of the expected range, H 2 predominantly decays into two scalars for κ ≥ 0.1, and has branching ratio to H 1 H 1 + A 1 A 1 close to unity for κ ≥ 0.25. In the computation of the next section we will take this branching ratio equal to unity. Note that the decays H 2 → H 1 h are suppressed by the square of the small mixing factors Z H . On the other hand, H 2 → hh is not suppressed, and is possible if kinematically open. When h decays into bb (which is the dominant decay mode of the 125 GeV Higgs boson) the decays H 2 → hh produce similar signals to H 2 → H 1 H 1 +A 1 A 1 , and for simplicity we consider benchmarks where H 2 → hh is kinematically forbidden. Higgs decays to final states other than bb also produce massive jets but a dedicated study is beyond the scope of this work.
The heavier CP-odd scalar A 2 is expected to decay into H 1 A 1 for analogous reasons. The partial widths for the possible modes are where |Z A A 2 A d | 2 is the small A d component of the SM singlet-like pseudoscalar A 2 . As it can be seen from Fig. 4 (right), the decay A 2 → H 1 A 1 is dominant for κ ≥ 0.05, and has branching ratio near unity for κ ≥ 0.1. Overall, in the benchmark points we have considered the decay of the Z boson into the scalars H 2 , A 2 produces the cascade decays depicted in Fig. 3 (right). Note, however, that longer cascades would in principle be possible. for example if we had considered a large mass differences between H i and A i (we have taken m H i = m A i ). In addition, cascade decays may involve the third pair H 3 , A 3 that we have assumed much heavier. In any case, these two simplified benchmarks depicted in Fig. 3 illustrate the two mechanisms of production of two multi-pronged large radius jets from the Z decay into new scalars.

Search strategies
The strategy to search for the dijet signals from multiple cascade decays of the Z boson into scalars takes advantage of the common features expected: (1) sizeable jet masses; (2) multi-pronged jet structure; (3) multiple b quarks that can be tagged.
Because we are not looking for a specific signal, a narrow jet mass window cannot be a priori imposed. Neither a dedicated tagger can be used; rather, a generic tagger for multipronged jets [27][28][29], or generic anomaly detection methods [30][31][32][33][34][35][36][37][38][39][40][41][42] must be used. In order to take advantage of the presence of b quarks, the sample is divided into different categories corresponding to the number of b tags, and the results are subsequently combined.
We perform a simulation including the backgrounds from dijet production (jj), bb, and tt. Event samples are generated with MadGraph [43] in 100 GeV bins starting at [300, 400] GeV up to [2.1, 2.2] TeV, plus a bin with p T ≥ 2.2 TeV. Each dijet sample contains 6 × 10 5 events, and each bb and tt sample 10 5 events. For each background process (jj, bb, tt) there are in total 20 samples in different p T bins which are combined with weight proportional to the cross section, so as to have good statistics across a wide range of transverse momentum. For the signal processes the relevant Lagrangian is implemented in Feynrules [44] and interfaced to MadGraph5 using the universal Feynrules output [45]. Hadronisation and parton showering performed with Pythia 8 [46] and detector simulation with Delphes 3.4 [47] using the configuration for the CMS detector. The reconstruction of jets and the analysis of their substructure is done using Fast-Jet [48].
In the analysis we use a collection of 'fat' jets with radius R = 0.8 reconstructed with the anti-k T algorithm [49] and groomed with Recursive Soft Drop [50]. In addition, for b tagging we use a collection of 'track jets' of radius R = 0.2, reconstructed using only tracks. This procedure is similar to the one applied for b-tagging of boosted H → bb by the ATLAS Collaboration [51]. A fat jet is considered to have a single b tag if there is one (and only one) b-tagged track jet within a distance ∆R = 0.8. A double b tag on the fat jet is considered if there are two or more b-tagged track jets within a distance ∆R = 0.8. The fat jet is untagged if there are no b-tagged tagged track jets within ∆R = 0.8. For the tagging of track jets we use an updated Run 2 parametric efficiency formula [52] in the 50% efficiency working point. Despite the fact that this parameterisation was not explicitly derived for R = 0.8 jets, we believe it still a sufficiently good approximation for our purpose, with misidentification probabilities around 0.14 for charm jets and 0.05 for light jets at the p T range of interest.
As event pre-selection we require the presence of two jets with mass m J ≥ 50 GeV, transverse momentum p T ≥ 400 GeV, rapidity difference |∆η| ≤ 1.3, and invariant mass m JJ ≥ 1 TeV. The presence of a signal can be detected as a bump in the dijet invariant mass distribution. We present in Fig. 5 this distribution for the relevant decays of the Z boson, in the benchmark with m Z = 2 TeV. Notice that the Z → tt contribution has a large tail at lower m JJ , due to the semilptonic top decays that are included in our simulation. For the rest of signals, the peak is slightly shifted due to the jet grooming. Although the different signal components do not peak exactly at the same mass, the combination of the several contributions is still possible. The background composition depends on m JJ , with tt and bb decreasing faster than jj because the latter includes processes with two valence quarks. At this level of event selection dijet production is overwhelmingly dominant. As a general rule, the higher the number of b tags, the smaller this background is because of the mistag suppression.
We use the MUST inspired tagger GenT in Ref. [29] that classifies any type of multipronged jets as signal, and QCD jets as background. The NN score X is shown in Fig. 6 for the QCD dijet background and the several signals arising from Z decays. (For each event with two jets, we use both jets in the plot.) We omit Z → jj which is backgroundlike and obviously suppressed by pre-selection cuts m J ≥ 50 GeV, and further by event selection. A dijet mass window m JJ ∈ [1.7, 2.3] TeV is required, and the 2 TeV Z benchmark is used. As it is expected, jets stemming from the five Z decay modes into massive particles are correctly classified as signal with a distribution leaning to the r.h.s. of the plot, while QCD jets are classified as background. Note that with this tagger the continuum tt production is not rejected as background, as it produces multi-pronged jets like the Z → tt signal ones.
For the final event selection we require, in addition, a lower cut on the NN score for the two jets that reduces the huge dijet background. The precise value of the cut is not optimised. Rather, we set a generic cut for each sample, such that the background is reduced to O(1) event per 100 GeV bin at m JJ 3 TeV. Specifically, we require • Sample 2DT: X ≥ 0.5; • Sample 1D1ST: X ≥ 0.6; • Samples 1DT, and 2ST: X ≥ 0.8; • Samples 1ST and 0T: X ≥ 0.9.
In the 2DT and 1D1ST samples the background is already small, therefore the cuts are rather mild. On the other hand, stringent cuts are required to reduce the background in the 1ST and 0T samples. The total background plus the Z signal is shown in Fig. 7 for the six event classes, before and after GenT tagging. The luminosity is taken as L = 139 fb −1 . Only in the 2DT class the presence of the signal is noticeable by eye, as a bump on the falling distribution after GenT tagging. For illustration, we collect in Table 2 the background composition in the region m JJ ∈ [1.7, 2.3] TeV in the six event classes.
Note that in some classes the background is dominated by tt, meaning that the b quark and GenT tagging strongly suppress dijet production. We have not attempted further suppression of the tt continuum background since there is a Z → tt signal too, especially contributing in the 2ST class. We prefer to keep an inclusive event selection to just reject dijet and bb production.
The expected significance of the Z signal in the different searches is computed by using the Monte Carlo predictions for signal plus background as pseudo-data, performing likelihood tests for the presence of narrow resonances over the expected background, using   the CL s method [53] with the asymptotic approximation of Ref. [54], and computing the p-value corresponding to each hypothesis for the resonance mass. The probability density functions of the potential narrow resonance signals are Gaussians with centre M (i.e. the resonance mass probed) and standard deviation of 0.05M (of the order of the experimental resolution). The results are given in figure 8, for a luminosity of 139 fb −1 . In the left panel we take M Z = 2 TeV with a coupling g Z = 0.3, and in the right panel we take M Z = 3 TeV, g Z = 0.5. These couplings are smaller than the 95% CL upper limits from searches, which are g Z ≤ 0.38 and g Z ≤ 0.69 for these two Z masses, respectively.
The potential significance of the Z signals in this generic search is remarkable, bearing in mind the values of the coupling we have used as benchmark. 2 For M Z = 2 TeV, g Z = 0.3, we have σ(pp → Z → Zh) = 0.6 fb, well below the current upper limit of 0.95 fb (see Fig. 2). For M Z = 3 TeV, g Z = 0.5 we have σ(pp → Z → W W ) = 0.16 fb, to be compared to the limit of 0.30 fb in Fig. 2. The sensitivity to the Z signals is driven by the 2DT class where the background can be made quite small with little reduction on the signal. To this class, Z → H 2 A 2 is the main contributing decay mode. In other event classes there are important contributions from Z → H 2 A 2 , Z → tt and Z → Zh. The decay Z → W W is only relevant for the 0T class.
Finally, let us comment that one could worry about the extra Z decay modes present in this model that might contribute to the signal in the searches targetting SM decay modes [3,4], which yield the strongent constraints on the Z coupling. For Z → Zh, which gives the most stringent limit for M Z = 2 TeV, the CMS Collaboration uses leptonic Z boson decays in Ref. [4], so the new modes Z → H 2 A 2 , Z → H 1 A 1 do not contribute. For Z → W W one expects some contribution from Z → H 1 A 1 in the signal region. Nevertheless, the limit on g Z from this channel is much looser for this Z mass. On the other hand, for M Z = 3 TeV the most stringent direct limit is from Z → W W . With the masses M H 1 ,A 1 = 96 GeV used, one expects some overlap with the jet mass window used in Ref. [3]. Nevertheless, the two-pronged jet tagging efficiencies are generically smaller for jets with b quarks, and the coupling used gives a cross section that is substantially smaller than the upper limit from this search.

Discussion
In this paper we have used as benchmark a Z model in which, in addition to decays into SM particles, Z → Zh, Z → W W , Z → tt, Z → jj, the new boson decays into new scalars Z → H i A i , i = 1, 2. The latter decay produces multi-pronged jets containing several b quarks. The goal of this paper has been to outline a general search strategy that could enable to detect such decays.
The presence of Z → Zh, Z → W W decays strongly constrains the model, that is, for a given Z mass, the non-observation of any excess in current searches [3,4] sets strong constraints on the Z coupling, which in turn constrains the possible cross section for the new modes σ(pp → Z → H i A i ) to be around 1.5 fb for M Z = 2 TeV, and 0.5 fb for M Z = 3 TeV, at most. In this quite constrained parameter space, we have seen that a generic search using (1) b tagging of sub-jets to divide the event sample into different categories according to the number of b tags; (2) a generic jet substructure tagger that selects any type of multi-pronged jets, is very effective, with an additional requirement of sizeable jet masses. (Here we have taken m J ≥ 50 GeV.) The benefit of this strategy is to be sensitive to the exotic decay modes of the Z boson that also produce two large-radius jets. In this context, it is worth pointing out that consistent models can be built in which the exotic Z decays into scalars are dominant [6], with branching ratio around 50%. In such case, a general search strategy to encompass all possible jet topologies is absolutely necessary.
We remark that more sophisticated strategies can be pursued that might enhance the significance of a possible signal. In particular, we have selected a generic requirement m J ≥ 50 GeV while in several of the event classes the signal concentrates on a relatively narrow m J interval. For example, we show in Fig. 9 the distribution of signal events for the 2DT class (the most sensitive one) in the m JJ , m J plane, where m J is the average of the two jet masses m J 1 , m J 2 in the event. Clearly, in this specific case the background can be substantially reduced, without affecting the signal, by increasing the lower cut on m J , and setting an upper cut. In an actual search, the optimisation can be done by considering a two-dimensional grid (m JJ , m J ) or even a three-dimensional one (m JJ , m J 1 , m J 2 ). The CMS Collaboration has already performed a two-dimensional search [1] for new resonances decaying into a pair of weak W, Z bosons, using a tagging variable for two-pronged jets. Therefore, generic searches such as the one proposed in this paper, with a classification by categories based on b tagging and a generic multi-pronged jet tagger, seem quite feasible and sensitive to new signals from heavy resonances.