Rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2014_I1307243

Measurements of jet vetoes and azimuthal decorrelations in dijet events produced in pp collisions at $\sqrt{s}$ = 7 TeV using the ATLAS detector
Experiment: ATLAS (LHC)
Inspire ID: 1307243
Status: VALIDATED
Authors:
  • James Robinson
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • pp -> jet jet x. sqrt(s) = 7000 GeV. Jets with pT > 60, 50 GeV

Additional jet activity in dijet events is measured using $pp$ collisions at ATLAS at a centre-of-mass energy of 7 TeV, for jets reconstructed using the anti-kt algorithm with radius parameter R=0.6. This is done using variables such as the fraction of dijet events without an additional jet in the rapidity interval bounded by the dijet subsystem and correlations between the azimuthal angles of the dijets. They are presented, both with and without a veto on additional jet activity in the rapidity interval, as a function of the mean transverse momentum of the dijets and of the rapidity interval size. The double differential dijet cross section is also measured as a function of the interval size and the azimuthal angle between the dijets. These variables probe differences in the approach to resummation of large logarithms when performing QCD calculations. The data are compared to POWHEG, interfaced to the PYTHIA 8 and HERWIG parton shower generators, as well as to HEJ with and without interfacing it to the ARIADNE parton shower generator. None of the theoretical predictions agree with the data across the full phase-space considered; however, POWHEG+PYTHIA 8 and HEJ+ARIADNE are found to provide the best agreement with the data. These measurements use the full data sample collected with the ATLAS detector in 7 TeV $pp$ collisions at the LHC and correspond to integrated luminosities of 36.1 pb$^-1$ and 4.5 fb$^-1$ for data collected during 2010 and 2011 respectively.

Source code: ATLAS_2014_I1307243.cc
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Tools/Logging.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Tools/BinnedHistogram.hh"

namespace Rivet {


  /// @brief ATLAS azimuthal decorrelation with jet veto analysis
  /// @author James Robinson <james.robinson@cern.ch>
  class ATLAS_2014_I1307243 : public Analysis {
  public:

    /// Constructor
    ATLAS_2014_I1307243()
      : Analysis("ATLAS_2014_I1307243"),
        _fiducialRegions{2010, 2011},
        _vetoScale{20*GeV, 30*GeV},
        _yFiducial{4.4, 2.4},
        _gapCategories{"inclusive", "gap"},
        _dy_max(8),
        _nEventsInAcceptance(0),
        _sumOfAcceptedWeights(0.)
    {   }


    /// Book histograms and initialise projections before the run
    void init() {

      /// Initialise and register projections here
      FinalState fs;
      FastJets fastJets(fs, FastJets::ANTIKT, 0.6);
      fastJets.useInvisibles(true);
      declare(fastJets, "AntiKt6JetsWithInvisibles");


      /// Book histograms
      foreach (const string& gapCategory, _gapCategories ) {
        const int gapCategoryOffset = (gapCategory == "inclusive") ? 0 : 1;

        // Temporary inclusive and gap histograms
        _h_tmp_events_dy[gapCategory] = bookHisto1D(1, 1, 1);
        _h_tmp_events_dy[gapCategory]->setPath("/TMP/" + toString(gapCategory) + "_events_dy");
        _h_tmp_events_pTbar[gapCategory] = bookHisto1D(2, 1, 1);
        _h_tmp_events_pTbar[gapCategory]->setPath("/TMP/" + toString(gapCategory) + "_events_pTbar");

        // Azimuthal moment histograms
        _h_profiled_cosDeltaPhi_dy[gapCategory]       = bookProfile1D( 5+4*gapCategoryOffset, 1, 1);
        _h_profiled_cosDeltaPhi_pTbar[gapCategory]    = bookProfile1D( 6+4*gapCategoryOffset, 1, 1);
        _h_C2C1_dy[gapCategory]                       = bookScatter2D( 7+4*gapCategoryOffset, 1, 1, false);
        _h_C2C1_pTbar[gapCategory]                    = bookScatter2D( 8+4*gapCategoryOffset, 1, 1, false);
        _h_profiled_cosTwoDeltaPhi_dy[gapCategory]    = bookProfile1D(37+2*gapCategoryOffset, 1, 1);
        _h_profiled_cosTwoDeltaPhi_pTbar[gapCategory] = bookProfile1D(38+2*gapCategoryOffset, 1, 1);

        // Gap fraction vs. Q0 and cross-section in dy slices
        for (size_t dyLow = 0; dyLow < _dy_max; ++dyLow ) {
          Histo1DPtr _h_tmp_events_Q0_single_dySlice = bookHisto1D( 29+dyLow, 1, 1);
          _h_tmp_events_Q0_single_dySlice->setPath("/TMP/" + toString(gapCategory) + "_events_dySlice_" + toString(dyLow) + "_" + toString(dyLow+1) + "_Q0");
          _h_tmp_events_Q0_dySlices[gapCategory].addHistogram( dyLow, dyLow+1, _h_tmp_events_Q0_single_dySlice );
          _h_crossSection_dphi_dySlices[gapCategory].addHistogram( dyLow, dyLow+1, bookHisto1D( 13+(_dy_max*gapCategoryOffset)+dyLow, 1, 1));
        }

      }

      // Number of jets in rapidity interval
      _h_profiled_nJets_rapidity_interval_dy    = bookProfile1D( 3, 1, 1);
      _h_profiled_nJets_rapidity_interval_pTbar = bookProfile1D( 4, 1, 1);
    }


