Pursuit of paired dijet resonances in the Run 2 dataset with ATLAS

New particles with large masses that decay into hadronically interacting particles are predicted by many models of physics beyond the Standard Model. A search for a massive resonance that decays into pairs of dijet resonances is performed using 140 fb − 1 of proton-proton collisions at ﬃﬃﬃ s p ¼ 13 TeV recorded by the ATLAS detector during Run 2 of the Large Hadron Collider. Resonances are searched for in the invariant mass of the tetrajet system, and in the average invariant mass of the pair of dijet systems. A data-driven background estimate is obtained by fitting the tetrajet and dijet invariant mass distributions with a four-parameter dijet function and a search for local excesses from resonant production of dijet pairs is performed. No significant excess of events beyond the Standard Model expectation is observed, and upper limits are set on the production cross sections of new physics scenarios.


Introduction
New massive particles that decay into hadronically interacting quarks and gluons are predicted in many scenarios of physics beyond the Standard Model (BSM) accessible at the Large Hadron Collider (LHC), including well-motivated models of particle dark matter [1][2][3][4][5][6][7][8] and models with large extra spatial dimensions [9][10][11][12][13][14][15].Quarks and gluons produced at high energies fragment and hadronize into collimated jets of particles [16], observable in particle detectors like ATLAS [17].The majority of Standard Model (SM) multĳet event production occurs via nonresonant quantum chromodynamics (QCD) processes, resulting in multĳet systems with smoothly falling invariant mass distributions.The large production crosssection of multĳet processes can make searching for fully hadronic BSM signatures challenging, especially without the presence of other distinguishing feaures like leptons and/or missing transverse momentum [18].However, when massive particles decay into pairs of jets ('dĳets') via -channel interactions, the invariant mass spectrum of the dĳet system exhibits the signature of the massive particle as a resonance around its mass value.While the rate of new particle production may be too low that no resonance is obvious, such models may be detected using data-driven techniques analyzing the smoothly-falling invariant mass distribution of the SM background.Searches for dĳet resonances have been a cornerstone of collider physics at the LHC [19][20][21][22][23][24][25][26][27][28][29][30][31][32][33][34][35][36] and at earlier colliders [37][38][39][40].This paper presents a search for a generic massive resonance  that decays into two pairs of intermediate resonances  with the same mass, each decaying into two partons and so typically producing a pair of dĳet systems.This decay structure is represented schematically in Figure 1.Examples of exotic physics models that could produce such a final state topology include scalar diquark [41][42][43][44] and coloron states [45][46][47][48], and additional new particle content such as a vector-like quarks that interact in pairs with the massive diquark or coloron and decay hadronically [49][50][51].The analysis is performed on the Run 2 proton-proton ( ) collision data recorded by the ATLAS experiment.Previous searches for signals with this resonance structure have been performed by the ATLAS [52][53][54][55] and CMS [56][57][58] collaborations at the LHC.Most recently, in Ref. [58], the CMS Collaboration studied final states where pairs of dĳet resonances are collimated but insufficiently boosted to be reconstructed in a single large-radius jet (e.g. a jet reconstructed with radius parameter =0.8-1.0 [59, 60]).A small, locally significant excess of events (3.9 standard deviations from two events, corresponding to a 1.6 standard deviation global significance) was observed with tetrajet resonance masses around  4j ∼ 8 TeV producing dĳet resonances with average masses around  2j ∼ 2 TeV.This prompted an investigation of such final state configurations using ATLAS data.As there are two resonances with different masses involved in the final state topology ( and ), both the tetrajet system and the average dĳet system invariant masses are separately studied using the BumpHunter algorithm [61][62][63][64][65].A data-driven background estimate is obtained by fitting these invariant mass distributions with a functional form.
An outline of the remainder of this paper is as follows.Section 2 provides overviews of the ATLAS detector, the Run 2   data sample and the signal and background Monte Carlo (MC) simulations used in this search.This is followed in Section 3 by a description of the analysis methodology including jet reconstruction, event selection, the data-driven background estimation procedure and the systematic uncertainties considered.The main results are presented in Section 4, interpreting the observed data in terms of upper limits on the production cross-sections of new physics scenarios.Concluding remarks are made in Section 5.

