Identifying a $Z'$ behind $b \to s \ell \ell$ anomalies at the LHC

Recent $b\to s\ell\ell$ anomalies may imply the existence of a new $Z'$ boson with left-handed $Z'bs$ and $Z'\mu\mu$ couplings. Such a $Z'$ may be directly observed at LHC via $b \bar s \to Z' \to \mu^+\mu^-$, and its relevance to $b\to s\ell\ell$ may be studied by searching for the process $gs \to Z'b \to \mu^+\mu^- b$. In this paper, we analyze the capability of the 14 TeV LHC to observe the $Z'$ in the $\mu^+ \mu^-$ and $\mu^+\mu^- b$ modes based on an effective model with major phenomenological constraints imposed. We find that both modes can be discovered with 3000 fb$^{-1}$ data if the $Z'bs$ coupling saturates the latest $B_s-\bar B_s$ mixing limit from UTfit at around $2\sigma$. Besides, a tiny right-handed $Z'bs$ coupling, if it exists, opens up the possibility of a relatively large left-handed counterpart, due to cancellation in the $B_s-\bar B_s$ mixing amplitude. In this case, we show that even a data sample of $\mathcal{O}(100)$ fb$^{-1}$ would enable discovery of both modes. We further study the impact of a $Z'bb$ coupling as large as the $Z'bs$ coupling. This scenario enables discovery of the $Z'$ in both modes with milder effects on the $B_s-\bar B_s$ mixing, but obscures the relevance of the $Z'$ to $b \to s\ell\ell$. Discrimination between the $Z'bs$ and $Z'bb$ couplings may come from the production cross section for the $Z'b\bar{b}$ final state. However, we do not find the prospect for this to be promising.


I. INTRODUCTION
Bottom-quark transitions of b → s have been of great interest as a means for studying physics beyond the standard model (SM) since the observation of the decay B → K * γ by the CLEO collaboration [1]. Involving a flavor-changing neutral current (FCNC), such processes are possible in the SM only at loop level, providing unique sensitivity to new physics (NP).
While the statistical significance of each discrepancy is not large enough, there is excitement about the possibility that their combination might suggest the presence of NP. To investigate this possibility, global-fit analyses based on the effective Hamiltonian formalism have been performed by several groups (see Refs. [15][16][17][18][19][20][21][22] for studies that include the recent R K * [6] result). These fits find that the tensions in the b → s observables can be simultaneously alleviated to a great extent by a NP contribution in a single Wilson coefficient. Most authors find this to be the coupling between the left-handed (LH) b → s current and either the LH or vector muon current. 1 In particular, the findings of the global-fit analyses motivate studying a new gauge boson, Z , with FCNC interactions. Many phenomenological studies of such a Z have been performed (see, e.g., ).
Produced at LHC via bs → Z and undergoing the decay Z → µ + µ − , the Z may be discovered in dimuon resonance searches. Such a discovery by itself, however, would not reveal the relevance of the observed resonance to the b → s anomalies. Comparison with searches in the electron mode can test lepton-flavor universality in the Z couplings, but one should also establish the coupling to the b → s current. In principle, this can be done with the decay Z → bs. However, this decay is suppressed relative to Z → µ + µ − due to the B s −B s mixing constraint, as we discuss below, and its detection suffers from overwhelming QCD background. On the other hand, the Z bs coupling can induce the process gs → Z b. Therefore, this coupling can be explored through production modes of the Z at LHC.
In this paper, we investigate the prospects for direct observation of the Z , as well as determination of the flavor structure of its couplings at LHC, with √ s = 14 TeV. We argue that achieving this dual goal requires measuring the cross sections of both pp → Z + X and pp → Z b + X, where X refers to additional activity in the pp collision. We employ an effective-model description of the Z couplings, and assume that it is the source of all the NP required to alleviate the tensions in b → sµ + µ − . We focus mainly on the role of the Z bs coupling in the Z production processes pp → Z + X and pp → Z b + X. The constraints on the leptonic coupling of the Z are much weaker than those on the Z bs coupling. Therefore, we use Z → µ + µ − as the main discovery mode. As the LH Z bs coupling is accompanied by a LH Z bb coupling in many UV complete models (see, e.g. Refs. [31,51,58,61]), we also study scenarios in which the Z bb coupling is of the same order as the Z bs coupling. We note that a larger Z bb coupling implies larger cross sections and easier discovery of the Z at LHC, but obscures the role of the Z in b → s . Therefore, this case is not the focus of our work. We also consider the process pp → Z bb + X for discrimination between the Z bs and Z bb couplings, where the latter coupling uniquely contributes via gg → Z bb.
We emphasize that our purpose is different from that of existing studies, which focus on discovery and/or constraint of the Z rather than testing its role in b → s . Recent studies on the impact of existing dimuon resonance searches on a Z motivated by the b → s anomalies can be found, e.g., in Refs. [64,71]. Ref. [71] also studies the future sensitivity at LHC, exploiting the use of additional b jets for background suppression. However, it targets a Z bb coupling that is much larger than the Z bs coupling. The future sensitivities at LHC and a 100 TeV pp collider are studied in Ref. [72] with an extrapolation of existing ATLAS limit.
We will show that the cross sections for bs → Z and gs → Z b are limited by a rather tight constraint on the Z bs coupling, which originates from the B s −B s mixing. This severely restricts the discovery potential of the Z , unless the Z bb coupling is comparable or larger than the Z bs coupling. As the B s −B s mixing constraint is indirect, the actual limit on the Z bs coupling depends on the details of the UV-complete model. In particular, a nonzero right-handed (RH) Z bs coupling can drastically change the constraint, due to the large (V − A) ⊗ (V + A) term [79] in the B s −B s mixing amplitude. Although there is no strong indication of a RH b → s current in the majority of the global-fit analyses, even a tiny RH Z bs coupling would allow for a large LH Z bs coupling due to the cancellation. This would significantly boost the Z production cross sections. We investigate the discovery potential in this case as well.
This paper is organized as follows. In Sec. II, we introduce the effective model for our collider study. In Sec. III, we evaluate existing phenomenological constraints on the relevant couplings of the Z boson to quarks and leptons. This is carried out for two representative Z mass values, m Z = 200 and 500 GeV. In Sec. IV, we study the signal and background cross sections for the three processes of interest, pp → Z + X, Z b + X and Z bb + X, given the coupling constraints. We then proceed to estimate the signal significances for the full integrated luminosity of the High-Luminosity LHC (HL-LHC) program, L = 3000 fb −1 . In Sec. V, we discuss the impact of a tiny but nonzero RH Z bs coupling, which allows discovery with smaller integrated luminosities. Summary and discussions are given in Sec. VI.

