Rivet analyses referenceCLEOC_2008_I777917Measurement of charm final-states, the total hadronic cross section and $R$ for energies between 3.92 and 4.26 GeVExperiment: CLEOC (CESR) Inspire ID: 777917 Status: VALIDATED Authors:
Beam energies: ANY Run details:
Measurement of charm final-states, the total hadronic cross section and $R$ for energies between 3.92 and 4.26 GeV The muonic cross section is also outputted to the yoda file so that ratio $R$ can be recalcuated if runs are combined. Source code: CLEOC_2008_I777917.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/UnstableParticles.hh"
5
6namespace Rivet {
7
8
9 /// @brief Add a short analysis description here
10 class CLEOC_2008_I777917 : public Analysis {
11 public:
12
13 /// Constructor
14 RIVET_DEFAULT_ANALYSIS_CTOR(CLEOC_2008_I777917);
15
16
17 /// @name Analysis methods
18 /// @{
19
20 /// Book histograms and initialise projections before the run
21 void init() {
22
23 // Initialise and register projections
24 declare(FinalState(), "FS");
25 declare(UnstableParticles(), "UFS");
26 book(_c_hadrons, "/TMP/sigma_hadrons");
27 book(_c_muons, "/TMP/sigma_muons");
28 book(_c_D0D0, "/TMP/sigma_D0D0");
29 book(_c_DpDm, "/TMP/sigma_DpDm");
30 book(_c_DsDs, "/TMP/sigma_DsDs");
31 book(_c_D0D0S, "/TMP/sigma_D0D0S");
32 book(_c_DpDmS, "/TMP/sigma_DpDmS");
33 book(_c_DsDsS, "/TMP/sigma_DsDsS");
34 book(_c_D0SD0S, "/TMP/sigma_D0SD0S");
35 book(_c_DpSDmS, "/TMP/sigma_DpSDmS");
36 book(_c_DsSDsS, "/TMP/sigma_DsSDsS");
37 book(_c_DD, "/TMP/sigma_DD");
38 book(_c_DDX, "/TMP/sigma_DDX");
39 book(_c_DSDpi, "/TMP/sigma_DSDpi");
40 book(_c_DSDSpi, "/TMP/sigma_DSDSpi");
41 }
42
43 void findChildren(const Particle & p,map<long,int> & nRes, int &ncount) {
44 for (const Particle &child : p.children()) {
45 if(child.children().empty()) {
46 nRes[child.pid()]-=1;
47 --ncount;
48 }
49 else
50 findChildren(child,nRes,ncount);
51 }
52 }
53
54
55 /// Perform the per-event analysis
56 void analyze(const Event& event) {
57 const FinalState& fs = apply<FinalState>(event, "FS");
58 // total hadronic and muonic cross sections
59 map<long,int> nCount;
60 int ntotal(0);
61 for (const Particle& p : fs.particles()) {
62 nCount[p.pid()] += 1;
63 ++ntotal;
64 }
65 // mu+mu- + photons
66 if(nCount[-13]==1 and nCount[13]==1 &&
67 ntotal==2+nCount[22])
68 _c_muons->fill();
69 // everything else
70 else
71 _c_hadrons->fill();
72 // identified final state with D mesons
73 const FinalState& ufs = apply<UnstableParticles>(event, "UFS");
74 for(unsigned int ix=0;ix<ufs.particles().size();++ix) {
75 const Particle& p1 = ufs.particles()[ix];
76 int id1 = abs(p1.pid());
77 if(id1 != 411 && id1 != 413 && id1 != 421 && id1 != 423 &&
78 id1 != 431 && id1 != 433)
79 continue;
80 // check fs
81 bool fs = true;
82 for (const Particle & child : p1.children()) {
83 if(child.pid()==p1.pid()) {
84 fs = false;
85 break;
86 }
87 }
88 if(!fs) continue;
89 // find the children
90 map<long,int> nRes = nCount;
91 int ncount = ntotal;
92 findChildren(p1,nRes,ncount);
93 bool matched=false;
94 int sign = p1.pid()/id1;
95 // loop over the other fs particles
96 for(unsigned int iy=ix+1;iy<ufs.particles().size();++iy) {
97 const Particle& p2 = ufs.particles()[iy];
98 fs = true;
99 for (const Particle & child : p2.children()) {
100 if(child.pid()==p2.pid()) {
101 fs = false;
102 break;
103 }
104 }
105 if(!fs) continue;
106 if(p2.pid()/abs(p2.pid())==sign) continue;
107 int id2 = abs(p2.pid());
108 if(id2 != 411 && id2 != 413 && id2 != 421 && id2 != 423 &&
109 id2 != 431 && id2 != 433)
110 continue;
111 if(!p2.parents().empty() && p2.parents()[0].pid()==p1.pid())
112 continue;
113 if((id1==411 || id1==421 || id1==431) && (id2==411 || id2==421 || id2==431 ))
114 _c_DDX->fill();
115 map<long,int> nRes2 = nRes;
116 int ncount2 = ncount;
117 findChildren(p2,nRes2,ncount2);
118 if(ncount2==0) {
119 matched=true;
120 for(auto const & val : nRes2) {
121 if(val.second!=0) {
122 matched = false;
123 break;
124 }
125 }
126 if(matched) {
127 if(id1==411 && id2==411) {
128 _c_DpDm->fill();
129 _c_DD ->fill();
130 }
131 else if(id1==421&& id2==421) {
132 _c_D0D0->fill();
133 _c_DD ->fill();
134 }
135 else if(id1==431&& id2==431) {
136 _c_DsDs->fill();
137 }
138 else if(id1==413 && id2==413) {
139 _c_DpSDmS->fill();
140 }
141 else if(id1==423&& id2==423) {
142 _c_D0SD0S->fill();
143 }
144 else if(id1==433&& id2==433) {
145 _c_DsSDsS->fill();
146 }
147 else if((id1==421 && id2==423) ||
148 (id1==423 && id2==421)) {
149 _c_D0D0S->fill();
150 }
151 else if((id1==411 && id2==413) ||
152 (id1==413 && id2==411)) {
153 _c_DpDmS->fill();
154 }
155 else if((id1==431 && id2==433) ||
156 (id1==433 && id2==431)) {
157 _c_DsDsS->fill();
158 }
159 }
160 }
161 else if(ncount2==1) {
162 int ipi=0;
163 if(nRes2[111]==1 && nRes2[211]==0 && nRes[-211]==0 )
164 ipi = 111;
165 else if(nRes2[111]==0 && nRes2[211]==1 && nRes[-211]==0 )
166 ipi = 211;
167 else if(nRes2[111]==0 && nRes2[211]==0 && nRes[-211]==1 )
168 ipi =-211;
169 if(ipi==0) continue;
170 matched=true;
171 for(auto const & val : nRes2) {
172 if(val.first==ipi)
173 continue;
174 else if(val.second!=0) {
175 matched = false;
176 break;
177 }
178 }
179 if(matched) {
180 bool Ddecay = false;
181 Particle mother = p1;
182 while (!mother.parents().empty()) {
183 mother = mother.parents()[0];
184 if(PID::isCharmMeson(mother.pid()) && mother.pid()!=p1.pid()) {
185 Ddecay = true;
186 break;
187 }
188 }
189 mother = p2;
190 while (!mother.parents().empty()) {
191 mother = mother.parents()[0];
192 if(PID::isCharmMeson(mother.pid()) && mother.pid()!=p1.pid()) {
193 Ddecay = true;
194 break;
195 }
196 }
197 if(Ddecay) continue;
198 if((id1==413 || id1==423 ) &&
199 (id2==413 || id2==423 )) {
200 _c_DSDSpi->fill();
201 }
202 else if((id1==411 || id1==421 ) &&
203 (id2==413 || id2==423 )) {
204 _c_DSDpi->fill();
205 }
206 else if((id1==413 || id1==423 ) &&
207 (id2==411 || id2==421 )) {
208 _c_DSDpi->fill();
209 }
210 }
211 }
212 }
213 }
214 }
215
216 /// Normalise histograms etc., after the run
217 void finalize() {
218 // R
219 Estimate0D R = *_c_hadrons/ *_c_muons;
220 double fact = crossSection()/picobarn/ sumOfWeights() /nanobarn;
221 double sig_h = _c_hadrons->val()*fact;
222 double err_h = _c_hadrons->err()*fact;
223 double sig_c = _c_DDX->val()*fact;
224 double err_c = _c_DDX->err()*fact;
225 double sig_m = _c_muons ->val()*fact;
226 double err_m = _c_muons ->err()*fact;
227 Estimate1DPtr charm, hadrons, muons, mult;
228 book(charm, 6,1,1);
229 book(hadrons, "sigma_hadrons");
230 book(muons, "sigma_muons" );
231 book(mult, 6,1,2);
232 for (auto& b : mult->bins()) {
233 if (inRange(sqrtS()/MeV, b.xMin(), b.xMax())) {
234 b.set(R.val(), R.errPos());
235 hadrons->bin(b.index()).set(sig_h, err_h);
236 charm ->bin(b.index()).set(sig_c, err_c);
237 muons ->bin(b.index()).set(sig_m, err_m);
238 }
239 }
240 for (unsigned int ix=1;ix<6;++ix) {
241 unsigned int imax(0);
242 if (ix<=3) imax = 4;
243 else imax = 3;
244 for (unsigned int iy=1;iy<imax;++iy) {
245 double sigma(0),error(0);
246 if(ix==1) {
247 if(iy==1) {
248 sigma = _c_D0D0->val()/picobarn;
249 error = _c_D0D0->err()/picobarn;
250 }
251 else if(iy==2) {
252 sigma = _c_D0D0S->val()/picobarn;
253 error = _c_D0D0S->err()/picobarn;
254 }
255 else if(iy==3) {
256 sigma = _c_D0SD0S->val()/picobarn;
257 error = _c_D0SD0S->err()/picobarn;
258 }
259 }
260 else if(ix==2) {
261 if(iy==1) {
262 sigma = _c_DpDm->val()/picobarn;
263 error = _c_DpDm->err()/picobarn;
264 }
265 else if(iy==2) {
266 sigma = _c_DpDmS->val()/picobarn;
267 error = _c_DpDmS->err()/picobarn;
268 }
269 else if(iy==3) {
270 sigma = _c_DpSDmS->val()/picobarn;
271 error = _c_DpSDmS->err()/picobarn;
272 }
273 }
274 else if(ix==3) {
275 if(iy==1) {
276 sigma = _c_DsDs->val()/picobarn;
277 error = _c_DsDs->err()/picobarn;
278 }
279 else if(iy==2) {
280 sigma = _c_DsDsS->val()/picobarn;
281 error = _c_DsDsS->err()/picobarn;
282 }
283 else if(iy==3) {
284 sigma = _c_DsSDsS->val()/picobarn;
285 error = _c_DsSDsS->err()/picobarn;
286 }
287 }
288 else if(ix==4) {
289 if(iy==1) {
290 sigma = _c_DSDpi->val()/picobarn;
291 error = _c_DSDpi->err()/picobarn;
292 }
293 else if(iy==2) {
294 sigma = _c_DSDSpi->val()/picobarn;
295 error = _c_DSDSpi->err()/picobarn;
296 }
297 }
298 else if(ix==5) {
299 if(iy==1) {
300 sigma = _c_DD->val()/nanobarn;
301 error = _c_DD->err()/nanobarn;
302 }
303 else if(iy==2) {
304 sigma = _c_DDX->val()/nanobarn;
305 error = _c_DDX->err()/nanobarn;
306 }
307 }
308 sigma *= crossSection()/picobarn/ sumOfWeights();
309 error *= crossSection()/picobarn/ sumOfWeights();
310 Estimate1DPtr mult;
311 book(mult, ix,1,iy);
312 for (auto& b : mult->bins()) {
313 if (inRange(sqrtS()/MeV, b.xMin(), b.xMax())) {
314 b.set(sigma, error);
315 }
316 }
317 }
318 }
319 }
320
321 /// @}
322
323
324 /// @name Histograms
325 /// @{
326 CounterPtr _c_D0D0, _c_DpDm,_c_DsDs;
327 CounterPtr _c_D0D0S, _c_DpDmS,_c_DsDsS;
328 CounterPtr _c_D0SD0S, _c_DpSDmS,_c_DsSDsS;
329 CounterPtr _c_DD, _c_DDX;
330 CounterPtr _c_DSDpi, _c_DSDSpi;
331 CounterPtr _c_hadrons, _c_muons;
332 /// @}
333
334
335 };
336
337
338 RIVET_DECLARE_PLUGIN(CLEOC_2008_I777917);
339
340
341}
|