1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/Beam.hh"
 4#include "Rivet/Projections/FinalState.hh"
 5
 6namespace Rivet {
 7
 8
 9  /// @brief e+e- to Upsilon to mu+mu-
10  class CRYSTAL_BALL_1992_I306832 : public Analysis {
11  public:
12
13    /// Constructor
14    RIVET_DEFAULT_ANALYSIS_CTOR(CRYSTAL_BALL_1992_I306832);
15
16
17    /// @name Analysis methods
18    /// @{
19
20    /// Book histograms and initialise projections before the run
21    void init() {
22      // Initialise and register projections
23      declare(Beam(), "Beams");
24      declare(FinalState(), "FS");
25      // histograms
26      for(unsigned int ix=0;ix<2;++ix)
27        book(_c[ix], 1+ix, 1, 1);
28      _eCent = getOption<string>("ECENT", std::to_string(sqrtS()/MeV));
29    }
30
31
32    /// Perform the per-event analysis
33    void analyze(const Event& event) {
34      // loop over FS particles
35      Particles fs = apply<FinalState>(event,"FS").particles();
36      Particles mm,mp;
37      for (const Particle& p : fs) {
38        if (p.abspid()==PID::MUON) {
39          if (p.pid()>0) mm.push_back(p);
40          else           mp.push_back(p);
41        }
42        else if(p.pid()!=PID::GAMMA) {
43          vetoEvent;
44        }
45      }
46      if (mm.size()==1 && mp.size()==1) {
47        _c[0]->fill(_eCent);
48        _c[1]->fill(_eCent);
49      }
50    }
51
52
53    /// Normalise histograms etc., after the run
54    void finalize() {
55      double fact = crossSection()/ sumOfWeights() /picobarn;
56      for (unsigned int ix=0; ix<2; ++ix)
57        scale(_c[ix],fact);
58    }
59
60    /// @}
61
62
63    /// @name Histograms
64    /// @{
65    BinnedHistoPtr<string> _c[2];
66    string _eCent;
67    /// @}
68
69
70  };
71
72
73  RIVET_DECLARE_PLUGIN(CRYSTAL_BALL_1992_I306832);
74
75}