II. EFFECTIVE MODEL
We describe the Z couplings to the SM fermions with the effective Lagrangian where P L,R = (1 ∓ γ 5 )/2, and g L bb , g L bs and g L,R µµ are coupling constants. For simplicity, since no significant CP violation has been observed in the relevant observables, we take the couplings to be real. In addition to the LH Z bs coupling and LH and RH Z µµ couplings motivated by the b → s + − global fits, we introduce the LH Z bb coupling predicted in many UV complete models. We take g L µµ to also be the coupling to the muon neutrino, as required by the SU(2) L gauge symmetry, and since observables containing neutrinos give meaningful constraints on the Z parameters, as shown below. The RH Z bs coupling is not included at this stage, as it is discussed only in Sec. V.
Since we take m Z m b , we integrate out the Z to obtain its contributions to the effective Hamiltonian for b → sµ + µ − transitions, given by (2) where N = − αG F √ 2π V tb V * ts , and are the Wilson coefficients, with the vector and axialvector muon couplings defined by With the parametrization of the Cabibbo-Kobayashi-Maskawa (CKM) matrix in Ref. [80], C NP 9,10 are treated as real to a good approximation.
Motivated by the global-fit analyses, we consider two possibilities for the chiral structure of the muon couplings: a vector coupling and a LH coupling. For each of these, we extract constraints on the couplings from global-fit analyses presented in Ref. [15]. The first analysis uses all available b → s data from LHCb, ATLAS, CMS and Belle. This yields the following two scenarios: (i) Vector coupling (g L µµ = g R µµ ) and, hence, C NP 10 = 0. In this case, the best fit value is Re C NP 9 = −1.11, with a 2-standard-deviation (2σ) range of (ii) LH coupling (g R µµ = 0), so that C NP 9 = −C NP 10 . For this scenario, the best fit value is Re C NP The authors of Ref. [15] also present results when taking into account only lepton-flavor universality observables, such as R K ( * ) measured by LHCb [5,6] and the differences Q 4 and Q 5 [81] between angular observables in B → K * µ + µ − and B → K * e + e − , measured by Belle [13]. This leads to the following two scenarios: In the rest of the paper we explore scenario (i) in detail, and occasionally comment on differences with respect to the other scenarios. Generally, these differences are small and do not significantly affect our main results in Section IV.

