rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2021_I1932460

Measurement of double-parton scattering in inclusive production of four jets with low transverse momentum in proton-proton collisions at $\sqrt{s}$ = 13 TeV.
Experiment: CMS (LHC)
Inspire ID: 1932460
Status: VALIDATED
Authors:
  • cms-pag-conveners-smp@cern.ch
  • Maxim Pieters
  • Hans Van Haevermaet
  • Pierre Van Mechelen
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • Event types Hard QCD; Cuts jet $p_\perp > 20$ GeV, $|\eta|$ < 4.7.

A measurement of inclusive four-jet production in proton-proton collisions at a center-of-mass energy of 13 TeV is presented. The transverse momenta of jets within $| \eta | < 4.7$ reach down to 35, 30, 25, and 20 GeV for the first-, second-, third-, and fourth- leading jet, respectively. Differential cross sections are measured as functions of the jet transverse momentum, jet pseudorapidity, and several other observables that describe the angular correlations between the jets. The measured distributions show sensitivity to different aspects of the underlying event, parton shower, and matrix element calculations. In particular, the interplay between angular correlations caused by parton shower and double-parton scattering contributions is shown to be important. The double-parton scattering contribution is extracted by means of a template fit to the data, using distributions for single-parton scattering obtained from Monte Carlo event generators and a double-parton scattering distribution constructed from inclusive single-jet events in data. The effective double-parton scattering cross section is calculated and discussed in view of previous measurements and of its dependence on the models used to describe the single-parton scattering background.

