10 #include <boost/algorithm/string.hpp>
11 #include <boost/lexical_cast.hpp>
24 std::string s_pythia_stream{
"PYTHIA8_INIT"};
29 #define PYTHIA8_NWEIGHTS nWeights
30 #define PYTHIA8_WEIGHT weight
31 #define PYTHIA8_WLABEL weightLabel
32 #define PYTHIA8_CONVERSION 1.0
34 #ifdef PYTHIA_VERSION_INTEGER
35 #undef PYTHIA8_NWEIGHTS
38 #if PYTHIA_VERSION_INTEGER > 8303
39 #define PYTHIA8_NWEIGHTS nWeightGroups
41 #define PYTHIA8_NWEIGHTS nVariationGroups
43 #define PYTHIA8_WEIGHT getGroupWeight
44 #define PYTHIA8_WLABEL getGroupName
53 std::string py8version()
55 static const std::string incdir (PY8INCLUDE_DIR);
56 std::string::size_type
pos = incdir.find (
"/pythia8/");
57 if (
pos == std::string::npos)
return "";
59 std::string::size_type pos2 = incdir.find (
"/",
pos);
60 if (pos2 == std::string::npos) pos2 = incdir.size();
61 return incdir.substr (
pos, pos2-
pos);
88 #ifndef PYTHIA8_3SERIES
105 m_runinfo = std::make_shared<HepMC3::GenRunInfo>();
107 struct HepMC3::GenRunInfo::ToolInfo
generator={std::string(
"Pythia8"),py8version(),std::string(
"Used generator")};
115 s_pythia_stream =
"PYTHIA8_INIT";
123 m_pythia->readString(
"Tune:pp = 5");
127 m_pythia->readString(
"PDF:pSet= LHAPDF6:cteq6l1");
130 m_pythia->readString(
"Next:numberShowEvent = 0");
133 m_pythia->settings.addFlag(
"AthenaPythia8ToHepMC:print_inconsistency",
true);
145 bool canSetHook =
true;
146 if (hook ==
"SuppressSmallPT") {
161 for(
const std::pair<const std::string, double> ¶m : Pythia8_UserHooks::UserHooksFactory::userSettings<double>()){
162 m_pythia->settings.addParm(param.first, param.second,
false,
false, 0., 0.);
165 for(
const std::pair<const std::string, int> ¶m : Pythia8_UserHooks::UserHooksFactory::userSettings<int>()){
166 m_pythia->settings.addMode(param.first, param.second,
false,
false, 0., 0.);
169 for(
const std::pair<const std::string, bool> ¶m : Pythia8_UserHooks::UserHooksFactory::userSettings<bool>()){
170 m_pythia->settings.addFlag(param.first, param.second);
173 for(
const std::pair<const std::string, std::string> ¶m : Pythia8_UserHooks::UserHooksFactory::userSettings<std::string>()){
174 m_pythia->settings.addWord(param.first, param.second);
179 ATH_MSG_ERROR(
"Unable to retrieve Athena Tool for custom Pythia processing");
180 return StatusCode::FAILURE;
191 if(
cmd.compare(
"")==0)
continue;
192 else if (
cmd.find(
"Beams:id") != std::string::npos ) {
193 ATH_MSG_ERROR(
"With command '" <<
cmd <<
"' you are trying to set a beam different from p: please use the Beam1/Beam2 properties instead:");
196 return StatusCode::FAILURE;
198 else if (
cmd.find(
"Beams:frameType") != std::string::npos ) {
200 ATH_MSG_WARNING(
" Found an explicit 'Beams:frameType' command: this will override transform beams momenta/energy parameters, regardless of its requested value. ");
207 return StatusCode::FAILURE;
219 return StatusCode::FAILURE;
226 ATH_MSG_INFO(
" !!!!!!!!!!!! WARNING ON PYTHIA RANDOM NUMBERS !!!!!!!!!! ");
227 ATH_MSG_INFO(
" THE ATHENA SERVICE AthRNGenSvc IS USED.");
228 ATH_MSG_INFO(
" !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ");
233 #if PYTHIA_VERSION_INTEGER >= 8310
238 s_pythia_stream =
"PYTHIA8";
240 ATH_MSG_INFO(
" !!!!!!!!!!!! WARNING ON PYTHIA RANDOM NUMBERS !!!!!!!!!! ");
241 ATH_MSG_INFO(
" THE STANDARD PYTHIA8 RANDOM NUMBER SERVICE IS USED.");
242 ATH_MSG_INFO(
" THE ATHENA SERVICE AthRNGSvc IS ***NOT*** USED.");
243 ATH_MSG_INFO(
" !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ");
254 std::vector<std::string> resonanceArgs;
257 if(resonanceArgs.size() != 2){
258 ATH_MSG_ERROR(
"Cannot Understand UserResonance job option!");
259 ATH_MSG_ERROR(
"You should specify it as a 'name:id1,id2,id3...'");
260 ATH_MSG_ERROR(
"Where name is the name of your UserResonance, and id1,id2,id3 are a comma separated list of PDG IDs to which it is applied");
263 std::vector<std::string> resonanceIds;
264 boost::split(resonanceIds, resonanceArgs.back(), boost::is_any_of(
","));
265 if(resonanceIds.size()==0){
266 ATH_MSG_ERROR(
"You did not specifiy any PDG ids to which your user resonance width should be applied!");
267 ATH_MSG_ERROR(
"You should specify a list as 'name:id1,id2,id3...'");
268 ATH_MSG_ERROR(
"Where name is the name of your UserResonance, and id1,id2,id3 are a comma separated list of PDG IDs to which it is applied");
272 for(std::vector<std::string>::const_iterator sId = resonanceIds.begin();
273 sId != resonanceIds.end(); ++sId){
274 int idResIn = boost::lexical_cast<int>(*sId);
279 #if PYTHIA_VERSION_INTEGER >= 8310
280 m_pythia->setResonancePtr(resonance);
282 m_pythia->setResonancePtr(resonance.get());
298 ATH_MSG_ERROR(
"Both LHE file and user process have been specified");
299 ATH_MSG_ERROR(
"LHE input does not make sense with a user process!");
303 canInit = canInit &&
m_pythia->readString(
"Beams:frameType = 4");
311 canInit = canInit &&
m_pythia->readString(
"Beams:frameType = 1");
319 #if PYTHIA_VERSION_INTEGER >= 8310
334 m_pythia->settings.writeFile(
"Settings_before.log",
true);
341 returnCode = StatusCode::FAILURE;
346 m_pythia->settings.writeFile(
"Settings_after.log",
true);
354 m_pythiaToHepMC.set_print_inconsistency(
m_pythia->settings.flag(
"AthenaPythia8ToHepMC:print_inconsistency") );
369 const EventContext& ctx = Gaudi::Hive::currentContext();
384 ATH_MSG_ERROR(
"Exceeded the max number of consecutive event failures.");
385 returnCode = StatusCode::FAILURE;
395 if(stat != StatusCode::SUCCESS) returnCode =
stat;
404 double eventWeight =
m_pythia->info.mergingWeight()*
m_pythia->info.weight();
407 if(returnCode != StatusCode::FAILURE &&
415 ATH_MSG_WARNING(
"Found event weight " << eventWeight <<
" between the float and double precision limits. Rejecting event.");
442 ATH_MSG_ERROR(
"Something wrong with this event - it contains fewer than 2 particles!");
444 return StatusCode::FAILURE;
452 auto evtlhe = std::make_shared<HepMC::GenEvent>();
454 auto extra = std::make_shared<HepMC::ShortEventAttribute>(evtlhe.get());
455 evt->add_attribute(
"LHERecord",
extra);
457 ATH_MSG_WARNING(
"LHE record saving requested, but not implemented with HEPMC2 functionality.");
482 ATH_MSG_DEBUG(
"No PDF information available in HepMC::GenEvent!");
490 return StatusCode::SUCCESS;
500 #ifndef PYTHIA8_304SERIES
503 if (hook->canEnhanceEmission()) {
507 #endif // not PYTHIA8_304SERIES
512 std::map<std::string,double> fWeights;
516 size_t atlas_specific_weights = 1;
517 fWeights[
"Default"]=eventWeight;
541 std::vector<std::string>::const_iterator
id =
m_weightIDs.begin()+atlas_specific_weights;
544 if(
m_pythia->info.getWeightsDetailedSize() != 0){
545 for(std::map<std::string, Pythia8::LHAwgt>::const_iterator wgt =
m_pythia->info.rwgt->wgts.begin();
546 wgt !=
m_pythia->info.rwgt->wgts.end(); ++wgt){
551 if(*
id != wgt->first){
552 ATH_MSG_ERROR(
"Mismatch in LHE3 weight id. Found "<<wgt->first<<
", expected "<<*
id);
553 return StatusCode::FAILURE;
558 std::map<std::string, Pythia8::LHAweight>::const_iterator weightName =
m_pythia->info.init_weights->find(wgt->first);
559 if(weightName !=
m_pythia->info.init_weights->end()){
572 for(
int iw = 1; iw <
m_pythia->info.PYTHIA8_NWEIGHTS(); ++iw){
588 fWeights[
"EXTRA_bare_LHE_weight"]=(-10.0)*
m_pythia->info.eventWeightLHEF;
597 std::map<std::string, Pythia8::LHAweight>::const_iterator
weight =
m_pythia->info.init_weights->find(
id);
606 ATH_MSG_ERROR(
"Something wrong when building list of weight names: " <<
m_weightNames.size() <<
" vs "<< fWeights.size() <<
", exiting ...");
607 return StatusCode::FAILURE;
613 if (!
evt->run_info()) {
614 evt->set_run_info(m_runinfo);
620 evt->weights().push_back(1.0);
626 evt->weights().resize(fWeights.size(), 1.0);
627 for (
auto w: fWeights) {
628 evt->weight(
w.first)=
w.second;
631 auto beams=
evt->beams();
639 evt->weights().clear();
641 auto beams=
evt->beam_particles();
642 ATH_MSG_DEBUG(
" Energy of the beams " << beams.first->momentum().e() );
649 return StatusCode::SUCCESS;
659 double xs =
info.sigmaGen();
663 ATH_MSG_DEBUG(
"Multiplying cross-section by CKKWL merging acceptance of "<<
m_nMerged <<
"/" <<
info.nAccepted() <<
" = " << accfactor
664 <<
": " << xs <<
" -> " << xs*accfactor);
671 std::cout <<
"Using FxFx cross section recipe: xs = "<<
m_sigmaTotal <<
" / " << 1e9*
info.nTried() << std::endl;
676 ATH_MSG_DEBUG(
"Multiplying cross-section by Pythia Modifier tool factor " << xsmod );
682 std::cout <<
"MetaData: cross-section (nb)= " << xs <<std::endl;
683 std::cout <<
"MetaData: generator= Pythia 8." << std::string(py8version()) <<std::endl;
686 std::cout<<
"MetaData: weights = ";
688 std::cout<<std::endl;
692 if (
info.nTried()>0)
ATH_MSG_INFO(
"Pythia8 efficiency (nAccepted/nTried %) = " <<
info.nAccepted()*100./
info.nTried());
693 else ATH_MSG_INFO(
"Pythia8 efficiency cannot be computed, nTried <=0");
702 return StatusCode::SUCCESS;
712 return s_pythia_stream;
719 std::string foundpath =
"";
730 #ifdef PYTHIA8_NWEIGHTS
731 #undef PYTHIA8_NWEIGHTS
732 #undef PYTHIA8_WEIGHT
733 #undef PYTHIA8_WLABEL
734 #undef PYTHIA8_CONVERSION