rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

BESII_2004_I652399

Mass and angular distributions in $J/\psi\to\omega \pi^+\pi^-$
Experiment: BESII (BEPC)
Inspire ID: 652399
Status: VALIDATED NOHEPDATA
Authors:
  • Peter Richardson
References:
  • Phys.Lett.B 598 (2004) 149-158
Beams: e- e+
Beam energies: (1.6, 1.6) GeV
Run details:
  • e+e- > J/psi

Mass and angular distributions in $J/\psi\to\omega \pi^+\pi^-$. The data were read from the figures in the paper and are not corrected for acceptance.

Source code: BESII_2004_I652399.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/Beam.hh"
  4#include "Rivet/Projections/UnstableParticles.hh"
  5#include "Rivet/Projections/DecayedParticles.hh"
  6
  7namespace Rivet {
  8
  9
 10  /// @brief J/psi -> omega pi+ pi-
 11  class BESII_2004_I652399 : public Analysis {
 12  public:
 13
 14    /// Constructor
 15    RIVET_DEFAULT_ANALYSIS_CTOR(BESII_2004_I652399);
 16
 17
 18    /// @name Analysis methods
 19    /// @{
 20
 21    /// Book histograms and initialise projections before the run
 22    void init() {
 23      // Initialise and register projections
 24      UnstableParticles ufs = UnstableParticles(Cuts::abspid==443);
 25      declare(ufs, "UFS");
 26      DecayedParticles PSI(ufs);
 27      PSI.addStable(223);
 28      declare(PSI, "PSI");
 29      declare(Beam(), "Beams");
 30      // histogram
 31      for (unsigned int ix=0; ix<2; ++ix) {
 32        book(_h_mass[ix],1,1,1+ix);
 33        for (unsigned int iy=0; iy<4; ++iy) {
 34          book(_h_angle[iy][ix],2,1+iy,1+ix);
 35        }
 36      }
 37    }
 38
 39    void findChildren(const Particle& p, Particles& pim, Particles& pip,
 40                      Particles& pi0, unsigned int &ncount) {
 41      for (const Particle& child : p.children()) {
 42        if (child.pid()==PID::PIPLUS) {
 43          pip.push_back(child);
 44          ncount+=1;
 45        }
 46        else if (child.pid()==PID::PIMINUS) {
 47          pim.push_back(child);
 48          ncount+=1;
 49        }
 50        else if (child.pid()==PID::PI0) {
 51          pi0.push_back(child);
 52          ncount+=1;
 53        }
 54        else if (child.children().empty()) {
 55          ncount+=1;
 56        }
 57        else {
 58          findChildren(child,pim,pip,pi0,ncount);
 59        }
 60      }
 61    }
 62
 63    /// Perform the per-event analysis
 64    void analyze(const Event& event) {
 65      // get the axis, direction of incoming electron
 66      const ParticlePair& beams = apply<Beam>(event, "Beams").beams();
 67      Vector3 axis;
 68      if (beams.first.pid()>0) {
 69        axis = beams.first .mom().p3().unit();
 70      }
 71      else {
 72        axis = beams.second.mom().p3().unit();
 73      }
 74      // J/psi decay products
 75      DecayedParticles PSI = apply<DecayedParticles>(event, "PSI");
 76      if ( PSI.decaying().size()!=1) vetoEvent;
 77      if (!PSI.modeMatches(0,3,mode)) vetoEvent;
 78      const Particle& Pip   = PSI.decayProducts()[0].at( 211)[0];
 79      const Particle& Pim   = PSI.decayProducts()[0].at(-211)[0];
 80      const Particle& omega = PSI.decayProducts()[0].at( 223)[0];
 81      // mass histograms
 82      FourMomentum pPiPi = Pip.mom()+Pim   .mom();
 83      const double mPiPi = pPiPi.mass();
 84      _h_mass[0]->fill(mPiPi/GeV);
 85      _h_mass[1]->fill((Pip.mom()+omega.mom()).mass()/GeV);
 86      _h_mass[1]->fill((Pim.mom()+omega.mom()).mass()/GeV);
 87      unsigned int iplot=0;
 88      if      (mPiPi>.3 && mPiPi<.5) iplot=0;
 89      else if (mPiPi>.5 && mPiPi<.7) iplot=1;
 90      else  return;
 91      if (abs(Pip.p3().unit().dot(axis))>0.84) return;
 92      if (abs(Pim.p3().unit().dot(axis))>0.84) return;
 93      LorentzTransform boost;
 94      if (PSI.decaying()[0].mom().p3().mod() > 1*MeV) {
 95        boost = LorentzTransform::mkFrameTransformFromBeta(PSI.decaying()[0].mom().betaVec());
 96      }
 97      pPiPi = boost.transform(pPiPi);
 98      // omega decay
 99      unsigned int ncount=0;
100      Particles pip, pim, pi0;
101      findChildren(omega,pim,pip,pi0,ncount);
102      if ( ncount!=3 || !(pim.size()==1 && pip.size()==1 && pi0.size()==1))  return;
103      if (abs(pip[0].p3().unit().dot(axis))>0.84)  return;
104      if (abs(pim[0].p3().unit().dot(axis))>0.84)  return;
105      FourMomentum pOmega = boost.transform(omega.mom());
106      Vector3 e1Z = pOmega.p3().unit();
107      Vector3 e1Y = e1Z.cross(axis).unit();
108      Vector3 e1X = e1Y.cross(e1Z).unit();
109      const double cOmega = e1Z.dot(axis);
110      _h_angle[1][iplot]->fill(cOmega);
111      // pion angles
112      LorentzTransform boost2 = LorentzTransform::mkFrameTransformFromBeta(pPiPi.betaVec());
113      FourMomentum pPi = boost2.transform(boost.transform(Pip.mom()));
114      const double cPi = pPi.p3().unit().dot(e1Z);
115      _h_angle[2][iplot]->fill(cPi);
116      Vector3 trans1 = pPi.p3().unit()-cPi*e1Z;
117      // omega angles
118      LorentzTransform boost3 = LorentzTransform::mkFrameTransformFromBeta(pOmega.betaVec());
119      FourMomentum ppip = boost3.transform(boost.transform(pip[0].mom()));
120      FourMomentum ppim = boost3.transform(boost.transform(pim[0].mom()));
121      Vector3 nW = ppip.p3().cross(ppim.p3()).unit();
122      const double bW = nW.dot(e1X);
123      _h_angle[3][iplot]->fill(bW);
124      Vector3 trans2 = ppip.p3()-ppip.p3().dot(e1Z)*e1Z;
125      const double chi = abs(atan2(trans1.cross(trans2).dot(e1Z),trans1.dot(trans2)));
126      _h_angle[0][iplot]->fill(chi/M_PI*180.);
127    }
128
129
130    /// Normalise histograms etc., after the run
131    void finalize() {
132      normalize(_h_mass, 1.0, false);
133      for (unsigned int iy=0; iy<4; ++iy) {
134        normalize(_h_angle[iy], 1.0, false);
135      }
136    }
137
138    /// @}
139
140
141    /// @name Histograms
142    /// @{
143    Histo1DPtr _h_mass[2],_h_angle[4][2];
144    const map<PdgId,unsigned int> mode = { { 223,1}, { 211,1},{-211,1}};
145    /// @}
146
147
148  };
149
150
151  RIVET_DECLARE_PLUGIN(BESII_2004_I652399);
152
153}