rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2011_S9120041

Traditional leading jet UE measurement at $\sqrt{s} = 0.9$ and 7 TeV
Experiment: CMS (LHC)
Inspire ID: 916908
Status: VALIDATED
Authors:
  • Mohammed Zakaria (mzakaria@ufl.edu)
References:
  • J. High Energy Phys 09 (2011) 109
Beams: p+ p+
Beam energies: (450.0, 450.0); (3500.0, 3500.0) GeV
Run details:
  • Requires inclusive inelastic events (non-diffractive and inelastic diffractive). The profile plots require large statistics.

A measurement of the underlying activity in scattering processes with a hard scale in the several-GeV region is performed in proton-proton collisions at Energies of 0.9 and 7\;TeV, using data collected by the CMS experiment at the LHC. The production of charged particles with pseudorapidity |eta| < 2 and transverse momentum $p_\perp > 0.5$\;GeV/$c$ is studied in the azimuthal region transverse to that of the leading set of charged particles forming a track-jet. Various comparisons are made between the two different energies and also beteen two sets of cuts on p_\perp for leading track jet p_\perp-leading $> 3$\;GeV and pT-leading $> 20$\;GeV. The activity is studied using 5 types of plots. Two profile plots for the multiplicity of charged particles and the scalar sum of p_\perp. and three distributions for the two previous quantities as well as p_\perp for all the particles in the transverse region.

Source code: CMS_2011_S9120041.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
144
145
// -*- C++ -*-

#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/ChargedFinalState.hh"
#include "Rivet/Projections/FastJets.hh"

using namespace std;

namespace Rivet {


  // UE charged particles vs. leading jet
  class CMS_2011_S9120041 : public Analysis {
  public:

    /// Constructor
    CMS_2011_S9120041() : Analysis("CMS_2011_S9120041") {}


    void init() {
      const ChargedFinalState cfs(-2.0, 2.0, 500*MeV);
      declare(cfs, "CFS");

      const ChargedFinalState cfsforjet(-2.5, 2.5, 500*MeV);
      const FastJets jetpro(cfsforjet, FastJets::SISCONE, 0.5);
      declare(jetpro, "Jets");

      if (fuzzyEquals(sqrtS(), 7.0*TeV)) {
        _h_Nch_vs_pT = bookProfile1D(1, 1, 1); // Nch vs. pT_max
        _h_Sum_vs_pT = bookProfile1D(2, 1, 1); // sum(pT) vs. pT_max
        _h_pT3_Nch   = bookHisto1D(5, 1, 1);   // transverse Nch,     pT_max > 3GeV
        _h_pT3_Sum   = bookHisto1D(6, 1, 1);   // transverse sum(pT), pT_max > 3GeV
        _h_pT3_pT    = bookHisto1D(7, 1, 1);   // transverse pT,      pT_max > 3GeV
        _h_pT20_Nch  = bookHisto1D(8, 1, 1);   // transverse Nch,     pT_max > 20GeV
        _h_pT20_Sum  = bookHisto1D(9, 1, 1);   // transverse sum(pT), pT_max > 20GeV
        _h_pT20_pT   = bookHisto1D(10, 1, 1);  // transverse pT,      pT_max > 20GeV
      }

      if (fuzzyEquals(sqrtS(), 0.9*TeV)) {
        _h_Nch_vs_pT = bookProfile1D(3, 1, 1); // Nch vs. pT_max
        _h_Sum_vs_pT = bookProfile1D(4, 1, 1); // sum(pT) vs. pT_max
        _h_pT3_Nch   = bookHisto1D(11, 1, 1);  // transverse Nch,     pT_max > 3GeV
        _h_pT3_Sum   = bookHisto1D(12, 1, 1);  // transverse sum(pT), pT_max > 3GeV
        _h_pT3_pT    = bookHisto1D(13, 1, 1);  // transverse pT,      pT_max > 3GeV
      }

      sumOfWeights3  = 0.0;
      sumOfWeights20 = 0.0;

      _nch_tot_pT3  = 0.0;
      _nch_tot_pT20 = 0.0;
    }


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

      // Find the lead jet, applying a restriction that the jets must be within |eta| < 2.
      FourMomentum p_lead;
      foreach (const Jet& j, apply<FastJets>(event, "Jets").jetsByPt(1.0*GeV)) {
        if (j.abseta() < 2.0) {
          p_lead = j.momentum();
          break;
        }
      }
      if (p_lead.isZero()) vetoEvent;
      const double philead = p_lead.phi();
      const double pTlead  = p_lead.pT();

      Particles particles = apply<ChargedFinalState>(event, "CFS").particlesByPt();

      int nTransverse = 0;
      double ptSumTransverse = 0.;
      foreach (const Particle& p, particles) {
        double dphi = fabs(deltaPhi(philead, p.phi()));
        if (dphi>PI/3. && dphi<PI*2./3.) {   // Transverse region
          nTransverse++;

          const double pT = p.pT()/GeV;
          ptSumTransverse += pT;

          if (pTlead > 3.0*GeV) _h_pT3_pT->fill(pT, weight);
          if (fuzzyEquals(sqrtS(), 7.0*TeV) && pTlead > 20.0*GeV) _h_pT20_pT->fill(pT, weight);
        }
      }

      const double area = 8./3. * PI;
      _h_Nch_vs_pT->fill(pTlead/GeV, 1./area*nTransverse, weight);
      _h_Sum_vs_pT->fill(pTlead/GeV, 1./area*ptSumTransverse, weight);
      if(pTlead > 3.0*GeV) {
        _h_pT3_Nch->fill(nTransverse, weight);
        _h_pT3_Sum->fill(ptSumTransverse, weight);
        sumOfWeights3 += weight;
        _nch_tot_pT3  += weight*nTransverse;
      }
      if (fuzzyEquals(sqrtS(), 7.0*TeV) && pTlead > 20.0*GeV) {
        _h_pT20_Nch->fill(nTransverse, weight);
        _h_pT20_Sum->fill(ptSumTransverse, weight);
        sumOfWeights20 += weight;
        _nch_tot_pT20  += weight*nTransverse;
      }
    }



    /// Normalise histograms etc., after the run
    void finalize() {
      normalize(_h_pT3_Nch);
      normalize(_h_pT3_Sum);
      if (sumOfWeights3 != 0.0) normalize(_h_pT3_pT, _nch_tot_pT3 / sumOfWeights3);

      if (fuzzyEquals(sqrtS(), 7.0*TeV)) {
        normalize(_h_pT20_Nch);
        normalize(_h_pT20_Sum);
        if (sumOfWeights20 != 0.0) normalize(_h_pT20_pT, _nch_tot_pT20 / sumOfWeights20);
      }
    }



  private:

    double sumOfWeights3;
    double sumOfWeights20;

    double _nch_tot_pT3;
    double _nch_tot_pT20;

    Profile1DPtr _h_Nch_vs_pT;
    Profile1DPtr _h_Sum_vs_pT;
    Histo1DPtr _h_pT3_Nch;
    Histo1DPtr _h_pT3_Sum;
    Histo1DPtr _h_pT3_pT;
    Histo1DPtr _h_pT20_Nch;
    Histo1DPtr _h_pT20_Sum;
    Histo1DPtr _h_pT20_pT;

  };


  // This global object acts as a hook for the plugin system
  DECLARE_RIVET_PLUGIN(CMS_2011_S9120041);
}