Simulation of Sextet Diquark Production

We present a method for simulating the production and decay of particles in the sextet representation of $SU(3)_C$ including the simulation of QCD radiation. First results from the Monte Carlo simulation of sextet diquark production at the LHC including both resonant and pair production are presented. We include limits on resonant diquark production from recent ATLAS results and perform the first simulation studies of the less model dependent pair production mechanism.


Introduction
Many models of Beyond the Standard Model (BSM) physics require the inclusion of diquarks. For example, diquarks appear in a number of Grand Unified Theories (GUTs) and have even been postulated as a form of dynamical symmetry breaking, giving rise to the masses of particles, [1,2]. The colour sextet diquark is, in group theory language, a rank 2 symmetric tensor formed from the direct product of two fundamental representations 3 ⊗ 3 = 6 ⊕3. As such it is the lowest colour representation which has not been observed and therefore investigation of sextet diquark production at the CERN Large Hadron Collider (LHC) is interesting in its own right.
The LHC experiments have started data taking at the high energy frontier ( √ s = 7 TeV), allowing probes of energy scales not previously seen. At the LHC because the fundamental collisions are between the quarks and gluons inside the colliding protons the strong force is the dominant interaction allowing Quantum Chromodynamics (QCD) to be studied at these new high energies. As a diquark is produced via strong interactions, and with the potential of a relatively low mass, diquarks may be seen in the early stages of LHC data taking. The LHC also favours the formation of diquarks from the valence quarks as it is a protonproton collider as opposed to a proton-antiproton collider, such as the Tevatron.
Due to their exotic colour structure and SU (3) C quantum numbers, diquarks will give rise to jets in the detector. The expected signals will be either a resonance in the invariant dijet mass distribution or the production of two equal mass dijet systems in four jet events in the case of pair production. In order to study the experimental signatures of diquark production, a Monte Carlo simulation is required that includes the production of sextet particles, their perturbative decays and the full Monte Carlo machinery of showering (including the exotic colour structure) and hadronization.
Although significant efforts have been made to study the resonant production of diquarks [3][4][5][6][7][8][9] and also pair production [5,10,11], a full study of experimental signatures including Monte Carlo simulations has not been performed. In this paper, we discuss: the implementation of the diquark model in general purpose Monte Carlo event generators; place constraints on the coupling as a function of mass based on the latest ATLAS results [12,13]; present some results of invariant mass distributions both for resonant and pair production.

