Rivet analyses
Analysis of Kinematic distributions in η and η′ meson decays
Experiment: ()
Status: VALIDATED
Authors: - Peter Richardson
References: none listed
Beams: * *
Beam energies: ANY
Run details: - Any type of process producing eta/eta prime mesons
Analysis of Kinematic distributions in η and η′ meson decays. The mass distributions in the decays η, η′ → γπ+π−, γγπ0, π0π0π0, π+π−π0 and η′ → π+π−η, π0π0η are produced. Based on an +old Herwig++ internal analysis.
Source
code:MC_DECAY_ETA.cc
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/UnstableParticles.hh"
namespace Rivet {
class MC_DECAY_ETA : public Analysis {
public:
/// Constructor
RIVET_DEFAULT_ANALYSIS_CTOR(MC_DECAY_ETA);
/// @name Analysis methods
/// @{
/// Book histograms and initialise projections before the run
void init() {
// Initialise and register projections
declare(UnstableParticles(), "UFS");
// Book histograms
double meta[2] = {547.45, 957.78};
for (unsigned int ix = 0; ix < 2; ++ix) {
std::ostringstream title;
title << "_" << ix;
_mgammagamma.push_back(Histo1DPtr());
book(_mgammagamma.back(), "mgammagamma" + title.str(), 200, 0., meta[ix]);
_mpi0gamma.push_back(Histo1DPtr());
book(_mpi0gamma.back(), "mpi0gamma" + title.str(), 200, 0., meta[ix]);
_mpipgamma.push_back(Histo1DPtr());
book(_mpipgamma.back(), "mpipgamma" + title.str(), 200, 0., meta[ix]);
_mpimgamma.push_back(Histo1DPtr());
book(_mpimgamma.back(), "mpimgamma" + title.str(), 200, 0., meta[ix]);
_photonenergy.push_back(Histo1DPtr());
book(_photonenergy.back(), "photonenergy" + title.str(), 200, 0., meta[ix]);
_mpippim.push_back(Histo1DPtr());
book(_mpippim.back(), "mpippim" + title.str(), 200, 0., meta[ix]);
_dpippim.push_back(Histo1DPtr());
book(_dpippim.back(), "dpippim" + title.str(), 200, 200., meta[ix]);
_dpi0pi0.push_back(Histo1DPtr());
book(_dpi0pi0.back(), "dpi0pi0" + title.str(), 200, 200., meta[ix]);
_dpi0pip.push_back(Histo1DPtr());
book(_dpi0pip.back(), "dpi0pip" + title.str(), 200, 200., meta[ix]);
_dpi0pim.push_back(Histo1DPtr());
book(_dpi0pim.back(), "dpi0pim" + title.str(), 200, 200., meta[ix]);
}
_dpi0pi0.push_back(Histo1DPtr());
book(_dpi0pi0.back(), "dpi0pi0_2", 200, 200., 500.);
_dpippim.push_back(Histo1DPtr());
book(_dpippim.back(), "dpippim_2", 200, 200., 500.);
book(_dpipeta, "dpipeta", 200, 500., meta[1]);
book(_dpimeta, "dpimeta", 200, 500., meta[1]);
book(_dpi0eta, "dpi0eta", 200, 500., meta[1]);
}
void findDecayProducts(const Particle& mother,
unsigned int& nstable,
Particles& pip,
Particles& pim,
Particles& pi0,
Particles& eta,
Particles& gamma) {
for (const Particle& p : mother.children()) {
int id = p.pid();
if (id == PID::ETA) {
eta.push_back(p);
++nstable;
}
else if (id == PID::PHOTON) {
gamma.push_back(p);
++nstable;
}
else if (id == PID::PIPLUS) {
pip.push_back(p);
++nstable;
}
else if (id == PID::PIMINUS) {
pim.push_back(p);
++nstable;
}
else if (id == PID::PI0) {
pi0.push_back(p);
++nstable;
}
else if (id == PID::K0S || id == PID::K0L || id == PID::KPLUS || id == PID::KMINUS)
++nstable;
else if (!p.children().empty()) {
findDecayProducts(p, nstable, pip, pim, pi0, eta, gamma);
}
else
++nstable;
}
}
/// Perform the per-event analysis
void analyze(const Event& event) {
// Loop over f_1 mesons
for (const Particle& meson :
apply<UnstableParticles>(event, "UFS").particles(Cuts::pid == 221 || Cuts::pid == 331)) {
unsigned int nstable(0);
Particles pip, pim, pi0, eta, gamma;
findDecayProducts(meson, nstable, pip, pim, pi0, eta, gamma);
unsigned int imeson = meson.pid() == 221 ? 0 : 1;
// pi0 gamma gamma
if (nstable == 3 && pi0.size() == 1 && gamma.size() == 2) {
_mgammagamma[imeson]->fill((gamma[0].momentum() + gamma[1].momentum()).mass() / MeV);
_mpi0gamma[imeson]->fill((pi0[0].momentum() + gamma[0].momentum()).mass() / MeV);
_mpi0gamma[imeson]->fill((pi0[0].momentum() + gamma[1].momentum()).mass() / MeV);
} // pi+pi-gamma analysis
else if (nstable == 3 && pip.size() == 1 && pim.size() == 1 && gamma.size() == 1) {
FourMomentum ptemp = pip[0].momentum() + pim[0].momentum();
double mpipi = ptemp.mass();
_mpippim[imeson]->fill(mpipi / MeV);
double egamma = 0.5 * (meson.mass() * meson.mass() - mpipi * mpipi) / meson.mass();
_photonenergy[imeson]->fill(egamma / MeV);
_mpipgamma[imeson]->fill((pip[0].momentum() + gamma[0].momentum()).mass() / MeV);
_mpimgamma[imeson]->fill((pim[0].momentum() + gamma[0].momentum()).mass() / MeV);
}
else if (nstable == 3 && pi0.size() == 3) {
_dpi0pi0[imeson]->fill((pi0[0].momentum() + pi0[1].momentum()).mass() / MeV);
_dpi0pi0[imeson]->fill((pi0[0].momentum() + pi0[2].momentum()).mass() / MeV);
_dpi0pi0[imeson]->fill((pi0[1].momentum() + pi0[2].momentum()).mass() / MeV);
}
else if (nstable == 3 && pip.size() == 1 && pim.size() == 1 && pi0.size() == 1) {
_dpi0pip[imeson]->fill((pi0[0].momentum() + pip[0].momentum()).mass() / MeV);
_dpi0pim[imeson]->fill((pi0[0].momentum() + pim[0].momentum()).mass() / MeV);
_dpippim[imeson]->fill((pip[0].momentum() + pim[0].momentum()).mass() / MeV);
}
else if (nstable == 3 && pi0.size() == 2 && eta.size() == 1) {
_dpi0pi0[2]->fill((pi0[0].momentum() + pi0[1].momentum()).mass() / MeV);
_dpi0eta->fill((pi0[0].momentum() + eta[0].momentum()).mass() / MeV);
_dpi0eta->fill((pi0[1].momentum() + eta[0].momentum()).mass() / MeV);
}
else if (nstable == 3 && pip.size() == 1 && pim.size() == 1 && eta.size() == 1) {
_dpippim[2]->fill((pip[0].momentum() + pim[0].momentum()).mass() / MeV);
_dpipeta->fill((pip[0].momentum() + eta[0].momentum()).mass() / MeV);
_dpimeta->fill((pim[0].momentum() + eta[0].momentum()).mass() / MeV);
}
}
}
/// Normalise histograms etc., after the run
void finalize() {
// normalize to unity
for (unsigned int ix = 0; ix < 2; ++ix) {
normalize(_mgammagamma[ix]);
normalize(_mpi0gamma[ix]);
normalize(_mpipgamma[ix]);
normalize(_mpimgamma[ix]);
normalize(_mpippim[ix]);
normalize(_photonenergy[ix]);
normalize(_dpippim[ix]);
normalize(_dpi0pi0[ix]);
normalize(_dpi0pip[ix]);
normalize(_dpi0pim[ix]);
}
normalize(_dpi0pi0[2]);
normalize(_dpippim[2]);
normalize(_dpipeta);
normalize(_dpimeta);
normalize(_dpi0eta);
}
/// @}
/**
* Histograms for the decay \f$\eta\to\pi^0\gamma\gamma\f$
*/
/// @{
/**
* Histogram for the mass of \f$\gamma\gamma\f$
*/
vector<Histo1DPtr> _mgammagamma;
/**
* Histogrma for the mass of \f$\pi^0\gamma\f$
*/
vector<Histo1DPtr> _mpi0gamma;
/// @}
/**
* Histograms for the decay \f$\eta\to\pi^+\pi^-\gamma\f$
*/
/// @{
/**
* Histogram for the mass of \f$\pi^+\gamma\f$
*/
vector<Histo1DPtr> _mpipgamma;
/**
* Histogram for the mass of \f$\pi^-\gamma\f$
*/
vector<Histo1DPtr> _mpimgamma;
/**
* Histogram for the mass of \f$\pi^+\pi^-\f$
*/
vector<Histo1DPtr> _mpippim;
/**
* Histogram for the photon energy
*/
vector<Histo1DPtr> _photonenergy;
/// @}
/**
* Histograms for the decay \f$\eta\pi\pi\pi\f$ and \f$\eta'\to\eta\pi\pi\f$.
*/
/// @{
/**
* Histogram for the mass of \f$\pi^+\pi^-\f$
*/
vector<Histo1DPtr> _dpippim;
/**
* Histogram for the mass of \f$\pi^0\pi^0\f$
*/
vector<Histo1DPtr> _dpi0pi0;
/**
* Histogram for the mass of \f$\pi^0\pi^+\f$
*/
vector<Histo1DPtr> _dpi0pip;
/**
* Histogram for the mass of \f$\pi^0\pi^-\f$
*/
vector<Histo1DPtr> _dpi0pim;
/**
* Histogram for the mass of \f$\pi^+\eta\f$
*/
Histo1DPtr _dpipeta;
/**
* Histogram for the mass of \f$\pi^-\eta\f$
*/
Histo1DPtr _dpimeta;
/**
* Histogram for the mass of \f$\pi^0\eta\f$
*/
Histo1DPtr _dpi0eta;
/// @}
};
RIVET_DECLARE_PLUGIN(MC_DECAY_ETA);
}