rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CDF_2005_S6217184

CDF Run II jet shape analysis
Experiment: CDF (Tevatron Run 2)
Inspire ID: 682179
Status: VALIDATED
Authors:
  • Lars Sonnenschein
  • Andy Buckley
References: Beams: p- p+
Beam energies: (980.0, 980.0) GeV
Run details:
  • QCD events at $\sqrt{s} = 1960$ GeV. Jet $p_\perp^\text{min}$ in plots is 37 GeV/$c$ -- choose a generator min $p_\perp$ well below this.

Measurement of jet shapes in inclusive jet production in $p \bar{p}$ collisions at center-of-mass energy $\sqrt{s} = 1.96$ TeV. The data cover jet transverse momenta from 37--380 GeV and absolute jet rapidities in the range 0.1--0.7.

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

namespace Rivet {


  /// @brief CDF Run II jet shape analysis
  /// @author Andy Buckley
  class CDF_2005_S6217184 : public Analysis {
  public:

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


    /// @name Analysis methods
    //@{

    void init() {
      // Set up projections
      const FinalState fs(-2.0, 2.0);
      declare(fs, "FS");
      FastJets fj(fs, FastJets::CDFMIDPOINT, 0.7);
      fj.useInvisibles();
      declare(fj, "Jets");

      // Specify pT bins
      _ptedges = {{ 37.0, 45.0, 55.0, 63.0, 73.0, 84.0, 97.0, 112.0, 128.0, 148.0,
                    166.0, 186.0, 208.0, 229.0, 250.0, 277.0, 304.0, 340.0, 380.0 }};

      // Register a jet shape projection and histogram for each pT bin
      for (size_t i = 0; i < 6; ++i) {
        for (size_t j = 0; j < 3; ++j) {
          const size_t k = i*3 + j;
          stringstream ss; ss << "JetShape" << k;
          const string pname = ss.str();
          _jsnames_pT[k] = pname;
          const JetShape jsp(fj, 0.0, 0.7, 7, _ptedges[k], _ptedges[k+1], 0.1, 0.7, RAPIDITY);
          declare(jsp, pname);
          _profhistRho_pT[k] = bookProfile1D(i+1, 1, j+1);
          _profhistPsi_pT[k] = bookProfile1D(6+i+1, 1, j+1);
        }
      }

      // Final histo
      _profhistPsi_vs_pT = bookScatter2D(13, 1, 1, true);
    }



    /// Do the analysis
    void analyze(const Event& evt) {

      // Get jets and require at least one to pass pT and y cuts
      const Jets jets = apply<FastJets>(evt, "Jets")
        .jetsByPt(Cuts::ptIn(_ptedges.front()*GeV, _ptedges.back()*GeV) && Cuts::absrap < 0.7);
      MSG_DEBUG("Jet multiplicity before cuts = " << jets.size());
      if (jets.size() == 0) {
        MSG_DEBUG("No jets found in required pT and rapidity range");
        vetoEvent;
      }

      // Calculate and histogram jet shapes
      const double weight = evt.weight();
      for (size_t ipt = 0; ipt < 18; ++ipt) {
        const JetShape& jsipt = apply<JetShape>(evt, _jsnames_pT[ipt]);
        for (size_t ijet = 0; ijet < jsipt.numJets(); ++ijet) {
          for (size_t rbin = 0; rbin < jsipt.numBins(); ++rbin) {
            const double r_rho = jsipt.rBinMid(rbin);
            MSG_DEBUG(ipt << " " << rbin << " (" << r_rho << ") " << jsipt.diffJetShape(ijet, rbin));
            /// @note Bin width Jacobian factor of 0.7/0.1 = 7 in the differential shapes plot
            _profhistRho_pT[ipt]->fill(r_rho/0.7, (0.7/0.1)*jsipt.diffJetShape(ijet, rbin), weight);
            const double r_Psi = jsipt.rBinMax(rbin);
            _profhistPsi_pT[ipt]->fill(r_Psi/0.7, jsipt.intJetShape(ijet, rbin), weight);
          }
        }
      }

    }


    // Finalize
    void finalize() {
      // Construct final 1-Psi(0.3/0.7) profile from Psi profiles
      for (size_t i = 0; i < _ptedges.size()-1; ++i) {
        // Get entry for rad_Psi = 0.2 bin
        /// @note Not a great handling of empty bins!
        Profile1DPtr ph_i = _profhistPsi_pT[i];
        const double y  = (ph_i->bin(2).effNumEntries() > 0) ? ph_i->bin(2).mean() : 0;
        const double ey = (ph_i->bin(2).effNumEntries() > 1) ? ph_i->bin(2).stdErr() : 0;
        _profhistPsi_vs_pT->point(i).setY(y, ey);
      }
    }

    //@}


  private:

    /// @name Analysis data
    //@{

    /// Jet \f$ p_\perp\f$ bins.
    vector<double> _ptedges; // This can't be a raw array if we want to initialise it non-painfully

    /// JetShape projection name for each \f$p_\perp\f$ bin.
    string _jsnames_pT[18];

    //@}


    /// @name Histograms
    //@{
    Profile1DPtr _profhistRho_pT[18];
    Profile1DPtr _profhistPsi_pT[18];
    Scatter2DPtr _profhistPsi_vs_pT;
    //@}

  };



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

}