Rivet analyses referenceATLAS_2018_I1667046Search for R-parity-violating SUSY in multi-jet final states at 13 TeVExperiment: ATLAS (LHC) Inspire ID: 1667046 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
This search uses 36.1/fb of data collected by the ATLAS detector in proton-proton collisions with a center-of-mass energy of $\sqrt{s}=13$ TeV at the LHC. The analysis is performed using requirements on the number of jets and the number of jets tagged as containing a $b$-hadron, as well as a topological observable formed by the scalar sum of masses of large-radius jets in the event. No significant excess above the expected Standard Model background was observed. Source code: ATLAS_2018_I1667046.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Tools/Cutflow.hh"
6
7namespace Rivet {
8
9
10 /// Search for R-parity-violating SUSY in multi-jet final states at 13 TeV
11 class ATLAS_2018_I1667046 : public Analysis {
12 public:
13
14 /// Constructor
15 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2018_I1667046);
16
17
18 /// @name Analysis methods
19 //@{
20
21 /// Book histograms and initialise projections before the run
22 void init() {
23
24 // Projections
25 const FinalState fs (Cuts::abseta < 4.9);
26 declare("SmallRJ", FastJets(fs, FastJets::ANTIKT, 0.4));
27 declare("LargeRJ", FastJets(fs, FastJets::ANTIKT, 1.0));
28
29 // Book histograms
30 book(_h_sigmaM, "sigmaM", 50, 200, 2000);
31 book(_h_modeta, "ModEta12", 42, 0, 4.2);
32
33 // Cutflows
34 _flows.addCutflow("CutFlow1",
35 {"NJet >= 4 ", "Delta12 < 1.4", "PJet1 > 400 GeV", "M SumJ > 1.0 ",
36 "NbJet > 0", "M SumJ > 1.0 & NbJet > 0"});
37 _flows.addCutflow("CutFlow2",
38 {"NJet >= 4 ", "Delta12 < 1.4", "NJet >= 5 ", "M SumJ > 0.8 ",
39 "NbJet > 0", "M SumJ > 0.8 & NbJet > 0"});
40 }
41
42
43 /// Perform the per-event analysis
44 void analyze(const Event& event) {
45 _flows.fillinit();
46
47 // Trim large-R jets and apply cuts
48 Jets LRJ_old = apply<FastJets>(event, "LargeRJ").jetsByPt(Cuts::abseta < 4.9);
49 Jets LRJJ;
50 fastjet::Filter trimmer(fastjet::JetDefinition(fastjet::kt_algorithm,0.2), fastjet::SelectorPtFractionMin(0.05));
51 for (Jet& j: LRJ_old) LRJJ.push_back(trimmer(j));
52 Jets LRJ = ifilter_select(LRJJ, Cuts::abseta < 2.0 && Cuts::pT > 200*GeV);
53 if (LRJ.size() < 4) vetoEvent;
54 LRJ = sortByPt(LRJ); // sorting for constructing sigmaM
55
56 // Small R jets need to pass some cuts + need to be BTagged.
57 const Jets SRJ = apply<FastJets>(event, "SmallRJ").jetsByPt(Cuts::abseta < 2.5 && Cuts::pT > 50*GeV);
58 Jets BT_SRJ = filter_select(SRJ, hasBTag(Cuts::pT > 5*GeV));
59 const int tagg = (BT_SRJ.size() == 0) ? 1 : 0; //if there are no B-TAGGED jets are present, tagg = 1.
60
61 // Now to find B-MATCHED Large R Jets: extra step, not useful for SR regions!
62 const Jets BM_LRJ = selectIfAnyDeltaRLess(LRJ, BT_SRJ, 1.0);
63
64
65 // Now to build observables SigmaM and deltaEta
66 // Add mass of leading four large R jets,
67 double sigmaM = 0.0;
68 for (const Jet& j : head(LRJ, 4)) sigmaM += j.mass();
69 _h_sigmaM->fill(sigmaM);
70
71 // Build deta between two leading large-R jets
72 const double delta_eta = fabs(deltaEta(LRJ[0], LRJ[1]));
73 _h_modeta->fill(delta_eta);
74
75 // CutFlow1
76 if (LRJ.size() >= 4) {
77 _flows["CutFlow1"].fill(1);
78 if (delta_eta < 1.4) {
79 _flows["CutFlow1"].fill(2);
80 if (LRJ[0].pT() > 400*GeV) {
81 _flows["CutFlow1"].fill(3);
82 if (sigmaM > 1000*GeV)
83 _flows["CutFlow1"].fill(4);
84 if (tagg == 0) {
85 _flows["CutFlow1"].fill(5);
86 if (sigmaM > 1000*GeV){
87 _flows["CutFlow1"].fill(6);
88 }
89 } //end of btagg loop
90 } //end of pT loop
91 } //end of delta loop
92 } //end of Njet>4 loop
93
94 // CutFlow2
95 if (LRJ.size() >= 4) {
96 _flows["CutFlow2"].fill(1);
97 if (delta_eta < 1.4) {
98 _flows["CutFlow2"].fill(2);
99 if (LRJ.size() >= 5) {
100 _flows["CutFlow2"].fill(3);
101 if (sigmaM > 800*GeV)
102 _flows["CutFlow2"].fill(4);
103 if (tagg == 0) {
104 _flows["CutFlow2"].fill(5);
105 if (sigmaM > 800*GeV){_flows["CutFlow2"].fill(6);
106 }
107 } //end of btagg loop
108 } //Njet>5 loop
109 } //end of delta loop
110 } //Njet>4 loop
111
112 }
113
114
115 /// Normalise histograms etc., after the run
116 void finalize() {
117 const double expected = 36.1*crossSection()/femtobarn;
118 normalize(_h_sigmaM, expected/sumOfWeights());
119 normalize(_h_modeta, expected/sumOfWeights());
120 // _flows.scale(99.7/numEvents());
121 // @todo suppress output in reentrant mode?
122 MSG_INFO(_flows);
123 }
124
125 //@}
126
127
128 private:
129
130 /// @name Histograms
131 Histo1DPtr _h_sigmaM, _h_modeta;
132
133 // Cutflows
134 Cutflows _flows;
135
136 };
137
138
139 // The hook for the plugin system
140 RIVET_DECLARE_PLUGIN(ATLAS_2018_I1667046);
141
142
143}
|