Event shape variables measured using multijet final states in proton-proton collisions at $\sqrt{s} =$ 13 TeV

The study of global event shape variables can provide sensitive tests of predictions for multijet production in proton-proton collisions. This paper presents a study of several event shape variables calculated using jet four momenta in proton-proton collisions at a centre-of-mass energy of 13 TeV and uses data recorded with the CMS detector at the LHC corresponding to an integrated luminosity of 2.2 fb$^{-1}$. After correcting for detector effects, the resulting distributions are compared with several theoretical predictions. The agreement generally improves as the energy, represented by the average transverse momentum of the two leading jets, increases.


Introduction
The production of quarks and gluons in hadron collisions and the process of hadron formation are subject to in-depth theoretical and experimental studies. The experiments at the CERN LHC have studied production of hadronic jets by measuring differential cross-sections, ratios of numbers of jets, angular distributions, etc., to deepen the understanding of quantum chromodynamics (QCD). While the production of quarks and gluons with large transverse momentum (p T ) is well described by calculations based on perturbative QCD, the hadronization process probes energy scales where perturbative calculations are not applicable. Instead, phenomenological models inspired by QCD are used to predict the experimental results.
Event shape variables (ESVs) are sensitive to the flow of energy in hadronic final states. These variables are safe from collinear and infrared divergences and have reduced experimental uncertainties [1]. Some distributions of ESVs are sensitive to the details of the hadronization process [2][3][4], so they can be used to tune parameters of Monte Carlo (MC) event generators, determine the strong coupling α S [5][6][7], and to search for new physics phenomena [8][9][10].
This paper reports a measurement of ESVs by the CMS Collaboration using hadronic jets in pp collisions at √ s = 13 TeV corresponding to an integrated luminosity (L int ) of 2.2 fb −1 . The following variables are studied: the complement of transverse thrust, total jet broadening, total jet mass, and total transverse jet mass. The theoretical uncertainties in the predictions of these ESVs can be reduced by careful choice of the quantity used to classify the energy scale of the events. Following Ref. [4], we use H T,2 = (p T,jet1 + p T,jet2 )/2, where p T,jet1 and p T,jet2 refer to the transverse momenta of the highest and second highest p T jets.The measured distributions are corrected for detector effects and compared with the predictions of QCD models implemented in the PYTHIA8 [27], MADGRAPH5 aMC@NLO+PYTHIA8 [28], and HERWIG++ [29] event generators.
The paper is organized as follows. The ESVs are discussed in Section 2. After briefly describing the elements of the CMS detector in Section 3, the jet reconstruction relevant to this analysis is described in Section 4. The data sample and event selection criteria are described in Section 5. Sections 6 and 7 present the unfolding technique and the systematic uncertainties, respectively. Section 8 contains comparisons between CMS data and theoretical predictions, and the results are summarized in Section 9.

Event shape variables
The four ESVs studied in this analysis are defined using the four-momenta of hadronic jets.
The complement of transverse thrust: The complement of thrust is defined as: where the thrust in the transverse plane is: Here, p T,i is the component of momentum of the i th jet perpendicular to the beam direction and thrust directionn T is the unit vector that maximizes the projection and defines the transverse thrust axis. The τ ⊥ is zero for a perfectly balanced two-jet event and is 1 − 2/π for an isotropic multijet event.
Total jet broadening: For each event, the transverse thrust axis is used to divide the event into upper (U) and lower (L) regions. The jets in U satisfy p T,i .n T > 0 and those in L have p T,i .n T < 0. For these two regions, the p T -weighted pseudorapidities and azimuthal angles are where X refers to the U or L regions. The jet broadening variable in each region is defined as where P T is the scalar p T sum of all the jets in the event. The total jet broadening is then defined as Total jet mass: The normalized squared invariant mass of the jets in the U and L regions of the event is defined by where M X is the invariant mass of the jets in the region X, and P is the scalar sum of the momenta of all central jets. The total jet mass is defined as the sum of the masses in the U and L regions, ρ Tot ≡ ρ U + ρ L .
Total transverse jet mass: The quantity corresponding to ρ Tot in the transverse plane, the total transverse jet mass (ρ T Tot ), is similarly calculated using p T,i of jets. These four ESVs probe different aspects of QCD [2] and are designed to have higher values for multijet, spherical events and lower values for back-to-back dijet events. While τ ⊥ is sensitive to the hard-scattering process, the jet masses and jet broadening depend more on the nonperturbative aspects of QCD, responsible for hadronisation process.

