Rivet analyses referenceALEPH_2000_I507531Measurement of $\pi^0$, $\eta$, $\eta^\prime$, $K^0_S$ and $\Lambda^0$ spectra in two and three jet eventsExperiment: ALEPH (LEP) Inspire ID: 507531 Status: VALIDATED Authors:
Beam energies: (45.6, 45.6) GeV Run details:
Measurement of $\pi^0$, $\eta$, $\eta^\prime$, $K^0_S$ and $\Lambda^0$ spectra in two and three jet events. In addition to the normal inclusive spectra the spectra in individual jets are measured for three jet events. Source code: ALEPH_2000_I507531.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/UnstableParticles.hh"
4#include "Rivet/Projections/FinalState.hh"
5#include "Rivet/Projections/FastJets.hh"
6#include "Rivet/Projections/Beam.hh"
7
8namespace Rivet {
9
10
11 /// @brief pi, eta, eta', K0, lambda spectra
12 class ALEPH_2000_I507531 : public Analysis {
13 public:
14
15 /// Constructor
16 RIVET_DEFAULT_ANALYSIS_CTOR(ALEPH_2000_I507531);
17
18
19 /// @name Analysis methods
20 /// @{
21
22 /// Book histograms and initialise projections before the run
23 void init() {
24
25 // Projections
26 declare(Beam(), "Beams");
27 declare(UnstableParticles(), "UFS");
28 declare(FinalState() , "FS");
29
30 // Histograms
31 // incl
32 book(_h_pi0 , 1,1,1);
33 book(_h_eta , 2,1,1);
34 book(_h_etaP, 3,1,1);
35 book(_h_K0 , 16,1,1);
36 book(_h_lam , 17,1,1);
37 // two jet
38 book(_h_2_pi0 , 4,1,1);
39 book(_h_2_eta , 5,1,1);
40 book(_h_2_etaP, 6,1,1);
41 book(_h_2_K0 , 18,1,1);
42 book(_h_2_lam , 19,1,1);
43 // three jet
44 book(_h_3_pi0 [0], 7,1,1);
45 book(_h_3_pi0 [1], 8,1,1);
46 book(_h_3_pi0 [2], 9,1,1);
47 book(_h_3_eta [0], 10,1,1);
48 book(_h_3_eta [1], 11,1,1);
49 book(_h_3_eta [2], 12,1,1);
50 book(_h_3_etaP[0], 13,1,1);
51 book(_h_3_etaP[1], 14,1,1);
52 book(_h_3_etaP[2], 15,1,1);
53 book(_h_3_K0 [0], 20,1,1);
54 book(_h_3_K0 [1], 21,1,1);
55 book(_h_3_K0 [2], 22,1,1);
56 book(_h_3_lam [0], 23,1,1);
57 book(_h_3_lam [1], 24,1,1);
58 book(_h_3_lam [2], 25,1,1);
59 book(_w2,"/TMP/W2");
60 book(_w3,"/TMP/W3");
61 }
62
63 void findDecayProducts(const Particle & parent, Particles & decay) {
64 for(const Particle & child : parent.children()) {
65 if(child.children().empty()) {
66 decay.push_back(child);
67 }
68 else {
69 findDecayProducts(child,decay);
70 }
71 }
72 }
73
74 /// Perform the per-event analysis
75 void analyze(const Event& event) {
76 // Get beams and average beam momentum
77 const ParticlePair& beams = apply<Beam>(event, "Beams").beams();
78 const double meanBeamMom = ( beams.first.p3().mod() +
79 beams.second.p3().mod() ) / 2.0;
80
81 Particles decay,fs;
82 // unstable particles
83 const UnstableParticles ufs = apply<UnstableParticles>(event, "UFS");
84 for(const Particle & part : ufs.particles(Cuts::pid==111 or Cuts::pid==221 or Cuts::pid==331 or
85 Cuts::pid==310 or Cuts::abspid==3122)) {
86 fs.push_back(part);
87 findDecayProducts(part,decay);
88 }
89 // FS particles
90 for(const Particle & part : apply<FinalState>(event, "FS").particles()) {
91 bool skip=false;
92 for(const Particle & dec :decay) {
93 if(dec.genParticle()==part.genParticle()) {
94 skip=true;
95 break;
96 }
97 }
98 if(skip) continue;
99 fs.push_back(part);
100 }
101 // Definition of the Durham algorithm
102 fastjet::JetDefinition durham_def(fastjet::ee_kt_algorithm, fastjet::E_scheme, fastjet::Best);
103 // pseudojets
104 vector<fastjet::PseudoJet> input_particles;
105 // Pseudo-jets from the non photons
106 unsigned int ix=0;
107 for (const Particle& p : fs ) {
108 const FourMomentum p4 = p.momentum();
109 input_particles.push_back(fastjet::PseudoJet(p4.px(), p4.py(), p4.pz(), p4.E()));
110 input_particles.back().set_user_index(ix);
111 ++ix;
112 }
113 // cluster the jets
114 fastjet::ClusterSequence clust_seq(input_particles, durham_def);
115 PseudoJets jets = fastjet::sorted_by_E(clust_seq.exclusive_jets_ycut(0.01));
116 if(jets.size()==2) _w2->fill();
117 else if(jets.size()==3) _w3->fill();
118 ix=0;
119 for(const Particle & part : fs) {
120 double xE = part.momentum().E()/meanBeamMom;
121 double xP = part.momentum().p3().mod()/meanBeamMom;
122 int ijet = jets.size()!=3 ? -1 : findJet(ix,jets);
123 if(part.pid()==111) {
124 _h_pi0->fill(xE);
125 if(jets.size()==2) {
126 _h_2_pi0->fill(xE);
127 }
128 else if(jets.size()==3) {
129 _h_3_pi0[ijet]->fill(xE);
130 }
131 }
132 else if(part.pid()==221) {
133 _h_eta->fill(xE);
134 if(jets.size()==2) {
135 _h_2_eta->fill(xE);
136 }
137 else if(jets.size()==3) {
138 _h_3_eta[ijet]->fill(xE);
139 }
140 }
141 else if(part.pid()==331) {
142 _h_etaP->fill(xE);
143 if(jets.size()==2) {
144 _h_2_etaP->fill(xE);
145 }
146 else if(jets.size()==3) {
147 _h_3_etaP[ijet]->fill(xE);
148 }
149 }
150 else if(part.pid()==310) {
151 double xi=-log(xP);
152 _h_K0->fill(xi);
153 if(jets.size()==2) {
154 _h_2_K0->fill(xi);
155 }
156 else if(jets.size()==3) {
157 _h_3_K0[ijet]->fill(xi);
158 }
159 }
160 else if(part.abspid()==3122) {
161 double xi=-log(xP);
162 _h_lam->fill(xi);
163 if(jets.size()==2) {
164 _h_2_lam->fill(xi);
165 }
166 else if(jets.size()==3) {
167 _h_3_lam[ijet]->fill(xi);
168 }
169 }
170 else {
171 break;
172 }
173 ix+=1;
174 }
175 }
176
177 int findJet(int id, const PseudoJets & jets) {
178 for(unsigned int ijet=0;ijet<jets.size();++ijet) {
179 for(const PseudoJet & con : jets[ijet].constituents()) {
180 if(con.user_index()==id)
181 return ijet;
182 }
183 }
184 return -1;
185 }
186
187
188 /// Normalise histograms etc., after the run
189 void finalize() {
190 scale(_h_pi0 ,1./sumOfWeights());
191 scale(_h_eta ,1./sumOfWeights());
192 scale(_h_etaP ,1./sumOfWeights());
193 scale(_h_K0 ,1./sumOfWeights());
194 scale(_h_lam ,1./sumOfWeights());
195 scale(_h_2_pi0 ,1./ *_w2);
196 scale(_h_2_eta ,1./ *_w2);
197 scale(_h_2_etaP ,1./ *_w2);
198 scale(_h_2_K0 ,1./ *_w2);
199 scale(_h_2_lam ,1./ *_w2);
200 for(unsigned int ix=0;ix<3;++ix ) {
201 scale(_h_3_pi0[ix] ,1./ *_w3);
202 scale(_h_3_eta[ix] ,1./ *_w3);
203 scale(_h_3_etaP[ix],1./ *_w3);
204 scale(_h_3_K0[ix] ,1./ *_w3);
205 scale(_h_3_lam[ix] ,1./ *_w3);
206 }
207 }
208
209 /// @}
210
211
212 /// @name Histograms
213 /// @{
214 Histo1DPtr _h_pi0 , _h_eta , _h_etaP , _h_K0 , _h_lam ;
215 Histo1DPtr _h_2_pi0 , _h_2_eta , _h_2_etaP , _h_2_K0 , _h_2_lam ;
216 Histo1DPtr _h_3_pi0[3], _h_3_eta[3], _h_3_etaP[3], _h_3_K0[3], _h_3_lam[3];
217 CounterPtr _w2,_w3;
218 /// @}
219
220
221 };
222
223
224 RIVET_DECLARE_PLUGIN(ALEPH_2000_I507531);
225
226
227}
|