1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/UnstableParticles.hh"
 4
 5namespace Rivet {
 6
 7
 8  /// @brief pi0 production cross section measured at STAR
 9  class STAR_2014_I1253360 : public Analysis {
10  public:
11
12    /// Constructor
13    RIVET_DEFAULT_ANALYSIS_CTOR(STAR_2014_I1253360);
14
15    const double ETA_MIN = 0.8;
16    const double ETA_MAX = 2.0;
17    const double DELTA_ETA = ETA_MAX - ETA_MIN;
18
19    /// @name Analysis methods
20    ///@{
21
22    /// Book histograms and initialise projections before the run
23    void init() {
24
25      // Initialise and register projections
26
27      const UnstableParticles fs(Cuts::eta > ETA_MIN && Cuts::eta < ETA_MAX && Cuts::pid == 111);
28      declare(fs, "fs");
29
30      // Book histograms
31      book(_h_pi0_pt, 8, 1, 1);
32      _axis = YODA::Axis<double>({5., 6., 7., 8., 9., 10., 12., 16.});
33
34    }
35
36
37    /// Perform the per-event analysis
38    void analyze(const Event& event) {
39
40      if (_edges.empty())  _edges = _h_pi0_pt->xEdges();
41
42      const Particles& particles = apply<UnstableParticles>(event, "fs").particles();
43
44      for (const Particle& p : particles) {
45        discfill(p.pT() / GeV, 1 / (p.pT() / GeV));
46      }
47
48    }
49
50    void discfill(const double value, double weight) {
51      const size_t idx = _axis.index(value);
52      string edge("OTHER");
53      if (idx && idx <= _edges.size()) {
54        edge = _edges[idx-1];
55        weight /= _axis.width(idx);
56      }
57      _h_pi0_pt->fill(edge, weight);
58    }
59
60    /// Normalise histograms etc., after the run
61    void finalize() {
62
63      scale(_h_pi0_pt, crossSection()/millibarn/sumW() / (2 * pi) / DELTA_ETA);
64
65    }
66
67    ///@}
68
69
70    /// @name Histograms
71    ///@{
72    BinnedHistoPtr<string> _h_pi0_pt;
73    YODA::Axis<double> _axis;
74    vector<string> _edges;
75    ///@}
76
77
78  };
79
80
81  RIVET_DECLARE_PLUGIN(STAR_2014_I1253360);
82
83}