III. ALLOWED PARAMETER SPACE
In Fig. 1 we show the various constraints on g L bs vs. g V µµ in scenario (i), for the representative Z -mass values of m Z = 200 and 500 GeV. The relevant inputs and constraint calculation methods are described in the remainder of this section.
Before embarking on this detailed description, we note that, as seen in Fig. 1, the limits on g V µµ are much weaker than those on g L bs . Therefore, the Z is likely to decay primarily to leptons, so that its decays into quarks can be ignored. The leptonic-decay dominance simplifies the discussion. In fact, it is also essential for direct observation of the Z at LHC, since searches with Z → bb, bs suffer from overwhelming QCD background. Thus, in the scenario (i), where g L µµ = g R µµ = 0, the dominant branching ratios are In scenario (ii), each of these branching ratios is 50%. To incorporate the results of b → s global fits, we use Eq. (3) to convert the 2σ constraints of Eq. (5) to the space of g L bs vs. g V µµ . The result is given by the blue hyperbolae in Fig. 1 As an example of the impact of the choice of scenario among those listed in Sec. II, we note that in scenario (i'), the resulting values of g V µµ are generically higher than those shown in Fig. 1, and have a wider range. Although this somewhat changes the Z → µ + µ − branching ratio, Eq. (7) is still satisfied. Therefore, the results we obtain in Section IV are not affected.
The most important constraint on the Z bs coupling comes from the B s mixing. With the tree-level Z exchange contribution, the total B s −B s mixing amplitude relative to the SM one (see, e.g. Ref. [82]) is given by where S 0 2.3 is an Inami-Lim function, g 2 is the SU(2) L gauge coupling, v 246 GeV, and a common QCD correction factor is assumed for the SM and NP contributions. The mass difference ∆m B 0 s = 2|M 12 | is precisely measured, at the per mill level [80], while the calculation of M 12 suffers from various sources of uncertainty. One of the dominant uncertainties is the CKM factor, with an uncertainty of ∼ 5% [83,84]. The other is the hadronic matrix element, obtained from the lattice. The average of N f = 2 + 1 lattice results compiled by FLAG in 2016 [85] implied a ∼ 12% uncertainty in M SM 12 . Recently, the situation was greatly improved with the advent of the accurate estimate by the Fermilab Lattice and MILC collaborations [86], which pushes down the uncertainty of the FLAG average to ∼ 6%. (See December 2017 update on the FLAG website [85].) A global analysis of the CKM parameters by CKMfitter [83], published before the recent lattice result [86], gives a constraint on M 12 with NP as |M 12 /M SM 12 | = 1.05 +0.14 −0.13 . On the other hand, the Summer 2016 result [84] by UTfit, which includes the result of Ref. [86], constrains NP with a better precision: |M 12 /M SM 12 | = 1.070±0.088. As we assume g L bs to be real, the Z contribution always enhances |M 12 |. If one takes these uncertainties to be Gaussian, these results imply |M 12 /M SM 12 | < 1.32 [83] or 1.25 [84] at 2σ for the CKMfitter and UTfit results, respectively. In this paper, we explore NP contributions to |M 12 /M SM 12 | at the level of up to 30%. Excluding larger contributions leads to the gray-shaded regions in Fig. 1. Future improvements in lattice calculations and measurements of the CKM parameters would tighten the constraint [87]. In Fig. 1, we illustrate the impact of possible future improvements by the vertical dotted lines for deviation of M 12 from SM by 15% or 5%.
We note that while the CKMfitter [83] and UTfit [84] results are tolerant to a NP contribution that enhances |M 12 |, there are studies that find the SM prediction of ∆m B 0 s = 2|M 12 | to be larger than the measured value, slightly favoring NP that reduces |M 12 |. In particular, a recent study [77], which adopts the 2017 FLAG result [85], finds the SM prediction to be 1.8σ above the measured value. Their result can be read as |M 12 /M SM 12 | 0.89 ± 0.06, which allows an enhancement by NP only up to ∼ 1% at 2σ. The rather small uncertainty is in part due to a smaller uncertainty of 2.1% assigned to the CKM factor. Addressing the discrepancy among the theoretical calculations is beyond the scope of this paper. Instead, the 1% vertical dotted lines are also shown in Fig. 1 for illustrating the impact of the result by Ref. [77].
The nonzero LH Z µµ coupling implies also the existence of a Z νν coupling, due to SU(2) L . Therefore, constraints are also set by B → K ( * ) νν. The effective Hamiltonian for b → sνν is [88] H where are lepton-flavor dependent Wilson coefficients. The SM contribution is lepton-flavor universal and is given by [89]. The Z contribution is given by Normalizing B → K ( * ) νν branching ratios by the SM ones and defining R ν Combining the charged and neutral B meson decays, the tightest limits [90] are set by Belle [91], who find at the 90% confidence level (C.L.). The tighter constraint comes from R ν K * , and is shown by the cyan lines in Fig. 1. The allowed region fully contains the blue hyperbolae favored by b → s The Z bs and Z bb couplings induce Z production at LHC via bs → Z and bb → Z . Hence, with Z → µ + µ − , these couplings are constrained by dimuon resonance searches at LHC.
We use the results from ATLAS, performed with 36.1 fb −1 at 13 TeV [92] and extract [93] the 95% credibility level limit: σ(pp → Z + X)B(Z → µ + µ − ) < 42 fb (9 fb) for m Z = 200 GeV (500 GeV). We calculate the Z production cross section at leading order (LO) using MadGraph5 aMC@NLO [94] with the NN23LO1 parton distribution function (PDF) set [95]. As the AT-LAS search does not veto additional activity in the event, we include also the processes gs → Z b, gb → Z b and gg → Z bb, Z bs in the cross-section calculation. We defer the more detailed discussion about Z production at LHC to Section IV. From the cross sections and the AT-LAS limits, we find (13) in scenario (i). In scenario (ii), where the Z couples to the LH muon current, the limits are weakened by an overall factor of 2/ √ 3 on the right-hand side due to the change in B(Z → µ + µ − ).
As long as |g L bb | |g L bs |, which is our scenario of interest, these limits are significantly weaker than those from the B s −B s mixing. Hence, they are not shown in Fig. 1. For flavor universal models, Z masses below m Z 3−4.5 TeV are ruled out by Ref. [92] with 95% CL. However, a very weakly coupled and flavor non-universal light Z such as described by Eq.(1), escapes the detection and could emerge in the future runs of LHC.
Muon pair production in the scattering of a muon neutrino and a nucleus N , known as neutrino trident production, tightly constrains Z µµ and Z ν µ ν µ couplings [96]. The ratio between the total ν µ N → ν µ N µ + µ − cross section and its SM prediction is given by [31] in scenario (i) with m Z 10 GeV. The measurement [97] by the CCFR collaboration is in a good agreement with SM, and implies σ/σ SM = 0.82 ± 0.28. We show the resulting 2σ upper limits on |g V µµ | by the horizontal solid red lines in Fig. 1.
The couplings of the Z boson with the muon and muon neutrino are modified by Z -loop contributions, which can lead to violation of the lepton-flavor universality in Z decays. In scenario (i), the vector and axial-vector Zµµ couplings relative to the SM-like Zee couplings are given by [31,50] where K(m 2 Z /m 2 Z ) is a loop function given in Ref. [98], and its real part is taken to match the convention of Ref. [99]. Here, the lepton-flavor universality in the SM case is exploited. Similarly, normalized Zνν couplings are given by where the factor of 1/3 effectively takes into account the fact that only Z → ν µνµ is affected by the Z among the three neutrino modes. The Z couplings were very precisely measured at SLC and LEP. Relevant results from the average of 14 electroweak measurements are g V e = −0.03816 ± 0.00047, g Ae = −0.50111 ± 0.00035, g V µ = −0.0367 ± 0.0023, g Aµ = −0.50120 ± 0.00054 and g V ν = g Aν = 0.5003 ± 0.0012 [99]. Of the four possible Z coupling ratios, we take only the one which is most sensitive to the effect of the Z , i.e. g Aµ /g Ae = 1.00018 ± 0.00128, where the uncertainties are added in quadrature. The resulting 2σ upper limits on |g V µµ | are shown by the horizontal red dashed lines in Fig. 1.
Nonzero values of g L bs or g L bb can alter the Zbb and Zss couplings at one loop. Taking the b and s quarks to be massless, we find that the Z loop with a nonzero g L bs modifies the LH Zbb and Zss couplings g Lb and g Ls relative to their SM values g SM Lb and g SM Ls in the same way as the LH Zµµ coupling, but with the replacement g L µµ → g L bs : The RH counterparts remain unchanged. The effect of the Z loop can be constrained by comparing the measured value g Lb = −0.4182 ± 0.0015 (g Ls = −0.423 ± 0.012) [99] to the corresponding SM prediction g SM [99]. Since g Lb is more precisely measured than g Ls , 2 we use g Lb to extract the limit on g L bs . Adding the errors in g Lb and g SM Lb in quadrature after symmetrizing the g SM Lb errors, we find the 2σ upper limit |g L bs | 0.34 (0.67) for m Z = 200 (500) GeV. These limits are much weaker than the ones obtained from the B s mixing, and we do not display them in Fig. 1. A similar conclusion can be made for the g L bb coupling as well. If both g L bs and g L bb are nonzero, an FCNC decay Z → bs, which is absent in the SM at tree-level, is induced by the one-loop Z contribution. A preliminary result by DELPHI [100] sets the 90% CL upper limit R b = q=d,s σ(e + e − → bq +bq)/σ(e + e − → hadrons) ≤ 2.6 × 10 −3 at the energy scale of the Z mass. Using B(Z → hadrons) 70% [80], one may rewrite the limit as q=d,s B(Z → bq +bq) 1.8 × 10 −3 . Since the Zloop-induced LH Zbs coupling is suppressed by the factor g L bs g L bb /(16π 2 ), the DELPHI limit is relevant only if both g L bs and g L bb are O(1) and the Z mass is not far from the Z mass. Since the B s -mixing constraint on g L bs is much tighter, and we concentrate on the case |g L bb | |g L bs |, the impact on B(Z → bs +bs) is generically far below the DELPHI limit in the scenarios considered in this paper.
The Z one-loop contribution to the muon anomalous magnetic moment a µ = (g µ − 2)/2 is [101] where scenario (i) and m Z m µ are assumed. The difference between the measured value of a µ and its SM prediction is (2.9 ± 0.9) × 10 −9 [102]. Assigning this difference to Eq. (18) yields the dark yellow 2σ regions in Fig. 1. Since the 2σ constraints from the neutrino trident cross section are tighter, the Z does not solve the tension in the muon g − 2. At the 3σ level, the g µ − 2 regions become compatible with the neutrino trident production constraints. Therefore, we ignore g µ − 2 in the rest of this paper.
We have discussed so far the constraints in scenario (i), summarized in Fig. 1 on the g V µµ vs. g L bs plane. From Fig. 1, we find the constraint on g L bs is very stringent due to the B s −B s mixing, while the constraint on g V µµ is much weaker. With the 30% NP effects allowed in the B s −B s mixing amplitude, the former constraint is |g L bs | 1.4 (3.6) × 10 −3 for m Z = 200 (500) GeV, imposing the lower limit on the Z µµ coupling |g V µµ | 0.031 (0.078) on the b → s hyperbolae. This validates the numerical values of the branching ratios in Eq. (7) to a good approximation, if |g L bb | is not too large compared to |g L bs |. The qualitative feature is the same in scenario (ii), where the Z µµ coupling is of LH, as the B s −B s mixing constraint does not change. Some of the other observables would give slightly different constraints on the Z µµ coupling, but the effect is minor to our study. The only notable difference from scenario (i) is the slight change in the value of B(Z → µ + µ − ) from Eq. (7). A nonzero axial-vector Z µµ coupling can make B s → µ + µ − relevant, but the latest LHCb result [103] does not exclude the allowed regions for the b → s anomalies.

IV. DISCOVERY AND IDENTIFICATION OF THE Z AT LHC
Having determined the constraints on the Z couplings, we proceed to study signatures for direct production of the on-shell Z in pp collisions with a center-of-mass energy of √ s = 14 TeV. The goal of this study is to ascertain the LHC potential for both discovery of the Z and determination of the flavor structure of its couplings. Therefore, motivated by the tensions in b → s , we fo-cus on the role of the Z bs coupling g L bs . If this is the dominant coupling to the quark sector, the Z will be primarily produced via the parton-level process bs → Z shown in Fig. 2. With the decay Z → µ + µ − , the Z may be discovered in the conventional dimuon resonance searches.
Such a discovery of a dimuon resonance, however, does not necessarily imply the existence of the Z bs coupling. In general, the process pp → Z + X may be facilitated by coupling to other quarks, particularly flavor-diagonal couplings. To test for dominance of the Z bs coupling, we propose to also search for pp → Z b + X (see Fig. 3 for typical parton-level processes). A Z bb coupling, predicted in many models (e.g., Refs [31,51,58,61]), motivated by the b → s anomalies, also contribute to pp → Z + X and pp → Z b + X via the parton level processes bb → Z and gb → Z b. Since a Z bb coupling also leads to pp → bbZ +X, via, for example, gg → Z bb, measuring the cross section for pp → Z bb + X may facilitate to discriminate the Z bb coupling from Z bs. Similarly, the production process pp → Z bs + X, occurring due to gg → Z bs, can in principle help probe the Z bs coupling.
We explore these signatures using the effective Lagrangian of Eq. (1) with scenario (i), i.e., with a vector Z µµ coupling. For each of the two Z mass values studied, we fix the couplings to the benchmark points shown by the red dots in Fig. 1: g L bs = 0.001, g V µµ = 0.04 (m Z = 200 GeV), g L bs = 0.0025, g V µµ = 0.1 (m Z = 500 GeV). (19) These values are selected such that g L bs leads to a 15% enhancement in the B s −B s -mixing amplitude M 12 . The value of g V µµ is then chosen so as to lie in the range given by Eqs. (5) and (3). We note that one may take a larger |g L bs | (with a smaller |g V µµ |), which would enlarge the pp → Z + X and pp → Z b + X cross sections by up to a factor of two, with the B s −B s -mixing constraint saturated at 2σ, i.e. a ∼30% enhancement in M 12 .
For the Z bb coupling, we study three cases for each benchmark point. The baseline case is g L bb = 0, which restricts assumptions about the Z couplings to the minimum needed to explain the b → s anomalies. In addition, we also explore the cases g L bb = g L bs and g L bb = 2g L bs , to study the impact of a nonzero g L bb . These choices of quark couplings satisfy the dimuon resonance search limits in Eq. (13) and maintain the Z branching ratios in Eq. (7).
In the following subsections, we mainly focus on the discovery potential of the Z in the production processes pp → Z + X, pp → Z b + X, and pp → Z bb + X, with the Z always decaying to µ + µ − . We use Monte Carlo event generator MadGraph5 aMC@NLO [94] to generate signal and background samples at LO with the NN23LO1 PDF set [95]. The effective Lagrangian of Eq. (1) is implemented in the FeynRules 2.0 [104] framework. The matrix elements for signal and background are generated with up to two additional jets and interfaced with PYTHIA 6.4 [105] for parton showering and hadronization. Matching is performed with the MLM prescription [106]. The generated events are passed into the Delphes 3.3.3 [107] fast detector simulation to incorporate detector effects based on ATLAS.
A. Observation of pp → Z + X Several SM processes constitute background for pp → Z + X, where we remind the reader that the Z decays into µ + µ − . The dominant background is due to the Drell-Yan (DY) events, pp → Z/γ * + X. The pp → tt events with semileptonic decay of both top quarks is the next largest background. Smaller backgrounds arise from pp → W t and V V , where V ≡ W, Z. Background may also arise from leptons produced in heavy-flavor decays or from jets faking leptons. These background sources are not well modeled by the simulation tools, and we ignore them, assuming that they can be reduced to subdominant level with lepton quality cuts without drastically impacting the results of our analysis.
We scale the LO cross sections obtained by Mad-Graph5 aMC@NLO as follows. The DY cross section is normalized to a NNLO QCD+NLO EW cross section by a factor of 1.27, obtained with FEWZ 3.1 [108] in the dimuon-invariant mass range m µµ > 106 GeV. We normalize the LO pp → tt and pp → W t cross sections to NNLO+NNLL cross sections by 1.84 [109] and 1.35 [110], respectively. The pp → W W , pp → W Z and pp → ZZ cross sections are normalized to NNLO QCD by 1.98 [111], 2.07 [112] and 1.74 [113], respectively. We do not apply correction factors to the signal cross sections throughout this paper.
We select events that contain at least two oppositely charged muons. The transverse momentum of each muon is required to satisfy p T µ > 50 GeV, and its pseudorapidity must be in the range |η µ | < 2.5. The two muons must satisfy ∆R µµ = ∆φ 2 µµ + ∆η 2 µµ > 0.4, where ∆φ µµ and ∆η µµ are the separations in azimuthal angle and pseudorapidity between the muons. Finally, we require the invariant mass of the two muons to satisfy |m µµ − m Z | < m cut , where m cut = 4 GeV and 16 GeV for m Z = 200 GeV and m Z = 500 GeV, respectively. These values are chosen so as to maximize the naive local significance of the no-signal hypothesis, S l = N S / √ N B , where N S and N B are the expected signal and background yields.
The invariant mass cut |m µµ −m Z | < m cut is not realistic for a discovery scenario, in which one does not know the true mass m Z . However, the value of S l thus obtained is a rough estimate of the one that will be found by the more sophisticated analysis that will eventually be performed with the full LHC data. One is actually interested in the global significance S g , which accounts for the probability to obtain the given value of S l anywhere in the dimuon-invariant mass range. Rigorous methods for estimating S g exist [114]. However, at this level of approximation, it is sufficient to use the crude estimate where P g and P l are the χ 2 probabilities corresponding to S g and S l , respectively, and m range is the size of the range of m µµ values explored in the analysis. Since cross sections drop to negligible levels at high m µµ , it is reasonable to take m µµ ∼ 2-3 TeV for this estimate.
The cross sections for signal and backgrounds after the cuts are listed in Table I for the benchmark points defined in Eq. (19) with the three choices of g L bb . The corresponding values of the local significance S l for an integrated luminosity L = 3000 fb −1 are summarized in Table II. Inserting the values of S l into Eq. (20), we conclude that the global significance will likely be greater than 5σ for the case g L bb = 2g L bs , allowing separate discovery by AT-LAS and CMS. For |g L bb | ≤ |g L bs |, the global significance will be under 5σ. Whether the 5σ mark will be passed by combining ATLAS and CMS results is beyond the precision of our rough estimate.
A larger |g L bs | can enhance the significance in each benchmark scenario. For the scenario of m Z = 200 (500) GeV with g L bb = 0, taking |g L bs | = 0.0013 (0.0034) with |g V µµ | = 0.031 (0.074) pushes the local significance slightly above 6σ, which may imply a global significance of 5σ. In this case, the B s −B s mixing amplitude |M 12 | is also enhanced from the SM one by 25% (28%), but still within the nominal 2σ allowed range, as discussed in Sec. III. NB for discovery of the process pp → Z + X with Z → µ + µ − , with an integrated luminosity of 3000 fb −1 , given the signal and background cross sections shown in Table I. m Z (GeV) Local significance g L bb = 0 g L bb = g L bs g L bb = 2g L bs 200 3.0 3.9 6.6 500 3.0 3.4 7. 2   TABLE IV. Local significance S l = NS/ √ NB for discovery of the process pp → Z b + X with Z → µ + µ − , with an integrated luminosity of 3000 fb −1 , given the signal and background cross sections shown in Table III.
The main SM background for pp → Z b + X is pp → tt. The second-largest background is Drell-Yan with at least one additional b jet, labeled as DY + b. Smaller contributions arise from DY + c, pp → W t, pp → V V , and DY + j, where j stands for a jet from a gluon or a u, d, or s quark. We normalize the cross sections for DY +b, DY + c, and DY + j to NNLO QCD by 1.83 [115]. The correction factors for pp → tt, pp → W t, and pp → V V are taken to be the same as in Section IV A.
We select simulated events that contain at least two opposite-charge muons. The muons are required to be in the pseudorapidity range |η µ | < 2.5, have minimal transverse momenta of p T µ > 50 (60) GeV for m Z = 200 (500) GeV, and be separated by ∆R µµ > 0.4. Jets are reconstructed using the anti-k T algorithm with radius parameter R = 0.5. It is assumed that a b-tagging algorithm reduces the efficiency for c jets and light jets by factors of 5 and 137, respectively [116]. Its efficiency for b jets is calculated in Delphes, accounting for the p T and η dependence. The leading b jet is required to have transverse momentum p T b > 30 GeV with |η b | < 2.5, and its separation from each of the two leading muons must sat-isfy ∆R bµ > 0.4. We reject events that have a second b-tagged jet with p T b > 30 GeV, slightly increasing the local significance. The missing transverse energy must be less than 40 GeV, in order to reduce pp → tt and pp → W t backgrounds. Finally, we apply the optimized dimuon-invariant mass cut |m µµ −m Z | < 5 (15) GeV for m Z = 200 (500) GeV.
The resulting cross sections are shown in Table III, and the corresponding local signal significances with 3000 fb −1 are summarized in Table IV. The local significances are slightly smaller than the corresponding ones in Table II, except in the case of m Z = 500 GeV and g bb = 2g L bs . Thus, we conclude that, like pp → Z +X, the process pp → Z b + X is likely be discovered at S g > 5 if g L bb ≥ 2g L bs in our benchmark points. By scaling the values in Table IV, we observe that a local significance of 6σ can be attained with |g L bs | 0.0014 (0.0036) for m Z = 200 (500) GeV, even if g L bb = 0, at the cost of a ∼ 30% enhancement in |M 12 /M SM 12 |.
C. Observation of pp → Z bb + X and pp → Z bs + X The dominant SM backgrounds for pp → Z bb + X are pp → tt, pp → W t and DY + b or c jets. pp → m Z (GeV) σ signal (fb) σ background (fb) g L bb = 0 g L bb = g L bs g L bb = 2g L bs DY + h. Local significance S l = NS/ √ NB for discovery of the process pp → Z bb + X with Z → µ + µ − , with an integrated luminosity of 3000 fb −1 , given the signal and background cross sections shown in Table V. V V gives a negligible contribution. We adopt the same correction factors for the background cross sections and follow the same event selection criteria as in Section IV C, and in addition require the subleading b jet to have p T b > 30 GeV, |η b | < 2.5, and to be separated from the leading b jet and each of the two leading muons by ∆R > 0.4, for both the Z masses.
The resulting cross sections are shown in Table V. The g L bb = 0 cases have tiny but nonzero cross sections, due to production via bs → Z g * , with g * → bb. By contrast, a nonzero g L bb induces the less suppressed process gg → Z bb. The corresponding local signal significances for 3000 fb −1 are given in Table VI. As the local significances are much less than 1, we conclude that observation of this process is not possible at LHC within the range of couplings explored here.
In general, |g L bb | may take a larger value, up to the limit of Eq. (13), namely, |g L bb | = 0.007 (0.024) for m Z = 200 (500) GeV with g L bs ∼ 0. With these values, we estimate the cross section of pp → Z bb + X to be 0.098 fb (0.055 fb) after the event selection cuts. This corresponds to a local significance of around 3.6σ (4.1σ) for an integrated luminosity of 3000 fb −1 . Thus, a global significance of 5σ is not likely. We note, however, that since we used the same QCD correction factors for the background cross sections as in the pp → Z b + X case, there is a greater uncertainty on these cross sections.
Generally, the cross section for pp → Z bb + X is strongly suppressed by the 3-body phase space. Since the same suppression applies for pp → Z bs + X, one expects the cross section for this process to be small as well. Moreover, the process pp → Z bs + X would also suffer from light-jet backgrounds which make the discovery not possible, given the B s −B s mixing constraint.

V. IMPACT OF THE RIGHT-HANDED Z bs COUPLING
In this section, we study an impact of a tiny but nonzero RH Z bs coupling by adding the following terms to the effective Lagrangian in Eq. (1): The resulting additional contributions to b → sµ + µ − are described by effective operators as in Eq. (2), with P L replaced by P R , and with C NP 9 and C NP 10 replaced by the Wilson coefficients There is no significant indication for nonzero C 9,10 in the majority of the b → s global fit analyses. However, even a tiny g R bs can drastically affect the B s −B s mixing, which is now given by [31] M 12 calculated with the hadronic matrix elements in Ref. [117]. The large negative coefficient of the g L bs g R bs term, which is partly due to renormalization group effects [79], means that a small value of g R bs allows for a large g L bs , due to cancellation between the terms. In Fig. 4 we show the B s −B s -mixing constraint on g R bs vs. g L bs , when |M 12 | is allowed to change by up to 30% of its SM value. Reducing the allowed NP contribution to |M 12 |, say, to 15%, would narrow the width of the tilted-crossshaped allowed region in Fig. 4. However, it would not change the conclusion, namely, that a large value of |g L bs | is allowed. What now becomes the most significant limit on g L bs is the ATLAS dimuon resonance search [92], shown by the solid lines, assuming Eq. (7).
The cancellation in M 12 requires g L bs g R bs > 0. This implies (Re C NP 9 )(Re C 9 ) > 0, contrary to the best-fit values for C NP 9 and C 9 , e.g., Refs. [15,16]. However, the cancellation requires only g R bs ∼ 0.1g L bs or C 9 ∼ 0.1C NP 9 . While the fits favor C NP 9 ∼ −1 and are consistent with , when allowing the Bs −Bs-mixing amplitude to deviate by up to 30% from its SM prediction. The gray regions are excluded. The solid lines show the exclusion by the ATLAS dimuon resonance search [92] for B(Z → µ + µ − ) 2/3. The red dots show our benchmark points. C 9 = 0, they cannot exclude a small negative C 9 . Indeed, the point (C NP 9 , C 9 ) = (−1.11, −0.1) is at the border of the 1σ ellipse in Ref. [15] with the assumption C NP 10 = C 10 = 0.
To illustrate the impact of such a possibly large Z bs coupling on the Z discovery potential, we consider scenario (i) with the following benchmark points for m Z = 200 and 500 GeV respectively (corresponding to the dots in Fig. 4): g L bs = 0.0032, g R bs = 0.00036, g V µµ = 0.012, g L bs = 0.0075, g R bs = 0.00083, g V µµ = 0.032.
Both points correspond to (C NP 9 , C 9 ) = (−1.11, −0.1). The effects of the tiny g R bs on other constraints, e.g. B → K ( * ) νν, are negligible. Taking g L bb = 0, we find that the pp → Z + X and pp → Z b + X cross sections are highly enhanced compared to the ones in Sec. IV, while the relatively large g L bs values lead to a non-negligible Z → bs branching ratio, slightly reducing B(Z → µ + µ − ) from ∼ 66% to ∼ 52% and ∼ 55% respectively for 200 and 500 GeV Z .
9.3 7.6 500 7.7 6.9 TABLE VII. Local significance NS/ √ NB for discovery of the process pp → Z + X → µ + µ − + X and pp → Z b + X → µ + µ − b + X with the integrated luminosity of 300 fb −1 , obtained by rescaling the results of Table II and IV to the benchmark points defined in Eq. (24) for g L bb = 0.
Taking account of these effects and rescaling the significances obtained in Sec. IV, we show in Table VII the local significances for discovery of pp → Z + X and pp → Z b+X with the integrated luminosity of 300 fb −1 , for the two benchmark points defined above. The results suggest that both processes can be discovered with the Run-3 dataset, even if the Z bb coupling vanishes. A larger g L bs in this scenario enhances the cross section for pp → Z bs + X process, but the discovery would still be beyond the reach of the HL-LHC.

VI. SUMMARY AND DISCUSSIONS
Observed tensions in b → s measurements can be explained by a new Z boson that couples to the lefthanded b → s current as well as to muons. In this paper, we have studied the collider phenomenology of such a Z based on an effective model introduced in Eq. (1). For this purpose, we first estimated phenomenological constraints on the Z bs and Z µµ couplings for the representative masses m Z = 200 and 500 GeV. The most important constraint is the B s mixing, which tightly constrains the LH Z bs coupling g L bs .
For fixed values of g L bs and m Z , the allowed coupling to muons is determined by global fits to b → s data, up to the value allowed by the constraint from the neutrino trident production, where the Z ν µ ν µ coupling is related to Z µµ coupling by the SU(2) L symmetry. We also introduced the Z bb coupling g L bb , which is mildly constrained by the dimuon resonance search at LHC. The resulting couplings are such that the Z decays mostly to µ + µ − and ν µνµ , with the two branching ratio values mildly depending on whether the muon coupling is vector-like or left-handed.
Given the coupling constraints, we explored the capability of the 14 TeV LHC to discover the Z and to determine the flavor structure of its couplings. For the sake of this dual goal, we studied the two processes pp → Z + X → µ + µ − + X and pp → Z b + X → µ + µ − b + X, where the former may be induced by bs → Z and/or bb → Z and the latter by gs → Z b and/or gb → Z b. We considered two representative Z masses of 200 and 500 GeV with three scenarios for the Z bb coupling: g L bb = 0, g L bs or 2g L bs . For g L bb = 0, we found that discovery of pp → Z + X and pp → Z b + X (with about 5σ global significance) with 3000 fb −1 data requires a large Z bs coupling, so that the B s mixing amplitude M 12 is enhanced by ∼ 30% or more relative to the SM expectation. This corresponds roughly to the 2σ upper limits of the global analyses [83,84] for the CKM parameters. With a nonzero g L bb , discovery in both the modes is possible without such a drastic effect on the B s mixing; in particular, for g L bb = 2g L bs , discovery is possible with a ∼ 15% deviation in the M 12 . For further discrimination between the Z bs and Z bb couplings, we also studied the process pp → Z bb + X → µ + µ − bb + X, predominantly arising from Z bb coupling. However, we found it to be not promising even with 3000 fb −1 integrated luminosity, due primarily to three-body phase-space suppression. The same conclusion applies to pp → Z bs + X, which gives direct access to the Z bs coupling.
The discovery potential of the Z is rather limited due to the B s mixing constraint. The B s mixing constraint, however, is only indirect and is susceptible to the details of the UV completion of the effective model. In particular, we illustrated that the existence of a tiny but nonzero right-handed Z bs coupling g R bs accommodates a large LH Z bs coupling due to the cancellation in the B s mixing amplitude, without conflicting with the b → s global fits. In this case we found that discovery in both the pp → Z + X and pp → Z b + X processes may occur even with O(100) fb −1 integrated luminosity.
Comments on the subtlety of the implementation of the B s mixing constraint are in order (see also Sec. III). As mentioned above, a 30% enhancement in the B s mixing amplitude M 12 by NP roughly corresponds to the 2σ upper limits by the latest global analyses of CKMfitter [83] and UTfit [84]. This may look rather tolerant, in view of the recent progress [86] in the estimation of the hadronic matrix element by lattice, which lead to ∼ 6% uncertainty in M SM 12 [85]. This is because the central values of |M 12 /M SM 12 | are greater than unity in Ref. [83] and Ref. [84], while the Z contribution always enhances |M 12 | relative to the SM value under the assumption of a real-valued g L bs with g R bs = 0. On the other hand, a recent study [77] finds the SM prediction of ∆m B 0 s = 2|M 12 | to be 1.8σ above the measured value, favoring |M 12 /M SM 12 | smaller than unity. This is opposite to the results by CKMfitter [83] and UTfit [84], although both UTfit (Summer 2016 result) and Ref. [77] take into account the recent lattice result [86]. If the result of Ref. [77] is the case, the Z contribution may enhance |M 12 | only up to ∼ 1% so that g L bs is strongly constrained. In this case the estimated signal significances at LHC would shrink down to insignificant values for |g L bb | |g L bs |, unless a tiny RH Z bs coupling exists for the cancellation in M 12 and/or g L bs is close to pure imaginary so that it gives a negative contribution in M 12 . The latter implies a nearly imaginary C NP 9 and would need a dedicated global analysis of b → s observables, as discussed in Ref. [77]. In any case, a consensus among the different groups seems to be still missing for the prediction of M 12 in the SM, and a better understanding would be required for its calculation. At the same time, improvements in lattice calculations and determinations of CKM parameters will also facilitate a more precise SM prediction for M 12 .
Although we considered Z bb and Z bs couplings to be the only couplings to the quark sector, the Z may also couple to other quarks in general. For instance, if a non zero Z cc coupling exists, the process pp → Z c + X can be induced at LHC. Such a process can mimic the pp → Z b + X signature if the final state c-jet gets misidentified as b-jet. This possibility can not be excluded yet as pointed out in Ref. [118], where a procedure to disentangle pp → Z c + X and pp → Z b + X is discussed with the simultaneous application of both cand b-tagging. We also remark that our estimation of the signal significances ignored various experimental uncertainties and the QCD corrections to the signal cross sections.
For illustration, we focused on m Z = 200 and 500 GeV. In general, heavier Z are possible. However, due to the fall in the parton luminosity with the resonance mass, the achievable significances are lower than those of 200 GeV and 500 GeV in both the pp → Z and pp → Z b processes.
Our results illustrate three possible scenarios for the LHC discovery and identification of a Z that might be behind the b → s anomalies. The first one is the case with the minimal assumption, where the LH Z bs coupling is the only coupling to the quark sector. In this case, the discovery of the pp → Z + X → µ + µ − + X and pp → Z b + X → µ + µ − b + X processes may occur with the full HL-LHC data, but should be accompanied by a ∼ 30% or larger enhancement in the B s mixing, which can be tested following future improvements in the estimation of the B s mixing. The second one is the case with a tiny but nonzero RH Z bs coupling such that the B s mixing remains SM-like due to the cancellation of the Z effects. In this case, the discovery of the two modes may occur with Run-3 data (or perhaps even Run-2 data); this scenario predicts a nonzero RH b → s current, with C 9 ∼ 0.1C NP 9 , which can be tested with improvements in b → s measurements by ATLAS, CMS, LHCb and Belle II. In particular, precise measurements of R K and R K * by LHCb with Run-2 or further dataset may pin down the chiral structure of the b → s current. The third scenario is the case with a flavor-conserving Z bb coupling much larger than Z bs. In this case, the two modes may be discovered with Run-3 data without a significant effect in the B s mixing and RH b → s current, but the role of the observed resonance in b → s is obscured.
Note added: While revising the manuscript we noticed that the CMS 13 TeV 36 fb −1 result [119] for a heavy resonance search in the dilepton final state is now available. We find that the extracted upper limits [120] on g L bb from Ref. [119] are comparable to those from AT-LAS [92] and do not change the conclusion of our results.