QCD next-to-leading-order predictions matched to parton showers for vector-like quark models

Vector-like quarks are featured by a wealth of beyond the Standard Model theories and are consequently an important goal of many LHC searches for new physics. Those searches, as well as most related phenomenological studies, however, rely on predictions evaluated at the leading-order accuracy in QCD and consider well-defined simplified benchmark scenarios. Adopting an effective bottom-up approach, we compute next-to-leading-order predictions for vector-like-quark pair production and single production in association with jets, with a weak or with a Higgs boson in a general new physics setup. We additionally compute vector-like-quark contributions to the production of a pair of Standard Model bosons at the same level of accuracy. For all processes under consideration, we focus both on total cross sections and on differential distributions, most these calculations being performed for the first time in our field. As a result, our work paves the way to precise extraction of experimental limits on vector-like quarks thanks to an accurate control of the shapes of the relevant observables and emphasise the extra handles that could be provided by novel vector-like-quark probes never envisaged so far.


Introduction
The Standard Model of particle physics is a successful theory of nature, although it exhibits many conceptual issues, like the hierarchy problem or the strong C P problem, and practical limitations like the absence of a viable candidate for explaining the dark matter pervading our universe. As a result, it is commonly acknowledged as an effective theory stemming from a more fundamental setup that has still to be observed and confirmed experimentally. This effective a e-mail: fuks@lpthe.jussieu.fr description has been recently strengthened with the discovery of a Higgs boson with properties very similar to those expected from the Standard Model in 2012. The null results of all collider searches for new particles predicted by most of beyond the Standard Model theories are, however, at the same time pushing the limits on the mass of these potential new particles to higher and higher energy scales.
Among the viable options for new physics, many extensions of the Standard Model predict the existence of additional quark species that should be observable during the next runs of the Large Hadron Collider (LHC) at CERN. One of the common feature of such theories consists of the predicted vector-like nature of the additional quarks, i.e. their left-handed and right-handed components lie in the same representation of the electroweak symmetry group. These quarks appear for instance in models with extra space-time dimensions that exhibit an extended gauge symmetry or a new strong dynamics giving rise to massive composite states [1][2][3][4][5]. As a result, vector-like quark searches play an important role in the ATLAS and CMS experimental program.
Current searches, relying on signatures induced by both the vector-like quark pair-production and single-production modes, impose strong constraints on the masses of the heavy quarks that are now bounded to be above about 750-1500 GeV [6][7][8][9][10][11][12], this wide range reflecting the wealth of options for describing how the new states decay into a pair of Standard Model particles. Care must, however, be taken with the interpretation of these limits as they are extracted once various simplifying assumptions are accounted for the new physics signal. Most of the bounds indeed assume that the quark partners decay, with a branching ratio of 100%, into third-generation top or bottom quarks, while more general situations where decays into second or first generation quarks are possible are less explored yet. Sizeable couplings to light quarks are still allowed by indirect constraints [13][14][15], which could have a severe impact on electroweak vectorlike quark processes at the LHC [16]. The latter, which induce vector-like quark decays, are driven by the couplings of the extra quarks to the weak and Higgs bosons. They admit a simple model-independent parameterisation [17] regardless of the representation of the new physics particles under the electroweak symmetry group, and this bottom-up approach is now adopted by the experimental collaborations. Most LHC searches for heavy quarks hence turn out to be agnostic of the ultraviolet completion of the model and can rather easily be reinterpreted in any framework, and options for combining different searches can also be considered.
Current vector-like quark searches, as well as most associated theoretical work, are, however, based on Monte Carlo simulations of the new physics signals where hard-scattering matrix-elements are evaluated at the leading-order (LO) accuracy in QCD. In addition, event samples featuring different final-state jet multiplicities are sometimes also merged in order to get a better control on the shapes of the key differential distributions. The formal precision of these calculations is nonetheless rather limited, which directly impacts the extraction of any limit on the properties of the new particle properties or the corresponding measurements in the case of a discovery. In this work, we build up a new procedure that allows us to make use of the MadGraph5_aMC@NLO framework [18] to compute total rates and differential distributions at the next-to-leading-order (NLO) accuracy in QCD for processes involving vector-like quarks. This more precisely concerns vector-like quark production on the one hand, and the production of Standard Model particles when vectorlike quark diagrams contribute, regardless of the strong or electroweak nature of the Born process.
Our methodology relies on the joint use of the Feyn-Rules [19] and NLOCT [20] packages, the latter making use of FeynArts [21], for automatically generating a UFO library [22] that contains all tree-level vertices and counterterms necessary for NLO QCD computations. This UFO library can then be further used by Mad-Graph5_aMC@NLO for event generation, at the LO and NLO accuracies in QCD as well as for loop-induced processes. Virtual loop contributions are numerically evaluated through the MadLoop module [23] and combined with the real-emission diagrams following the FKS subtraction method as implemented in MadFKS [24,25], and the matching to parton showers is finally achieved according to the MC@NLO prescription [26].
In Sect. 2, we detail how we have modified the modelindependent parameterisation of Ref. [17] to make it suitable for NLO calculations in QCD matched to parton showers for vector-like quark processes. We also define a series of benchmark scenarios for the phenomenological studies performed in Sect. 3. We investigate vector-like quark pair production and single production (in association with either a jet or a weak gauge or Higgs boson), as well as diboson production. We summarise our work in Sect. 4.

