1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/FinalState.hh"
 4#include "Rivet/Projections/UnstableParticles.hh"
 5
 6namespace Rivet {
 7
 8
 9  /// @brief e+ e- > c cbar at 10.52 GeV
10  class CLEO_1988_I23509 : public Analysis {
11  public:
12
13    /// Constructor
14    RIVET_DEFAULT_ANALYSIS_CTOR(CLEO_1988_I23509);
15
16
17    /// @name Analysis methods
18    /// @{
19
20    /// Book histograms and initialise projections before the run
21    void init() {
22      // projections
23      declare(FinalState(), "FS");
24      declare(UnstableParticles(Cuts::abspid==411  ||
25                                Cuts::abspid==421  ||
26                                Cuts::abspid==431  ||
27                                Cuts::abspid==4122 ||
28                                Cuts::abspid==4232 ||
29                                Cuts::abspid==4132 ||
30                                Cuts::abspid==4332), "UFS");
31      // histos
32      for (unsigned int ix=0; ix<2; ++ix) {
33        book(_h[ix],2+ix,1,1);
34      }
35      book(_h[2],"TMP/h_mu",refData<YODA::BinnedEstimate<string>>(2,1,2));
36    }
37
38
39    /// Perform the per-event analysis
40    void analyze(const Event& event) {
41      const FinalState& fs = apply<FinalState>(event, "FS");
42
43      map<long,int> nCount;
44      int ntotal(0);
45      for (const Particle& p : fs.particles()) {
46        nCount[p.pid()] += 1;
47        ++ntotal;
48      }
49      // mu+mu- + photons
50      if (nCount[-13]==1 and nCount[13]==1 && ntotal==2+nCount[22]) {
51        _h[2]->fill("10.5"s);
52      }
53      else { // everything else
54        const FinalState& ufs = apply<FinalState>(event, "UFS");
55        if(!ufs.particles().empty()) {
56          _h[0]->fill("10.5"s);
57          _h[1]->fill("10.5"s);
58        }
59      }
60    }
61
62
63    /// Normalise histograms etc., after the run
64    void finalize() {
65      scale(_h, crossSection()/sumOfWeights()/nanobarn);
66      BinnedEstimatePtr<string> tmp;
67      book(tmp, 2, 1, 2);
68      divide(_h[0], _h[2],tmp);
69      book(tmp, 3, 1, 2);
70      divide(_h[1], _h[2],tmp);
71    }
72
73    /// @}
74
75
76    /// @name Histograms
77    /// @{
78    BinnedHistoPtr<string> _h[3];
79    /// @}
80
81
82  };
83
84
85  RIVET_DECLARE_PLUGIN(CLEO_1988_I23509);
86
87}