Rivet analyses referenceMC_WJETSMonte Carlo validation observables for $W[e \, \nu]$ + jets productionExperiment: () Status: VALIDATED Authors:
Beams: * * Beam energies: ANY Run details:
Monte Carlo validation observables for $W[\ell \, \nu]$ + jets production 60 GeV $<m_W<$ 100 GeV cut applied. Source code: MC_WJETS.cc 1// -*- C++ -*-
2#include "Rivet/Analyses/MC_JETS_BASE.hh"
3#include "Rivet/Projections/PromptFinalState.hh"
4#include "Rivet/Projections/VetoedFinalState.hh"
5#include "Rivet/Projections/LeptonFinder.hh"
6#include "Rivet/Projections/MissingMomentum.hh"
7#include "Rivet/Projections/FastJets.hh"
8
9namespace Rivet {
10
11
12 /// @brief MC validation analysis for W + jets events
13 class MC_WJETS : public MC_JETS_BASE {
14 public:
15
16 /// Default constructor
17 MC_WJETS(string name="MC_WJETS")
18 : MC_JETS_BASE(name, 4, "Jets")
19 { }
20
21
22 /// @name Analysis methods
23 /// @{
24
25 /// Book histograms
26 void init() {
27
28 // Use analysis options
29 _dR = (getOption("SCHEME") == "BARE") ? 0.0 : 0.2;
30 _lepton = (getOption("LMODE") == "MU") ? PID::MUON : PID::ELECTRON;
31 const double ETACUT = getOption<double>("ABSETALMAX", 3.5);
32 const double PTCUT = getOption<double>("PTLMIN", 25.);
33 const Cut cut = Cuts::abseta < ETACUT && Cuts::pT > PTCUT*GeV;
34
35 // Define projections
36 declare("MET", MissingMomentum());
37 LeptonFinder lf(_dR, cut && Cuts::abspid == _lepton);
38 declare(lf, "Leptons");
39
40 // Set pT cut, jet alg, and clustering radius from input options
41 _jetptcut = getOption<double>("PTJMIN", 20.0) * GeV;
42 const double R = getOption<double>("R", 0.4);
43 JetAlg clusterAlgo;
44 const string algoopt = getOption("ALGO", "ANTIKT");
45 if ( algoopt == "KT" ) {
46 clusterAlgo = JetAlg::KT;
47 } else if ( algoopt == "CA" ) {
48 clusterAlgo = JetAlg::CA;
49 } else if ( algoopt == "ANTIKT" ) {
50 clusterAlgo = JetAlg::ANTIKT;
51 } else {
52 MSG_WARNING("Unknown jet clustering algorithm option " + algoopt + ". Defaulting to anti-kT");
53 clusterAlgo = JetAlg::ANTIKT;
54 }
55
56 // Find jets with clustering radius from input option
57 VetoedFinalState jetinput;
58 jetinput.addVetoOnThisFinalState(lf);
59 FastJets fj(jetinput, clusterAlgo, R);
60 declare(fj, "Jets");
61
62 book(_h_W_jet1_deta, "W_jet1_deta", 50, -5.0, 5.0);
63 book(_h_W_jet1_dR, "W_jet1_dR", 25, 0.5, 7.0);
64
65 MC_JETS_BASE::init();
66 }
67
68
69 /// Do the analysis
70 void analyze(const Event& event) {
71
72 // MET cut
73 const P4& pmiss = apply<MissingMom>(event, "MET").missingMom();
74 if (pmiss.pT() < 25*GeV) vetoEvent;
75
76 // Identify the closest-matching l+MET to m == mW
77 const Particles& ls = apply<LeptonFinder>(event, "Leptons").particles();
78 const int ifound = closestMatchIndex(ls, pmiss, Kin::mass, 80.4*GeV, 60*GeV, 100*GeV);
79 if (ifound < 0) vetoEvent;
80 const Particle& l = ls[ifound];
81 FourMomentum wmom = l.mom() + pmiss;
82
83 const Jets& jets = apply<FastJets>(event, "Jets").jetsByPt(Cuts::pT > _jetptcut);
84 if (jets.size() > 0) {
85 _h_W_jet1_deta->fill(wmom.eta()-jets[0].eta());
86 _h_W_jet1_dR->fill(deltaR(wmom, jets[0].momentum()));
87 }
88
89 MC_JETS_BASE::analyze(event);
90 }
91
92
93 /// Finalize
94 void finalize() {
95 scale(_h_W_jet1_deta, crossSection()/picobarn/sumOfWeights());
96 scale(_h_W_jet1_dR, crossSection()/picobarn/sumOfWeights());
97 MC_JETS_BASE::finalize();
98 }
99
100 /// @}
101
102
103 protected:
104
105 /// @name Parameters for specialised e/mu and dressed/bare subclassing
106 /// @{
107 double _dR;
108 PdgId _lepton;
109 /// @}
110
111
112 private:
113
114 /// @name Histograms
115 /// @{
116 Histo1DPtr _h_W_jet1_deta;
117 Histo1DPtr _h_W_jet1_dR;
118 /// @}
119
120 };
121
122
123 RIVET_DECLARE_PLUGIN(MC_WJETS);
124
125}
|