A model-independent parameterisation
for vector-like quark models

Model description
Vector-like quarks appear in many extensions of the Standard Model. They are usually included as fields lying in the fundamental representation of SU (3) c and carry colour charges similar to those of their Standard Model counterparts. However, they can lie in various representations of the weak interaction symmetry group SU (2) L and be assigned different hypercharge U (1) Y quantum numbers. Focusing on phenomenologically viable minimal models that comprise a single Standard Model Higgs field Φ, only weak triplets, doublets and singlets of vector-like quarks are allowed [27]. Consequently, the particle content of the theory can solely include four species of extra quarks, which we denote by X , T , B and Y , their respective electric charges being Q = 5/3, 2/3, −1/3 and −4/3. Although vector-like and Standard Model quarks with the same electric charge mix, the mixing pattern and the resulting phenomenology can be simplified when minimality requirements are imposed. Since we consider that the Higgs sector contains a single Standard-Model-like scalar Higgs field Φ, quark mixings are solely generated by its Yukawa interactions. The mass splitting between the vectorlike quarks is consequently also constrained to be small and connected to the Higgs vacuum expectation value v, so that the extra quarks will always directly decay into a gauge or Higgs boson and one of the Standard Model quarks [27]. A model-independent effective parameterisation apt to describe the phenomenology of this vector-like quark setup has been recently proposed [17], but it is not suitable for higher-order QCD calculations. The reason is that in the latter parameterisation, the strengths of the interactions of the vector-like and Standard Model quarks with a single Higgs boson depend on the masses of the model particles. As a result, the renormalisation of the quark masses and the one of the couplings are related, which prevents all ultraviolet divergences that arise at the NLO from cancelling. We therefore modify the modelling of Ref. [17] so that all the couplings of the vector-like quarks to a gauge or a Higgs boson are free parameters, and we have the following effective Lagrangian: being supplemented to the Standard Model Lagrangian L SM . The terms in the first and second lines consist of gaugeinvariant kinetic and mass terms for the vector-like quark fields (taken in the mass eigenbasis) after restricting the covariant derivatives to their QCD component, The coupling parameter g s denotes the strong coupling constant, G μ the gluon field and T (and f for further references) the fundamental representation matrices (the structure constants) of SU (3). Although the electroweak pieces of the covariant derivative could have been included, they have been omitted in order to simplify our model description since they are model-dependent and are expected to yield a negligible effect with respect to their strong interaction part. The third and fourth lines of Eq. (1) collects the effective interactions of the physical Higgs boson h with one Standard Model quark and its vector-like partner, generation indices being understood for clarity. As mentioned above, such interactions are yielded by the (flavour-changing) Yukawa couplings of the Higgs doublet Φ with the up-type (q u ), downtype (q d ) and vector-like quarks that induce a mixing of the Standard Model and the new physics quark sectors. The relevant elements of the mixing matrices have been included in the strengths of the effective interactionsκ.
In the last six lines of Eq. (1), we include the weak interactions of the Z -boson and W -boson with one Standard Model quark and one vector-like quark. In our conventions, we have factorised out the weak coupling g which represents the overall interaction strength, and the κ andκ parameters include the relevant elements of the quark mixing matrices, as for the Higgs interactions. Moreover, the c W parameter stands for the cosine of the electroweak mixing angle.
The L VLQ Lagrangian above is equivalent, at the tree level, to the one of Ref. [17] once we impose (3) In this notation, f stands for a generation index and Γ Q X denotes the kinematic factor of the partial decay width of the extra quark Q into a final state containing an X boson. The κ Q parameter encodes the magnitude of the coupling of the extra quark Q to the different electroweak bosons, while the ζ parameters refer to the mixing of the vector-like quarks with the Standard Model quarks, In this expression, we have represented the left-handed and right-handed 4 × 4 mixing matrices between the new quarks and the three Standard Model quarks by V L ,R . Finally, the ξ parameters of Eq. (3) determine the relative importance of the various decay modes of the vector-like quarks, their sum being equal to one. The calculation of differential and total cross sections for LHC processes at the NLO accuracy in QCD necessitates to evaluate, on the one hand, real-emission squared amplitudes and on the other hand, interferences of tree-level with virtual one-loop diagrams. The ultraviolet divergences that arise in the latter case are absorbed through the renormalisation of the fields and parameters appearing in L VLQ . This is achieved by replacing all fermionic and non-fermionic bare fields Ψ and Φ and bare parameters y by the corresponding renormalised quantities, where we truncate the renormalisation constants δ Z and δy at the first order in the strong coupling α s = g 2 s /(4π). While the wave-function renormalisation constants of the Standard Model quarks are not modified by the presence of the vectorlike quarks, the one of the gluon field is given, when the on-shell renormalisation scheme is adopted and when we include n f = 5 massless flavours of quarks, by where the B 0,1 and B 0,1 functions stand for standard twopoint Passarino-Veltman integrals and their derivatives [28]. The left-handed and right-handed wave function renormalisation constants δ Z L ,R Q and the mass renormalisation constants δm Q of a vector-like quark Q (with Q = T , B, X , Y ) are similar to the top-quark ones and read In these two expressions, C F = (n 2 c − 1)/(2n c ) is the quadratic Casimir invariant associated with the fundamental representation of SU (3), with n c = 3. Finally, in order to fix the renormalisation group running of α s so that it originates from gluons and the n f = 5 active light quark flavours, we renormalise α s by subtracting, from the gluon self-energy, the contributions of all massive particles evaluated at zeromomentum transfer, where the ultraviolet-divergent pieces are written in terms of 1 = 1 − γ E + log 4π with γ E being the Euler-Mascheroni constant and being connected to the number of space-time dimensions D = 4 − 2 . The first term in the right-hand side of Eq. (8) results from the Standard Model massless parton contributions, while the second term is connected to the massive states, namely the top quark and the four considered vector-like quark species. In Sect. 3, we will compute predictions at the NLO accuracy in QCD for processes involving vector-like quarks. We will rely on a numerical evaluation of the loop integrals in four dimensions, which necessitates the calculation of rational terms associated with the -dimensional pieces of the loop integrals. There exist two sets of such rational terms that are, respectively, connected to the loop-integral denominators (R 1 ) and numerators (R 2 ). While the former are universal, the latter are model-dependent and can be seen as a finite number of counterterm Feynman rules derived from the bare Lagrangian [29]. Starting from the L VLQ Lagrangian of Eq. (1), several R 2 counterterms with external gauge bosons are modified with respect to the Standard Model case, with In addition, the new R 2 counterterms involving external vector-like quarks are given by In the conventions of Eqs. (9) and (11), c i , μ i , and p i indicate the colour index (which can be associated either with the adjoint or the fundamental representation of SU (3)), the Lorentz index, and the four-momentum of the ith particle incoming to the R ...i... 2 vertex, respectively. Moreover, an explicit summation upon q, Q and f implies a summation over all quark species, the extra quark species and the Standard Model quark species, respectively.
In the phenomenological study undertaken in Sect. 3, the virtual one-loop contributions to the NLO predictions are evaluated with the MadLoop module [23], and then combined with the real contributions by means of the FKS subtraction method [24] as implemented in the MadFKS package [25]. Both MadLoop and MadFKS being part of Mad-Graph5_aMC@NLO [18], the entire calculation is entirely automated from the knowledge of the bare Lagrangian of Eq. (1) and the specification of the process of interest [30]. Technically, the translation of the model Lagrangian into a UFO library [22] that contains ultraviolet and R 2 counterterms and that could be used by MadGraph5_aMC@NLO is automatically performed with the FeynRules [19] and NLOCT [20] packages, the latter program taking care of the calculation of the one-loop ingredients of the model files. The corresponding FeynRules and UFO models have been made publicly available and can be downloaded from the webpage http://feynrules.irmp.ucl.ac.be/wiki/NLOModels.

