rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2019_I1772062

Soft-drop observables
Experiment: ATLAS (LHC)
Inspire ID: 1772062
Status: VALIDATED
Authors:
  • Jennifer Roloff
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • jet substructure observables for soft drop jets produced at 13 TeV in dijet events

Jet substructure quantities are measured using jets groomed with the soft-drop grooming procedure in dijet events from 32.9 fb$^{-1}$ of $pp$ collisions collected with the ATLAS detector at $\sqrt{s} = 13$ TeV. These observables are sensitive to a wide range of QCD phenomena. Some observables, such as the jet mass and opening angle between the two subjets which pass the soft-drop condition, can be described by a high-order (resummed) series in the strong coupling constant $\alpha_s$. Other observables, such as the momentum sharing between the two subjets, are nearly independent of $\alpha_s$. These observables can be constructed using all interacting particles or using only charged particles reconstructed in the inner tracking detectors. Track-based versions of these observables are not collinear safe, but are measured more precisely, and universal nonperturbative functions can absorb the collinear singularities. The unfolded data are directly compared with QCD calculations and hadron-level Monte Carlo simulations. The measurements are performed in different pseudorapidity regions, which are then used to extract quark and gluon jet shapes using the predicted quark and gluon fractions in each region. All of the parton shower and analytical calculations provide an excellent description of the data in most regions of phase space.

