rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2016_I1444991

Higgs-to-WW differential cross sections at 8 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1444991
Status: VALIDATED
Authors:
  • Kathrin Becker
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • gg -> H -> W W* -> enu munu production at 8 TeV

This paper describes a measurement of fiducial and differential cross sections of gluon-fusion Higgs boson production in the $H\rightarrow W W^\ast \rightarrow e\nu\mu\nu$ channel, using 20.3 fb$^{-1}$ of proton-proton collision data. The data were produced at a centre-of-mass energy of $\sqrt{s} = 8$ TeV at the CERN Large Hadron Collider and recorded by the ATLAS detector in 2012. Cross sections are measured from the observed $H\rightarrow W W^\ast \rightarrow e\nu\mu\nu$ signal yield in categories distinguished by the number of associated jets. The total cross section is measured in a fiducial region defined by the kinematic properties of the charged leptons and neutrinos. Differential cross sections are reported as a function of the number of jets, the Higgs boson transverse momentum, the dilepton rapidity, and the transverse momentum of the leading jet. The jet-veto efficiency, or fraction of events with no jets above a given transverse momentum threshold, is also reported. All measurements are compared to QCD predictions from Monte Carlo generators and fixed-order calculations, and are in agreement with the Standard Model predictions.

Source code: ATLAS_2016_I1444991.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/DressedLeptons.hh"
  5#include "Rivet/Projections/IdentifiedFinalState.hh"
  6#include "Rivet/Projections/PromptFinalState.hh"
  7#include "Rivet/Projections/VetoedFinalState.hh"
  8#include "Rivet/Projections/FastJets.hh"
  9#include "Rivet/Projections/VisibleFinalState.hh"
 10
 11namespace Rivet {
 12
 13
 14  class ATLAS_2016_I1444991 : public Analysis {
 15  public:
 16
 17    /// Constructor
 18    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_I1444991);
 19
 20
 21  public:
 22
 23    /// Book histograms and initialise projections before the run
 24    void init() {
 25
 26      // All particles within |eta| < 5.0
 27      const FinalState FS(Cuts::abseta < 5.0);
 28
 29      // Project photons for dressing
 30      IdentifiedFinalState photon_id(FS);
 31      photon_id.acceptIdPair(PID::PHOTON);
 32
 33      // Project dressed electrons with pT > 15 GeV and |eta| < 2.47
 34      IdentifiedFinalState el_id(FS);
 35      el_id.acceptIdPair(PID::ELECTRON);
 36      PromptFinalState el_bare(el_id);
 37      Cut cuts = (Cuts::abseta < 2.47) && ( (Cuts::abseta <= 1.37) || (Cuts::abseta >= 1.52) ) && (Cuts::pT > 15*GeV);
 38      DressedLeptons el_dressed_FS(photon_id, el_bare, 0.1, cuts, true);
 39      declare(el_dressed_FS,"EL_DRESSED_FS");
 40
 41      // Project dressed muons with pT > 15 GeV and |eta| < 2.5
 42      IdentifiedFinalState mu_id(FS);
 43      mu_id.acceptIdPair(PID::MUON);
 44      PromptFinalState mu_bare(mu_id);
 45      DressedLeptons mu_dressed_FS(photon_id, mu_bare, 0.1, Cuts::abseta < 2.5 && Cuts::pT > 15*GeV, true);
 46      declare(mu_dressed_FS,"MU_DRESSED_FS");
 47
 48      // get MET from generic invisibles
 49      VetoedFinalState inv_fs(FS);
 50      inv_fs.addVetoOnThisFinalState(VisibleFinalState(FS));
 51      declare(inv_fs, "InvisibleFS");
 52
 53      // Project jets
 54      FastJets jets(FS, FastJets::ANTIKT, 0.4);
 55      jets.useInvisibles(JetAlg::Invisibles::NONE);
 56      jets.useMuons(JetAlg::Muons::NONE);
 57      declare(jets, "jets");
 58
 59      // Book histograms
 60      book(_h_Njets        , 2,1,1);
 61      book(_h_PtllMET      , 3,1,1);
 62      book(_h_Yll          , 4,1,1);
 63      book(_h_PtLead       , 5,1,1);
 64      book(_h_Njets_norm   , 6,1,1);
 65      book(_h_PtllMET_norm , 7,1,1);
 66      book(_h_Yll_norm     , 8,1,1);
 67      book(_h_PtLead_norm  , 9,1,1);
 68      book(_h_JetVeto      , 10, 1, 1, true);
 69
 70      //histos for jetveto
 71      std::vector<double> ptlead25_bins = { 0., 25., 300. };
 72      std::vector<double> ptlead40_bins = { 0., 40., 300. };
 73      book(_h_pTj1_sel25 , "pTj1_sel25", ptlead25_bins);
 74      book(_h_pTj1_sel40 , "pTj1_sel40", ptlead40_bins);
 75    }
 76
 77
 78    /// Perform the per-event analysis
 79    void analyze(const Event& event) {
 80
 81      // Get final state particles
 82      const FinalState& ifs = applyProjection<FinalState>(event, "InvisibleFS");
 83      const vector<DressedLepton>& good_mu = applyProjection<DressedLeptons>(event, "MU_DRESSED_FS").dressedLeptons();
 84      const vector<DressedLepton>& el_dressed = applyProjection<DressedLeptons>(event, "EL_DRESSED_FS").dressedLeptons();
 85      const Jets& jets = applyProjection<FastJets>(event, "jets").jetsByPt(Cuts::pT>25*GeV && Cuts::abseta < 4.5);
 86
 87      //find good electrons
 88      vector<DressedLepton> good_el;
 89      for (const DressedLepton& el : el_dressed){
 90        bool keep = true;
 91        for (const DressedLepton& mu : good_mu) {
 92          keep &= deltaR(el, mu) >= 0.1;
 93        }
 94        if (keep)  good_el += el;
 95      }
 96
 97      // select only emu events
 98      if ((good_el.size() != 1) || good_mu.size() != 1)  vetoEvent;
 99
100      //built dilepton
101      FourMomentum dilep = good_el[0].momentum() + good_mu[0].momentum();
102      double Mll = dilep.mass();
103      double Yll = dilep.rapidity();
104      double DPhill = fabs(deltaPhi(good_el[0], good_mu[0]));
105      double pTl1 = (good_el[0].pT() > good_mu[0].pT())? good_el[0].pT() : good_mu[0].pT();
106
107      //get MET
108      FourMomentum met;
109      for (const Particle& p : ifs.particles())  met += p.momentum();
110
111      // do a few cuts before looking at jets
112      if (pTl1 <= 22. || DPhill >= 1.8 || met.pT() <= 20.)  vetoEvent;
113      if (Mll <= 10. || Mll >= 55.)  vetoEvent;
114
115      Jets jets_selected;
116      for (const Jet &j : jets) {
117        if( j.abseta() > 2.4 && j.pT()<=30*GeV ) continue;
118        bool keep = true;
119        for(DressedLepton el : good_el) {
120          keep &= deltaR(j, el) >= 0.3;
121        }
122        if (keep)  jets_selected += j;
123      }
124
125      double PtllMET = (met + good_el[0].momentum() + good_mu[0].momentum()).pT();
126
127      double Njets = jets_selected.size() > 2 ? 2 : jets_selected.size();
128      double pTj1 = jets_selected.size()? jets_selected[0].pT() : 0.1;
129
130      // Fill histograms
131      _h_Njets->fill(Njets);
132      _h_PtllMET->fill(PtllMET);
133      _h_Yll->fill(fabs(Yll));
134      _h_PtLead->fill(pTj1);
135      _h_Njets_norm->fill(Njets);
136      _h_PtllMET_norm->fill(PtllMET);
137      _h_Yll_norm->fill(fabs(Yll));
138      _h_PtLead_norm->fill(pTj1);
139      _h_pTj1_sel25->fill(pTj1);
140      _h_pTj1_sel40->fill(pTj1);
141    }
142
143
144    /// Normalise histograms etc., after the run
145    void finalize() {
146
147      const double xs = crossSectionPerEvent()/femtobarn;
148
149      /// @todo Normalise, scale and otherwise manipulate histograms here
150      scale(_h_Njets, xs);
151      scale(_h_PtllMET, xs);
152      scale(_h_Yll, xs);
153      scale(_h_PtLead, xs);
154      normalize(_h_Njets_norm);
155      normalize(_h_PtllMET_norm);
156      normalize(_h_Yll_norm);
157      normalize(_h_PtLead_norm);
158      scale(_h_pTj1_sel25, xs);
159      scale(_h_pTj1_sel40, xs);
160      normalize(_h_pTj1_sel25);
161      normalize(_h_pTj1_sel40);
162      // fill jet veto efficiency histogram
163      _h_JetVeto->point(0).setY(_h_pTj1_sel25->bin(0).sumW(), sqrt(_h_pTj1_sel25->bin(0).sumW2()));
164      _h_JetVeto->point(1).setY(_h_PtLead_norm->bin(0).sumW(), sqrt(_h_PtLead_norm->bin(0).sumW2()));
165      _h_JetVeto->point(2).setY(_h_pTj1_sel40->bin(0).sumW(), sqrt(_h_pTj1_sel25->bin(0).sumW2()));
166
167      scale(_h_PtLead_norm , 1000.); // curveball unit change in HepData, just for this one
168      scale(_h_PtllMET_norm, 1000.); // curveball unit change in HepData, and this one
169    }
170
171  private:
172
173    /// @name Histograms
174    //@{
175    Histo1DPtr _h_Njets;
176    Histo1DPtr _h_PtllMET;
177    Histo1DPtr _h_Yll;
178    Histo1DPtr _h_PtLead;
179    Histo1DPtr _h_Njets_norm;
180    Histo1DPtr _h_PtllMET_norm;
181    Histo1DPtr _h_Yll_norm;
182    Histo1DPtr _h_PtLead_norm;
183
184    Scatter2DPtr _h_JetVeto;
185
186    Histo1DPtr _h_pTj1_sel25;
187    Histo1DPtr _h_pTj1_sel40;
188
189  };
190
191
192  // The hook for the plugin system
193  RIVET_DECLARE_PLUGIN(ATLAS_2016_I1444991);
194
195}