rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2012_I1107658

Measurement of the underlying event activity in the Drell-Yan process at a centre-of-mass energy of 7 TeV
Experiment: CMS (LHC)
Inspire ID: 1107658
Status: VALIDATED
Authors:
  • Sunil Bansal
References:
  • CMS-QCD-11-012
  • CERN-PH-EP-2012-085
  • arXiv: 1204.1411
Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • Drell-Yan events with $Z/\gamma^* \to \mu\mu$. $m(\mu,\mu) > 20$ GeV

A measurement of the underlying event activity using Drell-Yan events using muonic final state. The production of charged particles with pseudorapidity $|\eta| < 2$ and transverse momentum $p_\perp > 0.5\,\text{GeV}/c$ is studied in towards, transverse and away region w.r.t. to the direction of di-muon system. The UE activity is measured in terms of of a particle density and an energy density. The particle density is computed as the average number of primary charged particles per unit pseudorapidity and per unit azimuth. The energy density is expressed in terms of the average of the scalar sum of the transverse momenta of primary charged particles per unit pseudorapidity and azimuth. The ratio of the energy and particle density is also reported in 3 regions. UE activity is studied as a function of invariant mass of muon pair ($M_{\mu\mu}$) by limiting the ISR contribution by requiring transverse momentum of muon pair $p_\perp(\mu\mu) < 5\,\text{GeV}/c$. The $p_\perp(\mu\mu)$ dependence is studied for the events having $M_{\mu\mu}$ in window of 81--101 GeV/$c$. The normalized charged particle multiplicity and $p_\perp$ spectrum of the charged particles in three regions also been reported for events having $M_{\mu\mu}$ in window of 81--101 GeV/$c$. Multiplicity and $p_\perp$ spectra in the transverse region are also reported, for events having $p_\perp(\mu\mu) < 5\,\text{GeV}/c$.

