1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/FinalState.hh"
 4#include "Rivet/Projections/ChargedFinalState.hh"
 5#include "Rivet/Projections/UnstableParticles.hh"
 6#include "Rivet/Projections/Beam.hh"
 7
 8namespace Rivet {
 9
10
11  /// @brief f'2 spectrum
12  class DELPHI_1996_I416741 : public Analysis {
13  public:
14
15    /// Constructor
16    RIVET_DEFAULT_ANALYSIS_CTOR(DELPHI_1996_I416741);
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_f2prime, 1, 1, 1);
30
31    }
32
33
34    /// Perform the per-event analysis
35    void analyze(const Event& event) {
36
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==335)) {
58        double xp = p.p3().mod()/meanBeamMom;
59	_h_f2prime->fill(xp);
60      }
61    }
62
63
64    /// Normalise histograms etc., after the run
65    void finalize() {
66      scale(_h_f2prime, 1./sumOfWeights());
67    }
68
69    /// @}
70
71
72    /// @name Histograms
73    /// @{
74    Histo1DPtr _h_f2prime;
75    /// @}
76
77
78  };
79
80
81  RIVET_DECLARE_PLUGIN(DELPHI_1996_I416741);
82
83
84}