    /// Perform the per-event analysis
    void analyze(const Event& event) {

      // Get the event weight
      const double weight( event.weight() );
      bool eventAccepted( false );

      for (int iFiducialRegion = 0; iFiducialRegion < 2; ++iFiducialRegion ) {

        // Retrieve all anti-kt R=0.6 jets above _pTMin and inside |_yFiducial|
        const Jets akt6Jets = apply<JetAlg>(event, "AntiKt6JetsWithInvisibles").jetsByPt( Cuts::absrap < _yFiducial.at(iFiducialRegion) );
        // If there are fewer than 2 jets then bail
        if ( akt6Jets.size() < 2 ) { vetoEvent; }

        // Require jets to be above {60, 50} GeV
        if ( akt6Jets.at(0).momentum().pT() < 60.*GeV || akt6Jets.at(1).momentum().pT() < 50.*GeV ) { vetoEvent; }

        // Identify gap boundaries
        double yMin( std::min( akt6Jets.at(0).momentum().rapidity(), akt6Jets.at(1).momentum().rapidity() ) );
        double yMax( std::max( akt6Jets.at(0).momentum().rapidity(), akt6Jets.at(1).momentum().rapidity() ) );

        // Determine azimuthal decorrelation quantities
        const double dy( yMax - yMin );
        const double dphi( mapAngle0ToPi( akt6Jets.at(0).momentum().phi() - akt6Jets.at(1).momentum().phi() ) );
        const double pTbar( (akt6Jets.at(0).momentum().pT() + akt6Jets.at(1).momentum().pT())/2.0 );

        // Impose minimum dy for the 2011 phase space
        if ( _fiducialRegions.at(iFiducialRegion) == 2011 && dy < 1.0 ) { vetoEvent; }

        // Construct gap candidates sample
        Jets gapCandidates;
        foreach( const Jet &j, akt6Jets ) {
          if ( inRange( j.momentum().rapidity(), yMin, yMax, OPEN, OPEN ) ) {
            gapCandidates.push_back( j );
          }
        }

        // Determine gap quantities
        unsigned int nJets_rapidity_interval( 0 );
        double maximumGapQ0( 0. );
        foreach( const Jet &jet, gapCandidates ) {
          const double pT( jet.momentum().pT() );
          if ( pT > _vetoScale.at(iFiducialRegion) ) { ++nJets_rapidity_interval; }
          if ( pT > maximumGapQ0 ) { maximumGapQ0 = pT; }
        }

        // Fill histograms
        if ( weight < 0.0 ) {
          MSG_DEBUG( "Negative weight " << weight << "found!" );
        }
        fillHistograms( _fiducialRegions.at(iFiducialRegion), dy, pTbar, dphi, nJets_rapidity_interval, maximumGapQ0, weight );
        eventAccepted = true;
      }

      // Count number of accepted events
      if ( eventAccepted ) {
        _nEventsInAcceptance++;
        _sumOfAcceptedWeights += weight;
      }
      return;
    }

    void fillHistograms( const unsigned int &fiducialRegion, const double &dy, const double &pTbar, const double &dphi, const double &nJets_rapidity_interval, const double &maximumGapQ0, const double &weight) {
      // Determine gap category
      vector<string> eventGapCategories{{"inclusive"}};
      if ( nJets_rapidity_interval == 0 ) { eventGapCategories += string("gap"); }

      // Fill histograms relevant for comparison with 2010 data
      if ( fiducialRegion == _fiducialRegions.at(0) ) {
        // Fill inclusive and gap histograms
        foreach( const string &gapCategory, eventGapCategories ) {
          _h_tmp_events_dy[gapCategory]->fill( dy, weight );
          _h_crossSection_dphi_dySlices[gapCategory].fill( dy, dphi / M_PI, weight );
          _h_profiled_cosDeltaPhi_dy[gapCategory]->fill( dy, cos(M_PI - dphi), weight );
          _h_profiled_cosTwoDeltaPhi_dy[gapCategory]->fill( dy, cos(2.0 * dphi), weight );
        }
        // Fill profiled nJets_rapidity_interval
        _h_profiled_nJets_rapidity_interval_dy->fill( dy, nJets_rapidity_interval, weight );
        // Fill Q0 histograms - can fill multiple points per event
        foreach( const YODA::HistoBin1D Q0_bin, _h_tmp_events_Q0_dySlices["inclusive"].getHistograms().at(0)->bins() ) {
          const double Q0( Q0_bin.xMid() );
          _h_tmp_events_Q0_dySlices["inclusive"].fill(dy, Q0, weight);
          if ( maximumGapQ0 <= Q0 ) { _h_tmp_events_Q0_dySlices["gap"].fill(dy, Q0, weight); }
        }

      // Fill histograms relevant for comparison with 2011 data
      } else if ( fiducialRegion == _fiducialRegions.at(1) ) {
        // Fill inclusive and gap histograms
        foreach( const string &gapCategory, eventGapCategories ) {
          _h_tmp_events_pTbar[gapCategory]->fill( pTbar, weight );
          _h_profiled_cosDeltaPhi_pTbar[gapCategory]->fill( pTbar, cos(M_PI - dphi), weight );
          _h_profiled_cosTwoDeltaPhi_pTbar[gapCategory]->fill( pTbar, cos(2.0 * dphi), weight );
        }
        // Fill profiled nJets_rapidity_interval
        _h_profiled_nJets_rapidity_interval_pTbar->fill( pTbar, nJets_rapidity_interval, weight );
      }
      return;
    }

