Sensitivities of new resonance couplings to W-bosons at the LHC

The search for new particles beyond the Standard Model (SM) is a cornerstone of high-energy physics, driven by unresolved questions such as the nature of dark matter, origin of neutrino mass, baryogenesis, and dark energy. Among the various extensions to the SM, new particles coupling to W-bosons represent a particularly intriguing avenue for exploration. In physics beyond the SM (BSM), there exist various types of neutral particles X that couple to W-bosons. Examples include (i) extended gauge models with additional gauge bosons ($ Z' $) [1], which couple to $ W^+W^- $ via mixing with SM gauge bosons γ or Z; (ii) extended scalar sectors like the Two-Higgs-Doublet Model (2HDM) [2] and supersymmetric models [3], where the new CP-even scalar couples to $ W^+W^- $ via the mixing with the SM Higgs boson, which may also be related to some recent findings at the Large Hadron Collider (LHC) [4]; and (iii) pseudo-scalars, such as axions [5−9] or axion-like particles (ALPs) a, which couple to $ W^+W^- $ via a dimension-five operator $ a W^+_{\mu\nu}\widetilde{W}^{-,\mu\nu} $. Besides, new neutral particles coupling to W-bosons have also emerged in theoretical frameworks such as composite Higgs models [10], extra-dimensional scenarios [11], and dark matter models [12].

As the new neutral particle X couples to di-W-bosons, it can be produced and decay via the X-W-W vertex at the LHC. In this study, we consider X produced in association with a W-boson through the s-channel exchange of a $ W^{(*)} $ boson, followed by its decay into a $ W^+W^- $ pair, resulting in a tri-W-boson final state. The mass of X is assumed to be greater than 170 GeV, ensuring that both its di-W-boson decay products are on-shell. We choose the heavy photophobic ALP [13] as a benchmark model in this work.

Axion models were originally proposed to address the strong CP problem in QCD [9, 14]. In these models, the axion arises as a pseudo-Nambu-Goldstone boson from a spontaneously broken U(1) symmetry at a high energy scale, a mechanism known as the Peccei-Quinn mechanism [5, 6]. In such original models [5−8, 15−21], coupling of the axions is strictly related to their mass $ m_a $. Under the extended concept of ALPs, the mass and couplings of ALPs are treated as independent parameters. This flexibility significantly broadens the ALP parameter space, making them promising candidates for astrophysical and collider-based searches [22−24]. At colliders, ALPs are typically studied through their interactions with SM particles [25−35]. Among these interactions, previous experiments have primarily focused on ALPs' coupling to the diphoton $ g_{a\gamma\gamma} $. The constraints on $ g_{a\gamma\gamma} $ are stringent across most of the ALP mass range, with limits coming from astrophysical phenomena for $ m_a \lesssim {\cal{O}}(\text{MeV}) $ [36, 37], low-energy $ e^-e^+ $ collisions for $ 0.1 \, \text{ GeV} \lesssim m_a \lesssim 10 \, \text{ GeV} $ [38−40], Pb-Pb collisions at the LHC for $ 10 \, \text{ GeV} \lesssim m_a \lesssim 100 \, \text{ GeV} $ [41, 42], and p-p collisions at the LHC for $ 100 \, \text{ GeV} \lesssim m_a \lesssim 2 \, \text{ TeV} $ [29]. Consequently, the ALP-photon coupling is expected to be small, motivating the study of photophobic ALPs [13]. Photophobic ALPs are characterized by suppressed couplings to diphotons, with their primary interactions occurring with other SM electroweak bosons.

Previous studies [13, 43, 44] on heavy photophobic ALPs at the LHC have derived limits based on the tri-W-boson final state by reinterpreting experimental analyses at center-of-mass energy $ \sqrt{s} = 8 $ or 13 TeV, with limited luminosities and mass ranges. Ref. [13] derived constraints on ALP masses in the range of $ 40-500 $ GeV by reinterpreting triboson searches in the SM framework at $ \sqrt{s} = 8 \; \text{ TeV} $, where ATLAS analyses of $ pp \to W^\pm(\to \mu^\pm \nu) \, W^\pm(\to\mu^\pm \nu) \, W^\mp(\to jj) $ with $ {\cal{L}} =$ 20.3 fb^–1 [45] were utilized as an ALP-mediated signal process $ p p \to W^\pm \, a (\to W^\pm W^\mp ) $. Ref. [43] reinterpreted the CMS analyses of the SM process $pp \to jj\, W^\pm (\to \ell^\pm \nu)\, W^\pm (\to \ell^\pm \nu)$ at $\sqrt{s} = 13$ TeV [46] as a non-resonant ALP-mediated vector boson scattering process, with the ALP serving as an off-shell mediator; the results revealed sensitivities for a photophobic ALP with a mass below 100 GeV at the High-Luminosity LHC (HL-LHC) with $ \sqrt{s} = 14 $ TeV and $ {\cal{L}} = 3 $ ab^–1. Ref. [44] reinterpreted the CMS analyses for the SM process $ pp \to W^\pm W^\pm W^\mp $ at $ \sqrt{s} = 13 $ TeV with $ {\cal{L}} =$ 35.9 fb^–1 [47] as an ALP-mediated signal process $ pp \to W^\pm\, a(\to W^\pm W^\mp) $, setting limits on ALP mass up to 1 TeV.

More recently, some of our authors conducted detailed analyses of the $ \gamma Z $ decay mode for heavy photophobic ALPs at the HL-LHC with $ \sqrt{s} = 14 $ TeV [48]. In that study, ALPs were produced with two jets through s-channel vector boson exchange and fusion, resulting in the process $ pp \to jj\, a (\to \gamma\, Z (\to \ell^+ \ell^-)) $, where $ \ell = e, \mu $. Machine-learning-based multivariate analyses were used to optimize background rejection. Discovery sensitivities for the ALP coupling to di-W-bosons, $ g_{aWW} $, were evaluated over masses from 100 to 4000 GeV with $ {\cal{L}} = 3 $ ab^–1 and 140 fb^–1. Other studies on heavy photophobic ALPs at the LHC were also reviewed in Ref. [48]. Besides, heavy photophobic ALPs have also been studied through a global fit analysis of electroweak precision observables [49].

We note that because the abovementioned experimental analyses were designed to target SM production processes rather than heavy photophobic ALPs, discrepancies in signal kinematics and background characteristics frequently result in conservative or potentially unreliable discovery sensitivities in existing reinterpretation studies. In this study, we conducted detailed analyses and derived discovery sensitivities at the HL-LHC for a new heavy neutral particle X that couples to W-bosons, using the process $ p p \to W^{\pm} X (\to W^{+} W^{-}) $ with a tri-W-boson final state. Compared with previous reinterpreting studies based on tri-W-boson final state, key advancements include (i) performing simulations for the signal process of $ p p \to W^{\pm} X (\to W^{+} W^{-}) $ and related background processes to fully capture signal and background kinematics at the detector-level; (ii) using machine-learning-based multivariate analysis (MVA) for optimal signal-to-background discrimination; (iii) presenting results over a wide mass range of $ 170-3000 $ GeV at $ \sqrt{s} = 14 \, \text{ TeV} $ with $ {\cal{L}} = 3 \, \text{ab}^{-1} $; and (iv) providing the model-independent discovery sensitivities for the production cross-section time branching ratio, $ \sigma(p p \to W^{\pm} X) \times \text{Br}(X \to W^{+} W^{-}) $, which can be applied to general models involving a new heavy neutral particle X that couples to W-bosons.

The remainder of this paper is organized as follows. In Sec. II, we introduce the signal process under study. Sec. III outlines the main SM backgrounds relevant to our signal. In Secs. IV and V, we describe the simulation setup and search strategy, respectively. Our results are presented in Sec. VI, and we conclude in Sec. VII. Additional details supporting the main text are provided in the appendices.

As shown in Fig. 1, we consider the production of X in association with a W-boson at the HL-LHC with center-of-mass energy $ \sqrt{s} $ = 14 TeV. X then decays to the di-W boson. The corresponding signal process is $ p p \to {W^{\pm}} \, X (\to W^{+} W^{-}) $, leading to the tri-W-boson final state. To eliminate the background effectively, two same-charged W-bosons are required to decay leptonically into muons, leading to the existence of same-sign muons in the final state. The remanent W-boson decays hadronically to di-jet, which has larger branching ratio and can increase the signal production cross-section.

Figure 1.  Production and decay of new particle X coupling to di-W boson at $ pp $ colliders, leading to the tri-W final state. The same-charged $ W^{\pm} $ bosons are considered to decay into $ \mu^{\pm} $ and $ \nu_{\mu} $ (or $ \overline{\nu}_{\mu} $), while the remaining $ W^\mp $ boson decays into di-jet.

For a concrete study, we take the ALP theory model as an example. The neutral particle X is assumed to be a heavy photophobic ALP that couples to electroweak gauge bosons only. The effective linear ALP Lagrangian is [50]

where $ W^b_{\mu\nu} $ and $ B_{\mu\nu} $ represent the field strength tensors of the $ {S U}(2)_L $ and $ {U}(1)_Y $ gauge groups, respectively, and $ \widetilde{O}^{\mu\nu} \equiv \dfrac{1}{2}\epsilon^{\mu\nu\alpha\beta} O_{\alpha\beta} $ $ (O = W,B) $ is the dual-field strength tensor. $ m_a $ and $ f_a $ denote the mass of ALP and its decay constant ¹, respectively, and they are assumed to be independent parameters in this work.

After electroweak symmetry breaking, the interactions between ALP and gauge bosons are generally expressed as

where $ F_{\mu\nu} $, $ Z_{\mu\nu} $, and $ W^+_{\mu\nu} $ are the field strength tensors of the electromagnetic, Z-, and W-fields. In addition, the coupling constants are expressed as

where $ c_{\tilde{W}} $ and $ c_{\tilde{B}} $ are the same constants as in Eq. (1), and $ c_{\theta} = \cos{\theta}_W $, $ s_{\theta} = \sin{\theta}_W $, and $ s_{2\theta} = \sin2{\theta}_W $ are all trigonometric functions of the Weinberg angle $ {\theta}_W $. In this work, focusing on the photophobic ALP scenario, we assume the coupling between ALP and di-photon $ g_{a\gamma\gamma} = \dfrac{4}{f_a}\,(s^2_{\theta}\,c_{\tilde{W}}\,+ c^2_{\theta}\,c_{\tilde{B}}) = 0 $. As a consequence, $ s^2_{\theta}\,c_{\tilde{W}}\,+c^2_{\theta}\,c_{\tilde{B}} = 0 $, making a constant proportional relationship between $ c_{\tilde{W}} $ and $ c_{\tilde{B}} $. Thus, $ g_{a{\gamma}Z} $ and $ g_{aZZ} $ can be expressed as functions of $ g_{_{aWW}} $:

where $ t_\theta = \tan \theta_W $.

Because the signal final state contains two same-sign muons plus di-jet and moderate missing energy, relevant SM processes that can mimic the signal are listed as follows:

(i) one-boson production with two jets: $ W^\pm jj $ and $ Z jj $;

(ii) two-boson production with two jets: $ W^\pm W^\pm jj $, $ W^+ W^-jj $, $ W^\pm Zjj $, and $ ZZjj $;

(iii) three-boson production: $ W^\pm W^\pm W^\mp $;

(iv) top quark pair production: $ t\bar{t} $.

After W-boson decays, the processes $ W^\pm W^\pm jj $ and $ W^\pm W^\pm W^\mp $ can result in the same final state as the signal and serve as irreducible backgrounds, despite their relatively low production cross-sections. The process $ W^{\pm}jj $, which has a significantly higher production cross-section, also contributes to the background due to the potential misidentification of muons. Additionally, charge mis-measurement of final-state muons can lead to backgrounds from processes like $ Zjj $, $ W^+ W^- jj $, and $ t \bar{t} $. Furthermore, the $ W^{\pm} Z jj $ process becomes relevant when one of the muons in the final state remains undetected, and the $ ZZjj $ process can contribute both through missed muons and mis-measured muon charges.

We first use the MadGraph5_aMC@NLO program (version 2.6.7) [51] to simulate the proton-proton collision events at the parton level, where the NNPDF23 parton distribution function (PDF) [52] of the proton is utilized. At the parton level, to produce data as close to the experimental results as possible, the following relatively loose thresholds are applied for both the signal and background processes: (a) the minimal transverse momentum of jets, photons, and charged leptons is set to 0.5 GeV, i.e., $ p_{T}(j/\gamma/l^{\pm})> $ 0.5 GeV; (b) the maximal pseudorapidity of jets is set to 10, i.e., $ |\eta(j)| < 10 $, and it is set to 5 for the photons and charged leptons, i.e., $ |\eta(\gamma/l^{\pm})| < 5 $; and (c) the solid angular difference $ \Delta R=\sqrt{{\Delta\eta}^2+{\Delta\phi}^2} $ between objects is set to 0.1. These threshold criteria are set in the ''run_card.dat'' file of the MadGraph program.

All signal and background events are then passed through the PYTHIA program (version 8.2) [53] for parton showering, hadronization, and decay of unstable particles. We use the Delphes program (version 3.4.2) [54] to perform the detector simulation, where the ATLAS detector configuration is adopted and the minimal transverse momenta of leptons and photons accepted by the detector is set to 2 GeV. Due to limited computing resources, when simulating background processes of $ W^\pm jj $, $ Z jj $, and $ W^\pm W^\mp jj $, the decay mode of one gauge boson is fixed to decay to a muon, i.e., $ W^\pm (\to \mu^\pm \nu_\mu)\, jj $, $ Z(\to \mu^+\mu^-)\, jj $, and $ W^\pm(\to \mu^\pm \nu_\mu)\, W^\mp jj $, respectively, so that these processes can still have sufficient event numbers after all selection criteria.

For the signal, we apply the ALP model file with the linear Lagrangian [50] in the Universal FeynRules Output (UFO) format [55] in the MadGraph5 program and generate the $ p \, p \,\rightarrow {W^{\pm}} \, a(\rightarrow W^+ W^-) $ events. A scan of the ALP mass, $ m_{a} $, is performed in the following manner: individual mass 170; 25 GeV increments in the mass range 200−300; 100 GeV increments in the mass range 300− 600 GeV; individual masses 750, 800, and 900; and 400 GeV increments in the mass range 1000−3000 GeV. For each ALP mass, the coupling $ g_{_{aWW}} $ is set to 1 TeV^–1. For each iteration of ALP mass $ m_a $, at least $ 10^6 $ events are generated to ensure that statistical uncertainties are minimized as much as possible within the bounds of computational resources.

To ensure consistency in our analysis, we use the production cross-sections calculated by the MadGraph5 program to estimate the event yields for both signal and background processes. Figure 2 illustrates the signal production cross-section $\sigma(p \, p \,\rightarrow W^{\pm} \, a)$, multiplied by the branching ratio ${\rm{Br}}(a \to W^+W^-)$, as a function of the ALP mass $m_a$ in the range of 170 GeV to 3000 GeV at the HL-LHC with $\sqrt{s} = 14~ \mathrm{TeV}$, where the coupling $g_{aWW}$ is fixed to $1 ~\mathrm{TeV}^{-1}$. As $m_a$ increases beyond 160 GeV, the phase space for the decay mode $a \to W^+W^-$ expands rapidly, leading to a sharp rise in the branching ratio, which reaches approximately 50% at $ \sim190 \ \mathrm{GeV} $ [48]. Beyond this point, the branching ratio increases more gradually, with only minor changes at higher masses. At the same time, the production cross-section for process $p \, p \,\rightarrow W^{\pm} \, a$ decreases steadily with increasing $m_a$. Consequently, the production cross-section times branching ratio peaks when $m_a \sim 190~ \mathrm{GeV}$, just above the $WW$ decay threshold of 160 GeV.

Figure 2.  (color online) Signal production cross-section $ \sigma(p \, p \,\rightarrow $$ {W^{\pm}} \, a) $ multiplied by branching ratio $ {\rm{Br}}(a \to W^+W^-) $ as a function of ALP mass $ m_{_{a}} $ from 170 to 3000 GeV at the HL-LHC with $ \sqrt{s} = 14 $ TeV, assuming the coupling $ g_{_{aWW}} = $ 1 TeV^–1.

The following preselection criteria were employed before a multivariate analysis was performed, where final-state muons and jets are ordered based on their transverse momenta and labeled as $O_i$ ($i = 1, 2, \ldots$), with $O = \mu, \, j$, respectively:

(1) Events are required to have exactly two muons, i.e., $ N (\mu) = $ 2.

(2) The minimal $ p_{_{\rm{T}}} $ of the two muons is 10 GeV, i.e., $ p_{_{\rm{T}}} (\mu) > 10 $ GeV.

(3) The two muons have the same charge, i.e., $ \mu^\pm \mu^\pm $.

(4) The minimal $ p_{_{\rm{T}}} $ of all jets is 30 GeV, i.e., $ p_{_{\rm{T}}}(j) > 30 $ GeV.

(5) Event are required to have at least two jets, i.e., $ N(j) \geq $ 2, and the number of $ b- $tagged jets must be zero, i.e., $ N(j_b) = 0 $.

Table 1 presents the expected numbers of events, $ N_{\rm{exp}} $, for the signal, with a benchmark mass of $m_a = $ 500 GeV and coupling $g_{aWW} = 1\ \mathrm{TeV}^{-1}$, as well as for the background processes. Here, $ N_{\rm{exp}} $ is calculated by the formula

	initial	(1)-(2)	(3)	(4)-(5)
$ W^\pm a(\rightarrow W^\pm W^\mp) $	$ 5.51 {\mkern-2mu\times\mkern-2mu} 10^{4} $	$ 1.27 {\mkern-2mu\times\mkern-2mu} 10^{3} $	$ 4.26 {\mkern-2mu\times\mkern-2mu} 10^{2} $	$ 2.54 {\mkern-2mu\times\mkern-2mu} 10^{2} $
$ W^{\pm}(\rightarrow \mu^{\pm}\nu_{\mu}) jj $	$ 3.62 {\mkern-2mu\times\mkern-2mu} 10^{11} $	$ 7.24 {\mkern-2mu\times\mkern-2mu} 10^{6} $	$ 2.69 {\mkern-2mu\times\mkern-2mu} 10^{6} $	$ 4.89 {\mkern-2mu\times\mkern-2mu} 10^{5} $
$ Z(\rightarrow \mu^{\pm}\mu^{\mp})jj $	$ 3.31 {\mkern-2mu\times\mkern-2mu} 10^{10} $	$ 1.13 {\mkern-2mu\times\mkern-2mu} 10^{10} $	$ 3.33 {\mkern-2mu\times\mkern-2mu} 10^{5} $	$ 6.67 {\mkern-2mu\times\mkern-2mu} 10^{4} $
$ W^\pm W^\pm jj $	$ 1.10 {\mkern-2mu\times\mkern-2mu} 10^{6} $	$ 5.22 {\mkern-2mu\times\mkern-2mu} 10^{3} $	$ 5.19 {\mkern-2mu\times\mkern-2mu} 10^{3} $	$ 3.88 {\mkern-2mu\times\mkern-2mu} 10^{3} $
$ W^\pm(\rightarrow \mu^{\pm}\nu_{\mu}) W^\mp jj $	$ 5.05 {\mkern-2mu\times\mkern-2mu} 10^{8} $	$ 2.43 {\mkern-2mu\times\mkern-2mu} 10^{7} $	$ 8.81 {\mkern-2mu\times\mkern-2mu} 10^{3} $	$ 4.53 {\mkern-2mu\times\mkern-2mu} 10^{3} $
$ W^\pm Zjj $	$ 1.03 {\mkern-2mu\times\mkern-2mu} 10^{9} $	$ 1.42 {\mkern-2mu\times\mkern-2mu} 10^{7} $	$ 6.17 {\mkern-2mu\times\mkern-2mu} 10^{5} $	$ 1.99 {\mkern-2mu\times\mkern-2mu} 10^{5} $
$ ZZjj $	$ 3.10 {\mkern-2mu\times\mkern-2mu} 10^{8} $	$ 8.42 {\mkern-2mu\times\mkern-2mu} 10^{6} $	$ 4.32 {\mkern-2mu\times\mkern-2mu} 10^{4} $	$ 1.55 {\mkern-2mu\times\mkern-2mu} 10^{4} $
$ W^\pm W^\pm W^\mp $	$ 3.79 {\mkern-2mu\times\mkern-2mu} 10^{5} $	$ 5.92 {\mkern-2mu\times\mkern-2mu} 10^{3} $	$ 1.73 {\mkern-2mu\times\mkern-2mu} 10^{3} $	$ 8.55 {\mkern-2mu\times\mkern-2mu} 10^{2} $
$ t\bar{t} $	$ 1.79 {\mkern-2mu\times\mkern-2mu} 10^{9} $	$ 1.22 {\mkern-2mu\times\mkern-2mu} 10^{7} $	$ 8.71 {\mkern-2mu\times\mkern-2mu} 10^{4} $	$ 2.91 {\mkern-2mu\times\mkern-2mu} 10^{4} $

Table 1.  Expected numbers of events for the signal with benchmark mass $ m_a $ = 500 GeV and coupling $ g_{_{aWW}} = $ 1 TeV^–1 and background processes after applying preselection criteria (1)−(5) sequentially at the HL-LHC with $ \sqrt{s} = $ 14 TeV and $ {\cal{L}} = $ 3 $ {\rm{ab}}^{-1} $.

where $ \sigma_{\rm{pro}} $ is the production cross-section of the signal or background process; $ {\cal{L}} $ is the integrated luminosity; and $ \epsilon_{\rm{pre}} $ is the preselection efficiency, which is evaluated based on our analyses. The numbers are shown after applying preselection criteria (1)−(5) sequentially at the HL-LHC, with $\sqrt{s} = 14\ \mathrm{TeV}$ and ${\cal{L}} = 3\ \mathrm{ab}^{-1}$. As shown in Table 1, after applying all preselection criteria, the number of background events for most processes is reduced by at least $10^4$-fold. Nevertheless, the total background rate remains much larger than the expected number of signal events.

After applying the preselection criteria, the Toolkit for Multivariate Analysis (TMVA) package [56] is employed to perform a multivariate analysis (MVA), which enhances the ability to distinguish between signal and background events. The events that pass the preselection steps are then subjected to analysis using the Boosted Decision Trees (BDT) algorithm within the TMVA framework. In this process, the discrimination between signal and background is achieved by utilizing the following set of kinematic observables:

(a) Missing transverse energy and its azimuth angle: $ \not {E}_{{\rm{T}}} $, $ \phi ({ \not {E}_{{\rm{T}}}}) $.

(b) The x, y, and z-components of momentum ($ p_x $, $ p_y $, $ p_z $) and energy (E) of the first two leading jets and muons: $ p_x (j_1) $, $ p_y (j_1) $, $ p_z (j_1) $, $ E (j_1) $; $ p_x (j_2) $, $ p_y (j_2) $, $ p_z (j_2) $, $ E (j_2) $; $ p_x (\mu_1) $, $ p_y (\mu_1) $, $ p_z (\mu_1) $, $ E (\mu_1) $; $ p_x (\mu_2) $, $ p_y (\mu_2) $, $ p_z (\mu_2) $, $ E (\mu_2) $.

(c) The number of charged tracks ($ N_{\rm{track}} $) and ratio ($ R_{\rm{E}} $) of the hadronic versus electromagnetic energy deposited in the calorimeter cells for the first two leading jets: $ N_{\rm{track}}(j_1) $, $ N_{\rm{track}}(j_2) $, $ R_{\rm{E}}(j_1) $, and $ R_{\rm{E}}(j_2) $. $ R_{\rm{E}} $ is typically greater than one for a jet.

(d) The solid angular separation $ \Delta R $ between $ j_1 $ and $ j_2 $, = invariant mass (m) of system of ($ j_1+j_2 $), and number (N) of jets: $ \Delta R(j_1,\, j_2) $, $ m(j_{1}+j_{2}) $, and $ N(j) $.

(e) Observables related to the isolation quality of muons. (i) The summed $ p_{\rm{T}} $ of other objects excluding the muon in an $ R = 0.4 $ cone around the muon: $ p^{\rm{iso}}_{\rm{T}}(\mu_1) $, $ p^{\rm{iso}}_{\rm{T}}(\mu_2) $, and the larger of the two $ p^{\rm{iso}}_{ {\rm{T}},\,{\rm{max}}}(\mu) $. (ii) The ratio of the transverse energy in a $ 3 \times 3 $ grid surrounding the muon to the $ p_{_{\rm{T}}} $ of the muon (''etrat'' in Ref. [57]), which is a percentage between 0 and 0.99: $ R_{\rm{grid}}(\mu_1) $, $ R_{\rm{grid}}(\mu_2) $. (iii) We make all possible combinations between $ \mu_1 / \mu_2 $ and every jet, compare the $ \Delta R $, and find the minimal value: $ \Delta R_{\rm{min}}(\mu,\, j) $. For well-isolated muons, $ p^{\rm{iso}}_{\rm{T}} $and $ R_{\rm{grid}} $ should be small, while $ \Delta R_{\rm{min}}(\mu,\, j) $ is large.

(f) To reconstruct the hadronically decaying W-boson ($ jj_{_W} $), we make all di-jet combinations and select the pair with the invariant mass closest to 80 GeV. The corresponding two jets are labeled as $ j_{_{W1}} $ and $ j_{_{W2}} $ sorted by their $ p_{_{\rm{T}}} $. We input the following observables of these two jets: $ p_x (j_{_{W1}}) $, $ p_y(j_{_{W1}}) $, $ p_z(j_{_{W1}}) $, $ E(j_{_{W1}}) $; $ p_x (j_{_{W2}}) $, $ p_y(j_{_{W2}}) $, $ p_z(j_{_{W2}}) $, $ E(j_{_{W2}}) $; $ N_{\rm{track}}(j_{_{W1}}) $, $ N_{\rm{track}}(j_{_{W2}}) $; $ R_{\rm{E}}(j_{_{W1}}) $, $ R_{\rm{E}}(j_{_{W2}}) $; $ \Delta R(j_{_{W1}},\, j_{_{W2}}) $, $ m(j_{_{W1}}+j_{_{W2}}) $.

(g) The transverse mass $ m_{\rm{T}} $ [58] of the system, which includes one muon and the missing transverse energy: $ m_{\rm{T}}(\mu_{1}+ \not {E}_{{\rm{T}}}) $, $ m_{\rm{T}}(\mu_{2}+ \not {E}_{{\rm{T}}}) $. Here, the transverse mass $ m_{\rm{T}} \equiv \sqrt{ (E_{\rm{T}}^{\rm{vis.}}+ \not {E}_{{\rm{T}}})^2 - ( \vec{p}_{\rm{T}}^{\rm{\,\,vis.}} + \vec{\not {p}}_{\rm{T}} )^2 } $, where $ E_{\rm{T}}^{\rm{vis.}} $ ($ \vec{p}_{\rm{T}}^{\rm{\,\,vis.}} $) is the transverse energy (momentum) of the visible object or system, while $ \not {E}_{{\rm{T}}} $ ($ \vec{\not {p}}_{\rm{T}} $) is the missing transverse energy (momentum). $ E_{\rm{T}}^{\rm{vis.}} = \sqrt{ (\vec{p}_{\rm{T}}^{\rm{\,\,vis.}})^2 + (m^{\rm{vis.}})^2 } $, where $ m^{\rm{vis.}} $ is the invariant mass of the visible object or system. $ \not {E}_{{\rm{T}}} = |\vec{\not {p}}_{\rm{T}}| $, assuming the invariant mass of the invisible object is zero.

(h) We calculate and compare $ \Delta R $ values of two combinations between $ jj_{_W} $ and $ \mu_1 / \mu_2 $. The muon with smaller (larger) $ \Delta R $ value is labeled as $ \mu_a $ ($ \mu_b $). We input the minimal value $ \Delta R(j_{_{W1}}+j_{_{W2}},\, \mu_a) $ and maximal value $ \Delta R(j_{_{W1}}+j_{_{W2}},\, \mu_b) $.

(i) The observables related to the reconstruction of the ALP mass are $ m(j_{_{W1}}+j_{_{W2}}+\mu_a) $, $ m(j_{_{W1}}+j_{_{W2}}+\mu_b) $, $ m_{\rm{T}}(j_{_{W1}}+ j_{_{W2}}+ \mu_a+ \not {E}_{{\rm{T}}}) $, and$ m_{\rm{T}}(j_{_{W1}}+j_{_{W2}}+\mu_b+ \not {E}_{{\rm{T}}}) $.

(j) The ovservables related to the reconstruction of the off-shell $ W^{(*)\pm} $ boson: $ m_{\rm{T}}(j_{_{W1}}+j_{_{W2}}+\mu_a+\mu_b+ \not {E}_{{\rm{T}}}) $.

Details of the observables $ N_{\rm{track}}(j_{i}) $, $ R_E(j_{i}) $, $ R_{\rm{grid}}(\mu_i) $, and $ p^{\rm{iso}}_{\rm{T}}(\mu) $ can be found in [57]. The hadronic W-boson ($jj_{_W}$) is reconstructed using the method described above. Appendix A presents distributions of kinematical observables after applying preselection criteria for signal and background processes at the HL-LHC with $\sqrt{s} = 14 \ {\rm{TeV}} $. Among them, for $ m_a $ values below 500 GeV, the endpoint positions of the distributions for both $ m(j_{_{W1}} + j_{_{W2}} + \mu_a) $ and $ m_{\rm{T}}(j_{_{W1}} + j_{_{W2}} + \mu_a + \not {E}_{{{\rm{T}}}}) $ align closely with $ m_a $; when $ m_a $ exceeds 500 GeV, the endpoint positions of the $ m_{\rm{T}}(j_{_{W1}} + j_{_{W2}} + \mu_a + \not {E}_{{{\rm{T}}}}) $ distribution are still related to $ m_a $, but this phenomenon becomes less obvious. Appendix B provides detailed information about the MVA analyses performed using the TMVA package [56].

Figure 3 presents the BDT response distributions for the total background and signal at the HL-LHC, with benchmark masses of $m_a = 400 \text{ GeV}$ and $m_a = 900 \text{ GeV}$. The left and right plots show similar trends in the distributions for both mass points. The clear separation between the signal and background distributions indicates that the BDT criteria are effective for background rejection at the HL-LHC. Additionally, the separation is more pronounced for $m_a = 900 \text{ GeV}$ compared to $m_a = $ 400 GeV, suggesting improved discriminating power at higher masses. Appendix C presents distributions of BDT responses after applying preselection criteria for the signal and background processes at the HL-LHC with $ \sqrt{s} = $ 14 TeV, assuming various $ m_a $ cases.

Figure 3.  (color online) BDT response distributions for total SM background and signal when ALP mass $ m_a $ = 400 GeV (left) and 900 GeV (right) at the HL-LHC with $ \sqrt{s} = $ 14 TeV.

After the preselection, the BDT cut is optimized according to the signal statistical significance calculated by Eq. (6) for each mass case:

where $ N_s $ and $ N_b $ are the expected numbers of events for signal and total background, respectively, after applying both the preselection and BDT criteria. Appendix D shows selection efficiencies of preselection and BDT criteria for signal and background processes at the HL-LHC with $ \sqrt{s} = $ 14 TeV, assuming different ALP masses, where ''−'' means the number of events can be reduced to negligible with $ {\cal{L}} = $ 3 ab⁻¹.

Using our search strategy, Fig. 4 presents the discovery sensitivities for the production cross-section $ \sigma(p p \to W^{\pm} X) $ multiplied by the branching ratio $ \text{Br}(X \to W^{+} W^{-}) $ as a function of $ m_X $ in the mass range of 170–3000 GeV at the HL-LHC, with $ \sqrt{s} = 14 \text{ TeV} $ and an integrated luminosity of $ {\cal{L}} = 3 \text{ab}^{-1} $. The red and green curves correspond to the 2-σ and 5-σ significances, respectively. As shown in the BDT distributions in Appendix C, the separation between the signal and background improves significantly for larger masses, leading to more effective background rejection. Consequently, the sensitivities for heavier masses do not decrease very rapidly, demonstrating strong discovery potential across a wide mass range.

Figure 4.  (color online) Discovery sensitivities on the production cross-section $ \sigma(p p \to W^{\pm} \, X) $ times branching ratio Br$ (X \to W^{+} W^{-}) $ in the mass range of 170−3000 GeV at the HL-LHC, with $ \sqrt{s} = $ 14 TeV and $ {\cal{L}} = $ 3 $ {\rm{ab}}^{-1} $. Red and green curves correspond to 2-σ and 5-σ significances, respectively.

In Fig. 5, for the concrete case of heavy photophobic ALPs, we show the discovery sensitivities on the coupling $ g_{aWW} $ with 2-σ and 5-σ significances as a function of the ALP mass $ m_a $ in the range of 170–3000 GeV at the HL-LHC, with $ \sqrt{s} = 14 \text{ TeV} $ and integrated luminosities of $3 \text{ab}^{-1} $ and 35.9 fb^–1. For comparison, we also display the 95% C.L. limit (blue curve) from Ref. [44], which is derived from reinterpreting the CMS analyses of the SM production $ pp \to W^\pm W^\pm W^\mp $ at $ \sqrt{s} = 13 \text{ TeV} $ with $ {\cal{L}} = 35.9 \text{ fb}^{-1} $. Our results demonstrate that at the HL-LHC with $ {\cal{L}} = 3 \text{ab}^{-1} $, the 2-σ sensitivity for $ g_{aWW} $ decreases from $ 1.06 \text{ TeV}^{-1} $ to $ 0.75 \text{ TeV}^{-1} $ as $ m_a $ increases from 170 GeV to 225 GeV, remains relatively stable for $ m_a $ in the range of 225–600 GeV, and then rises sharply to $ 9.42 \text{ TeV}^{-1} $ at $ m_a = 3000 \text{ GeV} $. This pattern reflects the varying sensitivity of the search strategy across different mass ranges.

Figure 5.  (color online) Discovery sensitivities with 2-σ and 5-σ significances on the coupling $ g_{_{aWW}} $ in the mass range of 170−3000 GeV at the HL-LHC, with $ \sqrt{s} = $ 14 TeV and $ {\cal{L}} = $ 3 $ {\rm{ab}}^{-1} $ and 35.9 fb^–1. The 95% C.L. limit (blue curve) [44] derived from reinterpreting CMS analyses of SM production $ pp \to W^\pm W^\pm W^\mp $ at $ \sqrt{s} = $ 13 TeV and $ {\cal{L}} = $ 35.9 fb^–1 is displayed for comparison.

One observes that with the same luminosity of 35.9 $ {\rm{fb}}^{-1} $, our strategy achieves better sensitivities for $ m_a $ above 300 GeV compared to the referenced approach [44]. This improvement is primarily due to differences in the optimization of the analysis strategies. The referenced studies focused on the SM process $ pp \to W^\pm W^\pm W^\mp $, while our strategy is specifically tailored to the tri-W final state produced by a heavy resonance. For larger masses, the kinematics of these two processes differ significantly, allowing our approach to reject background more effectively, as demonstrated by the BDT distributions in Appendix C.

In this paper, we proposed a search strategy at the HL-LHC for a new neutral particle X that couples to a pair of W-bosons. The particle X is produced with a W-boson and decays into two W-bosons, resulting in the process $ p p \to W^{\pm} X (\to W^{+} W^{-}) $ and a tri-W-boson final state. To suppress background, we focused on events where two same-charge W-bosons decay leptonically into muons, producing same-sign muons in the final state. The third W-boson decays hadronically into jets, leveraging its higher branching ratio to enhance the signal. As a case study, we used the heavy photophobic ALP as an example.

Signal and background events were simulated at the detector level. The signal production cross-section times branching ratio, $ \sigma(p \, p \rightarrow W^{\pm} \, a) \times \text{Br}(a \to W^+W^-)$, was evaluated as a function of the ALP mass $m_a$ at the HL-LHC with $\sqrt{s} = 14 \text{ TeV}$. Events were selected with exactly two same-sign di-muons and at least two non-b-tagged jets. A machine-learning-based MVA wasapplied to enhance signal-background discrimination. Distributions of key input variables and the corresponding BDT responses were presented, and selection efficiencies for both preselection and BDT criteria were provided for various ALP masses.

We presented the model-independent discovery sensitivities for the production cross-section times branching ratio $ \sigma(p p \to W^{\pm} X) \times \text{Br}(X \to W^{+} W^{-}) $ as a function of the particle mass $ m_X $ in the range of 170–3000 GeV at the HL-LHC, with $ \sqrt{s} = 14 \text{ TeV} $ and $ {\cal{L}} = 3 \text{ab}^{-1} $, which can be applied to general models involving a new neutral particle X that couples to W-bosons. The sensitivities for the coupling $ g_{aWW} $ of heavy photophobic ALPs were also shown for 2-σ and 5-σ significances in the same mass range under luminosities of $ 3 \text{ab}^{-1} $ and 35.9 $\text{fb}^{-1}$. At $ {\cal{L}} = 3 \text{ab}^{-1} $, the 2-σ sensitivity for $ g_{aWW} $ decreased from $ 1.06 \text{ TeV}^{-1} $ to $ 0.75 \text{ TeV}^{-1} $ as $ m_a $ increased from 170 GeV to 225 GeV, remained almost constant for $ m_a $ between 225 and 600 GeV, and then sharply rose to $ 9.42 \text{ TeV}^{-1} $ at $ m_a = 3000 \text{ GeV} $. Compared to previous limits derived from reinterpreting CMS analyses of SM production $ pp \to W^\pm W^\pm W^\mp $ at $ \sqrt{s} = 13 \text{ TeV} $, our strategy provides improved sensitivity for $ m_a $ above 300 GeV.

In this appendix, we provide detailed information about the MVA analyses performed using the TMVA package [56]. Figure B1 illustrates the importance values of input observables as a function of $ m_a $, evaluated during the MVA training. The importance values for BDT input variables are quantified by counting how often the variables are used to split decision tree nodes and by weighting each split occurrence by the separation gain-squared it has achieved and by the number of events in the node [56, 59]. The left plot highlights observables with the highest importance values, while those with the lowest importance values are displayed in the right panel.

Figure B1.  (color online) Importance values of input observables as a function of $ m_a $, evaluated during the MVA training.

Table B1 presents the correlations between two observables for the total background and signal with benchmark mass $ m_a = 400 $ GeV. The listed 12 pairs have strong correlations and remain stable under variations of $ m_a $. Figure B2 displays the correlations between pairs of observables for the signal as a function of $ m_a $. The correlations for the selected pairs exhibit significant dependence on $ m_a $.

observables	signal	background
$ p^{\rm{iso}}_{\rm{T}}(\mu_1) $, $ p^{\rm{iso}}_{{\rm{T}},\,{\rm{max}}}(\mu) $	$ 1.000 $	$ 0.998 $
$ m_{\rm{T}}(j_{_{W1}}+j_{_{W2}}+\mu_b+ \not {E}_{{\rm{T}}}) $, $ m_{\rm{T}}(j_{_{W1}}+j_{_{W2}}+\mu_a+\mu_b+ \not {E}_{{\rm{T}}}) $	$ 0.919 $	$ 0.974 $
$ m(j_{_{W1}}+j_{_{W2}}+\mu_b) $, $ m_{\rm{T}}(j_{_{W1}}+j_{_{W2}}+\mu_a+\mu_b+ \not {E}_{{\rm{T}}}) $	$ 0.797 $	$ 0.955 $
$ m(j_{_{W1}}+j_{_{W2}}+\mu_a) $, $ m(j_{_{W1}}+j_{_{W2}}) $	$ 0.845 $	$ 0.865 $
$ R_{\rm{E}}(j_1) $, $ R_{\rm{E}}(j_{_{W1}}) $	$ 0.762 $	$ 0.896 $
$ E (j_1) $, $ E (j_{_{W1}}) $	$ 0.747 $	$ 0.850 $
$ m(j_{_{W1}}+j_{_{W2}}+\mu_b) $, $ m_{\rm{T}}(j_{_{W1}}+j_{_{W2}}+\mu_b+ \not {E}_{{\rm{T}}}) $	$ 0.724 $	$ 0.926 $
$ m(j_{_{W1}}+j_{_{W2}}+\mu_b) $, $ m(j_{_{W1}}+j_{_{W2}}) $	$ 0.710 $	$ 0.873 $
$ m_{\rm{T}}(j_{_{W1}}+j_{_{W2}}+\mu_b+ \not {E}_{{\rm{T}}}) $, $ m(j_{_{W1}}+j_{_{W2}}) $	$ 0.715 $	$ 0.861 $
$ p_z(j_1) $, $ p_z(j_{_{W1}}) $	$ 0.791 $	$ 0.838 $
$ m(j_{_{W1}}+j_{_{W2}}+\mu_b) $, $ m(j_{_{W1}}+j_{_{W2}}+\mu_a) $	$ 0.740 $	$ 0.861 $
$ N_{\rm{track}}(j_1) $, $ N_{\rm{track}}(j_{_{W1}}) $	$ 0.761 $	$ 0.772 $

Table B1.  Correlations between two observables for the total background and signal with benchmark mass $ m_a = 400 $ GeV.

Figure B2.  (color online) Correlations between pairs of observables for the signal as a function of $ m_a $.

We thank Bin Diao and Zilong Ding for helpful discussions.