Rivet analyses referenceATLAS_2012_I1126136SUSY top partner search in jets with missing transverse momentumExperiment: ATLAS (LHC) Inspire ID: 1126136 Status: UNVALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
Search for direct pair production of supersymmetric top squarks, assuming the stop_1 decays into a top quark and the lightest supersymmetric particle, and that both top quarks decay to purely hadronic final states. This search has an integrated luminosity of 4.7 fb$^{-1}$ at $\sqrt{s} = 8$\,TeV. Source code: ATLAS_2012_I1126136.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/ChargedFinalState.hh"
5#include "Rivet/Projections/VisibleFinalState.hh"
6#include "Rivet/Projections/IdentifiedFinalState.hh"
7#include "Rivet/Projections/VetoedFinalState.hh"
8#include "Rivet/Projections/FastJets.hh"
9
10namespace Rivet {
11
12
13 class ATLAS_2012_I1126136 : public Analysis {
14 public:
15
16 /// Constructor
17 ATLAS_2012_I1126136()
18 : Analysis("ATLAS_2012_I1126136")
19 { }
20
21
22 /// @name Analysis methods
23 /// @{
24
25 /// Book histograms and initialize projections before the run
26 void init() {
27
28 // projection to find the electrons
29 IdentifiedFinalState elecs(Cuts::abseta < 2.47 && Cuts::pT > 20*GeV);
30 elecs.acceptIdPair(PID::ELECTRON);
31 declare(elecs, "elecs");
32
33 // projection to find the muons
34 IdentifiedFinalState muons(Cuts::abseta < 2.4 && Cuts::pT > 10*GeV);
35 muons.acceptIdPair(PID::MUON);
36 declare(muons, "muons");
37
38 // Jet finder
39 VetoedFinalState vfs;
40 vfs.addVetoPairId(PID::MUON);
41 declare(FastJets(vfs, JetAlg::ANTIKT, 0.4), "AntiKtJets04");
42
43 // for pTmiss
44 declare(VisibleFinalState(Cuts::abseta < 4.9),"vfs");
45
46 // Book histograms
47 book(_count_SR_A ,"count_SR_A" , 1, 0., 1.);
48 book(_count_SR_B ,"count_SR_B" , 1, 0., 1.);
49
50 book(_hist_mjjj1 ,"hist_mjjj1" , 30 , 0. , 600. );
51 book(_hist_mjjj2 ,"hist_mjjj2" , 30 , 0. , 600. );
52 book(_hist_ETmiss ,"hist_ETmiss", 20 , 100. , 600. );
53 book(_hist_mT2 ,"hist_mT2" , 200, 0. , 1000. );
54
55 }
56
57 /// Perform the per-event analysis
58 void analyze(const Event& event) {
59
60 // pTmiss
61 FourMomentum pTmiss;
62 for (const Particle& p : apply<VisibleFinalState>(event, "vfs").particles() ) {
63 pTmiss -= p.momentum();
64 }
65 double ETmiss = pTmiss.pT();
66
67 // require eTmiss > 150
68 if (ETmiss < 150*GeV) vetoEvent;
69
70 // get the candiate jets
71 Jets cand_jets = apply<FastJets>(event, "AntiKtJets04").jetsByPt(Cuts::pT > 20*GeV && Cuts::abseta < 4.5);
72
73 // find the electrons
74 Particles cand_e;
75 for( const Particle& e : apply<IdentifiedFinalState>(event, "elecs").particlesByPt()) {
76 // remove any leptons within 0.4 of any candidate jets
77 bool e_near_jet = false;
78 for ( const Jet& jet : cand_jets ) {
79 double dR = deltaR(e, jet);
80 if (inRange(dR, 0.2, 0.4)) {
81 e_near_jet = true;
82 break;
83 }
84 }
85 if ( e_near_jet ) continue;
86 cand_e.push_back(e);
87 }
88
89 // find the muons
90 Particles cand_mu;
91 for( const Particle& mu : apply<IdentifiedFinalState>(event, "muons").particlesByPt()) {
92 // remove any leptons within 0.4 of any candidate jets
93 bool mu_near_jet = false;
94 for ( const Jet& jet : cand_jets ) {
95 if ( deltaR(mu, jet) < 0.4 ) {
96 mu_near_jet = true;
97 break;
98 }
99 }
100 if ( mu_near_jet ) continue;
101 cand_mu.push_back(mu);
102 }
103
104 // veto events with leptons
105 if( ! cand_e.empty() || ! cand_mu.empty() )
106 vetoEvent;
107
108 // discard jets that overlap with electrons
109 Jets recon_jets;
110 for ( const Jet& jet : cand_jets ) {
111 if (jet.abseta() > 2.8 || jet.pT() < 30*GeV) continue;
112 bool away_from_e = true;
113 for (const Particle& e : cand_e ) {
114 if ( deltaR(e, jet) < 0.2 ) {
115 away_from_e = false;
116 break;
117 }
118 }
119 if ( away_from_e ) recon_jets.push_back( jet );
120 }
121
122 // find b jets
123 Jets tight_bjets,loose_bjets;
124 for(const Jet& jet : recon_jets) {
125 /// @todo Should be abseta?
126 if (!jet.bTagged() && jet.eta()>2.5) continue;
127 double prob = rand()/static_cast<double>(RAND_MAX);
128 if (prob <= 0.60) tight_bjets.push_back(jet);
129 if (prob <= 0.75) loose_bjets.push_back(jet);
130 }
131
132 // require >=1 tight or >=2 loose b-jets
133 if (! ( !tight_bjets.empty() || loose_bjets.size()>=2) )
134 vetoEvent;
135
136 // must be at least 6 jets with pT>30
137 if (recon_jets.size()<6 ) vetoEvent;
138
139 // hardest > 130
140 if (recon_jets[0].perp() < 130. ) vetoEvent;
141
142 // three hardest jets must be separated from etmiss
143 for (unsigned int ix=0;ix<3;++ix) {
144 if (deltaPhi(recon_jets[ix].momentum(),pTmiss)<0.2*PI)
145 vetoEvent;
146 }
147
148 // remove events with tau like jets
149 for (unsigned int ix=3;ix<recon_jets.size();++ix) {
150 // skip jets seperated from eTmiss
151 if (deltaPhi(recon_jets[ix].momentum(),pTmiss)>=0.2*PI)
152 continue;
153 // check the number of tracks between 1 and 4
154 unsigned int ncharged=0;
155 for ( const Particle & particle : recon_jets[ix].particles()) {
156 if (PID::charge3(particle.pid())!=0) ++ncharged;
157 }
158 if (ncharged==0 || ncharged>4) continue;
159 // calculate transverse mass and reject if < 100
160 double mT = 2.*recon_jets[ix].perp()*ETmiss
161 -recon_jets[ix].px()*pTmiss.px()
162 -recon_jets[ix].py()*pTmiss.py();
163 if (mT<100.) vetoEvent;
164 }
165
166 // if 2 loose b-jets apply mT cut
167 if (loose_bjets.size()>=2) {
168 // find b-jet closest to eTmiss
169 double minR(1e30);
170 unsigned int ijet(0);
171 for(unsigned int ix=0;ix<loose_bjets.size();++ix) {
172 double dR = deltaR(loose_bjets[ix].momentum(),pTmiss);
173 if(dR<minR) {
174 minR=dR;
175 ijet = ix;
176 }
177 }
178 double mT = 2.*loose_bjets[ijet].perp()*ETmiss
179 -loose_bjets[ijet].px()*pTmiss.px()
180 -loose_bjets[ijet].py()*pTmiss.py();
181 if(mT<170.) vetoEvent;
182 }
183
184 // 1 tight b-jet apply mT cut
185 if(tight_bjets.size()==1) {
186 for(unsigned int ix=0;ix<4;++ix) {
187 double mT = 2.*recon_jets[ix].perp()*ETmiss
188 -recon_jets[ix].px()*pTmiss.px()
189 -recon_jets[ix].py()*pTmiss.py();
190 if(mT<175.) vetoEvent;
191 }
192 }
193
194 // find the closest triplet of jets in (eta,phi)
195 unsigned int j1(0),j2(0),j3(0);
196 double minR2(1e30);
197 for(unsigned int i1=0;i1<recon_jets.size();++i1) {
198 for(unsigned int i2=i1+1;i2<recon_jets.size();++i2) {
199 for(unsigned int i3=i2+1;i3<recon_jets.size();++i3) {
200 double delR2 =
201 sqr(deltaR(recon_jets[i1].momentum(),recon_jets[i2].momentum())) +
202 sqr(deltaR(recon_jets[i1].momentum(),recon_jets[i3].momentum())) +
203 sqr(deltaR(recon_jets[i2].momentum(),recon_jets[i3].momentum()));
204 if(delR2<minR2) {
205 minR2=delR2;
206 j1=i1;
207 j2=i2;
208 j3=i3;
209 }
210 }
211 }
212 }
213 // 4-momentum and mass of first triplet
214 FourMomentum pjjj1 = recon_jets[j1].momentum() +
215 recon_jets[j2].momentum()+ recon_jets[j3].momentum();
216 double mjjj1 = pjjj1.mass();
217
218 // find the second triplet
219 unsigned int j4(0),j5(0),j6(0);
220 minR2=0.;
221 for(unsigned int i1=0;i1<recon_jets.size();++i1) {
222 if(i1==j1||i1==j2||i1==j3) continue;
223 for(unsigned int i2=i1+1;i2<recon_jets.size();++i2) {
224 if(i2==j1||i2==j2||i2==j3) continue;
225 for(unsigned int i3=i2+1;i3<recon_jets.size();++i3) {
226 if(i3==j1||i3==j2||i3==j3) continue;
227 double delR2 =
228 sqr(deltaR(recon_jets[i1].momentum(),recon_jets[i2].momentum())) +
229 sqr(deltaR(recon_jets[i1].momentum(),recon_jets[i3].momentum())) +
230 sqr(deltaR(recon_jets[i2].momentum(),recon_jets[i3].momentum()));
231 if(delR2<minR2) {
232 minR2=delR2;
233 j4=i1;
234 j5=i2;
235 j6=i3;
236 }
237 }
238 }
239 }
240
241 // 4-momentum and mass of first triplet
242 FourMomentum pjjj2 = recon_jets[j4].momentum() +
243 recon_jets[j5].momentum()+ recon_jets[j6].momentum();
244 double mjjj2 = pjjj2.mass();
245
246 _hist_mjjj1->fill(mjjj1);
247 _hist_mjjj2->fill(mjjj2);
248 // require triplets in 80<mjjj<270
249 if(mjjj1<80.||mjjj1>270.||mjjj2<80.||mjjj2>270.)
250 vetoEvent;
251
252 // counts in signal regions
253 _count_SR_A->fill(0.5);
254 if(ETmiss>260.) _count_SR_B->fill(0.5);
255
256 _hist_ETmiss->fill(ETmiss);
257 const double m_T2 = mT2(pjjj1, pjjj2, pTmiss, 0.0); // zero mass invisibles
258 _hist_mT2->fill(m_T2);
259 }
260 /// @}
261
262
263 void finalize() {
264
265 double norm = 4.7* crossSection()/sumOfWeights()/femtobarn;
266 scale(_count_SR_A , norm );
267 scale(_count_SR_B , norm );
268 scale(_hist_mjjj1 , 20.*norm );
269 scale(_hist_ETmiss, 50.*norm );
270 scale(_hist_mjjj2 , 20.*norm );
271 scale(_hist_mT2 , norm );
272
273 }
274
275 private:
276
277 /// @name Histograms
278 /// @{
279 Histo1DPtr _count_SR_A;
280 Histo1DPtr _count_SR_B;
281
282 Histo1DPtr _hist_mjjj1;
283 Histo1DPtr _hist_mjjj2;
284 Histo1DPtr _hist_ETmiss;
285 Histo1DPtr _hist_mT2;
286 /// @}
287
288 };
289
290 RIVET_DECLARE_PLUGIN(ATLAS_2012_I1126136);
291
292}
|