1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/Beam.hh"
 4#include "Rivet/Projections/ChargedFinalState.hh"
 5#include "Rivet/Projections/UnstableParticles.hh"
 6
 7namespace Rivet {
 8
 9
10  /// @brief Add a short analysis description here
11  class TOPAZ_1995_I381900 : public Analysis {
12  public:
13
14    /// Constructor
15    RIVET_DEFAULT_ANALYSIS_CTOR(TOPAZ_1995_I381900);
16
17
18    /// @name Analysis methods
19    //@{
20
21    /// Book histograms and initialise projections before the run
22    void init() {
23
24      // Initialise and register projections
25      declare(Beam(), "Beams");
26      declare(UnstableParticles(), "UFS");
27      const ChargedFinalState cfs;
28      declare(cfs, "CFS");
29
30      // Book histograms
31      book(_h_charged,1, 1, 1);
32      book(_h_pi     ,2, 1, 1);
33      book(_h_Kp     ,2, 1, 2);
34      book(_h_proton ,2, 1, 3);
35      book(_h_K0     ,3, 1, 1);
36      book(_wSum,"TMP/wSum");
37    }
38
39
40    /// Perform the per-event analysis
41    void analyze(const Event& event) {
42      const ChargedFinalState& cfs = apply<ChargedFinalState>(event, "CFS");
43      int nch = cfs.particles().size();
44      if(nch<5) vetoEvent;
45      // Get beams and average beam momentum
46      const ParticlePair& beams = apply<Beam>(event, "Beams").beams();
47      const double meanBeamMom = ( beams.first.p3().mod() +
48                                   beams.second.p3().mod() ) / 2.0;
49      MSG_DEBUG("Avg beam momentum = " << meanBeamMom);
50      _wSum->fill();
51
52      // neutral kaons
53      const UnstableParticles& ufs = apply<UnstableParticles>(event, "UFS");
54      for(const Particle & p : ufs.particles(Cuts::pid==130 || Cuts::pid==310)) {
55	double xi = -log(p.p3().mod()/meanBeamMom);
56	_h_K0->fill(xi);
57      }
58      // charged particles
59      for(const Particle & p : cfs.particles()) {
60	double xi = -log(p.p3().mod()/meanBeamMom);
61	_h_charged->fill(xi);
62	int id = abs(p.pid());
63	if(id==211)
64	  _h_pi->fill(xi);
65	else if(id==321)
66	  _h_Kp->fill(xi);
67	else if(id==2212)
68	  _h_proton->fill(xi);
69      }
70    }
71
72
73    /// Normalise histograms etc., after the run
74    void finalize() {
75      scale(_h_charged,1./ *_wSum);
76      scale(_h_pi     ,1./ *_wSum);
77      scale(_h_Kp     ,1./ *_wSum);
78      scale(_h_proton ,1./ *_wSum);
79      scale(_h_K0     ,1./ *_wSum);
80    }
81
82    //@}
83
84
85    /// @name Histograms
86    //@{
87    Histo1DPtr _h_charged,_h_pi,_h_Kp,_h_proton,_h_K0;
88    CounterPtr _wSum;
89    //@}
90
91
92  };
93
94
95  // The hook for the plugin system
96  RIVET_DECLARE_PLUGIN(TOPAZ_1995_I381900);
97
98
99}