Rivet analyses referenceBESIII_2018_I1689296Radiative $J/\psi$ decays to $K^0_SK^0_S$Experiment: BESIII (BEPC) Inspire ID: 1689296 Status: VALIDATED NOHEPDATA Authors:
Beam energies: (1.6, 1.6) GeV Run details:
Measurement of mass and angular distributions in the decay $J/\psi\to\gamma K^0_SK^0S$. Plots were read from the paper and are not corrected for efficiency/acceptance, although the paper states the backgrounds are small. Source code: BESIII_2018_I1689296.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 -> gamma KS0 KS0
11 class BESIII_2018_I1689296 : public Analysis {
12 public:
13
14 /// Constructor
15 RIVET_DEFAULT_ANALYSIS_CTOR(BESIII_2018_I1689296);
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(PID::K0S);
28 declare(PSI, "PSI");
29 declare(Beam(), "Beams");
30 // hisotgrams
31 for(unsigned int ix=0;ix<3;++ix) {
32 if(ix<2) book(_h_mass[ix],1,1,1+ix);
33 book(_h_angle[ix],2,1,1+ix);
34 }
35 }
36
37 // angle cuts due regions of BES calorimeter
38 bool vetoPhoton(const double & cTheta) {
39 return cTheta>0.92 || (cTheta>0.8 && cTheta<0.86);
40 }
41
42 /// Perform the per-event analysis
43 void analyze(const Event& event) {
44 // get the axis, direction of incoming electron
45 const ParticlePair& beams = apply<Beam>(event, "Beams").beams();
46 Vector3 axis;
47 if(beams.first.pid()>0)
48 axis = beams.first .momentum().p3().unit();
49 else
50 axis = beams.second.momentum().p3().unit();
51 // find the J/psi decays
52 static const map<PdgId,unsigned int> & mode = { { 310,2},{ 22,1}};
53 DecayedParticles PSI = apply<DecayedParticles>(event, "PSI");
54 if( PSI.decaying().size()!=1) vetoEvent;
55 if(!PSI.modeMatches(0,3,mode)) vetoEvent;
56 // particles
57 const Particles & K0 = PSI.decayProducts()[0].at(310);
58 const Particle & gam = PSI.decayProducts()[0].at( 22)[0];
59 double mKK = (K0[0].momentum()+K0[1].momentum()).mass();
60 _h_mass[0]->fill(mKK);
61 for(unsigned int ix=0;ix<2;++ix)
62 _h_mass[1]->fill((gam.momentum()+K0[ix].momentum()).mass());
63 double cTheta = axis.dot(gam.p3().unit());
64 if(vetoPhoton(abs(cTheta))) vetoEvent;
65 _h_angle[0]->fill(cTheta);
66 // remaining angles
67 LorentzTransform boost1 = LorentzTransform::mkFrameTransformFromBeta(PSI.decaying()[0].momentum().betaVec());
68 FourMomentum pGamma = boost1.transform(gam.momentum());
69 FourMomentum pKK = boost1.transform(K0[0].momentum()+K0[1].momentum());
70 Vector3 e1z = pGamma.p3().unit();
71 Vector3 e1y = e1z.cross(axis).unit();
72 Vector3 e1x = e1y.cross(e1z).unit();
73 LorentzTransform boost2 = LorentzTransform::mkFrameTransformFromBeta(pKK.betaVec());
74 Vector3 axis2 = boost2.transform(boost1.transform(K0[0].momentum())).p3().unit();
75 _h_angle[1]->fill(e1z.dot(axis2));
76 double phi = atan2(e1y.dot(axis2),e1x.dot(axis2));
77 _h_angle[2]->fill(phi);
78 }
79
80
81 /// Normalise histograms etc., after the run
82 void finalize() {
83 for(unsigned int ix=0;ix<3;++ix) {
84 if(ix<2) normalize(_h_mass[ix],1.,false);
85 normalize(_h_angle[ix],1.,false);
86 }
87 }
88
89 /// @}
90
91
92 /// @name Histograms
93 /// @{
94 Histo1DPtr _h_mass[2], _h_angle[3];
95 /// @}
96
97
98 };
99
100
101 RIVET_DECLARE_PLUGIN(BESIII_2018_I1689296);
102
103}
|