rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2011_I919017

Measurement of ATLAS track jet properties at 7 TeV
Experiment: ATLAS (LHC)
Inspire ID: 919017
Status: VALIDATED
Authors:
  • Seth Zenz
  • Andy Buckley
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • Min bias QCD at 7 TeV.

ATLAS measurement of track jet pT, multiplicity per jet, longitudinal fragmentation, transverse momentum, radius w.r.t jet axis distributions, with jets constructed from charged tracks with $pT > 300$ MeV, using the anti-$k_T$ jet algorithm with $R = 0.4, 0.6$.

Source code: ATLAS_2011_I919017.cc
   1// -*- C++ -*-
   2#include "Rivet/Analysis.hh"
   3#include "Rivet/Projections/ChargedFinalState.hh"
   4#include "Rivet/Projections/FastJets.hh"
   5
   6namespace Rivet {
   7
   8
   9  namespace {
  10
  11    inline double calcz(const Jet& j, const Particle& p) {
  12      const double num = j.p3().dot(p.p3());
  13      const double den = j.p3().mod2();
  14      return num/den;
  15    }
  16
  17    inline double calcptrel(const Jet& j, const Particle& p) {
  18      const double num = j.p3().cross(p.p3()).mod();
  19      const double den = j.p3().mod();
  20      return num/den;
  21    }
  22
  23    inline double calcr(const Jet& j, const Particle& p) {
  24      return deltaR(j.rapidity(), j.phi(), p.rapidity(), p.phi());
  25    }
  26
  27    // For annulus area kludge
  28    /// @todo Improve somehow... need normalisation *without* bin width factors!
  29    inline double calcrweight(const Jet& j, const Particle& p) {
  30      size_t nBins_r = 26;
  31      double bins_r[] = { 0.00, 0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.10,
  32                          0.12, 0.14, 0.16, 0.18, 0.20, 0.22, 0.24, 0.26, 0.28, 0.30,
  33                          0.35, 0.40, 0.45, 0.50, 0.55, 0.60 };
  34      double r = calcr(j,p);
  35      for (size_t bin = 0 ; bin < nBins_r ; bin++) {
  36        if (r < bins_r[bin+1]) {
  37          double up = bins_r[bin+1];
  38          double down = bins_r[bin];
  39          return ((up-down)/(M_PI*(up*up-down*down)));
  40        }
  41      }
  42      return 1.0;
  43    }
  44  }
  45
  46
  47  class ATLAS_2011_I919017 : public Analysis {
  48  public:
  49
  50    /// @name Constructors etc.
  51    //@{
  52
  53    /// Constructor
  54    ATLAS_2011_I919017()
  55      : Analysis("ATLAS_2011_I919017")
  56    {    }
  57
  58    //@}
  59
  60
  61  public:
  62
  63    /// @name Analysis methods
  64    //@{
  65
  66    /// Book histograms and initialise projections before the run
  67    void init() {
  68      ChargedFinalState cfs((Cuts::etaIn(-2.5, 2.5) && Cuts::pT >=  0.3*GeV));
  69      FastJets trkjets04(cfs, FastJets::ANTIKT, 0.4);
  70      FastJets trkjets06(cfs, FastJets::ANTIKT, 0.6);
  71      declare(trkjets04, "Jets04");
  72      declare(trkjets06, "Jets06");
  73
  74      // Book histograms
  75      book(_h_pt04_00_05 ,1, 1, 1);
  76      book(_h_pt06_00_05 ,2, 1, 1);
  77      book(_h_N04_00_05_04_06 ,1, 2, 1+5);
  78      book(_h_N06_00_05_04_06 ,2, 2, 1+5);
  79      book(_h_N04_00_05_06_10 ,1, 2, 2+5);
  80      book(_h_N06_00_05_06_10 ,2, 2, 2+5);
  81      book(_h_N04_00_05_10_15 ,1, 2, 3+5);
  82      book(_h_N06_00_05_10_15 ,2, 2, 3+5);
  83      book(_h_N04_00_05_15_24 ,1, 2, 4+5);
  84      book(_h_N06_00_05_15_24 ,2, 2, 4+5);
  85      book(_h_N04_00_05_24_40 ,1, 2, 5+5);
  86      book(_h_N06_00_05_24_40 ,2, 2, 5+5);
  87      book(_h_z04_00_05_04_06 ,1, 3, 1+5);
  88      book(_h_z06_00_05_04_06 ,2, 3, 1+5);
  89      book(_h_z04_00_05_06_10 ,1, 3, 2+5);
  90      book(_h_z06_00_05_06_10 ,2, 3, 2+5);
  91      book(_h_z04_00_05_10_15 ,1, 3, 3+5);
  92      book(_h_z06_00_05_10_15 ,2, 3, 3+5);
  93      book(_h_z04_00_05_15_24 ,1, 3, 4+5);
  94      book(_h_z06_00_05_15_24 ,2, 3, 4+5);
  95      book(_h_z04_00_05_24_40 ,1, 3, 5+5);
  96      book(_h_z06_00_05_24_40 ,2, 3, 5+5);
  97      book(_h_ptrel04_00_05_04_06 ,1, 4, 1+5);
  98      book(_h_ptrel06_00_05_04_06 ,2, 4, 1+5);
  99      book(_h_ptrel04_00_05_06_10 ,1, 4, 2+5);
 100      book(_h_ptrel06_00_05_06_10 ,2, 4, 2+5);
 101      book(_h_ptrel04_00_05_10_15 ,1, 4, 3+5);
 102      book(_h_ptrel06_00_05_10_15 ,2, 4, 3+5);
 103      book(_h_ptrel04_00_05_15_24 ,1, 4, 4+5);
 104      book(_h_ptrel06_00_05_15_24 ,2, 4, 4+5);
 105      book(_h_ptrel04_00_05_24_40 ,1, 4, 5+5);
 106      book(_h_ptrel06_00_05_24_40 ,2, 4, 5+5);
 107      book(_h_rdA04_00_05_04_06 ,1, 5, 1+5);
 108      book(_h_rdA06_00_05_04_06 ,2, 5, 1+5);
 109      book(_h_rdA04_00_05_06_10 ,1, 5, 2+5);
 110      book(_h_rdA06_00_05_06_10 ,2, 5, 2+5);
 111      book(_h_rdA04_00_05_10_15 ,1, 5, 3+5);
 112      book(_h_rdA06_00_05_10_15 ,2, 5, 3+5);
 113      book(_h_rdA04_00_05_15_24 ,1, 5, 4+5);
 114      book(_h_rdA06_00_05_15_24 ,2, 5, 4+5);
 115      book(_h_rdA04_00_05_24_40 ,1, 5, 5+5);
 116      book(_h_rdA06_00_05_24_40 ,2, 5, 5+5);
 117
 118      book(_h_pt04_05_10 ,1, 1, 2);
 119      book(_h_pt06_05_10 ,2, 1, 2);
 120      book(_h_N04_05_10_04_06 ,1, 2, 1+10);
 121      book(_h_N06_05_10_04_06 ,2, 2, 1+10);
 122      book(_h_N04_05_10_06_10 ,1, 2, 2+10);
 123      book(_h_N06_05_10_06_10 ,2, 2, 2+10);
 124      book(_h_N04_05_10_10_15 ,1, 2, 3+10);
 125      book(_h_N06_05_10_10_15 ,2, 2, 3+10);
 126      book(_h_N04_05_10_15_24 ,1, 2, 4+10);
 127      book(_h_N06_05_10_15_24 ,2, 2, 4+10);
 128      book(_h_N04_05_10_24_40 ,1, 2, 5+10);
 129      book(_h_N06_05_10_24_40 ,2, 2, 5+10);
 130      book(_h_z04_05_10_04_06 ,1, 3, 1+10);
 131      book(_h_z06_05_10_04_06 ,2, 3, 1+10);
 132      book(_h_z04_05_10_06_10 ,1, 3, 2+10);
 133      book(_h_z06_05_10_06_10 ,2, 3, 2+10);
 134      book(_h_z04_05_10_10_15 ,1, 3, 3+10);
 135      book(_h_z06_05_10_10_15 ,2, 3, 3+10);
 136      book(_h_z04_05_10_15_24 ,1, 3, 4+10);
 137      book(_h_z06_05_10_15_24 ,2, 3, 4+10);
 138      book(_h_z04_05_10_24_40 ,1, 3, 5+10);
 139      book(_h_z06_05_10_24_40 ,2, 3, 5+10);
 140      book(_h_ptrel04_05_10_04_06 ,1, 4, 1+10);
 141      book(_h_ptrel06_05_10_04_06 ,2, 4, 1+10);
 142      book(_h_ptrel04_05_10_06_10 ,1, 4, 2+10);
 143      book(_h_ptrel06_05_10_06_10 ,2, 4, 2+10);
 144      book(_h_ptrel04_05_10_10_15 ,1, 4, 3+10);
 145      book(_h_ptrel06_05_10_10_15 ,2, 4, 3+10);
 146      book(_h_ptrel04_05_10_15_24 ,1, 4, 4+10);
 147      book(_h_ptrel06_05_10_15_24 ,2, 4, 4+10);
 148      book(_h_ptrel04_05_10_24_40 ,1, 4, 5+10);
 149      book(_h_ptrel06_05_10_24_40 ,2, 4, 5+10);
 150      book(_h_rdA04_05_10_04_06 ,1, 5, 1+10);
 151      book(_h_rdA06_05_10_04_06 ,2, 5, 1+10);
 152      book(_h_rdA04_05_10_06_10 ,1, 5, 2+10);
 153      book(_h_rdA06_05_10_06_10 ,2, 5, 2+10);
 154      book(_h_rdA04_05_10_10_15 ,1, 5, 3+10);
 155      book(_h_rdA06_05_10_10_15 ,2, 5, 3+10);
 156      book(_h_rdA04_05_10_15_24 ,1, 5, 4+10);
 157      book(_h_rdA06_05_10_15_24 ,2, 5, 4+10);
 158      book(_h_rdA04_05_10_24_40 ,1, 5, 5+10);
 159      book(_h_rdA06_05_10_24_40 ,2, 5, 5+10);
 160
 161      book(_h_pt04_10_15 ,1, 1, 3);
 162      book(_h_pt06_10_15 ,2, 1, 3);
 163      book(_h_N04_10_15_04_06 ,1, 2, 1+15);
 164      book(_h_N06_10_15_04_06 ,2, 2, 1+15);
 165      book(_h_N04_10_15_06_10 ,1, 2, 2+15);
 166      book(_h_N06_10_15_06_10 ,2, 2, 2+15);
 167      book(_h_N04_10_15_10_15 ,1, 2, 3+15);
 168      book(_h_N06_10_15_10_15 ,2, 2, 3+15);
 169      book(_h_N04_10_15_15_24 ,1, 2, 4+15);
 170      book(_h_N06_10_15_15_24 ,2, 2, 4+15);
 171      book(_h_N04_10_15_24_40 ,1, 2, 5+15);
 172      book(_h_N06_10_15_24_40 ,2, 2, 5+15);
 173      book(_h_z04_10_15_04_06 ,1, 3, 1+15);
 174      book(_h_z06_10_15_04_06 ,2, 3, 1+15);
 175      book(_h_z04_10_15_06_10 ,1, 3, 2+15);
 176      book(_h_z06_10_15_06_10 ,2, 3, 2+15);
 177      book(_h_z04_10_15_10_15 ,1, 3, 3+15);
 178      book(_h_z06_10_15_10_15 ,2, 3, 3+15);
 179      book(_h_z04_10_15_15_24 ,1, 3, 4+15);
 180      book(_h_z06_10_15_15_24 ,2, 3, 4+15);
 181      book(_h_z04_10_15_24_40 ,1, 3, 5+15);
 182      book(_h_z06_10_15_24_40 ,2, 3, 5+15);
 183      book(_h_ptrel04_10_15_04_06 ,1, 4, 1+15);
 184      book(_h_ptrel06_10_15_04_06 ,2, 4, 1+15);
 185      book(_h_ptrel04_10_15_06_10 ,1, 4, 2+15);
 186      book(_h_ptrel06_10_15_06_10 ,2, 4, 2+15);
 187      book(_h_ptrel04_10_15_10_15 ,1, 4, 3+15);
 188      book(_h_ptrel06_10_15_10_15 ,2, 4, 3+15);
 189      book(_h_ptrel04_10_15_15_24 ,1, 4, 4+15);
 190      book(_h_ptrel06_10_15_15_24 ,2, 4, 4+15);
 191      book(_h_ptrel04_10_15_24_40 ,1, 4, 5+15);
 192      book(_h_ptrel06_10_15_24_40 ,2, 4, 5+15);
 193      book(_h_rdA04_10_15_04_06 ,1, 5, 1+15);
 194      book(_h_rdA06_10_15_04_06 ,2, 5, 1+15);
 195      book(_h_rdA04_10_15_06_10 ,1, 5, 2+15);
 196      book(_h_rdA06_10_15_06_10 ,2, 5, 2+15);
 197      book(_h_rdA04_10_15_10_15 ,1, 5, 3+15);
 198      book(_h_rdA06_10_15_10_15 ,2, 5, 3+15);
 199      book(_h_rdA04_10_15_15_24 ,1, 5, 4+15);
 200      book(_h_rdA06_10_15_15_24 ,2, 5, 4+15);
 201      book(_h_rdA04_10_15_24_40 ,1, 5, 5+15);
 202      book(_h_rdA06_10_15_24_40 ,2, 5, 5+15);
 203
 204      book(_h_pt04_15_19 ,1, 1, 4);
 205      book(_h_pt06_15_19 ,2, 1, 4);
 206      book(_h_N04_15_19_04_06 ,1, 2, 1+20);
 207      book(_h_N06_15_19_04_06 ,2, 2, 1+20);
 208      book(_h_N04_15_19_06_10 ,1, 2, 2+20);
 209      book(_h_N06_15_19_06_10 ,2, 2, 2+20);
 210      book(_h_N04_15_19_10_15 ,1, 2, 3+20);
 211      book(_h_N06_15_19_10_15 ,2, 2, 3+20);
 212      book(_h_N04_15_19_15_24 ,1, 2, 4+20);
 213      book(_h_N06_15_19_15_24 ,2, 2, 4+20);
 214      book(_h_N04_15_19_24_40 ,1, 2, 5+20);
 215      book(_h_N06_15_19_24_40 ,2, 2, 5+20);
 216      book(_h_z04_15_19_04_06 ,1, 3, 1+20);
 217      book(_h_z06_15_19_04_06 ,2, 3, 1+20);
 218      book(_h_z04_15_19_06_10 ,1, 3, 2+20);
 219      book(_h_z06_15_19_06_10 ,2, 3, 2+20);
 220      book(_h_z04_15_19_10_15 ,1, 3, 3+20);
 221      book(_h_z06_15_19_10_15 ,2, 3, 3+20);
 222      book(_h_z04_15_19_15_24 ,1, 3, 4+20);
 223      book(_h_z06_15_19_15_24 ,2, 3, 4+20);
 224      book(_h_z04_15_19_24_40 ,1, 3, 5+20);
 225      book(_h_z06_15_19_24_40 ,2, 3, 5+20);
 226      book(_h_ptrel04_15_19_04_06 ,1, 4, 1+20);
 227      book(_h_ptrel06_15_19_04_06 ,2, 4, 1+20);
 228      book(_h_ptrel04_15_19_06_10 ,1, 4, 2+20);
 229      book(_h_ptrel06_15_19_06_10 ,2, 4, 2+20);
 230      book(_h_ptrel04_15_19_10_15 ,1, 4, 3+20);
 231      book(_h_ptrel06_15_19_10_15 ,2, 4, 3+20);
 232      book(_h_ptrel04_15_19_15_24 ,1, 4, 4+20);
 233      book(_h_ptrel06_15_19_15_24 ,2, 4, 4+20);
 234      book(_h_ptrel04_15_19_24_40 ,1, 4, 5+20);
 235      book(_h_ptrel06_15_19_24_40 ,2, 4, 5+20);
 236      book(_h_rdA04_15_19_04_06 ,1, 5, 1+20);
 237      book(_h_rdA06_15_19_04_06 ,2, 5, 1+20);
 238      book(_h_rdA04_15_19_06_10 ,1, 5, 2+20);
 239      book(_h_rdA06_15_19_06_10 ,2, 5, 2+20);
 240      book(_h_rdA04_15_19_10_15 ,1, 5, 3+20);
 241      book(_h_rdA06_15_19_10_15 ,2, 5, 3+20);
 242      book(_h_rdA04_15_19_15_24 ,1, 5, 4+20);
 243      book(_h_rdA06_15_19_15_24 ,2, 5, 4+20);
 244      book(_h_rdA04_15_19_24_40 ,1, 5, 5+20);
 245      book(_h_rdA06_15_19_24_40 ,2, 5, 5+20);
 246
 247      book(_h_N04_00_19_04_06 ,1, 2, 1+0);
 248      book(_h_N06_00_19_04_06 ,2, 2, 1+0);
 249      book(_h_N04_00_19_06_10 ,1, 2, 2+0);
 250      book(_h_N06_00_19_06_10 ,2, 2, 2+0);
 251      book(_h_N04_00_19_10_15 ,1, 2, 3+0);
 252      book(_h_N06_00_19_10_15 ,2, 2, 3+0);
 253      book(_h_N04_00_19_15_24 ,1, 2, 4+0);
 254      book(_h_N06_00_19_15_24 ,2, 2, 4+0);
 255      book(_h_N04_00_19_24_40 ,1, 2, 5+0);
 256      book(_h_N06_00_19_24_40 ,2, 2, 5+0);
 257      book(_h_z04_00_19_04_06 ,1, 3, 1+0);
 258      book(_h_z06_00_19_04_06 ,2, 3, 1+0);
 259      book(_h_z04_00_19_06_10 ,1, 3, 2+0);
 260      book(_h_z06_00_19_06_10 ,2, 3, 2+0);
 261      book(_h_z04_00_19_10_15 ,1, 3, 3+0);
 262      book(_h_z06_00_19_10_15 ,2, 3, 3+0);
 263      book(_h_z04_00_19_15_24 ,1, 3, 4+0);
 264      book(_h_z06_00_19_15_24 ,2, 3, 4+0);
 265      book(_h_z04_00_19_24_40 ,1, 3, 5+0);
 266      book(_h_z06_00_19_24_40 ,2, 3, 5+0);
 267      book(_h_ptrel04_00_19_04_06 ,1, 4, 1+0);
 268      book(_h_ptrel06_00_19_04_06 ,2, 4, 1+0);
 269      book(_h_ptrel04_00_19_06_10 ,1, 4, 2+0);
 270      book(_h_ptrel06_00_19_06_10 ,2, 4, 2+0);
 271      book(_h_ptrel04_00_19_10_15 ,1, 4, 3+0);
 272      book(_h_ptrel06_00_19_10_15 ,2, 4, 3+0);
 273      book(_h_ptrel04_00_19_15_24 ,1, 4, 4+0);
 274      book(_h_ptrel06_00_19_15_24 ,2, 4, 4+0);
 275      book(_h_ptrel04_00_19_24_40 ,1, 4, 5+0);
 276      book(_h_ptrel06_00_19_24_40 ,2, 4, 5+0);
 277      book(_h_rdA04_00_19_04_06 ,1, 5, 1+0);
 278      book(_h_rdA06_00_19_04_06 ,2, 5, 1+0);
 279      book(_h_rdA04_00_19_06_10 ,1, 5, 2+0);
 280      book(_h_rdA06_00_19_06_10 ,2, 5, 2+0);
 281      book(_h_rdA04_00_19_10_15 ,1, 5, 3+0);
 282      book(_h_rdA06_00_19_10_15 ,2, 5, 3+0);
 283      book(_h_rdA04_00_19_15_24 ,1, 5, 4+0);
 284      book(_h_rdA06_00_19_15_24 ,2, 5, 4+0);
 285      book(_h_rdA04_00_19_24_40 ,1, 5, 5+0);
 286      book(_h_rdA06_00_19_24_40 ,2, 5, 5+0);
 287
 288      book(_sumofweights04, "_sumofweights04");
 289      book(_sumofweights06, "_sumofweights06");
 290      book(_numjets04_00_05_04_06, "_numjets04_00_05_04_06");
 291      book(_numjets04_00_05_06_10, "_numjets04_00_05_06_10");
 292      book(_numjets04_00_05_10_15, "_numjets04_00_05_10_15");
 293      book(_numjets04_00_05_15_24, "_numjets04_00_05_15_24");
 294      book(_numjets04_00_05_24_40, "_numjets04_00_05_24_40");
 295      book(_numjets06_00_05_04_06, "_numjets06_00_05_04_06");
 296      book(_numjets06_00_05_06_10, "_numjets06_00_05_06_10");
 297      book(_numjets06_00_05_10_15, "_numjets06_00_05_10_15");
 298      book(_numjets06_00_05_15_24, "_numjets06_00_05_15_24");
 299      book(_numjets06_00_05_24_40, "_numjets06_00_05_24_40");
 300      book(_numjets04_05_10_04_06, "_numjets04_05_10_04_06");
 301      book(_numjets04_05_10_06_10, "_numjets04_05_10_06_10");
 302      book(_numjets04_05_10_10_15, "_numjets04_05_10_10_15");
 303      book(_numjets04_05_10_15_24, "_numjets04_05_10_15_24");
 304      book(_numjets04_05_10_24_40, "_numjets04_05_10_24_40");
 305      book(_numjets06_05_10_04_06, "_numjets06_05_10_04_06");
 306      book(_numjets06_05_10_06_10, "_numjets06_05_10_06_10");
 307      book(_numjets06_05_10_10_15, "_numjets06_05_10_10_15");
 308      book(_numjets06_05_10_15_24, "_numjets06_05_10_15_24");
 309      book(_numjets06_05_10_24_40, "_numjets06_05_10_24_40");
 310      book(_numjets04_10_15_04_06, "_numjets04_10_15_04_06");
 311      book(_numjets04_10_15_06_10, "_numjets04_10_15_06_10");
 312      book(_numjets04_10_15_10_15, "_numjets04_10_15_10_15");
 313      book(_numjets04_10_15_15_24, "_numjets04_10_15_15_24");
 314      book(_numjets04_10_15_24_40, "_numjets04_10_15_24_40");
 315      book(_numjets06_10_15_04_06, "_numjets06_10_15_04_06");
 316      book(_numjets06_10_15_06_10, "_numjets06_10_15_06_10");
 317      book(_numjets06_10_15_10_15, "_numjets06_10_15_10_15");
 318      book(_numjets06_10_15_15_24, "_numjets06_10_15_15_24");
 319      book(_numjets06_10_15_24_40, "_numjets06_10_15_24_40");
 320      book(_numjets04_15_19_04_06, "_numjets04_15_19_04_06");
 321      book(_numjets04_15_19_06_10, "_numjets04_15_19_06_10");
 322      book(_numjets04_15_19_10_15, "_numjets04_15_19_10_15");
 323      book(_numjets04_15_19_15_24, "_numjets04_15_19_15_24");
 324      book(_numjets04_15_19_24_40, "_numjets04_15_19_24_40");
 325      book(_numjets06_15_19_04_06, "_numjets06_15_19_04_06");
 326      book(_numjets06_15_19_06_10, "_numjets06_15_19_06_10");
 327      book(_numjets06_15_19_10_15, "_numjets06_15_19_10_15");
 328      book(_numjets06_15_19_15_24, "_numjets06_15_19_15_24");
 329      book(_numjets06_15_19_24_40, "_numjets06_15_19_24_40");
 330      book(_numjets04_00_19_04_06, "_numjets04_00_19_04_06");
 331      book(_numjets04_00_19_06_10, "_numjets04_00_19_06_10");
 332      book(_numjets04_00_19_10_15, "_numjets04_00_19_10_15");
 333      book(_numjets04_00_19_15_24, "_numjets04_00_19_15_24");
 334      book(_numjets04_00_19_24_40, "_numjets04_00_19_24_40");
 335      book(_numjets06_00_19_04_06, "_numjets06_00_19_04_06");
 336      book(_numjets06_00_19_06_10, "_numjets06_00_19_06_10");
 337      book(_numjets06_00_19_10_15, "_numjets06_00_19_10_15");
 338      book(_numjets06_00_19_15_24, "_numjets06_00_19_15_24");
 339      book(_numjets06_00_19_24_40, "_numjets06_00_19_24_40");
 340    }
 341
 342
 343    /// Perform the per-event analysis
 344    void analyze(const Event& event) {
 345      const Jets& jets04 = apply<JetAlg>(event, "Jets04").jets();
 346      if (!jets04.empty()) {
 347        _sumofweights04->fill();
 348        for (const Jet& j : jets04) {
 349          const double jetpt = j.pT();
 350          if (j.absrap() < 0.5) {
 351            _h_pt04_00_05->fill(jetpt/GeV);
 352            if (inRange(jetpt/GeV, 4., 6.)) {
 353              _numjets04_00_05_04_06->fill();
 354              _h_N04_00_05_04_06->fill(j.particles().size());
 355              for (const Particle& p : j.particles()) {
 356                _h_z04_00_05_04_06->fill(calcz(j,p));
 357                _h_ptrel04_00_05_04_06->fill(calcptrel(j,p));
 358                _h_rdA04_00_05_04_06->fill(calcr(j,p),calcrweight(j,p));
 359              }
 360            }
 361            if (inRange(jetpt/GeV, 6., 10.)) {
 362              _numjets04_00_05_06_10->fill();
 363              _h_N04_00_05_06_10->fill(j.particles().size());
 364              for (const Particle& p : j.particles()) {
 365                _h_z04_00_05_06_10->fill(calcz(j,p));
 366                _h_ptrel04_00_05_06_10->fill(calcptrel(j,p));
 367                _h_rdA04_00_05_06_10->fill(calcr(j,p),calcrweight(j,p));
 368              }
 369            }
 370            if (inRange(jetpt/GeV, 10., 15.)) {
 371              _numjets04_00_05_10_15->fill();
 372              _h_N04_00_05_10_15->fill(j.particles().size());
 373              for (const Particle& p : j.particles()) {
 374                _h_z04_00_05_10_15->fill(calcz(j,p));
 375                _h_ptrel04_00_05_10_15->fill(calcptrel(j,p));
 376                _h_rdA04_00_05_10_15->fill(calcr(j,p),calcrweight(j,p));
 377              }
 378            }
 379            if (inRange(jetpt/GeV, 15., 24.)) {
 380              _numjets04_00_05_15_24->fill();
 381              _h_N04_00_05_15_24->fill(j.particles().size());
 382              for (const Particle& p : j.particles()) {
 383                _h_z04_00_05_15_24->fill(calcz(j,p));
 384                _h_ptrel04_00_05_15_24->fill(calcptrel(j,p));
 385                _h_rdA04_00_05_15_24->fill(calcr(j,p),calcrweight(j,p));
 386              }
 387            }
 388            if (inRange(jetpt/GeV, 24., 40.)) {
 389              _numjets04_00_05_24_40->fill();
 390              _h_N04_00_05_24_40->fill(j.particles().size());
 391              for (const Particle& p : j.particles()) {
 392                _h_z04_00_05_24_40->fill(calcz(j,p));
 393                _h_ptrel04_00_05_24_40->fill(calcptrel(j,p));
 394                _h_rdA04_00_05_24_40->fill(calcr(j,p),calcrweight(j,p));
 395              }
 396            }
 397          }
 398          if (j.absrap() > 0.5 && j.absrap() < 1.0) {
 399            _h_pt04_05_10->fill(jetpt/GeV);
 400            if (inRange(jetpt/GeV, 4., 6.)) {
 401              _numjets04_05_10_04_06->fill();
 402              _h_N04_05_10_04_06->fill(j.particles().size());
 403              for (const Particle& p : j.particles()) {
 404                _h_z04_05_10_04_06->fill(calcz(j,p));
 405                _h_ptrel04_05_10_04_06->fill(calcptrel(j,p));
 406                _h_rdA04_05_10_04_06->fill(calcr(j,p),calcrweight(j,p));
 407              }
 408            }
 409            if (inRange(jetpt/GeV, 6., 10.)) {
 410              _numjets04_05_10_06_10->fill();
 411              _h_N04_05_10_06_10->fill(j.particles().size());
 412              for (const Particle& p : j.particles()) {
 413                _h_z04_05_10_06_10->fill(calcz(j,p));
 414                _h_ptrel04_05_10_06_10->fill(calcptrel(j,p));
 415                _h_rdA04_05_10_06_10->fill(calcr(j,p),calcrweight(j,p));
 416              }
 417            }
 418            if (inRange(jetpt/GeV, 10., 15.)) {
 419              _numjets04_05_10_10_15->fill();
 420              _h_N04_05_10_10_15->fill(j.particles().size());
 421              for (const Particle& p : j.particles()) {
 422                _h_z04_05_10_10_15->fill(calcz(j,p));
 423                _h_ptrel04_05_10_10_15->fill(calcptrel(j,p));
 424                _h_rdA04_05_10_10_15->fill(calcr(j,p),calcrweight(j,p));
 425              }
 426            }
 427            if (inRange(jetpt/GeV, 15., 24.)) {
 428              _numjets04_05_10_15_24->fill();
 429              _h_N04_05_10_15_24->fill(j.particles().size());
 430              for (const Particle& p : j.particles()) {
 431                _h_z04_05_10_15_24->fill(calcz(j,p));
 432                _h_ptrel04_05_10_15_24->fill(calcptrel(j,p));
 433                _h_rdA04_05_10_15_24->fill(calcr(j,p),calcrweight(j,p));
 434              }
 435            }
 436            if (inRange(jetpt/GeV, 24., 40.)) {
 437              _numjets04_05_10_24_40->fill();
 438              _h_N04_05_10_24_40->fill(j.particles().size());
 439              for (const Particle& p : j.particles()) {
 440                _h_z04_05_10_24_40->fill(calcz(j,p));
 441                _h_ptrel04_05_10_24_40->fill(calcptrel(j,p));
 442                _h_rdA04_05_10_24_40->fill(calcr(j,p),calcrweight(j,p));
 443              }
 444            }
 445          }
 446          if (j.absrap() > 1.0 && j.absrap() < 1.5) {
 447            _h_pt04_10_15->fill(jetpt/GeV);
 448            if (inRange(jetpt/GeV, 4., 6.)) {
 449              _numjets04_10_15_04_06->fill();
 450              _h_N04_10_15_04_06->fill(j.particles().size());
 451              for (const Particle& p : j.particles()) {
 452                _h_z04_10_15_04_06->fill(calcz(j,p));
 453                _h_ptrel04_10_15_04_06->fill(calcptrel(j,p));
 454                _h_rdA04_10_15_04_06->fill(calcr(j,p),calcrweight(j,p));
 455              }
 456            }
 457            if (inRange(jetpt/GeV, 6., 10.)) {
 458              _numjets04_10_15_06_10->fill();
 459              _h_N04_10_15_06_10->fill(j.particles().size());
 460              for (const Particle& p : j.particles()) {
 461                _h_z04_10_15_06_10->fill(calcz(j,p));
 462                _h_ptrel04_10_15_06_10->fill(calcptrel(j,p));
 463                _h_rdA04_10_15_06_10->fill(calcr(j,p),calcrweight(j,p));
 464              }
 465            }
 466            if (inRange(jetpt/GeV, 10., 15.)) {
 467              _numjets04_10_15_10_15->fill();
 468              _h_N04_10_15_10_15->fill(j.particles().size());
 469              for (const Particle& p : j.particles()) {
 470                _h_z04_10_15_10_15->fill(calcz(j,p));
 471                _h_ptrel04_10_15_10_15->fill(calcptrel(j,p));
 472                _h_rdA04_10_15_10_15->fill(calcr(j,p),calcrweight(j,p));
 473              }
 474            }
 475            if (inRange(jetpt/GeV, 15., 24.)) {
 476              _numjets04_10_15_15_24->fill();
 477              _h_N04_10_15_15_24->fill(j.particles().size());
 478              for (const Particle& p : j.particles()) {
 479                _h_z04_10_15_15_24->fill(calcz(j,p));
 480                _h_ptrel04_10_15_15_24->fill(calcptrel(j,p));
 481                _h_rdA04_10_15_15_24->fill(calcr(j,p),calcrweight(j,p));
 482              }
 483            }
 484            if (inRange(jetpt/GeV, 24., 40.)) {
 485              _numjets04_10_15_24_40->fill();
 486              _h_N04_10_15_24_40->fill(j.particles().size());
 487              for (const Particle& p : j.particles()) {
 488                _h_z04_10_15_24_40->fill(calcz(j,p));
 489                _h_ptrel04_10_15_24_40->fill(calcptrel(j,p));
 490                _h_rdA04_10_15_24_40->fill(calcr(j,p),calcrweight(j,p));
 491              }
 492            }
 493          }
 494          if (j.absrap() > 1.5 && j.absrap() < 1.9) {
 495            _h_pt04_15_19->fill(jetpt/GeV);
 496            if (inRange(jetpt/GeV, 4., 6.)) {
 497              _numjets04_15_19_04_06->fill();
 498              _h_N04_15_19_04_06->fill(j.particles().size());
 499              for (const Particle& p : j.particles()) {
 500                _h_z04_15_19_04_06->fill(calcz(j,p));
 501                _h_ptrel04_15_19_04_06->fill(calcptrel(j,p));
 502                _h_rdA04_15_19_04_06->fill(calcr(j,p),calcrweight(j,p));
 503              }
 504            }
 505            if (inRange(jetpt/GeV, 6., 10.)) {
 506              _numjets04_15_19_06_10->fill();
 507              _h_N04_15_19_06_10->fill(j.particles().size());
 508              for (const Particle& p : j.particles()) {
 509                _h_z04_15_19_06_10->fill(calcz(j,p));
 510                _h_ptrel04_15_19_06_10->fill(calcptrel(j,p));
 511                _h_rdA04_15_19_06_10->fill(calcr(j,p),calcrweight(j,p));
 512              }
 513            }
 514            if (inRange(jetpt/GeV, 10., 15.)) {
 515              _numjets04_15_19_10_15->fill();
 516              _h_N04_15_19_10_15->fill(j.particles().size());
 517              for (const Particle& p : j.particles()) {
 518                _h_z04_15_19_10_15->fill(calcz(j,p));
 519                _h_ptrel04_15_19_10_15->fill(calcptrel(j,p));
 520                _h_rdA04_15_19_10_15->fill(calcr(j,p),calcrweight(j,p));
 521              }
 522            }
 523            if (inRange(jetpt/GeV, 15., 24.)) {
 524              _numjets04_15_19_15_24->fill();
 525              _h_N04_15_19_15_24->fill(j.particles().size());
 526              for (const Particle& p : j.particles()) {
 527                _h_z04_15_19_15_24->fill(calcz(j,p));
 528                _h_ptrel04_15_19_15_24->fill(calcptrel(j,p));
 529                _h_rdA04_15_19_15_24->fill(calcr(j,p),calcrweight(j,p));
 530              }
 531            }
 532            if (inRange(jetpt/GeV, 24., 40.)) {
 533              _numjets04_15_19_24_40->fill();
 534              _h_N04_15_19_24_40->fill(j.particles().size());
 535              for (const Particle& p : j.particles()) {
 536                _h_z04_15_19_24_40->fill(calcz(j,p));
 537                _h_ptrel04_15_19_24_40->fill(calcptrel(j,p));
 538                _h_rdA04_15_19_24_40->fill(calcr(j,p),calcrweight(j,p));
 539              }
 540            }
 541          } // 1.5 < rapidity < 1.9
 542          if (j.absrap() < 1.9) {
 543            if (inRange(jetpt/GeV, 4., 6.)) {
 544              _numjets04_00_19_04_06->fill();
 545              _h_N04_00_19_04_06->fill(j.particles().size());
 546              for (const Particle& p : j.particles()) {
 547                _h_z04_00_19_04_06->fill(calcz(j,p));
 548                _h_ptrel04_00_19_04_06->fill(calcptrel(j,p));
 549                _h_rdA04_00_19_04_06->fill(calcr(j,p),calcrweight(j,p));
 550              }
 551            }
 552            if (inRange(jetpt/GeV, 6., 10.)) {
 553              _numjets04_00_19_06_10->fill();
 554              _h_N04_00_19_06_10->fill(j.particles().size());
 555              for (const Particle& p : j.particles()) {
 556                _h_z04_00_19_06_10->fill(calcz(j,p));
 557                _h_ptrel04_00_19_06_10->fill(calcptrel(j,p));
 558                _h_rdA04_00_19_06_10->fill(calcr(j,p),calcrweight(j,p));
 559              }
 560            }
 561            if (inRange(jetpt/GeV, 10., 15.)) {
 562              _numjets04_00_19_10_15->fill();
 563              _h_N04_00_19_10_15->fill(j.particles().size());
 564              for (const Particle& p : j.particles()) {
 565                _h_z04_00_19_10_15->fill(calcz(j,p));
 566                _h_ptrel04_00_19_10_15->fill(calcptrel(j,p));
 567                _h_rdA04_00_19_10_15->fill(calcr(j,p),calcrweight(j,p));
 568              }
 569            }
 570            if (inRange(jetpt/GeV, 15., 24.)) {
 571              _numjets04_00_19_15_24->fill();
 572              _h_N04_00_19_15_24->fill(j.particles().size());
 573              for (const Particle& p : j.particles()) {
 574                _h_z04_00_19_15_24->fill(calcz(j,p));
 575                _h_ptrel04_00_19_15_24->fill(calcptrel(j,p));
 576                _h_rdA04_00_19_15_24->fill(calcr(j,p),calcrweight(j,p));
 577              }
 578            }
 579            if (inRange(jetpt/GeV, 24., 40.)) {
 580              _numjets04_00_19_24_40->fill();
 581              _h_N04_00_19_24_40->fill(j.particles().size());
 582              for (const Particle& p : j.particles()) {
 583                _h_z04_00_19_24_40->fill(calcz(j,p));
 584                _h_ptrel04_00_19_24_40->fill(calcptrel(j,p));
 585                _h_rdA04_00_19_24_40->fill(calcr(j,p),calcrweight(j,p));
 586              }
 587            }
 588          } // 0.0 < rapidity < 1.9
 589        } // each jet
 590      } // jets04 not empty
 591
 592      const Jets& jets06 = apply<JetAlg>(event, "Jets06").jets();
 593      if (!jets06.empty()) {
 594        _sumofweights06->fill();
 595        for (const Jet& j : jets06) {
 596          const double jetpt = j.pT();
 597          if (j.absrap() < 0.5) {
 598            _h_pt06_00_05->fill(jetpt/GeV);
 599            if (inRange(jetpt/GeV, 4., 6.)) {
 600              _numjets06_00_05_04_06->fill();
 601              _h_N06_00_05_04_06->fill(j.particles().size());
 602              for (const Particle& p : j.particles()) {
 603                _h_z06_00_05_04_06->fill(calcz(j,p));
 604                _h_ptrel06_00_05_04_06->fill(calcptrel(j,p));
 605                _h_rdA06_00_05_04_06->fill(calcr(j,p),calcrweight(j,p));
 606              }
 607            }
 608            if (inRange(jetpt/GeV, 6., 10.)) {
 609              _numjets06_00_05_06_10->fill();
 610              _h_N06_00_05_06_10->fill(j.particles().size());
 611              for (const Particle& p : j.particles()) {
 612                _h_z06_00_05_06_10->fill(calcz(j,p));
 613                _h_ptrel06_00_05_06_10->fill(calcptrel(j,p));
 614                _h_rdA06_00_05_06_10->fill(calcr(j,p),calcrweight(j,p));
 615              }
 616            }
 617            if (inRange(jetpt/GeV, 10., 15.)) {
 618              _numjets06_00_05_10_15->fill();
 619              _h_N06_00_05_10_15->fill(j.particles().size());
 620              for (const Particle& p : j.particles()) {
 621                _h_z06_00_05_10_15->fill(calcz(j,p));
 622                _h_ptrel06_00_05_10_15->fill(calcptrel(j,p));
 623                _h_rdA06_00_05_10_15->fill(calcr(j,p),calcrweight(j,p));
 624              }
 625            }
 626            if (inRange(jetpt/GeV, 15., 24.)) {
 627              _numjets06_00_05_15_24->fill();
 628              _h_N06_00_05_15_24->fill(j.particles().size());
 629              for (const Particle& p : j.particles()) {
 630                _h_z06_00_05_15_24->fill(calcz(j,p));
 631                _h_ptrel06_00_05_15_24->fill(calcptrel(j,p));
 632                _h_rdA06_00_05_15_24->fill(calcr(j,p),calcrweight(j,p));
 633              }
 634            }
 635            if (inRange(jetpt/GeV, 24., 40.)) {
 636              _numjets06_00_05_24_40->fill();
 637              _h_N06_00_05_24_40->fill(j.particles().size());
 638              for (const Particle& p : j.particles()) {
 639                _h_z06_00_05_24_40->fill(calcz(j,p));
 640                _h_ptrel06_00_05_24_40->fill(calcptrel(j,p));
 641                _h_rdA06_00_05_24_40->fill(calcr(j,p),calcrweight(j,p));
 642              }
 643            }
 644          }
 645          if (j.absrap() > 0.5 && j.absrap() < 1.0) {
 646            _h_pt06_05_10->fill(jetpt/GeV);
 647            if (inRange(jetpt/GeV, 4., 6.)) {
 648              _numjets06_05_10_04_06->fill();
 649              _h_N06_05_10_04_06->fill(j.particles().size());
 650              for (const Particle& p : j.particles()) {
 651                _h_z06_05_10_04_06->fill(calcz(j,p));
 652                _h_ptrel06_05_10_04_06->fill(calcptrel(j,p));
 653                _h_rdA06_05_10_04_06->fill(calcr(j,p),calcrweight(j,p));
 654              }
 655            }
 656            if (inRange(jetpt/GeV, 6., 10.)) {
 657              _numjets06_05_10_06_10->fill();
 658              _h_N06_05_10_06_10->fill(j.particles().size());
 659              for (const Particle& p : j.particles()) {
 660                _h_z06_05_10_06_10->fill(calcz(j,p));
 661                _h_ptrel06_05_10_06_10->fill(calcptrel(j,p));
 662                _h_rdA06_05_10_06_10->fill(calcr(j,p),calcrweight(j,p));
 663              }
 664            }
 665            if (inRange(jetpt/GeV, 10., 15.)) {
 666              _numjets06_05_10_10_15->fill();
 667              _h_N06_05_10_10_15->fill(j.particles().size());
 668              for (const Particle& p : j.particles()) {
 669                _h_z06_05_10_10_15->fill(calcz(j,p));
 670                _h_ptrel06_05_10_10_15->fill(calcptrel(j,p));
 671                _h_rdA06_05_10_10_15->fill(calcr(j,p),calcrweight(j,p));
 672              }
 673            }
 674            if (inRange(jetpt/GeV, 15., 24.)) {
 675              _numjets06_05_10_15_24->fill();
 676              _h_N06_05_10_15_24->fill(j.particles().size());
 677              for (const Particle& p : j.particles()) {
 678                _h_z06_05_10_15_24->fill(calcz(j,p));
 679                _h_ptrel06_05_10_15_24->fill(calcptrel(j,p));
 680                _h_rdA06_05_10_15_24->fill(calcr(j,p),calcrweight(j,p));
 681              }
 682            }
 683            if (inRange(jetpt/GeV, 24., 40.)) {
 684              _numjets06_05_10_24_40->fill();
 685              _h_N06_05_10_24_40->fill(j.particles().size());
 686              for (const Particle& p : j.particles()) {
 687                _h_z06_05_10_24_40->fill(calcz(j,p));
 688                _h_ptrel06_05_10_24_40->fill(calcptrel(j,p));
 689                _h_rdA06_05_10_24_40->fill(calcr(j,p),calcrweight(j,p));
 690              }
 691            }
 692          }
 693          if (j.absrap() > 1.0 && j.absrap() < 1.5) {
 694            _h_pt06_10_15->fill(jetpt/GeV);
 695            if (inRange(jetpt/GeV, 4., 6.)) {
 696              _numjets06_10_15_04_06->fill();
 697              _h_N06_10_15_04_06->fill(j.particles().size());
 698              for (const Particle& p : j.particles()) {
 699                _h_z06_10_15_04_06->fill(calcz(j,p));
 700                _h_ptrel06_10_15_04_06->fill(calcptrel(j,p));
 701                _h_rdA06_10_15_04_06->fill(calcr(j,p),calcrweight(j,p));
 702              }
 703            }
 704            if (inRange(jetpt/GeV, 6., 10.)) {
 705              _numjets06_10_15_06_10->fill();
 706              _h_N06_10_15_06_10->fill(j.particles().size());
 707              for (const Particle& p : j.particles()) {
 708                _h_z06_10_15_06_10->fill(calcz(j,p));
 709                _h_ptrel06_10_15_06_10->fill(calcptrel(j,p));
 710                _h_rdA06_10_15_06_10->fill(calcr(j,p),calcrweight(j,p));
 711              }
 712            }
 713            if (inRange(jetpt/GeV, 10., 15.)) {
 714              _numjets06_10_15_10_15->fill();
 715              _h_N06_10_15_10_15->fill(j.particles().size());
 716              for (const Particle& p : j.particles()) {
 717                _h_z06_10_15_10_15->fill(calcz(j,p));
 718                _h_ptrel06_10_15_10_15->fill(calcptrel(j,p));
 719                _h_rdA06_10_15_10_15->fill(calcr(j,p),calcrweight(j,p));
 720              }
 721            }
 722            if (inRange(jetpt/GeV, 15., 24.)) {
 723              _numjets06_10_15_15_24->fill();
 724              _h_N06_10_15_15_24->fill(j.particles().size());
 725              for (const Particle& p : j.particles()) {
 726                _h_z06_10_15_15_24->fill(calcz(j,p));
 727                _h_ptrel06_10_15_15_24->fill(calcptrel(j,p));
 728                _h_rdA06_10_15_15_24->fill(calcr(j,p),calcrweight(j,p));
 729              }
 730            }
 731            if (inRange(jetpt/GeV, 24., 40.)) {
 732              _numjets06_10_15_24_40->fill();
 733              _h_N06_10_15_24_40->fill(j.particles().size());
 734              for (const Particle& p : j.particles()) {
 735                _h_z06_10_15_24_40->fill(calcz(j,p));
 736                _h_ptrel06_10_15_24_40->fill(calcptrel(j,p));
 737                _h_rdA06_10_15_24_40->fill(calcr(j,p),calcrweight(j,p));
 738              }
 739            }
 740          }
 741          if (j.absrap() > 1.5 && j.absrap() < 1.9) {
 742            _h_pt06_15_19->fill(jetpt/GeV);
 743            if (inRange(jetpt/GeV, 4., 6.)) {
 744              _numjets06_15_19_04_06->fill();
 745              _h_N06_15_19_04_06->fill(j.particles().size());
 746              for (const Particle& p : j.particles()) {
 747                _h_z06_15_19_04_06->fill(calcz(j,p));
 748                _h_ptrel06_15_19_04_06->fill(calcptrel(j,p));
 749                _h_rdA06_15_19_04_06->fill(calcr(j,p),calcrweight(j,p));
 750              }
 751            }
 752            if (inRange(jetpt/GeV, 6., 10.)) {
 753              _numjets06_15_19_06_10->fill();
 754              _h_N06_15_19_06_10->fill(j.particles().size());
 755              for (const Particle& p : j.particles()) {
 756                _h_z06_15_19_06_10->fill(calcz(j,p));
 757                _h_ptrel06_15_19_06_10->fill(calcptrel(j,p));
 758                _h_rdA06_15_19_06_10->fill(calcr(j,p),calcrweight(j,p));
 759              }
 760            }
 761            if (inRange(jetpt/GeV, 10., 15.)) {
 762              _numjets06_15_19_10_15->fill();
 763              _h_N06_15_19_10_15->fill(j.particles().size());
 764              for (const Particle& p : j.particles()) {
 765                _h_z06_15_19_10_15->fill(calcz(j,p));
 766                _h_ptrel06_15_19_10_15->fill(calcptrel(j,p));
 767                _h_rdA06_15_19_10_15->fill(calcr(j,p),calcrweight(j,p));
 768              }
 769            }
 770            if (inRange(jetpt/GeV, 15., 24.)) {
 771              _numjets06_15_19_15_24->fill();
 772              _h_N06_15_19_15_24->fill(j.particles().size());
 773              for (const Particle& p : j.particles()) {
 774                _h_z06_15_19_15_24->fill(calcz(j,p));
 775                _h_ptrel06_15_19_15_24->fill(calcptrel(j,p));
 776                _h_rdA06_15_19_15_24->fill(calcr(j,p),calcrweight(j,p));
 777              }
 778            }
 779            if (inRange(jetpt/GeV, 24., 40.)) {
 780              _numjets06_15_19_24_40->fill();
 781              _h_N06_15_19_24_40->fill(j.particles().size());
 782              for (const Particle& p : j.particles()) {
 783                _h_z06_15_19_24_40->fill(calcz(j,p));
 784                _h_ptrel06_15_19_24_40->fill(calcptrel(j,p));
 785                _h_rdA06_15_19_24_40->fill(calcr(j,p),calcrweight(j,p));
 786              }
 787            }
 788          } // 1.5 < rapidity < 1.9
 789          if (j.absrap() < 1.9) {
 790            if (inRange(jetpt/GeV, 4., 6.)) {
 791              _numjets06_00_19_04_06->fill();
 792              _h_N06_00_19_04_06->fill(j.particles().size());
 793              for (const Particle& p : j.particles()) {
 794                _h_z06_00_19_04_06->fill(calcz(j,p));
 795                _h_ptrel06_00_19_04_06->fill(calcptrel(j,p));
 796                _h_rdA06_00_19_04_06->fill(calcr(j,p),calcrweight(j,p));
 797              }
 798            }
 799            if (inRange(jetpt/GeV, 6., 10.)) {
 800              _numjets06_00_19_06_10->fill();
 801              _h_N06_00_19_06_10->fill(j.particles().size());
 802              for (const Particle& p : j.particles()) {
 803                _h_z06_00_19_06_10->fill(calcz(j,p));
 804                _h_ptrel06_00_19_06_10->fill(calcptrel(j,p));
 805                _h_rdA06_00_19_06_10->fill(calcr(j,p),calcrweight(j,p));
 806              }
 807            }
 808            if (inRange(jetpt/GeV, 10., 15.)) {
 809              _numjets06_00_19_10_15->fill();
 810              _h_N06_00_19_10_15->fill(j.particles().size());
 811              for (const Particle& p : j.particles()) {
 812                _h_z06_00_19_10_15->fill(calcz(j,p));
 813                _h_ptrel06_00_19_10_15->fill(calcptrel(j,p));
 814                _h_rdA06_00_19_10_15->fill(calcr(j,p),calcrweight(j,p));
 815              }
 816            }
 817            if (inRange(jetpt/GeV, 15., 24.)) {
 818              _numjets06_00_19_15_24->fill();
 819              _h_N06_00_19_15_24->fill(j.particles().size());
 820              for (const Particle& p : j.particles()) {
 821                _h_z06_00_19_15_24->fill(calcz(j,p));
 822                _h_ptrel06_00_19_15_24->fill(calcptrel(j,p));
 823                _h_rdA06_00_19_15_24->fill(calcr(j,p),calcrweight(j,p));
 824              }
 825            }
 826            if (inRange(jetpt/GeV, 24., 40.)) {
 827              _numjets06_00_19_24_40->fill();
 828              _h_N06_00_19_24_40->fill(j.particles().size());
 829              for (const Particle& p : j.particles()) {
 830                _h_z06_00_19_24_40->fill(calcz(j,p));
 831                _h_ptrel06_00_19_24_40->fill(calcptrel(j,p));
 832                _h_rdA06_00_19_24_40->fill(calcr(j,p),calcrweight(j,p));
 833              }
 834            }
 835          }
 836        } // each jet
 837      } // jets06 not empty
 838
 839    } // end of event
 840
 841
 842    /// Normalise histograms etc., after the run
 843    void finalize() {
 844
 845      // pT histos: d2sigma_jet / deta dpT
 846      const double xsec = crossSection()/microbarn;
 847      safeinvscale(_h_pt04_00_05, _sumofweights04->val()*(2*0.5)/xsec);
 848      safeinvscale(_h_pt06_00_05, _sumofweights06->val()*(2*0.5)/xsec);
 849      safeinvscale(_h_pt04_05_10, _sumofweights04->val()*(2*0.5)/xsec);
 850      safeinvscale(_h_pt06_05_10, _sumofweights06->val()*(2*0.5)/xsec);
 851      safeinvscale(_h_pt04_10_15, _sumofweights04->val()*(2*0.5)/xsec);
 852      safeinvscale(_h_pt06_10_15, _sumofweights06->val()*(2*0.5)/xsec);
 853      safeinvscale(_h_pt04_15_19, _sumofweights04->val()*(2*0.4)/xsec);
 854      safeinvscale(_h_pt06_15_19, _sumofweights06->val()*(2*0.4)/xsec);
 855
 856      // N histos: 1/N_jet dN_jet / dN^{ch}_jet
 857      safeinvscale(_h_N04_00_05_04_06, _numjets04_00_05_04_06->val());
 858      safeinvscale(_h_N06_00_05_04_06, _numjets06_00_05_04_06->val());
 859      safeinvscale(_h_N04_00_05_06_10, _numjets04_00_05_06_10->val());
 860      safeinvscale(_h_N06_00_05_06_10, _numjets06_00_05_06_10->val());
 861      safeinvscale(_h_N04_00_05_10_15, _numjets04_00_05_10_15->val());
 862      safeinvscale(_h_N06_00_05_10_15, _numjets06_00_05_10_15->val());
 863      safeinvscale(_h_N04_00_05_15_24, _numjets04_00_05_15_24->val());
 864      safeinvscale(_h_N06_00_05_15_24, _numjets06_00_05_15_24->val());
 865      safeinvscale(_h_N04_00_05_24_40, _numjets04_00_05_24_40->val());
 866      safeinvscale(_h_N06_00_05_24_40, _numjets06_00_05_24_40->val());
 867      safeinvscale(_h_N04_05_10_04_06, _numjets04_05_10_04_06->val());
 868      safeinvscale(_h_N06_05_10_04_06, _numjets06_05_10_04_06->val());
 869      safeinvscale(_h_N04_05_10_06_10, _numjets04_05_10_06_10->val());
 870      safeinvscale(_h_N06_05_10_06_10, _numjets06_05_10_06_10->val());
 871      safeinvscale(_h_N04_05_10_10_15, _numjets04_05_10_10_15->val());
 872      safeinvscale(_h_N06_05_10_10_15, _numjets06_05_10_10_15->val());
 873      safeinvscale(_h_N04_05_10_15_24, _numjets04_05_10_15_24->val());
 874      safeinvscale(_h_N06_05_10_15_24, _numjets06_05_10_15_24->val());
 875      safeinvscale(_h_N04_05_10_24_40, _numjets04_05_10_24_40->val());
 876      safeinvscale(_h_N06_05_10_24_40, _numjets06_05_10_24_40->val());
 877      safeinvscale(_h_N04_10_15_04_06, _numjets04_10_15_04_06->val());
 878      safeinvscale(_h_N06_10_15_04_06, _numjets06_10_15_04_06->val());
 879      safeinvscale(_h_N04_10_15_06_10, _numjets04_10_15_06_10->val());
 880      safeinvscale(_h_N06_10_15_06_10, _numjets06_10_15_06_10->val());
 881      safeinvscale(_h_N04_10_15_10_15, _numjets04_10_15_10_15->val());
 882      safeinvscale(_h_N06_10_15_10_15, _numjets06_10_15_10_15->val());
 883      safeinvscale(_h_N04_10_15_15_24, _numjets04_10_15_15_24->val());
 884      safeinvscale(_h_N06_10_15_15_24, _numjets06_10_15_15_24->val());
 885      safeinvscale(_h_N04_10_15_24_40, _numjets04_10_15_24_40->val());
 886      safeinvscale(_h_N06_10_15_24_40, _numjets06_10_15_24_40->val());
 887      safeinvscale(_h_N04_15_19_04_06, _numjets04_15_19_04_06->val());
 888      safeinvscale(_h_N06_15_19_04_06, _numjets06_15_19_04_06->val());
 889      safeinvscale(_h_N04_15_19_06_10, _numjets04_15_19_06_10->val());
 890      safeinvscale(_h_N06_15_19_06_10, _numjets06_15_19_06_10->val());
 891      safeinvscale(_h_N04_15_19_10_15, _numjets04_15_19_10_15->val());
 892      safeinvscale(_h_N06_15_19_10_15, _numjets06_15_19_10_15->val());
 893      safeinvscale(_h_N04_15_19_15_24, _numjets04_15_19_15_24->val());
 894      safeinvscale(_h_N06_15_19_15_24, _numjets06_15_19_15_24->val());
 895      safeinvscale(_h_N04_15_19_24_40, _numjets04_15_19_24_40->val());
 896      safeinvscale(_h_N06_15_19_24_40, _numjets06_15_19_24_40->val());
 897      safeinvscale(_h_N04_00_19_04_06, _numjets04_00_19_04_06->val());
 898      safeinvscale(_h_N06_00_19_04_06, _numjets06_00_19_04_06->val());
 899      safeinvscale(_h_N04_00_19_06_10, _numjets04_00_19_06_10->val());
 900      safeinvscale(_h_N06_00_19_06_10, _numjets06_00_19_06_10->val());
 901      safeinvscale(_h_N04_00_19_10_15, _numjets04_00_19_10_15->val());
 902      safeinvscale(_h_N06_00_19_10_15, _numjets06_00_19_10_15->val());
 903      safeinvscale(_h_N04_00_19_15_24, _numjets04_00_19_15_24->val());
 904      safeinvscale(_h_N06_00_19_15_24, _numjets06_00_19_15_24->val());
 905      safeinvscale(_h_N04_00_19_24_40, _numjets04_00_19_24_40->val());
 906      safeinvscale(_h_N06_00_19_24_40, _numjets06_00_19_24_40->val());
 907      
 908      // z histos: 1/N_jet dN_track / dz_track->val()
 909      safeinvscale(_h_z04_00_05_04_06, _numjets04_00_05_04_06->val());
 910      safeinvscale(_h_z06_00_05_04_06, _numjets06_00_05_04_06->val());
 911      safeinvscale(_h_z04_00_05_06_10, _numjets04_00_05_06_10->val());
 912      safeinvscale(_h_z06_00_05_06_10, _numjets06_00_05_06_10->val());
 913      safeinvscale(_h_z04_00_05_10_15, _numjets04_00_05_10_15->val());
 914      safeinvscale(_h_z06_00_05_10_15, _numjets06_00_05_10_15->val());
 915      safeinvscale(_h_z04_00_05_15_24, _numjets04_00_05_15_24->val());
 916      safeinvscale(_h_z06_00_05_15_24, _numjets06_00_05_15_24->val());
 917      safeinvscale(_h_z04_00_05_24_40, _numjets04_00_05_24_40->val());
 918      safeinvscale(_h_z06_00_05_24_40, _numjets06_00_05_24_40->val());
 919      safeinvscale(_h_z04_05_10_04_06, _numjets04_05_10_04_06->val());
 920      safeinvscale(_h_z06_05_10_04_06, _numjets06_05_10_04_06->val());
 921      safeinvscale(_h_z04_05_10_06_10, _numjets04_05_10_06_10->val());
 922      safeinvscale(_h_z06_05_10_06_10, _numjets06_05_10_06_10->val());
 923      safeinvscale(_h_z04_05_10_10_15, _numjets04_05_10_10_15->val());
 924      safeinvscale(_h_z06_05_10_10_15, _numjets06_05_10_10_15->val());
 925      safeinvscale(_h_z04_05_10_15_24, _numjets04_05_10_15_24->val());
 926      safeinvscale(_h_z06_05_10_15_24, _numjets06_05_10_15_24->val());
 927      safeinvscale(_h_z04_05_10_24_40, _numjets04_05_10_24_40->val());
 928      safeinvscale(_h_z06_05_10_24_40, _numjets06_05_10_24_40->val());
 929      safeinvscale(_h_z04_10_15_04_06, _numjets04_10_15_04_06->val());
 930      safeinvscale(_h_z06_10_15_04_06, _numjets06_10_15_04_06->val());
 931      safeinvscale(_h_z04_10_15_06_10, _numjets04_10_15_06_10->val());
 932      safeinvscale(_h_z06_10_15_06_10, _numjets06_10_15_06_10->val());
 933      safeinvscale(_h_z04_10_15_10_15, _numjets04_10_15_10_15->val());
 934      safeinvscale(_h_z06_10_15_10_15, _numjets06_10_15_10_15->val());
 935      safeinvscale(_h_z04_10_15_15_24, _numjets04_10_15_15_24->val());
 936      safeinvscale(_h_z06_10_15_15_24, _numjets06_10_15_15_24->val());
 937      safeinvscale(_h_z04_10_15_24_40, _numjets04_10_15_24_40->val());
 938      safeinvscale(_h_z06_10_15_24_40, _numjets06_10_15_24_40->val());
 939      safeinvscale(_h_z04_15_19_04_06, _numjets04_15_19_04_06->val());
 940      safeinvscale(_h_z06_15_19_04_06, _numjets06_15_19_04_06->val());
 941      safeinvscale(_h_z04_15_19_06_10, _numjets04_15_19_06_10->val());
 942      safeinvscale(_h_z06_15_19_06_10, _numjets06_15_19_06_10->val());
 943      safeinvscale(_h_z04_15_19_10_15, _numjets04_15_19_10_15->val());
 944      safeinvscale(_h_z06_15_19_10_15, _numjets06_15_19_10_15->val());
 945      safeinvscale(_h_z04_15_19_15_24, _numjets04_15_19_15_24->val());
 946      safeinvscale(_h_z06_15_19_15_24, _numjets06_15_19_15_24->val());
 947      safeinvscale(_h_z04_15_19_24_40, _numjets04_15_19_24_40->val());
 948      safeinvscale(_h_z06_15_19_24_40, _numjets06_15_19_24_40->val());
 949      safeinvscale(_h_z04_00_19_04_06, _numjets04_00_19_04_06->val());
 950      safeinvscale(_h_z06_00_19_04_06, _numjets06_00_19_04_06->val());
 951      safeinvscale(_h_z04_00_19_06_10, _numjets04_00_19_06_10->val());
 952      safeinvscale(_h_z06_00_19_06_10, _numjets06_00_19_06_10->val());
 953      safeinvscale(_h_z04_00_19_10_15, _numjets04_00_19_10_15->val());
 954      safeinvscale(_h_z06_00_19_10_15, _numjets06_00_19_10_15->val());
 955      safeinvscale(_h_z04_00_19_15_24, _numjets04_00_19_15_24->val());
 956      safeinvscale(_h_z06_00_19_15_24, _numjets06_00_19_15_24->val());
 957      safeinvscale(_h_z04_00_19_24_40, _numjets04_00_19_24_40->val());
 958      safeinvscale(_h_z06_00_19_24_40, _numjets06_00_19_24_40->val());
 959
 960      // ptrel histos: 1/N_jet dN_track / dptrel_track
 961      safeinvscale(_h_ptrel04_00_05_04_06, _numjets04_00_05_04_06->val());
 962      safeinvscale(_h_ptrel06_00_05_04_06, _numjets06_00_05_04_06->val());
 963      safeinvscale(_h_ptrel04_00_05_06_10, _numjets04_00_05_06_10->val());
 964      safeinvscale(_h_ptrel06_00_05_06_10, _numjets06_00_05_06_10->val());
 965      safeinvscale(_h_ptrel04_00_05_10_15, _numjets04_00_05_10_15->val());
 966      safeinvscale(_h_ptrel06_00_05_10_15, _numjets06_00_05_10_15->val());
 967      safeinvscale(_h_ptrel04_00_05_15_24, _numjets04_00_05_15_24->val());
 968      safeinvscale(_h_ptrel06_00_05_15_24, _numjets06_00_05_15_24->val());
 969      safeinvscale(_h_ptrel04_00_05_24_40, _numjets04_00_05_24_40->val());
 970      safeinvscale(_h_ptrel06_00_05_24_40, _numjets06_00_05_24_40->val());
 971      safeinvscale(_h_ptrel04_05_10_04_06, _numjets04_05_10_04_06->val());
 972      safeinvscale(_h_ptrel06_05_10_04_06, _numjets06_05_10_04_06->val());
 973      safeinvscale(_h_ptrel04_05_10_06_10, _numjets04_05_10_06_10->val());
 974      safeinvscale(_h_ptrel06_05_10_06_10, _numjets06_05_10_06_10->val());
 975      safeinvscale(_h_ptrel04_05_10_10_15, _numjets04_05_10_10_15->val());
 976      safeinvscale(_h_ptrel06_05_10_10_15, _numjets06_05_10_10_15->val());
 977      safeinvscale(_h_ptrel04_05_10_15_24, _numjets04_05_10_15_24->val());
 978      safeinvscale(_h_ptrel06_05_10_15_24, _numjets06_05_10_15_24->val());
 979      safeinvscale(_h_ptrel04_05_10_24_40, _numjets04_05_10_24_40->val());
 980      safeinvscale(_h_ptrel06_05_10_24_40, _numjets06_05_10_24_40->val());
 981      safeinvscale(_h_ptrel04_10_15_04_06, _numjets04_10_15_04_06->val());
 982      safeinvscale(_h_ptrel06_10_15_04_06, _numjets06_10_15_04_06->val());
 983      safeinvscale(_h_ptrel04_10_15_06_10, _numjets04_10_15_06_10->val());
 984      safeinvscale(_h_ptrel06_10_15_06_10, _numjets06_10_15_06_10->val());
 985      safeinvscale(_h_ptrel04_10_15_10_15, _numjets04_10_15_10_15->val());
 986      safeinvscale(_h_ptrel06_10_15_10_15, _numjets06_10_15_10_15->val());
 987      safeinvscale(_h_ptrel04_10_15_15_24, _numjets04_10_15_15_24->val());
 988      safeinvscale(_h_ptrel06_10_15_15_24, _numjets06_10_15_15_24->val());
 989      safeinvscale(_h_ptrel04_10_15_24_40, _numjets04_10_15_24_40->val());
 990      safeinvscale(_h_ptrel06_10_15_24_40, _numjets06_10_15_24_40->val());
 991      safeinvscale(_h_ptrel04_15_19_04_06, _numjets04_15_19_04_06->val());
 992      safeinvscale(_h_ptrel06_15_19_04_06, _numjets06_15_19_04_06->val());
 993      safeinvscale(_h_ptrel04_15_19_06_10, _numjets04_15_19_06_10->val());
 994      safeinvscale(_h_ptrel06_15_19_06_10, _numjets06_15_19_06_10->val());
 995      safeinvscale(_h_ptrel04_15_19_10_15, _numjets04_15_19_10_15->val());
 996      safeinvscale(_h_ptrel06_15_19_10_15, _numjets06_15_19_10_15->val());
 997      safeinvscale(_h_ptrel04_15_19_15_24, _numjets04_15_19_15_24->val());
 998      safeinvscale(_h_ptrel06_15_19_15_24, _numjets06_15_19_15_24->val());
 999      safeinvscale(_h_ptrel04_15_19_24_40, _numjets04_15_19_24_40->val());
1000      safeinvscale(_h_ptrel06_15_19_24_40, _numjets06_15_19_24_40->val());
1001
1002      safeinvscale(_h_ptrel04_00_19_04_06, _numjets04_00_19_04_06->val());
1003      safeinvscale(_h_ptrel06_00_19_04_06, _numjets06_00_19_04_06->val());
1004      safeinvscale(_h_ptrel04_00_19_06_10, _numjets04_00_19_06_10->val());
1005      safeinvscale(_h_ptrel06_00_19_06_10, _numjets06_00_19_06_10->val());
1006      safeinvscale(_h_ptrel04_00_19_10_15, _numjets04_00_19_10_15->val());
1007      safeinvscale(_h_ptrel06_00_19_10_15, _numjets06_00_19_10_15->val());
1008      safeinvscale(_h_ptrel04_00_19_15_24, _numjets04_00_19_15_24->val());
1009      safeinvscale(_h_ptrel06_00_19_15_24, _numjets06_00_19_15_24->val());
1010      safeinvscale(_h_ptrel04_00_19_24_40, _numjets04_00_19_24_40->val());
1011      safeinvscale(_h_ptrel06_00_19_24_40, _numjets06_00_19_24_40->val());
1012
1013      // r histos: 1/N_jet dN_track / dA
1014      safeinvscale(_h_rdA04_00_05_04_06, _numjets04_00_05_04_06->val());
1015      safeinvscale(_h_rdA06_00_05_04_06, _numjets06_00_05_04_06->val());
1016      safeinvscale(_h_rdA04_00_05_06_10, _numjets04_00_05_06_10->val());
1017      safeinvscale(_h_rdA06_00_05_06_10, _numjets06_00_05_06_10->val());
1018      safeinvscale(_h_rdA04_00_05_10_15, _numjets04_00_05_10_15->val());
1019      safeinvscale(_h_rdA06_00_05_10_15, _numjets06_00_05_10_15->val());
1020      safeinvscale(_h_rdA04_00_05_15_24, _numjets04_00_05_15_24->val());
1021      safeinvscale(_h_rdA06_00_05_15_24, _numjets06_00_05_15_24->val());
1022      safeinvscale(_h_rdA04_00_05_24_40, _numjets04_00_05_24_40->val());
1023      safeinvscale(_h_rdA06_00_05_24_40, _numjets06_00_05_24_40->val());
1024      safeinvscale(_h_rdA04_05_10_04_06, _numjets04_05_10_04_06->val());
1025      safeinvscale(_h_rdA06_05_10_04_06, _numjets06_05_10_04_06->val());
1026      safeinvscale(_h_rdA04_05_10_06_10, _numjets04_05_10_06_10->val());
1027      safeinvscale(_h_rdA06_05_10_06_10, _numjets06_05_10_06_10->val());
1028      safeinvscale(_h_rdA04_05_10_10_15, _numjets04_05_10_10_15->val());
1029      safeinvscale(_h_rdA06_05_10_10_15, _numjets06_05_10_10_15->val());
1030      safeinvscale(_h_rdA04_05_10_15_24, _numjets04_05_10_15_24->val());
1031      safeinvscale(_h_rdA06_05_10_15_24, _numjets06_05_10_15_24->val());
1032      safeinvscale(_h_rdA04_05_10_24_40, _numjets04_05_10_24_40->val());
1033      safeinvscale(_h_rdA06_05_10_24_40, _numjets06_05_10_24_40->val());
1034      safeinvscale(_h_rdA04_10_15_04_06, _numjets04_10_15_04_06->val());
1035      safeinvscale(_h_rdA06_10_15_04_06, _numjets06_10_15_04_06->val());
1036      safeinvscale(_h_rdA04_10_15_06_10, _numjets04_10_15_06_10->val());
1037      safeinvscale(_h_rdA06_10_15_06_10, _numjets06_10_15_06_10->val());
1038      safeinvscale(_h_rdA04_10_15_10_15, _numjets04_10_15_10_15->val());
1039      safeinvscale(_h_rdA06_10_15_10_15, _numjets06_10_15_10_15->val());
1040      safeinvscale(_h_rdA04_10_15_15_24, _numjets04_10_15_15_24->val());
1041      safeinvscale(_h_rdA06_10_15_15_24, _numjets06_10_15_15_24->val());
1042      safeinvscale(_h_rdA04_10_15_24_40, _numjets04_10_15_24_40->val());
1043      safeinvscale(_h_rdA06_10_15_24_40, _numjets06_10_15_24_40->val());
1044      safeinvscale(_h_rdA04_15_19_04_06, _numjets04_15_19_04_06->val());
1045      safeinvscale(_h_rdA06_15_19_04_06, _numjets06_15_19_04_06->val());
1046      safeinvscale(_h_rdA04_15_19_06_10, _numjets04_15_19_06_10->val());
1047      safeinvscale(_h_rdA06_15_19_06_10, _numjets06_15_19_06_10->val());
1048      safeinvscale(_h_rdA04_15_19_10_15, _numjets04_15_19_10_15->val());
1049      safeinvscale(_h_rdA06_15_19_10_15, _numjets06_15_19_10_15->val());
1050      safeinvscale(_h_rdA04_15_19_15_24, _numjets04_15_19_15_24->val());
1051      safeinvscale(_h_rdA06_15_19_15_24, _numjets06_15_19_15_24->val());
1052      safeinvscale(_h_rdA04_15_19_24_40, _numjets04_15_19_24_40->val());
1053      safeinvscale(_h_rdA06_15_19_24_40, _numjets06_15_19_24_40->val());
1054
1055      safeinvscale(_h_rdA04_00_19_04_06, _numjets04_00_19_04_06->val());
1056      safeinvscale(_h_rdA06_00_19_04_06, _numjets06_00_19_04_06->val());
1057      safeinvscale(_h_rdA04_00_19_06_10, _numjets04_00_19_06_10->val());
1058      safeinvscale(_h_rdA06_00_19_06_10, _numjets06_00_19_06_10->val());
1059      safeinvscale(_h_rdA04_00_19_10_15, _numjets04_00_19_10_15->val());
1060      safeinvscale(_h_rdA06_00_19_10_15, _numjets06_00_19_10_15->val());
1061      safeinvscale(_h_rdA04_00_19_15_24, _numjets04_00_19_15_24->val());
1062      safeinvscale(_h_rdA06_00_19_15_24, _numjets06_00_19_15_24->val());
1063      safeinvscale(_h_rdA04_00_19_24_40, _numjets04_00_19_24_40->val());
1064      safeinvscale(_h_rdA06_00_19_24_40, _numjets06_00_19_24_40->val());
1065    }
1066
1067    //@}
1068
1069
1070  private:
1071
1072    void safeinvscale(Histo1DPtr h, double denom) {
1073      if (denom != 0) {
1074        scale(h, 1.0/denom);
1075      } else {
1076        normalize(h, 0);
1077      }
1078    }
1079
1080
1081    /// Event weights
1082    CounterPtr _sumofweights04, _sumofweights06;
1083
1084
1085    /// Jet counters
1086    CounterPtr _numjets04_00_05_04_06, _numjets04_00_05_06_10, _numjets04_00_05_10_15, _numjets04_00_05_15_24, _numjets04_00_05_24_40;
1087    CounterPtr _numjets06_00_05_04_06, _numjets06_00_05_06_10, _numjets06_00_05_10_15, _numjets06_00_05_15_24, _numjets06_00_05_24_40;
1088    CounterPtr _numjets04_05_10_04_06, _numjets04_05_10_06_10, _numjets04_05_10_10_15, _numjets04_05_10_15_24, _numjets04_05_10_24_40;
1089    CounterPtr _numjets06_05_10_04_06, _numjets06_05_10_06_10, _numjets06_05_10_10_15, _numjets06_05_10_15_24, _numjets06_05_10_24_40;
1090    CounterPtr _numjets04_10_15_04_06, _numjets04_10_15_06_10, _numjets04_10_15_10_15, _numjets04_10_15_15_24, _numjets04_10_15_24_40;
1091    CounterPtr _numjets06_10_15_04_06, _numjets06_10_15_06_10, _numjets06_10_15_10_15, _numjets06_10_15_15_24, _numjets06_10_15_24_40;
1092    CounterPtr _numjets04_15_19_04_06, _numjets04_15_19_06_10, _numjets04_15_19_10_15, _numjets04_15_19_15_24, _numjets04_15_19_24_40;
1093    CounterPtr _numjets06_15_19_04_06, _numjets06_15_19_06_10, _numjets06_15_19_10_15, _numjets06_15_19_15_24, _numjets06_15_19_24_40;
1094    CounterPtr _numjets04_00_19_04_06, _numjets04_00_19_06_10, _numjets04_00_19_10_15, _numjets04_00_19_15_24, _numjets04_00_19_24_40;
1095    CounterPtr _numjets06_00_19_04_06, _numjets06_00_19_06_10, _numjets06_00_19_10_15, _numjets06_00_19_15_24, _numjets06_00_19_24_40;
1096
1097
1098  private:
1099
1100    /// @name Histograms
1101    //@{
1102
1103    Histo1DPtr _h_pt04_00_05, _h_pt06_00_05;
1104    Histo1DPtr _h_N04_00_05_04_06, _h_N06_00_05_04_06;
1105    Histo1DPtr _h_N04_00_05_06_10, _h_N06_00_05_06_10;
1106    Histo1DPtr _h_N04_00_05_10_15, _h_N06_00_05_10_15;
1107    Histo1DPtr _h_N04_00_05_15_24, _h_N06_00_05_15_24;
1108    Histo1DPtr _h_N04_00_05_24_40, _h_N06_00_05_24_40;
1109    Histo1DPtr _h_z04_00_05_04_06, _h_z06_00_05_04_06;
1110    Histo1DPtr _h_z04_00_05_06_10, _h_z06_00_05_06_10;
1111    Histo1DPtr _h_z04_00_05_10_15, _h_z06_00_05_10_15;
1112    Histo1DPtr _h_z04_00_05_15_24, _h_z06_00_05_15_24;
1113    Histo1DPtr _h_z04_00_05_24_40, _h_z06_00_05_24_40;
1114    Histo1DPtr _h_ptrel04_00_05_04_06, _h_ptrel06_00_05_04_06;
1115    Histo1DPtr _h_ptrel04_00_05_06_10, _h_ptrel06_00_05_06_10;
1116    Histo1DPtr _h_ptrel04_00_05_10_15, _h_ptrel06_00_05_10_15;
1117    Histo1DPtr _h_ptrel04_00_05_15_24, _h_ptrel06_00_05_15_24;
1118    Histo1DPtr _h_ptrel04_00_05_24_40, _h_ptrel06_00_05_24_40;
1119    Histo1DPtr _h_rdA04_00_05_04_06, _h_rdA06_00_05_04_06;
1120    Histo1DPtr _h_rdA04_00_05_06_10, _h_rdA06_00_05_06_10;
1121    Histo1DPtr _h_rdA04_00_05_10_15, _h_rdA06_00_05_10_15;
1122    Histo1DPtr _h_rdA04_00_05_15_24, _h_rdA06_00_05_15_24;
1123    Histo1DPtr _h_rdA04_00_05_24_40, _h_rdA06_00_05_24_40;
1124
1125    Histo1DPtr _h_pt04_05_10, _h_pt06_05_10;
1126    Histo1DPtr _h_N04_05_10_04_06, _h_N06_05_10_04_06;
1127    Histo1DPtr _h_N04_05_10_06_10, _h_N06_05_10_06_10;
1128    Histo1DPtr _h_N04_05_10_10_15, _h_N06_05_10_10_15;
1129    Histo1DPtr _h_N04_05_10_15_24, _h_N06_05_10_15_24;
1130    Histo1DPtr _h_N04_05_10_24_40, _h_N06_05_10_24_40;
1131    Histo1DPtr _h_z04_05_10_04_06, _h_z06_05_10_04_06;
1132    Histo1DPtr _h_z04_05_10_06_10, _h_z06_05_10_06_10;
1133    Histo1DPtr _h_z04_05_10_10_15, _h_z06_05_10_10_15;
1134    Histo1DPtr _h_z04_05_10_15_24, _h_z06_05_10_15_24;
1135    Histo1DPtr _h_z04_05_10_24_40, _h_z06_05_10_24_40;
1136    Histo1DPtr _h_ptrel04_05_10_04_06, _h_ptrel06_05_10_04_06;
1137    Histo1DPtr _h_ptrel04_05_10_06_10, _h_ptrel06_05_10_06_10;
1138    Histo1DPtr _h_ptrel04_05_10_10_15, _h_ptrel06_05_10_10_15;
1139    Histo1DPtr _h_ptrel04_05_10_15_24, _h_ptrel06_05_10_15_24;
1140    Histo1DPtr _h_ptrel04_05_10_24_40, _h_ptrel06_05_10_24_40;
1141    Histo1DPtr _h_rdA04_05_10_04_06, _h_rdA06_05_10_04_06;
1142    Histo1DPtr _h_rdA04_05_10_06_10, _h_rdA06_05_10_06_10;
1143    Histo1DPtr _h_rdA04_05_10_10_15, _h_rdA06_05_10_10_15;
1144    Histo1DPtr _h_rdA04_05_10_15_24, _h_rdA06_05_10_15_24;
1145    Histo1DPtr _h_rdA04_05_10_24_40, _h_rdA06_05_10_24_40;
1146
1147    Histo1DPtr _h_pt04_10_15, _h_pt06_10_15;
1148    Histo1DPtr _h_N04_10_15_04_06, _h_N06_10_15_04_06;
1149    Histo1DPtr _h_N04_10_15_06_10, _h_N06_10_15_06_10;
1150    Histo1DPtr _h_N04_10_15_10_15, _h_N06_10_15_10_15;
1151    Histo1DPtr _h_N04_10_15_15_24, _h_N06_10_15_15_24;
1152    Histo1DPtr _h_N04_10_15_24_40, _h_N06_10_15_24_40;
1153    Histo1DPtr _h_z04_10_15_04_06, _h_z06_10_15_04_06;
1154    Histo1DPtr _h_z04_10_15_06_10, _h_z06_10_15_06_10;
1155    Histo1DPtr _h_z04_10_15_10_15, _h_z06_10_15_10_15;
1156    Histo1DPtr _h_z04_10_15_15_24, _h_z06_10_15_15_24;
1157    Histo1DPtr _h_z04_10_15_24_40, _h_z06_10_15_24_40;
1158    Histo1DPtr _h_ptrel04_10_15_04_06, _h_ptrel06_10_15_04_06;
1159    Histo1DPtr _h_ptrel04_10_15_06_10, _h_ptrel06_10_15_06_10;
1160    Histo1DPtr _h_ptrel04_10_15_10_15, _h_ptrel06_10_15_10_15;
1161    Histo1DPtr _h_ptrel04_10_15_15_24, _h_ptrel06_10_15_15_24;
1162    Histo1DPtr _h_ptrel04_10_15_24_40, _h_ptrel06_10_15_24_40;
1163    Histo1DPtr _h_rdA04_10_15_04_06, _h_rdA06_10_15_04_06;
1164    Histo1DPtr _h_rdA04_10_15_06_10, _h_rdA06_10_15_06_10;
1165    Histo1DPtr _h_rdA04_10_15_10_15, _h_rdA06_10_15_10_15;
1166    Histo1DPtr _h_rdA04_10_15_15_24, _h_rdA06_10_15_15_24;
1167    Histo1DPtr _h_rdA04_10_15_24_40, _h_rdA06_10_15_24_40;
1168
1169    Histo1DPtr _h_pt04_15_19, _h_pt06_15_19;
1170    Histo1DPtr _h_N04_15_19_04_06, _h_N06_15_19_04_06;
1171    Histo1DPtr _h_N04_15_19_06_10, _h_N06_15_19_06_10;
1172    Histo1DPtr _h_N04_15_19_10_15, _h_N06_15_19_10_15;
1173    Histo1DPtr _h_N04_15_19_15_24, _h_N06_15_19_15_24;
1174    Histo1DPtr _h_N04_15_19_24_40, _h_N06_15_19_24_40;
1175    Histo1DPtr _h_z04_15_19_04_06, _h_z06_15_19_04_06;
1176    Histo1DPtr _h_z04_15_19_06_10, _h_z06_15_19_06_10;
1177    Histo1DPtr _h_z04_15_19_10_15, _h_z06_15_19_10_15;
1178    Histo1DPtr _h_z04_15_19_15_24, _h_z06_15_19_15_24;
1179    Histo1DPtr _h_z04_15_19_24_40, _h_z06_15_19_24_40;
1180    Histo1DPtr _h_ptrel04_15_19_04_06, _h_ptrel06_15_19_04_06;
1181    Histo1DPtr _h_ptrel04_15_19_06_10, _h_ptrel06_15_19_06_10;
1182    Histo1DPtr _h_ptrel04_15_19_10_15, _h_ptrel06_15_19_10_15;
1183    Histo1DPtr _h_ptrel04_15_19_15_24, _h_ptrel06_15_19_15_24;
1184    Histo1DPtr _h_ptrel04_15_19_24_40, _h_ptrel06_15_19_24_40;
1185    Histo1DPtr _h_rdA04_15_19_04_06, _h_rdA06_15_19_04_06;
1186    Histo1DPtr _h_rdA04_15_19_06_10, _h_rdA06_15_19_06_10;
1187    Histo1DPtr _h_rdA04_15_19_10_15, _h_rdA06_15_19_10_15;
1188    Histo1DPtr _h_rdA04_15_19_15_24, _h_rdA06_15_19_15_24;
1189    Histo1DPtr _h_rdA04_15_19_24_40, _h_rdA06_15_19_24_40;
1190
1191    Histo1DPtr _h_N04_00_19_04_06, _h_N06_00_19_04_06;
1192    Histo1DPtr _h_N04_00_19_06_10, _h_N06_00_19_06_10;
1193    Histo1DPtr _h_N04_00_19_10_15, _h_N06_00_19_10_15;
1194    Histo1DPtr _h_N04_00_19_15_24, _h_N06_00_19_15_24;
1195    Histo1DPtr _h_N04_00_19_24_40, _h_N06_00_19_24_40;
1196    Histo1DPtr _h_z04_00_19_04_06, _h_z06_00_19_04_06;
1197    Histo1DPtr _h_z04_00_19_06_10, _h_z06_00_19_06_10;
1198    Histo1DPtr _h_z04_00_19_10_15, _h_z06_00_19_10_15;
1199    Histo1DPtr _h_z04_00_19_15_24, _h_z06_00_19_15_24;
1200    Histo1DPtr _h_z04_00_19_24_40, _h_z06_00_19_24_40;
1201    Histo1DPtr _h_ptrel04_00_19_04_06, _h_ptrel06_00_19_04_06;
1202    Histo1DPtr _h_ptrel04_00_19_06_10, _h_ptrel06_00_19_06_10;
1203    Histo1DPtr _h_ptrel04_00_19_10_15, _h_ptrel06_00_19_10_15;
1204    Histo1DPtr _h_ptrel04_00_19_15_24, _h_ptrel06_00_19_15_24;
1205    Histo1DPtr _h_ptrel04_00_19_24_40, _h_ptrel06_00_19_24_40;
1206    Histo1DPtr _h_rdA04_00_19_04_06, _h_rdA06_00_19_04_06;
1207    Histo1DPtr _h_rdA04_00_19_06_10, _h_rdA06_00_19_06_10;
1208    Histo1DPtr _h_rdA04_00_19_10_15, _h_rdA06_00_19_10_15;
1209    Histo1DPtr _h_rdA04_00_19_15_24, _h_rdA06_00_19_15_24;
1210    Histo1DPtr _h_rdA04_00_19_24_40, _h_rdA06_00_19_24_40;
1211
1212    //@}
1213
1214  };
1215
1216
1217  // This global object acts as a hook for the plugin system
1218  RIVET_DECLARE_PLUGIN(ATLAS_2011_I919017);
1219
1220}