The CMS detector
The central feature of the CMS apparatus is a superconducting solenoid of 6 m internal diameter, providing a magnetic field of 3.8 T. The solenoid volume holds a silicon pixel and strip tracker, a lead tungstate crystal electromagnetic calorimeter (ECAL), and a brass and scintillator hadron calorimeter (HCAL), each composed of a barrel and two endcap sections. Steel and quartz-fibre Cherenkov hadron forward calorimeters extend the pseudorapidity (η) coverage provided by the barrel and endcap detectors to the region 3.0 < |η| < 5.2. Muons are measured in gas-ionisation detectors embedded in the steel flux-return yoke outside the solenoid. In the region |η| < 1.74, the HCAL cells have widths of 0.087 in η and 0.087 radians in azimuthal angle (φ). For |η| < 1.48, the HCAL cells map onto 5×5 ECAL crystals arrays in the η-φ plane to form calorimeter towers projecting radially outwards from close to the nominal interaction point. At larger values of η, the size in η of the towers increases and the matching ECAL arrays contain fewer crystals. CMS uses a two stage online trigger to select events for offline analysis. In the first stage, a hardware-based level-1 (L1) trigger uses information from calorimeter and muon subsystems and selects event at a rate of about 100 kHz. In the second stage, a softwarebased high-level trigger (HLT), running on computer farms, uses full event information and reduces the event rate to about 1 KHz before data storage. A more detailed description of the CMS detector can be found in Ref. [30].

Jet reconstruction
The particle-flow (PF) event algorithm [31] reconstructs photons, electrons, charged and neutral hadrons, and muons with an optimised combination of information from the various elements of the CMS detector. The energy of a photon is directly obtained from the ECAL measurement. The energy of an electron is determined from a combination of the electron momentum at the primary interaction vertex as determined by the tracker, the energy of the corresponding ECAL cluster, and the energy sum of all bremsstrahlung photons spatially compatible with originating from the electron track. The momentum of a muon is obtained from the curvature of the corresponding track. The energy of a charged hadron is determined from a combination of its momentum measured in the tracker and the matching ECAL and HCAL energy deposits, corrected for zero-suppression effects and for the response function of the calorimeters to hadronic showers. Finally, the energy of a neutral hadron is obtained from the corresponding energy deposits in ECAL and HCAL.
Jets are reconstructed from photons, electrons, charged and neutral hadrons, and muons using the anti-k T clustering algorithm [32,33] with a distance parameter R = 0.4. Measurement of jet energy is affected by contamination from additional pp interactions in the same bunch crossing (pileup), as well as by the nonuniform and nonlinear response of the CMS calorimeters. The technique of charged-hadron subtraction [31] is used to reduce the contribution of particles that originate from pileup interactions to the jet energy measurement. The jet four-momentum is corrected for the difference observed in simulation between jets built from reconstructed particles and generator-level particles. The jet mass and direction are kept constant for the corrections, which are functions of the η and p T of the jet, as well as the energy density and jet area quantities defined in Ref. [34]. The latter are used to correct the energy offset introduced by the pileup interactions. The energy of the jets is further corrected using dijet, Z+jet, and γ+jet events, where the p T -balance of the event is exploited. The jet energy resolution typically amounts to 15% at 10 GeV, 8% at 100 GeV, and 4% at 1 TeV.

Collision data
This analysis uses pp collision data collected in 2015 at √ s = 13 TeV, corresponding to L int = 2.2 fb −1 . Events are selected at L1 and HLT that have jet p T or H T,2 thresholds, respectively, as shown in Table 1. The turn-on point for each trigger, offline H T,2 at which the trigger is 99% efficient, is used to define the H T,2 ranges for events.
Collision and simulated events are required to have at least three jets with p T > 30 GeV within the coverage of the tracker |η| < 2.4. For each event, three jets are used for the calculation of the ESVs. The jets with the highest and the second-highest p T are selected. From the remaining jets, the one with the highest recoil term is selected as the third jet. The recoil term for jet k is The data sample is divided into eight H T,2 ranges such that the uncertainty due to the trigger inefficiency is negligible. The ranges (in GeV) are: 73-93, 93-165, 165-225, 225-298, 298-365, 365-452, 452-557 and >557, as shown in Table 1, with the number of events in each range.

