rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

D0_2008_S7719523

Isolated $\gamma$ + jet cross-sections, differential in pT($\gamma$) for various $y$ bins
Experiment: D0 (Tevatron Run 2)
Inspire ID: 782968
Status: VALIDATED
Authors:
  • Andy Buckley
  • Gavin Hesketh
  • Frank Siegert
References: Beams: p- p+
Beam energies: (980.0, 980.0) GeV
Run details:
  • Produce only gamma + jet ($q,\bar{q},g$) hard processes (for Pythia 6, this means MSEL=10 and MSUB indices 14, 29 & 115 enabled). The lowest bin edge is at 30 GeV, so a kinematic pTmin cut is probably required to fill the histograms.

The process $p \bar{p}$ -> photon + jet + X as studied by the D0 detector at the Fermilab Tevatron collider at center-of-mass energy $\sqrt{s} = 1.96\;\text{TeV}$. Photons are reconstructed in the central rapidity region $|y_\gamma| < 1.0$ with transverse momenta in the range 30--400\;GeV, while jets are reconstructed in either the central $|y_\text{jet}| < 0.8$ or forward $1.5 < |y_\text{jet}| < 2.5$ rapidity intervals with $p_\perp^\text{jet} > 15\;\text{GeV}$. The differential cross section $\mathrm{d}^3 \sigma / \mathrm{d}{p_\perp^\gamma} \mathrm{d}{y_\gamma} \mathrm{d}{y_\text{jet}}$ is measured as a function of $p_\perp^\gamma$ in four regions, differing by the relative orientations of the photon and the jet. MC predictions have trouble with simultaneously describing the measured normalization and $p_\perp^\gamma$ dependence of the cross section in any of the four measured regions.

Source code: D0_2008_S7719523.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
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/LeadingParticlesFinalState.hh"
#include "Rivet/Projections/VetoedFinalState.hh"
#include "Rivet/Projections/FastJets.hh"

namespace Rivet {


  // A local scope function for division, handling the div-by-zero case
  /// @todo Why isn't the math divide() function being found?
  namespace {
    inline double _safediv(double a, double b, double result_if_err) {
      return (b != 0) ? a/b : result_if_err;
    }
  }


  /// @brief Measurement of isolated gamma + jet + X differential cross-sections
  ///
  /// Inclusive isolated gamma + jet cross-sections, differential in pT(gamma), for
  /// various photon and jet rapidity bins.
  ///
  /// @author Andy Buckley
  /// @author Gavin Hesketh
  class D0_2008_S7719523 : public Analysis {
  public:

    /// @name Constructors etc.
    //@{

    /// Constructor
    D0_2008_S7719523()
      : Analysis("D0_2008_S7719523")
    {    }

    //@}


    /// @name Analysis methods
    //@{

    /// Set up projections and book histograms
    void init() {
      // General FS
      FinalState fs;
      declare(fs, "FS");

      // Get leading photon
      LeadingParticlesFinalState photonfs(FinalState(-1.0, 1.0, 30.0*GeV));
      photonfs.addParticleId(PID::PHOTON);
      declare(photonfs, "LeadingPhoton");

      // FS excluding the leading photon
      VetoedFinalState vfs(fs);
      vfs.addVetoOnThisFinalState(photonfs);
      declare(vfs, "JetFS");

      // Jets
      FastJets jetpro(vfs, FastJets::D0ILCONE, 0.7);
      declare(jetpro, "Jets");

      // Histograms
      _h_central_same_cross_section = bookHisto1D(1, 1, 1);
      _h_central_opp_cross_section  = bookHisto1D(2, 1, 1);
      _h_forward_same_cross_section = bookHisto1D(3, 1, 1);
      _h_forward_opp_cross_section  = bookHisto1D(4, 1, 1);

      // Ratio histos to be filled by divide()
      _h_cen_opp_same = bookScatter2D(5, 1, 1);
      _h_fwd_opp_same = bookScatter2D(8, 1, 1);
      // Ratio histos to be filled manually, since the num/denom inputs don't match
      _h_cen_same_fwd_same = bookScatter2D(6, 1, 1, true);
      _h_cen_opp_fwd_same = bookScatter2D(7, 1, 1, true);
      _h_cen_same_fwd_opp = bookScatter2D(9, 1, 1, true);
      _h_cen_opp_fwd_opp = bookScatter2D(10, 1, 1, true);
    }



