rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2013_I1224539_DIJET

CMS jet mass measurement in dijet events
Experiment: CMS (LHC)
Inspire ID: 1224539
Status: VALIDATED
Authors:
  • salvatore.rappoccio@cern.ch
  • ntran@fnal.gov
  • kalanand@fnal.gov
  • dlopes@mail.cern.ch
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • 7 TeV $pp$ collisions. Events are required to have at least 2 jets with $(p_\perp_1 + p_\perp_2) / 2.0 > 220$ GeV with both jets satisfying $|y| < 2.4$.

Measurements of the mass spectra of the jets in dijet and $W/Z$+jet events from proton--proton collisions at a centre-of-mass energy of 7 TeV using a data sample of integrated luminosity of 5 fb$^-1$. The jets are reconstructed using the both the anti-$k_T$ algorithm with $R=0.7$ (AK7) and the Cambridge-Aachen algorithm with $R=0.8$ (CA8) and $R=1.2$ (CA12) with several grooming techniques applied (ungroomed, filtered, pruned and trimmed). See the text of the paper for more details. For the dijet events the distributions are presented as a function of the mean mass of the two leading jets in bins of the mean p_\perp of the two jets.

Source code: CMS_2013_I1224539_DIJET.cc
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/WFinder.hh"
#include "Rivet/Projections/ZFinder.hh"
#include "fastjet/tools/Filter.hh"
#include "fastjet/tools/Pruner.hh"

namespace Rivet {


  class CMS_2013_I1224539_DIJET : public Analysis {
  public:

    /// @name Constructors etc.
    //@{

    /// Constructor
    CMS_2013_I1224539_DIJET()
      : Analysis("CMS_2013_I1224539_DIJET"),
        _filter(fastjet::Filter(fastjet::JetDefinition(fastjet::cambridge_algorithm, 0.3), fastjet::SelectorNHardest(3))),
        _trimmer(fastjet::Filter(fastjet::JetDefinition(fastjet::kt_algorithm, 0.2), fastjet::SelectorPtFractionMin(0.03))),
        _pruner(fastjet::Pruner(fastjet::cambridge_algorithm, 0.1, 0.5))
    {    }

    //@}


  public:

    /// @name Analysis methods
    //@{

    /// Book histograms and initialise projections before the run
    void init() {
      FinalState fs(-2.4, 2.4, 0*GeV);
      declare(fs, "FS");

      // Jet collections
      declare(FastJets(fs, FastJets::ANTIKT, 0.7), "JetsAK7");
      declare(FastJets(fs, FastJets::CAM, 0.8), "JetsCA8");
      declare(FastJets(fs, FastJets::CAM, 1.2), "JetsCA12");

      // Histograms
      for (size_t i = 0; i < N_PT_BINS_dj; ++i ) {
        _h_ungroomedAvgJetMass_dj[i] = bookHisto1D(i+1+0*N_PT_BINS_dj, 1, 1);
        _h_filteredAvgJetMass_dj[i] = bookHisto1D(i+1+1*N_PT_BINS_dj, 1, 1);
        _h_trimmedAvgJetMass_dj[i] = bookHisto1D(i+1+2*N_PT_BINS_dj, 1, 1);
        _h_prunedAvgJetMass_dj[i] = bookHisto1D(i+1+3*N_PT_BINS_dj, 1, 1);
      }
    }


    // Find the pT histogram bin index for value pt (in GeV), to hack a 2D histogram equivalent
    /// @todo Use a YODA axis/finder alg when available
    size_t findPtBin(double ptJ) {
      const double ptBins_dj[N_PT_BINS_dj+1] = { 220.0, 300.0, 450.0, 500.0, 600.0, 800.0, 1000.0, 1500.0};
      for (size_t ibin = 0; ibin < N_PT_BINS_dj; ++ibin) {
        if (inRange(ptJ, ptBins_dj[ibin], ptBins_dj[ibin+1])) return ibin;
      }
      return N_PT_BINS_dj;
    }


    /// Perform the per-event analysis
    void analyze(const Event& event) {
      const double weight = event.weight();

      // Look at events with >= 2 jets
      const PseudoJets& psjetsAK7 = apply<FastJets>(event, "JetsAK7").pseudoJetsByPt( 50.0*GeV );
      if (psjetsAK7.size() < 2) vetoEvent;

      // Get the leading two jets and find their average pT
      const fastjet::PseudoJet& j0 = psjetsAK7[0];
      const fastjet::PseudoJet& j1 = psjetsAK7[1];
      double ptAvg = 0.5 * (j0.pt() + j1.pt());

      // Find the appropriate mean pT bin and escape if needed
      const size_t njetBin = findPtBin(ptAvg/GeV);
      if (njetBin >= N_PT_BINS_dj) vetoEvent;

      // Now run the substructure algs...
      fastjet::PseudoJet filtered0 = _filter(j0);
      fastjet::PseudoJet filtered1 = _filter(j1);
      fastjet::PseudoJet trimmed0 = _trimmer(j0);
      fastjet::PseudoJet trimmed1 = _trimmer(j1);
      fastjet::PseudoJet pruned0 = _pruner(j0);
      fastjet::PseudoJet pruned1 = _pruner(j1);

      // ... and fill the histograms
      _h_ungroomedAvgJetMass_dj[njetBin]->fill(0.5*(j0.m() + j1.m())/GeV, weight);
      _h_filteredAvgJetMass_dj[njetBin]->fill(0.5*(filtered0.m() + filtered1.m())/GeV, weight);
      _h_trimmedAvgJetMass_dj[njetBin]->fill(0.5*(trimmed0.m() + trimmed1.m())/GeV, weight);
      _h_prunedAvgJetMass_dj[njetBin]->fill(0.5*(pruned0.m() + pruned1.m())/GeV, weight);
    }


    /// Normalise histograms etc., after the run
    void finalize() {
      const double normalizationVal = 1000;
      for (size_t i = 0; i < N_PT_BINS_dj; ++i) {
        normalize(_h_ungroomedAvgJetMass_dj[i], normalizationVal);
        normalize(_h_filteredAvgJetMass_dj[i], normalizationVal);
        normalize(_h_trimmedAvgJetMass_dj[i], normalizationVal);
        normalize(_h_prunedAvgJetMass_dj[i], normalizationVal);
      }
    }

    //@}


  private:

    /// @name FastJet grooming tools (configured in constructor init list)
    //@{
    const fastjet::Filter _filter;
    const fastjet::Filter _trimmer;
    const fastjet::Pruner _pruner;
    //@}


    /// @name Histograms
    //@{
    enum BINS_dj { PT_220_300_dj=0, PT_300_450_dj, PT_450_500_dj, PT_500_600_dj,
                   PT_600_800_dj, PT_800_1000_dj, PT_1000_1500_dj, N_PT_BINS_dj };
    Histo1DPtr _h_ungroomedJet0pt, _h_ungroomedJet1pt;
    Histo1DPtr _h_ungroomedAvgJetMass_dj[N_PT_BINS_dj];
    Histo1DPtr _h_filteredAvgJetMass_dj[N_PT_BINS_dj];
    Histo1DPtr _h_trimmedAvgJetMass_dj[N_PT_BINS_dj];
    Histo1DPtr _h_prunedAvgJetMass_dj[N_PT_BINS_dj];
    //@}


  };



  // The hook for the plugin system
  DECLARE_RIVET_PLUGIN(CMS_2013_I1224539_DIJET);


}