rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ALEPH_2000_I507531

Measurement of $\pi^0$, $\eta$, $\eta^\prime$, $K^0_S$ and $\Lambda^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: ANY
Run details:
  • e+e- -> hadrons

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  // The hook for the plugin system
225  RIVET_DECLARE_PLUGIN(ALEPH_2000_I507531);
226
227
228}