Rivet analyses referenceATLAS_2014_I1288706Measurement of the low-mass Drell-Yan differential cross section at 7 TeVExperiment: ATLAS (LHC) Inspire ID: 1288706 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
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/DileptonFinder.hh"
5#include "Rivet/Particle.fhh"
6
7namespace Rivet {
8
9
10 /// Low-mass Drell-Yan differential cross section at 7 TeV
11 class ATLAS_2014_I1288706 : public Analysis {
12 public:
13
14 /// Constructor
15 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2014_I1288706);
16
17
18 /// @name Analysis methods
19 /// @{
20
21 /// Book histograms and initialise projections before the run
22 void init() {
23
24 // Set up projections
25 DileptonFinder zfinder_ext_dressed_mu(91.2*GeV, 0.1, Cuts::abseta < 2.4 && Cuts::pT > 6.0*GeV && Cuts::abspid == PID::MUON, Cuts::massIn(12.0*GeV, 66.0*GeV));
26 declare(zfinder_ext_dressed_mu, "DileptonFinder_ext_dressed_mu");
27
28 DileptonFinder zfinder_dressed_mu(91.2*GeV, 0.1, Cuts::abseta < 2.4 && Cuts::pT > 12*GeV && Cuts::abspid == PID::MUON, Cuts::massIn(26.0*GeV, 66.0*GeV));
29 declare(zfinder_dressed_mu, "DileptonFinder_dressed_mu");
30
31 DileptonFinder zfinder_dressed_el(91.2*GeV, 0.1, Cuts::abseta < 2.4 && Cuts::pT > 12*GeV && Cuts::abspid == PID::ELECTRON, Cuts::massIn(26.0*GeV, 66.0*GeV));
32 declare(zfinder_dressed_el, "DileptonFinder_dressed_el");
33
34 book(_hist_ext_mu_dressed, 1, 1, 1); // muon, dressed level, extended phase space
35 book(_hist_mu_dressed, 2, 1, 1); // muon, dressed level, nominal phase space
36 book(_hist_el_dressed, 2, 1, 2); // electron, dressed level, nominal phase space
37 }
38
39
40 /// Perform the per-event analysis
41 void analyze(const Event& event) {
42
43 const DileptonFinder& zfinder_ext_dressed_mu = apply<DileptonFinder>(event, "DileptonFinder_ext_dressed_mu");
44 const DileptonFinder& zfinder_dressed_mu = apply<DileptonFinder>(event, "DileptonFinder_dressed_mu" );
45 const DileptonFinder& zfinder_dressed_el = apply<DileptonFinder>(event, "DileptonFinder_dressed_el" );
46
47 fillPlots(zfinder_ext_dressed_mu, _hist_ext_mu_dressed, 9*GeV);
48 fillPlots(zfinder_dressed_mu, _hist_mu_dressed, 15*GeV);
49 fillPlots(zfinder_dressed_el, _hist_el_dressed, 15*GeV);
50 }
51
52
53 /// Helper function to fill the histogram if a Z is found with the required lepton cuts
54 void fillPlots(const DileptonFinder& zfinder, Histo1DPtr hist, double leading_pT) {
55 if (zfinder.bosons().size() != 1) return;
56 const FourMomentum el1 = zfinder.leptons()[0];
57 const FourMomentum el2 = zfinder.leptons()[1];
58 if (el1.pT() < leading_pT && el2.pT() < leading_pT) return;
59 const double mass = zfinder.bosons()[0].mass();
60 hist->fill(mass/GeV);
61 }
62
63
64 /// Normalise histograms etc., after the run
65 void finalize() {
66 scale(_hist_ext_mu_dressed, crossSection()/picobarn/sumOfWeights());
67 scale(_hist_mu_dressed, crossSection()/picobarn/sumOfWeights());
68 scale(_hist_el_dressed, crossSection()/picobarn/sumOfWeights());
69 }
70
71 /// @}
72
73
74 private:
75
76 /// Histograms
77 Histo1DPtr _hist_ext_mu_dressed, _hist_mu_dressed, _hist_el_dressed;
78
79 };
80
81
82 RIVET_DECLARE_PLUGIN(ATLAS_2014_I1288706);
83
84}
|