rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

OPAL_1997_S3608263

$K^{*0}$ meson production measured by OPAL at LEP 1.
Experiment: OPAL (LEP 1)
Inspire ID: 447146
Status: VALIDATED
Authors:
  • Peter Richardson
References:
  • Phys.Lett.B412:210-224,1997
  • hep-ex/9708022
Beams: e+ e-
Beam energies: (45.6, 45.6) GeV
Run details:
  • Hadronic Z decay events generated on the Z pole ($\sqrt{s} = 91.2$ GeV)

The $K^{*0}$ fragmentation function has been measured in hadronic $Z^0$ decays. In addition the helicity density matrix elements for inclusive $K^*(892)^0$ mesons from hadronic $Z^0$ decays have been measured over the full range of $K^{*0}$ momentum using data taken with the OPAL experiment at LEP. Both the fragmentation function and polarization measurements are implemented.

Source code: OPAL_1997_S3608263.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/Beam.hh"
  4#include "Rivet/Projections/FinalState.hh"
  5#include "Rivet/Projections/ChargedFinalState.hh"
  6#include "Rivet/Projections/UnstableParticles.hh"
  7#include "Rivet/Tools/BinnedHistogram.hh"
  8
  9namespace Rivet {
 10
 11
 12  /// @brief OPAL K*0 fragmentation function paper
 13  ///
 14  /// @author Peter Richardson
 15  class OPAL_1997_S3608263 : public Analysis {
 16  public:
 17
 18    RIVET_DEFAULT_ANALYSIS_CTOR(OPAL_1997_S3608263);
 19
 20
 21    /// @name Analysis methods
 22    /// @{
 23
 24    void init() {
 25      declare(Beam(), "Beams");
 26      declare(ChargedFinalState(), "FS");
 27      declare(UnstableParticles(), "UFS");
 28      book(_histXeK0   , 1, 1, 1);
 29      {Histo1DPtr temp; _h_ctheta.add(0.   ,0.01 ,book(temp, "ctheta_00",20,-1.,1.));}
 30      {Histo1DPtr temp; _h_ctheta.add(0.01 ,0.03 ,book(temp, "ctheta_01",20,-1.,1.));}
 31      {Histo1DPtr temp; _h_ctheta.add(0.03 ,0.10 ,book(temp, "ctheta_02",20,-1.,1.));}
 32      {Histo1DPtr temp; _h_ctheta.add(0.10 ,0.125,book(temp, "ctheta_03",20,-1.,1.));}
 33      {Histo1DPtr temp; _h_ctheta.add(0.125,0.14 ,book(temp, "ctheta_04",20,-1.,1.));}
 34      {Histo1DPtr temp; _h_ctheta.add(0.14 ,0.16 ,book(temp, "ctheta_05",20,-1.,1.));}
 35      {Histo1DPtr temp; _h_ctheta.add(0.16 ,0.20 ,book(temp, "ctheta_06",20,-1.,1.));}
 36      {Histo1DPtr temp; _h_ctheta.add(0.20 ,0.30 ,book(temp, "ctheta_07",20,-1.,1.));}
 37      {Histo1DPtr temp; _h_ctheta.add(0.30 ,0.40 ,book(temp, "ctheta_08",20,-1.,1.));}
 38      {Histo1DPtr temp; _h_ctheta.add(0.40 ,0.50 ,book(temp, "ctheta_09",20,-1.,1.));}
 39      {Histo1DPtr temp; _h_ctheta.add(0.50 ,0.70 ,book(temp, "ctheta_10",20,-1.,1.));}
 40      {Histo1DPtr temp; _h_ctheta.add(0.70 ,1.00 ,book(temp, "ctheta_11",20,-1.,1.));}
 41
 42      book(_h_ctheta_large,"ctheta_large",20,-1.,1.);
 43
 44      {Histo1DPtr temp; _h_alpha.add(0.   ,0.01 ,book(temp, "alpha_00",20,0.,0.5*M_PI));}
 45      {Histo1DPtr temp; _h_alpha.add(0.01 ,0.03 ,book(temp, "alpha_01",20,0.,0.5*M_PI));}
 46      {Histo1DPtr temp; _h_alpha.add(0.03 ,0.10 ,book(temp, "alpha_02",20,0.,0.5*M_PI));}
 47      {Histo1DPtr temp; _h_alpha.add(0.10 ,0.125,book(temp, "alpha_03",20,0.,0.5*M_PI));}
 48      {Histo1DPtr temp; _h_alpha.add(0.125,0.14 ,book(temp, "alpha_04",20,0.,0.5*M_PI));}
 49      {Histo1DPtr temp; _h_alpha.add(0.14 ,0.16 ,book(temp, "alpha_05",20,0.,0.5*M_PI));}
 50      {Histo1DPtr temp; _h_alpha.add(0.16 ,0.20 ,book(temp, "alpha_06",20,0.,0.5*M_PI));}
 51      {Histo1DPtr temp; _h_alpha.add(0.20 ,0.30 ,book(temp, "alpha_07",20,0.,0.5*M_PI));}
 52      {Histo1DPtr temp; _h_alpha.add(0.30 ,0.40 ,book(temp, "alpha_08",20,0.,0.5*M_PI));}
 53      {Histo1DPtr temp; _h_alpha.add(0.40 ,0.50 ,book(temp, "alpha_09",20,0.,0.5*M_PI));}
 54      {Histo1DPtr temp; _h_alpha.add(0.50 ,0.70 ,book(temp, "alpha_10",20,0.,0.5*M_PI));}
 55      {Histo1DPtr temp; _h_alpha.add(0.70 ,1.00 ,book(temp, "alpha_11",20,0.,0.5*M_PI));}
 56
 57      {Histo1DPtr temp; _h_alpha_low.add(0.30 ,0.40 ,book(temp, "alpha_low_00",20,0.,0.5*M_PI));}
 58      {Histo1DPtr temp; _h_alpha_low.add(0.40 ,0.50 ,book(temp, "alpha_low_01",20,0.,0.5*M_PI));}
 59      {Histo1DPtr temp; _h_alpha_low.add(0.50 ,0.70 ,book(temp, "alpha_low_02",20,0.,0.5*M_PI));}
 60      {Histo1DPtr temp; _h_alpha_low.add(0.70 ,1.00 ,book(temp, "alpha_low_03",20,0.,0.5*M_PI));}
 61
 62      {Histo1DPtr temp; _h_alpha_high.add(0.30 ,0.40 ,book(temp, "alpha_high_00",20,0.,0.5*M_PI));}
 63      {Histo1DPtr temp; _h_alpha_high.add(0.40 ,0.50 ,book(temp, "alpha_high_01",20,0.,0.5*M_PI));}
 64      {Histo1DPtr temp; _h_alpha_high.add(0.50 ,0.70 ,book(temp, "alpha_high_02",20,0.,0.5*M_PI));}
 65      {Histo1DPtr temp; _h_alpha_high.add(0.70 ,1.00 ,book(temp, "alpha_high_03",20,0.,0.5*M_PI));}
 66
 67
 68      book(_h_alpha_large     ,"_h_alpha_large"     ,20,0.,0.5*M_PI);
 69      book(_h_alpha_large_low ,"_h_alpha_large_low" ,20,0.,0.5*M_PI);
 70      book(_h_alpha_large_high,"_h_alpha_large_high",20,0.,0.5*M_PI);
 71    }
 72
 73    pair<double,double> calcRho(Histo1DPtr hist,unsigned int imode) {
 74      if(hist->numEntries()==0.) return make_pair(0.,0.);
 75      double sum1(0.),sum2(0.);
 76      for (auto bin : hist->bins() ) {
 77        double Oi = bin.area();
 78        if(Oi==0.) continue;
 79        double ai,bi;
 80        if(imode==0) {
 81          ai = 0.25*(bin.xMax()*(3.-sqr(bin.xMax())) - bin.xMin()*(3.-sqr(bin.xMin())));
 82          bi = 0.75*(bin.xMin()*(1.-sqr(bin.xMin())) - bin.xMax()*(1.-sqr(bin.xMax())));
 83        }
 84        else {
 85          ai = 2.*(bin.xMax()-bin.xMin())/M_PI;
 86          bi = 2.*(sin(2.*bin.xMax())-sin(2.*bin.xMin()))/M_PI;
 87        }
 88        double Ei = bin.areaErr();
 89        sum1 += sqr(bi/Ei);
 90        sum2 += bi/sqr(Ei)*(Oi-ai);
 91      }
 92      return make_pair(sum2/sum1,sqrt(1./sum1));
 93    }
 94
 95
 96    void analyze(const Event& e) {
 97      // First, veto on leptonic events by requiring at least 4 charged FS particles
 98      const FinalState& fs = apply<FinalState>(e, "FS");
 99      const size_t numParticles = fs.particles().size();
100
101      // Even if we only generate hadronic events, we still need a cut on numCharged >= 2.
102      if (numParticles < 2) {
103        MSG_DEBUG("Failed leptonic event cut");
104        vetoEvent;
105      }
106      MSG_DEBUG("Passed leptonic event cut");
107
108      // Get beams and average beam momentum
109      const ParticlePair& beams = apply<Beam>(e, "Beams").beams();
110      const double meanBeamMom = ( beams.first.p3().mod() +
111                                   beams.second.p3().mod() ) / 2.0;
112      MSG_DEBUG("Avg beam momentum = " << meanBeamMom);
113      Vector3 axis;
114      if(beams.first.pid()>0)
115        axis = beams.first .momentum().p3().unit();
116      else
117        axis = beams.second.momentum().p3().unit();
118
119      // Final state of unstable particles to get particle spectra
120      const UnstableParticles& ufs = apply<UnstableParticles>(e, "UFS");
121
122      for (const Particle& p : ufs.particles(Cuts::abspid==313)) {
123        double xp = p.p3().mod()/meanBeamMom;
124        _histXeK0->fill(xp);
125        int sign = p.pid()/313;
126        if(p.children().size()!=2) continue;
127        Particle kaon;
128        if(p.children()[0].pid()==sign*321 && p.children()[1].pid()==-sign*211) {
129          kaon = p.children()[0];
130        }
131        else if(p.children()[1].pid()==sign*321 && p.children()[0].pid()==-sign*211) {
132          kaon = p.children()[1];
133        }
134        else
135          continue;
136        // spin axes
137        Vector3 e1z = p.momentum().p3().unit();
138        Vector3 e1y = e1z.cross(axis).unit();
139        Vector3 e1x = e1y.cross(e1z).unit();
140        LorentzTransform boost = LorentzTransform::mkFrameTransformFromBeta(p.momentum().betaVec());
141        Vector3 axis1 = boost.transform(kaon.momentum()).p3().unit();
142        double ctheta = e1z.dot(axis1);
143        double phi = atan2(e1y.dot(axis1),e1x.dot(axis1));
144        double alpha = abs(abs(phi)-0.5*M_PI);
145        _h_ctheta.fill(xp,ctheta);
146        _h_alpha .fill(xp,alpha );
147        double cBeam = axis.dot(e1z);
148        (cBeam < 0.5 ? _h_alpha_low : _h_alpha_high).fill(xp,alpha);
149
150        if (xp > 0.3) {
151          _h_ctheta_large->fill(ctheta);
152          _h_alpha_large->fill(alpha);
153          (cBeam < 0.5 ? _h_alpha_large_low : _h_alpha_large_high)->fill(alpha);
154        }
155      }
156    }
157
158
159    /// Finalize
160    void finalize() {
161      scale(_histXeK0, 1./sumOfWeights());
162      const vector<double> x = {0., 0.01 ,0.03 ,0.10 ,0.125,0.14 ,0.16 ,0.20 ,0.30 ,0.40 ,0.50 ,0.70 ,1.00};
163      Scatter2DPtr h_rho00;
164      book(h_rho00,2,1,1);
165      Scatter2DPtr h_rho_off;
166      book(h_rho_off,2,1,2);
167      Scatter2DPtr h_ratio;
168      book(h_ratio,3,1,1);
169      Scatter2DPtr h_off_low;
170      book(h_off_low,4,1,1);
171      Scatter2DPtr h_off_high;
172      book(h_off_high,4,1,2);
173      for (unsigned int ix=0;ix<_h_ctheta.histos().size();++ix) {
174        // extract the rho00 component
175        normalize(_h_ctheta.histos()[ix]);
176        pair<double,double> rho00 = calcRho(_h_ctheta.histos()[ix],0);
177        h_rho00->addPoint(0.5*(x[ix]+x[ix+1]), rho00.first, make_pair(0.5*(x[ix+1]-x[ix]),0.5*(x[ix+1]-x[ix])),
178                          make_pair(rho00.second,rho00.second) );
179        // extract the rho 1 -1 component
180        normalize(_h_alpha.histos()[ix]);
181        pair<double,double> rho_off = calcRho(_h_alpha.histos()[ix],1);
182        h_rho_off->addPoint(0.5*(x[ix]+x[ix+1]), rho_off.first, make_pair(0.5*(x[ix+1]-x[ix]),0.5*(x[ix+1]-x[ix])),
183                            make_pair(rho_off.second,rho_off.second) );
184        if (ix<8) continue;
185        // at large xp also the ratio
186        double ratio = rho_off.first/(1.-rho00.first);
187        double dr    = ((rho_off.second - rho00.first*rho_off.second + rho_off.first*rho00.second))/sqr(1.-rho00.first);
188        h_ratio->addPoint(0.5*(x[ix]+x[ix+1]), ratio, make_pair(0.5*(x[ix+1]-x[ix]),0.5*(x[ix+1]-x[ix])),
189                          make_pair(dr,dr) );
190        unsigned int iy=ix-8;
191        // rho 1 -1 for cos theta <0.5
192        normalize(_h_alpha_low .histos()[iy]);
193        rho_off = calcRho(_h_alpha_low.histos()[iy],1);
194        h_off_low->addPoint(0.5*(x[ix]+x[ix+1]), rho_off.first, make_pair(0.5*(x[ix+1]-x[ix]),0.5*(x[ix+1]-x[ix])),
195                            make_pair(rho_off.second,rho_off.second) );
196        // rho 1 -1 for cos theta >0.5
197        normalize(_h_alpha_high.histos()[iy]);
198        rho_off = calcRho(_h_alpha_high.histos()[iy],1);
199        h_off_high->addPoint(0.5*(x[ix]+x[ix+1]), rho_off.first, make_pair(0.5*(x[ix+1]-x[ix]),0.5*(x[ix+1]-x[ix])),
200                             make_pair(rho_off.second,rho_off.second) );
201
202      }
203      // ratio for xp 0.3
204      normalize(_h_ctheta_large);
205      pair<double,double> rho00 = calcRho(_h_ctheta_large,0);
206      normalize(_h_alpha_large);
207      pair<double,double> rho_off = calcRho(_h_alpha_large,1);
208      double ratio = rho_off.first/(1.-rho00.first);
209      double dr    = ((rho_off.second - rho00.first*rho_off.second + rho_off.first*rho00.second))/sqr(1.-rho00.first);
210      Scatter2DPtr h_ratio_large;
211      book(h_ratio_large,3,2,1);
212      h_ratio_large->addPoint(0.65, ratio, make_pair(0.35,0.35),make_pair(dr,dr) );
213      // rho 1 -1 for xp >0.3 and cos theta < 0.5
214      normalize(_h_alpha_large_low );
215      rho_off = calcRho(_h_alpha_large_low,1);
216      book(h_off_low,4,2,1);
217      h_off_low->addPoint(0.65, rho_off.first, make_pair(0.35,0.35),make_pair(rho_off.second,rho_off.second) );
218      // rho 1 -1 for xp >0.3 and cos theta > 0.5
219      normalize(_h_alpha_large_high);
220      rho_off = calcRho(_h_alpha_large_high,1);
221      book(h_off_high,4,2,2);
222      h_off_high->addPoint(0.65, rho_off.first, make_pair(0.35,0.35),make_pair(rho_off.second,rho_off.second) );
223    }
224
225    /// @}
226
227
228  private:
229
230    /// @{
231    Histo1DPtr _histXeK0;
232    BinnedHistogram _h_ctheta,_h_alpha,_h_alpha_low,_h_alpha_high;
233    Histo1DPtr  _h_ctheta_large,_h_alpha_large, _h_alpha_large_low,_h_alpha_large_high;
234    /// @}
235
236  };
237
238
239
240  RIVET_DECLARE_ALIASED_PLUGIN(OPAL_1997_S3608263, OPAL_1997_I447146);
241
242}