Simulation
Monte Carlo simulations describe high energy collisions using: 1. a hard perturbative, either leading-or next-to-leading-order, matrix element to simulate the fundamental hard collision process; 2. the parton shower algorithm which evolves from the scale of the hard process to a cut-off scale, O(1 GeV), via the successive radiation of soft and collinear quarks and gluons; 3. the generation of multiple perturbative scattering processes to simulate the underlying event; 4. the perturbative decay of any fundamental particles, with lifetimes shorter than the timescale for hadron formation, followed again by the simulation of QCD radiation from the coloured decay products using the parton shower formalism; 5. a multiple partonic collision model is used to simulate the underlying event; 6. a hadronization model which describes the formation of hadrons at the cutoff scale from the quarks and gluons produced during the parton shower; 7. the decays of the unstable hadrons produced by hadronization.
The various calculations, approximations and models used in these simulations are reviewed in Ref. [14]. Simulating most models of BSM physics only requires the implementation of the various hard production and decay processes with the simulation of perturbative QCD radiation and hadronization proceeding in the same way as in the Standard Model (SM), provided the new particles decay before forming hadrons. However, in models including the production and decay of sextet particles or R-parity violating SUSY 1 the new colour structures require changes to both the simulation of QCD radiation and the subsequent hadronization.
As with all models the simulation starts with the calculation of the hard production and decay processes using the most general Lagrangian for the coupling of the sextet particles to the quarks [4,6,9,16,17] where q L is the left-handed quark doublet, u R and d R are the right-handed quark singlet fields, and q c ≡ Cq T is the charge conjugate quark field. The colour and generation indices are omitted to give a more compact notation and the subscripts on the scalar, Φ, and vector, V µ , fields denote the SM electroweak gauge quantum numbers: (SU (2) L , U (1) Y ). The Lagrangian is assumed to be flavour diagonal to avoid any flavour changing currents arising from the new interactions. The kinetic and QCD terms in the Lagrangian are for scalar diquarks, where Φ is the scalar diquark field and for vector diquarks, where V µ is the vector diquark field. The covariant derivative D µ has the standard form for Quantum Chromodynamics. The simulation of perturbative QCD radiation, relies on the large number of colours, N C , limit for both the treatment of perturbative QCD radiation and the subsequent hadronization. In this approach particles in the fundamental representation of SU (N C ) carry a colour, those in the antifundamental representation an anticolour and those in the adjoint representation both a colour and an anticolour. This allows us to consider the flow of colour, via colour lines, in hard interactions which is determined by the colour structure of the hard perturbative matrix elements.
Using this colour flow both: • soft gluon radiation, where in a given hard process there is a maximum angle for the radiation of gluon related to the colour flow of the process; • and hadronization, where the first step of both the string [18,19] and cluster [20] models is the formation of colour singlet systems; can be simulated. This is complicated in models involving sextet particles where in the large-N C limit the sextet particles possess two fundamental colours, appropriately symmetrized. This cannot be handled by conventional Monte Carlo simulations which require all the colours of the particles to have fundamental colours and/or anticolours. In order to simulate these particles we choose to represent (anti)sextet particles as having two (anti)colours.
Consider the production and subsequent decay of a scalar sextet particle. In order to simulate QCD radiation from the intermediate sextet resonance we have to simulate the production and decay separately. The matrix element for the process is where i, j are colour indices of the sextet particle, M iprod is the matrix element for the production of a scalar sextet particle with colour i, four-momentum p and mass m, and M j decay is the matrix element for the decay of a scalar sextet particle with colour j.
This can be rewritten using δ i j = K i abK ba j where K andK are the Clebsch-Gordan coefficients in the sextet and antisextet representations, respectively. Hence In order to consider the intermediate sextet particle as having two fundamental colours we have absorbed the Clebsch-Gordan into the redefined production, M abprod , and decay matrix elements, M ba decay . From this, the colour partners of the decay products can be determined and the usual angular ordering procedure applied [14,21]. The radiation pattern of gluons from the quarks for the resonant production of diquarks, The dipole radiation function for two massless particles i and j radiating a gluon is where p i,j are the 4-momenta of the radiating particles and, k and ω are the 4-momentum of energy of the gluon, respectively. In terms of the angle between the gluon and particle i, θ i , and the angle between the particles i and j, θ ij ,   The last term in Eq. 5 can be neglected, as usual, due to both the 1 N C suppression, compared to the leading term, and the dynamical suppression in the massless limit because there is no collinear singularity in this term. The radiation pattern is shown in Figure 1, where the outgoing quarks were held at 45 • and 225 • with respect to the incoming beam direction. The full radiation pattern, the result after neglecting the subleading terms and azimuthally averaging, and the improved angular ordered result, where the full result is used instead of the azimuthal average inside the angular-ordered region, are shown. Improved angular ordering, as implemented in Herwig++ performs well in the collinear limit.
The details of the colour decomposition for both the resonant and pair production of sextet particles and the radiation of gluons from sextet particles are presented in Appendix A.3. The same approach was recently used in Mad-Graph 5 [22] including the ability to automatically decompose the colours for higher representations of the gauge group. They also study the impact of matching the hardest perturbative emission in resonance production which we have not considered but will not effect on our results.

Phenomenology
In order to study the phenomenology simulations were performed for the scalar Φ 1,4/3 and the vector V µ+ 2,5/6 diquarks which were chosen as they can be produced as s-channel resonances from the partonic collision of the valence up quarks.
In all our analyses, jets were clustered using the anti-k T algorithm [23], as  implemented in the FastJet package [24], using a radius parameter of R = 0.6. This choice is typical for the ATLAS experiment at the LHC [12]. The LO * * PDFs of Ref. [25], which are the default choice in Herwig++, were used. There are phenomenological constraints on the diquark couplings from D 0 −D 0 mixing and non-strange pion decays [7]. For the up-type quarks there are constraints require that g uu R 0.1 and g cc R ∼ 0.
The g L couplings have to be constrained due to minimal flavour violation as the left-handed CKM matrix is well known [3]. It was therefore decided to take the couplings where the numbered indices refer to the generation. The value of the coupling will affect any studies involving jets as the width of the particle varies as a function of the couplings. As seen in Figure 2, a larger coupling produces a larger width, contributing to the smearing of a peak in the invariant mass spectrum of the jets. If the coupling is less than the value chosen above, any peak maybe enhanced compared to what is presented in the following sections. The width as a function of the diquark mass is shown in Figure 3.

