Rivet analyses referenceCMS_2013_I1218372Study of the underlying event at forward rapidity in proton--proton collisions at the LHCExperiment: CMS (LHC) Inspire ID: 1218372 Status: VALIDATED Authors:
Beam energies: (450.0, 450.0); (1380.0, 1380.0); (3500.0, 3500.0) GeV Run details:
Ratio of the energy deposited in the pseudorapidity range $-6.6 < \eta < -5.2$ for events with a charged particle jet with $|\eta|<2$ with respect to the energy in inclusive events, as a function of charged particle jet transverse momentum for $\sqrt{s}=$0.9, 2.76 and 7 TeV. Beam energy must be specified (in GeV) as analysis option "ENERGY" when rivet-merging samples. Source code: CMS_2013_I1218372.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/ChargedFinalState.hh"
5#include "Rivet/Projections/FastJets.hh"
6#include "Rivet/Projections/Beam.hh"
7#include "Rivet/Projections/VetoedFinalState.hh"
8
9namespace Rivet {
10
11
12 /// Ratio of energy in -6.6 < eta < -5.2 for events with a charged-particle jet
13 class CMS_2013_I1218372 : public Analysis {
14 public:
15
16 /// Constructor
17 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2013_I1218372);
18
19
20 void init() {
21
22 // gives the range of eta and min pT for the final state from which I get the jets
23 FastJets jetpro (ChargedFinalState((Cuts::etaIn(-2.5, 2.5) && Cuts::pT >= 0.3*GeV)), JetAlg::ANTIKT, 0.5);
24 declare(jetpro, "Jets");
25
26 // skip Neutrinos and Muons
27 VetoedFinalState fsv(FinalState((Cuts::etaIn(-7.0, -4.0))));
28 fsv.vetoNeutrinos();
29 fsv.addVetoPairId(PID::MUON);
30 declare(fsv, "fsv");
31
32 // for the hadron level selection
33 VetoedFinalState sfsv;
34 sfsv.vetoNeutrinos();
35 sfsv.addVetoPairId(PID::MUON);
36 declare(sfsv, "sfsv");
37
38 //counters
39 book(passedSumOfWeights, "passedSumOfWeights");
40 book(inclEflow, "inclEflow");
41
42 // Temporary histograms to fill the energy flow for leading jet events.
43 // Ratios are calculated in finalyze().
44 int id = 0;
45 if (isCompatibleWithSqrtS( 900)) id=1;
46 if (isCompatibleWithSqrtS(2760*GeV)) id=2;
47 if (isCompatibleWithSqrtS(7000*GeV)) id=3;
48 book(_h_ratio, id, 1, 1);
49 book(_tmp_jet , "TMP/eflow_jet" ,refData(id, 1, 1)); // Leading jet energy flow in pt
50 book(_tmp_njet, "TMP/number_jet" ,refData(id, 1, 1)); // Number of events in pt
51 }
52
53
54 /// Perform the per-event analysis
55 void analyze(const Event& event) {
56
57 // Skip if the event is empty
58 const FinalState& fsv = apply<FinalState>(event, "fsv");
59 if (fsv.empty()) vetoEvent;
60
61 // ====================== Minimum Bias selection
62
63 const FinalState& sfsv = apply<FinalState>(event, "sfsv");
64 Particles parts = sfsv.particles(cmpMomByRap);
65 if (parts.empty()) vetoEvent;
66
67 // find dymax
68 double dymax = 0;
69 int gap_pos = -1;
70 for (size_t i = 0; i < parts.size()-1; ++i) {
71 double dy = parts[i+1].rapidity() - parts[i].rapidity();
72 if (dy > dymax) {
73 dymax = dy;
74 gap_pos = i;
75 }
76 }
77
78 // calculate mx2 and my2
79 FourMomentum xmom;
80 for (int i=0; i<=gap_pos; ++i) {
81 xmom += parts[i].momentum();
82 }
83 double mx2 = xmom.mass2();
84 if (mx2<0) vetoEvent;
85
86 FourMomentum ymom;
87 for (size_t i=gap_pos+1; i<parts.size(); ++i) {
88 ymom += parts[i].momentum();
89 }
90 double my2 = ymom.mass2();
91 if (my2<0) vetoEvent;
92
93 // calculate xix and xiy and xidd
94 double xix = mx2 / sqr(sqrtS());
95 double xiy = my2 / sqr(sqrtS());
96 double xidd = mx2*my2 / sqr(sqrtS()*0.938*GeV);
97
98 // combine the selection: xi cuts
99 bool passedHadronCuts = false;
100 if (isCompatibleWithSqrtS( 900) && (xix > 0.1 || xiy > 0.4 || xidd > 0.5)) passedHadronCuts = true;
101 if (isCompatibleWithSqrtS(2760*GeV) && (xix > 0.07 || xiy > 0.2 || xidd > 0.5)) passedHadronCuts = true;
102 if (isCompatibleWithSqrtS(7000*GeV) && (xix > 0.04 || xiy > 0.1 || xidd > 0.5)) passedHadronCuts = true;
103 if (!passedHadronCuts) vetoEvent;
104
105 // ============================== MINIMUM BIAS EVENTS
106
107 // loop over particles to calculate the energy
108 passedSumOfWeights->fill();
109
110 for (const Particle& p : fsv.particles()) {
111 if (-5.2 > p.eta() && p.eta() > -6.6) inclEflow->fill(p.E()/GeV);
112 }
113
114 // ============================== JET EVENTS
115
116 const FastJets& jetpro = apply<FastJets>(event, "Jets");
117 const Jets& jets = jetpro.jetsByPt(Cuts::pT > 1.0*GeV);
118 if (jets.size()<1) vetoEvent;
119
120 if (fabs(jets[0].eta()) < 2.0) {
121 _tmp_njet->fill(jets[0].pT()/GeV);
122
123 // energy flow
124 for (const Particle& p : fsv.particles()) {
125 if (p.eta() > -6.6 && p.eta() < -5.2) { // ask for the CASTOR region
126 _tmp_jet->fill(jets[0].pT()/GeV, p.E()/GeV);
127 }
128 }
129 }
130
131 }// analysis
132
133 void finalize() {
134 scale(_tmp_jet, *passedSumOfWeights / *inclEflow);
135 divide(_tmp_jet, _tmp_njet, _h_ratio);
136 }
137
138 private:
139 // counters
140 CounterPtr passedSumOfWeights;
141 CounterPtr inclEflow;
142
143 // histograms
144 Estimate1DPtr _h_ratio;
145 Histo1DPtr _tmp_jet;
146 Histo1DPtr _tmp_njet;
147 };
148
149
150 RIVET_DECLARE_PLUGIN(CMS_2013_I1218372);
151
152}
|