ATLAS, the Run data, and simulation 2.1 The ATLAS detector
The ATLAS detector [17] at the LHC covers nearly the entire solid angle around the collision point. 1 It consists of an inner tracking detector surrounded by a thin superconducting solenoid, electromagnetic and 1 ATLAS uses a right-handed coordinate system with its origin at the nominal interaction point (IP) in the center of the detector and the -axis along the beam pipe.The -axis points from the IP to the center of the LHC ring, and the -axis points upwards.Cylindrical coordinates (, ) are used in the transverse plane,  being the azimuthal angle around the -axis.
hadron calorimeters, and a muon spectrometer incorporating three large superconducting air-core toroidal magnets.
The inner-detector system is immersed in a 2 T axial magnetic field and provides charged-particle tracking in the range of || < 2.5.The high-granularity silicon pixel detector covers the vertex region and typically provides four measurements per track, the first hit normally being in the insertable B-layer installed before Run 2 [66,67].It is followed by the silicon microstrip tracker , which usually provides eight measurements per track.These silicon detectors are complemented by the transition radiation tracker (TRT), which enables radially extended track reconstruction up to || = 2.0.The TRT also provides electron identification information based on the fraction of hits (typically 30 in total) above a higher energy-deposit threshold corresponding to transition radiation.
The calorimeter system covers the pseudorapidity range of || < 4.9.Within the region || < 3.2, electromagnetic calorimetry is provided by barrel and endcap high-granularity lead/liquid-argon (LAr) calorimeters, with an additional thin LAr presampler covering || < 1.8 to correct for energy loss in material upstream of the calorimeters.Hadron calorimetry is provided by the steel/scintillator-tile calorimeter, segmented into three barrel structures with || < 1.7, and two copper/LAr hadron endcap calorimeters.The solid angle coverage is completed with forward copper/LAr and tungsten/LAr calorimeter modules optimized for electromagnetic and hadronic energy measurements respectively.
The muon spectrometer comprises separate trigger and high-precision tracking chambers measuring the deflection of muons in a magnetic field generated by the superconducting air-core toroidal magnets.The field integral of the toroids ranges between 2.0 and 6.0 T m across most of the detector.Three layers of precision chambers, each consisting of layers of monitored drift tubes, cover the region || < 2.7, complemented by cathode-strip chambers in the forward region, where the background is highest.The muon trigger system covers the range of || < 2.4 with resistive-plate chambers in the barrel, and thin-gap chambers in the endcap regions.
Interesting events are selected by the first-level trigger system implemented in custom hardware, followed by selections made by algorithms implemented in software in the high-level trigger [68].The first-level trigger accepts events from the 40 MHz bunch crossings at a rate below 100 kHz, which the high-level trigger further reduces to record events to disk at about 1 kHz.
An extensive software suite [69] is used in data simulation, in the reconstruction and analysis of real and simulated data, in detector operations, and in the trigger and data acquisition systems of the experiment.

The Run 2 data sample
This analysis is performed using data from LHC   collisions with √  = 13 TeV, collected during 2015-2018 with the ATLAS detector.The total integrated luminosity of this data sample is 140 fb −1 .The uncertainty in the combined 2015-2018 integrated luminosity is 0.83% [70], obtained using the LUCID-2 detector [71] for the primary luminosity measurements.Due to the high instantaneous luminosity and the large total inelastic   cross-section, there are, on average, 33.7 simultaneous collisions ('pileup') in each bunch crossing.Data are required to satisfy certain quality requirements [72] to be included in the analysis.
During certain data-taking periods, modules of the tile calorimeter were disabled.A study of the impact of vetoing these disabled modules in MC and data was performed, and found to have a negligible impact on the background shape modeling and expected limits.