Benchmark scenarios
Throughout our phenomenological analysis, we adopt several series of benchmark scenarios in which one single vectorlike quark is light enough so that it could be reachable at the LHC. Moreover, for the sake of simplicity, we enforce its decay to proceed via a single channel. We denote each class of scenarios by the acronym QVi where the symbol Q can be either T , B, X or Y and refers to the nature of the relevant extra quark, the symbol V refers to the nature of the gauge boson which the vector-like quark Q decays into and i is a generation number related to the family which the quark Q mixes with. For instance, the scenario TW2 would correspond to a setup in which the Standard Model is supplemented by an extra up-type quark T that decays into a final state made of a W -boson and a strange quark with a branching ratio equal to 1.
These types of scenarios are motivated by several considerations. The mixings of the extra quark with the Standard Model sector are severely constrained by flavour-changing neutral current probes [17,31,32], LEP data [33,34] and atomic parity violation measurements [35]. Taking the parameterisation of the κ,κ andκ parameters of Eq. (3), sizeable mixings with all three generations are only allowed when the κ Q parameters are below 10 −2 -10 −3 [17]. Those bounds can, however, be relaxed when the mixing pattern is restricted to involve one or two quark generations. In our study, we enforce the vector-like quark mixing to only deal with one specific generation of Standard Model quarks, and we fix the values of the κ Q parameters to their current experimental limits of 0.07, 0.2 and 0.1 for mixings involving the first, second and third generation, respectively.