Source code: CMS_2021_I1932460.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Projections/PromptFinalState.hh"
  6
  7namespace Rivet {
  8
  9
 10  /// @brief Measurement of double-parton scattering in inclusive production of four jets with low transverse momentum in proton-proton collisions at $\sqrt{s}$ = 13 TeV.
 11  
 12  class CMS_2021_I1932460 : public Analysis {
 13  public:
 14
 15    /// Constructor
 16    DEFAULT_RIVET_ANALYSIS_CTOR(CMS_2021_I1932460);
 17
 18
 19    /// @name Analysis methods
 20    //@{
 21
 22    /// Book histograms and initialise projections before the run
 23    void init() {
 24
 25      // Initialise and register projections
 26
 27      // the basic final-state projection: 
 28      // all final-state particles within 
 29      // the given eta acceptance: 4.7 (jet range) + 0.4 (cone size)
 30      const FinalState fs(Cuts::abseta < 5.1); 
 31
 32      // the final-state particles declared above are clustered using FastJet with
 33      // the anti-kT algorithm and a jet-radius parameter 0.4
 34      // muons and neutrinos are excluded from the clustering
 35      FastJets jetfs(fs, FastJets::ANTIKT, 0.4, JetAlg::Muons::NONE, JetAlg::Invisibles::NONE);
 36      declare(jetfs, "jets");
 37
 38      // Book histograms
 39      book(_h["JetPt1"], 1, 1, 1);
 40      book(_h["JetPt2"], 2, 1, 1);
 41      book(_h["JetPt3"], 3, 1, 1);
 42      book(_h["JetPt4"], 4, 1, 1);
 43      book(_h["JetEta1"], 5, 1, 1);
 44      book(_h["JetEta2"], 6, 1, 1);
 45      book(_h["JetEta3"], 7, 1, 1);
 46      book(_h["JetEta4"], 8, 1, 1);
 47      
 48      book(_h["DeltaPhiSoft_binNorm"], 9, 1, 1);
 49      book(_h["DeltaPhi3_binNorm"], 10, 1, 1);
 50      book(_h["DeltaY_binNorm"], 11, 1, 1);
 51      book(_h["DeltaPhiY_binNorm"], 12, 1, 1);
 52      book(_h["DeltaPtSoft_binNorm"], 13, 1, 1);
 53      book(_h["DeltaS_binNorm"], 14, 1, 1);
 54      
 55      book(_h["DeltaPhiSoft"], 45, 1, 1);
 56      book(_h["DeltaPhi3"], 46, 1, 1);
 57      book(_h["DeltaY"], 47, 1, 1);
 58      book(_h["DeltaPhiY"], 48, 1, 1);
 59      book(_h["DeltaPtSoft"], 49, 1, 1);
 60      book(_h["DeltaS"], 50, 1, 1);
 61
 62    }
 63
 64
 65    /// Perform the per-event analysis
 66    void analyze(const Event& event) {
 67
 68      // retrieve clustered jets, sorted by pT, with a minimum pT cut 10 GeV and eta range 4.7 (similar to PFJet collection)
 69      Jets jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::abseta < 4.7 && Cuts::pT > 10*GeV);
 70      
 71      // fill only if there are at least 4 jets
 72      if (jets.size() < 4) vetoEvent;
 73      
 74      double pt0 = jets[0].pt();
 75      double pt1 = jets[1].pt();
 76      double pt2 = jets[2].pt();
 77      double pt3 = jets[3].pt();
 78
 79      // pt selection
 80      if (pt0 < 35.0 || pt1 < 30.0 || pt2 < 25.0 || pt3 < 20.0) vetoEvent;
 81
 82      double phi0 = jets[0].phi();
 83      double phi1 = jets[1].phi();
 84      double phi2 = jets[2].phi();
 85      double phi3 = jets[3].phi();
 86
 87      _h["JetPt1"]->fill(pt0);
 88      _h["JetPt2"]->fill(pt1);
 89      _h["JetPt3"]->fill(pt2);
 90      _h["JetPt4"]->fill(pt3);
 91
 92      _h["JetEta1"]->fill(jets[0].eta());
 93      _h["JetEta2"]->fill(jets[1].eta());
 94      _h["JetEta3"]->fill(jets[2].eta());
 95      _h["JetEta4"]->fill(jets[3].eta());
 96
 97      // delta phi and eta of the 2 soft jets
 98      _h["DeltaPhiSoft"]->fill(abs(deltaPhi(phi2, phi3)));
 99      _h["DeltaPhiSoft_binNorm"]->fill(abs(deltaPhi(phi2, phi3)));
100
101      // delta pt between soft jets
102      double DptSoft = sqrt(pow(pt2*cos(phi2) + pt3*cos(phi3), 2) + pow(pt2*sin(phi2) + pt3*sin(phi3), 2))/(pt2 + pt3);
103      _h["DeltaPtSoft"]->fill(DptSoft);
104      _h["DeltaPtSoft_binNorm"]->fill(DptSoft);
105
106      // delta S
107      if (pt0 > 50.0 && pt1 > 30.0 && pt2 > 30.0 && pt3 > 30.0) {
108        double phiH = atan2(pt0*sin(phi0) + pt1*sin(phi1) , pt0*cos(phi0) + pt1*cos(phi1));
109        double phiS = atan2(pt2*sin(phi2) + pt3*sin(phi3) , pt2*cos(phi2) + pt3*cos(phi3));
110        double DS = abs(deltaPhi(phiH, phiS));
111        _h["DeltaS"]->fill(DS);
112        _h["DeltaS_binNorm"]->fill(DS);   
113      }
114
115      // delta Y: most remote jets in rapidity, find min & max eta
116      double mineta = 99999;
117      double maxeta = -99999;
118      int minetapos = -1;
119      int maxetapos = -1;
120  
121      for (int i = 0; i < 4; ++i) {
122        if (jets[i].eta() < mineta) {
123          mineta = jets[i].eta();
124          minetapos = i;
125        }
126        if (jets[i].eta() > maxeta) {
127          maxeta = jets[i].eta();
128          maxetapos = i;
129        }  
130      }
131  
132      _h["DeltaY"]->fill(abs(jets[minetapos].eta() - jets[maxetapos].eta()));
133      _h["DeltaY_binNorm"]->fill(abs(jets[minetapos].eta() - jets[maxetapos].eta()));
134
135      // Delta phi Y: azimuthal angle between most remote jets in eta
136      _h["DeltaPhiY"]->fill(abs(deltaPhi(jets[minetapos].phi(), jets[maxetapos].phi())));
137      _h["DeltaPhiY_binNorm"]->fill(abs(deltaPhi(jets[minetapos].phi(), jets[maxetapos].phi())));
138
139      // delta phi3
140      double minphi3 = 999;
141      for (int iphi1 = 0; iphi1 < 4; ++iphi1) {
142        for (int iphi2 = 0; iphi2 < 4; ++iphi2) {
143          for (int iphi3 = 0; iphi3 < 4; ++iphi3) {
144            if ( !(iphi1 == iphi2 || iphi2 == iphi3 || iphi1 == iphi3) ) {
145              double temp_phi1= jets[iphi1].phi();
146              double temp_phi2= jets[iphi2].phi();
147              double temp_phi3= jets[iphi3].phi();
148              double temp_minphi3 = abs(deltaPhi(temp_phi1, temp_phi2)) + abs(deltaPhi(temp_phi2, temp_phi3));
149              if (temp_minphi3 < minphi3) minphi3 = temp_minphi3;
150            }
151          }
152        }
153      }
154  
155      _h["DeltaPhi3"]->fill(minphi3);
156      _h["DeltaPhi3_binNorm"]->fill(minphi3);
157      
158    }
159
160
161    /// Normalise histograms etc., after the run
162    void finalize() {
163
164      scale(_h["JetPt1"], crossSection()/picobarn/sumOfWeights()); // norm to cross section
165      scale(_h["JetPt2"], crossSection()/picobarn/sumOfWeights());
166      scale(_h["JetPt3"], crossSection()/picobarn/sumOfWeights());
167      scale(_h["JetPt4"], crossSection()/picobarn/sumOfWeights());
168      scale(_h["JetEta1"], crossSection()/picobarn/sumOfWeights());
169      scale(_h["JetEta2"], crossSection()/picobarn/sumOfWeights());
170      scale(_h["JetEta3"], crossSection()/picobarn/sumOfWeights());
171      scale(_h["JetEta4"], crossSection()/picobarn/sumOfWeights());
172      scale(_h["DeltaPhiSoft"], crossSection()/picobarn/sumOfWeights());
173      scale(_h["DeltaPhi3"], crossSection()/picobarn/sumOfWeights());
174      scale(_h["DeltaY"], crossSection()/picobarn/sumOfWeights());
175      scale(_h["DeltaPhiY"], crossSection()/picobarn/sumOfWeights());
176      scale(_h["DeltaPtSoft"], crossSection()/picobarn/sumOfWeights());
177      scale(_h["DeltaS"], crossSection()/picobarn/sumOfWeights());
178      
179      // create bin normalised histograms
180      
181      // Correct for binwidths: Rivet automatically normalises histograms to binwidth when plotting AFTER the normalisation executed here. 
182      // So we must calculate an extra correction here so that finally our bin-normalised histograms end up around 1 as in the paper.
183      // in YODA bin index starts from 0
184      
185      // For DeltaY, which has variable binwidths we need to do following steps
186      // divide histograms by binwidth
187      for (unsigned int i = 0; i < _h["DeltaY_binNorm"]->numBins(); ++i) {
188        _h["DeltaY_binNorm"]->bin(i).scaleW(1.0/_h["DeltaY_binNorm"]->bin(i).xWidth());
189      }
190      
191      // normalise to average of first 4 bins
192      scale(_h["DeltaY_binNorm"], 1.0/(_h["DeltaY_binNorm"]->integralRange(0,3)/4.0));
193      
194      // multiply again with binwidth
195      for (unsigned int i = 0; i < _h["DeltaY_binNorm"]->numBins(); ++i) {
196        _h["DeltaY_binNorm"]->bin(i).scaleW(_h["DeltaY_binNorm"]->bin(i).xWidth());
197      }
198      
199      // DeltaPhiSoft and DeltaPhi3 histograms have uniform binwidths, so multiply with first binwidth is sufficient
200      scale(_h["DeltaPhiSoft_binNorm"], _h["DeltaPhiSoft_binNorm"]->bin(0).xWidth()/(_h["DeltaPhiSoft_binNorm"]->integralRange(0,4)/5.0));
201      scale(_h["DeltaPhi3_binNorm"], _h["DeltaPhi3_binNorm"]->bin(0).xWidth()/(_h["DeltaPhi3_binNorm"]->integralRange(0,3)/4.0));
202      
203      // DeltaPhiY, DeltaPtSoft and DeltaS are normalised to last bin
204      scale(_h["DeltaPhiY_binNorm"], _h["DeltaPhiY_binNorm"]->bin(11).xWidth()/_h["DeltaPhiY_binNorm"]->bin(11).sumW() );
205      scale(_h["DeltaPtSoft_binNorm"], _h["DeltaPtSoft_binNorm"]->bin(7).xWidth()/_h["DeltaPtSoft_binNorm"]->bin(7).sumW() );
206      scale(_h["DeltaS_binNorm"], _h["DeltaS_binNorm"]->bin(6).xWidth()/_h["DeltaS_binNorm"]->bin(6).sumW() );
207
208    }
209    
210    map<string, Histo1DPtr> _h;
211
212  };
213
214
215  // The hook for the plugin system
216  DECLARE_RIVET_PLUGIN(CMS_2021_I1932460);
217
218
219}