Rivet analyses referenceATLAS_2012_I1082009$D^{*\pm}$ production in jetsExperiment: ATLAS (LHC) Inspire ID: 1082009 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
Measurement of $D^{*\pm}$ meson production in jets from proton-proton collisions at a centre-of-mass energy of $\sqrt{s}=7$ TeV at the LHC. The measurement is based on a data sample recorded with the ATLAS detector with an integrated luminosity of $0.30\,\text{pb}^{-1}$ for jets with transverse momentum between 25 and 70 GeV in the pseudorapidity range $|eta| < 2.5$. Source code: ATLAS_2012_I1082009.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/IdentifiedFinalState.hh"
4#include "Rivet/Projections/VetoedFinalState.hh"
5#include "Rivet/Projections/MissingMomentum.hh"
6#include "Rivet/Projections/FastJets.hh"
7#include "Rivet/Projections/LeadingParticlesFinalState.hh"
8#include "Rivet/Projections/UnstableParticles.hh"
9
10namespace Rivet {
11
12
13 class ATLAS_2012_I1082009 : public Analysis {
14 public:
15
16 /// @name Constructors etc.
17 //@{
18
19 /// Constructor
20 ATLAS_2012_I1082009()
21 : Analysis("ATLAS_2012_I1082009")
22 { }
23
24 //@}
25
26
27 public:
28
29 /// @name Analysis methods
30 //@{
31
32 /// Book histograms and initialise projections before the run
33 void init() {
34
35 // Input for the jets: No neutrinos, no muons
36 VetoedFinalState veto;
37 veto.addVetoPairId(PID::MUON);
38 veto.vetoNeutrinos();
39 FastJets jets(veto, FastJets::ANTIKT, 0.6);
40 declare(jets, "jets");
41 // unstable final-state for D*
42 declare(UnstableParticles(), "UFS");
43
44 book(_weight25_30, "_weight_25_30");
45 book(_weight30_40, "_weight_30_40");
46 book(_weight40_50, "_weight_40_50");
47 book(_weight50_60, "_weight_50_60");
48 book(_weight60_70, "_weight_60_70");
49 book(_weight25_70, "_weight_25_70");
50
51 book(_h_pt25_30 , 8,1,1);
52 book(_h_pt30_40 , 9,1,1);
53 book(_h_pt40_50 ,10,1,1);
54 book(_h_pt50_60 ,11,1,1);
55 book(_h_pt60_70 ,12,1,1);
56 book(_h_pt25_70 ,13,1,1);
57 }
58
59
60 /// Perform the per-event analysis
61 void analyze(const Event& event) {
62 // get the jets
63 Jets jets;
64 for (const Jet& jet : apply<FastJets>(event, "jets").jetsByPt(25.0*GeV)) {
65 if ( jet.abseta() < 2.5 ) jets.push_back(jet);
66 }
67 // get the D* mesons
68 const UnstableParticles& ufs = apply<UnstableParticles>(event, "UFS");
69 Particles Dstar;
70 for (const Particle& p : ufs.particles()) {
71 const int id = p.abspid();
72 if(id==413) Dstar.push_back(p);
73 }
74
75 // loop over the jobs
76 for (const Jet& jet : jets ) {
77 double perp = jet.perp();
78 bool found = false;
79 double z(0.);
80 if(perp<25.||perp>70.) continue;
81 for(const Particle & p : Dstar) {
82 if(p.perp()<7.5) continue;
83 if(deltaR(p, jet.momentum())<0.6) {
84 Vector3 axis = jet.p3().unit();
85 z = axis.dot(p.p3())/jet.E();
86 if(z<0.3) continue;
87 found = true;
88 break;
89 }
90 }
91 _weight25_70->fill();
92 if(found) _h_pt25_70->fill(z);
93 if(perp>=25.&&perp<30.) {
94 _weight25_30->fill();
95 if(found) _h_pt25_30->fill(z);
96 }
97 else if(perp>=30.&&perp<40.) {
98 _weight30_40->fill();
99 if(found) _h_pt30_40->fill(z);
100 }
101 else if(perp>=40.&&perp<50.) {
102 _weight40_50->fill();
103 if(found) _h_pt40_50->fill(z);
104 }
105 else if(perp>=50.&&perp<60.) {
106 _weight50_60->fill();
107 if(found) _h_pt50_60->fill(z);
108 }
109 else if(perp>=60.&&perp<70.) {
110 _weight60_70->fill();
111 if(found) _h_pt60_70->fill(z);
112 }
113 }
114 }
115
116
117 /// Normalise histograms etc., after the run
118 void finalize() {
119 scale(_h_pt25_30,1./ *_weight25_30);
120 scale(_h_pt30_40,1./ *_weight30_40);
121 scale(_h_pt40_50,1./ *_weight40_50);
122 scale(_h_pt50_60,1./ *_weight50_60);
123 scale(_h_pt60_70,1./ *_weight60_70);
124 scale(_h_pt25_70,1./ *_weight25_70);
125 }
126
127 //@}
128
129
130 private:
131
132 /// @name Histograms
133 //@{
134 CounterPtr _weight25_30,_weight30_40,_weight40_50;
135 CounterPtr _weight50_60,_weight60_70,_weight25_70;
136
137 Histo1DPtr _h_pt25_30;
138 Histo1DPtr _h_pt30_40;
139 Histo1DPtr _h_pt40_50;
140 Histo1DPtr _h_pt50_60;
141 Histo1DPtr _h_pt60_70;
142 Histo1DPtr _h_pt25_70;
143 //@}
144
145 };
146
147 // The hook for the plugin system
148 RIVET_DECLARE_PLUGIN(ATLAS_2012_I1082009);
149
150}
|