rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

PLUTO_1983_I191161

Measurement of average transverse momentum w.r.t thurst and Sphericity Axes for $E_{\text{CMS}}=7.7\to31.6$ GeV
Experiment: PLUTO (DORIS/Petra)
Inspire ID: 191161
Status: VALIDATED
Authors:
  • Peter Richardson
No references listed
Beams: e- e+
Beam energies: (3.9, 3.9); (4.7, 4.7); (6.0, 6.0); (6.5, 6.5); (8.5, 8.5); (11.0, 11.0); (13.8, 13.8); (15.4, 15.4) GeV
Run details:
  • $e^+e^-\to$ hadrons. Beam energy must be specified as analysis option "ENERGY" when rivet-merging samples.

Measurements of the average transverse momentum w.r.t thurst and Sphericity Axes for $E_{\text{CMS}}=7.7\to31.6$ GeV. Beam energy must be specified as analysis option "ENERGY" when rivet-merging samples.

Source code: PLUTO_1983_I191161.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/Beam.hh"
  4#include "Rivet/Projections/Thrust.hh"
  5#include "Rivet/Projections/Sphericity.hh"
  6#include "Rivet/Projections/FinalState.hh"
  7#include "Rivet/Projections/ChargedFinalState.hh"
  8
  9namespace Rivet {
 10
 11
 12  /// @brief average pT w.r.t. thrust and sphericity axes
 13  class PLUTO_1983_I191161 : public Analysis {
 14  public:
 15
 16    /// Constructor
 17    RIVET_DEFAULT_ANALYSIS_CTOR(PLUTO_1983_I191161);
 18
 19
 20    /// @name Analysis methods
 21    /// @{
 22
 23    /// Book histograms and initialise projections before the run
 24    void init() {
 25
 26      // Initialise and register projections
 27      const FinalState fs;
 28      declare(fs, "FS");
 29      declare(ChargedFinalState(), "CFS");
 30      // Thrust
 31      const Thrust thrust(fs);
 32      declare(thrust, "Thrust");
 33      const Sphericity sphericity(fs);
 34      declare(sphericity, "Sphericity");
 35
 36      sqs = 1.0;
 37      if (isCompatibleWithSqrtS(7.7)) sqs = 7.7;
 38      else if (isCompatibleWithSqrtS(9.4)) sqs = 9.4;
 39      else if (isCompatibleWithSqrtS(12.)) sqs = 12.;
 40      else if (isCompatibleWithSqrtS(13.)) sqs = 13.;
 41      else if (isCompatibleWithSqrtS(17.)) sqs = 17.;
 42      else if (isCompatibleWithSqrtS(22.)) sqs = 22.;
 43      else if (isCompatibleWithSqrtS(27.6)) sqs = 27.6;
 44      else if (isCompatibleWithSqrtS(30.8)) sqs = 30.8;
 45      else MSG_ERROR("Beam energy " << sqrtS() << " not supported!");
 46    }
 47
 48
 49    /// Perform the per-event analysis
 50    void analyze(const Event& event) {
 51      // at least 4 charged particles
 52      if (apply<ChargedFinalState>(event, "CFS").particles().size()<4) vetoEvent;
 53      // Sphericities
 54      MSG_DEBUG("Calculating sphericity");
 55      const Sphericity& sphericity = apply<Sphericity>(event, "Sphericity");
 56      MSG_DEBUG("Calculating thrust");
 57      const Thrust& thrust = apply<Thrust>(event, "Thrust");
 58      const FinalState & fs = apply<FinalState>(event, "FS");
 59      // remove pi0->gammagamma decay products and replace with the pi0s
 60      // needed to get the average pTs right
 61      Particles fsParticles;
 62      set<ConstGenParticlePtr> pi0;
 63      for(const Particle & p : fs.particles()) {
 64        if ((p.pid()!=PID::PHOTON && p.abspid()!=PID::ELECTRON)|| p.parents().empty() || p.parents()[0].pid()!=PID::PI0)
 65          fsParticles.push_back(p);
 66        else {
 67          if (pi0.find(p.parents()[0].genParticle())==pi0.end())
 68            fsParticles.push_back(p.parents()[0]);
 69        }
 70      }
 71      double pT_T_sum(0.),pT2_T_sum(0.);
 72      double pT_S_sum(0.),pT2_S_sum(0.);
 73      for(const Particle & p : fsParticles) {
 74        const Vector3 mom3 = p.p3();
 75        const double pTinT = dot(mom3, thrust.thrustMajorAxis());
 76        const double pToutT = dot(mom3, thrust.thrustMinorAxis());
 77        const double pTinS = dot(mom3, sphericity.sphericityMajorAxis());
 78        const double pToutS = dot(mom3, sphericity.sphericityMinorAxis());
 79        const double pT2_T = sqr(pTinT) + sqr(pToutT);
 80        const double pT2_S = sqr(pTinS) + sqr(pToutS);
 81        const double pT_T  = sqrt(pT2_T);
 82        const double pT_S  = sqrt(pT2_S);
 83        pT_T_sum  += sqrt(pT2_T);
 84        pT2_T_sum +=      pT2_T ;
 85        pT_S_sum  += sqrt(pT2_S);
 86        pT2_S_sum +=      pT2_S ;
 87        _p_thrust_pt      .fill(pT_T /MeV         );
 88        _p_thrust_pt2     .fill(pT2_T/1e3/sqr(MeV));
 89        _p_sphere_pt      .fill(pT_S /MeV         );
 90        _p_sphere_pt2     .fill(pT2_S/1e3/sqr(MeV));
 91      }
 92      _p_thrust_sum_pt  .fill(pT_T_sum /GeV               );
 93      _p_thrust_sum_pt2 .fill(pT2_T_sum/GeV               );
 94      _p_sphere_sum_pt  .fill(pT_S_sum /GeV               );
 95      _p_sphere_sum_pt2 .fill(pT2_S_sum/GeV               );
 96    }
 97
 98
 99    /// Normalise histograms etc., after the run
100    void finalize() {
101      for(unsigned int ix=1;ix<3;++ix) {
102        for(unsigned int iy=1;iy<5;++iy) {
103          double value = 0.0, error = 0.0;
104          if (ix==1) {
105            if (iy==1) {
106              value = _p_thrust_pt.xMean();
107              error = _p_thrust_pt.xStdErr();
108            }
109            else if (iy==2) {
110              value = _p_thrust_pt2.xMean();
111              error = _p_thrust_pt2.xStdErr();
112            }
113            else if (iy==3) {
114              value = _p_thrust_sum_pt.xMean();
115              error = _p_thrust_sum_pt.xStdErr();
116            }
117            else if (iy==4) {
118              value = _p_thrust_sum_pt2.xMean();
119              error = _p_thrust_sum_pt2.xStdErr();
120            }
121          }
122          else {
123            if (iy==1) {
124              value = _p_sphere_pt.xMean();
125              error = _p_sphere_pt.xStdErr();
126            }
127            else if (iy==2) {
128              value = _p_sphere_pt2.xMean();
129              error = _p_sphere_pt2.xStdErr();
130            }
131            else if (iy==3) {
132              value = _p_sphere_sum_pt.xMean();
133              error = _p_sphere_sum_pt.xStdErr();
134            }
135            else if (iy==4) {
136              value = _p_sphere_sum_pt2.xMean();
137              error = _p_sphere_sum_pt2.xStdErr();
138            }
139          }
140          Scatter2D temphisto(refData(ix, 1, iy));
141          Scatter2DPtr mult;
142          book(mult, ix, 1, iy);
143          for (size_t b = 0; b < temphisto.numPoints(); b++) {
144            const double x  = temphisto.point(b).x();
145            pair<double,double> ex = temphisto.point(b).xErrs();
146            pair<double,double> ex2 = ex;
147            if (ex2.first ==0.) ex2. first=0.0001;
148            if (ex2.second==0.) ex2.second=0.0001;
149            if (inRange(sqs, x-ex2.first, x+ex2.second)) {
150              mult->addPoint(x, value, ex, make_pair(error,error));
151            }
152            else {
153              mult->addPoint(x, 0., ex, make_pair(0.,.0));
154            }
155          }
156        }
157      }
158    }
159
160    /// @}
161
162
163    /// @name Histograms
164    //@{
165    YODA::Dbn1D _p_thrust_pt, _p_thrust_pt2, _p_thrust_sum_pt, _p_thrust_sum_pt2;
166    YODA::Dbn1D _p_sphere_pt, _p_sphere_pt2, _p_sphere_sum_pt, _p_sphere_sum_pt2;
167    double sqs;
168    //@}
169
170
171  };
172
173
174  // The hook for the plugin system
175  RIVET_DECLARE_PLUGIN(PLUTO_1983_I191161);
176
177
178}