1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/Beam.hh"
 4#include "Rivet/Projections/FinalState.hh"
 5#include "Rivet/Projections/ChargedFinalState.hh"
 6#include "Rivet/Projections/UnstableParticles.hh"
 7
 8namespace Rivet {
 9
10
11  /// @brief K*0, Phi and production
12  class DELPHI_1996_I420528 : public Analysis {
13  public:
14
15    /// Constructor
16    RIVET_DEFAULT_ANALYSIS_CTOR(DELPHI_1996_I420528);
17
18
19    /// @name Analysis methods
20    /// @{
21
22    /// Book histograms and initialise projections before the run
23    void init() {
24      declare(Beam(), "Beams");
25      declare(ChargedFinalState(), "FS");
26      declare(UnstableParticles(), "UFS");
27
28      // Book histograms
29      book(_h_Kstar,1, 1, 1);
30      book(_h_phi  ,3, 1, 1);
31
32    }
33
34
35    /// Perform the per-event analysis
36    void analyze(const Event& event) {
37      // First, veto on leptonic events by requiring at least 4 charged FS particles
38      const FinalState& fs = apply<FinalState>(event, "FS");
39      const size_t numParticles = fs.particles().size();
40
41      // Even if we only generate hadronic events, we still need a cut on numCharged >= 2.
42      if (numParticles < 2) {
43        MSG_DEBUG("Failed leptonic event cut");
44        vetoEvent;
45      }
46      MSG_DEBUG("Passed leptonic event cut");
47
48      // Get beams and average beam momentum
49      const ParticlePair& beams = apply<Beam>(event, "Beams").beams();
50      const double meanBeamMom = ( beams.first.p3().mod() +
51                                   beams.second.p3().mod() ) / 2.0;
52      MSG_DEBUG("Avg beam momentum = " << meanBeamMom);
53
54      // Final state of unstable particles to get particle spectra
55      const UnstableParticles& ufs = apply<UnstableParticles>(event, "UFS");
56
57      for (const Particle& p : ufs.particles(Cuts::abspid==313 or Cuts::pid==333)) {
58        const int id = p.abspid();
59        double xp = p.p3().mod()/meanBeamMom;
60        switch (id) {
61        case 333:
62          _h_phi->fill(xp);
63          break;
64        case 313:
65          _h_Kstar->fill(xp);
66          break;
67        }
68      }
69    }
70
71
72    /// Normalise histograms etc., after the run
73    void finalize() {
74      scale(_h_Kstar , 1./sumOfWeights());
75      scale(_h_phi , 1./sumOfWeights());
76    }
77
78    /// @}
79
80
81    /// @name Histograms
82    /// @{
83    Histo1DPtr _h_Kstar, _h_phi;
84    /// @}
85
86
87  };
88
89
90  RIVET_DECLARE_PLUGIN(DELPHI_1996_I420528);
91
92
93}