    /// Normalise histograms etc., after the run
    void finalize() {
      // Normalise cross-section plots to correct cross-section
      const double xs_pb( crossSection()/picobarn );
      const double nEventsGenerated( sumOfWeights() );
      double xs_norm_factor( xs_pb/nEventsGenerated );
      const double ySpan(1); // all dy spans are 1
      foreach (const string& gapCategory, _gapCategories ) {
        _h_crossSection_dphi_dySlices[gapCategory].scale(xs_norm_factor/ySpan/M_PI, this);
      }

      // Create C2/C1 scatter from profiles
      foreach (const string& gapCategory, _gapCategories ) {
        divide( _h_profiled_cosTwoDeltaPhi_dy[gapCategory], _h_profiled_cosDeltaPhi_dy[gapCategory], _h_C2C1_dy[gapCategory] );
        divide( _h_profiled_cosTwoDeltaPhi_pTbar[gapCategory], _h_profiled_cosDeltaPhi_pTbar[gapCategory], _h_C2C1_pTbar[gapCategory] );
      }

      // Fill simple gap fractions
      Scatter2DPtr h_gap_fraction_dy    = bookScatter2D( 1, 1, 1);
      Scatter2DPtr h_gap_fraction_pTbar = bookScatter2D( 2, 1, 1, false);
      Rivet::Analysis::efficiency( *_h_tmp_events_dy["gap"].get(), *_h_tmp_events_dy["inclusive"].get(), h_gap_fraction_dy );
      Rivet::Analysis::efficiency( *_h_tmp_events_pTbar["gap"].get(), *_h_tmp_events_pTbar["inclusive"].get(), h_gap_fraction_pTbar );

      // Register and fill Q0 gap fractions
      for (unsigned int dyLow = 0; dyLow < _dy_max; ++dyLow ) {
        Scatter2DPtr h_gap_fraction_Q0 = bookScatter2D( 29+dyLow, 1, 1, false);
        Rivet::Analysis::efficiency( *_h_tmp_events_Q0_dySlices["gap"].getHistograms().at(dyLow).get(), *_h_tmp_events_Q0_dySlices["inclusive"].getHistograms().at(dyLow).get(), h_gap_fraction_Q0 );
      }

      /// Print summary information
      MSG_INFO( "Cross-Section/pb     : " << xs_pb );
      MSG_INFO( "Sum of weights       : " << sumOfWeights() << "  (" << _sumOfAcceptedWeights << " accepted)" );
      MSG_INFO( "nEvents              : " << numEvents() << "  (" << _nEventsInAcceptance << " accepted)" );
      MSG_INFO( "Normalisation factor : " << xs_norm_factor );
    }


  private:

    /// Member variables
    vector<unsigned int> _fiducialRegions;
    vector<double> _vetoScale, _yFiducial;
    vector<string> _gapCategories;

    unsigned int _dy_max;
    unsigned int _nEventsInAcceptance;
    double _sumOfAcceptedWeights;

    /// Histograms
    // Number of events : gap and non-gap : necessary input for gap fractions but not saved as output
    map<string, Histo1DPtr> _h_tmp_events_dy;
    map<string, Histo1DPtr> _h_tmp_events_pTbar;
    map<string, BinnedHistogram<double> > _h_tmp_events_Q0_dySlices;

    // Number of jets in rapidity interval
    Profile1DPtr _h_profiled_nJets_rapidity_interval_dy;
    Profile1DPtr _h_profiled_nJets_rapidity_interval_pTbar;

    // Azimuthal moment histograms
    map<string, Profile1DPtr> _h_profiled_cosDeltaPhi_dy;
    map<string, Profile1DPtr> _h_profiled_cosDeltaPhi_pTbar;
    map<string, Profile1DPtr> _h_profiled_cosTwoDeltaPhi_dy;
    map<string, Profile1DPtr> _h_profiled_cosTwoDeltaPhi_pTbar;
    map<string, Scatter2DPtr> _h_C2C1_dy;
    map<string, Scatter2DPtr> _h_C2C1_pTbar;

    // Cross-section vs. #Delta#phi
    map<string, BinnedHistogram<double> > _h_crossSection_dphi_dySlices;
  };

  // The hook for the plugin system
  DECLARE_RIVET_PLUGIN(ATLAS_2014_I1307243);

}