Rivet analyses referenceBESIII_2016_I1384778Azimuthal asymmetries in inclusive charged pion-pair production at $\sqrt{s}=3.65$ GeVExperiment: BESIII (BEPC) Inspire ID: 1384778 Status: VALIDATED Authors:
Beam energies: (1.8, 1.8) GeV Run details:
Measurement of azimuthal asymmetries in inclusive charged pion-pair production at $\sqrt{s}=3.65$ GeV by the BESII experiment Source code: BESIII_2016_I1384778.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/Beam.hh"
5
6namespace Rivet {
7
8
9 /// @brief Collins assymmetry
10 class BESIII_2016_I1384778 : public Analysis {
11 public:
12
13 /// Constructor
14 RIVET_DEFAULT_ANALYSIS_CTOR(BESIII_2016_I1384778);
15
16
17 /// @name Analysis methods
18 /// @{
19
20 /// Book histograms and initialise projections before the run
21 void init() {
22 declare(Beam(), "Beams");
23 declare(FinalState(Cuts::abspid==PID::PIPLUS), "FS");
24 // book the histograms
25 _h_L.resize(6);
26 _h_U.resize(6);
27 _h_C.resize(6);
28 for(size_t ix=0; _h_L.size(); ++ix) {
29 const string pre = "/TMP/h_z1z2_"+to_string(ix+1);
30 book(_h_L[ix], pre+"_L", 20, 0., M_PI);
31 book(_h_U[ix], pre+"_U", 20, 0., M_PI);
32 book(_h_C[ix], pre+"_C", 20, 0., M_PI);
33 }
34 book(_h_pT_L, {0., 0.2, 0.3, 0.45, 0.8, 1.4});
35 book(_h_pT_U, {0., 0.2, 0.3, 0.45, 0.8, 1.4});
36 book(_h_pT_C, {0., 0.2, 0.3, 0.45, 0.8, 1.4});
37 for (size_t ix=1; ix < _h_pT_L->numBins()+1; ++ix) {
38 const string pre = "/TMP/h_pT_"+to_string(ix);
39 book(_h_pT_L->bin(ix), pre+"_L", 20, 0.0, M_PI);
40 book(_h_pT_U->bin(ix), pre+"_U", 20, 0.0, M_PI);
41 book(_h_pT_C->bin(ix), pre+"_C", 20, 0.0, M_PI);
42 }
43 }
44
45
46 /// Perform the per-event analysis
47 void analyze(const Event& event) {
48 // get the axis, direction of incoming electron
49 const ParticlePair& beams = apply<Beam>(event, "Beams").beams();
50 Vector3 axis;
51 if (beams.first.pid()>0) {
52 axis = beams.first .momentum().p3().unit();
53 }
54 else {
55 axis = beams.second.momentum().p3().unit();
56 }
57 // loop over pions pair, using index to avoid double counting
58 Particles pions = apply<FinalState>(event, "FS").particles();
59 for (unsigned int i1=0; i1<pions.size(); ++i1) {
60 const double x1=2.*pions[i1].momentum().t()/sqrtS();
61 // cut on z1
62 if (x1<0.2||x1>0.9) continue;
63 // cos theta cut
64 if (abs(cos(pions[i1].momentum().p3().polarAngle()))>0.93) continue;
65 for(unsigned int i2=i1+1;i2<pions.size();++i2) {
66 // cut on z2
67 const double x2=2.*pions[i2].momentum().t()/sqrtS();
68 if (x2<0.2||x2>0.9) continue;
69 // cos theta cut
70 if (abs(cos(pions[i2].momentum().p3().polarAngle()))>0.93) continue;
71 // cut on opening angle (>120 degrees)
72 if (pions[i1].momentum().p3().angle(pions[i2].momentum().p3())>2.*M_PI/3.) {
73 continue;
74 }
75 Particle p1=pions[i1], p2=pions[i2];
76 double z1(x1),z2(x2);
77 // randomly order the particles
78 if (rand()/static_cast<double>(RAND_MAX) < 0.5) {
79 swap(p1,p2);
80 swap(z1,z2);
81 }
82 // particle 2 defines the z axis
83 Vector3 ez = p2.momentum().p3().unit();
84 // beam and 2 define the plane (y is normal to plane)
85 Vector3 ey = ez.cross(axis).unit();
86 // x by cross product
87 Vector3 ex = ey.cross(ez).unit();
88 // phi
89 double phi = ex.angle(p1.momentum().p3());
90 // hists vs z1,z2
91 unsigned int ibin=0;
92 if (z1<=.3&&z2<=.3) {
93 ibin=0;
94 }
95 else if (z1>0.5&&z2>0.5) {
96 ibin=5;
97 }
98 else if (min(z1,z2)<=0.3) {
99 if (max(z1,z2)>0.5) {
100 ibin=2;
101 }
102 else {
103 ibin=1;
104 }
105 }
106 else {
107 if (max(z1,z2)>0.5) {
108 ibin=4;
109 }
110 else {
111 ibin=3;
112 }
113 }
114 _h_C[ibin]->fill(phi);
115 if (p1.pid()==p2.pid()) {
116 _h_L[ibin]->fill(phi);
117 }
118 else {
119 _h_U[ibin]->fill(phi);
120 }
121 // hists vs pT
122 double pPar2 = sqr(ez.dot(p1.momentum().p3()));
123 double pPerp = sqrt(p1.momentum().p3().mod2()-pPar2);
124 _h_pT_C->fill(pPerp,phi);
125 if(p1.pid()==p2.pid()) {
126 _h_pT_L->fill(pPerp,phi);
127 }
128 else {
129 _h_pT_U->fill(pPerp,phi);
130 }
131 }
132 }
133 }
134
135 pair<double,double> calcAsymmetry(Estimate1DPtr hist) {
136 double sum1(0.),sum2(0.);
137 for (const auto& b : hist->bins()) {
138 double Oi = b.val();
139 if(Oi==0. || std::isnan(Oi) ) continue;
140 double ai = 1.;
141 double bi = 0.5*(sin(2.*b.xMax())-sin(2.*b.xMin()))/b.xWidth();
142 double Ei = b.errAvg();
143 sum1 += sqr(bi/Ei);
144 sum2 += bi/sqr(Ei)*(Oi-ai);
145 }
146 return make_pair(sum2/sum1, sqrt(1./sum1));
147 }
148
149
150 /// Normalise histograms etc., after the run
151 void finalize() {
152 // ratios
153 Estimate1DPtr _h_z_UL,_h_z_UC;
154 book(_h_z_UL,1,1,5);
155 book(_h_z_UC,1,1,6);
156 for (size_t ix=0; ix < _h_L.size(); ++ix) {
157 normalize(_h_L[ix]);
158 normalize(_h_U[ix]);
159 normalize(_h_C[ix]);
160 const string pre = "/TMP/R_z1z2_"+to_string(ix+1);
161 Estimate1DPtr R1;
162 book(R1, pre+"_UL");
163 divide(_h_U[ix], _h_L[ix], R1);
164 Estimate1DPtr R2;
165 book(R2, pre+"_UC");
166 divide(_h_U[ix], _h_C[ix], R2);
167 pair<double,double> asym1 = calcAsymmetry(R1);
168 _h_z_UL->bin(ix+1).set(asym1.first, asym1.second);
169 pair<double,double> asym2 = calcAsymmetry(R2);
170 _h_z_UC->bin(ix+1).set(asym2.first, asym2.second);
171 }
172 Estimate1DPtr _h_pT_UL,_h_pT_UC;
173 book(_h_pT_UL,2,1,4);
174 book(_h_pT_UC,2,1,5);
175 for (size_t ix=1; ix < _h_pT_L->numBins()+1; ++ix) {
176 normalize(_h_pT_L->bin(ix));
177 normalize(_h_pT_U->bin(ix));
178 normalize(_h_pT_C->bin(ix));
179 const string pre = "/TMP/R_pT_"+to_string(ix);
180 Estimate1DPtr R1;
181 book(R1, pre+"_UL");
182 divide(_h_pT_U->bin(ix), _h_pT_L->bin(ix), R1);
183 Estimate1DPtr R2;
184 book(R2, pre+"_UC");
185 divide(_h_pT_U->bin(ix), _h_pT_C->bin(ix), R2);
186 pair<double,double> asym1 = calcAsymmetry(R1);
187 _h_pT_UL->bin(ix).set(asym1.first, asym1.second);
188 pair<double,double> asym2 = calcAsymmetry(R2);
189 _h_pT_UC->bin(ix).set(asym2.first, asym2.second);
190 }
191 }
192 /// @}
193
194
195 /// @name Histograms
196 /// @{
197 vector<Histo1DPtr> _h_L,_h_U,_h_C;
198 Histo1DGroupPtr _h_pT_L,_h_pT_U,_h_pT_C;
199 /// @}
200
201
202 };
203
204
205 RIVET_DECLARE_PLUGIN(BESIII_2016_I1384778);
206
207}
|