Rivet analyses referenceBELLE_2007_I723333Cross section for $e^+e^-\to D^{*+}D^{*-}$ and $D^\pm D^{*\mp}$ below 5 GeVExperiment: BELLE (KEKB) Inspire ID: 723333 Status: VALIDATED Authors:
Beam energies: ANY Run details:
Cross section for $e^+e^-\to D^{*+}D^{*-}$ and $D^\pm D^{*\mp}$ below 5 GeV Source code: BELLE_2007_I723333.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 BELLE_2007_I723333 : public Analysis {
11 public:
12
13 /// Constructor
14 RIVET_DEFAULT_ANALYSIS_CTOR(BELLE_2007_I723333);
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(_nDSS, "/TMP/nDSS");
27 book(_nDS, "/TMP/nDS");
28 }
29
30 void findChildren(const Particle & p,map<long,int> & nRes, int &ncount) {
31 for(const Particle &child : p.children()) {
32 if(child.children().empty()) {
33 nRes[child.pid()]-=1;
34 --ncount;
35 }
36 else
37 findChildren(child,nRes,ncount);
38 }
39 }
40
41 /// Perform the per-event analysis
42 void analyze(const Event& event) {
43 const FinalState& fs = apply<FinalState>(event, "FS");
44
45 map<long,int> nCount;
46 int ntotal(0);
47 for (const Particle& p : fs.particles()) {
48 nCount[p.pid()] += 1;
49 ++ntotal;
50 }
51 const FinalState& ufs = apply<FinalState>(event, "UFS");
52
53 for(unsigned int ix=0;ix<ufs.particles().size();++ix) {
54 const Particle& p1 = ufs.particles()[ix];
55 if(abs(p1.pid())!=413) continue;
56 map<long,int> nRes = nCount;
57 int ncount = ntotal;
58 findChildren(p1,nRes,ncount);
59 bool matched=false;
60 int sign = -p1.pid()/abs(p1.pid());
61 for(unsigned int iy=0;iy<ufs.particles().size();++iy) {
62 if(ix==iy) continue;
63 const Particle& p2 = ufs.particles()[iy];
64 if(!p2.parents().empty() && p2.parents()[0].pid()==p1.pid())
65 continue;
66 if(p2.pid()!=sign*413 && p2.pid()!=sign*411) continue;
67 map<long,int> nRes2 = nRes;
68 int ncount2 = ncount;
69 findChildren(p2,nRes2,ncount2);
70 if(ncount2!=0) continue;
71 matched=true;
72 for(auto const & val : nRes2) {
73 if(val.second!=0) {
74 matched = false;
75 break;
76 }
77 }
78 if(matched) {
79 sign = abs(p2.pid());
80 break;
81 }
82 }
83 if(matched) {
84 if(sign==411)
85 _nDS->fill();
86 else if(sign==413)
87 _nDSS->fill();
88 }
89 }
90 }
91
92
93 /// Normalise histograms etc., after the run
94 void finalize() {
95 for(unsigned int ix=1;ix<3;++ix) {
96 double sigma,error;
97 if(ix==1) {
98 sigma = _nDSS->val();
99 error = _nDSS->err();
100 }
101 else {
102 sigma = _nDS->val();
103 error = _nDS->err();
104 }
105 sigma *= crossSection()/ sumOfWeights() /nanobarn;
106 error *= crossSection()/ sumOfWeights() /nanobarn;
107 Scatter2D temphisto(refData(ix, 1, 1));
108 Scatter2DPtr mult;
109 book(mult, ix, 1, 1);
110 for (size_t b = 0; b < temphisto.numPoints(); b++) {
111 const double x = temphisto.point(b).x();
112 pair<double,double> ex = temphisto.point(b).xErrs();
113 pair<double,double> ex2 = ex;
114 if(ex2.first ==0.) ex2. first=0.0001;
115 if(ex2.second==0.) ex2.second=0.0001;
116 if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) {
117 mult->addPoint(x, sigma, ex, make_pair(error,error));
118 }
119 else {
120 mult->addPoint(x, 0., ex, make_pair(0.,.0));
121 }
122 }
123 }
124 }
125
126 //@}
127
128
129 /// @name Histograms
130 //@{
131 CounterPtr _nDSS,_nDS;
132 //@}
133
134
135 };
136
137
138 // The hook for the plugin system
139 RIVET_DECLARE_PLUGIN(BELLE_2007_I723333);
140
141
142}
|