LHC phenomenology
In this section, we compute total cross sections and differential distributions both at the LO and NLO accuracy for several processes involving vector-like quarks. We study the genuine effects of the NLO corrections, as well as the induced reduction of the theoretical uncertainties. We then investigate the effects of matching fixed-order calculations to parton showers, both at LO and NLO. In Sects. 3.1 and 3.2, we, respectively, focus on vector-like quark pair and single For the considered processes, the central (total and differential) cross-section values are computed after setting the renormalisation and factorisation scales to the average transverse mass of the final-state particles and by using the NLO set of the NNPDF 3.0 parton density functions (PDF) [36] accessed via the LHAPDF 6 library [37]. Scale uncertainties are derived by varying both scales independently by a factor of two up and down, and the PDF uncertainties are extracted following the recommendations of Ref. [38], and both contributions to the theoretical uncertainties are added in quadrature.
3.1 Vector-like quark pair production at the LHC Vector-like quark pair production is in general dominated by QCD contributions, which has the advantage to be independent of the model details. Model-dependent electroweak diagrams induced by the last four lines of the Lagrangian of Eq. (1) may, however, be non-negligible, in particular when the final state of interest can be produced from the scattering of one or two valence quarks. We focus on the production of a pair of vector-like quarks including the cases where they have the same electric charge, with Q being either T , B, X or Y . While the first of these three subprocesses receives both strong (diagrams of the first line of Fig. 1) and electroweak (first diagram of the second line of Fig. 1) contributions, the latter two subprocesses can only be mediated by the t-channel exchange of a weak or Higgs boson (last two diagrams of the second line of Fig. 1). NLO corrections to the strong contributions to the production of a pair of vector-like quarks (the diagrams of the first line of Fig. 1) can be automatically calculated within the MadGraph5_aMC@NLO framework. Thanks to the upgrading of the model to NLO as explained in the previous section, it is now sufficient to type in the program shell, taking the example of TT production, import model VLQ_NLO_UFO generate p p > tp tp˜ [QCD] output launch With this set of instructions, we first import the UFO model associated with the Lagrangian of Eq. (1) and then start the calculation of the cross section and the generation of Monte Carlo events, at the NLO accuracy in QCD, for the production of a pair of TT quark and antiquark (whose UFO names are tp and tp˜). Other vector-like quark processes can be obtained by replacing the tp symbol by bp (for a B quark), x (for an X quark) and y (for a Y quark). LO event generation can be achieved in the same way once the [QCD] tag is omitted. We recall that the syntax is case insensitive.
For the electroweak channels (the diagrams of the second line of Fig. 1), the command to be typed in reads generate p p > tp tp˜QCD=0 [QCD] still for the same example of TT production. Other channels (including the production of a pair of heavy quarks carrying the same electric charge) can be addressed similarly. Care must, however, be taken as mixed electroweak and QCD loops appear at NLO. In our approach, we focus on NLO calculations in QCD, and not in QED or in the context of the electroweak theory. As a result, MadGraph5_aMC@NLO automatically discards Feynman diagrams tagged as an electroweak correction to a QCD graph, although in our case all diagrams should be kept for a proper cancellation of all divergences. One possible way to cure this issue would be to include in the UFO model library all ultraviolet and R 2 counterterms that would be necessary for undertaking mixed NLO calculations in QCD and QED and to implement in MadFKS the necessary subtraction terms. This, however, goes beyond the scope of this work, so that in order to maintain automation from the user standpoint, we have instead released a public script that should be called for diagram generation in order to prevent MadGraph5_aMC@NLO from discarding any loop diagram that would be tagged as an electroweak correction to a QCD Born diagram. More precisely, the script allows for the inclusion of all box diagrams containing at least two strong interaction vertices, but it removes weakboson or Higgs-boson loop contributions that consist of an electroweak correction to a QCD Born process. Additionally, diagrams exhibiting the t-channel exchange of a weak or Higgs boson but with an additional gluon are kept. Not using the script instead yields the removal of several necessary box diagram contributions.
Finally, QCD and electroweak diagrams can interfere. Because of the mixing of QCD and electroweak interaction orders at the one-loop level and the missing counterterms and subtraction terms that have been mentioned above, MadGraph5_aMC@NLO cannot currently be used for the calculation of these interferences beyond the LO accuracy. We therefore rely on LO simulations and reweight instead the results to include a K -factor approximating the effect of the QCD corrections, denoted by K (int) NLOQCD , taken as the ratio of the NLO to the LO (differential) results. We choose this K (int) NLOQCD factor to be the geometrical average of the K (QCD) NLOQCD and K (EW) NLOQCD K -factors obtained in the context of the QCD and electroweak diagrams taken independently, where our notations explicitly indicate the nature of the corresponding underlying Born process. We have additionally checked that, for the central value, our procedure yields numerical differences that are of at most 2%. The related MadGraph5_aMC@NLO command allowing for generating interference events is, again for the example of TT production, which is a standard command for LO event generation in MadGraph5_aMC@NLO.
All MadGraph5_aMC@NLO scripts necessary for differential and total cross-section calculations and event generation are available from the webpage http://feynrules. irmp.ucl.ac.be/wiki/NLOModels, together with the UFO and FeynRules models.
In Fig. 2 and Table 1, we present LO and NLO total cross sections for the production of a pair of vector-like T quarks for the different scenarios introduced in Sect. 2.2, and we depict the dependence of the cross sections on the vector-like quark mass. We first focus on the pure QCD contribution that is independent of the vector-like quark nature (diagrams of the upper line of Fig. 1). The genuine NLO contributions are found to be important as they first induce a shift in the cross section of about 50% within the entire probed m Q range and next reduce the dependence of the rate on the unphysical factorisation and renormalisation scales. The uncertainty band is indeed significantly reduced when NLO effects are accounted for, the scale dependence being reduced to the level of about ±10% over the entire mass range. At the NLO, the scale dependence induced by the virtual contributions indeed partially compensates for the one stemming from the Born and real-emission diagrams. Our results for the pure QCD case (strong production Born diagrams only) agree with the literature, and we recall that this corresponds to the production of a top-quark pair with a different top-quark mass [39]. It is in contrast the first calculation including the impact of the electroweak diagrams at the NLO accuracy in QCD.
In the considered scenarios, electroweak contributions to the production of a pair of vector-like quarks possibly carrying the same electric charge are found to be important only for the TH1 scenario, the associated results being shown on Fig. 2 as dashed and dotted bands. These bands are, respectively, related to the production of a pair of vector-like quarkantiquark (dashed) and to the production of a pair of heavy quark regardless of their electric charge (dotted). In this last case, the contributions from the three processes of Eq. (12) are summed over. The TH1 scenario features an extra quark that mixes with the first generation of Standard Model quarks so that parton density effects could lead to an enhancement of the production rate due to quark-antiquark, quark-quark and antiquark-antiquark (electroweak) scattering diagrams involving one or two initial valence quarks. This is particularly pronounced for setups featuring heavy vector-like quarks that require one to probe large Bjorken-x phase-space regions. As a consequence, the central cross-section values and the scale and parton density uncertainties are different from the pure QCD context since non-QCD diagrams (featuring a different initial state) dominate. This is illustrated in Table 1 for a few mass choices. Whereas the production rates are always larger, the uncertainties can be either smaller or larger than in the QCD case. We observe a huge gain in cross section for TH1 scenarios with a very heavy extra quark. This stems from Eq. (3) that shows that the coupling of the extra quark to the Higgs boson and a lighter quark has been taken proportional to the vectorlike quark mass, and is thus enhanced for large values of m Q . In principle, such a coupling should also be proportional to the related mixing matrix element ζ that compensates this enhancement, as shown in Eq. (3) and in Eq. (4). Setting ζ = 1, this feature is translated into the adopted value for the κ Q parameters.
Similar properties can be found in the context of the production of B, X and Y vector-like quarks, as illustrated in Appendix A.
Accurate differential distributions are often helpful for setting more precise exclusion limits and refine the experimental search strategies. Our implementation in the Mad-Graph5_aMC@NLO platform can be used to this aim, and we present in Fig. 3 several differential distributions in several observables, including NLO and parton-shower effects. We have chosen the BH2 class of benchmark scenarios with a vector-like quark mass set either to 500 GeV (upper series of curves on the figure) or 1500 GeV (lower series of curves on the figure). In our calculations, we have made use of the MadSpin [40] and MadWidth [41] programs for automatically taking care of the heavy quark decays in a way in which both off-shell and spin correlation effects are retained, matched the fixed-order calculation with the parton-shower and hadronisation infrastructure as modelled by the Pythia 8 package [42], and we have reconstructed all final-state jets by means of the anti-k T algorithm [43] with a radius parameter set to 0.5 as implemented in FastJet [44]. As shown on Figs. 3 and 4, our predictions confirm the total cross-section findings of Table 4 (see Appendix A). The contributions of the electroweak diagrams are, respectively, negligible and significant for light and heavy vector-like quarks. Moreover, the parton density uncertainties dominate for setups exhibiting a large M B value, rendering the theoretical predictions barely reliable.
In the two upper panels of Fig. 3, we study the properties of the first produced B-quark and show its transverse momentum p T (B 1 ) and pseudorapidity η(B 1 ) distributions. In the case of a light B-quark, the parton-shower effects (green and red solid lines) slightly affect the shapes of the spectra predicted by the fixed-order calculations (blue and brown solid lines), both at the LO and NLO accuracies. For heavier vector-like quarks, slight modifications can be observed in the small p T and large |η| region although the accurate modelling of the first extra jet does not yield any impact at the level of the individual B-quarks. These differences are largely covered by the theoretical uncertainties stemming from the poor knowledge of the parton densities in the relevant phase-space regions. These PDF uncertainties are dominant so that the reduction at the level of the scale uncertainties has a small impact. However, PDF uncertainties are expected to improve in the coming years thanks to new LHC data, so that it is mandatory to have NLO calculations available to get a better control on the predictions. In contrast, parton-shower effects are directly visible when distributions related to the two B-quarks considered as a pair are considered (Fig. 3). Focusing on the related transverse-momentum distributions, the fixed-order predictions (for which only the p T (B 1 + B 2 ) = 0 bin is populated at LO) diverge at small p T due to soft and collinear radiation giving rise to large logarithms that must be resummed to all orders to obtain reliable predictions. This resummation is effectively achieved by matching the fixed-order calculations to parton showers, and the resulting distribution exhibits a reliable behaviour with a peak for p T (B 1 + B 2 ) of about 10-20 GeV. Uncertainties originating from the choice of the shower algorithm and its inherent free parameters are, however, not estimated. The magnitude of the electroweak diagrams is also studied (dashed lines). In the case of lighter B-quarks, the theoretical predictions are essentially driven by the QCD contributions so that no differences can be noticed. This contrasts with the heavy B-quark case where electroweak diagrams enhance the total rate by about 30% (see Appendix A) and also impact the differential distributions both in terms of normalisation and shape. This originates from the different topologies of the electroweak diagrams that feature a t-channel colourless boson exchange. We now include the vector-like quark decays into a Higgs boson and a strange quark and reconstruct the final-state jets as detailed above. Considering only hard central jets for which |η| < 2.4 and p T > 30 GeV, we present the transverse-momentum distributions of the three leading jets as well as the spectrum of the H T variable defined as the scalar sum of the transverse momenta of all final-state jets and leptons in Fig. 4, the generated events being inclusive in the Higgs-boson decays. Leptons are included only if their transverse momentum is larger than 30 GeV, their pseudorapidity smaller than 2.4 in absolute value and if they are isolated from any jet by an angular distance in the transverse plane ΔR of at least 0.5. Focusing on fixed-order predictions matched to parton showers, we observe that the first two leading jets are in general hard, as they result from the decay of two heavy coloured particles. In contrast, the structure of the p T dependence of the third jet is more representative of the one expected from a radiation jet, this jet being most of the time arising from initial-state or final-state radiation. Turning to the H T distribution (lower right panel of the figure), one notices a peak for very small H T values when the B-quark mass is fixed at 500 GeV. This feature arises from events where the two jets originating from the heavy quark decays are mis-reconstructed, the leading jet being thus the radiation jet so that the associated jet activity in the events is not significantly large.

