rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

D0_1996_S3324664

Azimuthal decorrelation of jets widely separated in rapidity
Experiment: D0 (Tevatron Run 1)
Inspire ID: 416886
Status: VALIDATED
Authors:
  • Frank Siegert
References: Beams: p- p+
Beam energies: (900.0, 900.0) GeV
Run details:
  • $p \bar{p} \to \,\text{jets}$ at 1800 \text{GeV}

First measurement of the azimuthal decorrelation between jets with pseudorapidity separation up to five units. The data were accumulated using the D0 detector during Tevatron Run 1 at $\sqrt{s} = 1.8 \text{TeV}$. Requires di-jet events with at least one jet above 50 GeV and two jets above 20 GeV.

Source code: D0_1996_S3324664.cc
 1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Tools/BinnedHistogram.hh"
 4#include "Rivet/Projections/FastJets.hh"
 5#include "Rivet/Projections/FinalState.hh"
 6
 7namespace Rivet {
 8
 9
10  /// D0 azimuthal correlation of jets widely separated in rapidity
11  class D0_1996_S3324664 : public Analysis {
12  public:
13
14    RIVET_DEFAULT_ANALYSIS_CTOR(D0_1996_S3324664);
15
16
17    /// @name Analysis methods
18    /// @{
19
20    void init() {
21      const FinalState fs;
22      declare(fs, "FS");
23      /// @todo Use correct jet algorithm
24      declare(FastJets(fs, FastJets::D0ILCONE, 0.7), "ConeJets");
25
26      book(_h_deta ,1, 1, 1);
27      {Histo1DPtr tmp; _h_dphi.add(0.0, 2.0, book(tmp, 2, 1, 1));}
28      {Histo1DPtr tmp; _h_dphi.add(2.0, 4.0, book(tmp, 2, 1, 2));}
29      {Histo1DPtr tmp; _h_dphi.add(4.0, 6.0, book(tmp, 2, 1, 3));}
30      book(_h_cosdphi_deta ,3, 1, 1);
31    }
32
33
34    void analyze(const Event& event) {
35      const double weight = 1.0;
36      Jets jets = apply<FastJets>(event, "ConeJets").jets(Cuts::Et > 20*GeV && Cuts::abseta<3, cmpMomByEt);
37
38      if (jets.size() < 2) vetoEvent;
39
40      FourMomentum minjet = jets[0].momentum();
41      FourMomentum maxjet = jets[1].momentum();
42      double mineta = minjet.eta();
43      double maxeta = maxjet.eta();
44
45      for (const Jet& jet : jets) {
46        double eta = jet.eta();
47        if (eta < mineta) {
48          minjet = jet.momentum();
49          mineta = eta;
50        } else if (eta > maxeta) {
51          maxjet = jet.momentum();
52          maxeta = eta;
53        }
54      }
55
56      if (minjet.Et() < 50*GeV && maxjet.Et() < 50.0*GeV) vetoEvent;
57
58      double deta = maxjet.eta()-minjet.eta();
59      double dphi = mapAngle0To2Pi(maxjet.phi()-minjet.phi());
60
61      _h_deta->fill(deta, weight);
62      _h_dphi.fill(deta, 1.0-dphi/M_PI, weight);
63      _h_cosdphi_deta->fill(deta, cos(M_PI-dphi), weight);
64    }
65
66
67    void finalize() {
68      // Normalised to #events
69      normalize(_h_deta, 8830.); // fixed norm OK
70
71      // Normalied to 1/(4pi)
72      for (Histo1DPtr histo : _h_dphi.histos()) {
73        normalize(histo, 1./(4.*M_PI));
74      }
75
76    }
77
78    /// @}
79
80
81  private:
82
83    /// @name Histograms
84    /// @{
85    Histo1DPtr _h_deta;
86    BinnedHistogram _h_dphi;
87    Profile1DPtr _h_cosdphi_deta;
88    /// @}
89
90  };
91
92
93
94  RIVET_DECLARE_ALIASED_PLUGIN(D0_1996_S3324664, D0_1996_I416886);
95
96}