1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/Beam.hh"
 4#include "Rivet/Projections/UnstableParticles.hh"
 5
 6namespace Rivet {
 7
 8
 9  /// @brief Spectrum for Lambda_c(2625)
10  class ARGUS_1993_I357132 : public Analysis {
11  public:
12
13    /// Constructor
14    RIVET_DEFAULT_ANALYSIS_CTOR(ARGUS_1993_I357132);
15
16
17    /// @name Analysis methods
18    ///@{
19
20    /// Book histograms and initialise projections before the run
21    void init() {
22      // projections
23      declare(Beam(), "Beams");
24      declare(UnstableParticles(), "UFS");
25      // book histos
26      book(_h_rate,1,1,1);
27      book(_h_x,2,1,1);
28    }
29    
30    void findDecayProducts(Particle parent, Particles & Lambda_c, Particles & pions,unsigned int & nstable) {
31      for(const Particle & p : parent.children()) {
32	if(p.abspid()==4122) {
33	  Lambda_c.push_back(p);
34	  ++nstable;
35	}
36	else if(p.abspid()==PID::PIPLUS) {
37	  pions.push_back(p);
38	  ++nstable;
39	}
40	else if(!p.children().empty())
41	  findDecayProducts(p,Lambda_c,pions,nstable);
42	else
43	  ++nstable;
44      }
45    }
46
47
48    /// Perform the per-event analysis
49    void analyze(const Event& event) {
50      static const int id2625 = 4124;
51      // Get beams and average beam momentum
52      const ParticlePair& beams = apply<Beam>(event, "Beams").beams();
53      const double Emax = ( beams.first.p3().mod() + beams.second.p3().mod() ) / 2.0;
54      const double Pmax = sqrt(sqr(Emax)-sqr(2.625));
55      const UnstableParticles& ufs = apply<UnstableParticles>(event, "UFS");
56      for (const Particle& p : ufs.particles(Cuts::abspid==id2625)) {
57	double xp = p.momentum().p3().mod()/Pmax;
58	_h_x->fill(xp);
59	Particles Lambda_c,pions;
60	unsigned int nstable(0);
61	findDecayProducts(p,Lambda_c,pions,nstable);
62	if(nstable==3&&pions.size()==2&&Lambda_c.size()==1)
63	  _h_rate->fill(10.58);
64      }
65    }
66
67
68    /// Normalise histograms etc., after the run
69    void finalize() {
70      normalize(_h_x);
71      scale(_h_rate,crossSection()/picobarn/sumOfWeights());
72    }
73
74    ///@}
75
76
77    /// @name Histograms
78    ///@{
79    Histo1DPtr _h_x,_h_rate;
80    ///@}
81
82
83  };
84
85
86  RIVET_DECLARE_PLUGIN(ARGUS_1993_I357132);
87
88}