Rivet analyses referenceJADE_OPAL_2000_I513337Jet rates in $e^+e^-$ at JADE [35--44 GeV] and OPAL [91--189 GeV].Experiment: JADE_OPAL (PETRA and LEP) Inspire ID: 513337 Status: VALIDATED Authors:
Beam energies: (17.5, 17.5); (22.0, 22.0); (45.6, 45.6); (66.5, 66.5); (80.5, 80.5); (86.0, 86.0); (91.5, 91.5); (94.5, 94.5) GeV Run details:
Differential and integrated jet rates for Durham and JADE jet algorithms. The integration cut value used for the integrated rate observables is not well-defined in the paper: the midpoint of the differential bin has been used thanks to information from Stefan Kluth and Christoph Pahl. We anyway recommend that the differential plots be preferred over the integrated ones for MC generator validation and tuning, to minimise correlations. Source code: JADE_OPAL_2000_I513337.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FastJets.hh"
4#include "Rivet/Projections/FinalState.hh"
5
6namespace Rivet {
7
8
9 /// @brief Jet rates in \f$ e^+ e^- \f$ at OPAL and JADE
10 ///
11 /// @author Frank Siegert
12 class JADE_OPAL_2000_I513337 : public Analysis {
13 public:
14
15 /// Constructor
16 RIVET_DEFAULT_ANALYSIS_CTOR(JADE_OPAL_2000_I513337)
17
18
19 /// @name Analysis methods
20 /// @{
21
22 void init() {
23 // Projections
24 const FinalState fs;
25 declare(fs, "FS");
26 FastJets jadeJets = FastJets(fs, JetAlg::JADE, 0.7, JetMuons::ALL, JetInvisibles::DECAY);
27 FastJets durhamJets = FastJets(fs, JetAlg::DURHAM, 0.7, JetMuons::ALL, JetInvisibles::DECAY);
28 declare(jadeJets, "JadeJets");
29 declare(durhamJets, "DurhamJets");
30
31 // Histos
32 int offset = 0;
33 switch (int(sqrtS()/GeV + 0.5)) {
34 case 35: offset = 7; break;
35 case 44: offset = 8; break;
36 case 91: offset = 9; break;
37 case 133: offset = 10; break;
38 case 161: offset = 11; break;
39 case 172: offset = 12; break;
40 case 183: offset = 13; break;
41 case 189: offset = 14; break;
42 default: break;
43 }
44 for (size_t i = 0; i < 5; ++i) {
45 book(_h_R_Jade[i] ,offset, 1, i+1);
46 book(_h_R_Durham[i] ,offset+9, 1, i+1);
47 if (i < 4) book(_h_y_Durham[i], offset+17, 1, i+1);
48 }
49 }
50
51
52
53 void analyze(const Event& e) {
54 MSG_DEBUG("Num particles = " << apply<FinalState>(e, "FS").particles().size());
55
56 const FastJets& jadejet = apply<FastJets>(e, "JadeJets");
57 if (jadejet.clusterSeq()) {
58 const double y_23 = jadejet.clusterSeq()->exclusive_ymerge_max(2);
59 const double y_34 = jadejet.clusterSeq()->exclusive_ymerge_max(3);
60 const double y_45 = jadejet.clusterSeq()->exclusive_ymerge_max(4);
61 const double y_56 = jadejet.clusterSeq()->exclusive_ymerge_max(5);
62
63 for (const auto& b : _h_R_Jade[0]->bins()) {
64 const string xedge = b.xEdge();
65 const double ycut = std::stod(xedge);
66 if (y_23 < ycut) {
67 _h_R_Jade[0]->fill(xedge);
68 }
69 }
70 for (const auto& b : _h_R_Jade[1]->bins()) {
71 const string xedge = b.xEdge();
72 const double ycut = std::stod(xedge);
73 if (y_34 < ycut && y_23 > ycut) {
74 _h_R_Jade[1]->fill(xedge);
75 }
76 }
77 for (const auto& b : _h_R_Jade[2]->bins()) {
78 const string xedge = b.xEdge();
79 const double ycut = std::stod(xedge);
80 if (y_45 < ycut && y_34 > ycut) {
81 _h_R_Jade[2]->fill(xedge);
82 }
83 }
84 for (const auto& b : _h_R_Jade[3]->bins()) {
85 const string xedge = b.xEdge();
86 const double ycut = std::stod(xedge);
87 if (y_56 < ycut && y_45 > ycut) {
88 _h_R_Jade[3]->fill(xedge);
89 }
90 }
91 for (const auto& b : _h_R_Jade[4]->bins()) {
92 const string xedge = b.xEdge();
93 const double ycut = std::stod(xedge);
94 if (y_56 > ycut) {
95 _h_R_Jade[4]->fill(xedge);
96 }
97 }
98 }
99
100 const FastJets& durjet = apply<FastJets>(e, "DurhamJets");
101 if (durjet.clusterSeq()) {
102 const double y_23 = durjet.clusterSeq()->exclusive_ymerge_max(2);
103 const double y_34 = durjet.clusterSeq()->exclusive_ymerge_max(3);
104 const double y_45 = durjet.clusterSeq()->exclusive_ymerge_max(4);
105 const double y_56 = durjet.clusterSeq()->exclusive_ymerge_max(5);
106
107 _h_y_Durham[0]->fill(y_23);
108 _h_y_Durham[1]->fill(y_34);
109 _h_y_Durham[2]->fill(y_45);
110 _h_y_Durham[3]->fill(y_56);
111
112 for (const auto& b : _h_R_Durham[0]->bins()) {
113 const string xedge = b.xEdge();
114 const double ycut = std::stod(xedge);
115 if (y_23 < ycut) {
116 _h_R_Durham[0]->fill(xedge);
117 }
118 }
119 for (const auto& b : _h_R_Durham[1]->bins()) {
120 const string xedge = b.xEdge();
121 const double ycut = std::stod(xedge);
122 if (y_34 < ycut && y_23 > ycut) {
123 _h_R_Durham[1]->fill(xedge);
124 }
125 }
126 for (const auto& b : _h_R_Durham[2]->bins()) {
127 const string xedge = b.xEdge();
128 const double ycut = std::stod(xedge);
129 if (y_45 < ycut && y_34 > ycut) {
130 _h_R_Durham[2]->fill(xedge);
131 }
132 }
133 for (const auto& b : _h_R_Durham[2]->bins()) {
134 const string xedge = b.xEdge();
135 const double ycut = std::stod(xedge);
136 if (y_56 < ycut && y_45 > ycut) {
137 _h_R_Durham[3]->fill(xedge);
138 }
139 }
140 for (const auto& b : _h_R_Durham[4]->bins()) {
141 const string xedge = b.xEdge();
142 const double ycut = std::stod(xedge);
143 if (y_56 > ycut) {
144 _h_R_Durham[4]->fill(xedge);
145 }
146 }
147 }
148 }
149
150
151
152 /// Finalize
153 void finalize() {
154 normalize(_h_y_Durham);
155 scale(_h_R_Jade, 100/sumOfWeights());
156 scale(_h_R_Durham, 100/sumOfWeights());
157 }
158
159 /// @}
160
161
162 private:
163
164 /// @name Histograms
165 /// @{
166 BinnedHistoPtr<string> _h_R_Jade[5];
167 BinnedHistoPtr<string> _h_R_Durham[5];
168 Histo1DPtr _h_y_Durham[4];
169 /// @}
170
171 };
172
173
174
175 RIVET_DECLARE_ALIASED_PLUGIN(JADE_OPAL_2000_I513337, JADE_OPAL_2000_S4300807);
176
177}
|