Rivet analyses

Dijet cross-sections in photoproduction

Experiment: H1 (HERA)

Inspire ID: 581409

Status: VALIDATED

Authors: - Vithyaban Anjelo Narendra

References: - Eur.Phys.J.C 25 (2002) 13-23 - arXiv: hep-ex/0201006

Beams: p+ e+, e+ p+, p+ e-, e- p+

Beam energies: (820.0, 27.6); (27.6, 820.0); (820.0, 27.6); (27.6, 820.0)GeV

Run details: - Photoproduction in HERA I

Dijet cross sections as functions of several jet observables are measured in photoproduction using the H1 detector at HERA. The data sample comprises e+p data with an integrated luminosity of 34.9 pb−1. Jets are selected using the inclusive k algorithm with a minimum transverse energy of 25 GeV for the leading jet. The phase space covers longitudinal proton momentum fraction xp and photon longitudinal momentum fraction xγ in the ranges 0.05 < xp < 0.6 and 0.1 < xγ < 1. The predictions of next-to-leading order perturbative QCD, including recent photon and proton parton densities, are found to be compatible with the data in a wide kinematical range.

Source code:H1_2002_I581409.cc

// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/DISFinalState.hh"

namespace Rivet {

  /// @brief Dijet cross-sections in photoproduction
  class H1_2002_I581409 : public Analysis {
  public:

    /// Constructor
    RIVET_DEFAULT_ANALYSIS_CTOR(H1_2002_I581409);


    /// @name Analysis methods
    /// @{

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

      // DIS final state
      const DISFinalState dfs(DISFrame::LAB);
      declare(dfs, "dfs");

      declare(FastJets(dfs, JetAlg::KT, 1.0), "jets");

      // Book histograms
      // Table1, T2, F2
      book(_h["mjj"], 1, 1, 1);

      // Table2, T2, F2A
      book(_h["Et_mean"], 2, 1, 1);

      // Table3, T2, F2B
      book(_h["Et_max"], 3, 1, 1);

      // Table4 , T3 F3AB
      book(_h["eta_bar_A"], 4, 1, 1);
      book(_h["eta_bar_B"], 4, 1, 2);

      // Table5, T3, F3CD
      book(_h["eta_bar_C"], 5, 1, 1);
      book(_h["eta_bar_D"], 5, 1, 2);

      // Table6, T4, F4AB
      book(_h["x_gamma_A"], 6, 1, 1);
      book(_h["x_gamma_B"], 6, 1, 2);

      // Table7, T5, F4CD
      book(_h["x_p_C"], 7, 1, 1);
      book(_h["x_p_D"], 7, 1, 2);

      // Table8, T6, F5,6
      book(_h["x_gamma_bylowEt"],  8, 1, 1);
      book(_h["x_gamma_byhighEt"], 8, 1, 2);

      // Table9, T7, F8AB
      book(_h["cos_theta_A"], 9, 1, 1);
      book(_h["cos_theta_B"], 9, 1, 2);

      // Table10, T7, F8CD
      book(_h["cos_theta_C"], 10, 1, 1);
      book(_h["cos_theta_D"], 10, 1, 2);

    }


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

      // DIS Kinematics & Warnings
      const DISKinematics& kin = apply<DISFinalState>(event, "dfs").kinematics();
      if ( kin.failed() ) vetoEvent;

      // inelasticity cuts
      if (!inRange(kin.y(), 0.1, 0.9)) vetoEvent;

      // Q^2 cuts
      if (kin.Q2() > 1*GeV2) vetoEvent;

      const int orientation = kin.orientation();

      // Jets
      const Jets jets = apply<FastJets>(event, "jets").jets(cmpMomByEt);

      // Veto event by number of jets
      if (jets.size() < 2) vetoEvent;

      // Define jet numbers
      const Jet& j1 = jets[0];
      const Jet& j2 = jets[1];

      // Cuts by jet energy
      if (j1.Et() < 25*GeV) vetoEvent;
      if (j2.Et() < 15*GeV) vetoEvent;

      // Define eta (with orientations)
      const double eta1 = orientation*j1.eta(), eta2 = orientation*j2.eta();

      // Cuts by eta
      if (eta1 < -0.5 || eta1 > 2.5) vetoEvent;
      if (eta2 < -0.5 || eta2 > 2.5) vetoEvent;

      // Invariant Mass
      double mjj = (j1.mom() + j2.mom()).mass();

      //Define Et_mean
      const double Et_mean = 0.5*(j1.Et() + j2.Et());

      // Define eta_bar
      const double eta_bar = 0.5*(eta1 + eta2);

      // Define x_gamma
      const double x_gamma = (j1.Et() * exp(-eta1) + j2.Et() * exp(-eta2)) / (2*kin.y()*kin.beamLepton().E());

      // Define x_p
      const double x_p = (j1.Et() * exp(eta1) + j2.Et() * exp(eta2)) / (2*kin.beamHadron().E());

      // Define cos_theta_*
      const double cos_theta = abs(tanh(0.5*(eta1 - eta2)));

      //fill histos
      _h["mjj"]->fill(mjj/GeV);
      _h["Et_mean"]->fill(Et_mean/GeV);
      _h["Et_max"]->fill(j1.Et()/GeV);

      if (inRange(kin.y(), 0.1, 0.5)) {
        if (inRange(j1.Et(), 25*GeV, 35*GeV)) _h["eta_bar_A"]->fill(eta_bar);
        if (inRange(j1.Et(), 35*GeV, 80*GeV)) _h["eta_bar_B"]->fill(eta_bar);
      }

      if (inRange(kin.y(), 0.5, 0.9)) {
        if (inRange(j1.Et(), 25*GeV, 35*GeV)) _h["eta_bar_C"]->fill(eta_bar);
        if (inRange(j1.Et(), 35*GeV, 80*GeV)) _h["eta_bar_D"]->fill(eta_bar);
      }

      if (x_p < 0.1) _h["x_gamma_A"]->fill(x_gamma);
      if (x_p > 0.1) _h["x_gamma_B"]->fill(x_gamma);

      if (x_gamma < 0.8) _h["x_p_C"]->fill(x_p);
      if (x_gamma > 0.8) _h["x_p_D"]->fill(x_p);

      if (inRange(j1.Et(), 25*GeV, 35*GeV)) _h["x_gamma_bylowEt"]->fill(x_gamma);
      if (inRange(j1.Et(), 35*GeV, 80*GeV)) _h["x_gamma_byhighEt"]->fill(x_gamma);

      if (x_gamma < 0.8) _h["cos_theta_A"]->fill(cos_theta);
      if (x_gamma > 0.8) _h["cos_theta_B"]->fill(cos_theta);

      if (mjj > 65*GeV) {
        if (x_gamma < 0.8) _h["cos_theta_C"]->fill(cos_theta);
        if (x_gamma > 0.8) _h["cos_theta_D"]->fill(cos_theta);
      }

    }


    /// Normalise histograms etc., after the run
    void finalize() {
      scale(_h, crossSection()/picobarn/sumW());
    }

    /// @}


    /// @name Histograms
    /// @{
    map<string, Histo1DPtr> _h;
    /// @}

  };

  RIVET_DECLARE_PLUGIN(H1_2002_I581409);
}