1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
| // -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/Beam.hh"
#include "Rivet/Projections/UnstableParticles.hh"
namespace Rivet {
/// @brief Add a short analysis description here
class CLEO_2000_I526554 : public Analysis {
public:
/// Constructor
RIVET_DEFAULT_ANALYSIS_CTOR(CLEO_2000_I526554);
/// @name Analysis methods
//@{
/// Book histograms and initialise projections before the run
void init() {
declare(Beam(), "Beams");
declare(UnstableParticles(), "UFS");
// Book histograms
book(_h_Ds_star1 , 1, 1, 1);
book(_h_Ds , 2, 1, 1);
book(_h_Ds_star2 , 3, 1, 1);
book(_h_Ds_primary, 4, 1, 1);
}
/// Perform the per-event analysis
void analyze(const Event& event) {
// Loop through unstable FS particles and look for charmed mesons
const UnstableParticles& ufs = apply<UnstableParticles>(event, "UFS");
const Beam beamproj = apply<Beam>(event, "Beams");
const ParticlePair& beams = beamproj.beams();
const FourMomentum mom_tot = beams.first.momentum() + beams.second.momentum();
LorentzTransform cms_boost;
if (mom_tot.p3().mod() > 1*MeV)
cms_boost = LorentzTransform::mkFrameTransformFromBeta(mom_tot.betaVec());
const double s = sqr(beamproj.sqrtS());
for (const Particle& p : ufs.particles(Cuts::abspid==431 or Cuts::abspid==433)) {
// 3-momentum in CMS frame
const double mom = cms_boost.transform(p.momentum()).vector3().mod();
const int pdgid = p.abspid();
double mH2(0.),xp(0.);
bool primary = true;
switch (pdgid) {
case 431:
// MSG_DEBUG("D_s found");
mH2 = sqr(1.96834);
xp = mom/sqrt(s/4.0 - mH2);
_h_Ds->fill(xp);
for(const Particle & mother : p.parents()) {
if(PID::isCharmMeson(mother.pid())) {
primary = false;
break;
}
}
if(primary)
_h_Ds_primary->fill(xp);
break;
case 433:
MSG_DEBUG("D_s* found");
mH2 = sqr(2.1122);
xp = mom/sqrt(s/4.0 - mH2);
_h_Ds_star1->fill(xp);
_h_Ds_star2->fill(xp);
break;
default:
break;
}
}
}
/// Normalise histograms etc., after the run
void finalize() {
// BR(D_s->phi pi) x BR(phi->K+K-)
double br2 = 0.0227;
// x D_s* -> gamma d_s
double br1 = br2*.935;
// cross section factor
double fact = crossSection()/picobarn/sumOfWeights();
// normalize the cross sections (bin width)
scale(_h_Ds_star1 , br1*fact*0.03 );
scale(_h_Ds , br2*fact*0.03 );
scale(_h_Ds_star2 , br2*fact*0.03 );
scale(_h_Ds_primary, br2*fact*0.03 );
}
//@}
/// @name Histograms
//@{
Histo1DPtr _h_Ds,_h_Ds_primary,_h_Ds_star1,_h_Ds_star2;
//@}
};
// The hook for the plugin system
RIVET_DECLARE_PLUGIN(CLEO_2000_I526554);
}
|