1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/PromptFinalState.hh"
  4#include "Rivet/Projections/MissingMomentum.hh"
  5#include "Rivet/Projections/LeptonFinder.hh"
  6
  7namespace Rivet {
  8
  9
 10  /// ATLAS W pT analysis at 7 TeV
 11  class ATLAS_2011_I925932 : public Analysis {
 12  public:
 13
 14    /// Constructor
 15    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2011_I925932);
 16
 17
 18    /// @name Analysis methods
 19    /// @{
 20
 21    void init() {
 22      // Set up projections
 23      Cut cuts = Cuts::abseta < 2.4 && Cuts::pT > 20*GeV;
 24      declare("MET", MissingMomentum());
 25      LeptonFinder ef_dressed(0.2, cuts && Cuts::abspid == PID::ELECTRON);
 26      declare(ef_dressed, "Elecs_dressed");
 27      LeptonFinder ef_bare(0.0, cuts && Cuts::abspid == PID::ELECTRON);
 28      declare(ef_bare, "Elecs_bare");
 29      LeptonFinder mf_dressed(0.2, cuts && Cuts::abspid == PID::MUON);
 30      declare(mf_dressed, "Muons_dressed");
 31      LeptonFinder mf_bare(0.0, cuts && Cuts::abspid == PID::MUON);
 32      declare(mf_bare, "Muons_bare");
 33
 34      // Book histograms
 35      book(_hist_wpt_dressed_el, 1, 1, 2);
 36      book(_hist_wpt_bare_el,    1, 1, 3);
 37      book(_hist_wpt_dressed_mu, 2, 1, 2);
 38      book(_hist_wpt_bare_mu,    2, 1, 3);
 39    }
 40
 41
 42    /// Do the analysis
 43    void analyze(const Event& event) {
 44
 45      // W reco, starting with MET
 46      const P4& pmiss = apply<MissingMom>(event, "MET").missingMom();
 47      if (pmiss.Et() < 25*GeV) vetoEvent;
 48
 49      // Identify the closest-matching l+MET to m == mW
 50      const Particles& esd = apply<LeptonFinder>(event, "Elecs_dressed").particles();
 51      const int iedfound = closestMatchIndex(esd, pmiss, Kin::mass, 80.4*GeV, 0*GeV, 1000*GeV);
 52      const Particles& esb = apply<LeptonFinder>(event, "Elecs_bare").particles();
 53      const int iebfound = closestMatchIndex(esb, pmiss, Kin::mass, 80.4*GeV, 0*GeV, 1000*GeV);
 54      const Particles& musd = apply<LeptonFinder>(event, "Muons_dressed").particles();
 55      const int imdfound = closestMatchIndex(musd, pmiss, Kin::mass, 80.4*GeV, 0*GeV, 1000*GeV);
 56      const Particles& musb = apply<LeptonFinder>(event, "Muons_bare").particles();
 57      const int imbfound = closestMatchIndex(musb, pmiss, Kin::mass, 80.4*GeV, 0*GeV, 1000*GeV);
 58
 59      MSG_DEBUG("Found " << int(iedfound >= 0) + int(imdfound >= 0) << " dressed W -> e/mu nu");
 60      MSG_DEBUG("Found " << int(iebfound >= 0) + int(imbfound >= 0) << " bare W -> e/mu nu");
 61      if (iedfound < 0 && iebfound < 0 && imdfound < 0 && imbfound < 0) {
 62        MSG_DEBUG("No W bosons found");
 63        vetoEvent;
 64      }
 65
 66      // "Dressed" electron
 67      if (iedfound >= 0) {
 68        const Particle& l = esd[iedfound];
 69        if (mT(pmiss, l) > 40*GeV) _hist_wpt_dressed_el->fill((l.mom() + pmiss).pT()/GeV);
 70      }
 71
 72      // "Bare" electron
 73      if (iebfound >= 0) {
 74        const Particle& l = esb[iebfound];
 75        if (mT(pmiss, l) > 40*GeV) _hist_wpt_bare_el->fill((l.mom() + pmiss).pT()/GeV);
 76      }
 77
 78      // "Dressed" muon
 79      if (imdfound >= 0) {
 80        const Particle& l = musd[imdfound];
 81        if (mT(pmiss, l) > 40*GeV) _hist_wpt_dressed_mu->fill((l.mom() + pmiss).pT()/GeV);
 82      }
 83
 84      // "Bare" muon
 85      if (imbfound >= 0) {
 86        const Particle& l = musb[imbfound];
 87        if (mT(pmiss, l) > 40*GeV) _hist_wpt_bare_mu->fill((l.mom() + pmiss).pT()/GeV);
 88      }
 89
 90    }
 91
 92
 93    // Normalize histos
 94    void finalize() {
 95      normalize(_hist_wpt_dressed_el);
 96      normalize(_hist_wpt_bare_el);
 97      normalize(_hist_wpt_dressed_mu);
 98      normalize(_hist_wpt_bare_mu);
 99    }
100
101    /// @}
102
103
104  private:
105
106    Histo1DPtr _hist_wpt_dressed_el, _hist_wpt_bare_el, _hist_wpt_dressed_mu, _hist_wpt_bare_mu;
107
108  };
109
110
111  RIVET_DECLARE_PLUGIN(ATLAS_2011_I925932);
112
113}