Rivet analyses referenceMC_ZVBFMonte Carlo validation observables for VBF $Z[\ell^+ \, \ell^-]$ + 2 jet productionExperiment: () Status: VALIDATED Authors:
Beams: * * Beam energies: ANY Run details:
Available observables are the pT of jets 1-4, jet multiplicity, $\Delta\eta(Z, \text{jet1})$, $\Delta R(\text{jet2}, \text{jet3})$, $m_{jj}$, $H_\text{T}$ and third-jet centrality. Source code: MC_ZVBF.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/ZFinder.hh"
4#include "Rivet/Projections/FastJets.hh"
5
6namespace Rivet {
7
8
9 /// @brief MC validation analysis for Zjj events
10 class MC_ZVBF : public Analysis {
11 public:
12
13 /// Default constructor
14 RIVET_DEFAULT_ANALYSIS_CTOR(MC_ZVBF);
15
16 /// @name Analysis methods
17 //@{
18
19 /// Initialize
20 void init() {
21 _dR=0.1;
22 if (getOption("SCHEME") == "BARE") _dR = 0.0;
23 _lepton=PID::ELECTRON;
24 if (getOption("LMODE") == "MU") _lepton = PID::MUON;
25
26 FinalState fs;
27 Cut cut = Cuts::abseta < 3.5 && Cuts::pT > 25*GeV;
28 ZFinder zfinder(fs, cut, _lepton, 65*GeV, 115*GeV, _dR, ZFinder::ClusterPhotons::NODECAY, ZFinder::AddPhotons::YES);
29 declare(zfinder, "ZFinder");
30 FastJets jetpro(zfinder.remainingFinalState(), FastJets::ANTIKT, 0.4);
31 declare(jetpro, "Jets");
32
33 const double sqrts = sqrtS() ? sqrtS() : 14*TeV;
34 book(_h["gap_inc"], "N_gapjets_inclusive", 8, -0.5, 7.5);
35 book(_h["gap_exc"], "N_gapjets_exclusive", 8, -0.5, 7.5);
36 book(_h["Z_jet1_deta"], "Z_jet1_deta", 50, -5, 5);
37 book(_h["Z_jet1_dR"], "Z_jet1_dR", 25, 0.5, 7.0);
38 book(_h["HT"], "jets_HT", logspace(40, 50, sqrts/GeV/2.0));
39 book(_h["mjj"], "m_jj", 40, 200.0, sqrts/GeV/2.0);
40 book(_h["jve_mjj"], "_jve_mjj", 40, 200.0, sqrts/GeV/2.0);
41 book(_h["pTV"] ,"Z_pT", logspace(100, 1.0, 0.5*sqrts/GeV));
42 book(_h["dphi"], "dphi_jj", 20., -1., 1.);
43 book(_h["drap"], "drap_jj", 20., -10., 10.);
44 book(_h["3JC"], "jet_3_centrality", 25., -2.5, 2.5);
45
46 book(_s["jve_mjj"], "jet_veto_efficiency_mjj");
47
48 for (size_t i = 0; i < 4; ++i) {
49 const string pTname = "jet_pT_" + to_str(i+1);
50 const double pTmax = 1.0/(double(i)+2.0) * sqrts/GeV/2.0;
51 const int nbins_pT = 100/(i+1);
52 if (pTmax > 10) { // Protection aginst logspace exception, needed for LEP
53 book(_h[pTname], pTname, logspace(nbins_pT, 10.0, pTmax));
54 }
55 const string etaname = "jet_eta_" + to_str(i+1);
56 book(_h[etaname], etaname, (i > 1 ? 25 : 50), -5.0, 5.0);
57 const string rapname = "jet_y_" + to_str(i+1);
58 book(_h[rapname], rapname, (i > 1 ? 25 : 50), -5.0, 5.0);
59 const string phiname = "jet_phi_" + to_str(i+1);
60 book(_h[phiname], phiname, (i > 1 ? 25 : 50), -1.0, 1.0);
61 }
62 }
63
64
65 /// Do the analysis
66 void analyze(const Event & e) {
67 MSG_TRACE("MC_ZVBF: running ZFinder");
68 const ZFinder& zfinder = apply<ZFinder>(e, "ZFinder");
69 if (zfinder.bosons().size() != 1) vetoEvent;
70 const FourMomentum& zmom = zfinder.bosons()[0].momentum();
71 MSG_TRACE("MC_ZVBF: have exactly one Z boson candidate");
72
73 const Jets& jets = apply<FastJets>(e, "Jets").jetsByPt(Cuts::absrap < 5 && Cuts::pT > 30*GeV);
74 if (jets.size() < 2) {
75 MSG_TRACE("MC_ZVBF: does not have at least two valid jets");
76 vetoEvent;
77 }
78
79 Jet tag1 = jets.at(0);
80 Jet tag2 = jets.at(1);
81 const double mjj = (tag1.mom() + tag2.mom()).mass()/GeV;
82 if (mjj < 200.) {
83 MSG_TRACE("MC_ZVBF: should have at least 200 GeV in Mjj");
84 vetoEvent;
85 }
86
87 // jet kinematics
88 for (size_t i = 0; i < min(4u, jets.size()); ++i) {
89 const string pTname = "jet_pT_" + to_str(i+1);
90 const string etaname = "jet_eta_" + to_str(i+1);
91 const string rapname = "jet_y_" + to_str(i+1);
92 const string phiname = "jet_phi_" + to_str(i+1);
93 _h[pTname]->fill(jets[i].pT()/GeV);
94 _h[etaname]->fill(jets[i].eta());
95 _h[rapname]->fill(jets[i].rap());
96 _h[phiname]->fill(mapAngleMPiToPi(jets[i].phi())/M_PI);
97 }
98
99 size_t n_gap = 0;
100 // start loop at the 3rd hardest pT jet
101 for (size_t i = 2; i < jets.size(); ++i) {
102 const Jet j = jets.at(i);
103 if (isBetween(j, tag1, tag2)) ++n_gap;
104 }
105
106 // gap-jet multiplicities
107 _h["gap_exc"]->fill(n_gap);
108 for (size_t i = 0; i <= 7; ++i) {
109 if (n_gap >= i) {
110 _h["gap_inc"]->fill(i);
111 }
112 }
113
114 _h["jve_mjj"]->fill(mjj);
115 if (n_gap) {
116 // third-jet centrality
117 const double rap1 = jets[0].rap();
118 const double rap2 = jets[1].rap();
119 const double rap3 = jets[2].rap();
120 const double JC = (rap3 - 0.5*(rap1 + rap2))/(rap1 - rap2);
121 _h["3JC"]->fill(JC);
122 }
123 else {
124 MSG_TRACE("MC_ZVBF: should satisfy a CJV");
125 const double HT = sum(jets, Kin::pT, 0.0)/GeV;
126 _h["HT"]->fill(HT);
127 _h["mjj"]->fill(mjj);
128 _h["pTV"]->fill(zmom.pT()/GeV);
129 _h["dphi"]->fill(signedDeltaPhi(tag1, tag2));
130 _h["drap"]->fill(tag1.rap() - tag2.rap());
131
132 // Z tagging jet correlations
133 _h["Z_jet1_deta"]->fill(zmom.eta()-jets[0].eta());
134 _h["Z_jet1_dR"]->fill(deltaR(zmom, jets[0].momentum()));
135 }
136 }
137
138
139 /// Finalize
140 void finalize() {
141 scale(_h, crossSection()/femtobarn/sumOfWeights());
142 efficiency(_h["mjj"], _h["jve_mjj"], _s["jve_mjj"]);
143 }
144
145 //@}
146
147 // check if jet is between tagging jets
148 bool isBetween(const Jet &probe, const Jet &boundary1, const Jet &boundary2) {
149 double y_p = probe.rapidity();
150 double y_b1 = boundary1.rapidity();
151 double y_b2 = boundary2.rapidity();
152
153 double y_min = std::min(y_b1, y_b2);
154 double y_max = std::max(y_b1, y_b2);
155
156 return (y_p > y_min && y_p < y_max);
157 }
158
159 double signedDeltaPhi(Jet &j1, Jet &j2) {
160 double dphijj = 0.;
161 if (j1.rap() > j2.rap()) dphijj = j1.phi() - j2.phi();
162 else dphijj = j2.phi() - j1.phi();
163 return mapAngleMPiToPi(dphijj)/M_PI;
164 }
165
166 private:
167
168 /// @name Parameters for specialised e/mu and dressed/bare subclassing
169 //@{
170 double _dR;
171 PdgId _lepton;
172 //@}
173
174 /// @name Histograms
175 //@{
176 map<string,Histo1DPtr> _h;
177 map<string,Scatter2DPtr> _s;
178 //@}
179
180 };
181
182 RIVET_DECLARE_PLUGIN(MC_ZVBF);
183}
|