Simulated event samples
Samples of MC simulated signal and background (multĳet) events are used for optimisation, estimation of possible signal contributions, and validation of background estimation strategies.Pythia 8.230 [73,74] was used as the nominal MC generator for both the signal and the background events.Events were simulated using the A14 set of tuned parameters ('tune') [75], the Lund string hadronization model and the NNPDF2.3lo[76] leading-order (LO) parton distribution function (PDF) set.The Pythia parton shower algorithm uses a dipole-style transverse momentum ( T ) ordered evolution, and its renormalization and factorization scales were set to the geometric mean of the squared transverse masses of the outgoing particles.EvtGen [77] was used to model decays of heavy flavor hadrons.
Signal samples were generated with Pythia 8.230 using the process  ′ →  , where the mass of the  ′ corresponds to the  mass, and the / masses were both set to be equal to the  mass.The  and  full widths at half maximum of a relativistic Breit-Wigner were set to 0.1 GeV, so that the width of the resonance is determined by the detector resolution (typically ranging between 1-4% for both  4j and ⟨ 2j ⟩).These exotic  bosons were forced to decay into quark-antiquark pairs, and decays into top-antitop quark pairs were disabled.Representative  4j and ⟨ 2j ⟩ distributions are shown in Figure 2 for several different choices of  = ⟨ 2j ⟩/ 4j (see Section 3.2) with   = 6000 GeV.The signal distributions have clear peaks near the generated signal masses.Background samples of 'hard-QCD' multĳet processes were simulated using the same Pythia settings.These samples were used to optimize aspects of the analysis in early stages, although they are not used for the final background estimate.Additional background multĳet samples were simulated with Sherpa 2.2.5 [78] to test the robustness of the background estimation procedure, using the default AHADIC cluster hadronization model [79].This sample includes LO matrix element calculations for 2 → 2 processes, and used the Sherpa parton shower algorithm based on Catani-Seymour dipole subtraction [80].It used the CT14nnlo next-to-next-to-leading-order (NNLO) PDF [81] set for matrix element calculations and CT10 for multi-parton interactions [82].
Simulated background events were passed through a detailed detector simulation [83] based on Geant4 [84], while simulated signal events were reconstructed with a fast simulation that uses a parameterisation of the ATLAS calorimeter response [85].In both cases, the samples were overlayed with minimum-bias interactions simulated using Pythia 8 with the A3 tune [86] and the NNPDF2.3loPDF set to represent pileup interactions.The distribution of the average number of pileup interactions in simulation is reweighted during data analysis to match that observed in Run 2 data.
Additional details of the MC samples used in this measurement may be found in Ref. [87].

Particle flow jets
Jets are reconstructed from particle flow objects [88] using the anti-  algorithm [89] as implemented in FastJet [90], using a jet radius parameter  = 0.4.The ATLAS particle flow algorithm combines measurements from the ATLAS inner detector and calorimeter systems [91] to improve the jet energy resolution (JER), reduce sensitivity to pileup effects, and improve the jet reconstruction efficiency (especially at low jet  T ) relative to the jet reconstruction based on calorimeter signals alone.Jets are required to have a  T > 60 GeV and a rapidity  satisfying || < 2.4.The jet energy scale (JES) of particle flow jets is calibrated using a combination of simulation-based and in situ corrections [92].

