Rivet analyses referenceLENA_1981_I164397Charged particle multiplicities and thrust in $\Upsilon(1S,2S)$ decays and nearby continuumExperiment: LENA (DORIS) Inspire ID: 164397 Status: VALIDATED Authors:
Beams: * * Beam energies: ANY Run details:
Measurement of the charged particle multiplicities and a thrust-like variable in $\Upsilon(1S,2S)$ decays and nearby continuum. As LENA had no magnetic field the momenta of the charged particles were not measured and therefore a thrust-like variable $T^\prime$ was constructed using only the directions of the charged particles. Beam energy must be specified as analysis option "ENERGY" when rivet-merging samples. Source code: LENA_1981_I164397.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/UnstableParticles.hh"
4#include "Rivet/Projections/ChargedFinalState.hh"
5#include "Rivet/Projections/Thrust.hh"
6
7namespace Rivet {
8
9
10 /// @brief Thrust like variable at Upsilon(1s,2S)
11 class LENA_1981_I164397 : public Analysis {
12 public:
13
14 /// Constructor
15 RIVET_DEFAULT_ANALYSIS_CTOR(LENA_1981_I164397);
16
17
18 /// @name Analysis methods
19 /// @{
20
21 /// Book histograms and initialise projections before the run
22 void init() {
23 declare(UnstableParticles(), "UFS");
24 declare(ChargedFinalState(), "FS");
25
26 book(_weightSum_cont, "TMP/weightSum_cont");
27 book(_weightSum_Ups1, "TMP/weightSum_Ups1");
28 book(_weightSum_Ups2, "TMP/weightSum_Ups2");
29 book(_charge_cont, "TMP/charge_cont");
30 book(_charge_Ups1, "TMP/charge_Ups1");
31 book(_charge_Ups2, "TMP/charge_Ups2");
32
33 if(isCompatibleWithSqrtS(9.5149*GeV,1e-2)) {
34 book(_hist_T_cont ,4, 1, 1);
35 }
36 else if(isCompatibleWithSqrtS(9.9903*GeV,1e-2)) {
37 book(_hist_T_cont ,4, 1, 2);
38 }
39 book(_hist_T_Ups1 ,4, 1, 3);
40 book(_hist_T_Ups2 ,4, 1, 4);
41 }
42
43 /// Recursively walk the decay tree to find the charged decay products of @a p
44 void findDecayProducts(Particle mother, Particles& charged) {
45 for(const Particle & p: mother.children()) {
46 const int id = p.pid();
47 if(!p.children().empty())
48 findDecayProducts(p, charged);
49 else if(PID::isCharged(id))
50 charged.push_back(p);
51 }
52 }
53
54 // defn of thrust in paper used just the direction
55 double thrustPrime(const LorentzTransform & boost, const Particles & particles) {
56 vector<Vector3> vecs;
57 for(const Particle & p : particles) {
58 vecs.push_back(boost.transform(p.momentum()).p3().unit());
59 }
60 Thrust thrust;
61 thrust.calc(vecs);
62 return thrust.thrust();
63 }
64
65 /// Perform the per-event analysis
66 void analyze(const Event& event) {
67 // Find the Upsilons among the unstables
68 const UnstableParticles& ufs = apply<UnstableParticles>(event, "UFS");
69 Particles upsilons = ufs.particles(Cuts::pid==553 or Cuts::pid==100553);
70 if (upsilons.empty()) {
71 MSG_DEBUG("No Upsilons found => continuum event");
72 _weightSum_cont->fill();
73 Particles cfs = apply<ChargedFinalState>(event, "FS").particles();
74 _charge_cont->fill(cfs.size());
75 if(_hist_T_cont) {
76 LorentzTransform boost;
77 _hist_T_cont->fill(thrustPrime(boost,cfs));
78 }
79 }
80 // Upsilon(s) found
81 else {
82 for (const Particle& ups : upsilons) {
83 const int parentId = ups.pid();
84 Particles charged;
85 // boost to rest frame (if required)
86 LorentzTransform boost;
87 if (ups.p3().mod() > 1*MeV)
88 boost = LorentzTransform::mkFrameTransformFromBeta(ups.momentum().betaVec());
89 // Find the decay products we want
90 findDecayProducts(ups, charged);
91 if(parentId==553) {
92 _weightSum_Ups1->fill();
93 _charge_Ups1->fill(charged.size());
94 _hist_T_Ups1->fill(thrustPrime(boost,charged));
95 }
96 else {
97 _weightSum_Ups2->fill();
98 _charge_Ups2->fill(charged.size());
99 _hist_T_Ups2->fill(thrustPrime(boost,charged));
100 }
101 }
102 }
103
104 }
105
106
107 /// Normalise histograms etc., after the run
108 void finalize() {
109 // charged particle multiplicity
110 if(_weightSum_cont->val()>0. ) {
111 scale(_charge_cont,1./ *_weightSum_cont );
112 if(_hist_T_cont) scale(_hist_T_cont,1./ *_weightSum_cont );
113 }
114 if(_weightSum_Ups1->val()>0. ) {
115 scale(_charge_Ups1,1./ *_weightSum_Ups1 );
116 scale(_hist_T_Ups1,1./ *_weightSum_Ups1 );
117 }
118 if(_weightSum_Ups2->val()>0. ) {
119 scale(_charge_Ups2,1./ *_weightSum_Ups2 );
120 scale(_hist_T_Ups2,1./ *_weightSum_Ups2 );
121 }
122 Estimate1DPtr _mult;
123 book(_mult, 3, 1, 1);
124 for (auto& b : _mult->bins()) {
125 // Upsilon 1S
126 if(b.index()==4) {
127 if (_weightSum_Ups1->val()>0.) {
128 b.set(_charge_Ups1->val(), _charge_Ups1->err());
129 }
130 }
131 // Upsilon 2S
132 else if(b.index()==7) {
133 if (_weightSum_Ups2->val()>0.) {
134 b.set(_charge_Ups2->val(), _charge_Ups2->err());
135 }
136 }
137 else {
138 if (inRange(sqrtS()/GeV, b.xMin(), b.xMax()) && _weightSum_cont->val()>0.) {
139 b.set(_charge_cont->val(), _charge_cont->err());
140 }
141 }
142 }
143 }
144
145 /// @}
146
147
148 /// @name Histograms
149 /// @{
150 CounterPtr _weightSum_cont, _weightSum_Ups1, _weightSum_Ups2;
151 CounterPtr _charge_cont, _charge_Ups1, _charge_Ups2;
152 Histo1DPtr _hist_T_cont,_hist_T_Ups1,_hist_T_Ups2;
153 /// @}
154
155
156 };
157
158
159 RIVET_DECLARE_PLUGIN(LENA_1981_I164397);
160
161}
|