1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/UnstableParticles.hh"
 4
 5namespace Rivet {
 6
 7
 8  /// @brief Xi_c -> xi pi pi
 9  class BELLE_2018_I1698390 : public Analysis {
10  public:
11
12    /// Constructor
13    RIVET_DEFAULT_ANALYSIS_CTOR(BELLE_2018_I1698390);
14
15
16    /// @name Analysis methods
17    /// @{
18
19    /// Book histograms and initialise projections before the run
20    void init() {
21      // projection
22      declare(UnstableParticles(),"UFS");
23      // histo
24      book(_h,1,1,1);
25    }
26
27    void findDecay(int & sign, const Particle & parent, Particles & xi, Particles & pi,
28		   unsigned int & nstable) {
29      for(const Particle & p : parent.children()) {
30	if(p.pid()==sign*3312) {
31	  ++nstable;
32	  xi.push_back(p);
33	}
34	else if(p.pid()==sign*211) {
35	  ++nstable;
36	  pi.push_back(p);
37	}
38	else if(p.children().empty())
39	  ++nstable;
40	else
41	  findDecay(sign,p,xi,pi,nstable);
42      }
43    }
44
45
46    /// Perform the per-event analysis
47    void analyze(const Event& event) {
48      for(const Particle & p : apply<UnstableParticles>(event,"UFS").particles(Cuts::abspid==4232)) {
49	int sign = p.pid()/p.abspid();
50	Particles xi,pi;
51	unsigned int nstable(0);
52	findDecay(sign,p,xi,pi,nstable);
53	LorentzTransform boost = LorentzTransform::mkFrameTransformFromBeta(p.momentum().betaVec());
54	if(nstable==3&&xi.size()==1&&pi.size()==2) {
55	  double p1 = boost.transform(pi[0].momentum()).p3().mod();
56	  double p2 = boost.transform(pi[1].momentum()).p3().mod();
57	  double mass = p1>p2 ?
58	    (pi[1].momentum()+xi[0].momentum()).mass() :
59	    (pi[0].momentum()+xi[0].momentum()).mass();
60	  _h->fill(mass);
61	}
62      }
63    }
64
65
66    /// Normalise histograms etc., after the run
67    void finalize() {
68      normalize(_h);
69    }
70
71    /// @}
72
73
74    /// @name Histograms
75    /// @{
76    Histo1DPtr _h;
77    /// @}
78
79
80  };
81
82
83  RIVET_DECLARE_PLUGIN(BELLE_2018_I1698390);
84
85}