Resonance Production
If the diquark has an appropriate mass and coupling it may be resonantly produced at the LHC. The resonance production and subsequent decay of a general diquark (scalar or vector) is shown in Figure 4, where the decay of the diquark will depend on its mass and unknown couplings to the quarks, g. Figure 5 shows the cross section for the production of scalar and vector diquarks from incoming uu quarks (resonant production) and for incoming gluons (pair production). The resonant production cross section depends quadratically on the unknown diquark coupling to quarks, which as been assumed to be 0.1 in this plot, whereas the pair production cross section is independent of this coupling.
The diquark will decay into two quarks giving rise to at least two jets. The search for the production of a diquark via resonant production should therefore be in the dijet invariant mass spectrum, where a smeared peak is expected around the diquark mass. The primary background to this search channel is QCD 2 → 2 scattering processes.
The signal and background were simulated using Herwig++. The analysis and modelling of the backgrounds followed that suggested in Ref. [12]. The transverse momenta and pseudorapidities of the jets were required to be p 150 GeV, p 2 T > 60 GeV and |η 1,2 | < 2.5 where 1 is the hardest jet and 2 is the subleading jet. In addition we required that the dijet invariant mass, m jj satisfied m jj > 300 GeV and the rapidity difference between the leading and subleading jet was |∆η 12 | < 1.3. The dijet invariant mass spectrum after these cuts is shown in Figures 6 and 7 for √ s = 7 and 14 TeV, respectively. The diquarks were simulated at masses of 500 GeV, 800 GeV, 1200 GeV, 1600 GeV and 2000 GeV.
As simulating the QCD m jj spectrum at high masses is difficult, a functional form was fitted to the low masses, where the a i are fitted parameters and x = m jj / √ s, and extrapolated out into the high mass region.
The results of Ref. [12] can be used to impose constraints on the diquark coupling as a function of the diquark mass. The event selection from Ref. [12] was used to reproduce the correct acceptance. This requires that the event contains at least two jets with p T > 150 GeV and a subleading jet with p T > 30 GeV. Both the leading and subleading p T jets must satisfy |η j | < 2.5 with ∆η 12 < 1.3 and m jj > 150 GeV.
The signal, after the cuts, was fitted to a Gaussian distribution, with the mean fixed, m, at the simulated diquark mass to obtain the standard deviation of the distribution, σ, so that the results presented in Ref. [12] could be used to obtain the limits on the diquark coupling. As suggested in [12], long tails were removed by taking a window around the diquark mass of ±20% for the fit. If the σ/m value obtained was below the range of that given in the paper, then the number of events associated with the lowest σ/m for that mass was used. This allows a conservative estimate for the excluded coupling, as opposed to one which may be obtained by extrapolation into the unknown region. The limit on the diquark coupling is shown in Figure 8 where because the statistical errors were negligible, the bands shown come from varying the scale from 50% to 200% of the default scale choice, i.e. the diquark mass. Following the work of Ref. [12] the ATLAS collaboration has released an  [12,13]. The band reflects the uncertainty from varying the scale between 50% and 200% of the diquark mass.
updated analysis [13], including additional data corresponding to an integrated luminosity of 163 pb −1 . This analysis included slightly harder cuts requiring p T > 180 GeV and m jj > 170 GeV in addition to the cuts used in Ref. [12]. The limit obtained from this higher integrated luminosity analysis is also shown in Figure 8. We note that ATLAS performed better than the expected median limit in the Ref. [12] and worse than the expected median limit in [13] in the 1400 − 1600 GeV mass range, giving rise to the overlap in Figure 8.

