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
src
Py8InterfaceBase.cc
Go to the documentation of this file.
1
#include "
GeneratorInterface/Pythia8Interface/interface/Py8InterfaceBase.h
"
2
#include "
GeneratorInterface/Pythia8Interface/interface/RandomP8.h
"
3
#include "
GeneratorInterface/Core/interface/RNDMEngineAccess.h
"
4
5
#include "
FWCore/Utilities/interface/Exception.h
"
6
7
using namespace
Pythia8;
8
9
namespace
gen
{
10
11
Py8InterfaceBase::Py8InterfaceBase(
edm::ParameterSet
const
& ps )
12
{
13
14
randomEngine
= &
getEngineReference
();
15
16
fMasterGen.reset(
new
Pythia);
17
fDecayer.reset(
new
Pythia);
18
19
fMasterGen->readString(
"Next:numberShowEvent = 0"
);
20
fDecayer->readString(
"Next:numberShowEvent = 0"
);
21
22
// RandomP8* RP8 = new RandomP8();
23
fMasterGen->setRndmEnginePtr(
new
RandomP8
() );
24
fDecayer->setRndmEnginePtr(
new
RandomP8
() );
25
26
fParameters = ps.
getParameter
<
edm::ParameterSet
>(
"PythiaParameters"
);
27
28
pythiaPylistVerbosity
= ps.
getUntrackedParameter
<
int
>(
"pythiaPylistVerbosity"
, 0);
29
pythiaHepMCVerbosity
= ps.
getUntrackedParameter
<
bool
>(
"pythiaHepMCVerbosity"
,
false
);
30
maxEventsToPrint
= ps.
getUntrackedParameter
<
int
>(
"maxEventsToPrint"
, 0);
31
32
}
33
34
bool
Py8InterfaceBase::readSettings(
int
)
35
{
36
37
for
(
ParameterCollector::const_iterator
line
= fParameters.begin();
38
line
!= fParameters.end(); ++
line
)
39
{
40
if
(
line
->find(
"Random:"
) != std::string::npos)
41
throw
cms::Exception
(
"PythiaError"
) <<
"Attempted to set random number "
42
"using Pythia commands. Please use "
"the RandomNumberGeneratorService."
43
<< std::endl;
44
45
if
(!fMasterGen->readString(*
line
))
throw
cms::Exception
(
"PythiaError"
)
46
<<
"Pythia 8 did not accept \""
47
<< *
line
<<
"\"."
<< std::endl;
48
}
49
50
return
true
;
51
52
}
53
54
bool
Py8InterfaceBase::declareStableParticles(
const
std::vector<int>& pdgIds )
55
{
56
57
for
(
size_t
i
=0;
i
<pdgIds.size();
i
++ )
58
{
59
// FIXME: need to double check if PID's are the same in Py6 & Py8,
60
// because the HepPDT translation tool is actually for **Py6**
61
//
62
// well, actually it looks like Py8 operates in PDT id's rather than Py6's
63
//
64
// int PyID = HepPID::translatePDTtoPythia( pdgIds[i] );
65
int
PyID = pdgIds[
i
];
66
std::ostringstream pyCard ;
67
pyCard << PyID <<
":mayDecay=false"
;
68
fMasterGen->readString( pyCard.str() );
69
// alternative:
70
// set the 2nd input argument warn=false
71
// - this way Py8 will NOT print warnings about unknown particle code(s)
72
// fMasterPtr->readString( pyCard.str(), false )
73
}
74
75
return
true
;
76
77
}
78
79
bool
Py8InterfaceBase:: declareSpecialSettings(
const
std::vector<std::string>& settings )
80
{
81
82
for
(
unsigned
int
iss=0; iss<settings.size(); iss++ )
83
{
84
if
( settings[iss].
find
(
"QED-brem-off"
) == std::string::npos )
continue
;
85
fMasterGen->readString(
"TimeShower:QEDshowerByL=off"
);
86
}
87
88
return
true
;
89
}
90
91
void
Py8InterfaceBase::statistics
()
92
{
93
94
fMasterGen->statistics();
95
return
;
96
97
}
98
99
}
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
i
int i
Definition:
DBlmapReader.cc:9
Py8InterfaceBase.h
spr::find
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition:
FindCaloHit.cc:7
relval_steps.gen
def gen
Definition:
relval_steps.py:239
gen::getEngineReference
CLHEP::HepRandomEngine & getEngineReference()
Definition:
RNDMEngineAccess.cc:9
geometryCSVtoXML.line
tuple line
Definition:
geometryCSVtoXML.py:15
RNDMEngineAccess.h
step1_ZMM_7Tev.pythiaPylistVerbosity
tuple pythiaPylistVerbosity
Definition:
step1_ZMM_7Tev.py:63
benchmark_cfg.statistics
tuple statistics
Definition:
benchmark_cfg.py:95
gen::ParameterCollector::const_iterator
Definition:
ParameterCollector.h:35
Exception.h
RandomP8.h
randomEngine
CLHEP::HepRandomEngine * randomEngine
Definition:
RandomP8.cc:5
cms::Exception
Definition:
Exception.h:68
RandomP8
Definition:
RandomP8.h:6
step1_ZMM_7Tev.maxEventsToPrint
tuple maxEventsToPrint
Definition:
step1_ZMM_7Tev.py:67
edm::ParameterSet
Definition:
ParameterSet.h:35
step1_ZMM_7Tev.pythiaHepMCVerbosity
tuple pythiaHepMCVerbosity
Definition:
step1_ZMM_7Tev.py:65
Generated for CMSSW Reference Manual by
1.8.5