rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

OPAL_2000_I474010

Correlations for $\Lambda^0$ and $\bar\Lambda^0$ production at LEP1
Experiment: OPAL (LEP)
Inspire ID: 474010
Status: VALIDATED
Authors:
  • Peter Richardson
References:
  • Eur.Phys.J. C13 (2000) 185-195
Beams: e- e+
Beam energies: ANY
Run details:
  • e+e- to hadrons

Measurement of the correlations between the production of $\Lambda^0$ and $\bar\Lambda^0$ baryons at LEP1, which is useful for testing models of baryon production.

Source code: OPAL_2000_I474010.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/Thrust.hh"
  5#include "Rivet/Projections/FastJets.hh"
  6#include "Rivet/Projections/UnstableParticles.hh"
  7
  8
  9namespace {
 10
 11unsigned int factorial(unsigned int n) {
 12  if(n<=1) return 0;
 13  unsigned int out=1;
 14  for (unsigned int i=n;i>0;--i) out *= i;
 15  return out;
 16}
 17
 18double rapidityAxis(Rivet::Vector4 momentum, Rivet::Vector3 axis) {
 19  // assume axis is unit vector
 20  double plong = momentum.vector3().dot(axis);
 21  return 0.5*log((momentum.t()+plong)/(momentum.t()-plong));
 22}
 23
 24double cosThetaStar(const Rivet::Particle & p1, const Rivet::Particle & p2, Rivet::Vector3 axis) {
 25  Rivet::FourMomentum psum = p1.momentum()+p2.momentum();
 26  Rivet::LorentzTransform trans(Rivet::LorentzTransform::mkObjTransformFromBeta(-psum.betaVec()));
 27  Rivet::FourMomentum m1 = trans.transform(p1.momentum());
 28  Rivet::FourVector a4;
 29  a4.setX(axis.x());
 30  a4.setY(axis.y());
 31  a4.setZ(axis.z());
 32  a4.setT(1.);
 33  a4 = trans.transform(a4);
 34  return fabs(m1.vector3().dot(a4.vector3())/m1.vector3().mod()/a4.vector3().mod());
 35}
 36
 37}
 38
 39namespace Rivet {
 40
 41
 42  /// @brief lambda anti-lambda correlations
 43  class OPAL_2000_I474010 : public Analysis {
 44  public:
 45
 46    /// Constructor
 47    RIVET_DEFAULT_ANALYSIS_CTOR(OPAL_2000_I474010);
 48
 49
 50    /// @name Analysis methods
 51    //@{
 52
 53    /// Book histograms and initialise projections before the run
 54    void init() {
 55      // projections
 56      const FinalState fs;
 57      declare(fs, "FS");
 58      FastJets durhamjets(fs, FastJets::DURHAM, 0.7);
 59      durhamjets.useInvisibles(JetAlg::Invisibles::ALL);
 60      declare(durhamjets, "DurhamJets");
 61      const Thrust thrust(fs);
 62      declare(thrust, "Thrust");
 63      declare(UnstableParticles(), "UFS");
 64      // histograms
 65      book(_histLLBar    ,  1, 1, 1);
 66      book(_histLL       ,  1, 1, 2);
 67      book(_histLLBarCorr,  1, 1, 3);
 68      book(_histLLBar_2jet    ,  2, 1, 1);
 69      book(_histLL_2jet       ,  2, 1, 2);
 70      book(_histLLBarCorr_2jet,  2, 1, 3);
 71      book(_histLLBar_3jet    ,  3, 1, 1);
 72      book(_histLL_3jet       ,  3, 1, 2);
 73      book(_histLLBarCorr_3jet,  3, 1, 3);
 74      book(_histdcosTheta,  4, 1, 1);
 75      book(_histdy       ,  5, 1, 1);
 76      // weights
 77      book(_weight2Jet,"/TMP/W2Jet");
 78      book(_weight3Jet,"/TMP/W3Jet");
 79
 80    }
 81
 82
 83    /// Perform the per-event analysis
 84    void analyze(const Event& event) {
 85      // First, veto on leptonic events by requiring at least 4 charged FS particles
 86      const FinalState& fs = applyProjection<FinalState>(event, "FS");
 87      const size_t numParticles = fs.particles().size();
 88
 89      // Even if we only generate hadronic events, we still need a cut on numCharged >= 2.
 90      if (numParticles < 2) {
 91        MSG_DEBUG("Failed leptonic event cut");
 92        vetoEvent;
 93      }
 94      MSG_DEBUG("Passed leptonic event cut");
 95
 96      const Thrust& thrust = applyProjection<Thrust>(event, "Thrust");
 97      Vector3 axis = thrust.thrustAxis();
 98
 99      // Final state of unstable particles to get particle spectra
100      const UnstableParticles & ufs = applyProjection<UnstableParticles>(event, "UFS");
101
102      Particles lambda,lambdaBar;
103      for (const Particle& p : ufs.particles()) {
104        const int id = p.pid();
105	if(id==3122)
106	  lambda.push_back(p);
107	else if(id==-3122)
108	  lambdaBar.push_back(p);
109      }
110
111      // get the number of jets
112      const FastJets& durjet = applyProjection<FastJets>(event, "DurhamJets");
113
114      unsigned int njet = durjet.clusterSeq()->n_exclusive_jets_ycut(0.005);
115      if(njet==2)
116	_weight2Jet->fill();
117      else if(njet==3)
118	_weight3Jet->fill();
119
120      // need a pair of baryons
121      if(lambda.size()+lambdaBar.size()<2)
122	vetoEvent;
123
124      // multiplicities
125      int nsame = (factorial(lambda.size())+factorial(lambdaBar.size()))/2;
126      int ndiff = lambda.size()*lambdaBar.size();
127      _histLLBar    ->fill(_histLLBar    ->bin(0).xMid(), ndiff       );
128      _histLL       ->fill(_histLL       ->bin(0).xMid(),       nsame );
129      _histLLBarCorr->fill(_histLLBarCorr->bin(0).xMid(),(ndiff-nsame));
130      
131      // multiplicities for different numbers of jets
132      if(njet==2) {
133	_histLLBar_2jet    ->fill(_histLLBar_2jet    ->bin(0).xMid(), ndiff       );
134	_histLL_2jet       ->fill(_histLL_2jet       ->bin(0).xMid(),       nsame );
135	_histLLBarCorr_2jet->fill(_histLLBarCorr_2jet->bin(0).xMid(),(ndiff-nsame));
136      }
137      else if(njet==3) {
138	_histLLBar_3jet    ->fill(_histLLBar_3jet    ->bin(0).xMid(), ndiff       );
139	_histLL_3jet       ->fill(_histLL_3jet       ->bin(0).xMid(),       nsame );
140	_histLLBarCorr_3jet->fill(_histLLBarCorr_3jet->bin(0).xMid(),(ndiff-nsame));
141      }
142      // uncorrelated lambda
143      for(unsigned int ix=0;ix<lambda.size();++ix) {
144	double y1 = rapidityAxis(lambda[ix].momentum(),axis);
145	for(unsigned int iy=ix+1;iy<lambda.size();++iy) {
146	  double y2 = rapidityAxis(lambda[iy].momentum(),axis);
147	  if(njet==2) _histdy->fill(fabs(y1-y2),-1.0);
148	  double ctheta = cosThetaStar(lambda[ix],lambda[iy],axis);
149	  _histdcosTheta->fill(ctheta,-1.0);
150	}
151      }
152      // uncorrelated lambdabar
153      for(unsigned int ix=0;ix<lambdaBar.size();++ix) {
154	double y1 = rapidityAxis(lambdaBar[ix].momentum(),axis);
155	for(unsigned int iy=ix+1;iy<lambdaBar.size();++iy) {
156	  double y2 = rapidityAxis(lambdaBar[iy].momentum(),axis);
157	  if(njet==2) _histdy->fill(fabs(y1-y2),-1.0);
158	  double ctheta = cosThetaStar(lambdaBar[ix],lambdaBar[iy],axis);
159	  _histdcosTheta->fill(ctheta,-1.0);
160	}
161      }
162      // all lambda lambdabar
163      for(unsigned int ix=0;ix<lambda.size();++ix) {
164	double y1 = rapidityAxis(lambda[ix].momentum(),axis);
165	for(unsigned int iy=0;iy<lambdaBar.size();++iy) {
166	  double y2 = rapidityAxis(lambdaBar[iy].momentum(),axis);
167	  if(njet==2) _histdy->fill(fabs(y1-y2));
168	  double ctheta = cosThetaStar(lambda[ix],lambdaBar[iy],axis);
169	  _histdcosTheta->fill(ctheta);
170	}
171      }
172    }
173
174
175    /// Normalise histograms etc., after the run
176    void finalize() {
177      scale(_histLLBar,1./sumOfWeights());
178      scale(_histLL,1./sumOfWeights());
179      scale(_histLLBarCorr,1./sumOfWeights());
180      scale(_histLLBar_2jet,1./ *_weight2Jet);
181      scale(_histLL_2jet,1./ *_weight2Jet);
182      scale(_histLLBarCorr_2jet,1./ *_weight2Jet);
183      scale(_histLLBar_3jet,1./ *_weight3Jet);
184      scale(_histLL_3jet,1./ *_weight3Jet);
185      scale(_histLLBarCorr_3jet,1./ *_weight3Jet);
186      normalize(_histdcosTheta, 1.);
187      normalize(_histdy       , 1.);
188    }
189
190    //@}
191
192
193    /// @name Histograms
194    //@{
195
196    CounterPtr _weight2Jet;
197    CounterPtr _weight3Jet;
198
199    Histo1DPtr _histLLBar;
200    Histo1DPtr _histLL;
201    Histo1DPtr _histLLBarCorr;
202
203    Histo1DPtr _histLLBar_2jet;
204    Histo1DPtr _histLL_2jet;
205    Histo1DPtr _histLLBarCorr_2jet;
206
207    Histo1DPtr _histLLBar_3jet;
208    Histo1DPtr _histLL_3jet;
209    Histo1DPtr _histLLBarCorr_3jet;
210
211    Histo1DPtr _histdcosTheta;
212    Histo1DPtr _histdy;
213    //@}
214
215
216  };
217
218
219  // The hook for the plugin system
220  RIVET_DECLARE_PLUGIN(OPAL_2000_I474010);
221
222
223}