1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/FinalState.hh"
 4
 5namespace Rivet {
 6
 7
 8  /// @brief  R measurement
 9  class TASSO_1980_I143690 : public Analysis {
10  public:
11
12    /// Constructor
13    RIVET_DEFAULT_ANALYSIS_CTOR(TASSO_1980_I143690);
14
15
16    /// @name Analysis methods
17    /// @{
18
19    /// Book histograms and initialise projections before the run
20    void init() {
21      // Initialise and register projections
22      declare(FinalState(), "FS");
23      // counters for R
24      book(_c_hadrons, "/TMP/sigma_hadrons", refData<YODA::BinnedEstimate<string>>(1,1,1));
25      book(_c_muons,   "/TMP/sigma_muons"  , refData<YODA::BinnedEstimate<string>>(1,1,1));
26      for (const string& en : _c_hadrons.binning().edges<0>()) {
27        const size_t idx = en.find("-");
28        if(idx!=std::string::npos) {
29          const double emin = std::stod(en.substr(0,idx));
30          const double emax = std::stod(en.substr(idx+1,string::npos));
31          if(inRange(sqrtS()/GeV, emin, emax)) {
32            _ecms = en;
33            break;
34          }
35        }
36        else {
37          const double end = std::stod(en)*GeV;
38          if (isCompatibleWithSqrtS(end)) {
39            _ecms = en;
40            break;
41          }
42        }
43      }
44      if (_ecms.empty())
45        MSG_ERROR("Beam energy incompatible with analysis.");
46    }
47
48
49    /// Perform the per-event analysis
50    void analyze(const Event& event) {
51      const FinalState& fs = apply<FinalState>(event, "FS");
52
53      map<long,int> nCount;
54      int ntotal(0);
55      for (const Particle& p : fs.particles()) {
56        nCount[p.pid()] += 1;
57        ++ntotal;
58      }
59      if (nCount[-13]==1 and nCount[13]==1 && ntotal==2+nCount[22]) {
60        // mu+mu- + photons
61        _c_muons->fill(_ecms);
62      }
63      else {
64        // everything else
65        _c_hadrons->fill(_ecms);
66      }
67    }
68
69
70    /// Normalise histograms etc., after the run
71    void finalize() {
72      BinnedEstimatePtr<string> mult;
73      book(mult, 1, 1, 1);
74      divide(_c_hadrons,_c_muons,mult);
75    }
76
77    /// @}
78
79
80    /// @name Histograms
81    /// @{
82    BinnedHistoPtr<string> _c_hadrons, _c_muons;
83    string _ecms;
84    /// @}
85
86
87  };
88
89
90  RIVET_DECLARE_PLUGIN(TASSO_1980_I143690);
91
92
93}