Simulated events
Events are simulated using PYTHIA v8.212, MADGRAPH5 aMC@NLO V5 2.2.2+PYTHIA8, and HERWIG++ v2.7.1. The NNPDF3.0 [35] parton distribution function (PDF) set is used. The PYTHIA8 and HERWIG++ event generators use leading order 2→2 matrix element (ME) calculations and parton shower (PS) for generation of multijet topologies. The PYTHIA8 event generator uses a p T -ordered PS, and the underlying event description is based on the multiple parton interaction (MPI) model. Events are generated with two PYTHIA8 tunes: CUETP8M1 [36] and Monash [37]. Minimum bias data collected by the CMS experiment were used to derive the PYTHIA8 CUETP8M1 tune, which is based on the Monash tune. The MADGRAPH5 aMC@NLO generator uses ME calculations to generate hard-scattering events with two to four partons and PYTHIA8 CUETP8M1 for subsequent fragmentation and hadronization. The MLM [38] matching procedure is used to avoid double counting of jets between the ME calculation and the PS description. The HERWIG++ generator uses an angular-ordered PS. For simulated events, particle-level jets are obtained by applying the anti-k T clustering algorithm to all generated stable particles, excluding neutrinos, with R = 0.4.
The simulation events are passed through a complete and detailed reconstruction in the CMS detector using the same reconstruction as the collision events.

Unfolding of distributions
A reconstructed collision event differs from the true event because of finite resolution of the detector, detector acceptances, and uncertainties and efficiencies of measurement. Hence, the detector-level distributions obtained from data are unfolded to estimate the underlying particlelevel distributions, which can be compared with predictions from theoretical models as well as with results obtained by other experiments.
Simulated events passing through the complete detector simulation, event reconstruction, and selection chain are used to construct the response matrix for an ESV, which relates its particlelevel distribution with that at detector level. The response matrix incorporates all the experimental effects and is subsequently used as input for the unfolding of the observed distribution in data. Some events that satisfy the selection criteria at the particle level might not at the detector level, leading to an inefficiency. The reverse may also happen, leading to misidentification. Further, an event may migrate from one H T,2 range to another. The corresponding efficiency and misidentification rates are also incorporated in the unfolding process, and they contribute to the related uncertainty of the unfolding process.
To investigate possible bias due to the choice of an MC generator to construct the response matrices, we generate event samples from three different generators: PYTHIA8 CUETP8M1, MADGRAPH5 aMC@NLO, and HERWIG++. Each detector level distribution is unfolded using these three response matrices and the corresponding particle-level distributions are compared.
No evidence for significant bias is observed.
Two different methods, which are implemented in RooUnfold [39], are used for unfolding the observed distributions: D'Agostini iteration with early stopping [40], and Singular Value Decomposition (SVD) [41]. The difference between the unfolded distributions produced with these two methods is much smaller than 1%. Our unfolding is done using the D'Agostini iteration and PYTHIA8 CUETP8M1 is used for constructing the response matrix. The SVD method is used as a cross-check.

