1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/ChargedFinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5
  6namespace Rivet {
  7
  8
  9  /// @brief jet rates at 91 GeV
 10  class DELPHI_1990_I297698 : public Analysis {
 11  public:
 12
 13    /// Constructor
 14    DEFAULT_RIVET_ANALYSIS_CTOR(DELPHI_1990_I297698);
 15
 16
 17    /// @name Analysis methods
 18    ///@{
 19
 20    /// Book histograms and initialise projections before the run
 21    void init() {
 22      const ChargedFinalState cfs;
 23      declare(cfs, "FS");
 24      declare(FastJets(cfs, FastJets::JADE, 0.7), "JadeJets");
 25      // histos
 26      book(_h_2,1,1,1);
 27      book(_h_3,1,1,2);
 28      book(_h_4,1,1,3);
 29      book(_h_5,1,1,4);
 30    }
 31
 32
 33    /// Perform the per-event analysis
 34    void analyze(const Event& event) {
 35      const FinalState& fs = apply<FinalState>(event, "FS");
 36      const size_t numParticles = fs.particles().size();
 37      if (numParticles < 2) {
 38        MSG_DEBUG("Failed leptonic event cut");
 39        vetoEvent;
 40      }
 41      MSG_DEBUG("Passed leptonic event cut");
 42      const FastJets& jets = apply<FastJets>(event, "JadeJets");
 43      if (jets.clusterSeq()) {
 44        const double y_23 = jets.clusterSeq()->exclusive_ymerge_max(2);
 45        const double y_34 = jets.clusterSeq()->exclusive_ymerge_max(3);
 46        const double y_45 = jets.clusterSeq()->exclusive_ymerge_max(4);
 47        const double y_56 = jets.clusterSeq()->exclusive_ymerge_max(5);
 48        for (size_t i = 0; i < _h_2->numBins(); ++i) {
 49          double ycut = _h_2->bin(i).xMid();
 50          double width = _h_2->bin(i).xWidth();
 51          if (y_23 < ycut) {
 52            _h_2->fillBin(i, width);
 53          }
 54        }
 55        for (size_t i = 0; i < _h_3->numBins(); ++i) {
 56          double ycut = _h_3->bin(i).xMid();
 57          double width = _h_3->bin(i).xWidth();
 58          if (y_34 < ycut && y_23 > ycut) {
 59            _h_3->fillBin(i, width);
 60          }
 61        }
 62        for (size_t i = 0; i < _h_4->numBins(); ++i) {
 63          double ycut = _h_4->bin(i).xMid();
 64          double width = _h_4->bin(i).xWidth();
 65          if (y_45 < ycut && y_34 > ycut) {
 66            _h_4->fillBin(i, width);
 67          }
 68        }
 69        for (size_t i = 0; i < _h_5->numBins(); ++i) {
 70          double ycut = _h_5->bin(i).xMid();
 71          double width = _h_5->bin(i).xWidth();
 72          if (y_56 < ycut && y_45 > ycut) {
 73            _h_5->fillBin(i, width);
 74          }
 75        }
 76      }
 77    }
 78
 79
 80    /// Normalise histograms etc., after the run
 81    void finalize() {
 82      scale(_h_2, 100/sumOfWeights());
 83      scale(_h_3, 100/sumOfWeights());
 84      scale(_h_4, 100/sumOfWeights());
 85      scale(_h_5, 100/sumOfWeights());
 86    }
 87
 88    ///@}
 89
 90
 91    /// @name Histograms
 92    ///@{
 93    Histo1DPtr _h_2,_h_3,_h_4,_h_5;
 94    ///@}
 95
 96
 97  };
 98
 99
100  DECLARE_RIVET_PLUGIN(DELPHI_1990_I297698);
101
102}