Source code: ATLAS_2019_I1772062.cc
  1#include "Rivet/Analysis.hh"
  2#include "Rivet/Projections/FinalState.hh"
  3#include "Rivet/Projections/FastJets.hh"
  4#include "Rivet/Particle.hh"
  5#include "Rivet/Projections/ChargedFinalState.hh"
  6#include "Rivet/Projections/VetoedFinalState.hh"
  7
  8#include "fastjet/contrib/SoftDrop.hh"
  9#include "fastjet/Selector.hh"
 10
 11
 12namespace Rivet {
 13
 14
 15  /// @brief Soft-drop mass at 13 TeV
 16  class ATLAS_2019_I1772062: public Analysis {
 17  public:
 18
 19    /// Constructor
 20    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2019_I1772062);
 21
 22
 23    void getQuarkGluon(Rivet::Histo1DPtr hForward, Rivet::Histo1DPtr hCentral, Rivet::Histo1DPtr hQuark, Rivet::Histo1DPtr hGluon, int ptbin, string parName, size_t beta) {
 24
 25      int nBins = rhoBins.size() - 1;
 26      if (parName == "rg" || parName == "trg") nBins = rgBins.size()-1;
 27      if ((parName == "zg" || parName == "tzg") && beta == 0) nBins = zgBinsBeta0.size()-1;
 28      if ((parName == "zg" || parName == "tzg") && beta != 0) nBins = zgBins.size()-1;
 29
 30      double FGC = gluonFractionCentral[ptbin];
 31      double FGF = gluonFractionForward[ptbin];
 32      double FQC = 1.-FGC;
 33      double FQF = 1.-FGF;
 34
 35      for (size_t i=0; i<hGluon->numBins(); i++) {
 36        double binCenter = hGluon->bin(i).midpoint();
 37        double gVal = 0., qVal = 0.;
 38        if ((FQF -  FQC) != 0.) {
 39          gVal = (FQF * hCentral->bin(ptbin*(nBins) + i).height() - FQC * hForward->bin(ptbin*(nBins) + i).height()) / (FQF - FQC);
 40          qVal = (FGF * hCentral->bin(ptbin*(nBins) + i).height() - FGC * hForward->bin(ptbin*(nBins) + i).height()) / (FQC - FQF);
 41        }
 42        hGluon->fill(binCenter, gVal);
 43        hQuark->fill(binCenter, qVal);
 44      }
 45
 46      histNorm(hQuark, parName);
 47      histNorm(hGluon, parName);
 48    }
 49
 50    void ptNorm(Rivet::Histo1DPtr ptBinnedHist, std::string var, size_t beta) {
 51      size_t varNormBin1 = 0;
 52      size_t varNormBin2 = 0;
 53      size_t nBins = 10;
 54
 55      if (var=="m" || var=="tm") {
 56        varNormBin1 = normBin1;
 57        varNormBin2 = normBin2;
 58      }
 59      if (var=="zg" || var=="tzg") {
 60        if(beta==0){
 61          varNormBin2 = zgBinsBeta0.size()-1;
 62          nBins = zgBinsBeta0.size()-1;
 63        }
 64        else {
 65          varNormBin2 = zgBins.size()-1;
 66          nBins = zgBins.size()-1;
 67        }
 68      }
 69      if (var=="rg" || var=="trg") {
 70        varNormBin2 = rgBins.size()-1;
 71        nBins = rgBins.size()-1;
 72      }
 73
 74      for (size_t k=0; k< ptBins.size()-1; ++k) {
 75        double normalization = 0;
 76
 77        for (size_t j=varNormBin1; j<varNormBin2; ++j) {
 78          double binWidth = 1.;
 79          if (var=="m" || var=="tm") {
 80            binWidth = rhoBins[j+1] - rhoBins[j];
 81          }
 82          if (var=="zg" || var=="tzg") {
 83            if(beta==0){
 84              binWidth = zgBinsBeta0[j+1]- zgBinsBeta0[j];
 85            }
 86            else {
 87              binWidth = zgBins[j+1]- zgBins[j];
 88            }
 89          }
 90
 91          if (var=="rg" || var=="trg") {
 92            binWidth = rgBins[j+1]- rgBins[j];
 93            if (j==nBins-1) ptBinnedHist->bin( k*nBins+j ).scaleW(2.);
 94            //normalization += ptBinnedHist->bin(k*nBins+j).height()*binWidth;
 95            normalization += ptBinnedHist->bin(k*(nBins) + j).height()*binWidth;
 96          }
 97          else{
 98            normalization += ptBinnedHist->bin(k*(nBins) + j).height()*binWidth;
 99          }
100        }
101
102        if (normalization == 0) continue;
103
104        for (unsigned int j=0; j<nBins; j++) {
105          if (var=="rg" || var=="trg") {
106            ptBinnedHist->bin(k*(nBins) + j ).scaleW(1. / (normalization) );
107          }
108          else{
109            ptBinnedHist->bin(k*(nBins) + j ).scaleW(1. / (normalization));
110          }
111
112        }
113      }
114
115      return;
116    }
117
118    void histNorm(Rivet::Histo1DPtr hist, std::string var) {
119      if (var=="m" || var=="tm") {
120        double norm = 0.;
121        for (size_t i = normBin1; i < normBin2; i++) { //only normalize in the resummation region.
122          norm+=hist->bin(i).area();
123        }
124        if (norm > 0.) {
125          hist->scaleW(1.0/(norm));
126        }
127      }
128      else if ( var=="zg" || var=="tzg") {
129        normalize(hist);
130      }
131      else{
132        normalize(hist);
133      }
134    }
135
136
137    int return_bin(float pT, float rho, std::string whichvar, size_t beta) {
138      // First thing's first
139      if (pT < ptBins[0]) return -100;
140
141      if (whichvar=="m" && rho < pow(10,-4.5)) return -100;
142      if (whichvar=="tm" && rho < pow(10,-4.5)) return -100;
143
144      if (whichvar=="zg" && rho <= 0) return -100;
145      if (whichvar=="tzg" && rho <= 0) return -100;
146
147      if (whichvar=="rg" && rho <= -1.2) return -100;
148      if (whichvar=="trg" && rho <= -1.2) return -100;
149
150      if (whichvar=="id" && rho <= 1) return -100;
151
152      int pTbin = 1;
153      if (pT < ptBins[0]) pTbin = 0; //should not happen
154      else if (pT < ptBins[1]) pTbin = 1;
155      else if (pT < ptBins[2]) pTbin = 2;
156      else if (pT < ptBins[3]) pTbin = 3;
157      else if (pT < ptBins[4]) pTbin = 4;
158      else pTbin = 5;
159      if (pTbin == 0) return -1;
160
161      int rhobin = 1.;
162      if ((whichvar=="m") || (whichvar=="tm"))
163        {
164          if (rho < pow(10,-4.5)) rhobin = 0; //this should not happen.
165          else if (rho < pow(10,-4.1)) rhobin = 1;
166          else if (rho < pow(10,-3.7)) rhobin = 2;
167          else if (rho < pow(10,-3.3)) rhobin = 3;
168          else if (rho < pow(10,-2.9)) rhobin = 4;
169          else if (rho < pow(10,-2.5)) rhobin = 5;
170          else if (rho < pow(10,-2.1)) rhobin = 6;
171          else if (rho < pow(10,-1.7)) rhobin = 7;
172          else if (rho < pow(10,-1.3)) rhobin = 8;
173          else if (rho < pow(10,-0.9)) rhobin = 9;
174          else if (rho < pow(10,-0.5)) rhobin = 10;
175          else rhobin = 10;
176          return rhobin*1. + (pTbin*1.-1.)*10.-1;
177        }
178
179      // zg
180      else if (((whichvar=="zg")||(whichvar=="tzg")) && beta == 0)
181        {
182          if (rho < 0.10) return -10;
183          else if (rho < 0.15) rhobin = 1;
184          else if (rho < 0.20) rhobin = 2;
185          else if (rho < 0.25) rhobin = 3;
186          else if (rho < 0.30) rhobin = 4;
187          else if (rho < 0.35) rhobin = 5;
188          else if (rho < 0.40) rhobin = 6;
189          else if (rho < 0.45) rhobin = 7;
190          else if (rho < 0.50) rhobin = 8;
191          else rhobin = 8;
192          return rhobin*1. + (pTbin*1.-1.)*8.-1;
193        }
194      else if (((whichvar=="zg")||(whichvar=="tzg")) && beta != 0)
195        {
196          if (rho < 0.00) return -10;
197          else if (rho < 0.05) rhobin = 1;
198          else if (rho < 0.10) rhobin = 2;
199          else if (rho < 0.15) rhobin = 3;
200          else if (rho < 0.20) rhobin = 4;
201          else if (rho < 0.25) rhobin = 5;
202          else if (rho < 0.30) rhobin = 6;
203          else if (rho < 0.35) rhobin = 7;
204          else if (rho < 0.40) rhobin = 8;
205          else if (rho < 0.45) rhobin = 9;
206          else if (rho < 0.50) rhobin = 10;
207          else rhobin = 10;
208          return rhobin*1. + (pTbin*1.-1.)*10.-1;
209        }
210
211      //rg
212      else if ((whichvar=="rg")||(whichvar=="trg"))
213        {
214          if (rho < -1.2) return -10;
215          else if (rho < -1.0) rhobin = 1;
216          else if (rho < -0.8) rhobin = 2;
217          else if (rho < -0.6) rhobin = 3;
218          else if (rho < -0.4) rhobin = 4;
219          else if (rho < -0.2) rhobin = 5;
220          else if (rho < -0.1) rhobin = 6;
221          else rhobin = 6;
222          return rhobin*1. + (pTbin*1.-1.)*6.-1;
223        }
224
225      return -100;
226    }
227
228
229
230    /// Book cuts and projections
231    void init() {
232      // All final state particles
233      const FinalState fs(Cuts::abseta < 5.0);
234
235      FastJets jets(fs, FastJets::ANTIKT, 0.8, JetAlg::Muons::NONE, JetAlg::Invisibles::NONE);
236      declare(jets, "jets");
237
238      ChargedFinalState tracks(Cuts::pT > 0.5*GeV && Cuts::abseta < 2.5);
239      declare(tracks, "tracks");
240
241      normBin1 = 2; normBin2 = 7;
242      gluonFractionCentral = {0.75, 0.72, 0.66, 0.61, 0.54};
243      gluonFractionForward = {0.70, 0.64, 0.57, 0.51, 0.43};
244      rgBins = {-1.2, -1.0, -0.8, -0.6, -0.4, -0.2, -0.1};
245      zgBins = {0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5};
246      zgBinsBeta0 = {0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5};
247      rhoBins = {-4.5, -4.1, -3.7, -3.3, -2.9, -2.5, -2.1, -1.7, -1.3, -0.9, -0.5};
248      ptBins = {300, 400, 600, 800, 1000, 2000};
249
250
251
252      book(_h["_h_Table1"], 1,1,1);   // Cluster, Rho, beta=0, inclusive pT, both jets
253      book(_h["_h_Table2"], 2,1,1);   // Track, Rho, beta=0,   inclusive pT, both jets
254      book(_h["_h_Table3"], 3,1,1);   // Cluster, Rho, beta=1, inclusive pT, both jets
255      book(_h["_h_Table4"], 4,1,1);   // Track, Rho, beta=1,   inclusive pT, both jets
256      book(_h["_h_Table5"], 5,1,1);   // Cluster, Rho, beta=2, inclusive pT, both jets
257      book(_h["_h_Table6"], 6,1,1);   // Track, Rho, beta=2,   inclusive pT, both jets
258
259      book(_h["_h_Table7"], 7,1,1);   // Cluster, zg, beta=0,  inclusive pT, both jets
260      book(_h["_h_Table8"], 8,1,1);   // Track, zg, beta=0,    inclusive pT, both jets
261      book(_h["_h_Table9"], 9,1,1);   // Cluster, zg, beta=1,  inclusive pT, both jets
262      book(_h["_h_Table10"], 10,1,1);  // Track, zg, beta=1,    inclusive pT, both jets
263      book(_h["_h_Table11"], 11,1,1);  // Cluster, zg, beta=2,  inclusive pT, both jets
264      book(_h["_h_Table12"], 12,1,1);  // Track, zg, beta=2,    inclusive pT, both jets
265
266      book(_h["_h_Table13"], 13,1,1);  // Cluster, rg, beta=0,  inclusive pT, both jets
267      book(_h["_h_Table14"], 14,1,1);  // Track, rg, beta=0,    inclusive pT, both jets
268      book(_h["_h_Table15"], 15,1,1);  // Cluster, rg, beta=1,  inclusive pT, both jets
269      book(_h["_h_Table16"], 16,1,1);  // Track, rg, beta=1,    inclusive pT, both jets
270      book(_h["_h_Table17"], 17,1,1);  // Cluster, rg, beta=2,  inclusive pT, both jets
271      book(_h["_h_Table18"], 18,1,1);  // Track, rg, beta=2,    inclusive pT, both jets
272
273
274      book(_h["_h_Table19"], 19,1,1);  // Cluster, Rho, beta=0, inclusive pT, central jet
275      book(_h["_h_Table20"], 20,1,1);  // Track, Rho, beta=0,   inclusive pT, central jet
276      book(_h["_h_Table21"], 21,1,1);  // Cluster, Rho, beta=1, inclusive pT, central jet
277      book(_h["_h_Table22"], 22,1,1);  // Track, Rho, beta=1,   inclusive pT, central jet
278      book(_h["_h_Table23"], 23,1,1);  // Cluster, Rho, beta=2, inclusive pT, central jet
279      book(_h["_h_Table24"], 24,1,1);  // Track, Rho, beta=2,   inclusive pT, central jet
280
281      book(_h["_h_Table25"], 25,1,1);  // Cluster, zg, beta=0,  inclusive pT, central jet
282      book(_h["_h_Table26"], 26,1,1);  // Track, zg, beta=0,    inclusive pT, central jet
283      book(_h["_h_Table27"], 27,1,1);  // Cluster, zg, beta=1,  inclusive pT, central jet
284      book(_h["_h_Table28"], 28,1,1);  // Track, zg, beta=1,    inclusive pT, central jet
285      book(_h["_h_Table29"], 29,1,1);  // Cluster, zg, beta=2,  inclusive pT, central jet
286      book(_h["_h_Table30"], 30,1,1);  // Track, zg, beta=2,    inclusive pT, central jet
287
288      book(_h["_h_Table31"], 31,1,1);  // Cluster, rg, beta=0,  inclusive pT, central jet
289      book(_h["_h_Table32"], 32,1,1);  // Track, rg, beta=0,    inclusive pT, central jet
290      book(_h["_h_Table33"], 33,1,1);  // Cluster, rg, beta=1,  inclusive pT, central jet
291      book(_h["_h_Table34"], 34,1,1);  // Track, rg, beta=1,    inclusive pT, central jet
292      book(_h["_h_Table35"], 35,1,1);  // Cluster, rg, beta=2,  inclusive pT, central jet
293      book(_h["_h_Table36"], 36,1,1);  // Track, rg, beta=2,    inclusive pT, central jet
294
295
296      book(_h["_h_Table37"], 37,1,1);  // Cluster, Rho, beta=0, inclusive pT, forward jet
297      book(_h["_h_Table38"], 38,1,1);  // Track, Rho, beta=0,   inclusive pT, forward jet
298      book(_h["_h_Table39"], 39,1,1);  // Cluster, Rho, beta=1, inclusive pT, forward jet
299      book(_h["_h_Table40"], 40,1,1);  // Track, Rho, beta=1,   inclusive pT, forward jet
300      book(_h["_h_Table41"], 41,1,1);  // Cluster, Rho, beta=2, inclusive pT, forward jet
301      book(_h["_h_Table42"], 42,1,1);  // Track, Rho, beta=2,   inclusive pT, forward jet
302
303      book(_h["_h_Table43"], 43,1,1);  // Cluster, zg, beta=0,  inclusive pT, forward jet
304      book(_h["_h_Table44"], 44,1,1);  // Track, zg, beta=0,    inclusive pT, forward jet
305      book(_h["_h_Table45"], 45,1,1);  // Cluster, zg, beta=1,  inclusive pT, forward jet
306      book(_h["_h_Table46"], 46,1,1);  // Track, zg, beta=1,    inclusive pT, forward jet
307      book(_h["_h_Table47"], 47,1,1);  // Cluster, zg, beta=2,  inclusive pT, forward jet
308      book(_h["_h_Table48"], 48,1,1);  // Track, zg, beta=2,    inclusive pT, forward jet
309
310      book(_h["_h_Table49"], 49,1,1);  // Cluster, rg, beta=0,  inclusive pT, forward jet
311      book(_h["_h_Table50"], 50,1,1);  // Track, rg, beta=0,    inclusive pT, forward jet
312      book(_h["_h_Table51"], 51,1,1);  // Cluster, rg, beta=1,  inclusive pT, forward jet
313      book(_h["_h_Table52"], 52,1,1);  // Track, rg, beta=1,    inclusive pT, forward jet
314      book(_h["_h_Table53"], 53,1,1);  // Cluster, rg, beta=2,  inclusive pT, forward jet
315      book(_h["_h_Table54"], 54,1,1);  // Track, rg, beta=2,    inclusive pT, forward jet
316
317
318
319
320      book(_h["_h_Table55"], 55,1,1);  // Cluster, Rho, beta=0, pT binned, both jets
321      book(_h["_h_Table56"], 56,1,1);  // Track, Rho, beta=0,   pT binned, both jets
322      book(_h["_h_Table57"], 57,1,1);  // Cluster, Rho, beta=1, pT binned, both jets
323      book(_h["_h_Table58"], 58,1,1);  // Track, Rho, beta=1,   pT binned, both jets
324      book(_h["_h_Table59"], 59,1,1);  // Cluster, Rho, beta=2, pT binned, both jets
325      book(_h["_h_Table60"], 60,1,1);  // Track, Rho, beta=2,   pT binned, both jets
326
327      book(_h["_h_Table61"], 61,1,1);  // Cluster, zg, beta=0,  pT binned, both jets
328      book(_h["_h_Table62"], 62,1,1);  // Track, zg, beta=0,    pT binned, both jets
329      book(_h["_h_Table63"], 63,1,1);  // Cluster, zg, beta=1,  pT binned, both jets
330      book(_h["_h_Table64"], 64,1,1);  // Track, zg, beta=1,    pT binned, both jets
331      book(_h["_h_Table65"], 65,1,1);  // Cluster, zg, beta=2,  pT binned, both jets
332      book(_h["_h_Table66"], 66,1,1);  // Track, zg, beta=2,    pT binned, both jets
333
334      book(_h["_h_Table67"], 67,1,1);  // Cluster, rg, beta=0,  pT binned, both jets
335      book(_h["_h_Table68"], 68,1,1);  // Track, rg, beta=0,    pT binned, both jets
336      book(_h["_h_Table69"], 69,1,1);  // Cluster, rg, beta=1,  pT binned, both jets
337      book(_h["_h_Table70"], 70,1,1);  // Track, rg, beta=1,    pT binned, both jets
338      book(_h["_h_Table71"], 71,1,1);  // Cluster, rg, beta=2,  pT binned, both jets
339      book(_h["_h_Table72"], 72,1,1);  // Track, rg, beta=2,    pT binned, both jets
340
341
342      book(_h["_h_Table73"], 73,1,1);  // Cluster, Rho, beta=0, pT binned, central jet
343      book(_h["_h_Table74"], 74,1,1);  // Track, Rho, beta=0,   pT binned, central jet
344      book(_h["_h_Table75"], 75,1,1);  // Cluster, Rho, beta=1, pT binned, central jet
345      book(_h["_h_Table76"], 76,1,1);  // Track, Rho, beta=1,   pT binned, central jet
346      book(_h["_h_Table77"], 77,1,1);  // Cluster, Rho, beta=2, pT binned, central jet
347      book(_h["_h_Table78"], 78,1,1);  // Track, Rho, beta=2,   pT binned, central jet
348
349      book(_h["_h_Table79"], 79,1,1);  // Cluster, zg, beta=0,  pT binned, central jet
350      book(_h["_h_Table80"], 80,1,1);  // Track, zg, beta=0,    pT binned, central jet
351      book(_h["_h_Table81"], 81,1,1);  // Cluster, zg, beta=1,  pT binned, central jet
352      book(_h["_h_Table82"], 82,1,1); // Track, zg, beta=1,    pT binned, central jet
353      book(_h["_h_Table83"], 83,1,1); // Cluster, zg, beta=2,  pT binned, central jet
354      book(_h["_h_Table84"], 84,1,1); // Track, zg, beta=2,    pT binned, central jet
355
356      book(_h["_h_Table85"], 85,1,1); // Cluster, rg, beta=0,  pT binned, central jet
357      book(_h["_h_Table86"], 86,1,1); // Track, rg, beta=0,    pT binned, central jet
358      book(_h["_h_Table87"], 87,1,1); // Cluster, rg, beta=1,  pT binned, central jet
359      book(_h["_h_Table88"], 88,1,1); // Track, rg, beta=1,    pT binned, central jet
360      book(_h["_h_Table89"], 89,1,1); // Cluster, rg, beta=2,  pT binned, central jet
361      book(_h["_h_Table90"], 90,1,1); // Track, rg, beta=2,    pT binned, central jet
362
363
364      book(_h["_h_Table91"], 91,1,1); // Cluster, Rho, beta=0, pT binned, forward jet
365      book(_h["_h_Table92"], 92,1,1); // Track, Rho, beta=0,   pT binned, forward jet
366      book(_h["_h_Table93"], 93,1,1); // Cluster, Rho, beta=1, pT binned, forward jet
367      book(_h["_h_Table94"], 94,1,1); // Track, Rho, beta=1,   pT binned, forward jet
368      book(_h["_h_Table95"], 95,1,1); // Cluster, Rho, beta=2, pT binned, forward jet
369      book(_h["_h_Table96"], 96,1,1); // Track, Rho, beta=2,   pT binned, forward jet
370
371      book(_h["_h_Table97"], 97,1,1); // Cluster, zg, beta=0,  pT binned, forward jet
372      book(_h["_h_Table98"], 98,1,1); // Track, zg, beta=0,    pT binned, forward jet
373      book(_h["_h_Table99"], 99,1,1); // Cluster, zg, beta=1,  pT binned, forward jet
374      book(_h["_h_Table100"], 100,1,1); // Track, zg, beta=1,    pT binned, forward jet
375      book(_h["_h_Table101"], 101,1,1); // Cluster, zg, beta=2,  pT binned, forward jet
376      book(_h["_h_Table102"], 102,1,1); // Track, zg, beta=2,    pT binned, forward jet
377
378      book(_h["_h_Table103"], 103,1,1); // Cluster, rg, beta=0,  pT binned, forward jet
379      book(_h["_h_Table104"], 104,1,1); // Track, rg, beta=0,    pT binned, forward jet
380      book(_h["_h_Table105"], 105,1,1); // Cluster, rg, beta=1,  pT binned, forward jet
381      book(_h["_h_Table106"], 106,1,1); // Track, rg, beta=1,    pT binned, forward jet
382      book(_h["_h_Table107"], 107,1,1); // Cluster, rg, beta=2,  pT binned, forward jet
383      book(_h["_h_Table108"], 108,1,1); // Track, rg, beta=2,    pT binned, forward jet
384
385
386
387
388      book(_h["_h_Table109"], 109,1,1);   // Cluster, Rho, beta=0, quark jet
389      book(_h["_h_Table110"], 110,1,1);   // Track, Rho, beta=0,   quark jet
390      book(_h["_h_Table111"], 111,1,1);   // Cluster, Rho, beta=1, quark jet
391      book(_h["_h_Table112"], 112,1,1);   // Track, Rho, beta=1,   quark jet
392      book(_h["_h_Table113"], 113,1,1);   // Cluster, Rho, beta=2, quark jet
393      book(_h["_h_Table114"], 114,1,1);   // Track, Rho, beta=2,   quark jet
394
395      book(_h["_h_Table115"], 115,1,1);   // Cluster, zg, beta=0,  quark jet
396      book(_h["_h_Table116"], 116,1,1);   // Track, zg, beta=0,    quark jet
397      book(_h["_h_Table117"], 117,1,1);   // Cluster, zg, beta=1,  quark jet
398      book(_h["_h_Table118"], 118,1,1);  // Track, zg, beta=1,    quark jet
399      book(_h["_h_Table119"], 119,1,1);  // Cluster, zg, beta=2,  quark jet
400      book(_h["_h_Table120"], 120,1,1);  // Track, zg, beta=2,    quark jet
401
402      book(_h["_h_Table121"], 121,1,1);  // Cluster, rg, beta=0,  quark jet
403      book(_h["_h_Table122"], 122,1,1);  // Track, rg, beta=0,    quark jet
404      book(_h["_h_Table123"], 123,1,1);  // Cluster, rg, beta=1,  quark jet
405      book(_h["_h_Table124"], 124,1,1);  // Track, rg, beta=1,    quark jet
406      book(_h["_h_Table125"], 125,1,1);  // Cluster, rg, beta=2,  quark jet
407      book(_h["_h_Table126"], 126,1,1);  // Track, rg, beta=2,    quark jet
408
409
410      book(_h["_h_Table127"], 127,1,1);   // Cluster, Rho, beta=0, gluon jet
411      book(_h["_h_Table128"], 128,1,1);   // Track, Rho, beta=0,   gluon jet
412      book(_h["_h_Table129"], 129,1,1);   // Cluster, Rho, beta=1, gluon jet
413      book(_h["_h_Table130"], 130,1,1);   // Track, Rho, beta=1,   gluon jet
414      book(_h["_h_Table131"], 131,1,1);   // Cluster, Rho, beta=2, gluon jet
415      book(_h["_h_Table132"], 132,1,1);   // Track, Rho, beta=2,   gluon jet
416
417      book(_h["_h_Table133"], 133,1,1);   // Cluster, zg, beta=0,  gluon jet
418      book(_h["_h_Table134"], 134,1,1);   // Track, zg, beta=0,    gluon jet
419      book(_h["_h_Table135"], 135,1,1);   // Cluster, zg, beta=1,  gluon jet
420      book(_h["_h_Table136"], 136,1,1);  // Track, zg, beta=1,    gluon jet
421      book(_h["_h_Table137"], 137,1,1);  // Cluster, zg, beta=2,  gluon jet
422      book(_h["_h_Table138"], 138,1,1);  // Track, zg, beta=2,    gluon jet
423
424      book(_h["_h_Table139"], 139,1,1);  // Cluster, rg, beta=0,  gluon jet
425      book(_h["_h_Table140"], 140,1,1);  // Track, rg, beta=0,    gluon jet
426      book(_h["_h_Table141"], 141,1,1);  // Cluster, rg, beta=1,  gluon jet
427      book(_h["_h_Table142"], 142,1,1);  // Track, rg, beta=1,    gluon jet
428      book(_h["_h_Table143"], 143,1,1);  // Cluster, rg, beta=2,  gluon jet
429      book(_h["_h_Table144"], 144,1,1);  // Track, rg, beta=2,    gluon jet
430
431    }
432
433
434    void analyze(const Event& event) {
435
436      const Jets& myJets = apply<FastJets>(event, "jets").jetsByPt(200*GeV);
437      const Particles& tracks = apply<ChargedFinalState>(event, "tracks").particlesByPt();
438
439      if (myJets.size() < 2)  vetoEvent;
440      if (myJets[0].pT() > 1.5*myJets[1].pT())  vetoEvent;
441      if (myJets[0].abseta() > 1.5 || myJets[1].abseta() > 1.5) vetoEvent;
442
443      std::vector<bool> isCentral;
444      isCentral.push_back(true);
445      isCentral.push_back(false);
446      if (myJets[0].abseta() > myJets[1].abseta()) {
447        isCentral[0] = false;
448        isCentral[1] = true;
449      }
450
451      for (size_t i = 0; i < 2; ++i) {
452        if (myJets[i].pT() < 300*GeV) continue;
453
454        vector<fastjet::PseudoJet> charged_constituents;
455        for (const Particle& p : tracks) {
456          const double dr = deltaR(myJets[i], p, PSEUDORAPIDITY);
457          if (dr > 0.8) continue;
458          if (abs(p.pid()) == 13) continue;
459          charged_constituents.push_back(p);
460        }
461
462        fastjet::ClusterSequence cs_ca(myJets[i].constituents(), fastjet::JetDefinition(fastjet::cambridge_algorithm, 0.8));
463        vector<fastjet::PseudoJet> myJet_ca = fastjet::sorted_by_pt(cs_ca.inclusive_jets(10.0));
464
465        fastjet::ClusterSequence cs_ca_charged(charged_constituents, fastjet::JetDefinition(fastjet::cambridge_algorithm, 0.8));
466        vector<fastjet::PseudoJet> myJet_ca_charged = fastjet::sorted_by_pt(cs_ca_charged.inclusive_jets(10.0));
467
468        if (myJet_ca.size()==0) continue;
469        if (myJet_ca_charged.size()==0) continue;
470
471        // grooming parameters that are scanned.
472        vector<size_t> betas = { 0, 1, 2 };
473        for (size_t ibeta : betas) {
474
475          fastjet::contrib::SoftDrop sd(double(ibeta), 0.1); //beta, zcut
476          fastjet::PseudoJet sdJet = sd(myJet_ca[0]);
477          fastjet::PseudoJet sdJet_charged = sd(myJet_ca_charged[0]);
478
479          double rho2              = pow(sdJet.m()/myJets[i].pT(),2);
480          double log10rho2         = log(rho2)/log(10.);
481          double zg                = sdJet.structure_of<fastjet::contrib::SoftDrop>().symmetry();
482          double rg                = sdJet.structure_of<fastjet::contrib::SoftDrop>().delta_R();
483          rg = (rg > 0) ? log(rg) / log(10.) : -100;
484
485          double rho2_charged      = pow(sdJet_charged.m()/myJet_ca_charged[0].pt(),2);
486          double log10rho2_charged = log(rho2_charged)/log(10.);
487          double zg_charged        = sdJet_charged.structure_of<fastjet::contrib::SoftDrop>().symmetry();
488          double rg_charged             = sdJet.structure_of<fastjet::contrib::SoftDrop>().delta_R();
489          rg_charged = (rg_charged > 0) ? log(rg_charged) / log(10.) : -100;
490
491          double pt_log10rho2 = return_bin(myJets[i].pT()/GeV, rho2, "m", ibeta);
492          double pt_zg = return_bin(myJets[i].pT()/GeV, zg, "zg", ibeta);
493          double pt_rg = return_bin(myJets[i].pT()/GeV, rg, "rg", ibeta);
494
495          double pt_log10rho2_charged = return_bin(myJets[i].pT()/GeV, rho2_charged, "tm", ibeta);
496          double pt_zg_charged = return_bin(myJets[i].pT()/GeV,        zg_charged,   "tzg", ibeta);
497          double pt_rg_charged = return_bin(myJets[i].pT()/GeV,        rg_charged,   "trg", ibeta);
498
499
500          if (ibeta==0)  {
501            _h["_h_Table1"]->fill(  log10rho2);
502            _h["_h_Table2"]->fill(  log10rho2_charged);
503            _h["_h_Table7"]->fill(  zg);
504            _h["_h_Table8"]->fill(  zg_charged);
505            if (rg > -1.2)
506              _h["_h_Table13"]->fill(  rg);
507            if (rg_charged > -1.2)
508              _h["_h_Table14"]->fill(  rg_charged);
509
510            _h["_h_Table55"]->fill(  pt_log10rho2);
511            _h["_h_Table56"]->fill(  pt_log10rho2_charged);
512            _h["_h_Table61"]->fill(  pt_zg);
513            _h["_h_Table62"]->fill(  pt_zg_charged);
514            _h["_h_Table67"]->fill(  pt_rg);
515            _h["_h_Table68"]->fill(  pt_rg_charged);
516
517            if (isCentral[i]) {
518              _h["_h_Table19"]->fill(  log10rho2);
519              _h["_h_Table20"]->fill(  log10rho2_charged);
520              _h["_h_Table25"]->fill(  zg);
521              _h["_h_Table26"]->fill(  zg_charged);
522              if (rg > -1.2)
523                _h["_h_Table31"]->fill(  rg);
524              if (rg_charged > -1.2)
525                _h["_h_Table32"]->fill(  rg_charged);
526
527              _h["_h_Table73"]->fill(  pt_log10rho2);
528              _h["_h_Table74"]->fill(  pt_log10rho2_charged);
529              _h["_h_Table79"]->fill(  pt_zg);
530              _h["_h_Table80"]->fill(  pt_zg_charged);
531              _h["_h_Table85"]->fill(  pt_rg);
532              _h["_h_Table86"]->fill(  pt_rg_charged);
533            }
534
535            if (!isCentral[i]) {
536              _h["_h_Table37"]->fill(  log10rho2);
537              _h["_h_Table38"]->fill(  log10rho2_charged);
538              _h["_h_Table43"]->fill(  zg);
539              _h["_h_Table44"]->fill(  zg_charged);
540              if (rg > -1.2)
541                _h["_h_Table49"]->fill(  rg);
542              if (rg_charged > -1.2)
543                _h["_h_Table50"]->fill(  rg_charged);
544
545              _h["_h_Table91"]->fill(  pt_log10rho2);
546              _h["_h_Table92"]->fill(  pt_log10rho2_charged);
547              _h["_h_Table97"]->fill(  pt_zg);
548              _h["_h_Table98"]->fill(  pt_zg_charged);
549              _h["_h_Table103"]->fill(  pt_rg);
550              _h["_h_Table104"]->fill(  pt_rg_charged);
551            }
552          }
553          if (ibeta==1)  {
554            _h["_h_Table3"]->fill(  log10rho2);
555            _h["_h_Table4"]->fill(  log10rho2_charged);
556            _h["_h_Table9"]->fill(  zg);
557            _h["_h_Table10"]->fill(  zg_charged);
558            if (rg > -1.2)
559              _h["_h_Table15"]->fill(  rg);
560            if (rg_charged > -1.2)
561              _h["_h_Table16"]->fill(  rg_charged);
562
563            _h["_h_Table57"]->fill(  pt_log10rho2);
564            _h["_h_Table58"]->fill(  pt_log10rho2_charged);
565            _h["_h_Table63"]->fill(  pt_zg);
566            _h["_h_Table64"]->fill(  pt_zg_charged);
567            _h["_h_Table69"]->fill(  pt_rg);
568            _h["_h_Table70"]->fill(  pt_rg_charged);
569
570            if (isCentral[i]) {
571              _h["_h_Table21"]->fill(  log10rho2);
572              _h["_h_Table22"]->fill(  log10rho2_charged);
573              _h["_h_Table27"]->fill(  zg);
574              _h["_h_Table28"]->fill(  zg_charged);
575              if (rg > -1.2)
576                _h["_h_Table33"]->fill(  rg);
577              if (rg_charged > -1.2)
578                _h["_h_Table34"]->fill(  rg_charged);
579
580              _h["_h_Table75"]->fill(  pt_log10rho2);
581              _h["_h_Table76"]->fill(  pt_log10rho2_charged);
582              _h["_h_Table81"]->fill(  pt_zg);
583              _h["_h_Table82"]->fill(  pt_zg_charged);
584              _h["_h_Table87"]->fill(  pt_rg);
585              _h["_h_Table88"]->fill(  pt_rg_charged);
586            }
587
588            if (!isCentral[i]) {
589              _h["_h_Table39"]->fill(  log10rho2);
590              _h["_h_Table40"]->fill(  log10rho2_charged);
591              _h["_h_Table45"]->fill(  zg);
592              _h["_h_Table46"]->fill(  zg_charged);
593              if (rg > -1.2)
594                _h["_h_Table51"]->fill(  rg);
595              if (rg_charged > -1.2)
596                _h["_h_Table52"]->fill(  rg_charged);
597
598              _h["_h_Table93"]->fill(  pt_log10rho2);
599              _h["_h_Table94"]->fill(  pt_log10rho2_charged);
600              _h["_h_Table99"]->fill(  pt_zg);
601              _h["_h_Table100"]->fill(  pt_zg_charged);
602              _h["_h_Table105"]->fill(  pt_rg);
603              _h["_h_Table106"]->fill(  pt_rg_charged);
604            }
605          }
606          if (ibeta==2)  {
607            _h["_h_Table5"]->fill(  log10rho2);
608            _h["_h_Table6"]->fill(  log10rho2_charged);
609            _h["_h_Table11"]->fill(  zg);
610            _h["_h_Table12"]->fill(  zg_charged);
611            if (rg > -1.2)
612              _h["_h_Table17"]->fill(  rg);
613            if (rg_charged > -1.2)
614              _h["_h_Table18"]->fill(  rg_charged);
615
616            _h["_h_Table59"]->fill(  pt_log10rho2);
617            _h["_h_Table60"]->fill(  pt_log10rho2_charged);
618            _h["_h_Table65"]->fill(  pt_zg);
619            _h["_h_Table66"]->fill(  pt_zg_charged);
620            _h["_h_Table71"]->fill(  pt_rg);
621            _h["_h_Table72"]->fill(  pt_rg_charged);
622
623            if (isCentral[i]) {
624              _h["_h_Table23"]->fill(  log10rho2);
625              _h["_h_Table24"]->fill(  log10rho2_charged);
626              _h["_h_Table29"]->fill(  zg);
627              _h["_h_Table30"]->fill(  zg_charged);
628              if (rg > -1.2)
629                _h["_h_Table35"]->fill(  rg);
630              if (rg_charged > -1.2)
631                _h["_h_Table36"]->fill(  rg_charged);
632
633              _h["_h_Table77"]->fill(  pt_log10rho2);
634              _h["_h_Table78"]->fill(  pt_log10rho2_charged);
635              _h["_h_Table83"]->fill(  pt_zg);
636              _h["_h_Table84"]->fill(  pt_zg_charged);
637              _h["_h_Table89"]->fill(  pt_rg);
638              _h["_h_Table90"]->fill(  pt_rg_charged);
639            }
640
641            if (!isCentral[i]) {
642              _h["_h_Table41"]->fill(  log10rho2);
643              _h["_h_Table42"]->fill(  log10rho2_charged);
644              _h["_h_Table47"]->fill(  zg);
645              _h["_h_Table48"]->fill(  zg_charged);
646              if (rg > -1.2)
647                _h["_h_Table53"]->fill(  rg);
648              if (rg_charged > -1.2)
649                _h["_h_Table54"]->fill(  rg_charged);
650
651              _h["_h_Table95"]->fill(  pt_log10rho2);
652              _h["_h_Table96"]->fill(  pt_log10rho2_charged);
653              _h["_h_Table101"]->fill(  pt_zg);
654              _h["_h_Table102"]->fill(  pt_zg_charged);
655              _h["_h_Table107"]->fill(  pt_rg);
656              _h["_h_Table108"]->fill(  pt_rg_charged);
657            }
658          }
659
660        }
661      }
662    }
663
664    void finalize() {
665      //Normalization comes here.
666      histNorm(_h["_h_Table1"], "m");
667      histNorm(_h["_h_Table2"], "m");
668      histNorm(_h["_h_Table3"], "m");
669      histNorm(_h["_h_Table4"], "m");
670      histNorm(_h["_h_Table5"], "m");
671      histNorm(_h["_h_Table6"], "m");
672      histNorm(_h["_h_Table7"], "zg");
673      histNorm(_h["_h_Table8"], "zg");
674      histNorm(_h["_h_Table9"], "zg");
675      histNorm(_h["_h_Table10"], "zg");
676      histNorm(_h["_h_Table11"], "zg");
677      histNorm(_h["_h_Table12"], "zg");
678      histNorm(_h["_h_Table13"], "rg");
679      histNorm(_h["_h_Table14"], "rg");
680      histNorm(_h["_h_Table15"], "rg");
681      histNorm(_h["_h_Table16"], "rg");
682      histNorm(_h["_h_Table17"], "rg");
683      histNorm(_h["_h_Table18"], "rg");
684
685      histNorm(_h["_h_Table19"], "m");
686      histNorm(_h["_h_Table20"], "m");
687      histNorm(_h["_h_Table21"], "m");
688      histNorm(_h["_h_Table22"], "m");
689      histNorm(_h["_h_Table23"], "m");
690      histNorm(_h["_h_Table24"], "m");
691      histNorm(_h["_h_Table25"], "zg");
692      histNorm(_h["_h_Table26"], "zg");
693      histNorm(_h["_h_Table27"], "zg");
694      histNorm(_h["_h_Table28"], "zg");
695      histNorm(_h["_h_Table29"], "zg");
696      histNorm(_h["_h_Table30"], "zg");
697      histNorm(_h["_h_Table31"], "rg");
698      histNorm(_h["_h_Table32"], "rg");
699      histNorm(_h["_h_Table33"], "rg");
700      histNorm(_h["_h_Table34"], "rg");
701      histNorm(_h["_h_Table35"], "rg");
702      histNorm(_h["_h_Table36"], "rg");
703
704
705      histNorm(_h["_h_Table37"], "m");
706      histNorm(_h["_h_Table38"], "m");
707      histNorm(_h["_h_Table39"], "m");
708      histNorm(_h["_h_Table40"], "m");
709      histNorm(_h["_h_Table41"], "m");
710      histNorm(_h["_h_Table42"], "m");
711      histNorm(_h["_h_Table43"], "zg");
712      histNorm(_h["_h_Table44"], "zg");
713      histNorm(_h["_h_Table45"], "zg");
714      histNorm(_h["_h_Table46"], "zg");
715      histNorm(_h["_h_Table47"], "zg");
716      histNorm(_h["_h_Table48"], "zg");
717      histNorm(_h["_h_Table49"], "rg");
718      histNorm(_h["_h_Table50"], "rg");
719      histNorm(_h["_h_Table51"], "rg");
720      histNorm(_h["_h_Table52"], "rg");
721      histNorm(_h["_h_Table53"], "rg");
722      histNorm(_h["_h_Table54"], "rg");
723
724
725      ptNorm(_h["_h_Table55"], "m", 0);
726      ptNorm(_h["_h_Table56"], "m", 0);
727      ptNorm(_h["_h_Table57"], "m", 1);
728      ptNorm(_h["_h_Table58"], "m", 1);
729      ptNorm(_h["_h_Table59"], "m", 2);
730      ptNorm(_h["_h_Table60"], "m", 2);
731      ptNorm(_h["_h_Table61"], "zg", 0);
732      ptNorm(_h["_h_Table62"], "zg", 0);
733      ptNorm(_h["_h_Table63"], "zg", 1);
734      ptNorm(_h["_h_Table64"], "zg", 1);
735      ptNorm(_h["_h_Table65"], "zg", 2);
736      ptNorm(_h["_h_Table66"], "zg", 2);
737      ptNorm(_h["_h_Table67"], "rg", 0);
738      ptNorm(_h["_h_Table68"], "rg", 0);
739      ptNorm(_h["_h_Table69"], "rg", 1);
740      ptNorm(_h["_h_Table70"], "rg", 1);
741      ptNorm(_h["_h_Table71"], "rg", 2);
742      ptNorm(_h["_h_Table72"], "rg", 2);
743
744      ptNorm(_h["_h_Table73"], "m", 0);
745      ptNorm(_h["_h_Table74"], "m", 0);
746      ptNorm(_h["_h_Table75"], "m", 1);
747      ptNorm(_h["_h_Table76"], "m", 1);
748      ptNorm(_h["_h_Table77"], "m", 2);
749      ptNorm(_h["_h_Table78"], "m", 2);
750      ptNorm(_h["_h_Table79"], "zg", 0);
751      ptNorm(_h["_h_Table80"], "zg", 0);
752      ptNorm(_h["_h_Table81"], "zg", 1);
753      ptNorm(_h["_h_Table82"], "zg", 1);
754      ptNorm(_h["_h_Table83"], "zg", 2);
755      ptNorm(_h["_h_Table84"], "zg", 2);
756      ptNorm(_h["_h_Table85"], "rg", 0);
757      ptNorm(_h["_h_Table86"], "rg", 0);
758      ptNorm(_h["_h_Table87"], "rg", 1);
759      ptNorm(_h["_h_Table88"], "rg", 1);
760      ptNorm(_h["_h_Table89"], "rg", 2);
761      ptNorm(_h["_h_Table90"], "rg", 2);
762
763
764      ptNorm(_h["_h_Table91"], "m", 0);
765      ptNorm(_h["_h_Table92"], "m", 0);
766      ptNorm(_h["_h_Table93"], "m", 1);
767      ptNorm(_h["_h_Table94"], "m", 1);
768      ptNorm(_h["_h_Table95"], "m", 2);
769      ptNorm(_h["_h_Table96"], "m", 2);
770      ptNorm(_h["_h_Table97"], "zg", 0);
771      ptNorm(_h["_h_Table98"], "zg", 0);
772      ptNorm(_h["_h_Table99"], "zg", 1);
773      ptNorm(_h["_h_Table100"], "zg", 1);
774      ptNorm(_h["_h_Table101"], "zg", 2);
775      ptNorm(_h["_h_Table102"], "zg", 2);
776      ptNorm(_h["_h_Table103"], "rg", 0);
777      ptNorm(_h["_h_Table104"], "rg", 0);
778      ptNorm(_h["_h_Table105"], "rg", 1);
779      ptNorm(_h["_h_Table106"], "rg", 1);
780      ptNorm(_h["_h_Table107"], "rg", 2);
781      ptNorm(_h["_h_Table108"], "rg", 2);
782
783
784
785      getQuarkGluon(_h["_h_Table91"], _h["_h_Table73"], _h["_h_Table109"], _h["_h_Table127"], 2, "m", 0);
786      getQuarkGluon(_h["_h_Table92"], _h["_h_Table74"], _h["_h_Table110"], _h["_h_Table128"], 2, "m", 0);
787      getQuarkGluon(_h["_h_Table93"], _h["_h_Table75"], _h["_h_Table111"], _h["_h_Table129"], 2, "m", 1);
788      getQuarkGluon(_h["_h_Table94"], _h["_h_Table76"], _h["_h_Table112"], _h["_h_Table130"], 2, "m", 1);
789      getQuarkGluon(_h["_h_Table95"], _h["_h_Table77"], _h["_h_Table113"], _h["_h_Table131"], 2, "m", 2);
790      getQuarkGluon(_h["_h_Table96"], _h["_h_Table78"], _h["_h_Table114"], _h["_h_Table132"], 2, "m", 2);
791      getQuarkGluon(_h["_h_Table97"], _h["_h_Table79"], _h["_h_Table115"], _h["_h_Table133"], 2, "zg", 0);
792      getQuarkGluon(_h["_h_Table98"], _h["_h_Table80"], _h["_h_Table116"], _h["_h_Table134"], 2, "zg", 0);
793      getQuarkGluon(_h["_h_Table99"], _h["_h_Table81"], _h["_h_Table117"], _h["_h_Table135"], 2, "zg", 1);
794      getQuarkGluon(_h["_h_Table100"], _h["_h_Table82"], _h["_h_Table118"], _h["_h_Table136"], 2, "zg", 1);
795      getQuarkGluon(_h["_h_Table101"], _h["_h_Table83"], _h["_h_Table119"], _h["_h_Table137"], 2, "zg", 2);
796      getQuarkGluon(_h["_h_Table102"], _h["_h_Table84"], _h["_h_Table120"], _h["_h_Table138"], 2, "zg", 2);
797      getQuarkGluon(_h["_h_Table103"], _h["_h_Table85"], _h["_h_Table121"], _h["_h_Table139"], 2, "rg", 0);
798      getQuarkGluon(_h["_h_Table104"], _h["_h_Table86"], _h["_h_Table122"], _h["_h_Table140"], 2, "rg", 0);
799      getQuarkGluon(_h["_h_Table105"], _h["_h_Table87"], _h["_h_Table123"], _h["_h_Table141"], 2, "rg", 1);
800      getQuarkGluon(_h["_h_Table106"], _h["_h_Table88"], _h["_h_Table124"], _h["_h_Table142"], 2, "rg", 1);
801      getQuarkGluon(_h["_h_Table107"], _h["_h_Table89"], _h["_h_Table125"], _h["_h_Table143"], 2, "rg", 2);
802      getQuarkGluon(_h["_h_Table108"], _h["_h_Table90"], _h["_h_Table126"], _h["_h_Table144"], 2, "rg", 2);
803    }
804
805  protected:
806
807    size_t normBin1, normBin2;
808    vector<double> gluonFractionCentral, gluonFractionForward;
809    vector<double> rgBins, zgBins, rhoBins, ptBins, zgBinsBeta0;
810
811  private:
812    map<string, Histo1DPtr> _h;
813
814  };
815
816  RIVET_DECLARE_PLUGIN(ATLAS_2019_I1772062);
817}