rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

L3_1998_I467929

$\tau$ polarization at LEP1
Experiment: L3 (LEP)
Inspire ID: 467929
Status: VALIDATED
Authors:
  • Peter Richardson
References:
  • Phys.Lett.B 429 (1998) 387-398
Beams: e+ e-
Beam energies: (45.6, 45.6) GeV
Run details:
  • e+ e- > tau+ tau-

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