Rivet analyses referenceATLAS_2014_I1282447W + charm production at 7 TeVExperiment: ATLAS (LHC) Inspire ID: 1282447 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
This routine implements the measurement of $W$ boson production in association with a single charm quark with $4.6 \text{fb}^{-1}$ of $pp$ collision data at $\sqrt{s} = 7$ TeV collected with the ATLAS detector at the Large Hadron Collider. Results are quoted for the $W$ boson production decaying into either an electron or muon and the charm quark being identified as any charmed hadron above 5 GeV inside $\Delta R = 0.4$ of a jet with more than 25 GeV. Alternatively the presence of the charm quark is indicated by the presence of a $D$ or a $D^*$ meson with a $p_\text{T}$ above 8 GeV. The cross sections are quoted for the number of opposite sign minus same sign events, where the signs under consideration are the charge of the $W$ boson and the charmed hadrons tagging the event. Given are the integrated and differential cross sections as a function of the pseudorapidity of the lepton from the $W$-boson decay are measured. Additionally, the cross-section ratios $\sigma(W^+ + bar{c}) / \sigma(W^- + c)$ as well as the $p_\text{T}$ dependent cross sections of the $D$/$D^*$ mesons normalized to the $W$ inclusive cross sections are published. The measured data is unfolded to the Born level. One should therefore take care to run on samples without QED radiation off of the electrons. IMPORTANT NOTICE --- For the MC predictions in the paper, the branching fractions to $D$ and $D^*$ mesons have been adapted to be 0.2256 ($D$) and 0.2287 ($D^*$) (LEP/HERA combination). Some suggestion on how to do post-processing -- in case separate samples for $W^+ + c$ and $W^- + c$ and inclusive $W$ are used -- are included in the cc-file. This post-processing is needed to properly display the histograms for the cross section ratio plots $\sigma(W^+ + bar{c}) / \sigma(W^- + c)$ as well as for the cross section ratios of $W + D^{(*)}$ production over inclusive $W$ production ($\sigma(W^{+/-} D^{(*)}) / \sigma(W^{+/-})$) as a function of the $D^{(*)}$ meson transverse momentum. Source code: ATLAS_2014_I1282447.cc 1// -*- C++ -*-
2// ATLAS W+c analysis
3
4//////////////////////////////////////////////////////////////////////////
5/*
6 Description of rivet analysis ATLAS_2014_I1282447 W+c production
7
8 This rivet routine implements the ATLAS W+c analysis.
9 Apart from those histograms, described and published on HEP Data, here
10 are some helper histograms defined, these are:
11
12 d02-x01-y01, d02-x01-y02 and d08-x01-y01 are ratios, the nominator ("_plus")
13 and denominator ("_minus") histograms are also given, so that the ratios can
14 be reconstructed if need be (e.g. when running on separate samples).
15
16 d05 and d06 are ratios over inclusive W production.
17 The routine has to be run on a sample for inclusive W production in order to
18 make sure the denominator ("_winc") is correctly filled.
19
20 The ratios can be constructed using the following sample code:
21 python divideWCharm.py
22
23 import yoda
24 hists_wc = yoda.read("Rivet_Wc.yoda")
25 hists_winc = yoda.read("Rivet_Winc.yoda")
26
27 ## division histograms --> ONLY for different plus minus runs
28 # (merge before using yodamerge Rivet_plus.yoda Rivet_minus.yoda > Rivet_Wc.yoda)
29
30 d02y01_plus = hists_wc["/ATLAS_2014_I1282447/d02-x01-y01_plus"]
31 d02y01_minus = hists_wc["/ATLAS_2014_I1282447/d02-x01-y01_minus"]
32 ratio_d02y01 = d02y01_plus.divide(d02y01_minus)
33 ratio_d02y01.path = "/ATLAS_2014_I1282447/d02-x01-y01"
34
35 d02y02_plus = hists_wc["/ATLAS_2014_I1282447/d02-x01-y02_plus"]
36 d02y02_minus = hists_wc["/ATLAS_2014_I1282447/d02-x01-y02_minus"]
37 ratio_d02y02= d02y02_plus.divide(d02y02_minus)
38 ratio_d02y02.path = "/ATLAS_2014_I1282447/d02-x01-y02"
39
40 d08y01_plus = hists_wc["/ATLAS_2014_I1282447/d08-x01-y01_plus"]
41 d08y01_minus = hists_wc["/ATLAS_2014_I1282447/d08-x01-y01_minus"]
42 ratio_d08y01= d08y01_plus.divide(d08y01_minus)
43 ratio_d08y01.path = "/ATLAS_2014_I1282447/d08-x01-y01"
44
45 # inclusive cross section
46 h_winc = hists_winc["/ATLAS_2014_I1282447/d05-x01-y01"]
47 h_d = hists_wc["/ATLAS_2014_I1282447/d01-x01-y02"]
48 h_dstar= hists_wc["/ATLAS_2014_I1282447/d01-x01-y03"]
49
50 ratio_wd = h_d.divide(h_winc)
51 ratio_wd.path = "/ATLAS_2014_I1282447/d05-x01-y02"
52
53 ratio_wdstar = h_d.divide(h_winc)
54 ratio_wdstar.path = "/ATLAS_2014_I1282447/d05-x01-y03"
55
56 # pT differential
57 h_winc_plus = hists_winc["/ATLAS_2014_I1282447/d06-x01-y01_winc"]
58 h_winc_minus = hists_winc["/ATLAS_2014_I1282447/d06-x01-y02_winc"]
59
60 h_wd_plus = hists_wc["/ATLAS_2014_I1282447/d06-x01-y01_wplus"]
61 h_wd_minus = hists_wc["/ATLAS_2014_I1282447/d06-x01-y02_wminus"]
62 h_wdstar_plus = hists_wc["/ATLAS_2014_I1282447/d06-x01-y03_wplus"]
63 h_wdstar_minus = hists_wc["/ATLAS_2014_I1282447/d06-x01-y04_wminus"]
64
65 ratio_wd_plus = h_wd_plus.divide(h_winc_plus)
66 ratio_wd_plus.path = "/ATLAS_2014_I1282447/d06-x01-y01"
67 ratio_wd_minus = h_wd_plus.divide(h_winc_minus)
68 ratio_wd_minus.path = "/ATLAS_2014_I1282447/d06-x01-y02"
69
70 ratio_wdstar_plus = h_wdstar_plus.divide(h_winc_plus)
71 ratio_wdstar_plus.path = "/ATLAS_2014_I1282447/d06-x01-y03"
72 ratio_wdstar_minus = h_wdstar_plus.divide(h_winc_minus)
73 ratio_wdstar_minus.path = "/ATLAS_2014_I1282447/d06-x01-y04"
74
75 ratio_wd_plus = h_wd_plus.divide(h_winc_plus)
76 ratio_wd_plus.path = "/ATLAS_2014_I1282447/d06-x01-y01"
77 ratio_wd_minus = h_wd_plus.divide(h_winc_minus)
78 ratio_wd_minus.path = "/ATLAS_2014_I1282447/d06-x01-y02"
79
80 h_winc_plus= hists_winc["/ATLAS_2014_I1282447/d06-x01-y01_winc"]
81 h_winc_minus= hists_winc["/ATLAS_2014_I1282447/d06-x01-y02_winc"]
82
83 ## copy other histograms for plotting
84
85 d01x01y01= hists_wc["/ATLAS_2014_I1282447/d01-x01-y01"]
86 d01x01y01.path = "/ATLAS_2014_I1282447/d01-x01-y01"
87
88 d01x01y02= hists_wc["/ATLAS_2014_I1282447/d01-x01-y02"]
89 d01x01y02.path = "/ATLAS_2014_I1282447/d01-x01-y02"
90
91 d01x01y03= hists_wc["/ATLAS_2014_I1282447/d01-x01-y03"]
92 d01x01y03.path = "/ATLAS_2014_I1282447/d01-x01-y03"
93
94 d03x01y01= hists_wc["/ATLAS_2014_I1282447/d03-x01-y01"]
95 d03x01y01.path = "/ATLAS_2014_I1282447/d03-x01-y01"
96
97 d03x01y02= hists_wc["/ATLAS_2014_I1282447/d03-x01-y02"]
98 d03x01y02.path = "/ATLAS_2014_I1282447/d03-x01-y02"
99
100 d04x01y01= hists_wc["/ATLAS_2014_I1282447/d04-x01-y01"]
101 d04x01y01.path = "/ATLAS_2014_I1282447/d04-x01-y01"
102
103 d04x01y02= hists_wc["/ATLAS_2014_I1282447/d04-x01-y02"]
104 d04x01y02.path = "/ATLAS_2014_I1282447/d04-x01-y02"
105
106 d04x01y03= hists_wc["/ATLAS_2014_I1282447/d04-x01-y03"]
107 d04x01y03.path = "/ATLAS_2014_I1282447/d04-x01-y03"
108
109 d04x01y04= hists_wc["/ATLAS_2014_I1282447/d04-x01-y04"]
110 d04x01y04.path = "/ATLAS_2014_I1282447/d04-x01-y04"
111
112 d07x01y01= hists_wc["/ATLAS_2014_I1282447/d07-x01-y01"]
113 d07x01y01.path = "/ATLAS_2014_I1282447/d07-x01-y01"
114
115 yoda.write([ratio_d02y01,ratio_d02y02,ratio_d08y01, ratio_wd ,ratio_wdstar,ratio_wd_plus,ratio_wd_minus ,ratio_wdstar_plus,ratio_wdstar_minus,d01x01y01,d01x01y02,d01x01y03,d03x01y01,d03x01y02,d04x01y01,d04x01y02,d04x01y03,d04x01y04,d07x01y01],"validation.yoda")
116
117*/
118//////////////////////////////////////////////////////////////////////////
119
120#include "Rivet/Analysis.hh"
121#include "Rivet/Projections/UnstableParticles.hh"
122#include "Rivet/Projections/WFinder.hh"
123#include "Rivet/Projections/FastJets.hh"
124#include "Rivet/Projections/ChargedFinalState.hh"
125#include "Rivet/Projections/VetoedFinalState.hh"
126
127namespace Rivet {
128
129
130
131 class ATLAS_2014_I1282447 : public Analysis {
132 public:
133
134 /// Constructor
135 ATLAS_2014_I1282447() : Analysis("ATLAS_2014_I1282447")
136 {
137
138 }
139
140
141 /// @name Analysis methods
142 //@{
143
144 /// Book histograms and initialise projections before the run
145 void init() {
146
147 /// @todo Initialise and register projections here
148 UnstableParticles fs;
149
150 Cut cuts = Cuts::etaIn(-2.5, 2.5) & (Cuts::pT > 20*GeV);
151
152 /// should use sample WITHOUT QED radiation off the electron
153 WFinder wfinder_born_el(fs, cuts, PID::ELECTRON, 25*GeV, 8000*GeV, 15*GeV, 0.1, WFinder::ChargedLeptons::PROMPT, WFinder::ClusterPhotons::ALL, WFinder::AddPhotons::YES);
154 declare(wfinder_born_el, "WFinder_born_el");
155
156 WFinder wfinder_born_mu(fs, cuts, PID::MUON , 25*GeV, 8000*GeV, 15*GeV, 0.1, WFinder::ChargedLeptons::PROMPT, WFinder::ClusterPhotons::ALL, WFinder::AddPhotons::YES);
157 declare(wfinder_born_mu, "WFinder_born_mu");
158
159 // all hadrons that could be coming from a charm decay --
160 // -- for safety, use region -3.5 - 3.5
161 declare(UnstableParticles(Cuts::abseta <3.5), "hadrons");
162
163 // Input for the jets: no neutrinos, no muons, and no electron which passed the electron cuts
164 // also: NO electron, muon or tau (needed due to ATLAS jet truth reconstruction feature)
165 VetoedFinalState veto;
166
167 veto.addVetoOnThisFinalState(wfinder_born_el);
168 veto.addVetoOnThisFinalState(wfinder_born_mu);
169 veto.addVetoPairId(PID::ELECTRON);
170 veto.addVetoPairId(PID::MUON);
171 veto.addVetoPairId(PID::TAU);
172
173 FastJets jets(veto, FastJets::ANTIKT, 0.4);
174 declare(jets, "jets");
175
176 // Book histograms
177
178 // charge separated integrated cross sections
179 book(_hist_wcjet_charge ,"d01-x01-y01");
180 book(_hist_wd_charge ,"d01-x01-y02");
181 book(_hist_wdstar_charge ,"d01-x01-y03");
182
183 // charge integrated total cross sections
184 book(_hist_wcjet_ratio,"d02-x01-y01");
185 book(_hist_wd_ratio ,"d02-x01-y02");
186
187 book(_hist_wcjet_minus ,"d02-x01-y01_minus");
188 book(_hist_wd_minus ,"d02-x01-y02_minus");
189
190 book(_hist_wcjet_plus ,"d02-x01-y01_plus");
191 book(_hist_wd_plus ,"d02-x01-y02_plus");
192
193 // eta distributions
194 book(_hist_wplus_wcjet_eta_lep ,"d03-x01-y01");
195 book(_hist_wminus_wcjet_eta_lep ,"d03-x01-y02");
196
197 book(_hist_wplus_wdminus_eta_lep ,"d04-x01-y01");
198 book(_hist_wminus_wdplus_eta_lep ,"d04-x01-y02");
199 book(_hist_wplus_wdstar_eta_lep ,"d04-x01-y03");
200 book(_hist_wminus_wdstar_eta_lep ,"d04-x01-y04");
201
202 // ratio of cross section (WD over W inclusive) // postprocess!
203 book(_hist_w_inc ,"d05-x01-y01");
204 book(_hist_wd_winc_ratio ,"d05-x01-y02");
205 book(_hist_wdstar_winc_ratio,"d05-x01-y03");
206
207 // ratio of cross section (WD over W inclusive -- function of pT of D meson)
208 book(_hist_wplusd_wplusinc_pt_ratio ,"d06-x01-y01");
209 book(_hist_wminusd_wminusinc_pt_ratio ,"d06-x01-y02");
210 book(_hist_wplusdstar_wplusinc_pt_ratio ,"d06-x01-y03");
211 book(_hist_wminusdstar_wminusinc_pt_ratio,"d06-x01-y04");
212
213 // could use for postprocessing!
214 book(_hist_wplusd_wplusinc_pt ,"d06-x01-y01_wplus");
215 book(_hist_wminusd_wminusinc_pt ,"d06-x01-y02_wminus");
216 book(_hist_wplusdstar_wplusinc_pt ,"d06-x01-y03_wplus");
217 book(_hist_wminusdstar_wminusinc_pt ,"d06-x01-y04_wminus");
218
219 book(_hist_wplus_winc ,"d06-x01-y01_winc");
220 book(_hist_wminus_winc ,"d06-x01-y02_winc");
221
222 // jet multiplicity of charge integrated W+cjet cross section (+0 or +1 jet in addition to the charm jet)
223 book(_hist_wcjet_jets ,"d07-x01-y01");
224
225 // jet multiplicity of W+cjet cross section ratio (+0 or +1 jet in addition to the charm jet)
226 book(_hist_wcjet_jets_ratio ,"d08-x01-y01");
227 book(_hist_wcjet_jets_plus ,"d08-x01-y01_plus");
228 book(_hist_wcjet_jets_minus ,"d08-x01-y01_minus");
229
230 }
231
232
233 /// Perform the per-event analysis
234 void analyze(const Event& event) {
235
236 double charge_weight = 0; // account for OS/SS events
237
238 int lepton_charge = 0;
239 double lepton_eta = 0.;
240
241 /// Find leptons
242 const WFinder& wfinder_born_el = apply<WFinder>(event, "WFinder_born_el");
243 const WFinder& wfinder_born_mu = apply<WFinder>(event, "WFinder_born_mu");
244
245 if (wfinder_born_el.empty() && wfinder_born_mu.empty()) {
246 MSG_DEBUG("No W bosons found");
247 vetoEvent;
248 }
249
250 bool keepevent = false;
251
252 //check electrons
253 if (!wfinder_born_el.empty()) {
254 const FourMomentum nu = wfinder_born_el.constituentNeutrinos()[0];
255 if (wfinder_born_el.mT() > 40*GeV && nu.pT() > 25*GeV) {
256 keepevent = true;
257 lepton_charge = wfinder_born_el.constituentLeptons()[0].charge();
258 lepton_eta = fabs(wfinder_born_el.constituentLeptons()[0].pseudorapidity());
259 }
260 }
261
262 //check muons
263 if (!wfinder_born_mu.empty()) {
264 const FourMomentum nu = wfinder_born_mu.constituentNeutrinos()[0];
265 if (wfinder_born_mu.mT() > 40*GeV && nu.pT() > 25*GeV) {
266 keepevent = true;
267 lepton_charge = wfinder_born_mu.constituentLeptons()[0].charge();
268 lepton_eta = fabs(wfinder_born_mu.constituentLeptons()[0].pseudorapidity());
269 }
270 }
271
272 if (!keepevent) {
273 MSG_DEBUG("Event does not pass mT and MET cuts");
274 vetoEvent;
275 }
276
277 if (lepton_charge > 0) {
278 _hist_wplus_winc->fill(10.);
279 _hist_wplus_winc->fill(16.);
280 _hist_wplus_winc->fill(30.);
281 _hist_wplus_winc->fill(60.);
282 _hist_w_inc->fill(+1);
283 }
284 else if (lepton_charge < 0) {
285 _hist_wminus_winc->fill(10.);
286 _hist_wminus_winc->fill(16.);
287 _hist_wminus_winc->fill(30.);
288 _hist_wminus_winc->fill(60.);
289 _hist_w_inc->fill(-1);
290 }
291
292 // Find hadrons in the event
293 const UnstableParticles& fs = apply<UnstableParticles>(event, "hadrons");
294
295 /// FIND Different channels
296 // 1: wcjet
297 // get jets
298 const Jets& jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT>25.0*GeV && Cuts::abseta<2.5);
299 // loop over jets to select jets used to match to charm
300 Jets js;
301 int matched_charmHadron = 0;
302 double charm_charge = 0.;
303 int njets = 0;
304 int nj = 0;
305 bool mat_jet = false;
306
307 double ptcharm = 0;
308 if (matched_charmHadron > -1) {
309 for (const Jet& j : jets) {
310 mat_jet = false;
311
312 njets += 1;
313 for (const Particle& p : fs.particles()) {
314 /// @todo Avoid touching HepMC!
315 ConstGenParticlePtr part = p.genParticle();
316 if (p.hasCharm()) {
317 //if(isFromBDecay(p)) continue;
318 if (p.fromBottom()) continue;
319 if (p.pT() < 5*GeV ) continue;
320 if (hasCharmedChildren(part)) continue;
321 if (deltaR(p, j) < 0.3) {
322 mat_jet = true;
323 if (p.pT() > ptcharm) {
324 charm_charge = part->pdg_id();
325 ptcharm = p.pT();
326 }
327 }
328 }
329 }
330 if (mat_jet) nj++;
331 }
332
333 if (charm_charge * lepton_charge > 0) charge_weight = -1;
334 else charge_weight = +1;
335
336 if (nj == 1) {
337 if (lepton_charge > 0) {
338 _hist_wcjet_charge ->fill( 1, charge_weight);
339 _hist_wcjet_plus ->fill( 0, charge_weight);
340 _hist_wplus_wcjet_eta_lep ->fill(lepton_eta, charge_weight);
341 _hist_wcjet_jets_plus ->fill(njets-1 , charge_weight);
342 }
343 else if (lepton_charge < 0) {
344 _hist_wcjet_charge ->fill( -1, charge_weight);
345 _hist_wcjet_minus ->fill( 0, charge_weight);
346 _hist_wminus_wcjet_eta_lep->fill(lepton_eta, charge_weight);
347 _hist_wcjet_jets_minus ->fill(njets-1 , charge_weight);
348 }
349
350 _hist_wcjet_jets->fill(njets-1, charge_weight);
351 }
352 }
353
354 // // 1/2: w+d(*) meson
355
356 for (const Particle& p : fs.particles()) {
357
358 /// @todo Avoid touching HepMC!
359 ConstGenParticlePtr part = p.genParticle();
360 if (p.pT() < 8*GeV) continue;
361 if (fabs(p.eta()) > 2.2) continue;
362
363 // W+D
364 if (abs(part->pdg_id()) == 411) {
365 if (lepton_charge * part->pdg_id() > 0) charge_weight = -1;
366 else charge_weight = +1;
367
368 // fill histos
369 if (lepton_charge > 0) {
370 _hist_wd_charge ->fill( 1, charge_weight);
371 _hist_wd_plus ->fill( 0, charge_weight);
372 _hist_wplus_wdminus_eta_lep->fill(lepton_eta, charge_weight);
373 _hist_wplusd_wplusinc_pt ->fill( p.pT(), charge_weight);
374 }
375 else if (lepton_charge < 0) {
376 _hist_wd_charge ->fill( -1, charge_weight);
377 _hist_wd_minus ->fill( 0, charge_weight);
378 _hist_wminus_wdplus_eta_lep->fill(lepton_eta, charge_weight);
379 _hist_wminusd_wminusinc_pt ->fill(p.pT() , charge_weight);
380 }
381 }
382
383 // W+Dstar
384 if ( abs(part->pdg_id()) == 413 ) {
385 if (lepton_charge*part->pdg_id() > 0) charge_weight = -1;
386 else charge_weight = +1;
387
388 if (lepton_charge > 0) {
389 _hist_wdstar_charge->fill(+1, charge_weight);
390 _hist_wd_plus->fill( 0, charge_weight);
391 _hist_wplus_wdstar_eta_lep->fill( lepton_eta, charge_weight);
392 _hist_wplusdstar_wplusinc_pt->fill( p.pT(), charge_weight);
393 }
394 else if (lepton_charge < 0) {
395 _hist_wdstar_charge->fill(-1, charge_weight);
396 _hist_wd_minus->fill(0, charge_weight);
397 _hist_wminus_wdstar_eta_lep->fill(lepton_eta, charge_weight);
398 _hist_wminusdstar_wminusinc_pt->fill(p.pT(), charge_weight);
399 }
400 }
401 }
402
403 }
404
405
406 /// Normalise histograms etc., after the run
407 void finalize() {
408
409 const double sf = crossSection() / sumOfWeights();
410
411 // norm to cross section
412 // d01
413 scale(_hist_wcjet_charge, sf);
414 scale(_hist_wd_charge, sf);
415 scale(_hist_wdstar_charge, sf);
416
417 //d02
418 scale(_hist_wcjet_plus, sf);
419 scale(_hist_wcjet_minus, sf);
420 scale(_hist_wd_plus, sf);
421 scale(_hist_wd_minus, sf);
422
423 divide(_hist_wcjet_plus, _hist_wcjet_minus, _hist_wcjet_ratio);
424 divide(_hist_wd_plus, _hist_wd_minus, _hist_wd_ratio );
425
426 //d03
427 scale(_hist_wplus_wcjet_eta_lep, sf);
428 scale(_hist_wminus_wcjet_eta_lep, sf);
429
430 //d04
431 scale(_hist_wplus_wdminus_eta_lep, crossSection()/sumOfWeights());
432 scale(_hist_wminus_wdplus_eta_lep, crossSection()/sumOfWeights());
433
434 scale(_hist_wplus_wdstar_eta_lep , crossSection()/sumOfWeights());
435 scale(_hist_wminus_wdstar_eta_lep, crossSection()/sumOfWeights());
436
437 //d05
438 scale(_hist_w_inc, 0.01 * sf); // in percent --> /100
439 divide(_hist_wd_charge, _hist_w_inc, _hist_wd_winc_ratio );
440 divide(_hist_wdstar_charge, _hist_w_inc, _hist_wdstar_winc_ratio);
441
442 //d06, in percentage!
443 scale(_hist_wplusd_wplusinc_pt, sf);
444 scale(_hist_wminusd_wminusinc_pt, sf);
445 scale(_hist_wplusdstar_wplusinc_pt, sf);
446 scale(_hist_wminusdstar_wminusinc_pt, sf);
447
448 scale(_hist_wplus_winc, 0.01 * sf); // in percent --> /100
449 scale(_hist_wminus_winc, 0.01 * sf); // in percent --> /100
450
451 divide(_hist_wplusd_wplusinc_pt, _hist_wplus_winc , _hist_wplusd_wplusinc_pt_ratio );
452 divide(_hist_wminusd_wminusinc_pt, _hist_wminus_winc, _hist_wminusd_wminusinc_pt_ratio );
453 divide(_hist_wplusdstar_wplusinc_pt, _hist_wplus_winc , _hist_wplusdstar_wplusinc_pt_ratio );
454 divide(_hist_wminusdstar_wminusinc_pt, _hist_wminus_winc, _hist_wminusdstar_wminusinc_pt_ratio);
455
456
457 //d07
458 scale(_hist_wcjet_jets, sf);
459
460 //d08
461 scale(_hist_wcjet_jets_minus, sf);
462 scale(_hist_wcjet_jets_plus, sf);
463 divide(_hist_wcjet_jets_plus, _hist_wcjet_jets_minus , _hist_wcjet_jets_ratio);
464 }
465
466 //@}
467
468
469 private:
470
471 // Data members like post-cuts event weight counters go here
472
473 // Check whether particle comes from b-decay
474 bool isFromBDecay(const Particle& p) {
475
476 /// @todo I think we can just replicated the original behaviour with this call
477 /// Note slight difference to Rivet's native Particle::fromBottom method!
478 return p.hasAncestorWith([](const Particle &p)->bool{return p.hasBottom();});
479 /*
480 bool isfromB = false;
481
482 if (p.genParticle() == nullptr) return false;
483
484 ConstGenParticlePtr part = p.genParticle();
485 ConstGenVertexPtr ivtx = part->production_vertex();
486 while (ivtx) {
487 if (ivtx->particles_in().size() < 1) {
488 isfromB = false;
489 break;
490 }
491 const HepMC::GenVertex::particles_in_const_iterator iPart_invtx = ivtx->particles_in_const_begin();
492 part = (*iPart_invtx);
493 if (!part) {
494 isfromB = false;
495 break;
496 }
497 isfromB = PID::hasBottom(part->pdg_id());
498
499 if (isfromB == true) break;
500 ivtx = part->production_vertex();
501 if ( part->pdg_id() == 2212 || !ivtx ) break; // reached beam
502 }
503 return isfromB;
504 */
505 }
506
507
508 // Check whether particle has charmed children
509 /// @todo Use built-in method and avoid HepMC!
510 bool hasCharmedChildren(ConstGenParticlePtr part) {
511
512 bool hasCharmedChild = false;
513 if (part == nullptr) return false;
514
515 ConstGenVertexPtr ivtx = part->end_vertex();
516 if (ivtx == nullptr) return false;
517
518 // if (ivtx->particles_out_size() < 2) return false;
519 //HepMC::GenVertex::particles_out_const_iterator iPart_invtx = ivtx->particles_out_const_begin();
520 //HepMC::GenVertex::particles_out_const_iterator end_invtx = ivtx->particles_out_const_end();
521
522 for(ConstGenParticlePtr p2: HepMCUtils::particles(ivtx, Relatives::CHILDREN)){
523 if (p2 == part) continue;
524 hasCharmedChild = PID::hasCharm(p2->pdg_id());
525 if (hasCharmedChild == true) break;
526 hasCharmedChild = hasCharmedChildren(p2);
527 if (hasCharmedChild == true) break;
528 }
529 return hasCharmedChild;
530 }
531
532
533 private:
534
535 /// @name Histograms
536 //@{
537
538 //d01-x01-
539 Histo1DPtr _hist_wcjet_charge;
540 Histo1DPtr _hist_wd_charge;
541 Histo1DPtr _hist_wdstar_charge;
542
543 //d02-x01-
544 Scatter2DPtr _hist_wcjet_ratio;
545 Scatter2DPtr _hist_wd_ratio;
546 Histo1DPtr _hist_wcjet_plus;
547 Histo1DPtr _hist_wd_plus;
548
549 Histo1DPtr _hist_wcjet_minus;
550 Histo1DPtr _hist_wd_minus;
551
552 //d03-x01-
553 Histo1DPtr _hist_wplus_wcjet_eta_lep;
554 Histo1DPtr _hist_wminus_wcjet_eta_lep;
555
556 //d04-x01-
557 Histo1DPtr _hist_wplus_wdminus_eta_lep;
558 Histo1DPtr _hist_wminus_wdplus_eta_lep;
559
560 //d05-x01-
561 Histo1DPtr _hist_wplus_wdstar_eta_lep;
562 Histo1DPtr _hist_wminus_wdstar_eta_lep;
563
564
565 // postprocessing histos
566 //d05-x01
567 Histo1DPtr _hist_w_inc;
568 Scatter2DPtr _hist_wd_winc_ratio;
569 Scatter2DPtr _hist_wdstar_winc_ratio;
570
571 //d06-x01
572 Histo1DPtr _hist_wplus_winc;
573 Histo1DPtr _hist_wminus_winc;
574
575 Scatter2DPtr _hist_wplusd_wplusinc_pt_ratio;
576 Scatter2DPtr _hist_wminusd_wminusinc_pt_ratio;
577 Scatter2DPtr _hist_wplusdstar_wplusinc_pt_ratio;
578 Scatter2DPtr _hist_wminusdstar_wminusinc_pt_ratio;
579
580 Histo1DPtr _hist_wplusd_wplusinc_pt ;
581 Histo1DPtr _hist_wminusd_wminusinc_pt;
582 Histo1DPtr _hist_wplusdstar_wplusinc_pt;
583 Histo1DPtr _hist_wminusdstar_wminusinc_pt;
584
585 // d07-x01
586 Histo1DPtr _hist_wcjet_jets ;
587
588 //d08-x01
589 Scatter2DPtr _hist_wcjet_jets_ratio ;
590 Histo1DPtr _hist_wcjet_jets_plus ;
591 Histo1DPtr _hist_wcjet_jets_minus;
592 //@}
593
594 };
595
596
597 // The hook for the plugin system
598 RIVET_DECLARE_PLUGIN(ATLAS_2014_I1282447);
599
600}
|