rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

MC_LEADJETUE

Underlying event in leading jet events, extended to LHC
Experiment: ()
Status: VALIDATED
Authors:
  • Andy Buckley
No references listed
Beams: * *
Beam energies: ANY
Run details:
  • LHC pp QCD interactions at 0.9, 10 or 14 TeV. Particles with $c \tau > 10$ mm should be set stable. Several $p_\perp^\text{min}$ cutoffs are probably required to fill the profile histograms.

Rick Field's measurement of the underlying event in leading jet events, extended to the LHC. As usual, the leading jet of the defines an azimuthal toward/transverse/away decomposition, in this case the event is accepted within $|\eta| < 2$, as in the CDF 2008 version of the analysis. Since this isn't the Tevatron, I've chosen to use $k_\perp$ rather than midpoint jets.

Source code: MC_LEADJETUE.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/ChargedFinalState.hh"
  5#include "Rivet/Projections/FastJets.hh"
  6
  7namespace Rivet {
  8
  9
 10  /// @brief MC validation analysis for underlying event in jet events
 11  /// @author Andy Buckley
 12  class MC_LEADJETUE : public Analysis {
 13  public:
 14
 15    /// Constructor
 16    MC_LEADJETUE()
 17      : Analysis("MC_LEADJETUE")
 18    {    }
 19
 20
 21    /// @name Analysis methods
 22    //@{
 23
 24    // Book histograms
 25    void init() {
 26      // Final state for the jet finding
 27      const FinalState fsj((Cuts::etaIn(-4.0, 4.0)));
 28      declare(fsj, "FSJ");
 29      declare(FastJets(fsj, FastJets::KT, 0.7), "Jets");
 30
 31      // Charged final state for the distributions
 32      const ChargedFinalState cfs((Cuts::etaIn(-1.0, 1.0) && Cuts::pT >=  0.5*GeV));
 33      declare(cfs, "CFS");
 34
 35      const double maxpt1 = 500.0;
 36      book(_hist_pnchg      ,"trans-nchg", 50, 0.0, maxpt1);
 37      book(_hist_pmaxnchg   ,"trans-maxnchg", 50, 0.0, maxpt1);
 38      book(_hist_pminnchg   ,"trans-minnchg", 50, 0.0, maxpt1);
 39      book(_hist_pcptsum    ,"trans-ptsum", 50, 0.0, maxpt1);
 40      book(_hist_pmaxcptsum ,"trans-maxptsum", 50, 0.0, maxpt1);
 41      book(_hist_pmincptsum ,"trans-minptsum", 50, 0.0, maxpt1);
 42      book(_hist_pcptave    ,"trans-ptavg", 50, 0.0, maxpt1);
 43    }
 44
 45
 46    // Do the analysis
 47    void analyze(const Event& e) {
 48
 49      const FinalState& fsj = apply<FinalState>(e, "FSJ");
 50      if (fsj.particles().empty()) {
 51        MSG_DEBUG("Failed multiplicity cut");
 52        vetoEvent;
 53      }
 54
 55      const FastJets& jetpro = apply<FastJets>(e, "Jets");
 56      const Jets jets = jetpro.jetsByPt();
 57      MSG_DEBUG("Jet multiplicity = " << jets.size());
 58
 59      // Require the leading jet to be within |eta| < 2
 60      if (jets.size() < 1 || fabs(jets[0].eta()) > 2) {
 61        MSG_DEBUG("Failed jet cut");
 62        vetoEvent;
 63      }
 64
 65      const double jetphi = jets[0].phi();
 66      const double jetpT  = jets[0].pT();
 67      MSG_DEBUG("Leading jet: pT = " << jetpT/GeV << " GeV"
 68                << ", eta = " << jets[0].eta()
 69                << ", phi = " << jetphi);
 70
 71      // Get the final states to work with for filling the distributions
 72      const FinalState& cfs = apply<ChargedFinalState>(e, "CFS");
 73
 74      size_t numOverall(0), numToward(0), numTrans1(0), numTrans2(0), numAway(0);
 75      double ptSumTrans1(0.0), ptSumTrans2(0.0);
 76      double ptMaxOverall(0.0), ptMaxToward(0.0), ptMaxTrans1(0.0), ptMaxTrans2(0.0), ptMaxAway(0.0);
 77
 78      // Calculate all the charged stuff
 79      for (const Particle& p : cfs.particles()) {
 80        const double dPhi = deltaPhi(p.phi(), jetphi);
 81        const double pT = p.pT();
 82        const double phi = p.phi();
 83        const double rotatedphi = phi - jetphi;
 84
 85        //ptSumOverall += pT;
 86        ++numOverall;
 87        if (pT > ptMaxOverall) ptMaxOverall = pT;
 88
 89        if (dPhi < PI/3.0) {
 90          //ptSumToward += pT;
 91          ++numToward;
 92          if (pT > ptMaxToward) ptMaxToward = pT;
 93        }
 94        else if (dPhi < 2*PI/3.0) {
 95          if (rotatedphi <= PI) {
 96            ptSumTrans1 += pT;
 97            ++numTrans1;
 98            if (pT > ptMaxTrans1) ptMaxTrans1 = pT;
 99          } else {
100            ptSumTrans2 += pT;
101            ++numTrans2;
102            if (pT > ptMaxTrans2) ptMaxTrans2 = pT;
103          }
104        }
105        else {
106          //ptSumAway += pT;
107          ++numAway;
108          if (pT > ptMaxAway) ptMaxAway = pT;
109        }
110      }
111
112
113      // Fill the histograms
114      //_hist_tnchg->fill(jetpT/GeV, numToward/(4*PI/3));
115      _hist_pnchg->fill(jetpT/GeV, (numTrans1+numTrans2)/(4*PI/3));
116      _hist_pmaxnchg->fill(jetpT/GeV, (numTrans1>numTrans2 ? numTrans1 : numTrans2)/(2*PI/3));
117      _hist_pminnchg->fill(jetpT/GeV, (numTrans1<numTrans2 ? numTrans1 : numTrans2)/(2*PI/3));
118      //_hist_pdifnchg->fill(jetpT/GeV, abs(numTrans1-numTrans2)/(2*PI/3));
119      //_hist_anchg->fill(jetpT/GeV, numAway/(4*PI/3));
120
121      //_hist_tcptsum->fill(jetpT/GeV, ptSumToward/GeV/(4*PI/3));
122      _hist_pcptsum->fill(jetpT/GeV, (ptSumTrans1+ptSumTrans2)/GeV/(4*PI/3));
123      _hist_pmaxcptsum->fill(jetpT/GeV, (ptSumTrans1>ptSumTrans2 ? ptSumTrans1 : ptSumTrans2)/GeV/(2*PI/3));
124      _hist_pmincptsum->fill(jetpT/GeV, (ptSumTrans1<ptSumTrans2 ? ptSumTrans1 : ptSumTrans2)/GeV/(2*PI/3));
125      //_hist_pdifcptsum->fill(jetpT/GeV, fabs(ptSumTrans1-ptSumTrans2)/GeV/(2*PI/3));
126      //_hist_acptsum->fill(jetpT/GeV, ptSumAway/GeV/(4*PI/3));
127
128      //if (numToward > 0) {
129      //  _hist_tcptave->fill(jetpT/GeV, ptSumToward/GeV/numToward);
130      //  _hist_tcptmax->fill(jetpT/GeV, ptMaxToward/GeV);
131      //}
132      if ((numTrans1+numTrans2) > 0) {
133        _hist_pcptave->fill(jetpT/GeV, (ptSumTrans1+ptSumTrans2)/GeV/(numTrans1+numTrans2));
134        //_hist_pcptmax->fill(jetpT/GeV, (ptMaxTrans1 > ptMaxTrans2 ? ptMaxTrans1 : ptMaxTrans2)/GeV);
135      }
136      //if (numAway > 0) {
137      //  _hist_acptave->fill(jetpT/GeV, ptSumAway/GeV/numAway);
138      //  _hist_acptmax->fill(jetpT/GeV, ptMaxAway/GeV);
139      //}
140    }
141
142
143    void finalize() {
144      //
145    }
146
147
148  private:
149
150    Profile1DPtr _hist_pnchg;
151    Profile1DPtr _hist_pmaxnchg;
152    Profile1DPtr _hist_pminnchg;
153    Profile1DPtr _hist_pcptsum;
154    Profile1DPtr _hist_pmaxcptsum;
155    Profile1DPtr _hist_pmincptsum;
156    Profile1DPtr _hist_pcptave;
157
158  };
159
160
161
162  // The hook for the plugin system
163  RIVET_DECLARE_PLUGIN(MC_LEADJETUE);
164
165}