rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2016_I1479760

Hard double-parton scattering in four-jet events at 7 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1479760
Status: VALIDATED
Authors:
  • Orel Gueta
  • Christian Gutschow
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • Four-jet events at 7 TeV

Inclusive four-jet events produced in $pp$ collisions at a centre-of-mass energy of $\sqrt{s} = 7$ TeV have been analysed for the presence of hard double-parton scattering using data corresponding to an integrated luminosity of $(37.3 \pm 1.3)$ $\mathrm{pb}^{-1}$, collected with the ATLAS detector at the LHC. The distributions of observables sensitive to the contribution of hard double-parton scattering are provided for events containing at least four jets with $p_{\mathrm{T}} \geq 20$ GeV and $\eta \leq 4.4$, and at least one of which has $p_{\text{T}} \geq 42.5$ GeV.

Source code: ATLAS_2016_I1479760.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5
  6namespace Rivet {
  7
  8
  9  /// Hard double-parton scattering in four-jet events at 7 TeV
 10  class ATLAS_2016_I1479760 : public Analysis {
 11  public:
 12
 13
 14    /// Constructor
 15    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_I1479760);
 16
 17
 18    /// Book histograms and initialise projections before the run
 19    void init() {
 20
 21      /// Declare AntiKt 0.6 jets without muons and neutrinos
 22      FastJets fastJets(FinalState(), JetAlg::ANTIKT, 0.6);
 23      fastJets.useInvisibles(JetInvisibles::NONE);
 24      fastJets.useMuons(JetMuons::NONE);
 25      declare(fastJets, "AntiKt6Jets");
 26
 27      book(_hists["deltaPt34"]       ,  1, 1, 1);
 28      book(_hists["deltaPhi34"]      ,  2, 1, 1);
 29      book(_hists["deltaPt12"]       ,  3, 1, 1);
 30      book(_hists["deltaPt13"]       ,  4, 1, 1);
 31      book(_hists["deltaPt23"]       ,  5, 1, 1);
 32      book(_hists["deltaPt14"]       ,  6, 1, 1);
 33      book(_hists["deltaPt24"]       ,  7, 1, 1);
 34      book(_hists["deltaPhi12"]      ,  8, 1, 1);
 35      book(_hists["deltaPhi13"]      ,  9, 1, 1);
 36      book(_hists["deltaPhi23"]      , 10, 1, 1);
 37      book(_hists["deltaPhi14"]      , 11, 1, 1);
 38      book(_hists["deltaPhi24"]      , 12, 1, 1);
 39      book(_hists["deltaY12"]        , 13, 1, 1);
 40      book(_hists["deltaY34"]        , 14, 1, 1);
 41      book(_hists["deltaY13"]        , 15, 1, 1);
 42      book(_hists["deltaY23"]        , 16, 1, 1);
 43      book(_hists["deltaY14"]        , 17, 1, 1);
 44      book(_hists["deltaY24"]        , 18, 1, 1);
 45      book(_hists["deltaPhiPlanes12"], 19, 1, 1);
 46      book(_hists["deltaPhiPlanes13"], 20, 1, 1);
 47      book(_hists["deltaPhiPlanes14"], 21, 1, 1);
 48    }
 49
 50
 51    /// Calculate the DeltaPt variable
 52    double calcDeltaPt(const Jet& j1, const Jet& j2) {
 53      return  (j1.momentum() + j2.momentum()).pT() / (j1.pT() + j2.pT());
 54    }
 55
 56    /// Calculate the DeltaPhi variable between event planes
 57    double calcDeltaPhiPlanes(const Jet& j1, const Jet& j2, const Jet& j3, const Jet& j4) {
 58      const FourMomentum sumVec1 = j1.momentum() + j2.momentum();
 59      const FourMomentum sumVec2 = j3.momentum() + j4.momentum();
 60      return  deltaPhi(sumVec1, sumVec2);
 61    }
 62
 63
 64    /// Perform the per-event analysis
 65    void analyze(const Event& event) {
 66
 67      // Retrieve all anti-kt R=0.6 jets with pT above 20 GeV and eta < 4.4
 68      const Jets jets = apply<JetFinder>(event, "AntiKt6Jets").jetsByPt(Cuts::pT >= 20*GeV && Cuts::abseta <= 4.4);
 69
 70      // Require at least 4 jets, with the leading jet pT above 42.5 GeV
 71      if (jets.size() < 4) vetoEvent;
 72      if (jets[0].pT() < 42.5*GeV) vetoEvent;
 73
 74      /// Fill histograms
 75      _hists["deltaPt12"]->fill( calcDeltaPt( jets[0], jets[1] ));
 76      _hists["deltaPt34"]->fill( calcDeltaPt( jets[2], jets[3] ));
 77      _hists["deltaPt13"]->fill( calcDeltaPt( jets[0], jets[2] ));
 78      _hists["deltaPt23"]->fill( calcDeltaPt( jets[1], jets[2] ));
 79      _hists["deltaPt14"]->fill( calcDeltaPt( jets[0], jets[3] ));
 80      _hists["deltaPt24"]->fill( calcDeltaPt( jets[1], jets[3] ));
 81      //
 82      _hists["deltaPhi12"]->fill( deltaPhi( jets[0],jets[1] ));
 83      _hists["deltaPhi34"]->fill( deltaPhi( jets[2],jets[3] ));
 84      _hists["deltaPhi13"]->fill( deltaPhi( jets[0],jets[2] ));
 85      _hists["deltaPhi23"]->fill( deltaPhi( jets[1],jets[2] ));
 86      _hists["deltaPhi14"]->fill( deltaPhi( jets[0],jets[3] ));
 87      _hists["deltaPhi24"]->fill( deltaPhi( jets[1],jets[3] ));
 88      //
 89      _hists["deltaY12"]->fill( deltaRap( jets[0], jets[1] ));
 90      _hists["deltaY34"]->fill( deltaRap( jets[2], jets[3] ));
 91      _hists["deltaY13"]->fill( deltaRap( jets[0], jets[2] ));
 92      _hists["deltaY23"]->fill( deltaRap( jets[1], jets[2] ));
 93      _hists["deltaY14"]->fill( deltaRap( jets[0], jets[3] ));
 94      _hists["deltaY24"]->fill( deltaRap( jets[1], jets[3] ));
 95      //
 96      _hists["deltaPhiPlanes12"]->fill( calcDeltaPhiPlanes(jets[0], jets[1], jets[2], jets[3] ));
 97      _hists["deltaPhiPlanes13"]->fill( calcDeltaPhiPlanes(jets[0], jets[2], jets[1], jets[3] ));
 98      _hists["deltaPhiPlanes14"]->fill( calcDeltaPhiPlanes(jets[0], jets[3], jets[1], jets[2] ));
 99    }
100
101
102    /// Post-run processing
103    void finalize() {
104      normalize(_hists);
105    }
106
107    /// Histograms
108    map<string, Histo1DPtr> _hists;
109
110  };
111
112
113  RIVET_DECLARE_PLUGIN(ATLAS_2016_I1479760);
114
115}