Rivet analyses referenceCDF_2008_I768451Measurement of differential Z/γ∗ + jet + X cross sectionsExperiment: CDF (Tevatron Run 2) Inspire ID: 768451 Status: VALIDATED Authors:
Beam energies: (980.0, 980.0) GeV Run details:
Cross sections as a function of jet transverse momentum in 1 and 2 jet events, and jet multiplicity in pˉp collisions at √s=1.96 TeV, based on an integrated luminosity of 1.7fb−1. The measurements cover the rapidity region |yjet|<2.1 and the transverse momentum range pjet⊥>30GeV/c. Source code: CDF_2008_I768451.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/IdentifiedFinalState.hh"
5#include "Rivet/Projections/VetoedFinalState.hh"
6#include "Rivet/Projections/FastJets.hh"
7
8namespace Rivet {
9
10
11 /// @brief Measurement differential Z/\f$ \gamma^* \f$ + jet + \f$ X \f$ cross sections
12 ///
13 /// @author Frank Siegert
14 class CDF_2008_I768451 : public Analysis {
15 public:
16
17 RIVET_DEFAULT_ANALYSIS_CTOR(CDF_2008_I768451);
18
19
20 /// @name Analysis methods
21 /// @{
22
23 /// Book histograms
24 void init() {
25 // Full final state
26 FinalState fs((Cuts::etaIn(-5.0, 5.0)));
27 declare(fs, "FS");
28
29 // Leading electrons in tracking acceptance
30 IdentifiedFinalState elfs(Cuts::abseta < 5 && Cuts::pT > 25*GeV);
31 elfs.acceptIdPair(PID::ELECTRON);
32 declare(elfs, "LeadingElectrons");
33
34 book(_h_jet_multiplicity ,1, 1, 1);
35 book(_h_jet_pT_cross_section_incl_1jet ,2, 1, 1);
36 book(_h_jet_pT_cross_section_incl_2jet ,3, 1, 1);
37 }
38
39
40 /// Do the analysis
41 void analyze(const Event & event) {
42 // Skip if the event is empty
43 const FinalState& fs = apply<FinalState>(event, "FS");
44 if (fs.empty()) {
45 MSG_DEBUG("Skipping event " << numEvents() << " because no final state pair found");
46 vetoEvent;
47 }
48
49 // Find the Z candidates
50 const FinalState & electronfs = apply<FinalState>(event, "LeadingElectrons");
51 std::vector<std::pair<Particle, Particle> > Z_candidates;
52 Particles all_els=electronfs.particles();
53 for (size_t i=0; i<all_els.size(); ++i) {
54 for (size_t j=i+1; j<all_els.size(); ++j) {
55 bool candidate=true;
56 double mZ = FourMomentum(all_els[i].momentum()+all_els[j].momentum()).mass()/GeV;
57 if (mZ < 66.0 || mZ > 116.0) {
58 candidate = false;
59 }
60 double abs_eta_0 = fabs(all_els[i].eta());
61 double abs_eta_1 = fabs(all_els[j].eta());
62 if (abs_eta_1 < abs_eta_0) {
63 double tmp = abs_eta_0;
64 abs_eta_0 = abs_eta_1;
65 abs_eta_1 = tmp;
66 }
67 if (abs_eta_0 > 1.0) {
68 candidate = false;
69 }
70 if (!(abs_eta_1 < 1.0 || (inRange(abs_eta_1, 1.2, 2.8)))) {
71 candidate = false;
72 }
73 if (candidate) {
74 Z_candidates.push_back(make_pair(all_els[i], all_els[j]));
75 }
76 }
77 }
78 if (Z_candidates.size() != 1) {
79 MSG_DEBUG("Skipping event " << numEvents() << " because no unique electron pair found ");
80 vetoEvent;
81 }
82
83 // Now build the jets on a FS without the electrons from the Z (including QED radiation)
84 Particles jetparts;
85 for (const Particle& p : fs.particles()) {
86 bool copy = true;
87 if (p.pid() == PID::PHOTON) {
88 FourMomentum p_e0 = Z_candidates[0].first.momentum();
89 FourMomentum p_e1 = Z_candidates[0].second.momentum();
90 FourMomentum p_P = p.momentum();
91 if (deltaR(p_e0, p_P) < 0.2) copy = false;
92 if (deltaR(p_e1, p_P) < 0.2) copy = false;
93 } else {
94 if (HepMCUtils::uniqueId(p.genParticle()) == HepMCUtils::uniqueId(Z_candidates[0].first.genParticle())) copy = false;
95 if (HepMCUtils::uniqueId(p.genParticle()) == HepMCUtils::uniqueId(Z_candidates[0].second.genParticle())) copy = false;
96 }
97 if (copy) jetparts.push_back(p);
98 }
99
100 // Proceed to lepton dressing
101 const PseudoJets pjs = mkPseudoJets(jetparts);
102 const auto jplugin = make_shared<fastjet::CDFMidPointPlugin>(0.7, 0.5, 1.0);
103 const Jets jets_all = mkJets(fastjet::ClusterSequence(pjs, jplugin.get()).inclusive_jets());
104 const Jets jets_cut = sortByPt(select(jets_all, Cuts::pT > 30*GeV && Cuts::abseta < 2.1));
105 // FastJets jetpro(JetAlg::CDFMIDPOINT, 0.7);
106 // jetpro.calc(jetparts);
107 // // Take jets with pt > 30, |eta| < 2.1:
108 // const Jets& jets = jetpro.jets();
109 // Jets jets_cut;
110 // for (const Jet& j, jets) {
111 // if (j.pT()/GeV > 30.0 && j.abseta() < 2.1) {
112 // jets_cut.push_back(j);
113 // }
114 // }
115 // // Sort by pT:
116 // sort(jets_cut.begin(), jets_cut.end(), cmpMomByPt);
117
118 // Return if there are no jets:
119 MSG_DEBUG("Num jets above 30 GeV = " << jets_cut.size());
120 if (jets_cut.empty()) {
121 MSG_DEBUG("No jets pass cuts ");
122 vetoEvent;
123 }
124
125 // Cut on Delta R between Z electrons and *all* jets
126 for (const Jet& j : jets_cut) {
127 if (deltaR(Z_candidates[0].first, j) < 0.7) vetoEvent;
128 if (deltaR(Z_candidates[0].second, j) < 0.7) vetoEvent;
129 }
130
131 // Fill histograms
132 for (size_t njet=1; njet<=jets_cut.size(); ++njet) {
133 _h_jet_multiplicity->fill(njet);
134 }
135 for (const Jet& j : jets_cut) {
136 if (jets_cut.size() > 0) {
137 _h_jet_pT_cross_section_incl_1jet->fill(j.pT());
138 }
139 if (jets_cut.size() > 1) {
140 _h_jet_pT_cross_section_incl_2jet->fill(j.pT());
141 }
142 }
143 }
144
145
146 /// Rescale histos
147 void finalize() {
148 const double invlumi = crossSection()/femtobarn/sumOfWeights();
149 scale(_h_jet_multiplicity, invlumi);
150 scale(_h_jet_pT_cross_section_incl_1jet, invlumi);
151 scale(_h_jet_pT_cross_section_incl_2jet, invlumi);
152 }
153
154 /// @}
155
156
157 private:
158
159 /// @name Histograms
160 /// @{
161 Histo1DPtr _h_jet_multiplicity;
162 Histo1DPtr _h_jet_pT_cross_section_incl_1jet;
163 Histo1DPtr _h_jet_pT_cross_section_incl_2jet;
164 /// @}
165
166 };
167
168
169
170 RIVET_DECLARE_ALIASED_PLUGIN(CDF_2008_I768451, CDF_2008_S7540469);
171
172}
|