Pair Production
The pair production of diquarks (scalar and vector) occurs via the Feynman diagrams shown in Figure 9. The pair production process has one main advantage over the resonant production, it does not depend on the unknown diquark coupling. Instead, the pair production cross section depends only on the SU (3) C representation, mass and spin of the particle. The pair production processes therefore has the potential to distinguish between whether a particle in the antitriplet or sextet representation was produced due to the dependency of the cross section on the colour representation.
To date there have been no studies of the experimental signals of diquark pair production. The cross section has been calculated [5,10,11] and some work towards a jet study, no Monte Carlo study has been performed.
The pair production and subsequent decay of diquarks is expected to give four jets, with two pairs of jets forming systems with the mass of the diquark. The Fig. 9: The four diagrams contributing to the pair production of a diquark.
backgrounds to the pair production of diquarks are: • vector boson W W , ZZ and ZW pair production; • vector boson, W and Z, production in association with additional jets; • top quark pair, tt, production; • QCD jet production.
The analysis proceeded by placing cuts on the four hardest jets: p 1 T > 150 GeV, p 2 T > 100 GeV, p 3 T > 60 GeV and p 4 T > 30 GeV, where the four jets i = 1, 4 are ordered in p T such that the first jet is the hardest. All four jets were required to have pseudorapidity |η i | < 3. Two pairs of jets were then formed with the pairing selected that minimized the mass difference between the two pairs of jets. If, after pairing, the two hardest jets are in the same pair of jets, the event was vetoed. The mass difference between the pairs was required to be less than 20 GeV.
The signal and backgrounds were simulated for the production of the Φ 1,1/3 and V µ+ 2,5/6 diquarks giving the results shown in Figures 10 and 11 for √ s = 7 and 14 TeV, respectively. As the backgrounds are dominated by QCD scattering, i.e. the contribution of the QCD scattering processes is approximately one hundred times that of all the other backgrounds combined, only the sum of the backgrounds is shown. As for the resonant production, the plots show the production of diquarks with masses of 500 GeV, 800 GeV and 1200 GeV, 1600 GeV and 2000 GeV. The low mass QCD background was fitted with Eq. 10 and extended out into the high mass region. A window of ±50 GeV was taken around the diquark mass and the S √ B was calculated for a number of luminosities the results of which are shown in Figure 12. We see that that the vector diquark manifests itself more prominently than the scalar, which is consistent with the increased cross section of the vector in the pair production process as seen in Figure 5. It will be hard to observe a scalar diquark using the pair production process at √ s = 7 TeV while with L = 5 fb −1 it should be possible to observe a vector diquark with mass less than 700 GeV.  There is a marked increase in discovery potential at the increased energy and luminosities running at √ s = 14 TeV brings. A vector diquark in the mass range presented here (< 2000 GeV) should be seen with L = 10 fb −1 , whereas even with L = 100 fb −1 only a low mass ( 1050 GeV) scalar diquark has potential for discovery.

Conclusions
In this paper we have presented a method for simulating the production and decay of particles in the sextet colour representation. This approach has been implemented in Herwig++ and will be available in a forthcoming release. This new simulation was used to simulate the production and decay of vector and scalar sextet diquarks at energies relevant to the LHC. Based on the findings and the latest ATLAS search for new particles in two-jet final states, new constraints have been put on the couplings of the diquarks to SM particles.
We have presented the first studies of the pair production mechanism, which is independent of the unknown coupling of the sextet diquark to quarks. This process has a promising search reach with the possibility of observing vector diquarks with masses less than 710 GeV at √ s = 7 TeV and both vector, for masses less than 2 TeV, and scalar, for masses less than 1 TeV, with the LHC running at design energy. Hopefully the availability of a Monte Carlo simulation of these processes will allow a more detailed experimental study.

Acknowledgements
We are grateful to all the other members of the Herwig++ collaboration for valuable discussions. We acknowledge the use of the UK Grid for Particle Physics in producing the results. This work was supported by the Science and Technology Facilities Council. DW acknowledges support by the STFC studentship ST/F007299/1.

A Colour Decomposition
Herwig++, as with all general purpose Monte Carlo generators, has all the machinery set up to work with the (anti)fundamental representation of SU (3) C . Therefore exotic color representations must be decomposed into a fundamental representation basis.
The results of the this decomposition for the resonant and pair production of sextet diquarks and the radiation of a gluon from a sextet diquark are explicitly outlined below.

A.1 Resonant Production
The diagram for resonant production and decay is shown in Figure 4. Following decomposition into the fundamental representation, there are two unique colour flows associated with this process as shown in Figure 13. The colour factor associated with these colour flows is N C (N C +1) 2 .
(a) (b) Fig. 13: Unique colour flows associated with the resonant production of a diquark.

A.2 Pair Production
The diagrams contributing to the pair production of diquarks are shown in Figure 9. Each of these diagrams is decomposed into a colour factor and colourless component. By taking the colour factors for the diagrams, adding and squaring, a table of colour factors can be produced for each term in the total matrix element squared.
The pair production process has twelve unique colour flows, as shown in Figure 14. The colours in Figure 14 have no physical meaning and are included as a visual aid.
From these colour flows, Table 1 is produced, where c 1 = Tab. 1: Associated colour factors for the diagrams shown in Figure 14. The values of c i are given in the text.

A.3 Shower
The splitting functions were decomposed in the same way as for the pair production. There are four unique colour flows associated with a diquark emitting a gluon, as shown in Figure 15, again where the colours are included as a visual aid.
Only the colour prefactor of the existing splitting functions is changed. The colour prefactor is given by 10 3 , i.e. the diquarks radiate 2 1 2 times more than a particle in the octet representation.
During the shower, it is assumed that gluons did not branch to form diquarks owing to the large diquark mass.