Systematic uncertainties
There are multiple sources of uncertainties in the unfolding process, and the contributions from each individual source are added in quadrature to obtain the total uncertainty. Figure 1 shows the total uncertainty and the contributions from various sources as a function of each ESV for the specific range 225 < H T,2 < 298 GeV.
To estimate the effect of each source, the four-momentum of each jet is scaled up and down by the corresponding uncertainty, the ESV is calculated, and the response matrix obtained with the nominal JES is used to unfold the distributions obtained with the nominal, scaled up, and scaled down JES values. For each bin of the unfolded distribution, the larger of the differences between the nominal, and the varied ones is taken as the systematic uncertainty. The systematic uncertainties due to different sources are then added in quadrature. For most bins in the distribution of an ESV, the uncertainty is 4-6%. However, it reaches about 12% for the highest and lowest bins of ρ Tot , lowest bins of ρ T Tot , and about 8% for the highest bins of B Tot . Typically JES is the largest source of systematic uncertainty in the ESVs.
• Jet energy resolution (JER): The JER is obtained from the ratio of p T of the two jets in dijet events as a function of p T and η [42]. It has been observed that the JER is worse in data compared to simulation. Hence, extra smearing is applied to the simulated events, and different response matrices are constructed. The detector-level distribution of an ESV is unfolded with the different response matrices incorporating the uncertainty due to JER. The estimated uncertainties in the ESVs are of the order of 1%.
• Unfolding: The detector-level distribution of an ESV obtained from simulated events of PYTHIA8 CUETP8M1 is unfolded with two response matrices derived from MAD-GRAPH5 aMC@NLO and HERWIG++, and compared with the corresponding particlelevel distribution in the same sample. Similar exercises are carried out for the MAD-GRAPH5 aMC@NLO sample using PYTHIA8 CUETP8M1 and HERWIG++ response matrices, and for the HERWIG++ sample using PYTHIA8 CUETP8M1 and MAD-GRAPH5 aMC@NLO response matrices. Out of these six differences for each bin, the largest is taken as the systematic uncertainty. In the closure tests of the individual response matrices, if, for a particular bin, the difference in the unfolded and generated values is larger than the uncertainty already assigned, the larger one is taken as the uncertainty due to the unfolding for that bin. The bias inherent in the DAgostini method is estimated by using different generators. The difference in the unfolded results is included as an unfolding uncertainty. The uncertainty due to unfolding is of the order of 2%, except for a few lowest, and highest bins where it dominates the total uncertainty.
• Parton distribution function: The uncertainty due to the PDFs in the particle-level distribution of an ESV is estimated using the 100 sets of NNPDF3.0 replicas. The standard deviation of the 100 values thus obtained for a bin is taken as the uncertainty due to PDFs for that bin. For most bins, the uncertainty due to the PDFs is less than 1%, but increases for higher values of the variables. For B Tot the uncertainty due to the PDFs increases very rapidly (>20%) and dominates for the last few bins.
The contribution of other sources of systematic uncertainty, i.e., pileup, and trigger efficiency are negligible.

