rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2014_I1282447

W + charm production at 7 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1282447
Status: VALIDATED
Authors:
  • Kristin Lohwasser
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • W+c and inclusive double production in pp collisions.

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#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/UnstableParticles.hh"
  4#include "Rivet/Projections/ChargedFinalState.hh"
  5#include "Rivet/Projections/PromptFinalState.hh"
  6#include "Rivet/Projections/VetoedFinalState.hh"
  7#include "Rivet/Projections/MissingMomentum.hh"
  8#include "Rivet/Projections/LeptonFinder.hh"
  9#include "Rivet/Projections/FastJets.hh"
 10
 11namespace Rivet {
 12
 13
 14  /// @brief ATLAS W+c at 7 TeV
 15  ///
 16  /// This routine implements the ATLAS 7 TeV pp W+c analysis.
 17  ///
 18  /// Apart from those histograms described and published on HepData, there
 19  /// are some helper histograms:
 20  ///
 21  /// - d02-x01-y01, d02-x01-y02 and d08-x01-y01 are ratios, the nominator ("_plus")
 22  ///   and denominator ("_minus") histograms are also given, so that the ratios can
 23  ///   be reconstructed if need be (e.g. when running on separate samples).
 24  ///
 25  /// - d05 and d06 are ratios over inclusive W production.
 26  ///   The routine has to be run on a sample for inclusive W production in order to
 27  ///   make sure the denominator ("_winc") is correctly filled.
 28  ///
 29  /// The ratios can be constructed using the following sample code:
 30  ///
 31  /// ```python
 32  /// ## python divideWCharm.py
 33  ///
 34  /// import yoda
 35  /// hists_wc   = yoda.read("Rivet_Wc.yoda")
 36  /// hists_winc = yoda.read("Rivet_Winc.yoda")
 37  ///
 38  /// ## division histograms --> ONLY for different plus minus runs
 39  /// # (merge before using yodamerge Rivet_plus.yoda Rivet_minus.yoda > Rivet_Wc.yoda)
 40  ///
 41  /// d02y01_plus = hists_wc["/ATLAS_2014_I1282447/d02-x01-y01_plus"]
 42  /// d02y01_minus = hists_wc["/ATLAS_2014_I1282447/d02-x01-y01_minus"]
 43  /// ratio_d02y01 =  d02y01_plus.divide(d02y01_minus)
 44  /// ratio_d02y01.path = "/ATLAS_2014_I1282447/d02-x01-y01"
 45  ///
 46  /// d02y02_plus = hists_wc["/ATLAS_2014_I1282447/d02-x01-y02_plus"]
 47  /// d02y02_minus = hists_wc["/ATLAS_2014_I1282447/d02-x01-y02_minus"]
 48  /// ratio_d02y02=  d02y02_plus.divide(d02y02_minus)
 49  /// ratio_d02y02.path = "/ATLAS_2014_I1282447/d02-x01-y02"
 50  ///
 51  /// d08y01_plus = hists_wc["/ATLAS_2014_I1282447/d08-x01-y01_plus"]
 52  /// d08y01_minus = hists_wc["/ATLAS_2014_I1282447/d08-x01-y01_minus"]
 53  /// ratio_d08y01=  d08y01_plus.divide(d08y01_minus)
 54  /// ratio_d08y01.path = "/ATLAS_2014_I1282447/d08-x01-y01"
 55  ///
 56  /// # inclusive cross section
 57  /// h_winc = hists_winc["/ATLAS_2014_I1282447/d05-x01-y01"]
 58  /// h_d    = hists_wc["/ATLAS_2014_I1282447/d01-x01-y02"]
 59  /// h_dstar= hists_wc["/ATLAS_2014_I1282447/d01-x01-y03"]
 60  ///
 61  /// ratio_wd      =  h_d.divide(h_winc)
 62  /// ratio_wd.path = "/ATLAS_2014_I1282447/d05-x01-y02"
 63  ///
 64  /// ratio_wdstar      =  h_d.divide(h_winc)
 65  /// ratio_wdstar.path = "/ATLAS_2014_I1282447/d05-x01-y03"
 66  ///
 67  /// # pT differential
 68  /// h_winc_plus  = hists_winc["/ATLAS_2014_I1282447/d06-x01-y01_winc"]
 69  /// h_winc_minus = hists_winc["/ATLAS_2014_I1282447/d06-x01-y02_winc"]
 70  ///
 71  /// h_wd_plus      = hists_wc["/ATLAS_2014_I1282447/d06-x01-y01_wplus"]
 72  /// h_wd_minus     = hists_wc["/ATLAS_2014_I1282447/d06-x01-y02_wminus"]
 73  /// h_wdstar_plus  = hists_wc["/ATLAS_2014_I1282447/d06-x01-y03_wplus"]
 74  /// h_wdstar_minus = hists_wc["/ATLAS_2014_I1282447/d06-x01-y04_wminus"]
 75  ///
 76  /// ratio_wd_plus       =  h_wd_plus.divide(h_winc_plus)
 77  /// ratio_wd_plus.path  = "/ATLAS_2014_I1282447/d06-x01-y01"
 78  /// ratio_wd_minus      =  h_wd_plus.divide(h_winc_minus)
 79  /// ratio_wd_minus.path = "/ATLAS_2014_I1282447/d06-x01-y02"
 80  ///
 81  /// ratio_wdstar_plus       =  h_wdstar_plus.divide(h_winc_plus)
 82  /// ratio_wdstar_plus.path  = "/ATLAS_2014_I1282447/d06-x01-y03"
 83  /// ratio_wdstar_minus      =  h_wdstar_plus.divide(h_winc_minus)
 84  /// ratio_wdstar_minus.path = "/ATLAS_2014_I1282447/d06-x01-y04"
 85  ///
 86  /// ratio_wd_plus =  h_wd_plus.divide(h_winc_plus)
 87  /// ratio_wd_plus.path = "/ATLAS_2014_I1282447/d06-x01-y01"
 88  /// ratio_wd_minus =  h_wd_plus.divide(h_winc_minus)
 89  /// ratio_wd_minus.path = "/ATLAS_2014_I1282447/d06-x01-y02"
 90  ///
 91  /// h_winc_plus= hists_winc["/ATLAS_2014_I1282447/d06-x01-y01_winc"]
 92  /// h_winc_minus= hists_winc["/ATLAS_2014_I1282447/d06-x01-y02_winc"]
 93  ///
 94  /// ## copy other histograms for plotting
 95  ///
 96  /// d01x01y01= hists_wc["/ATLAS_2014_I1282447/d01-x01-y01"]
 97  /// d01x01y01.path = "/ATLAS_2014_I1282447/d01-x01-y01"
 98  ///
 99  /// d01x01y02= hists_wc["/ATLAS_2014_I1282447/d01-x01-y02"]
100  /// d01x01y02.path = "/ATLAS_2014_I1282447/d01-x01-y02"
101  ///
102  /// d01x01y03= hists_wc["/ATLAS_2014_I1282447/d01-x01-y03"]
103  /// d01x01y03.path = "/ATLAS_2014_I1282447/d01-x01-y03"
104  ///
105  /// d03x01y01= hists_wc["/ATLAS_2014_I1282447/d03-x01-y01"]
106  /// d03x01y01.path = "/ATLAS_2014_I1282447/d03-x01-y01"
107  ///
108  /// d03x01y02= hists_wc["/ATLAS_2014_I1282447/d03-x01-y02"]
109  /// d03x01y02.path = "/ATLAS_2014_I1282447/d03-x01-y02"
110  ///
111  /// d04x01y01= hists_wc["/ATLAS_2014_I1282447/d04-x01-y01"]
112  /// d04x01y01.path = "/ATLAS_2014_I1282447/d04-x01-y01"
113  ///
114  /// d04x01y02= hists_wc["/ATLAS_2014_I1282447/d04-x01-y02"]
115  /// d04x01y02.path = "/ATLAS_2014_I1282447/d04-x01-y02"
116  ///
117  /// d04x01y03= hists_wc["/ATLAS_2014_I1282447/d04-x01-y03"]
118  /// d04x01y03.path = "/ATLAS_2014_I1282447/d04-x01-y03"
119  ///
120  /// d04x01y04= hists_wc["/ATLAS_2014_I1282447/d04-x01-y04"]
121  /// d04x01y04.path = "/ATLAS_2014_I1282447/d04-x01-y04"
122  ///
123  /// d07x01y01= hists_wc["/ATLAS_2014_I1282447/d07-x01-y01"]
124  /// d07x01y01.path = "/ATLAS_2014_I1282447/d07-x01-y01"
125  ///
126  /// 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")
127  /// ```
128  ///
129  class ATLAS_2014_I1282447 : public Analysis {
130  public:
131
132    /// Constructor
133    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2014_I1282447);
134
135
136    /// @name Analysis methods
137    /// @{
138
139    /// Book histograms and initialise projections before the run
140    void init() {
141
142      // Initialise and register projections here
143      Cut cuts = Cuts::abseta < 2.5 && Cuts::pT > 20*GeV;
144
145      /// @note Should use sample WITHOUT QED radiation off the electron
146      ///
147      /// @warning This plus PhotonOrigin::ALL would mean dressing only with pi0 -> y y photons!?!
148      declare("MET", MissingMomentum());
149      LeptonFinder ef(0.1, cuts && Cuts::abspid == PID::ELECTRON);
150      declare(ef, "Elecs");
151      LeptonFinder mf(0.1, cuts && Cuts::abspid == PID::MUON);
152      declare(mf, "Muons");
153
154      // all hadrons that could be coming from a charm decay --
155      // -- for safety, use region -3.5 - 3.5
156      declare(UnstableParticles(Cuts::abseta < 3.5), "hadrons");
157
158      // Input for the jets: no neutrinos, no muons, and no electron which passed the electron cuts
159      // also: NO electron, muon or tau (needed due to ATLAS jet truth reconstruction feature)
160      VetoedFinalState veto;
161      veto.addVetoOnThisFinalState(ef);
162      veto.addVetoOnThisFinalState(mf);
163      veto.addVetoPairId(PID::ELECTRON);
164      veto.addVetoPairId(PID::MUON);
165      veto.addVetoPairId(PID::TAU);
166
167      FastJets jets(veto, JetAlg::ANTIKT, 0.4);
168      declare(jets, "jets");
169
170      // Book histograms
171
172      // charge separated integrated cross sections
173      book(_hist_wcjet_charge  ,"d01-x01-y01");
174      book(_hist_wd_charge     ,"d01-x01-y02");
175      book(_hist_wdstar_charge ,"d01-x01-y03");
176
177      // charge integrated total cross sections
178      book(_hist_wcjet_ratio,"d02-x01-y01");
179      book(_hist_wd_ratio   ,"d02-x01-y02");
180
181      book(_hist_wcjet_minus ,"d02-x01-y01_minus");
182      book(_hist_wd_minus    ,"d02-x01-y02_minus");
183
184      book(_hist_wcjet_plus  ,"d02-x01-y01_plus");
185      book(_hist_wd_plus     ,"d02-x01-y02_plus");
186
187      // eta distributions
188      book(_hist_wplus_wcjet_eta_lep   ,"d03-x01-y01");
189      book(_hist_wminus_wcjet_eta_lep  ,"d03-x01-y02");
190
191      book(_hist_wplus_wdminus_eta_lep ,"d04-x01-y01");
192      book(_hist_wminus_wdplus_eta_lep ,"d04-x01-y02");
193      book(_hist_wplus_wdstar_eta_lep  ,"d04-x01-y03");
194      book(_hist_wminus_wdstar_eta_lep ,"d04-x01-y04");
195
196      // ratio of cross section (WD over W inclusive) // postprocess!
197      book(_hist_w_inc             ,"d05-x01-y01");
198      book(_hist_wd_winc_ratio    ,"d05-x01-y02");
199      book(_hist_wdstar_winc_ratio,"d05-x01-y03");
200
201      // ratio of cross section (WD over W inclusive -- function of pT of D meson)
202      book(_hist_wplusd_wplusinc_pt_ratio      ,"d06-x01-y01");
203      book(_hist_wminusd_wminusinc_pt_ratio    ,"d06-x01-y02");
204      book(_hist_wplusdstar_wplusinc_pt_ratio  ,"d06-x01-y03");
205      book(_hist_wminusdstar_wminusinc_pt_ratio,"d06-x01-y04");
206
207      // could use for postprocessing!
208      book(_hist_wplusd_wplusinc_pt       ,"d06-x01-y01_wplus");
209      book(_hist_wminusd_wminusinc_pt     ,"d06-x01-y02_wminus");
210      book(_hist_wplusdstar_wplusinc_pt   ,"d06-x01-y03_wplus");
211      book(_hist_wminusdstar_wminusinc_pt ,"d06-x01-y04_wminus");
212
213      book(_hist_wplus_winc  ,"d06-x01-y01_winc");
214      book(_hist_wminus_winc ,"d06-x01-y02_winc");
215
216      // jet multiplicity of charge integrated W+cjet cross section (+0 or +1 jet in addition to the charm jet)
217      book(_hist_wcjet_jets  ,"d07-x01-y01");
218
219      // jet multiplicity of W+cjet cross section ratio (+0 or +1 jet in addition to the charm jet)
220      book(_hist_wcjet_jets_ratio ,"d08-x01-y01");
221      book(_hist_wcjet_jets_plus   ,"d08-x01-y01_plus");
222      book(_hist_wcjet_jets_minus  ,"d08-x01-y01_minus");
223    }
224
225
226    /// Perform the per-event analysis
227    void analyze(const Event& event) {
228
229      double charge_weight = 0; // account for OS/SS events
230      int    lepton_charge = 0;
231      double lepton_eta    = 0.;
232
233      /// Find leptons
234      const P4& pmiss = apply<MissingMom>(event, "MET").missingMom();
235      if (pmiss.Et() < 25*GeV) vetoEvent;
236
237      // Identify the closest-matching l+MET to m == mW
238      const Particles& es = apply<LeptonFinder>(event, "Elecs").particles();
239      const int iefound = closestMatchIndex(es, pmiss, Kin::mass, 80.4*GeV, 15*GeV);
240      const Particles& mus = apply<LeptonFinder>(event, "Muons").particles();
241      const int imfound = closestMatchIndex(mus, pmiss, Kin::mass, 80.4*GeV, 15*GeV);
242      if (iefound < 0 && imfound < 0) {
243        MSG_DEBUG("No W bosons found");
244        vetoEvent;
245      }
246      if (iefound >= 0 && imfound >= 0) {
247        MSG_DEBUG("Multiple W bosons found");
248        vetoEvent;
249      }
250
251      bool keepevent = false;
252
253      //check electrons
254      if (iefound >= 0) {
255        const Particle& el = es[iefound];
256        if (mT(pmiss, el) > 40*GeV) {
257          keepevent = true;
258          lepton_charge = el.charge();
259          lepton_eta = el.abseta();
260        }
261      }
262
263      //check muons
264      if (imfound >= 0) {
265        const Particle& mu = mus[imfound];
266        if (mT(pmiss, mu) > 40*GeV) {
267          keepevent = true;
268          lepton_charge = mu.charge();
269          lepton_eta = mu.abseta();
270        }
271      }
272
273      if (!keepevent) {
274        MSG_DEBUG("Event does not pass mT cuts");
275        vetoEvent;
276      }
277
278      if (lepton_charge > 0) {
279        _hist_wplus_winc->fill(10.);
280        _hist_wplus_winc->fill(16.);
281        _hist_wplus_winc->fill(30.);
282        _hist_wplus_winc->fill(60.);
283        _hist_w_inc->fill(+1);
284      }
285      else if (lepton_charge < 0) {
286        _hist_wminus_winc->fill(10.);
287        _hist_wminus_winc->fill(16.);
288        _hist_wminus_winc->fill(30.);
289        _hist_wminus_winc->fill(60.);
290        _hist_w_inc->fill(-1);
291      }
292
293      // Find hadrons in the event
294      const UnstableParticles& fs = apply<UnstableParticles>(event, "hadrons");
295
296      /// FIND Different channels
297      // 1: wcjet
298      // get jets
299      const Jets& jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 25.0*GeV && Cuts::abseta < 2.5);
300      // loop over jets to select jets used to match to charm
301      Jets js;
302      int    matched_charmHadron = 0;
303      double charm_charge = 0.;
304      int    njets = 0;
305      int    nj = 0;
306      bool   mat_jet = false;
307
308      double ptcharm = 0;
309      if (matched_charmHadron > -1) {
310        for (const Jet& j : jets) {
311          mat_jet = false;
312
313          njets += 1;
314          for (const Particle& p : fs.particles()) {
315            /// @todo Avoid touching HepMC!
316            ConstGenParticlePtr part = p.genParticle();
317            if (p.hasCharm()) {
318              //if(isFromBDecay(p)) continue;
319              if (p.fromBottom()) continue;
320              if (p.pT() < 5*GeV ) continue;
321              if (hasCharmedChildren(part)) continue;
322              if (deltaR(p, j) < 0.3) {
323                mat_jet = true;
324                if (p.pT() > ptcharm) {
325                  charm_charge = part->pdg_id();
326                  ptcharm = p.pT();
327                }
328              }
329            }
330          }
331          if (mat_jet) nj++;
332        }
333
334        if (charm_charge * lepton_charge > 0)  charge_weight = -1;
335        else charge_weight = +1;
336
337        if (nj == 1)  {
338          if (lepton_charge > 0) {
339            _hist_wcjet_charge        ->fill(         1, charge_weight);
340            _hist_wcjet_plus          ->fill(         0, charge_weight);
341            _hist_wplus_wcjet_eta_lep ->fill(lepton_eta, charge_weight);
342            _hist_wcjet_jets_plus     ->fill(njets-1   , charge_weight);
343          }
344          else if (lepton_charge < 0) {
345            _hist_wcjet_charge        ->fill(        -1, charge_weight);
346            _hist_wcjet_minus         ->fill(         0, charge_weight);
347            _hist_wminus_wcjet_eta_lep->fill(lepton_eta, charge_weight);
348            _hist_wcjet_jets_minus    ->fill(njets-1   , charge_weight);
349          }
350
351          _hist_wcjet_jets->fill(njets-1, charge_weight);
352        }
353      }
354
355      // // 1/2: w+d(*) meson
356
357      for (const Particle& p : fs.particles()) {
358
359        /// @todo Avoid touching HepMC!
360        ConstGenParticlePtr part = p.genParticle();
361        if (p.pT() < 8*GeV)       continue;
362        if (fabs(p.eta()) > 2.2)  continue;
363
364        // W+D
365        if (abs(part->pdg_id()) == 411) {
366          if (lepton_charge * part->pdg_id() > 0)  charge_weight = -1;
367          else                                     charge_weight = +1;
368
369          // fill histos
370          if (lepton_charge > 0) {
371            _hist_wd_charge            ->fill(         1, charge_weight);
372            _hist_wd_plus              ->fill(         0, charge_weight);
373            _hist_wplus_wdminus_eta_lep->fill(lepton_eta, charge_weight);
374            _hist_wplusd_wplusinc_pt   ->fill(    p.pT(), charge_weight);
375          }
376          else if (lepton_charge < 0) {
377            _hist_wd_charge            ->fill(        -1, charge_weight);
378            _hist_wd_minus             ->fill(         0, charge_weight);
379            _hist_wminus_wdplus_eta_lep->fill(lepton_eta, charge_weight);
380            _hist_wminusd_wminusinc_pt ->fill(p.pT()    , charge_weight);
381          }
382        }
383
384        // W+Dstar
385        if ( abs(part->pdg_id()) == 413 ) {
386          if (lepton_charge*part->pdg_id() > 0) charge_weight = -1;
387          else charge_weight = +1;
388
389          if (lepton_charge > 0) {
390            _hist_wdstar_charge->fill(+1, charge_weight);
391            _hist_wd_plus->fill( 0, charge_weight);
392            _hist_wplus_wdstar_eta_lep->fill( lepton_eta, charge_weight);
393            _hist_wplusdstar_wplusinc_pt->fill(  p.pT(), charge_weight);
394          }
395          else if (lepton_charge < 0) {
396            _hist_wdstar_charge->fill(-1, charge_weight);
397            _hist_wd_minus->fill(0, charge_weight);
398            _hist_wminus_wdstar_eta_lep->fill(lepton_eta, charge_weight);
399            _hist_wminusdstar_wminusinc_pt->fill(p.pT(), charge_weight);
400          }
401        }
402      }
403
404    }
405
406
407    /// Normalise histograms etc., after the run
408    void finalize() {
409
410      const double sf = crossSection()/picobarn / sumOfWeights();
411
412      // norm to cross section
413      // d01
414      scale(_hist_wcjet_charge,  sf);
415      scale(_hist_wd_charge,     sf);
416      scale(_hist_wdstar_charge, sf);
417
418      //d02
419      scale(_hist_wcjet_plus,  sf);
420      scale(_hist_wcjet_minus, sf);
421      scale(_hist_wd_plus,     sf);
422      scale(_hist_wd_minus,    sf);
423
424      divide(_hist_wcjet_plus, _hist_wcjet_minus, _hist_wcjet_ratio);
425      divide(_hist_wd_plus,    _hist_wd_minus,    _hist_wd_ratio   );
426
427      //d03
428      scale(_hist_wplus_wcjet_eta_lep,  sf);
429      scale(_hist_wminus_wcjet_eta_lep, sf);
430
431      //d04
432      scale(_hist_wplus_wdminus_eta_lep, crossSection()/picobarn/sumOfWeights());
433      scale(_hist_wminus_wdplus_eta_lep, crossSection()/picobarn/sumOfWeights());
434
435      scale(_hist_wplus_wdstar_eta_lep , crossSection()/picobarn/sumOfWeights());
436      scale(_hist_wminus_wdstar_eta_lep, crossSection()/picobarn/sumOfWeights());
437
438      //d05
439      scale(_hist_w_inc, 0.01 * sf); // in percent --> /100
440      divide(_hist_wd_charge,     _hist_w_inc, _hist_wd_winc_ratio    );
441      divide(_hist_wdstar_charge, _hist_w_inc, _hist_wdstar_winc_ratio);
442
443      //d06, in percentage!
444      scale(_hist_wplusd_wplusinc_pt,       sf);
445      scale(_hist_wminusd_wminusinc_pt,     sf);
446      scale(_hist_wplusdstar_wplusinc_pt,   sf);
447      scale(_hist_wminusdstar_wminusinc_pt, sf);
448
449      scale(_hist_wplus_winc,  0.01 * sf); // in percent --> /100
450      scale(_hist_wminus_winc, 0.01 * sf); // in percent --> /100
451
452      divide(_hist_wplusd_wplusinc_pt,       _hist_wplus_winc , _hist_wplusd_wplusinc_pt_ratio      );
453      divide(_hist_wminusd_wminusinc_pt,     _hist_wminus_winc, _hist_wminusd_wminusinc_pt_ratio    );
454      divide(_hist_wplusdstar_wplusinc_pt,   _hist_wplus_winc , _hist_wplusdstar_wplusinc_pt_ratio  );
455      divide(_hist_wminusdstar_wminusinc_pt, _hist_wminus_winc, _hist_wminusdstar_wminusinc_pt_ratio);
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    ///
475    /// @note Slight difference wrt Rivet's native Particle::fromBottom method!
476    bool isFromBDecay(const Particle& p) {
477      return p.hasAncestorWith([](const Particle& p)->bool{return p.hasBottom();});
478    }
479
480    /// Check whether particle has charmed children
481    ///
482    /// @todo Use built-in method and avoid HepMC!
483    bool hasCharmedChildren(ConstGenParticlePtr part) {
484
485      bool hasCharmedChild = false;
486      if (part == nullptr)  return false;
487
488      ConstGenVertexPtr ivtx = part->end_vertex();
489      if (ivtx == nullptr)  return false;
490
491      for (ConstGenParticlePtr p2: HepMCUtils::particles(ivtx, Relatives::CHILDREN)) {
492        if (p2 == part)  continue;
493        hasCharmedChild = PID::hasCharm(p2->pdg_id());
494        if (hasCharmedChild == true)  break;
495        hasCharmedChild = hasCharmedChildren(p2);
496        if (hasCharmedChild == true)  break;
497      }
498      return hasCharmedChild;
499    }
500
501
502  private:
503
504    /// @name Histograms
505    /// @{
506
507    //d01-x01-
508    Histo1DPtr   _hist_wcjet_charge;
509    Histo1DPtr   _hist_wd_charge;
510    Histo1DPtr   _hist_wdstar_charge;
511
512    //d02-x01-
513    Estimate1DPtr _hist_wcjet_ratio;
514    Estimate1DPtr _hist_wd_ratio;
515    Histo1DPtr _hist_wcjet_plus;
516    Histo1DPtr _hist_wd_plus;
517
518    Histo1DPtr _hist_wcjet_minus;
519    Histo1DPtr _hist_wd_minus;
520
521    //d03-x01-
522    Histo1DPtr _hist_wplus_wcjet_eta_lep;
523    Histo1DPtr _hist_wminus_wcjet_eta_lep;
524
525    //d04-x01-
526    Histo1DPtr _hist_wplus_wdminus_eta_lep;
527    Histo1DPtr _hist_wminus_wdplus_eta_lep;
528
529    //d05-x01-
530    Histo1DPtr _hist_wplus_wdstar_eta_lep;
531    Histo1DPtr _hist_wminus_wdstar_eta_lep;
532
533    // postprocessing histos
534    //d05-x01
535    Histo1DPtr _hist_w_inc;
536    Estimate1DPtr _hist_wd_winc_ratio;
537    Estimate1DPtr _hist_wdstar_winc_ratio;
538
539    //d06-x01
540    Histo1DPtr _hist_wplus_winc;
541    Histo1DPtr _hist_wminus_winc;
542
543    Estimate1DPtr _hist_wplusd_wplusinc_pt_ratio;
544    Estimate1DPtr _hist_wminusd_wminusinc_pt_ratio;
545    Estimate1DPtr _hist_wplusdstar_wplusinc_pt_ratio;
546    Estimate1DPtr _hist_wminusdstar_wminusinc_pt_ratio;
547
548    Histo1DPtr _hist_wplusd_wplusinc_pt ;
549    Histo1DPtr _hist_wminusd_wminusinc_pt;
550    Histo1DPtr _hist_wplusdstar_wplusinc_pt;
551    Histo1DPtr _hist_wminusdstar_wminusinc_pt;
552
553    // d07-x01
554    Histo1DPtr _hist_wcjet_jets ;
555
556    //d08-x01
557    Estimate1DPtr  _hist_wcjet_jets_ratio ;
558    Histo1DPtr    _hist_wcjet_jets_plus ;
559    Histo1DPtr    _hist_wcjet_jets_minus;
560    /// @}
561
562  };
563
564
565  RIVET_DECLARE_PLUGIN(ATLAS_2014_I1282447);
566
567}