rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2014_I1288706

Measurement of the low-mass Drell-Yan differential cross section at 7 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1288706
Status: VALIDATED
Authors:
  • Elena Yatsenko
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • Run with inclusive $Z$ events, with $Z/\gamma^*$ decays to electrons and/or muons at $\sqrt{s} = 7$ TeV.

Measurements of the differential cross section for the process $Z/\gamma^\ast \rightarrow l$ $(l = e, \mu)$ as a function of dilepton invariant mass in pp collisions at $\sqrt{s} = 7$ TeV. The measurement is performed in the $e$ and $\mu$ channels for invariant masses between 26 GeV and 66 GeV in the fiducial region $p_\text{T}^\text{leading} > 15$ GeV, $p_\text{T}^\text{subleading} > 12$ GeV, $|\eta| < 2.4$ using an integrated luminosity of 1.6 $\text{fb}^{-1}$. The analysis is extended to invariant masses as low as 12 GeV in the muon channel within a fiducal region of $p_\text{T}^\text{leading} > 9$ GeV, $p_\text{T}^\text{subleading} > 6$ GeV, $|\eta| < 2.4$ with 35 $\text{pb}^{-1}$.

Source code: ATLAS_2014_I1288706.cc
 1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/FinalState.hh"
 4#include "Rivet/Projections/ZFinder.hh"
 5#include "Rivet/Particle.fhh"
 6
 7namespace Rivet {
 8
 9
10  class ATLAS_2014_I1288706 : public Analysis {
11  public:
12
13    /// Constructor
14    ATLAS_2014_I1288706()
15      : Analysis("ATLAS_2014_I1288706")
16    {
17     }
18
19
20    /// @name Analysis methods
21    //@{
22
23    /// Book histograms and initialise projections before the run
24    void init() {
25
26       // Set up projections
27       FinalState fs;
28
29       ZFinder zfinder_ext_dressed_mu(fs, Cuts::abseta<2.4 && Cuts::pT>6.0*GeV, PID::MUON, 12.0*GeV, 66.0*GeV, 0.1);
30       declare(zfinder_ext_dressed_mu, "ZFinder_ext_dressed_mu");
31
32       ZFinder zfinder_dressed_mu(fs, Cuts::abseta<2.4 && Cuts::pT>12*GeV, PID::MUON, 26.0*GeV, 66.0*GeV, 0.1);
33       declare(zfinder_dressed_mu, "ZFinder_dressed_mu");
34
35       ZFinder zfinder_dressed_el(fs, Cuts::abseta<2.4 && Cuts::pT>12*GeV, PID::ELECTRON, 26.0*GeV, 66.0*GeV, 0.1);
36       declare(zfinder_dressed_el, "ZFinder_dressed_el");
37
38       book(_hist_ext_mu_dressed ,1, 1, 1); // muon, dressed level, extended phase space
39       book(_hist_mu_dressed	    ,2, 1, 1); // muon, dressed level, nominal phase space
40       book(_hist_el_dressed	    ,2, 1, 2); // electron, dressed level, nominal phase space
41    }
42
43
44    /// Perform the per-event analysis
45    void analyze(const Event& event) {
46
47
48      const ZFinder& zfinder_ext_dressed_mu = apply<ZFinder>(event, "ZFinder_ext_dressed_mu");
49      const ZFinder& zfinder_dressed_mu     = apply<ZFinder>(event, "ZFinder_dressed_mu"    );
50      const ZFinder& zfinder_dressed_el     = apply<ZFinder>(event, "ZFinder_dressed_el"    );
51
52      fillPlots(zfinder_ext_dressed_mu, _hist_ext_mu_dressed, 9*GeV);
53      fillPlots(zfinder_dressed_mu,     _hist_mu_dressed,    15*GeV);
54      fillPlots(zfinder_dressed_el,     _hist_el_dressed,    15*GeV);
55    }
56
57
58    /// Helper function to fill the histogram if a Z is found with the required lepton cuts
59    void fillPlots(const ZFinder& zfinder, Histo1DPtr hist, double leading_pT) {
60      if (zfinder.bosons().size() != 1) return;
61      const FourMomentum el1 = zfinder.rawParticles()[0];
62      const FourMomentum el2 = zfinder.rawParticles()[1];
63      if (el1.pT() < leading_pT && el2.pT() < leading_pT) return;
64      const double mass = zfinder.bosons()[0].mass();
65      hist->fill(mass/GeV);
66    }
67
68
69    /// Normalise histograms etc., after the run
70    void finalize() {
71      scale(_hist_ext_mu_dressed, crossSection()/sumOfWeights());
72      scale(_hist_mu_dressed,     crossSection()/sumOfWeights());
73      scale(_hist_el_dressed,     crossSection()/sumOfWeights());
74    }
75
76    //@}
77
78
79  private:
80
81    /// Histograms
82    Histo1DPtr _hist_ext_mu_dressed, _hist_mu_dressed, _hist_el_dressed;
83
84  };
85
86
87  // The hook for the plugin system
88  RIVET_DECLARE_PLUGIN(ATLAS_2014_I1288706);
89
90}