Processing math: 100%
rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ALEPH_2000_I507531

Measurement of π0, η, η, K0S and Λ0 spectra in two and three jet events
Experiment: ALEPH (LEP)
Inspire ID: 507531
Status: VALIDATED
Authors:
  • Peter Richardson
References:
  • Eur.Phys.J. C16 (2000) 613
Beams: e- e+
Beam energies: (45.6, 45.6) GeV
Run details:
  • e+e- -> hadrons

Measurement of π0, η, η, K0S and Λ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}