Rivet analyses referenceL3_1998_I467929$\tau$ polarization at LEP1Experiment: L3 (LEP) Inspire ID: 467929 Status: VALIDATED Authors:
Beam energies: (45.6, 45.6) GeV Run details:
Measurement of the $\tau$ lepton polarization in $e^+e^-\to\tau^+\tau^-$ at the $Z^0$ pole by the L3 experiment at LEP1. Source code: L3_1998_I467929.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/Beam.hh"
4#include "Rivet/Projections/ChargedFinalState.hh"
5#include "Rivet/Projections/UnstableParticles.hh"
6
7namespace Rivet {
8
9
10 /// @brief e+e- > tau+ tau-
11 class L3_1998_I467929 : public Analysis {
12 public:
13
14 /// Constructor
15 RIVET_DEFAULT_ANALYSIS_CTOR(L3_1998_I467929);
16
17
18 /// @name Analysis methods
19 /// @{
20
21 /// Book histograms and initialise projections before the run
22 void init() {
23 // Initialise and register projections
24 declare(Beam(), "Beams");
25 declare(ChargedFinalState(), "FS");
26 declare(UnstableParticles(), "UFS");
27 // book hists
28 const vector<double> edges = {-0.94, -0.83, -0.72, -0.55,-0.35,-0.12, 0.12, 0.35, 0.55, 0.72, 0.83, 0.94};
29 book(_h_e, edges); _h_e->maskBins({2, 10});
30 book(_h_mu, edges); _h_mu->maskBins({2, 10});
31 book(_h_pi, edges); _h_pi->maskBins({2, 10});
32 book(_h_rho, edges); _h_rho->maskBins({2, 10});
33 unsigned int iy=1;
34 for (size_t ix=0; ix < _h_e->numBins(); ++ix) {
35 const string suff = std::to_string(ix);
36 if (iy==2 || iy==10) ++iy;
37 book(_h_e->bin(iy), "_h_e_"+suff, 20, -1.0, 1.0);
38 book(_h_mu->bin(iy), "_h_mu_"+suff, 20, -1.0, 1.0);
39 book(_h_pi->bin(iy), "_h_pi_"+suff, 20, -1.0, 1.0);
40 book(_h_rho->bin(iy), "_h_rho_"+suff, 20, -1.0, 1.0);
41 iy+=1;
42 }
43 }
44
45 void findTau(const Particle & p, unsigned int & nprod, Particles& piP,
46 Particles& pi0, Particles& ell, Particles& nu_ell, Particles& nu_tau) {
47 for (const Particle& child : p.children()) {
48 if(child.pid()==PID::ELECTRON || child.pid()==PID::MUON) {
49 ++nprod;
50 ell.push_back(child);
51 }
52 else if(child.pid()==PID::NU_EBAR || child.pid()==PID::NU_MUBAR) {
53 ++nprod;
54 nu_ell.push_back(child);
55 }
56 else if(child.pid()==PID::PIMINUS) {
57 ++nprod;
58 piP.push_back(child);
59 }
60 else if(child.pid()==PID::PI0) {
61 ++nprod;
62 pi0.push_back(child);
63 }
64 else if(child.pid()==PID::NU_TAU) {
65 ++nprod;
66 nu_tau.push_back(child);
67 }
68 else if(child.pid()==PID::GAMMA) {
69 continue;
70 }
71 else if(child.children().empty() || child.pid()==221 || child.pid()==331) {
72 ++nprod;
73 }
74 else {
75 findTau(child,nprod,piP,pi0,ell,nu_ell,nu_tau);
76 }
77 }
78 }
79
80 /// Perform the per-event analysis
81 void analyze(const Event& event) {
82 // require 2 chanrged particles to veto hadronic events
83 if(apply<ChargedFinalState>(event, "FS").particles().size()!=2) vetoEvent;
84 // Get beams and average beam momentum
85 const ParticlePair& beams = apply<Beam>(event, "Beams").beams();
86 Vector3 axis;
87 if (beams.first.pid()>0) {
88 axis = beams.first .momentum().p3().unit();
89 }
90 else {
91 axis = beams.second.momentum().p3().unit();
92 }
93 // loop over tau leptons
94 for (const Particle& p : apply<UnstableParticles>(event, "UFS").particles(Cuts::pid==15)) {
95 unsigned int nprod(0);
96 Particles piP, pi0, ell, nu_ell, nu_tau;
97 findTau(p,nprod,piP, pi0, ell, nu_ell, nu_tau);
98 LorentzTransform boost1 = LorentzTransform::mkFrameTransformFromBeta(p.momentum().betaVec());
99 double cBeam = axis.dot(p.momentum().p3().unit());
100 if(nprod==2 && nu_tau.size()==1 && piP.size()==1) {
101 FourMomentum pPi = boost1.transform(piP[0].momentum());
102 double cTheta = pPi.p3().unit().dot(p.momentum().p3().unit());
103 _h_pi->fill(cBeam,cTheta);
104 }
105 else if(nprod==3 && nu_tau.size()==1 && ell.size()==1 && nu_ell.size()==1) {
106 if(ell[0].pid()==PID::ELECTRON) {
107 _h_e->fill(cBeam,2.*ell[0].momentum().t()/sqrtS());
108 }
109 else {
110 _h_mu->fill(cBeam,2.*ell[0].momentum().t()/sqrtS());
111 }
112 }
113 else if(nprod==3 && nu_tau.size()==1 && piP.size()==1&& pi0.size()==1) {
114 FourMomentum pRho = boost1.transform(piP[0].momentum()+pi0[0].momentum());
115 double cTheta = pRho.p3().unit().dot(p.momentum().p3().unit());
116 _h_rho->fill(cBeam,cTheta);
117 }
118 }
119 }
120
121 pair<double,double> calcP(Histo1DPtr hist,unsigned int imode) {
122 if(hist->numEntries()==0.) return make_pair(0.,0.);
123 double sum1(0.),sum2(0.);
124 for (const auto& bin : hist->bins() ) {
125 double Oi = bin.sumW();
126 if(Oi==0.) continue;
127 double ai(0.),bi(0.);
128 // tau -> pi/rho nu
129 if(imode==0) {
130 ai = 0.5*(bin.xMax()-bin.xMin());
131 bi = 0.5*ai*(bin.xMax()+bin.xMin());
132 }
133 // lepton mode
134 else {
135 ai = (-5*bin.xMin() + 3*pow(bin.xMin(),3) - pow(bin.xMin(),4) + 5*bin.xMax() - 3*pow(bin.xMax(),3) + pow(bin.xMax(),4))/3.;
136 bi = ( -bin.xMin() + 3*pow(bin.xMin(),3) - 2*pow(bin.xMin(),4) + bin.xMax() - 3*pow(bin.xMax(),3) + 2*pow(bin.xMax(),4))/3.;
137 }
138 double Ei = bin.errW();
139 sum1 += sqr(bi/Ei);
140 sum2 += bi/sqr(Ei)*(Oi-ai);
141 }
142 return make_pair(sum2/sum1,sqrt(1./sum1));
143 }
144
145 /// Normalise histograms etc., after the run
146 void finalize() {
147 Estimate1DPtr _h_P;
148 book(_h_P, 1, 1, 1);
149 unsigned int iy=1;
150 for (size_t ix=0; ix < _h_e->numBins(); ++ix) {
151 if (iy==2 || iy==10) ++iy;
152 normalize(_h_e->bin(iy));
153 pair<double,double> P_e = calcP(_h_e->bin(iy), 1);
154 double s1 = P_e.first/sqr(P_e.second);
155 double s2 = 1./sqr(P_e.second);
156 normalize(_h_mu->bin(iy));
157 pair<double,double> P_mu = calcP(_h_mu->bin(iy), 1);
158 s1 += P_mu.first/sqr(P_mu.second);
159 s2 += 1./sqr(P_mu.second);
160 normalize(_h_pi->bin(iy));
161 pair<double,double> P_pi = calcP(_h_pi->bin(iy), 0);
162 s1 += P_pi.first/sqr(P_pi.second);
163 s2 += 1./sqr(P_pi.second);
164 normalize(_h_rho->bin(iy));
165 pair<double,double> P_rho = calcP(_h_rho->bin(iy), 0);
166 P_rho.first /=0.46;
167 P_rho.second /=0.46;
168 s1 += P_rho.first/sqr(P_rho.second);
169 s2 += 1./sqr(P_rho.second);
170 // average
171 _h_P->bin(iy).set(s1/s2, sqrt(1./s2));
172 ++iy;
173 }
174 }
175
176 /// @}
177
178
179 /// @name Histograms
180 /// @{
181 Histo1DGroupPtr _h_e,_h_mu,_h_pi,_h_rho;
182 /// @}
183
184
185 };
186
187
188 RIVET_DECLARE_PLUGIN(L3_1998_I467929);
189
190}
|