Results
The modelling of initial-state radiation (ISR), final-state radiation (FSR) of gluons, and MPI in PYTHIA8 CUETP8M1 is tested by studying each aspect individually, via the comparison of simulated ESV distributions with data, as shown in Figure 2. This study shows that the effect of disabling ISR results in a very large shift of the ESVs to lower values, i.e., reducing the spherical nature of the multijet events. The effect of disabling the FSR is small compared to the ISR, and the effect of MPI is even smaller.
The unfolded distributions for the ESVs obtained from data are compared with the particlelevel predictions of various MC generators, as shown in The MPI parameters in the PYTHIA8 Monash and CUETP8M1 tunes are very similar. The predictions of these two tunes agree well for the four ESVs studied. In general, the agreement between them improves with increasing H T,2 . Both tunes show good agreement with data for the τ ⊥ and ρ T Tot variables, except for the two lowest ranges of H T,2 , and both overestimate the multijet contribution to ρ Tot and B Tot . We note that τ ⊥ and ρ T Tot variables are evaluated in the transverse plane, whereas B Tot and ρ Tot are evaluated using both longitudinal and transverse components of the jets. This indicates that the treatment of the energy flow in the transverse plane is modelled well in the Monash and CUETP8M1 tunes of PYTHIA8, whereas the energy  The HERWIG++ generator shows good agreement with data for all four ESVs studied, and it is better than the CUETP8M1 and Monash tunes of PYTHIA8 in predicting ρ Tot and B Tot . This implies its better treatment of energy flow out of the transverse plane. Although both PYTHIA8 and HERWIG++ use a PS approach to generate multijet events and hadronization, the former uses string fragmentation and a p T -ordered shower, whereas the latter uses cluster fragmentation and angular-ordered shower.
The MADGRAPH5 aMC@NLO generator shows good agreement with data for τ ⊥ and ρ T Tot and its agreement with data for ρ Tot and B Tot is much better compared to the CUETP8M1 and Monash tunes of PYTHIA8. The ME approach for generating multiparton hard scattering processes models the transverse as well as longitudinal flows of energy better than PYTHIA8.
The following features emerge from the comparison plots of the four ESVs.. Agreement between data and benchmark event generators improves with H T,2 . Figure 11 shows the evolution of the mean value of each ESV with H T,2 and confirms the above observations. With higher H T,2 , the initial partons are more boosted, and hence the event tends to be less spherical. Also, α S decreases with H T,2 , resulting in less emission of hard gluons, which further spoils the multijet, spherical nature of the event. Thus, the mean value of each ESV decreases with increasing  Figure 2: The effects of MPI, ISR, and FSR in PYTHIA8 CUETP8M1 on τ ⊥ (upper left), B Tot (upper right), ρ Tot (lower left) and ρ T Tot (lower right) for a typical range 225 < H T,2 < 298 GeV. The ratio plots for simulation (MC) with respect to data are shown in the lower panel of each plot. The inner gray band represents the statistical uncertainty and the yellow band represents the total uncertainty (systematic + statistical) in each plot.  , total jet mass (ρ Tot ) (lower left) and total transverse jet mass (ρ T Tot ) (lower right) for 73 < H T,2 < 93 GeV. In each ratio plot, the inner gray band represents statistical uncertainty and the yellow band represents the total uncertainty (systematic and statistical components added in quadrature) on data and the MC predictions include only statistical uncertainty.  Figure 4: Normalized differential distributions of unfolded data compared with theoretical (MC) predictions of PYTHIA8 CUETP8M1 (red line), PYTHIA8 Monash (blue dash-dotted line), MADGRAPH5 aMC@NLO (pink dash-dot-dotted line) and HERWIG++ (brown dash-dot-dotted line) as a function of ESV: complement of transverse thrust (τ ⊥ ) (upper left), total jet broadening (B Tot ) (upper right), total jet mass (ρ Tot ) (lower left) and total transverse jet mass (ρ T Tot ) (lower right) for 93 < H T,2 < 165 GeV. In each ratio plot, the inner gray band represents statistical uncertainty and the yellow band represents the total uncertainty (systematic and statistical components added in quadrature) on data and the MC predictions include only statistical uncertainty.  Figure 5: Normalized differential distributions of unfolded data compared with theoretical (MC) predictions of PYTHIA8 CUETP8M1 (red line), PYTHIA8 Monash (blue dash-dotted line), MADGRAPH5 aMC@NLO (pink dash-dot-dotted line) and HERWIG++ (brown dash-dot-dotted line) as a function of ESV: complement of transverse thrust (τ ⊥ ) (upper left), total jet broadening (B Tot ) (upper right), total jet mass (ρ Tot ) (lower left) and total transverse jet mass (ρ T Tot ) (lower right) for 165 < H T,2 < 225 GeV. In each ratio plot, the inner gray band represents statistical uncertainty and the yellow band represents the total uncertainty (systematic and statistical components added in quadrature) on data and the MC predictions include only statistical uncertainty.  Figure 6: Normalized differential distributions of unfolded data compared with theoretical (MC) predictions of PYTHIA8 CUETP8M1 (red line), PYTHIA8 Monash (blue dash-dotted line), MADGRAPH5 aMC@NLO (pink dash-dot-dotted line) and HERWIG++ (brown dash-dot-dotted line) as a function of ESV: complement of transverse thrust (τ ⊥ ) (upper left), total jet broadening (B Tot ) (upper right), total jet mass (ρ Tot ) (lower left) and total transverse jet mass (ρ T Tot ) (lower right) for 225 < H T,2 < 298 GeV. In each ratio plot, the inner gray band represents statistical uncertainty and the yellow band represents the total uncertainty (systematic and statistical components added in quadrature) on data and the MC predictions include only statistical uncertainty.  Figure 7: Normalized differential distributions of unfolded data compared with theoretical (MC) predictions of PYTHIA8 CUETP8M1 (red line), PYTHIA8 Monash (blue dash-dotted line), MADGRAPH5 aMC@NLO (pink dash-dot-dotted line) and HERWIG++ (brown dash-dot-dotted line) as a function of ESV: complement of transverse thrust (τ ⊥ ) (upper left), total jet broadening (B Tot ) (upper right), total jet mass (ρ Tot ) (lower left) and total transverse jet mass (ρ T Tot ) (lower right) for 298 < H T,2 < 365 GeV. In each ratio plot, the inner gray band represents statistical uncertainty and the yellow band represents the total uncertainty (systematic and statistical components added in quadrature) on data and the MC predictions include only statistical uncertainty.  Figure 8: Normalized differential distributions of unfolded data compared with theoretical (MC) predictions of PYTHIA8 CUETP8M1 (red line), PYTHIA8 Monash (blue dash-dotted line), MADGRAPH5 aMC@NLO (pink dash-dot-dotted line) and HERWIG++ (brown dash-dot-dotted line) as a function of ESV: complement of transverse thrust (τ ⊥ ) (upper left), total jet broadening (B Tot ) (upper right), total jet mass (ρ Tot ) (lower left) and total transverse jet mass (ρ T Tot ) (lower right) for 365 < H T,2 < 452 GeV. In each ratio plot, the inner gray band represents statistical uncertainty and the yellow band represents the total uncertainty (systematic and statistical components added in quadrature) on data and the MC predictions include only statistical uncertainty.  Figure 9: Normalized differential distributions of unfolded data compared with theoretical (MC) predictions of PYTHIA8 CUETP8M1 (red line), PYTHIA8 Monash (blue dash-dotted line), MADGRAPH5 aMC@NLO (pink dash-dot-dotted line) and HERWIG++ (brown dash-dot-dotted line) as a function of ESV: complement of transverse thrust (τ ⊥ ) (upper left), total jet broadening (B Tot ) (upper right), total jet mass (ρ Tot ) (lower left) and total transverse jet mass (ρ T Tot ) (lower right) for 452 < H T,2 < 557 GeV. In each ratio plot, the inner gray band represents statistical uncertainty and the yellow band represents the total uncertainty (systematic and statistical components added in quadrature) on data and the MC predictions include only statistical uncertainty.  Figure 10: Normalized differential distributions of unfolded data compared with theoretical (MC) predictions of PYTHIA8 CUETP8M1 (red line), PYTHIA8 Monash (blue dash-dotted line), MADGRAPH5 aMC@NLO (pink dash-dot-dotted line) and HERWIG++ (brown dash-dot-dotted line) as a function of ESV: complement of transverse thrust (τ ⊥ ) (upper left), total jet broadening (B Tot ) (upper right), total jet mass (ρ Tot ) transverse jet mass (ρ T Tot ) (lower left) and total transverse jet mass (ρ T Tot ) (lower right) for H T,2 > 557 GeV. In each ratio plot, the inner gray band represents statistical uncertainty and the yellow band represents the total uncertainty (systematic and statistical components added in quadrature) on data and the MC predictions include only statistical uncertainty.