Single vector-like quark production in association with jets
Single vector-like quark production mechanisms are of a pure electroweak nature. The associated predictions are thus model-dependent as the sizes of the electroweak vectorlike quark couplings are free parameters of the model. Comparing vector-like quark single and pair production, the latter gets an enhancement originating from the presence of strong diagram contributions (first line of Fig. 1) together with a phase-space suppression for large vector-like quark mass values. In contrast, electroweak single vectorlike quark production is less suppressed for large vectorlike quark masses, which could compensate the weakness of the involved interaction vertices and make this channel the main LHC discovery mode for a heavy vectorlike quark. As a results, several ATLAS and CMS vector-  like quark searches also target their single-production mode [10,[45][46][47][48].
In this section, we focus on vector-like quark single production in association with jets, Other single-production mechanisms exist, with, for instance, a final-state gauge or a Higgs boson, but we refer to Sect. 3.3 for the latter. A representative set of Feynman diagrams related to single vector-like quark production with jets is shown on Fig. 5, and NLO cross-section calculation and event generation can be achieved with MadGraph5_aMC@NLO by typing in the program interpreter the command for single T production. Other processes with a different final-state vector-like quark can be accounted for with a similar syntax, and LO event generation only necessitates to remove the [QCD] tag. As for electroweak contributions to vector-like quark pair production, mixed QCD and electroweak loop diagrams appear at the NLO level and must be treated consistently for getting ultraviolet-finite results. This step being not automated, we provide scripts that steer the event generation process on the UFO model webpage. In Fig. 6 and Table 2, we present LO and NLO total cross section for single B quark production for the BZ1, BZ2, BW1 and BW2 scenarios introduced in Sect. 2.2 for which single vector-like B quark production occurs via Z -boson or Wboson exchanges. This consists of the first NLO predictions for a single vector-like quark production process. Although all depicted total cross sections exhibit a similar order of magnitude, a small hierarchy is observed. It is driven by an interplay of the parton densities that enhance mechanisms involving first generation quarks and of the new physics coupling strengths that are much larger for vector-like quark mixings with second generation (κ Q = 0.2) than with first generation quarks (κ Q = 0.07). For vector-like B-quark masses smaller than 500 GeV, single-production cross sections are slightly suppressed by a factor of 2 or 3 with respect to the strong production of a pair of B quarks (see Appendix A), but feature a similar order of magnitude for m B ∈ [500, 800] GeV. In contrast, single-B production always dominates for heavier vector-like quarks by up to two orders of magnitude for larger values of m B . For scenarios with smaller κ Q values, the changes in the relative importance of the two production channels can, however, be shifted towards higher masses.
In Fig. 7, we turn to the study of differential distributions related to inclusive single T -quark production at the LHC, focusing on the TZ1 scenario as an illustrative benchmark point and for a vector-like quark mass m T fixed to 500 GeV or 1500 GeV. Event generation and reconstruction is performed following the guidelines mentioned in Sect  Fig. 7 Differential distributions depicting the properties of a singly produced vector-like T (orT ) quark. We present its transverse momentum (upper), rapidity (lower left) and pseudorapidity (lower right) spectrum and compare fixed-order predictions at the LO accuracy (purple curve) and NLO accuracy (blue curve), as well as prediction including the matching of these two calculations to parton showers (green and red curves for the LO and NLO cases, respectively). We have fixed the heavy quark mass either to 500 GeV (upper series of curves) or to 1500 GeV (lower series of curves) we present the transverse momentum p T (T ), rapidity y(T ) and pseudorapidity η(T ) spectrum of the vector-like quark. We show predictions both at the fixed-order (purple and blue curves for the LO and NLO accuracy, respectively) and after matching the results to parton showers (green and red curves for the LO and NLO accuracy, respectively). Theoretical uncertainties originating from scale variations and parton densities are included for the matched predictions. In general, NLO effects only moderately affect the shape of the different spectra but in contrast drastically reduce the theoretical scale uncertainties and allow for a better prediction of the spectrum normalisation. Similarly to the pair-production case, matching to parton showers only mildly impacts fixed-order predictions for the properties of the produced heavy quarks.
Regardless of m T , the transverse-momentum distribution of singly produced T quarks exhibits a typical steeply falling behaviour with increasing values of p T (T ) (with a peak at low p T , which is invisible due to the binning choice), and the vector-like quark is quite forward by virtue of the process topology, with |η| > 2 in average. In Fig. 8, we study the properties of the decay products of the heavy quarks that consist of a Z -boson and an up quark in the TZ1 scenario. We first present the transverse-momentum spectrum of the three leading jets with a pseudorapidity satisfying |η| < 2.5 and a transverse-momentum larger than 30 GeV. Being inclusive in the Z -boson decay, the jets can originate either from the heavy T -quark decay, from the Zboson decay or from initial-or final-state radiation. Focusing isolated lepton transverse momenta. We compare LO (green) and NLO (red) predictions after matching the fixed-order calculations to parton showers. We have fixed the heavy quark mass either to 500 GeV or to 1500 GeV on the leading jet p T distribution, we observe that it peaks at half the T -quark mass, which shows that the leading jet is often issued from the heavy quark decay. In contrast, the second jet transverse-momentum distribution exhibit a plateau extending up to half the T -mass, which we conclude that it could alternatively originate either directly from the T -quark decay, or from the Z -boson decay. The third jet finally shows a different behaviour, the distribution peaking this time at a lower p T value, so that it is likely to be connected to the hard process.
In addition to the leading jet transverse-momentum distributions, we also present the distribution in the H T variable defined as the scalar sum of the transverse momentum of the final-state jets and isolated leptons, the latter being only considered if their pseudorapidity fulfils η| < 2.4 and their transverse momentum p T > 30 GeV. Moreover, we require the leptons to be well separated from any jet, imposing the angular distance ΔR to be larger than 0.5. The H T variable exhibits a non-trivial structure which peaks both at the Tquark mass m T and at m T /2, the second peak arising from cases where the Z -boson decays invisibly.
In all cases, the shapes turn out to be only slightly affected by the NLO effects and the uncertainties are drastically reduced.  for di-Higgs and diboson production, respectively, using the $$ symbol to remove any possible intermediate resonance.

Other processes impacted by
In the last case, the pair of symbols v v stands either for w+ w-or z z according to the weak boson under consideration. LO event generation can finally be achieved by removing the [QCD] tag.
Additional single vector-like quark production processes where the heavy quark is produced in association with a Standard Model boson can be considered as extra mechanisms useful for seeking for vector-like quarks (diagrams shown in the rightmost part of Fig. 9). Such processes can be simulated with MadGraph5_aMC@NLO, by typing in the commands generate p p > tp h $$ tp˜[QCD] add process p p > tp˜h $$ tp [QCD] for T /T production, as an example. Once again, the $$ symbol is used in order to avoid intermediate resonances. V Q production can be undertaken similarly.
As an illustrative example, we show in Table 3 and Fig. 10 total rates at LO and NLO for Higgs-boson (upper) and Zboson (lower) pair production, respectively, for scenarios featuring a T quark interacting with the first generation of Standard Model quarks and either the Higgs boson or the Zboson. We observe gigantic K -factors in the case of di-Higgs production. This enhancement is connected to an interplay of two effects. Turning back to the observations made in Sect. 3.1, the vector-like quark coupling to light quarks and a Higgs boson has strength that is proportional to the heavy quark mass, and thus becomes stronger and stronger with increasing vector-like quark masses. In addition, a new channel where the initial state is comprised of a gluon and a quark opens up at NLO. This component of the NLO cross section turns out to dominate due to the large gluon density in the proton. As a result, the total cross section for vector-likequark-mediated di-Higgs production is more or less constant with the heavy quark mass at the NLO QCD accuracy, which contrasts with the LO case.

Conclusions
We have modified a previously introduced model-independent parameterisation suitable for the study of the vector-like quark phenomenology at the LHC so that it is now suitable for NLO calculations in QCD matched to parton showers within the MadGraph5_aMC@NLO framework. We have illustrated its usage in the context of vector-like quark pair production and vector-like quark single production in association either with a jet or with a weak or Higgs boson. For all showcased processes, we have considered QCD and electroweak diagram contributions and investigated NLO and parton-shower effects on the normalisation and shapes of the associated kinematical distributions.
We have found that NLO K -factors are important, both globally (at the total-rate level) and at the differential distribution level and could hence potentially impact limits currently extracted from vector-like quark search results of the ATLAS and CMS collaborations. We have in particular noticed the existence of potentially huge K -factors for new physics scenarios involving the coupling of a heavy vector-like quark to first generation quarks and a Higgs boson due to new production channels that open at the NLO accuracy. This motivates further investigations, in particular to assess how an experimental analysis including detector effects could benefit from the gain in cross section stemming from the new topologies that dominate at NLO and to determine the impact on the current vector-like quark limits and LHC discovery potential.  quarks for the different scenarios introduced in Sect. 2.2. We depict the dependence of the cross sections on the vector-like quark mass and recall that the pure QCD contributions can be found in Table 2 as they are independent of the vectorlike quark nature. Results for X and Y pair production are not shown as non-QCD contributions are negligible for all considered scenarios.
Appendix B: Total cross sections for the single production of a T , X or Y quark In Fig. 12, Table 5 and Table 6, we present LO and NLO total cross sections related to the production of a single vectorlike T , X and Y quark for different scenarios introduced in Sect. 2.2. In each case, we depict the dependence of the cross sections on the vector-like quark mass and study the uncertainties stemming from scale variation.