Rivet analyses referenceARGUS_1988_I266892Baryon-antibaryon correlations at the $\Upsilon(1S)$ and nearby continuumExperiment: ARGUS (DORIS) Inspire ID: 266892 Status: VALIDATED NOHEPDATA Authors:
Beam energies: (4.7, 4.7); (5.0, 5.0) GeV Run details:
Correlations in the production of $p\bar{p}$ and hyperon-antihyperon pairs. Source code: ARGUS_1988_I266892.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/Thrust.hh"
4#include "Rivet/Projections/FinalState.hh"
5#include "Rivet/Projections/UnstableParticles.hh"
6
7namespace Rivet {
8
9
10 /// @brief baryon correlations
11 class ARGUS_1988_I266892 : public Analysis {
12 public:
13
14 /// Constructor
15 RIVET_DEFAULT_ANALYSIS_CTOR(ARGUS_1988_I266892);
16
17
18 /// @name Analysis methods
19 /// @{
20
21 /// Book histograms and initialise projections before the run
22 void init() {
23 // projections
24 const FinalState fs;
25 declare(fs , "FS" );
26 declare(UnstableParticles(), "UFS" );
27 declare(Thrust(fs) , "Thrust");
28 // histos
29 for (unsigned int ix=0; ix<3; ++ix) {
30 book(_nBB2[ix], "TMP/nBB2_"+toString(ix+1),
31 refData<YODA::BinnedEstimate<string>>(4+2*ix, 1, 1));
32 book(_n[ix], "TMP/n_"+toString(ix+1),
33 refData<YODA::BinnedEstimate<string>>(4+2*ix, 1, 1));
34 book(_nBB [ix],3+2*ix,1,1);
35 book(_c[ix],"TMP/c_"+toString(ix+1));
36 if (ix==2) continue;
37 book(_h_pt[ix],1,1,1+ix);
38 }
39 book(_h_pt[2],2,1,1);
40 }
41
42 /// Recursively walk the decay tree to find the stable decay products of @a p
43 void findDecayProducts(const Particle& mother, Particles& final,
44 Particles& protons, Particles& lambda, Particles& xi,
45 Particles& l1520) {
46 for (const Particle & p: mother.children()) {
47 if (p.abspid()==3122 ) lambda.push_back(p);
48 else if (p.abspid()==3312 ) xi .push_back(p);
49 else if (p.abspid()==102134) l1520 .push_back(p);
50 if (!p.children().empty()) {
51 findDecayProducts(p, final, protons, lambda, xi, l1520);
52 }
53 else {
54 final.push_back(p);
55 if (p.abspid()==2212) protons.push_back(p);
56 }
57 }
58 }
59
60
61 /// Perform the per-event analysis
62 void analyze(const Event& event) {
63 if(_edges.empty()) _edges = _nBB2[1]->xEdges();
64 // Find the Upsilons among the unstables
65 const UnstableParticles& ufs = apply<UnstableParticles>(event, "UFS");
66 Particles upsilons = ufs.particles(Cuts::pid==553);
67 if (upsilons.empty()) {
68 _c[1]->fill();
69 _c[2]->fill();
70 Thrust thrust = apply<Thrust>(event, "Thrust");
71 Vector3 axis = thrust.thrustAxis();
72 Particles protons = apply<FinalState>(event, "FS").particles(Cuts::abspid==2212);
73 for (unsigned int ix=0;ix<protons.size();++ix) {
74 double d1 = protons[ix].mom().p3().dot(axis);
75 Vector3 pt1 = protons[ix].mom().p3()-d1*axis;
76 for (unsigned int iy=ix+1;iy<protons.size();++iy) {
77 if (protons[ix].pid()*protons[iy].pid()>0) continue;
78 double d2 = protons[iy].mom().p3().dot(axis);
79 Vector3 pt2 = protons[iy].mom().p3()-d2*axis;
80 double phi = pt1.angle(pt2);
81 if (phi<0.) phi +=M_PI;
82 phi *=180./M_PI;
83 if(phi==180.) phi=phi-1e-10;
84 if (d1*d2>0) _h_pt[0]->fill(phi);
85 else _h_pt[2]->fill(phi);
86 }
87 }
88 Particles lambda = ufs.particles(Cuts::abspid==3122);
89 for (unsigned int ix=0; ix<lambda.size(); ++ix) {
90 _n[0]->fill(_edges[1]);
91 for (unsigned int iy=ix+1; iy<lambda.size(); ++iy) {
92 if (lambda[ix].pid()*lambda[iy].pid()<0) {
93 _nBB [0]->fill(_edges[1]);
94 _nBB2[0]->fill(_edges[1]);
95 }
96 }
97 }
98 Particles xi = ufs.particles(Cuts::abspid==3312);
99 for (unsigned int ix=0;ix<xi.size();++ix) {
100 _n[1]->fill(_edges[1]);
101 _n[1]->fill(_edges[2]);
102 for (unsigned int iy=0;iy<lambda.size();++iy) {
103 if (xi[ix].pid()*lambda[iy].pid()<0) {
104 _nBB [1]->fill(_edges[1]);
105 _nBB2[1]->fill(_edges[1]);
106 _nBB2[1]->fill(_edges[2]);
107 }
108 }
109 }
110 Particles l1520 = ufs.particles(Cuts::abspid==102134);
111 for (unsigned int ix=0; ix<l1520.size(); ++ix) {
112 _n[2]->fill(_edges[1]);
113 _n[2]->fill(_edges[2]);
114 for (unsigned int iy=0; iy<lambda.size(); ++iy) {
115 if (l1520[ix].pid()*lambda[iy].pid()<0) {
116 _nBB [2]->fill(_edges[1]);
117 _nBB2[2]->fill(_edges[1]);
118 _nBB2[2]->fill(_edges[2]);
119 }
120 }
121 }
122 }
123 else {
124 for (const Particle& ups : upsilons) {
125 _c[0]->fill();
126 _c[2]->fill();
127 LorentzTransform boost;
128 if (ups.p3().mod() > 1*MeV) {
129 boost = LorentzTransform::mkFrameTransformFromBeta(ups.mom().betaVec());
130 }
131 Particles final, protons, lambda, xi, l1520;
132 findDecayProducts(ups, final, protons, lambda, xi, l1520);
133 vector<FourMomentum> mom; mom.reserve(final.size());
134 for (const Particle& p : final) {
135 mom.push_back(boost.transform(p.mom()));
136 }
137 Thrust thrust;
138 thrust.calc(mom);
139 Vector3 axis = thrust.thrustAxis();
140 for (unsigned int ix=0; ix<protons.size(); ++ix) {
141 double d1 = protons[ix].mom().p3().dot(axis);
142 Vector3 pt1 = protons[ix].mom().p3()-d1*axis;
143 for (unsigned int iy=ix+1;iy<protons.size();++iy) {
144 if (protons[ix].pid()*protons[iy].pid()>0) continue;
145 double d2 = protons[iy].mom().p3().dot(axis);
146 Vector3 pt2 = protons[iy].mom().p3()-d2*axis;
147 double phi = pt1.angle(pt2);
148 if (phi<0.) phi +=M_PI;
149 phi *=180./M_PI;
150 if(phi==180.) phi=phi-1e-10;
151 if (d1*d2>0) _h_pt[1]->fill(phi);
152 }
153 }
154 for (unsigned int ix=0; ix<lambda.size(); ++ix) {
155 _n[0]->fill(_edges[0]);
156 for (unsigned int iy=ix+1; iy<lambda.size(); ++iy) {
157 if (lambda[ix].pid()*lambda[iy].pid()<0) {
158 _nBB [0]->fill(_edges[0]);
159 _nBB2[0]->fill(_edges[0]);
160 }
161 }
162 }
163 for (unsigned int ix=0; ix<xi.size(); ++ix) {
164 _n[1]->fill(_edges[0]);
165 _n[1]->fill(_edges[2]);
166 for (unsigned int iy=0; iy<lambda.size(); ++iy) {
167 if(xi[ix].pid()*lambda[iy].pid()<0) {
168 _nBB [1]->fill(_edges[0]);
169 _nBB2[1]->fill(_edges[0]);
170 _nBB2[1]->fill(_edges[2]);
171 }
172 }
173 }
174 for (unsigned int ix=0; ix<l1520.size(); ++ix) {
175 _n[2]->fill(_edges[0]);
176 _n[2]->fill(_edges[2]);
177 for (unsigned int iy=0; iy<lambda.size(); ++iy) {
178 if (l1520[ix].pid()*lambda[iy].pid()<0) {
179 _nBB [2]->fill(_edges[0]);
180 _nBB2[2]->fill(_edges[0]);
181 _nBB2[2]->fill(_edges[2]);
182 }
183 }
184 }
185 }
186 }
187 }
188
189
190 /// Normalise histograms etc., after the run
191 void finalize() {
192 normalize(_h_pt, 1.0, false);
193 for (unsigned int ix=0; ix<3; ++ix) {
194 for(unsigned int iy=0;iy<_nBB[ix]->numBins();++iy)
195 _nBB[ix]->bin(iy+1).scaleW( 1./_c[iy]->val());
196 BinnedEstimatePtr<string> tmp;
197 book(tmp,4+2*ix, 1, 1);
198 divide(_nBB2[ix], _n[ix],tmp);
199 if(ix==0) scale(tmp,2.);
200 }
201 }
202
203 /// @}
204
205
206 /// @name Histograms
207 /// @{
208 Histo1DPtr _h_pt[3];
209 BinnedHistoPtr<string> _nBB[3],_nBB2[3],_n[3];
210 CounterPtr _c[3];
211 vector<string> _edges;
212 /// @}
213
214
215 };
216
217
218 RIVET_DECLARE_PLUGIN(ARGUS_1988_I266892);
219
220}
|
