1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/UnstableParticles.hh"
 4#include "Rivet/Projections/DecayedParticles.hh"
 5
 6namespace Rivet {
 7
 8
 9  /// @brief B0 -> rho+ rho-
10  class BABAR_2004_I632736 : public Analysis {
11  public:
12
13    /// Constructor
14    RIVET_DEFAULT_ANALYSIS_CTOR(BABAR_2004_I632736);
15
16
17    /// @name Analysis methods
18    /// @{
19
20    /// Book histograms and initialise projections before the run
21    void init() {
22      // projections
23      UnstableParticles ufs = UnstableParticles(Cuts::abspid==511);
24      declare(ufs, "UFS");
25      DecayedParticles B0(ufs);
26      B0.addStable( 213);
27      B0.addStable(-213);
28      declare(B0, "B0");
29      // histos
30      book(_p[0],1,1,1);
31      book(_p[1],"TMP/wgt");
32    }
33
34
35    /// Perform the per-event analysis
36    void analyze(const Event& event) {
37      static const map<PdgId,unsigned int> & mode   = { { 213,1}, {-213,1} };
38      DecayedParticles B0 = apply<DecayedParticles>(event, "B0");
39      for(unsigned int ix=0;ix<B0.decaying().size();++ix) {
40       	int sign = B0.decaying()[ix].pid()/B0.decaying()[ix].abspid();
41      	if(!B0.modeMatches(ix,2,mode)) continue;
42	const Particle & rhop = B0.decayProducts()[ix].at( sign*213)[0];
43	//const Particle & rhom = B0.decayProducts()[ix].at(-sign*213)[0];
44	// compute the helicity angle
45	if(rhop.children().size()!=2) continue;
46	Particle pip;
47	if(rhop.children()[0].pid()==sign*211 &&
48	   rhop.children()[1].pid()==111)
49	  pip = rhop.children()[0];
50	else if(rhop.children()[1].pid()==sign*211 &&
51		rhop.children()[0].pid()==111)
52	  pip = rhop.children()[1];
53	  else
54	    continue;
55	LorentzTransform boost1 =
56	  LorentzTransform::mkFrameTransformFromBeta(B0.decaying()[ix].momentum().betaVec());
57	FourMomentum prho = boost1.transform(rhop.momentum());
58	FourMomentum ppi  = boost1.transform(pip);
59	LorentzTransform boost2 =
60	  LorentzTransform::mkFrameTransformFromBeta(prho.betaVec());
61	ppi = boost2.transform(ppi);
62	double cTheta = ppi.p3().unit().dot(prho.p3().unit());
63	_p[0]->fill(-0.5*(1.-5.*sqr(cTheta)));
64        _p[1]->fill();
65      }
66    }
67
68
69    /// Normalise histograms etc., after the run
70    void finalize() {
71      scale(_p[0],1./ *_p[1]);
72    }
73
74    /// @}
75
76
77    /// @name Histograms
78    /// @{
79    CounterPtr _p[2];
80    /// @}
81
82
83  };
84
85
86  RIVET_DECLARE_PLUGIN(BABAR_2004_I632736);
87
88}