1// -*- C++ -*-
 2#include "Rivet/Analyses/MC_KTSPLITTINGS_BASE.hh"
 3#include "Rivet/Projections/FastJets.hh"
 4#include "Rivet/Projections/VetoedFinalState.hh"
 5#include "Rivet/Projections/LeptonFinder.hh"
 6#include "Rivet/Projections/MissingMomentum.hh"
 7
 8namespace Rivet {
 9
10
11  /// @brief MC validation analysis for e+mu W^+W^- + jets events
12  class MC_WWKTSPLITTINGS : public MC_KTSPLITTINGS_BASE {
13  public:
14
15    /// Default constructor
16    MC_WWKTSPLITTINGS()
17      : MC_KTSPLITTINGS_BASE("MC_WWKTSPLITTINGS", 4, "Jets")
18    {    }
19
20
21    /// @name Analysis methods
22    /// @{
23
24    /// Book histograms
25    void init() {
26      declare("MET", MissingMomentum());
27
28      // Find electrons with cuts from input options
29      const double ETAECUT = getOption<double>("ABSETAEMAX", 3.5);
30      const double PTECUT = getOption<double>("PTEMIN", 25.);
31      const Cut cut_e = Cuts::abseta < ETAECUT && Cuts::pT > PTECUT*GeV;
32      LeptonFinder ef(0.1, cut_e && Cuts::abspid == PID::ELECTRON);
33      declare(ef, "Elecs");
34
35      // Find muons with cuts from input options
36      const double ETAMUCUT = getOption<double>("ABSETAMUMAX", 3.5);
37      const double PTMUCUT = getOption<double>("PTMUMIN", 25.);
38      const Cut cut_m = Cuts::abseta < ETAMUCUT && Cuts::pT > PTMUCUT*GeV;
39      LeptonFinder mf(0.2, cut_m && Cuts::abspid == PID::MUON);
40      declare(mf, "Muons");
41
42      // Find jets with clustering radius from input option
43      VetoedFinalState jetinput;
44      jetinput
45	.addVetoOnThisFinalState(ef)
46	.addVetoOnThisFinalState(mf);
47      const double R = getOption<double>("R", 0.6);
48      FastJets fj(jetinput, JetAlg::KT, R);
49      declare(fj, "Jets");
50
51      MC_KTSPLITTINGS_BASE::init();
52    }
53
54
55
56    /// Do the analysis
57    void analyze(const Event& event) {
58
59      // MET cut
60      const P4& pmiss = apply<MissingMom>(event, "MET").missingMom();
61      if (pmiss.pT() < 25*GeV) vetoEvent;
62
63      // Identify the closest-matching l+MET to m == mW
64      /// @note Dubious strategy, given there are two neutrinos...
65      const Particles& es = apply<LeptonFinder>(event, "Elecs").particles();
66      const int iefound = closestMatchIndex(es, pmiss, Kin::mass, 80.4*GeV, 60*GeV, 100*GeV);
67      const Particles& mus = apply<LeptonFinder>(event, "Muons").particles();
68      const int imfound = closestMatchIndex(mus, pmiss, Kin::mass, 80.4*GeV, 60*GeV, 100*GeV);
69
70      // Require two valid W candidates
71      if (iefound < 0 || imfound < 0) vetoEvent;
72
73      MC_KTSPLITTINGS_BASE::analyze(event);
74    }
75
76
77    /// Finalize
78    void finalize() {
79      MC_KTSPLITTINGS_BASE::finalize();
80    }
81
82    /// @}
83
84  };
85
86
87  RIVET_DECLARE_PLUGIN(MC_WWKTSPLITTINGS);
88
89}