1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/ChargedFinalState.hh"
  4#include "Rivet/Projections/SmearedParticles.hh"
  5
  6namespace Rivet {
  7
  8
  9  /// Multiplicities and pT spectra from STAR for pp at 200 GeV
 10  class STAR_2008_I793126 : public Analysis {
 11  public:
 12
 13    /// Constructor
 14    RIVET_DEFAULT_ANALYSIS_CTOR(STAR_2008_I793126);
 15
 16
 17    /// @name Analysis methods
 18    /// @{
 19
 20    /// Book histograms and initialise projections before the run
 21    void init() {
 22      const ChargedFinalState cfs(Cuts::abseta < 0.5 && Cuts::pT >  0.2*GeV);
 23      const SmearedParticles lfs(cfs, [](const Particle& p) {
 24          // Track reconstruction efficiencies for tracks with pT from 0 to 600 MeV
 25          // in steps of 50 MeV. The efficiency is assumed to be 0.88 for pT >= 600 MeV
 26          const static vector<double> TRKEFF = {0,0,0.38,0.72,0.78,0.81,0.82,0.84,0.85,0.86,0.87,0.88};
 27          const size_t idx = size_t(min(floor(p.pT()/MeV/50), 11.));
 28          return TRKEFF[idx];
 29        });
 30      declare(lfs, "FS");
 31
 32      book(_h_dNch           ,1, 1, 1);
 33      book(_h_dpT_Pi         ,2, 1, 1);
 34      book(_h_dpT_Piplus     ,2, 1, 2);
 35      book(_h_dpT_Kaon       ,2, 1, 3);
 36      book(_h_dpT_Kaonplus   ,2, 1, 4);
 37      book(_h_dpT_AntiProton ,2, 1, 5);
 38      book(_h_dpT_Proton     ,2, 1, 6);
 39      // book(nCutsPassed, "nCutsPassed");
 40      // book(nPi, "nPi");
 41      // book(nPiPlus, "nPiPlus");
 42      // book(nKaon, "nKaon");
 43      // book(nKaonPlus, "nKaonPlus");
 44      // book(nProton, "nProton");
 45      // book(nAntiProton, "nAntiProton");
 46    }
 47
 48
 49    /// Perform the per-event analysis
 50    void analyze(const Event& event) {
 51      const ParticleFinder& charged = apply<ParticleFinder>(event, "FS");
 52
 53      // Vertex reconstruction efficiencies as a function of charged multiplicity.
 54      // For events with more than 23 reconstructed tracks the efficiency is 100%.
 55      double vtxeffs[24] = { 0.000000,0.512667,0.739365,0.847131,0.906946,0.940922,0.959328,0.96997,
 56                             0.975838,0.984432,0.988311,0.990327,0.990758,0.995767,0.99412,0.992271,
 57                             0.996631,0.994802,0.99635,0.997384,0.998986,0.996441,0.994513,1.000000 };
 58
 59      double vtxeff = 1.0;
 60      if (charged.particles().size() < 24) {
 61        vtxeff = vtxeffs[charged.particles().size()];
 62      }
 63
 64      const double weight = vtxeff;
 65
 66      for (const Particle& p : charged.particles()) {
 67        double pT = p.pT()/GeV;
 68        double y = p.rapidity();
 69        if (fabs(y) < 0.1) {
 70          // nCutsPassed->fill(weight);
 71          const PdgId id = p.pid();
 72          switch (id) {
 73          case -211:
 74            _h_dpT_Pi->fill(pT, weight/(TWOPI*pT*0.2));
 75            // nPi->fill(weight);
 76            break;
 77          case 211:
 78            _h_dpT_Piplus->fill(pT, weight/(TWOPI*pT*0.2));
 79            // nPiPlus->fill(weight);
 80            break;
 81          case -321:
 82            _h_dpT_Kaon->fill(pT, weight/(TWOPI*pT*0.2));
 83            // nKaon->fill(weight);
 84            break;
 85          case 321:
 86            _h_dpT_Kaonplus->fill(pT, weight/(TWOPI*pT*0.2));
 87            // nKaonPlus->fill(weight);
 88            break;
 89          case -2212:
 90            _h_dpT_AntiProton->fill(pT, weight/(TWOPI*pT*0.2));
 91            // nAntiProton->fill(weight);
 92            break;
 93          case 2212:
 94            _h_dpT_Proton->fill(pT, weight/(TWOPI*pT*0.2));
 95            // nProton->fill(weight);
 96            break;
 97          }
 98        }
 99        else {
100          continue;
101        }
102      }
103      _h_dNch->fill(charged.particles().size(),weight);
104    }
105
106
107    /// Normalise histograms etc., after the run
108    void finalize() {
109      //double nTot = nPi + nPiPlus + nKaon + nKaonPlus + nProton + nAntiProton;
110      normalize(_h_dNch);
111
112      /// @todo Norm to data!
113      normalize(_h_dpT_Pi        , 0.389825 );
114      normalize(_h_dpT_Piplus    , 0.396025 );
115      normalize(_h_dpT_Kaon      , 0.03897  );
116      normalize(_h_dpT_Kaonplus  , 0.04046  );
117      normalize(_h_dpT_AntiProton, 0.0187255);
118      normalize(_h_dpT_Proton    , 0.016511 );
119    }
120
121    /// @}
122
123
124    /// @name Histograms
125    /// @{
126    Histo1DPtr _h_dNch;
127    Histo1DPtr _h_dpT_Pi, _h_dpT_Piplus;
128    Histo1DPtr _h_dpT_Kaon, _h_dpT_Kaonplus;
129    Histo1DPtr _h_dpT_AntiProton, _h_dpT_Proton;
130    Profile1DPtr _h_pT_vs_Nch;
131    //CounterPtr nCutsPassed, nPi, nPiPlus, nKaon, nKaonPlus, nProton, nAntiProton;
132    ///@}
133
134  };
135
136
137
138  RIVET_DECLARE_ALIASED_PLUGIN(STAR_2008_I793126, STAR_2008_S7869363);
139
140}