rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

D0_2009_S8320160

Dijet angular distributions
Experiment: D0 (Tevatron Run 2)
Inspire ID: 824127
Status: VALIDATED
Authors:
  • Frank Siegert
References: Beams: p- p+
Beam energies: (980.0, 980.0) GeV
Run details:
  • $p \bar{p} \to$ jets at 1960 GeV

Dijet angular distributions in different bins of dijet mass from 0.25 TeV to above 1.1 TeV in $p \bar{p}$ collisions at $\sqrt{s} = 1.96$ TeV, based on an integrated luminosity of 0.7 fb$^{-1}$.

Source code: D0_2009_S8320160.cc
 1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Tools/BinnedHistogram.hh"
 4#include "Rivet/Projections/FinalState.hh"
 5#include "Rivet/Projections/ChargedFinalState.hh"
 6#include "Rivet/Projections/FastJets.hh"
 7
 8namespace Rivet {
 9
10
11  /// D0 dijet angular distributions
12  class D0_2009_S8320160 : public Analysis {
13  public:
14
15    RIVET_DEFAULT_ANALYSIS_CTOR(D0_2009_S8320160);
16
17
18    /// @name Analysis methods
19    /// @{
20
21    // Book histograms
22    void init() {
23      FinalState fs;
24      FastJets conefinder(fs, FastJets::D0ILCONE, 0.7);
25      declare(conefinder, "ConeFinder");
26
27      {Histo1DPtr tmp; _h_chi_dijet.add(250., 300., book(tmp, 1, 1, 1));}
28      {Histo1DPtr tmp; _h_chi_dijet.add(300., 400., book(tmp, 2, 1, 1));}
29      {Histo1DPtr tmp; _h_chi_dijet.add(400., 500., book(tmp, 3, 1, 1));}
30      {Histo1DPtr tmp; _h_chi_dijet.add(500., 600., book(tmp, 4, 1, 1));}
31      {Histo1DPtr tmp; _h_chi_dijet.add(600., 700., book(tmp, 5, 1, 1));}
32      {Histo1DPtr tmp; _h_chi_dijet.add(700., 800., book(tmp, 6, 1, 1));}
33      {Histo1DPtr tmp; _h_chi_dijet.add(800., 900., book(tmp, 7, 1, 1));}
34      {Histo1DPtr tmp; _h_chi_dijet.add(900.,1000., book(tmp, 8, 1, 1));}
35      {Histo1DPtr tmp; _h_chi_dijet.add(1000.,1100.,book(tmp, 9, 1, 1));}
36      {Histo1DPtr tmp; _h_chi_dijet.add(1100.,1960, book(tmp, 10, 1, 1));}
37    }
38
39
40    /// Do the analysis
41    void analyze(const Event & e) {
42      const double weight = 1.0;
43
44      const Jets& jets = apply<JetAlg>(e, "ConeFinder").jetsByPt();
45      if (jets.size() < 2) vetoEvent;
46
47      FourMomentum j0(jets[0].momentum());
48      FourMomentum j1(jets[1].momentum());
49      double y0 = j0.rapidity();
50      double y1 = j1.rapidity();
51
52      if (fabs(y0+y1)>2) vetoEvent;
53
54      double mjj = FourMomentum(j0+j1).mass();
55      double chi = exp(fabs(y0-y1));
56      if(chi<16.)  _h_chi_dijet.fill(mjj, chi, weight);
57    }
58
59
60    /// Finalize
61    void finalize() {
62      for (Histo1DPtr hist : _h_chi_dijet.histos()) {
63        normalize(hist);
64      }
65    }
66
67    /// @}
68
69
70  private:
71
72    /// @name Histograms
73    /// @{
74    BinnedHistogram _h_chi_dijet;
75    /// @}
76
77  };
78
79
80
81  RIVET_DECLARE_ALIASED_PLUGIN(D0_2009_S8320160, D0_2009_I824127);
82
83}