Main Page
Namespaces
Classes
Package Documentation
CVS Directory
WorkBook
Offline Guide
Release schedule
•
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Properties
Friends
Macros
Pages
src
GeneratorInterface
Pythia8Interface
plugins
Py8JetGun.cc
Go to the documentation of this file.
1
2
#include "
GeneratorInterface/Core/interface/GeneratorFilter.h
"
3
#include "
GeneratorInterface/ExternalDecays/interface/ExternalDecayDriver.h
"
4
5
#include "
GeneratorInterface/Pythia8Interface/interface/Py8GunBase.h
"
6
#include "
GeneratorInterface/Pythia8Interface/interface/RandomP8.h
"
7
8
namespace
gen
{
9
10
class
Py8JetGun
:
public
Py8GunBase
{
11
12
public
:
13
14
Py8JetGun
(
edm::ParameterSet
const
& );
15
~Py8JetGun
() {}
16
17
bool
generatePartonsAndHadronize
();
18
const
char
*
classname
()
const
;
19
20
private
:
21
22
// PtGun particle(s) characteristics
23
double
fMinEta
;
24
double
fMaxEta
;
25
double
fMinP
;
26
double
fMaxP
;
27
double
fMinE
;
28
double
fMaxE
;
29
30
};
31
32
// implementation
33
//
34
Py8JetGun::Py8JetGun
(
edm::ParameterSet
const
& ps )
35
:
Py8GunBase
(ps)
36
{
37
38
// ParameterSet defpset ;
39
edm::ParameterSet
pgun_params =
40
ps.
getParameter
<
edm::ParameterSet
>(
"PGunParameters"
);
// , defpset ) ;
41
fMinEta
= pgun_params.
getParameter
<
double
>(
"MinEta"
);
// ,-2.2);
42
fMaxEta
= pgun_params.
getParameter
<
double
>(
"MaxEta"
);
// , 2.2);
43
fMinP
= pgun_params.
getParameter
<
double
>(
"MinP"
);
// , 0.);
44
fMaxP
= pgun_params.
getParameter
<
double
>(
"MaxP"
);
// , 0.);
45
fMinE
= pgun_params.
getParameter
<
double
>(
"MinE"
);
// , 0.);
46
fMaxE
= pgun_params.
getParameter
<
double
>(
"MaxE"
);
// , 0.);
47
48
}
49
50
bool
Py8JetGun::generatePartonsAndHadronize
()
51
{
52
53
fMasterGen
->event.reset();
54
55
double
totPx = 0.;
56
double
totPy = 0.;
57
double
totPz = 0.;
58
double
totE = 0.;
59
double
totM = 0.;
60
double
phi
,
eta
, the, ee,
pp
;
61
62
for
(
size_t
i
=0;
i
<
fPartIDs
.size();
i
++ )
63
{
64
65
int
particleID =
fPartIDs
[
i
];
// this is PDG - need to convert to Py8 ???
66
67
// FIXME !!!
68
// Ouch, it's using bare randomEngine pointer - that's NOT safe.
69
// Need to hold a pointer somewhere properly !!!
70
//
71
phi = 2. *
M_PI
*
randomEngine
->flat() ;
72
the = acos( -1. + 2.*
randomEngine
->flat() );
73
74
// from input
75
//
76
ee = (
fMaxE
-
fMinE
)*
randomEngine
->flat() +
fMinE
;
77
78
double
mass = (
fMasterGen
->particleData).mass( particleID );
79
// double mass = (pythia->particleData).m0( particleID );
80
81
pp =
sqrt
( ee*ee - mass*mass );
82
83
double
px = pp *
sin
(the) *
cos
(phi);
84
double
py = pp *
sin
(the) *
sin
(phi);
85
double
pz = pp *
cos
(the);
86
87
if
( !((
fMasterGen
->particleData).isParticle( particleID )) )
88
{
89
particleID = std::fabs(particleID) ;
90
}
91
(
fMasterGen
->event).
append
( particleID, 1, 0, 0, px, py, pz, ee, mass );
92
93
// values for computing total mass
94
//
95
totPx += px;
96
totPy += py;
97
totPz += pz;
98
totE += ee;
99
100
}
101
102
totM =
sqrt
( totE*totE - (totPx*totPx+totPy*totPy+totPz*totPz) );
103
104
//now the boost (from input params)
105
//
106
pp = (
fMaxP
-
fMinP
)*
randomEngine
->flat() +
fMinP
;
107
ee =
sqrt
( totM*totM + pp*pp );
108
109
//the boost direction (from input params)
110
//
111
phi = (
fMaxPhi
-
fMinPhi
)*
randomEngine
->flat() +
fMinPhi
;
112
eta = (
fMaxEta
-
fMinEta
)*
randomEngine
->flat() +
fMinEta
;
113
the = 2.*atan(
exp
(-eta));
114
115
double
betaX = pp/ee *
std::sin
(the) *
std::cos
(phi);
116
double
betaY = pp/ee *
std::sin
(the) *
std::sin
(phi);
117
double
betaZ = pp/ee *
std::cos
(the);
118
119
// boost all particles
120
//
121
(
fMasterGen
->event).bst( betaX, betaY, betaZ );
122
123
if
( !
fMasterGen
->next() )
return
false
;
124
125
event
().reset(
new
HepMC::GenEvent);
126
toHepMC
.fill_next_event(
fMasterGen
->event,
event
().
get
() );
127
128
return
true
;
129
130
}
131
132
const
char
*
Py8JetGun::classname
()
const
133
{
134
return
"Py8JetGun"
;
135
}
136
137
typedef
edm::GeneratorFilter<gen::Py8JetGun, gen::ExternalDecayDriver>
Pythia8JetGun
;
138
139
}
// end namespace
140
141
using
gen::Pythia8JetGun
;
142
DEFINE_FWK_MODULE
(
Pythia8JetGun
);
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
i
int i
Definition:
DBlmapReader.cc:9
create_public_lumi_plots.exp
tuple exp
Definition:
create_public_lumi_plots.py:1104
gen::Py8InterfaceBase::fMasterGen
std::auto_ptr< Pythia8::Pythia > fMasterGen
Definition:
Py8InterfaceBase.h:33
gen::Py8JetGun::fMinEta
double fMinEta
Definition:
Py8JetGun.cc:23
gen::Py8JetGun::~Py8JetGun
~Py8JetGun()
Definition:
Py8JetGun.cc:15
gen::Py8JetGun::fMinE
double fMinE
Definition:
Py8JetGun.cc:27
createTree.pp
tuple pp
Definition:
createTree.py:15
gen::Py8JetGun::classname
const char * classname() const
Definition:
Py8JetGun.cc:132
DEFINE_FWK_MODULE
#define DEFINE_FWK_MODULE(type)
Definition:
MakerMacros.h:17
funct::sin
Sin< T >::type sin(const T &t)
Definition:
Sin.h:22
gen::Py8GunBase::fMinPhi
double fMinPhi
Definition:
Py8GunBase.h:56
gen::Py8JetGun::fMaxP
double fMaxP
Definition:
Py8JetGun.cc:26
gen::Py8InterfaceBase::toHepMC
HepMC::I_Pythia8 toHepMC
Definition:
Py8InterfaceBase.h:35
python.multivaluedict.append
def append
Definition:
multivaluedict.py:73
ExternalDecayDriver.h
eta
T eta() const
Definition:
Basic3DVectorLD.h:179
gen::Py8JetGun
Definition:
Py8JetGun.cc:10
gen::Py8JetGun::fMaxE
double fMaxE
Definition:
Py8JetGun.cc:28
relval_steps.gen
def gen
Definition:
relval_steps.py:239
gen::Py8JetGun::Py8JetGun
Py8JetGun(edm::ParameterSet const &)
Definition:
Py8JetGun.cc:34
Py8GunBase.h
mathSSE::sqrt
T sqrt(T t)
Definition:
SSEVec.h:48
gen::Py8GunBase
Definition:
Py8GunBase.h:25
funct::cos
Cos< T >::type cos(const T &t)
Definition:
Cos.h:22
GeneratorFilter.h
gen::Py8GunBase::fPartIDs
std::vector< int > fPartIDs
Definition:
Py8GunBase.h:55
gen::Py8GunBase::event
std::auto_ptr< HepMC::GenEvent > & event()
Definition:
Py8GunBase.h:50
edm::GeneratorFilter
Definition:
GeneratorFilter.h:34
M_PI
#define M_PI
Definition:
BFit3D.cc:3
RandomP8.h
randomEngine
CLHEP::HepRandomEngine * randomEngine
Definition:
RandomP8.cc:5
gen::Pythia8JetGun
edm::GeneratorFilter< gen::Py8JetGun, gen::ExternalDecayDriver > Pythia8JetGun
Definition:
Py8JetGun.cc:137
gen::Py8GunBase::fMaxPhi
double fMaxPhi
Definition:
Py8GunBase.h:57
edm::ParameterSet
Definition:
ParameterSet.h:35
gen::Py8JetGun::fMinP
double fMinP
Definition:
Py8JetGun.cc:25
gen::Py8JetGun::fMaxEta
double fMaxEta
Definition:
Py8JetGun.cc:24
gen::Py8JetGun::generatePartonsAndHadronize
bool generatePartonsAndHadronize()
Definition:
Py8JetGun.cc:50
phi
Definition:
DDAxes.h:10
Generated for CMSSW Reference Manual by
1.8.5