Event selection
To be considered for analysis, all detector-level events are required to have at least one primary vertex reconstructed from two or more inner-detector tracks with  T > 500 MeV.Events are required to have at least four jets, from which two dĳet pairs are reconstructed.In the selected events, an event selection similar to that of Ref. [58] is applied to allow direct comparisons of the two searches.The two dĳet pairs are determined by minimizing the Δ between the jets, defined as where A, B, C, and D are the ordering of the four highest- T jets in the event that minimizes Δ.The value of 0.8 ensures that the reconstructed pair of dĳet systems are collimated, but not so boosted that they will be reconstructed as a large- jet.Once the two dĳet pairs AB and CD are selected, the dĳet systems are required to satisfy angular requirements In addition, the mass asymmetry between the two dĳet pairs is required to satisfy After the event selection procedure is complete, the observables of interest are the invariant mass of the tetrajet system,  4j (a proxy for   ), and the average invariant mass of the two dĳet systems, ⟨ 2j ⟩ (a proxy  For the background distribution,  is less correlated with  4j than ⟨ 2j ⟩.The analysis is performed in regions of  rather than ⟨ 2j ⟩, to reduce background sculpting in the tetrajet and average dĳet invariant mass spectra due to the selection criteria.Twelve different  regions are used, evenly spaced to cover 0.10 <  < 0.34.For each  region, separate fits are performed for  4j and ⟨ 2j ⟩.
The combined acceptance times efficiency of all analysis selections is between 12%-45% for signal events as a function of the signal particle masses   and   for 2000 GeV <   < 10000 GeV, and 500 GeV <   < 3300 GeV.
Events in data are required to have been selected by one of several single-jet triggers, whose thresholds varied depending on the data-taking period during Run 2. The particular triggers used to select events in a given run period were always the lowest unprescaled triggers recording data during that time.The single jet triggers become fully efficient for different values of  4j and ⟨ 2j ⟩, depending on the  bin that is selected.To use fully efficient triggers while also retaining sensitivity to the widest range of values, different minimum  4j and ⟨ 2j ⟩ thresholds are imposed on the data and simulated signal events used in the interpretation of the search (see Section 4).The trigger thresholds were optimised such that they were at least 99.5% efficient for each trigger, and are listed in Table 1 for the various  regions used in the search.Over 220 000 events in the Run 2 data sample satisfy the analysis selections.

Signal templates
The results are interpreted using model-independent and model-dependent strategies.For the modelindependent results, the signal templates are modeled as Gaussian distributions, with a mean equal to   and   for the  4j and ⟨ 2j ⟩ distributions respectively.Studying a range of template widths is important, as the theoretical width of a signal can vary across possible signal models.In this interpretation, template widths ranging from 5% to 15% for both  4j and ⟨ 2j ⟩were used.The upper end of the template width is determined from the results of the spurious signal test described in Section 3.4.
For the model-dependent limits, the shape of the  4j and ⟨ 2j ⟩ distributions are parameterized using a Crystal Ball function [93], which provides a good description of the shape of mass distribution.Signal samples are produced with a limited set of signal masses, and these templates are used as inputs to interpolation between mass points to provide a finer signal grid.The interpolation is done separately for each  region by morphing the parameterized Crystal Ball fits of the signal shape.

Background estimation
Non-resonant QCD processes, which constitute the SM background for this search, result in multĳet systems with smoothly falling invariant mass distributions.To estimate this background in the search regions, a parametric function is fit to the observed data distributions in 1 GeV bins: where  1 ,  .For the background estimation, a 4-parameter fit is used, where  5 is set to zero, while the 5-parameter fit is used to produce pseudodata to validate the fit strategy.Three-parameter fit functions were also studied but did not have sufficient flexibility to describe the background.
The background distribution is fit using a binned, maximum-likelihood fit.In background-only fits, the signal strength is set to zero, while in the signal-plus-background fits, the signal strength is left as a free parameter.For the model-dependent interpretation, the signal probability density function is defined as the Crystal Ball function fit to the simulated signal events.For the model-independent interpretations, the signal is parameterized as a Gaussian distribution, where the signal width is set to be a fixed fraction of the signal peak.
The data-driven background fitting procedure was validated with MC simulation using several cross-checks, including 'spurious signal tests' and 'signal injection tests'.Tests involving signal-plus-background fits are performed using both templates derived from the simulated signal samples described in Section 2.3, and for model-independent Gaussian signal shapes, using signal widths of 5%, 10%, and 15% of the signal peak.
The spurious signal test evaluates whether the fitting procedure is biased in a manner that will produce a non-zero extracted signal when fitting a data sample with no true signal.This test is performed for a 4-parameter fit function by performing a signal plus background fit with a Gaussian signal hypothesis of a specified width to 100 pseudodata distributions that are generated from background-only fits to the data distribution with an 5-parameter function, to provide more flexibility to the background distribution than to the fit function.Signal widths for both the  and  ranging from 5% to 15% of the signal mean are tested, as well as using the signal templates directly.For each pseudodata distribution, the number of extracted signal events,   , is determined, and the median value and standard deviation of   across all pseudodata distributions are taken to be  spur and  spur respectively.To satisfy the spurious signal requirements,  spur / spur is required to be less than 0.5.For  4j and ⟨ 2j ⟩, all  regions satisfy this criterion for the signal templates, and for Gaussian signals with widths of 5%, 10%, and 15%.
The signal injection test is performed to ensure that the background fit is able to extract a signal component with the expected signal strength.Simulated signal models with the Gaussian templates with signal widths of 5% to 15% and signal templates are included in the fitted background distribution with a given signal cross-section selected to be in the range of 0-5, where  =   / √   and the number of signal and background events are determined using a 2 window around the injected signal peak in each test.For these studies, the extracted signal strength is scaled to be the number of signal events within a 2 window around the injected signal peak, to match the definition used for the injection.The injected signals in this study were extracted for 100 pseudodata distributions, with the requirement that the median extracted significance is within 0.5 of the injected significance for fits to both the  4j and ⟨ 2j ⟩ distributions in all analysis  regions.All signal templates and Gaussian signals that passed the spurious signal tests also passed the signal injection tests.

Systematic and statistical uncertainties
When interpreting the analysis in terms of candidate signal models, the impact of various experimental and theoretical sources of uncertainty is considered.The uncertainties in the luminosity and parton distribution functions are included as Gaussian constraints on the yield, while the uncertainties in the jet energy scale, jet energy resolution, tune, and theoretical renormalization and factorization scale are included as Gaussian constraints on the shape of the distribution.The uncertainty in the background modeling is implemented as an additional 'spurious' signal-like contribution.
Luminosity.The uncertainty in the combined 2015-2018 integrated luminosity is 0.83% [70], obtained using the LUCID-2 detector [71] for the primary luminosity measurements.It is treated as a single normalization uncertainty applied as a scale factor to the signal models.

Parton distribution functions.
The theoretical uncertainty envelope associated with the NNPDF 2.3 LO set of PDFs is propagated through the analysis, where their impact is primarily on the normalization of the signal events.The change in analysis selection efficiency is recalculated for each provided PDF variation, and the standard deviation of all such variations is taken as a measure of the systematic uncertainty due to the PDFs.This uncertainty is a sub-1% effect for all signal models considered.
Jet energy scale and resolution.Systematic uncertainties in the  = 0.4 JES and JER are evaluated using a series of in situ measurements and simulation-based techniques, documented in Ref. [92].Improvements have been made to the component of the jet energy scale uncertainty related to the extrapolation of single-hadron response measurements [95,96] and combined test-beam results [97,98] into jets.These uncertainties are reduced by roughly a factor of two compared to those reported in Ref. [92].Uncertainties due to differences between the gluon-initiated jet energy response of different MC generator setups have also been reduced ('jet flavor response' in Ref. [92]), by performing more granular comparisons of the effect of different parton shower and hadronization models on the jet response using the samples documented in Ref. [87].Following the improved procedure compared with that documented in Ref. [92], the most significant source of uncertainty in the JES now originates from the absolute in situ JES calibration.
Variation of initial-state   value in the A14 set of tuned parameters.The A14 set of tuned parameters used in the Pythia 8 signal simulation includes a pair of 'eigentune variations' that can be used to assess the sensitivity of an analysis to the value of the QCD coupling,   , in initial-state radiation (ISR) [75].The value of   was varied between 0.115 − 0.140 from its initial value of 0.127.The impact of this variation is negligible compared with other systematic variations.
Theoretical renormalization and factorization scale variations.The QCD renormalization and factorization scales (  ,   ) used in the parton shower of the Pythia 8 signal samples are each varied up and down by a factor of two, via weights provided by the Pythia event generator [99].These variations assess the sensitivity of the analysis to parton shower configurations that contain branchings that may compromise the Pythia parton shower's underlying assumptions.Scale variations for such configurations will result in a large variation for that shower.The theoretical uncertainties resulting from these scale variations in the mean of the signal distribution are typically less than 0.5%, and are smaller than the JES uncertainties.
Background modeling A systematic uncertainty to cover potential modeling biases is accounted for using the spurious signal  modeling .The value of  modeling is determined as the envelope of | spur | over  4j and ⟨ 2j ⟩ respectively.This is implemented as an additional signal contribution, such that where  signal ( , ) is the number of extracted signal events at a given   or   , L, , and  are the integrated luminosity, acceptance, and efficiency factors respectively and  modeling is a nuisance parameter associated with the modeling uncertainty.The acceptance is defined as the fraction of simulated events at generator level passing the analysis selection cuts, while the efficiency is the fraction of reconstructed events passing the selection.

Search results
Example tetrajet and average dĳet invariant mass distributions in data, together with the corresponding fitted background estimates, are shown in Figure 4 for two representative  regions.The example  regions are selected to show the highest tetrajet invariant mass, and the most significant localized excess observed in data.The data are well-described by the 4-parameter fit function in all  regions, and the global  2 -value ranges from 0.74 to 1.00 for the  4j spectra, and from 0.08 to 1.00 for the ⟨ 2j ⟩ spectra.The BumpHunter [62, 63] algorithm, as implemented in pyBumpHunter [64,65], is used to quantify the statistical significance of possible resonant signals that may be present in the  4j and ⟨ 2j ⟩ distributions.This is performed using mass bins where the bin width is determined by the mass resolution of  4j or ⟨ 2j ⟩ as a function of the mass, where the mass resolution is determined using a Gaussian fit to the mass response distribution.The width of the invariant mass window scanned by BumpHunter is varied between two and six resolution bins, and all possible windows of the  4j and ⟨ 2j ⟩ distributions are scanned, in each  region.For each scanned window, BumpHunter evaluates the statistical significance of the observed difference between the data distribution and the background fit.The BumpHunter -value is defined as the smallest observed probability for the data in a given window to deviate from the background prediction by the observed amount due to a Poissonian fluctuation of the background, using pseudo-experiments generated from the background prediction.The most significant localized excesses identified by the BumpHunter algorithm are found at 3200 GeV in the  region from 0.24 <  < 0.26 for the  4j spectra with a global significance of 0.53 standard deviations, and at 800 GeV in the  region from 0.26 <  < 0.28 for the ⟨ 2j ⟩ spectra with a global significance of 1.98 standard deviations.[GeV]

Cross-section limits
As no signal is observed, limits can be placed on the range of possible production cross-sections for the hypothetical  and  bosons.
The numbers of signal and background events are estimated from maximum-likelihood fits of the signalplus-background models to the corresponding  4j and ⟨ 2j ⟩ distributions.Systematic uncertainties described in Section 3.5 are included in the fits via nuisance parameters constrained by Gaussian penalty terms.The -value is determined from a profile-likelihood-ratio test statistic [100].The local -value for compatibility with the background-only hypothesis when testing a given signal hypothesis ( 0 ) is evaluated based on the asymptotic approximation.Global significance values are computed from background-only pseudo-experiments to account for the trial factors due to scanning both the signal mass and the width hypotheses.The expected and observed 95% confidence level (CL) exclusion limits on the product of the cross-section, branching ratio and acceptance are computed using a modified frequentist-approach CL s [101], in an asymptotic approximation to the test-statistic distribution.
Figure 5 shows the 95% CL upper limits on the allowed cross-sections of these particles as a function of their mass, derived using the signal templates used to optimize the analysis for two representative  regions.
Results are interpolated linearly in the logarithm of the cross-section.Similar results are shown in Figure 6 for the Gaussian signal templates with 5%, 10%, and 15% signal widths.A summary of the limits for all generated signal masses is shown in Figure 7 for   and   as a function of .Overall, the limits are smooth as a function of the mass and , and flatten out in the high-mass region where the background estimation predicts significantly less than one event.
The relative contribution of statistical and systematic uncertainties (see Section 3.5) on the final analysis sensitivity was assessed by repeating the limit-setting procedure while including only statistical sources of uncertainty.The analysis sensitivity was not observed to significantly differ during this test.
To better illustrate different types of events passing the event selection, two event displays are shown in Figure 8.The first shows the event with the highest four-jet mass, with a value of  4j = 6.6 TeV and corresponding ⟨ 2j ⟩ = 2.2 TeV, while the second shows the event with the highest- T fourth jet that is selected ( 4j = 5.2 TeV, ⟨ 2j ⟩ = 0.90 TeV).
2000 3000 4000 5000 6000 7000 8000 9000 10000 =13 TeV, 140 fb s < 0.26 α 0.24 < 95% CL upper limits Observed: Figure 5: The expected and observed 95% confidence exclusion limits on the signal cross-section times acceptance (), efficiency (), and branching ratio (BR) as a function of (a, c, e)   and (b, d, f)   using the signal templates and a 4-parameter fit function for (a, b) 0.24 <  < 0.26, (c, d) 0.26 <  < 0.28, and (e, f) 0.32 <  < 0.34.Observed and expected limits are indicated with markers or a dashed line, respectively.The shaded bands around the expected limit indicates the (darker band) 1 and (lighter band) 2 uncertainty range.Figure 6: The expected and observed limits on the signal cross-section times acceptance (), efficiency (), and branching ratio (BR) as a function of (a, c, e)   and (b, d, f)   for Gaussian signal templates using a 4-parameter fit function for (a, b) 0.24 <  < 0.26, (c, d) 0.26 <  < 0.28, and (e, f) 0.32 <  < 0.34.Observed and expected limits corresponding to different choices of template widths (5%, 10% and 15%) are indicated as different sets of markers or line styles, respectively.The shaded bands around the expected limit for templates with 5% width indicate the (darker band) 1 and (lighter band) 2 uncertainty range.

Conclusion
A search for the production of a generic massive resonance  that decays into two pairs of intermediate resonances , each decaying into two jets, is performed using 140 fb −1 of proton-proton collisions with √  = 13 TeV collected with the ATLAS Detector during Run 2 of the LHC.Such a resonant signal in multĳet events could be manifested in many models of physics beyond the Standard Model, including in well-motivated models of particle dark matter and models with large extra spatial dimensions.
A data-driven background estimate is obtained by fitting these invariant mass distributions with a functional form.The tetrajet system and average dĳet system invariant masses are then studied using the BumpHunter algorithm.No significant excess of events beyond the Standard Model expectation is observed.The most significant localized excesses are found at 3200 GeV in the  region from 0.24 <  < 0.26 for the  4j spectra (global significance of 0.53 standard deviations), and at 800 GeV in the  region from 0.26 <  < 0.28 for the ⟨ 2j ⟩ spectra (global significance of 1.98 standard deviations).The highest tetrajet invariant mass observed is  4j = 6.6 TeV, with a corresponding ⟨ 2j ⟩ value of 2.2 TeV.Using the observed data, upper limits are set on the production cross-sections of new physics scenarios as a function of the  and  masses in both the model-dependent and model-independent interpretations.Data distributions from this search are openly available on the HEPData platform [102] for use in future reinterpretations.
[18] ATLAS Collaboration, Performance of missing transverse momentum reconstruction with the ATLAS detector using proton-proton collisions at        [

Figure 1 :
Figure 1: A schematic representation of the signal topology studied in this analysis: a massive new particle  decays into two new particles with intermediate mass , each decaying into a dĳet system.
) and 3(b) after the event selection for the Run 2 data and for a simulated signal sample with   = 6 TeV and   = 2 TeV, respectively.As shown, requirements on ⟨ 2j ⟩ are correlated with  4j and therefore they can sculpt the background mass distribution.Figures3(c)-3(f) illustrate the correlation between  4j , ⟨ 2j ⟩ and .

Figure 4 :
Figure 4: The (a, c, e) tetrajet and (b, d, f) average dĳet invariant mass distributions in data are shown, along with the fitted background estimates for (a, b) 0.24 <  < 0.26, (c, d) 0.26 <  < 0.28, and (e, f) 0.32 <  < 0.34.The bottom panel of each figure illustrates the fit residuals in terms of standard deviations ().

Figure 7 :Figure 8 :
Figure 7: The observed 95% confidence exclusion limits on the signal cross-section times acceptance (), efficiency (), and branching ratio (BR) as a function of (a)   , and (b)   for signals templates using a 4-parameter fit function.

Table 1 :
Table of selections for the minimum ⟨ 2j ⟩ and  4j values considered when selecting events for a given  bin.These selections are based on a requirement that the single-jet triggers used in the search are at least 99.5% efficient in the selected region.
).The ratio of these quantities,  = ⟨ 2j ⟩/ 4j , is used to parameterize the kinematic space studied in this search in terms of the Lorentz boost of the  decay products.The correlations between  4j and ⟨ 2j ⟩ are shown in Figures3(a