Source code: CMS_2012_I1107658.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/ZFinder.hh"
  5#include "Rivet/Projections/FinalState.hh"
  6#include "Rivet/Projections/ChargedFinalState.hh"
  7#include "Rivet/Projections/VetoedFinalState.hh"
  8
  9namespace Rivet {
 10
 11
 12  /// Underlying event activity in the Drell-Yan process at 7 TeV
 13  class CMS_2012_I1107658 : public Analysis {
 14  public:
 15
 16    /// Constructor
 17    CMS_2012_I1107658()
 18      : Analysis("CMS_2012_I1107658")
 19    {   }
 20
 21
 22    /// Initialization
 23    void init() {
 24
 25      /// @note Using a bare muon Z (but with a clustering radius!?)
 26      Cut cut = Cuts::abseta < 2.4 && Cuts::pT > 20*GeV;
 27      ZFinder zfinder(FinalState(), cut, PID::MUON, 4*GeV, 140*GeV, 0.2, ZFinder::ClusterPhotons::NONE);
 28      declare(zfinder, "ZFinder");
 29
 30      ChargedFinalState cfs((Cuts::etaIn(-2, 2) && Cuts::pT >=  500*MeV));
 31      VetoedFinalState nonmuons(cfs);
 32      nonmuons.addVetoPairId(PID::MUON);
 33      declare(nonmuons, "nonmuons");
 34
 35      book(_h_Nchg_towards_pTmumu                 ,1, 1, 1);
 36      book(_h_Nchg_transverse_pTmumu              ,2, 1, 1);
 37      book(_h_Nchg_away_pTmumu                    ,3, 1, 1);
 38      book(_h_pTsum_towards_pTmumu                ,4, 1, 1);
 39      book(_h_pTsum_transverse_pTmumu             ,5, 1, 1);
 40      book(_h_pTsum_away_pTmumu                   ,6, 1, 1);
 41      book(_h_avgpT_towards_pTmumu                ,7, 1, 1);
 42      book(_h_avgpT_transverse_pTmumu             ,8, 1, 1);
 43      book(_h_avgpT_away_pTmumu                   ,9, 1, 1);
 44      book(_h_Nchg_towards_plus_transverse_Mmumu  ,10, 1, 1);
 45      book(_h_pTsum_towards_plus_transverse_Mmumu ,11, 1, 1);
 46      book(_h_avgpT_towards_plus_transverse_Mmumu ,12, 1, 1);
 47      book(_h_Nchg_towards_zmass_81_101           ,13, 1, 1);
 48      book(_h_Nchg_transverse_zmass_81_101        ,14, 1, 1);
 49      book(_h_Nchg_away_zmass_81_101              ,15, 1, 1);
 50      book(_h_pT_towards_zmass_81_101             ,16, 1, 1);
 51      book(_h_pT_transverse_zmass_81_101          ,17, 1, 1);
 52      book(_h_pT_away_zmass_81_101                ,18, 1, 1);
 53      book(_h_Nchg_transverse_zpt_5               ,19, 1, 1);
 54      book(_h_pT_transverse_zpt_5                 ,20, 1, 1);
 55    }
 56
 57
 58    /// Perform the per-event analysis
 59    void analyze(const Event& event) {
 60      const double weight = 1.0;
 61      const ZFinder& zfinder = apply<ZFinder>(event, "ZFinder");
 62
 63      if (zfinder.bosons().size() != 1) vetoEvent;
 64
 65      double Zpt = zfinder.bosons()[0].pT()/GeV;
 66      double Zphi = zfinder.bosons()[0].phi();
 67      double Zmass = zfinder.bosons()[0].mass()/GeV;
 68
 69      Particles particles = apply<VetoedFinalState>(event, "nonmuons").particles();
 70
 71      int nTowards = 0;
 72      int nTransverse = 0;
 73      int nAway = 0;
 74      double ptSumTowards = 0;
 75      double ptSumTransverse = 0;
 76      double ptSumAway = 0;
 77
 78      for (const Particle& p : particles) {
 79        double dphi = fabs(deltaPhi(Zphi, p.phi()));
 80        double pT = p.pT();
 81
 82        if ( dphi < M_PI/3 ) {
 83          nTowards++;
 84          ptSumTowards += pT;
 85          if (Zmass > 81. && Zmass < 101.) _h_pT_towards_zmass_81_101->fill(pT, weight);
 86        } else if ( dphi < 2.*M_PI/3 ) {
 87          nTransverse++;
 88          ptSumTransverse += pT;
 89          if (Zmass > 81. && Zmass < 101.) _h_pT_transverse_zmass_81_101->fill(pT, weight);
 90          if (Zpt < 5.) _h_pT_transverse_zpt_5->fill(pT, weight);
 91        } else {
 92          nAway++;
 93          ptSumAway += pT;
 94          if (Zmass > 81. && Zmass < 101.) _h_pT_away_zmass_81_101->fill(pT, weight);
 95        }
 96
 97      } // Loop over particles
 98
 99
100      const double area = 8./3.*M_PI;
101      if (Zmass > 81. && Zmass < 101.) {
102        _h_Nchg_towards_pTmumu->         fill(Zpt, 1./area * nTowards, weight);
103        _h_Nchg_transverse_pTmumu->      fill(Zpt, 1./area * nTransverse, weight);
104        _h_Nchg_away_pTmumu->            fill(Zpt, 1./area * nAway, weight);
105        _h_pTsum_towards_pTmumu->        fill(Zpt, 1./area * ptSumTowards, weight);
106        _h_pTsum_transverse_pTmumu->     fill(Zpt, 1./area * ptSumTransverse, weight);
107        _h_pTsum_away_pTmumu->           fill(Zpt, 1./area * ptSumAway, weight);
108        if (nTowards > 0)    _h_avgpT_towards_pTmumu->    fill(Zpt, ptSumTowards/nTowards, weight);
109        if (nTransverse > 0) _h_avgpT_transverse_pTmumu-> fill(Zpt, ptSumTransverse/nTransverse, weight);
110        if (nAway > 0)       _h_avgpT_away_pTmumu->       fill(Zpt, ptSumAway/nAway, weight);
111        _h_Nchg_towards_zmass_81_101->   fill(nTowards, weight);
112        _h_Nchg_transverse_zmass_81_101->fill(nTransverse, weight);
113        _h_Nchg_away_zmass_81_101->      fill(nAway, weight);
114      }
115
116      if (Zpt < 5.) {
117        _h_Nchg_towards_plus_transverse_Mmumu->fill(Zmass, (nTowards + nTransverse)/(2.*area), weight);
118        _h_pTsum_towards_plus_transverse_Mmumu->fill(Zmass, (ptSumTowards + ptSumTransverse)/(2.*area), weight);
119        if ((nTowards + nTransverse) > 0) _h_avgpT_towards_plus_transverse_Mmumu->fill(Zmass, (ptSumTowards + ptSumTransverse)/(nTowards + nTransverse), weight);
120        _h_Nchg_transverse_zpt_5->fill(nTransverse, weight);
121      }
122
123    }
124
125
126    /// Normalise histograms etc., after the run
127    void finalize() {
128      scale(_h_pT_towards_zmass_81_101,    safediv(1, _h_Nchg_towards_zmass_81_101->integral(), 0));
129      scale(_h_pT_transverse_zmass_81_101, safediv(1, _h_Nchg_transverse_zmass_81_101->integral(), 0));
130      scale(_h_pT_away_zmass_81_101,       safediv(1, _h_Nchg_away_zmass_81_101->integral(), 0));
131      scale(_h_pT_transverse_zpt_5,        safediv(1, _h_Nchg_transverse_zpt_5->integral(), 0));
132      normalize(_h_Nchg_towards_zmass_81_101);
133      normalize(_h_Nchg_transverse_zmass_81_101);
134      normalize(_h_Nchg_away_zmass_81_101);
135      normalize(_h_Nchg_transverse_zpt_5);
136    }
137
138
139  private:
140
141
142    /// @name Histogram objects
143    //@{
144    Profile1DPtr _h_Nchg_towards_pTmumu;
145    Profile1DPtr _h_Nchg_transverse_pTmumu;
146    Profile1DPtr _h_Nchg_away_pTmumu;
147    Profile1DPtr _h_pTsum_towards_pTmumu;
148    Profile1DPtr _h_pTsum_transverse_pTmumu;
149    Profile1DPtr _h_pTsum_away_pTmumu;
150    Profile1DPtr _h_avgpT_towards_pTmumu;
151    Profile1DPtr _h_avgpT_transverse_pTmumu;
152    Profile1DPtr _h_avgpT_away_pTmumu;
153    Profile1DPtr _h_Nchg_towards_plus_transverse_Mmumu;
154    Profile1DPtr _h_pTsum_towards_plus_transverse_Mmumu;
155    Profile1DPtr _h_avgpT_towards_plus_transverse_Mmumu;
156    Histo1DPtr _h_Nchg_towards_zmass_81_101;
157    Histo1DPtr _h_Nchg_transverse_zmass_81_101;
158    Histo1DPtr _h_Nchg_away_zmass_81_101;
159    Histo1DPtr _h_pT_towards_zmass_81_101;
160    Histo1DPtr _h_pT_transverse_zmass_81_101;
161    Histo1DPtr _h_pT_away_zmass_81_101;
162    Histo1DPtr _h_Nchg_transverse_zpt_5;
163    Histo1DPtr _h_pT_transverse_zpt_5;
164    //@}
165
166  };
167
168
169  // Hook for the plugin system
170  RIVET_DECLARE_PLUGIN(CMS_2012_I1107658);
171
172}