Summary
This paper presents the first measurement at √ s = 13 TeV of four event shape variables: complement of transverse thrust (τ ⊥ ), total jet broadening (B Tot ), total jet mass (ρ Tot ), and total transverse jet mass (ρ T Tot ) using proton-proton collision data. It also covers a wider range of energy than the analysis at √ s = 7 TeV [19,22]. Data are compared with theoretical predictions from event generators PYTHIA8, HERWIG++, and MADGRAPH5 aMC@NLO+PYTHIA8. The PYTHIA8 generator describes the flow of energy in the transverse plane well as seen in the τ ⊥ and ρ T Tot distributions. HERWIG++ and MADGRAPH5 aMC@NLO show good agreement with the data for all the four event shape variables and are better than PYTHIA8 in predicting ρ Tot and B Tot . A study of the effects of initial state radiation, final state radiation, and multiple parton interactions in PYTHIA8 is also presented.

Acknowledgments
We congratulate our colleagues in the CERN accelerator departments for the excellent performance of the LHC and thank the technical and administrative staffs at CERN and at other CMS institutes for their contributions to the success of the CMS effort. In addition, we gratefully acknowledge the computing centres and personnel of the Worldwide LHC Computing Grid for delivering so effectively the computing infrastructure essential to our analyses. [12] DELPHI Collaboration, "Tuning and test of fragmentation models based on identified particles and precision event shape data", Z. Phys. C 73 (1996) 11, doi:10.1007/s002880050295.