rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2012_I1126136

SUSY top partner search in jets with missing transverse momentum
Experiment: ATLAS (LHC)
Inspire ID: 1126136
Status: UNVALIDATED
Authors:
  • Peter Richardson
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • BSM signal events at 7000 GeV.

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, FastJets::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      const double weight = 1.0;
 60
 61      // pTmiss
 62      FourMomentum pTmiss;
 63      for (const Particle& p : apply<VisibleFinalState>(event, "vfs").particles() ) {
 64        pTmiss -= p.momentum();
 65      }
 66      double ETmiss = pTmiss.pT();
 67
 68      // require eTmiss > 150
 69      if (ETmiss < 150*GeV) vetoEvent;
 70
 71      // get the candiate jets
 72      Jets cand_jets;
 73      for ( const Jet& jet : apply<FastJets>(event, "AntiKtJets04").jetsByPt(20.0*GeV) ) {
 74        if (jet.abseta() < 4.5) cand_jets.push_back(jet);
 75      }
 76
 77      // find the electrons
 78      Particles cand_e;
 79      for( const Particle& e : apply<IdentifiedFinalState>(event, "elecs").particlesByPt()) {
 80        // remove any leptons within 0.4 of any candidate jets
 81        bool e_near_jet = false;
 82        for ( const Jet& jet : cand_jets ) {
 83          double dR = deltaR(e, jet);
 84          if (inRange(dR, 0.2, 0.4)) {
 85            e_near_jet = true;
 86            break;
 87          }
 88        }
 89        if ( e_near_jet ) continue;
 90        cand_e.push_back(e);
 91      }
 92
 93      // find the muons
 94      Particles cand_mu;
 95      for( const Particle& mu : apply<IdentifiedFinalState>(event, "muons").particlesByPt()) {
 96        // remove any leptons within 0.4 of any candidate jets
 97        bool mu_near_jet = false;
 98        for ( const Jet& jet : cand_jets ) {
 99          if ( deltaR(mu, jet) < 0.4 ) {
100            mu_near_jet = true;
101            break;
102          }
103        }
104        if ( mu_near_jet ) continue;
105        cand_mu.push_back(mu);
106      }
107
108      // veto events with leptons
109      if( ! cand_e.empty() || ! cand_mu.empty() )
110        vetoEvent;
111
112      // discard jets that overlap with electrons
113      Jets recon_jets;
114      for ( const Jet& jet : cand_jets ) {
115        if (jet.abseta() > 2.8 || jet.pT() < 30*GeV) continue;
116        bool away_from_e = true;
117        for (const Particle& e : cand_e ) {
118          if ( deltaR(e, jet) < 0.2 ) {
119            away_from_e = false;
120            break;
121          }
122        }
123        if ( away_from_e ) recon_jets.push_back( jet );
124      }
125
126      // find b jets
127      Jets tight_bjets,loose_bjets;
128      for(const Jet& jet : recon_jets) {
129        /// @todo Should be abseta?
130        if (!jet.bTagged() && jet.eta()>2.5) continue;
131        double prob = rand()/static_cast<double>(RAND_MAX);
132        if (prob <= 0.60) tight_bjets.push_back(jet);
133        if (prob <= 0.75) loose_bjets.push_back(jet);
134      }
135
136      // require >=1 tight or >=2 loose b-jets
137      if (! ( !tight_bjets.empty() || loose_bjets.size()>=2) )
138        vetoEvent;
139
140      // must be at least 6 jets with pT>30
141      if (recon_jets.size()<6 ) vetoEvent;
142
143      // hardest > 130
144      if (recon_jets[0].perp() < 130. ) vetoEvent;
145
146      // three hardest jets must be separated from etmiss
147      for (unsigned int ix=0;ix<3;++ix) {
148        if (deltaPhi(recon_jets[ix].momentum(),pTmiss)<0.2*PI)
149          vetoEvent;
150      }
151
152      // remove events with tau like jets
153      for (unsigned int ix=3;ix<recon_jets.size();++ix) {
154        // skip jets seperated from eTmiss
155        if (deltaPhi(recon_jets[ix].momentum(),pTmiss)>=0.2*PI)
156          continue;
157        // check the number of tracks between 1 and 4
158        unsigned int ncharged=0;
159        for ( const Particle & particle : recon_jets[ix].particles()) {
160          if (PID::charge3(particle.pid())!=0) ++ncharged;
161        }
162        if (ncharged==0 || ncharged>4) continue;
163        // calculate transverse mass and reject if < 100
164        double mT = 2.*recon_jets[ix].perp()*ETmiss
165          -recon_jets[ix].px()*pTmiss.px()
166          -recon_jets[ix].py()*pTmiss.py();
167        if (mT<100.) vetoEvent;
168      }
169
170      // if 2 loose b-jets apply mT cut
171      if (loose_bjets.size()>=2) {
172        // find b-jet closest to eTmiss
173        double minR(1e30);
174        unsigned int ijet(0);
175        for(unsigned int ix=0;ix<loose_bjets.size();++ix) {
176          double dR = deltaR(loose_bjets[ix].momentum(),pTmiss);
177          if(dR<minR) {
178            minR=dR;
179            ijet = ix;
180          }
181        }
182        double  mT = 2.*loose_bjets[ijet].perp()*ETmiss
183          -loose_bjets[ijet].px()*pTmiss.px()
184          -loose_bjets[ijet].py()*pTmiss.py();
185        if(mT<170.) vetoEvent;
186      }
187
188      // 1 tight b-jet apply mT cut
189      if(tight_bjets.size()==1) {
190        for(unsigned int ix=0;ix<4;++ix) {
191          double mT = 2.*recon_jets[ix].perp()*ETmiss
192            -recon_jets[ix].px()*pTmiss.px()
193            -recon_jets[ix].py()*pTmiss.py();
194          if(mT<175.) vetoEvent;
195        }
196      }
197
198      // find the closest triplet of jets in (eta,phi)
199      unsigned int j1(0),j2(0),j3(0);
200      double minR2(1e30);
201      for(unsigned int i1=0;i1<recon_jets.size();++i1) {
202        for(unsigned int i2=i1+1;i2<recon_jets.size();++i2) {
203          for(unsigned int i3=i2+1;i3<recon_jets.size();++i3) {
204            double delR2 =
205              sqr(deltaR(recon_jets[i1].momentum(),recon_jets[i2].momentum())) +
206              sqr(deltaR(recon_jets[i1].momentum(),recon_jets[i3].momentum())) +
207              sqr(deltaR(recon_jets[i2].momentum(),recon_jets[i3].momentum()));
208            if(delR2<minR2) {
209              minR2=delR2;
210              j1=i1;
211              j2=i2;
212              j3=i3;
213            }
214          }
215        }
216      }
217      // 4-momentum and mass of first triplet
218      FourMomentum pjjj1 = recon_jets[j1].momentum() +
219        recon_jets[j2].momentum()+ recon_jets[j3].momentum();
220      double mjjj1 = pjjj1.mass();
221
222      // find the second triplet
223      unsigned int j4(0),j5(0),j6(0);
224      minR2=0.;
225      for(unsigned int i1=0;i1<recon_jets.size();++i1) {
226        if(i1==j1||i1==j2||i1==j3) continue;
227        for(unsigned int i2=i1+1;i2<recon_jets.size();++i2) {
228          if(i2==j1||i2==j2||i2==j3) continue;
229          for(unsigned int i3=i2+1;i3<recon_jets.size();++i3) {
230            if(i3==j1||i3==j2||i3==j3) continue;
231            double delR2 =
232              sqr(deltaR(recon_jets[i1].momentum(),recon_jets[i2].momentum())) +
233              sqr(deltaR(recon_jets[i1].momentum(),recon_jets[i3].momentum())) +
234              sqr(deltaR(recon_jets[i2].momentum(),recon_jets[i3].momentum()));
235            if(delR2<minR2) {
236              minR2=delR2;
237              j4=i1;
238              j5=i2;
239              j6=i3;
240            }
241          }
242        }
243      }
244
245      // 4-momentum and mass of first triplet
246      FourMomentum pjjj2 = recon_jets[j4].momentum() +
247        recon_jets[j5].momentum()+ recon_jets[j6].momentum();
248      double mjjj2 = pjjj2.mass();
249
250      _hist_mjjj1->fill(mjjj1,weight);
251      _hist_mjjj2->fill(mjjj2,weight);
252      // require triplets in 80<mjjj<270
253      if(mjjj1<80.||mjjj1>270.||mjjj2<80.||mjjj2>270.)
254        vetoEvent;
255
256      // counts in signal regions
257      _count_SR_A->fill(0.5,weight);
258      if(ETmiss>260.) _count_SR_B->fill(0.5,weight);
259
260      _hist_ETmiss->fill(ETmiss,weight);
261      const double m_T2 = mT2(pjjj1, pjjj2, pTmiss, 0.0); // zero mass invisibles
262      _hist_mT2->fill(m_T2,weight);
263    }
264    //@}
265
266
267    void finalize() {
268
269      double norm = 4.7* crossSection()/sumOfWeights()/femtobarn;
270      scale(_count_SR_A ,     norm );
271      scale(_count_SR_B ,     norm );
272      scale(_hist_mjjj1 , 20.*norm );
273      scale(_hist_ETmiss, 50.*norm );
274      scale(_hist_mjjj2 , 20.*norm );
275      scale(_hist_mT2   ,     norm );
276
277    }
278
279  private:
280
281    /// @name Histograms
282    //@{
283    Histo1DPtr _count_SR_A;
284    Histo1DPtr _count_SR_B;
285
286    Histo1DPtr _hist_mjjj1;
287    Histo1DPtr _hist_mjjj2;
288    Histo1DPtr _hist_ETmiss;
289    Histo1DPtr _hist_mT2;
290    //@}
291
292  };
293
294  // The hook for the plugin system
295  RIVET_DECLARE_PLUGIN(ATLAS_2012_I1126136);
296
297}