Rivet analyses referenceBELLE_2017_I1607562Invariant-mass and fractional-energy dependence of inclusive production of di-hadrons at $\sqrt{s}=10.58$ GeVExperiment: BELLE (KEKB) Inspire ID: 1607562 Status: VALIDATED Authors:
Beam energies: (5.3, 5.3) GeV Run details:
Measurement of the double differential cross secrion for the production of hadron pairs by BELLE Source code: BELLE_2017_I1607562.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/Thrust.hh"
5#include "Rivet/Tools/BinnedHistogram.hh"
6
7namespace Rivet {
8
9
10 /// @brief BELLE double differential cross section
11 class BELLE_2017_I1607562 : public Analysis {
12 public:
13
14 /// Constructor
15 RIVET_DEFAULT_ANALYSIS_CTOR(BELLE_2017_I1607562);
16
17
18 /// @name Analysis methods
19 ///@{
20
21 /// Book histograms and initialise projections before the run
22 void init() {
23 // projections
24 FinalState fs;
25 declare(fs,"FS");
26 declare(Thrust(fs),"Thrust");
27 // histograms
28 double bins[17]={0.20,0.25,0.30,0.35,0.40,0.45,0.50,0.55,0.60,
29 0.65,0.70,0.75,0.80,0.85,0.90,0.95,1.00};
30 for(unsigned int ip=0;ip<6;++ip) {
31 for(unsigned int iy=0;iy<16;++iy) {
32 Histo1DPtr temp;
33 book(temp,1,ip+1,iy+1);
34 _h_all [ip].add(bins[iy],bins[iy+1],temp);
35 book(temp,2,ip+1,iy+1);
36 _h_strong[ip].add(bins[iy],bins[iy+1],temp);
37 }
38 }
39 }
40
41 bool isWeak(const Particle & p) {
42 bool weak = false;
43 if(p.parents().empty()) return weak;
44 Particle parent = p.parents()[0];
45 while (!parent.parents().empty()) {
46 if(parent.abspid()==411 || parent.abspid()==421 || parent.abspid()==431 ||
47 parent.abspid()==4122 || parent.abspid()==4232 || parent.abspid()==4132 ||
48 parent.abspid()==4332) {
49 weak=true;
50 break;
51 }
52 parent = parent.parents()[0];
53 }
54 return weak;
55 }
56
57 /// Perform the per-event analysis
58 void analyze(const Event& event) {
59 // get thrust and apply cut
60 const Thrust thrust = apply<Thrust>(event,"Thrust");
61 if(thrust.thrust()<0.8) vetoEvent;
62 // get thrust axis
63 Vector3 axis = thrust.thrustAxis();
64 Particles charged = apply<FinalState>(event,"FS").particles(Cuts::abspid==PID::KPLUS or
65 Cuts::abspid==PID::PIPLUS);
66 for(unsigned int ix=0;ix<charged.size();++ix) {
67 double dot1 = axis.dot(charged[ix].momentum().p3());
68 bool weak1 = isWeak(charged[ix]);
69 if(2.*charged[ix].momentum().t()/sqrtS()<0.1) continue;
70 for(unsigned int iy=ix+1;iy<charged.size();++iy) {
71 if(2.*charged[iy].momentum().t()/sqrtS()<0.1) continue;
72 double dot2 = axis.dot(charged[iy].momentum().p3());
73 bool weak2 = isWeak(charged[iy]);
74 if(dot1*dot2<0.) continue;
75 FourMomentum p = charged[ix].momentum()+charged[iy].momentum();
76 double z12 = 2.*p.t()/sqrtS();
77 double m12 = p.mass();
78 bool strong = !weak1 && !weak2;
79 if(charged[ix].pid()==PID::PIPLUS) {
80 if(charged[iy].pid()==PID::PIPLUS) {
81 _h_all[1].fill(z12,m12);
82 if(strong) _h_strong[1].fill(z12,m12);
83 }
84 else if(charged[iy].pid()==PID::PIMINUS) {
85 _h_all[0].fill(z12,m12);
86 if(strong) _h_strong[0].fill(z12,m12);
87 }
88 else if(charged[iy].pid()==PID::KPLUS) {
89 _h_all[3].fill(z12,m12);
90 if(strong) _h_strong[3].fill(z12,m12);
91 }
92 else if(charged[iy].pid()==PID::KMINUS) {
93 _h_all[2].fill(z12,m12);
94 if(strong) _h_strong[2].fill(z12,m12);
95 }
96 }
97 else if(charged[ix].pid()==PID::PIMINUS) {
98 if(charged[iy].pid()==PID::PIPLUS) {
99 _h_all[0].fill(z12,m12);
100 if(strong) _h_strong[0].fill(z12,m12);
101 }
102 else if(charged[iy].pid()==PID::PIMINUS) {
103 _h_all[1].fill(z12,m12);
104 if(strong) _h_strong[1].fill(z12,m12);
105 }
106 else if(charged[iy].pid()==PID::KPLUS) {
107 _h_all[2].fill(z12,m12);
108 if(strong) _h_strong[2].fill(z12,m12);
109 }
110 else if(charged[iy].pid()==PID::KMINUS) {
111 _h_all[3].fill(z12,m12);
112 if(strong) _h_strong[3].fill(z12,m12);
113 }
114 }
115 else if(charged[ix].pid()==PID::KPLUS) {
116 if(charged[iy].pid()==PID::PIPLUS) {
117 _h_all[3].fill(z12,m12);
118 if(strong) _h_strong[3].fill(z12,m12);
119 }
120 else if(charged[iy].pid()==PID::PIMINUS) {
121 _h_all[2].fill(z12,m12);
122 if(strong) _h_strong[2].fill(z12,m12);
123 }
124 else if(charged[iy].pid()==PID::KPLUS) {
125 _h_all[5].fill(z12,m12);
126 if(strong) _h_strong[5].fill(z12,m12);
127 }
128 else if(charged[iy].pid()==PID::KMINUS) {
129 _h_all[4].fill(z12,m12);
130 if(strong) _h_strong[4].fill(z12,m12);
131 }
132 }
133 else if(charged[ix].pid()==PID::KMINUS) {
134 if(charged[iy].pid()==PID::PIPLUS) {
135 _h_all[2].fill(z12,m12);
136 if(strong) _h_strong[2].fill(z12,m12);
137 }
138 else if(charged[iy].pid()==PID::PIMINUS) {
139 _h_all[3].fill(z12,m12);
140 if(strong) _h_strong[3].fill(z12,m12);
141 }
142 else if(charged[iy].pid()==PID::KPLUS) {
143 _h_all[4].fill(z12,m12);
144 if(strong) _h_strong[4].fill(z12,m12);
145 }
146 else if(charged[iy].pid()==PID::KMINUS) {
147 _h_all[5].fill(z12,m12);
148 if(strong) _h_strong[5].fill(z12,m12);
149 }
150 }
151 }
152 }
153 }
154
155
156 /// Normalise histograms etc., after the run
157 void finalize() {
158 double fact = crossSection()/nanobarn/sumOfWeights();
159 for(unsigned int ix=0;ix<6;++ix) {
160 _h_all [ix].scale(fact,this);
161 _h_strong[ix].scale(fact,this);
162 }
163 }
164
165 ///@}
166
167
168 /// @name Histograms
169 ///@{
170 BinnedHistogram _h_all[6],_h_strong[6];
171 ///@}
172
173
174 };
175
176
177 RIVET_DECLARE_PLUGIN(BELLE_2017_I1607562);
178
179}
|