rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2012_I1112263

3 lepton plus missing transverse energy SUSY search
Experiment: ATLAS (LHC)
Inspire ID: 1112263
Status: VALIDATED
Authors:
  • Peter Richardson
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • BSM signal events at 7000 GeV.

Search for SUSY using events with 3 leptons in association with missing transverse energy in proton-proton collisions at a centre-of-mass energy of 7 TeV. The data sample has a total integrated luminosity of 2.06 fb$^{-1}$. There is no reference data and in addition to the control plots from the paper the number of events in the two signal regions, correctly normalized to an integrated luminosity 2.06 fb$^{-1}$, are calculated.

Source code: ATLAS_2012_I1112263.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/VetoedFinalState.hh"
  7#include "Rivet/Projections/IdentifiedFinalState.hh"
  8#include "Rivet/Projections/FastJets.hh"
  9#include "Rivet/Tools/RivetMT2.hh"
 10
 11namespace Rivet {
 12
 13
 14  /// @author Peter Richardson
 15  class ATLAS_2012_I1112263 : public Analysis {
 16  public:
 17
 18    /// Constructor
 19    ATLAS_2012_I1112263()
 20      : Analysis("ATLAS_2012_I1112263")
 21    {    }
 22
 23
 24    /// @name Analysis methods
 25    //@{
 26
 27    /// Book histograms and initialise projections before the run
 28    void init() {
 29
 30      // projection to find the electrons
 31      IdentifiedFinalState elecs(Cuts::abseta < 2.47 && Cuts::pT > 10*GeV);
 32      elecs.acceptIdPair(PID::ELECTRON);
 33      declare(elecs, "elecs");
 34
 35      // projection to find the muons
 36      IdentifiedFinalState muons(Cuts::abseta < 2.4 && Cuts::pT > 10*GeV);
 37      muons.acceptIdPair(PID::MUON);
 38      declare(muons, "muons");
 39
 40      // for pTmiss
 41      declare(VisibleFinalState(Cuts::abseta < 4.9),"vfs");
 42
 43      VetoedFinalState vfs;
 44      vfs.addVetoPairId(PID::MUON);
 45
 46      /// Jet finder
 47      declare(FastJets(vfs, FastJets::ANTIKT, 0.4), "AntiKtJets04");
 48
 49      // all tracks (to do deltaR with leptons)
 50      declare(ChargedFinalState(Cuts::abseta < 3.0),"cfs");
 51
 52      // Book histograms
 53      {Histo1DPtr tmp; _hist_leptonpT_SR1.push_back(book(tmp,"hist_lepton_pT_1_SR1",11,0.,220.));}
 54      {Histo1DPtr tmp; _hist_leptonpT_SR1.push_back(book(tmp,"hist_lepton_pT_2_SR1", 7,0.,140.));}
 55      {Histo1DPtr tmp; _hist_leptonpT_SR1.push_back(book(tmp,"hist_lepton_pT_3_SR1", 8,0.,160.));}
 56      {Histo1DPtr tmp; _hist_leptonpT_SR2.push_back(book(tmp,"hist_lepton_pT_1_SR2",11,0.,220.));}
 57      {Histo1DPtr tmp; _hist_leptonpT_SR2.push_back(book(tmp,"hist_lepton_pT_2_SR2", 7,0.,140.));}
 58      {Histo1DPtr tmp; _hist_leptonpT_SR2.push_back(book(tmp,"hist_lepton_pT_3_SR2", 8,0.,160.));}
 59      book(_hist_etmiss_SR1_A ,"hist_etmiss_SR1_A",15,10.,310.);
 60      book(_hist_etmiss_SR1_B ,"hist_etmiss_SR1_B", 9,10.,190.);
 61      book(_hist_etmiss_SR2_A ,"hist_etmiss_SR2_A",15,10.,310.);
 62      book(_hist_etmiss_SR2_B ,"hist_etmiss_SR2_B", 9,10.,190.);
 63      book(_hist_mSFOS,"hist_mSFOF",9,0.,180.);
 64
 65      book(_count_SR1 ,"count_SR1", 1, 0., 1.);
 66      book(_count_SR2 ,"count_SR2", 1, 0., 1.);
 67    }
 68
 69
 70    /// Perform the per-event analysis
 71    void analyze(const Event& event) {
 72
 73      // Get the jet candidates
 74      Jets cand_jets;
 75      for (const Jet& jet : apply<FastJets>(event, "AntiKtJets04").jetsByPt(20.0*GeV) ) {
 76        if ( fabs( jet.eta() ) < 2.8 ) {
 77          cand_jets.push_back(jet);
 78        }
 79      }
 80
 81      // Candidate muons
 82      Particles cand_mu;
 83      Particles chg_tracks =
 84        apply<ChargedFinalState>(event, "cfs").particles();
 85      for ( const Particle & mu : apply<IdentifiedFinalState>(event, "muons").particlesByPt() ) {
 86        double pTinCone = -mu.pT();
 87        for ( const Particle & track : chg_tracks ) {
 88          if ( deltaR(mu.momentum(),track.momentum()) <= 0.2 )
 89            pTinCone += track.pT();
 90        }
 91        if ( pTinCone < 1.8*GeV )
 92          cand_mu.push_back(mu);
 93      }
 94
 95      // Candidate electrons
 96      Particles cand_e;
 97      for ( const Particle & e : apply<IdentifiedFinalState>(event, "elecs").particlesByPt() ) {
 98        double eta = e.eta();
 99        // Remove electrons with pT<15 in old veto region
100        // (NOT EXPLICIT IN THIS PAPER BUT IN SIMILAR 4 LEPTON PAPER and THIS DESCRPITION
101        //  IS MUCH WORSE SO ASSUME THIS IS DONE)
102        if ( fabs(eta)>1.37 && fabs(eta) < 1.52 && e.perp()< 15.*GeV)
103          continue;
104        double pTinCone = -e.perp();
105        for ( const Particle & track : chg_tracks ) {
106          if ( deltaR(e.momentum(),track.momentum()) <= 0.2 )
107            pTinCone += track.pT();
108        }
109        if (pTinCone/e.perp()<0.1) {
110          cand_e.push_back(e);
111        }
112      }
113
114      // Resolve jet/lepton ambiguity
115      // (NOT EXPLICIT IN THIS PAPER BUT IN SIMILAR 4 LEPTON PAPER and THIS DESCRPITION
116      //  IS MUCH WORSE SO ASSUME THIS IS DONE)
117      Jets recon_jets;
118      for ( const Jet& jet : cand_jets ) {
119        bool away_from_e = true;
120        for ( const Particle & e : cand_e ) {
121          if ( deltaR(e.momentum(),jet.momentum()) <= 0.2 ) {
122            away_from_e = false;
123            break;
124          }
125        }
126        if ( away_from_e )
127          recon_jets.push_back( jet );
128      }
129
130      // Only keep electrons more than R=0.4 from jets
131      Particles recon_e;
132      for ( const Particle & e : cand_e ) {
133        bool away = true;
134        for ( const Jet& jet : recon_jets ) {
135          if ( deltaR(e.momentum(),jet.momentum()) < 0.4 ) {
136            away = false;
137            break;
138          }
139        }
140        // ... and 0.1 from any muons
141        if ( ! away ) {
142          for ( const Particle & mu : cand_e ) {
143            if ( deltaR(mu.momentum(),e.momentum()) < 0.1 ) {
144              away = false;
145              break;
146            }
147          }
148        }
149        if ( away )
150          recon_e.push_back( e );
151      }
152      // Only keep muons more than R=0.4 from jets
153      Particles recon_mu;
154      for ( const Particle & mu : cand_mu ) {
155        bool away = true;
156        for ( const Jet& jet : recon_jets ) {
157          if ( deltaR(mu.momentum(),jet.momentum()) < 0.4 ) {
158            away = false;
159            break;
160          }
161        }
162        // ... and 0.1 from any electrona
163        if ( ! away ) {
164          for ( const Particle & e : cand_e ) {
165            if ( deltaR(mu.momentum(),e.momentum()) < 0.1 ) {
166              away = false;
167              break;
168            }
169          }
170        }
171        if ( away )
172          recon_mu.push_back( mu );
173      }
174
175      // pTmiss
176      Particles vfs_particles =
177        apply<VisibleFinalState>(event, "vfs").particles();
178      FourMomentum pTmiss;
179      for ( const Particle & p : vfs_particles ) {
180        pTmiss -= p.momentum();
181      }
182      double eTmiss = pTmiss.pT();
183
184      // Now only use recon_jets, recon_mu, recon_e
185
186      // Reject events with wrong number of leptons
187      if ( recon_mu.size() + recon_e.size() != 3 ) {
188        MSG_DEBUG("To few charged leptons left after selection");
189        vetoEvent;
190      }
191
192      // ATLAS calo problem
193      if (rand()/static_cast<double>(RAND_MAX) <= 0.42) {
194        for ( const Particle & e : recon_e ) {
195          double eta = e.eta();
196          double phi = e.azimuthalAngle(MINUSPI_PLUSPI);
197          if (inRange(eta, -0.1, 1.5) && inRange(phi, -0.9, -0.5)) vetoEvent;
198        }
199        for ( const Jet & jet : recon_jets ) {
200          const double eta = jet.rapidity();
201          const double phi = jet.azimuthalAngle(MINUSPI_PLUSPI);
202          if (jet.perp() > 40*GeV && inRange(eta, -0.1, 1.5) && inRange(phi, -0.9, -0.5)) vetoEvent;
203        }
204      }
205
206      if ( !( !recon_e .empty() && recon_e[0] .perp() > 25*GeV) &&
207           !( !recon_mu.empty() && recon_mu[0].perp() > 20*GeV) ) {
208        MSG_DEBUG("Hardest lepton fails trigger");
209        vetoEvent;
210      }
211
212      // eTmiss cut
213      if (eTmiss < 50*GeV) vetoEvent;
214
215      // Check at least 1 SFOS pair
216      double mSFOS=1e30, mdiff=1e30*GeV;
217      size_t nSFOS=0;
218      for (size_t ix = 0; ix < recon_e.size(); ++ix) {
219        for (size_t iy = ix+1; iy < recon_e.size(); ++iy) {
220          if (recon_e[ix].pid()*recon_e[iy].pid() > 0) continue;
221          ++nSFOS;
222          double mtest = (recon_e[ix].momentum() + recon_e[iy].momentum()).mass();
223          // Veto is mass<20
224          if (mtest < 20*GeV) vetoEvent;
225          if (fabs(mtest - 90*GeV) < mdiff) {
226            mSFOS = mtest;
227            mdiff = fabs(mtest - 90*GeV);
228          }
229        }
230      }
231      for (size_t ix = 0; ix < recon_mu.size(); ++ix) {
232        for (size_t iy = ix+1; iy < recon_mu.size(); ++iy) {
233          if (recon_mu[ix].pid()*recon_mu[iy].pid() > 0) continue;
234          ++nSFOS;
235          double mtest = (recon_mu[ix].momentum() + recon_mu[iy].momentum()).mass();
236          // Veto is mass < 20*GeV
237          if (mtest < 20*GeV) vetoEvent;
238          if (fabs(mtest - 90*GeV) < mdiff) {
239            mSFOS = mtest;
240            mdiff = fabs(mtest - 90*GeV);
241          }
242        }
243      }
244      // Require at least 1 SFOS pair
245      if (nSFOS == 0) vetoEvent;
246      // b-jet veto in SR!
247      if (mdiff > 10*GeV) {
248        for (const Jet & jet : recon_jets ) {
249          if (jet.bTagged() && rand()/static_cast<double>(RAND_MAX) <= 0.60) vetoEvent;
250        }
251      }
252
253	  // Histogram filling
254	  const double weight = 1.0;
255
256      // Region SR1, Z depleted
257      if (mdiff > 10*GeV) {
258        _count_SR1->fill(0.5, weight);
259        _hist_etmiss_SR1_A->fill(eTmiss, weight);
260        _hist_etmiss_SR1_B->fill(eTmiss, weight);
261        _hist_mSFOS->fill(mSFOS, weight);
262      }
263      // Region SR2, Z enriched
264      else {
265        _count_SR2->fill(0.5, weight);
266        _hist_etmiss_SR2_A->fill(eTmiss, weight);
267        _hist_etmiss_SR2_B->fill(eTmiss, weight);
268      }
269      // Make the control plots
270      // lepton pT
271      size_t ie=0, imu=0;
272      for (size_t ix = 0; ix < 3; ++ix) {
273        Histo1DPtr hist = (mdiff > 10*GeV) ? _hist_leptonpT_SR1[ix] :  _hist_leptonpT_SR2[ix];
274        double pTe  = (ie  < recon_e .size()) ? recon_e [ie ].perp() : -1*GeV;
275        double pTmu = (imu < recon_mu.size()) ? recon_mu[imu].perp() : -1*GeV;
276        if (pTe > pTmu) {
277          hist->fill(pTe, weight);
278          ++ie;
279        } else {
280          hist->fill(pTmu, weight);
281          ++imu;
282        }
283      }
284
285    }
286
287
288    void finalize() {
289      const double norm = crossSection()/femtobarn*2.06/sumOfWeights();
290      // These are number of events at 2.06fb^-1 per 20 GeV
291      for (size_t ix = 0; ix < 3; ++ix) {
292        scale(_hist_leptonpT_SR1[ix], norm*20.);
293        scale(_hist_leptonpT_SR2[ix], norm*20.);
294      }
295      scale(_hist_etmiss_SR1_A, norm*20.);
296      scale(_hist_etmiss_SR1_B, norm*20.);
297      scale(_hist_etmiss_SR2_A, norm*20.);
298      scale(_hist_etmiss_SR2_B, norm*20.);
299      scale(_hist_mSFOS, norm*20.);
300      // These are number of events at 2.06fb^-1
301      scale(_count_SR1, norm);
302      scale(_count_SR2, norm);
303    }
304
305    //@}
306
307
308  private:
309
310  /// @name Histograms
311  //@{
312  vector<Histo1DPtr> _hist_leptonpT_SR1, _hist_leptonpT_SR2;
313  Histo1DPtr _hist_etmiss_SR1_A, _hist_etmiss_SR1_B, _hist_etmiss_SR2_A, _hist_etmiss_SR2_B;
314  Histo1DPtr _hist_mSFOS, _count_SR1, _count_SR2;
315  //@}
316
317  };
318
319  // The hook for the plugin system
320  RIVET_DECLARE_PLUGIN(ATLAS_2012_I1112263);
321
322}