    /// Do the analysis
    void analyze(const Event& event) {
      const double weight = event.weight();

      // Get the photon
      const FinalState& photonfs = apply<FinalState>(event, "LeadingPhoton");
      if (photonfs.particles().size() != 1) {
        vetoEvent;
      }
      const FourMomentum photon = photonfs.particles().front().momentum();

      // Isolate photon by ensuring that a 0.4 cone around it contains less than 7% of the photon's energy
      double egamma = photon.E();
      double eta_P = photon.eta();
      double phi_P = photon.phi();
      double econe = 0.0;
      foreach (const Particle& p, apply<FinalState>(event, "JetFS").particles()) {
        if (deltaR(eta_P, phi_P, p.eta(), p.phi()) < 0.4) {
          econe += p.E();
          // Veto as soon as E_cone gets larger
          if (econe/egamma > 0.07) {
            MSG_DEBUG("Vetoing event because photon is insufficiently isolated");
            vetoEvent;
          }
        }
      }

      Jets jets = apply<FastJets>(event, "Jets").jetsByPt(15.0*GeV);
      if (jets.empty()) vetoEvent;

      FourMomentum leadingJet = jets[0].momentum();
      if (deltaR(eta_P, phi_P, leadingJet.eta(), leadingJet.phi()) < 0.7) {
        vetoEvent;
      }

      int photon_jet_sign = sign( leadingJet.rapidity() * photon.rapidity() );

      // Veto if leading jet is outside plotted rapidity regions
      const double abs_y1 = fabs(leadingJet.rapidity());
      if (inRange(abs_y1, 0.8, 1.5) || abs_y1 > 2.5) {
        MSG_DEBUG("Leading jet falls outside acceptance range; |y1| = " << abs_y1);
        vetoEvent;
      }

      // Fill histos
      if (fabs(leadingJet.rapidity()) < 0.8) {
        Histo1DPtr h = (photon_jet_sign >= 1) ? _h_central_same_cross_section : _h_central_opp_cross_section;
        h->fill(photon.pT(), weight);
      } else if (inRange( fabs(leadingJet.rapidity()), 1.5, 2.5)) {
        Histo1DPtr h = (photon_jet_sign >= 1) ? _h_forward_same_cross_section : _h_forward_opp_cross_section;
        h->fill(photon.pT(), weight);
      }

    }


    /// Finalize
    void finalize() {
      const double lumi_gen = sumOfWeights()/crossSection();
      const double dy_photon = 2.0;
      const double dy_jet_central = 1.6;
      const double dy_jet_forward = 2.0;

      // Cross-section ratios (6 plots)
      // Central/central and forward/forward ratios
      divide(_h_central_opp_cross_section, _h_central_same_cross_section, _h_cen_opp_same);
      divide(_h_forward_opp_cross_section, _h_forward_same_cross_section, _h_fwd_opp_same);
      // Central/forward ratio combinations
      /// @note The central/forward histo binnings are not the same! Hence the need to do these by hand :-(
      for (size_t i = 0; i < _h_cen_same_fwd_same->numPoints(); ++i) {
        const YODA::HistoBin1D& cen_same_bini = _h_central_same_cross_section->bin(i);
        const YODA::HistoBin1D& cen_opp_bini = _h_central_opp_cross_section->bin(i);
        const YODA::HistoBin1D& fwd_same_bini = _h_central_same_cross_section->bin(i);
        const YODA::HistoBin1D& fwd_opp_bini = _h_central_opp_cross_section->bin(i);
        _h_cen_same_fwd_same->point(i).setY(_safediv(cen_same_bini.sumW(), fwd_same_bini.sumW(), 0),
                                            add_quad(cen_same_bini.relErr(), fwd_same_bini.relErr()));
        _h_cen_opp_fwd_same->point(i).setY(_safediv(cen_opp_bini.sumW(), fwd_same_bini.sumW(), 0),
                                           add_quad(cen_opp_bini.relErr(), fwd_same_bini.relErr()));
        _h_cen_same_fwd_opp->point(i).setY(_safediv(cen_same_bini.sumW(), fwd_opp_bini.sumW(), 0),
                                           add_quad(cen_same_bini.relErr(), fwd_opp_bini.relErr()));
        _h_cen_opp_fwd_opp->point(i).setY(_safediv(cen_opp_bini.sumW(), fwd_opp_bini.sumW(), 0),
                                          add_quad(cen_opp_bini.relErr(), fwd_opp_bini.relErr()));
      }

      // Use generator cross section for remaining histograms
      // Each of these needs the additional factor 2 because the
      // y_photon * y_jet requirement reduces the corresponding 2D "bin width"
      // by a factor 1/2.
      scale(_h_central_same_cross_section, 2.0/lumi_gen * 1.0/dy_photon * 1.0/dy_jet_central);
      scale(_h_central_opp_cross_section, 2.0/lumi_gen * 1.0/dy_photon * 1.0/dy_jet_central);
      scale(_h_forward_same_cross_section, 2.0/lumi_gen * 1.0/dy_photon * 1.0/dy_jet_forward);
      scale(_h_forward_opp_cross_section, 2.0/lumi_gen * 1.0/dy_photon * 1.0/dy_jet_forward);
    }

    //@}

  private:

    /// @name Histograms
    //@{
    Histo1DPtr _h_central_same_cross_section;
    Histo1DPtr _h_central_opp_cross_section;
    Histo1DPtr _h_forward_same_cross_section;
    Histo1DPtr _h_forward_opp_cross_section;

    Scatter2DPtr _h_cen_opp_same;
    Scatter2DPtr _h_fwd_opp_same;
    Scatter2DPtr _h_cen_same_fwd_same;
    Scatter2DPtr _h_cen_opp_fwd_same;
    Scatter2DPtr _h_cen_same_fwd_opp;
    Scatter2DPtr _h_cen_opp_fwd_opp;
    //@}

  };



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

}