#include "SHERPA/Initialization/Initialization_Handler.H" #include "SHERPA/PerturbativePhysics/Hard_Decay_Handler.H" #include "SHERPA/PerturbativePhysics/Shower_Handler.H" #include "SHERPA/SoftPhysics/Beam_Remnant_Handler.H" #include "SHERPA/SoftPhysics/Colour_Reconnection_Handler.H" #include "SHERPA/SoftPhysics/Hadron_Decay_Handler.H" #include "SHERPA/SoftPhysics/Hadron_Init.H" #include "SHERPA/SoftPhysics/Soft_Collision_Handler.H" #include "SHERPA/PerturbativePhysics/MI_Handler.H" #include "SHERPA/SoftPhysics/Soft_Photon_Handler.H" #include "SHERPA/Tools/Event_Reader_Base.H" #include "SHERPA/Main/Filter.H" #include "PHASIC++/Scales/Core_Scale_Setter.H" #include "MODEL/Main/Model_Base.H" #include "MODEL/Main/Running_AlphaS.H" #include "METOOLS/Currents/C_Spinor.H" #include "PDF/Main/Structure_Function.H" #include "PDF/Main/Intact.H" #include "PDF/Main/PDF_Base.H" #include "REMNANTS/Tools/Remnants_Parameters.H" #include "ATOOLS/Org/Message.H" #include "ATOOLS/Org/CXXFLAGS_PACKAGES.H" #include "ATOOLS/Math/Scaling.H" #include "ATOOLS/Phys/Spinor.H" #include "ATOOLS/Phys/Variations.H" #include "ATOOLS/Phys/Fragmentation_Base.H" #include "ATOOLS/Org/Shell_Tools.H" #include "ATOOLS/Math/Variable.H" #include "ATOOLS/Org/Data_Writer.H" #include "SHERPA/Single_Events/Hadron_Decays.H" #include "ATOOLS/Org/Library_Loader.H" #include "PHASIC++/Main/Phase_Space_Handler.H" #include "PHASIC++/Selectors/Selector.H" #include "PHASIC++/Process/ME_Generator_Base.H" #include "PHASIC++/Channels/Channel_Generator.H" #include "PDF/Main/NLOMC_Base.H" #include "PDF/Main/Shower_Base.H" #include "ATOOLS/Org/MyStrStream.H" #include "ATOOLS/Math/Random.H" #include "ATOOLS/Org/My_MPI.H" #include "ATOOLS/Org/Scoped_Settings.H" #include "ATOOLS/Org/Run_Parameter.H" using namespace SHERPA; using namespace MODEL; using namespace BEAM; using namespace PDF; using namespace REMNANTS; using namespace ATOOLS; using namespace std; typedef void (*PDF_Init_Function)(); typedef void (*PDF_Exit_Function)(); Initialization_Handler::Initialization_Handler() : m_mode(eventtype::StandardPerturbative), m_savestatus(false), p_model(NULL), p_beamspectra(NULL), p_mehandler(NULL), p_harddecays(NULL), p_beamremnants(NULL), p_reconnections(NULL), p_fragmentation(NULL), p_hdhandler(NULL), p_softphotons(NULL), p_evtreader(NULL), p_variations(NULL), p_filter(NULL) { RegisterDefaults(); Settings& s = Settings::GetMainSettings(); // configure runtime parameters if (s["SAVE_STATUS"].Get() != "") { std::string savestatus(s["SAVE_STATUS"].Get()); if (savestatus[savestatus.size() - 1] != '/') savestatus += "/"; rpa->gen.SetVariable("SHERPA_STATUS_PATH", rpa->gen.Variable("SHERPA_RUN_PATH") + "/" + savestatus); m_savestatus=true; } else { rpa->gen.SetVariable("SHERPA_STATUS_PATH", ""); } m_evtform = s["EVENT_INPUT"].Get(); if (m_evtform != "") { m_mode=eventtype::EventReader; msg_Out()<<"Sherpa will read in events as "<SetCheck(s["CHECK_LIBLOCK"].Get()); rpa->Init(); CheckVersion(); LoadLibraries(); ShowParameterSyntax(); ran->InitExternal(); rpa->gen.SetSoftSC(s["HADRON_DECAYS"]["Spin_Correlations"].Get()); rpa->gen.SetHardSC(s["HARD_DECAYS"]["Spin_Correlations"].Get()); exh->AddTerminatorObject(this); } void Initialization_Handler::RegisterDefaults() { Settings& s = Settings::GetMainSettings(); s["BEAM_REMNANTS"].SetDefault(true); s["EVENT_GENERATION_MODE"].SetDefault("PartiallyUnweighted"); s["EVENT_TYPE"].SetDefault("StandardPerturbative"); s["SOFT_COLLISIONS"].UseNoneReplacements().SetDefault("None"); s["BEAM_RESCATTERING"].UseNoneReplacements().SetDefault("None"); s["EVT_FILE_PATH"].SetDefault("."); s["ANALYSIS_OUTPUT"].SetDefault("Analysis/"); s["RESULT_DIRECTORY"].SetDefault("Results"); s["CHECK_LIBLOCK"].SetDefault(0); s["OUTPUT_PRECISION"].SetDefault(12); s["FILE_SIZE"].SetDefault(std::numeric_limits::max()); s["WRITE_REFERENCES_FILE"].SetDefault(true); s["MODEL"].SetDefault("SM"); s["FRAGMENTATION"].SetDefault("Ahadic").UseNoneReplacements(); s["HARD_DECAYS"]["Enabled"].SetDefault(false); s["N_COLOR"].SetDefault(3.0); std::string frag{ s["FRAGMENTATION"].Get() }; s["HADRON_DECAYS"]["Model"].SetDefault("HADRONS++").UseNoneReplacements(); s["HADRON_DECAYS"]["Max_Proper_Lifetime"].SetDefault(10.0); s["HADRON_DECAYS"]["Spin_Correlations"].SetDefault(0); auto hdenabled = s["HARD_DECAYS"]["Enabled"].Get(); s["HARD_DECAYS"]["Spin_Correlations"].SetDefault(hdenabled); s["EVENT_INPUT"].SetDefault(""); s["STATUS_PATH"].SetDefault(""); s["PATH"].SetDefault(""); s["SAVE_STATUS"].SetDefault(""); s.DeclareVectorSettingsWithEmptyDefault({ "EVENT_OUTPUT", "ANALYSIS", "SHERPA_LDADD", "SHERPA_VERSION", "PDF_LIBRARY", "PDF_SET", "PDF_SET_VERSIONS", "MPI_PDF_LIBRARY", "MPI_PDF_SET", "MPI_PDF_SET_VERSIONS", "BBR_PDF_LIBRARY", "BBR_PDF_SET", "BBR_PDF_SET_VERSIONS", "VARIATIONS", "SCALE_VARIATIONS", "PDF_VARIATIONS", "QCUT_VARIATIONS", "BUNCHES", "NLO_CSS_DISALLOW_FLAVOUR", "MASSIVE_PS", "MASSLESS_PS" }); s.DeclareMatrixSettingsWithEmptyDefault({ "CSS_ENHANCE", "ASSOCIATED_CONTRIBUTIONS_VARIATIONS" }); s["EVENT_OUTPUT"].UseNoneReplacements(); s["VARIATIONS"].UseNoneReplacements(); s["SCALE_VARIATIONS"].UseNoneReplacements(); s["PDF_VARIATIONS"].UseNoneReplacements(); s["QCUT_VARIATIONS"].UseNoneReplacements().SetSynonyms({"CKKW_VARIATIONS"}); s["PDF_LIBRARY"].UseNoneReplacements(); s["MPI_PDF_LIBRARY"].UseNoneReplacements(); s["BBR_PDF_LIBRARY"].UseNoneReplacements(); s["ANALYSIS"].UseNoneReplacements(); s["SHOW_ME_GENERATORS"].SetDefault(0); s["SHOW_PS_GENERATORS"].SetDefault(0); s["SHOW_NLOMC_GENERATORS"].SetDefault(0); s["SHOW_SHOWER_GENERATORS"].SetDefault(0); s["SHOW_KFACTOR_SYNTAX"].SetDefault(0); s["SHOW_SCALE_SYNTAX"].SetDefault(0); s["SHOW_SELECTOR_SYNTAX"].SetDefault(0); s["SHOW_MODEL_SYNTAX"].SetDefault(0); s["SHOW_FILTER_SYNTAX"].SetDefault(0); s["SHOW_ANALYSIS_SYNTAX"].SetDefault(0); s["SHOW_VARIABLE_SYNTAX"].SetDefault(0); s["SHOW_PDF_SETS"].SetDefault(0); s["ISR_E_ORDER"].SetDefault(1); s["ISR_E_SCHEME"].SetDefault(2); s["KFACTOR"].SetDefault("None").UseNoneReplacements(); s["SCALES"].SetDefault("METS{MU_F2}{MU_R2}{MU_Q2}"); s["SCALE_FACTOR"].SetDefault(1.0); s["FACTORIZATION_SCALE_FACTOR"].SetDefault(1.0); s["RENORMALIZATION_SCALE_FACTOR"].SetDefault(1.0); s["RESUMMATION_SCALE_FACTOR"].SetDefault(1.0); s["USR_WGT_MODE"].SetDefault(true); Scoped_Settings metssettings{ Settings::GetMainSettings()["METS"] }; metssettings["CLUSTER_MODE"].SetDefault(0); s["NNLOqT_FOMODE"].SetDefault(0); // m_mtmode != 0 to reweight the whole cross section by full mt-dependent LO // gg->H cross section and add mt-dependence higher order correcions of the // Wilson coefficient for the ggH coupling s["HNNLO_MTOP_MODE"].SetDefault(0); // m_kfmode = [001]_2 to enable factorized matching of the Wilson coefficient for ggH coupling; // = [010]_2 to enable individual matching of the Wilson coefficient for ggH coupling; // = [100]_2 to remove delta(pT) part of NNLO K factor for a separate LO parton shower. s["HNNLO_KF_MODE"].SetDefault(0); // shower settings (shower classes are rarely singletons, so we either // register settings here or we prevent SetDefault... to called more than once // otherwise s["SHOWER_GENERATOR"].SetDefault("CSS").UseNoneReplacements(); std::string showergen{ s["SHOWER_GENERATOR"].Get() }; s["JET_CRITERION"].SetDefault(showergen); s["NLOMC_GENERATOR"].SetDefault(showergen); s["CSS_EVOLUTION_SCHEME"].SetDefault(30+30*100); s["CSS_KFACTOR_SCHEME"].SetDefault(1); s["CSS_SCALE_SCHEME"].SetDefault(0); s["CSS_SCALE_VARIATION_SCHEME"].SetDefault(1); // TODO: Should this be set to 3.0 for the new Dire default? See the manual // Sherpa section on master for details s["CSS_FS_PT2MIN"].SetDefault(1.0); s["CSS_IS_PT2MIN"].SetDefault(2.0); s["CSS_PT2MIN_GSPLIT_FACTOR"].SetDefault(1.0); s["CSS_FS_AS_FAC"].SetDefault(1.0); s["CSS_IS_AS_FAC"].SetDefault(0.5); s["CSS_PDF_FAC"].SetDefault(1.0); s["CSS_SCALE_FACTOR"].SetDefault(1.); s["CSS_MASS_THRESHOLD"].SetDefault(0.0); s["CSS_FORCED_DECAYS"].SetDefault(true); s["CSS_FORCED_GLUON_SCALING"].SetDefault(-3./2.); s["VIRTUAL_EVALUATION_FRACTION"].SetDefault(1.0); s["CSS_RECO_CHECK"].SetDefault(0); s["CSS_MAXEM"].SetDefault(std::numeric_limits::max()); s["CSS_REWEIGHT"].SetDefault(showergen != "None"); s["OUTPUT_ME_ONLY_VARIATIONS"].SetDefault(showergen != "None"); s["CSS_MAX_REWEIGHT_FACTOR"].SetDefault(1e3); s["REWEIGHT_MCATNLO_EM"].SetDefault(1); s["CSS_REWEIGHT_SCALE_CUTOFF"].SetDefault(5.0); s["CSS_KIN_SCHEME"].SetDefault(0); s["NLO_CSS_KIN_SCHEME"].SetDefault(0); s["CSS_OEF"].SetDefault(3.0); s["CSS_KMODE"].SetDefault(2); s["CSS_RESPECT_Q2"].SetDefault(false); s["CSS_CKFMODE"].SetDefault(1); s["CSS_PDFCHECK"].SetDefault(1); s["CSS_QCD_MODE"].SetDefault(1); s["CSS_EW_MODE"].SetDefault(false); s["CSS_RECO_DECAYS"].SetDefault(0); s["CSS_MAXPART"].SetDefault(std::numeric_limits::max()); s["CSS_PDF_MIN"].SetDefault(1.0e-4); s["CSS_PDF_MIN_X"].SetDefault(1.0e-2); s["CSS_WEIGHT_CHECK"].SetDefault(false); s["CSS_CMODE"].SetDefault(1); s["CSS_NCOL"].SetDefault(3); s["CSS_RECALC_FACTOR"].SetDefault(4.0); s["CSS_TC_ENHANCE"].SetDefault(1.0); s["CSS_COUPLING_SCHEME"].SetDefault(1); s["CSS_ME_CORRECTION"].SetDefault(0); s["CSS_KERNEL_TYPE"].SetDefault(15); s["NLO_CSS_RECALC_FACTOR"].SetDefault(2.0); s["NLO_CSS_PSMODE"].SetDefault(0); s["NLO_CSS_WEIGHT_CHECK"].SetDefault(0); s["NLO_CSS_MAXEM"].SetDefault(1); s["MI_CSS_KFACTOR_SCHEME"].SetDefault(0); s["MI_CSS_IS_PT2MIN"].SetDefault(4.0); s["MI_CSS_FS_PT2MIN"].SetDefault(1.0); s["MI_CSS_PT2MIN_GSPLIT_FACTOR"].SetDefault(1.0); s["MI_CSS_IS_AS_FAC"].SetDefault(0.66); s["MI_CSS_FS_AS_FAC"].SetDefault(0.66); s["MI_CSS_KIN_SCHEME"].SetDefault(1); s["COMIX_DEFAULT_GAUGE"].SetDefault(1); s["DIPOLES"]["SCHEME"].SetDefault(subscheme::CSS); s["DIPOLES"]["KAPPA"].SetDefault(2.0/3.0); s["COUPLINGS"].SetDefault("Alpha_QCD 1"); s["EXTRAXS_CSS_APPROX_ME"].SetDefault(false); s["RESPECT_MASSIVE_FLAG"].SetDefault(false); } Initialization_Handler::~Initialization_Handler() { if (m_savestatus) { msg_Error()<gen.Variable("SHERPA_STATUS_PATH")<<"'."<gen.Variable("SHERPA_STATUS_PATH"),493); exh->PrepareTerminate(); } if (p_evtreader) { delete p_evtreader; p_evtreader = NULL; } if (p_mehandler) { delete p_mehandler; p_mehandler = NULL; } if (p_reconnections) { delete p_reconnections; p_reconnections = NULL; } if (p_fragmentation) { delete p_fragmentation; p_fragmentation = NULL; } if (p_beamremnants) { delete p_beamremnants; p_beamremnants = NULL; } if (p_harddecays) { delete p_harddecays; p_harddecays = NULL; } if (p_hdhandler) { delete p_hdhandler; p_hdhandler = NULL; } if (p_softphotons) { delete p_softphotons; p_softphotons = NULL; } if (p_beamspectra) { delete p_beamspectra; p_beamspectra = NULL; } if (p_model) { delete p_model; p_model = NULL; } if (p_variations) { delete p_variations; p_variations = NULL; } if (p_filter) { delete p_filter; p_filter = NULL; } while (m_analyses.size()>0) { delete m_analyses.back(); m_analyses.pop_back(); } while (m_outputs.size()>0) { delete m_outputs.back(); m_outputs.pop_back(); } while (m_isrhandlers.size()>0) { delete m_isrhandlers.begin()->second; m_isrhandlers.erase(m_isrhandlers.begin()); } while (m_mihandlers.size()>0) { delete m_mihandlers.begin()->second; m_mihandlers.erase(m_mihandlers.begin()); } while (m_remnanthandlers.size()>0) { REMNANTS::Remnant_Handler * rh = m_remnanthandlers.begin()->second; for (REMNANTS::Remnant_Handler_Map::iterator rit=m_remnanthandlers.begin(); rit!=m_remnanthandlers.end();) { if (rit->second==rh) m_remnanthandlers.erase(rit++); else rit++; } delete rh; } while (m_schandlers.size()>0) { delete m_schandlers.begin()->second; m_schandlers.erase(m_schandlers.begin()); } while (m_showerhandlers.size()>0) { delete m_showerhandlers.begin()->second; m_showerhandlers.erase(m_showerhandlers.begin()); } PHASIC::Phase_Space_Handler::DeleteInfo(); exh->RemoveTerminatorObject(this); for (set::iterator pdflib=m_pdflibs.begin(); pdflib!=m_pdflibs.end(); ++pdflib) { if (*pdflib=="None") continue; void *exit(s_loader->GetLibraryFunction(*pdflib,"ExitPDFLib")); if (exit==NULL) PRINT_INFO("Error: Cannot unload PDF library "+*pdflib); else ((PDF_Exit_Function)exit)(); } //delete m_defsets[PDF::isr::hard_process]; // = new std::string[2]; //delete m_defsets[PDF::isr::hard_subprocess]; // = new std::string[2]; //delete m_defsets[PDF::isr::bunch_rescatter]; // = new std::string[2]; } void Initialization_Handler::CheckVersion() { std::vector versioninfo{ Settings::GetMainSettings()["SHERPA_VERSION"].GetVector() }; if (versioninfo.empty()) return; std::string currentversion(ToString(SHERPA_VERSION)+"." +ToString(SHERPA_SUBVERSION)); if (versioninfo.size()==1 && versioninfo[0]!=currentversion) { THROW(normal_exit,"Run card request Sherpa "+versioninfo[0] +". This is Sherpa "+currentversion); } else if (versioninfo.size()==2) { if (versioninfo[0]==currentversion || versioninfo[1]==currentversion) return; size_t min1(versioninfo[0].find(".",0)), min2(versioninfo[0].find(".",min1+1)), max1(versioninfo[1].find(".",0)), max2(versioninfo[1].find(".",max1+1)); size_t minmajvers(ToType(versioninfo[0].substr(0,min1))), minminvers(ToType(versioninfo[0].substr(min1+1,min2))), minbugvers(ToType(versioninfo[0].substr(min2+1))), maxmajvers(ToType(versioninfo[1].substr(0,max1))), maxminvers(ToType(versioninfo[1].substr(max1+1,max2))), maxbugvers(ToType(versioninfo[1].substr(max2+1))), curmajvers(ToType(currentversion.substr(0,max1))), curminvers(ToType(currentversion.substr(max1+1,max2))), curbugvers(ToType(currentversion.substr(max2+1))); if (!(CompareVersions(minmajvers,minminvers,minbugvers, curmajvers,curminvers,curbugvers) *CompareVersions(curmajvers,curminvers,curbugvers, maxmajvers,maxminvers,maxbugvers))) THROW(normal_exit,"Run card request Sherpa "+versioninfo[0] +"-"+versioninfo[1] +". This is Sherpa "+currentversion); } else THROW(not_implemented,"SHERPA_VERSION information not recognised."); } bool Initialization_Handler::CompareVersions (const size_t& a1,const size_t& b1,const size_t& c1, const size_t& a2,const size_t& b2,const size_t& c2) { if (a1 ldadd = Settings::GetMainSettings()["SHERPA_LDADD"].GetVector(); for (size_t i(0);iLoadLibrary(ldadd[i])) { THROW(fatal_error,"Cannot load extra library."); } else msg_Info()<()); if (helpi>0) { msg->SetLevel(2); PHASIC::ME_Generator_Base::ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } helpi = s["SHOW_PS_GENERATORS"].Get(); if (helpi>0) { msg->SetLevel(2); PHASIC::Channel_Generator::ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } helpi = s["SHOW_NLOMC_GENERATORS"].Get(); if (helpi>0) { msg->SetLevel(2); PDF::NLOMC_Base::ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } helpi = s["SHOW_SHOWER_GENERATORS"].Get(); if (helpi>0) { msg->SetLevel(2); PDF::Shower_Base::ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } helpi = s["SHOW_KFACTOR_SYNTAX"].Get(); if (helpi>0) { msg->SetLevel(2); PHASIC::KFactor_Setter_Base::ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } helpi = s["SHOW_SCALE_SYNTAX"].Get(); if (helpi>0) { msg->SetLevel(2); if (helpi&1) PHASIC::Scale_Setter_Base::ShowSyntax(helpi); if (helpi&2) PHASIC::Core_Scale_Setter::ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } helpi = s["SHOW_SELECTOR_SYNTAX"].Get(); if (helpi>0) { msg->SetLevel(2); PHASIC::Selector_Base::ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } helpi = s["SHOW_MODEL_SYNTAX"].Get(); if (helpi>0) { msg->SetLevel(2); MODEL::Model_Base::ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } helpi = s["SHOW_FILTER_SYNTAX"].Get(); if (helpi>0) { msg->SetLevel(2); Filter::ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } helpi = s["SHOW_ANALYSIS_SYNTAX"].Get(); if (helpi>0) { msg->SetLevel(2); InitializeTheAnalyses(); for (Analysis_Vector::iterator it=m_analyses.begin(); it!=m_analyses.end(); ++it) (*it)->ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } helpi = s["SHOW_VARIABLE_SYNTAX"].Get(); if (helpi>0) { msg->SetLevel(2); ATOOLS::Variable_Base::ShowVariables(helpi); THROW(normal_exit,"Syntax shown."); } } std::string StripSectionTags(const std::string &name) { if (name.find('|')!=std::string::npos) return name.substr(0,name.find('|')); return name; } void Initialization_Handler::PrepareTerminate() { Settings& s = Settings::GetMainSettings(); std::string path(rpa->gen.Variable("SHERPA_STATUS_PATH")+"/"); if (path=="/") return; Copy(s.GetPath(), path + s.GetPath()); Data_Writer writer; writer.SetOutputFile(path+"cmd"); writer.SetVectorType(vtc::vertical); std::vector lines = { "SHERPA_RUN_PATH = "+rpa->gen.Variable("SHERPA_RUN_PATH"), "SHERPA_CPP_PATH = "+rpa->gen.Variable("SHERPA_CPP_PATH"), "SHERPA_LIB_PATH = "+rpa->gen.Variable("SHERPA_LIB_PATH"), }; writer.VectorToFile(lines); } bool Initialization_Handler::InitializeTheFramework(int nr) { Settings& s = Settings::GetMainSettings(); bool okay = true; const int defgauge{ s["COMIX_DEFAULT_GAUGE"].Get() }; Spinor::SetDefaultGauge(defgauge); Spinor::SetDefaultGauge(defgauge); SetGlobalVariables(); std::string stag(rpa->gen.Variable("RNG_SEED")); while (stag.find(' ')!=std::string::npos) stag.replace(stag.find(' '),1,"-"); s.AddTag("RNG_SEED", stag); Hadron_Init().Init(); okay = okay && InitializeTheModel(); if (m_mode==eventtype::StandardPerturbative) { std::string eventtype{ s["EVENT_TYPE"].Get() }; if (eventtype=="StandardPerturbative") m_mode=eventtype::StandardPerturbative; else if (eventtype=="MinimumBias") { m_mode=eventtype::MinimumBias; if (s["SOFT_COLLISIONS"].Get()==string("Amisic")) s["MI_HANDLER"].OverrideScalar("Amisic"); else if (s["SOFT_COLLISIONS"].Get()==string("Shrimps")) s["MI_HANDLER"].OverrideScalar("None"); s["ME_GENERATORS"].OverrideScalar("None"); } else if (eventtype=="HadronDecay") { m_mode=eventtype::HadronDecay; s["MI_HANDLER"].OverrideScalar("None"); s["ME_GENERATORS"].OverrideScalar("None"); } else { THROW(not_implemented,"Unknown event type '"+eventtype+"'"); } } okay = okay && InitializeTheBeams(); okay = okay && InitializeThePDFs(); okay = okay && InitializeTheRemnants(); if (!p_model->ModelInit(m_isrhandlers)) THROW(critical_error,"Model cannot be initialized"); p_model->InitializeInteractionModel(); if (!CheckBeamISRConsistency()) return 0.; if (m_mode==eventtype::EventReader) { std::string infile; size_t bpos(m_evtform.find('[')), epos(m_evtform.rfind(']')); if (bpos!=std::string::npos && epos!=std::string::npos) { infile=m_evtform.substr(bpos+1,epos-bpos-1); m_evtform=m_evtform.substr(0,bpos); } std::string libname(m_evtform); if (libname.find('_')) libname=libname.substr(0,libname.find('_')); if (!s_loader->LoadLibrary("Sherpa"+libname+"Input")) THROW(missing_module,"Cannot load input library Sherpa"+libname+"Input."); p_evtreader = Event_Reader_Base::Getter_Function::GetObject (m_evtform,Input_Arguments(s.GetPath(), infile, p_model, m_isrhandlers[isr::hard_process])); if (p_evtreader==NULL) THROW(fatal_error,"Event reader not found"); msg_Events()<<"SHERPA will read in the events."<gen.NumberOfEvents()>0) { okay = okay && InitializeTheUnderlyingEvents(); okay = okay && InitializeTheSoftCollisions(); okay = okay && InitializeTheColourReconnections(); okay = okay && InitializeTheFragmentation(); okay = okay && InitializeTheHadronDecays(); okay = okay && InitializeTheSoftPhotons(); okay = okay && InitializeTheIO(); okay = okay && InitializeTheFilter(); okay = okay && InitializeTheReweighting(Variations_Mode::all); okay = okay && InitializeTheAnalyses(); } else { okay = okay && InitializeTheReweighting(Variations_Mode::nominal_only); } return okay; } bool Initialization_Handler::CheckBeamISRConsistency() { if (p_model->Name()==std::string("ADD")) { double ms = p_model->ScalarConstant("M_s"); if (msgen.Ecms()) { msg_Error()<<"WARNING in Initialization_Handler::CheckBeamISRConsistency :"<gen.Ecms()); string name=p_model->Name(); if (name==std::string("ADD")) { double mcut2 = sqr(p_model->ScalarConstant("M_cut")); // if ISR & beam -> apply mcut on ISR only // if beam only -> apply mcut on Beam smax = Min(smax,mcut2); for (size_t i=1;i<3;++i) { isr::id id=(isr::id)i; if (m_isrhandlers[id]->On()) { m_isrhandlers[id]->SetFixedSprimeMax(smax); m_isrhandlers[id]->SetFixedSprimeMin(smin); } else if (p_beamspectra->On()) { p_beamspectra->SetSprimeMax(smax); } } } if (!(p_beamspectra->CheckConsistency(m_bunch_particles))) { msg_Error()<<"Error in Initialization of the Sherpa framework : "<GetBeam(0)<<" -> " <Flav(0)<<" and " <GetBeam(1)<<" -> " <Flav(1)<(); std::string outpath=s["EVT_FILE_PATH"].Get(); for (size_t i=0; iLoadLibrary("Sherpa"+libname+"Output")) THROW(missing_module,"Cannot load output library Sherpa"+libname+"Output."); out=Output_Base::Getter_Function::GetObject (outputs[i], Output_Arguments(outpath, outfile)); } if (out==NULL) THROW(fatal_error,"Cannot initialize "+outputs[i]+" output"); m_outputs.push_back(out); } return true; } bool Initialization_Handler::InitializeTheModel() { Settings& s = Settings::GetMainSettings(); if (p_model) delete p_model; std::string name(s["MODEL"].Get()); p_model=Model_Base::Model_Getter_Function:: GetObject(name, Model_Arguments(true)); if (p_model==NULL) { if (!s_loader->LoadLibrary("Sherpa"+name)) THROW(missing_module,"Cannot load model library Sherpa"+name+"."); p_model=Model_Base::Model_Getter_Function:: GetObject(name, Model_Arguments(true)); } if (p_model==NULL) THROW(not_implemented,"Model not implemented"); MODEL::s_model=p_model; return 1; } bool Initialization_Handler::InitializeTheBeams() { if (p_beamspectra) { delete p_beamspectra; p_beamspectra = NULL; } p_beamspectra = new Beam_Spectra_Handler(); p_beamspectra->Output(); return 1; } bool Initialization_Handler::InitializeThePDFs() { Settings& settings = Settings::GetMainSettings(); // Load all necessary PDF libraries and check that they are available LoadPDFLibraries(settings); // Define bunch flavours DefineBunchFlavours(settings); // Initialisation of PDF sets for (size_t i=1;i<4;++i) InitISRHandler((isr::id)(i),settings); msg_Info()<<"ISR handling:"<()!=string("None"); for (size_t pid=1;pid<4;pid++) { PDF::isr::id pc; if (pid==1) { msg_Info()<<" PDFs for hard scattering: "; pc = PDF::isr::hard_process; } if (pid==2) { msg_Info()<<" PDFs for multiple parton interactions: "; pc = PDF::isr::hard_subprocess; } if (pid==3 && needs_resc) { msg_Info()<<" PDFs for beam re-scattering: "; pc = PDF::isr::bunch_rescatter; } if (pid!=3 || (pid==3 && needs_resc)) { for (size_t beam=0;beam<2;beam++) { if (m_isrhandlers[pc]->Type(beam)==PDF::isrtype::lepton || m_isrhandlers[pc]->Type(beam)==PDF::isrtype::hadron) { msg_Info() << m_isrhandlers[pc]->PDF(beam)->Set(); } else { msg_Info() << "None"; } if (beam==0 && (m_isrhandlers[pc]->Type(1)==PDF::isrtype::lepton || m_isrhandlers[pc]->Type(1)==PDF::isrtype::hadron)) msg_Info()<<" + "; } msg_Info()<<"\n"; } } return 1; } void Initialization_Handler::LoadPDFLibraries(Settings& settings) { std::vector pdflibs = settings["PDF_LIBRARY"].GetVector(); std::vector mpilibs = settings["MPI_PDF_LIBRARY"].GetVector(); std::vector bbrlibs = settings["BBR_PDF_LIBRARY"].GetVector(); m_defsets[PDF::isr::hard_process] = std::array(); m_defsets[PDF::isr::hard_subprocess] = std::array(); m_defsets[PDF::isr::bunch_rescatter] = std::array(); for (size_t beam=0;beam<2;++beam) { ///////////////////////////////////////////////////////// // define bunch particle-dependent PDF libraries and sets here ///////////////////////////////////////////////////////// std::string deflib("None"), defset; if (p_beamspectra->GetBeam(beam)->Bunch(0).Kfcode()==kf_p_plus) { deflib = PDF::pdfdefs->DefaultPDFLibrary(kf_p_plus); defset = PDF::pdfdefs->DefaultPDFSet(kf_p_plus); } else if (p_beamspectra->GetBeam(beam)->Bunch(0).Kfcode()==kf_e || p_beamspectra->GetBeam(beam)->Bunch(0).Kfcode()==kf_mu) { deflib = PDF::pdfdefs->DefaultPDFLibrary(kf_e); defset = PDF::pdfdefs->DefaultPDFSet(kf_e); } else if (p_beamspectra->GetBeam(beam)->Bunch(0).IsPhoton()) { deflib = PDF::pdfdefs->DefaultPDFLibrary(kf_photon); defset = PDF::pdfdefs->DefaultPDFSet(kf_photon); } // fix PDFs and default sets for the hard_process here if (pdflibs.empty()) m_pdflibs.insert(deflib); else m_pdflibs.insert(pdflibs[Min(beam,pdflibs.size()-1)]); m_defsets[PDF::isr::hard_process][beam] = defset; // fix PDFs and default sets for the MPI's / hard_subprocesses here // we may have to define defaults here. if (!mpilibs.empty()) { std::string libname = mpilibs[Min(beam,mpilibs.size()-1)]; if (m_pdflibs.find(libname)==m_pdflibs.end()) m_pdflibs.insert(libname); m_defsets[PDF::isr::hard_subprocess][beam] = defset; } else m_defsets[PDF::isr::hard_subprocess][beam] = defset; // fix PDFs and default sets for the beam rescattering here // EPA is the only configuration at the moment where we allow additional // scattering/interactions of the incoming beams if (m_mode==eventtype::StandardPerturbative && settings["BEAM_RESCATTERING"].Get()!=string("None") && p_beamspectra->GetBeam(beam)->Beam().IsHadron() && p_beamspectra->GetBeam(beam)->Bunch(0).Kfcode()==kf_photon && p_beamspectra->GetBeam(beam)->Bunch(1)==p_beamspectra->GetBeam(beam)->Beam()) { if (!bbrlibs.empty()) { std::string libname = bbrlibs[Min(beam,bbrlibs.size()-1)]; if (m_pdflibs.find(libname)==m_pdflibs.end()) m_pdflibs.insert(libname); m_defsets[PDF::isr::bunch_rescatter][beam] = std::string("None"); } else { m_pdflibs.insert(PDF::pdfdefs->DefaultPDFLibrary(kf_p_plus)); m_defsets[PDF::isr::bunch_rescatter][beam] = PDF::pdfdefs->DefaultPDFSet(kf_p_plus); } } else m_defsets[PDF::isr::bunch_rescatter][beam] = string(""); } // add LHAPDF if necessary and load the relevant libraries if (Variations::NeedsLHAPDF6Interface()) { m_pdflibs.insert("LHAPDFSherpa"); } for (set::iterator pdflib=m_pdflibs.begin(); pdflib!=m_pdflibs.end();++pdflib) { if (*pdflib=="None") continue; if (*pdflib=="LHAPDFSherpa") { #ifdef USING__LHAPDF s_loader->AddPath(std::string(LHAPDF_PATH)+"/lib"); s_loader->LoadLibrary("LHAPDF"); #else THROW(fatal_error, "LHAPDF has not been enabled during configuration.") #endif } void *init(s_loader->GetLibraryFunction(*pdflib,std::string("InitPDFLib"))); if (init==NULL) THROW(fatal_error,"Cannot load PDF library "+*pdflib); ((PDF_Init_Function)init)(); } // PDF set listing output int helpi{ settings["SHOW_PDF_SETS"].Get() }; if (helpi>0) { msg->SetLevel(2); PDF::PDF_Base::ShowSyntax(helpi); THROW(normal_exit,"Syntax shown."); } } void Initialization_Handler::InitISRHandler(const PDF::isr::id & pid,Settings& settings) { if (m_isrhandlers.find(pid)!=m_isrhandlers.end()) delete m_isrhandlers[pid]; bool needs_resc = settings["BEAM_RESCATTERING"].Get()!=string("None"); ///////////////////////////////////////////////////////////// // make sure rescatter ISR bases are only initialised if necessary ///////////////////////////////////////////////////////////// if (pid==PDF::isr::bunch_rescatter && !needs_resc) return; std::string tag = ( pid==PDF::isr::hard_process ? string("PDF_SET") : pid==PDF::isr::hard_subprocess ? string("MPI_PDF_SET") : string("BBR_PDF_SET") ); std::string vtag = tag+string("_VERSIONS"); ///////////////////////////////////////////////////////////// // read sets and versions for relevant part of event generation // (hard process, MPI, beam rescattering) ///////////////////////////////////////////////////////////// std::vector sets = settings[tag].GetVector(); std::vector versions = settings[vtag].GetVector(); if (sets.size() > 2) { THROW(fatal_error, "You can not specify more than two PDF sets."); } if (versions.size() > 2) { THROW(fatal_error, "You can not specify more than two PDF set versions."); } std::array isrbases = {}; for (size_t beam=0;beam<2;beam++) { isrbases[beam] = NULL; // fix actual set and version for beam number and part of // event generation here - assume central (i.e. version 0 by default) std::string set = ( (sets.size()==0 || sets[Min(beam,sets.size()-1)]=="Default") ? m_defsets[pid][beam] : sets[Min(beam,sets.size()-1)] ); int version = (versions.size()== 0 ? 0 : versions[Min(beam,versions.size()-1)] ); int order = -1, scheme = -1; // special treatment of electron PDF if ((pid==PDF::isr::hard_process || pid==PDF::isr::hard_subprocess) && set == "PDFe") { order = settings["ISR_E_ORDER"].Get(); scheme = settings["ISR_E_SCHEME"].Get(); } // Initialise the actual PDF and make sure it arrives at the right place. // Here we need to be careful about flavours .... (I think I have it ...) Flavour flav = ( (pid != PDF::isr::bunch_rescatter) ? m_bunch_particles[beam] : p_beamspectra->GetBeam(beam)->Beam() ); PDF_Arguments args = PDF_Arguments(flav, beam, set, version, order, scheme); if (pid!=PDF::isr::bunch_rescatter) { PDF_Base * pdfbase = PDF_Base::PDF_Getter_Function::GetObject(set,args); if (m_bunch_particles[beam].IsHadron() && pdfbase==NULL) THROW(critical_error,"PDF '"+set+"' does not exist in any of the loaded" +" libraries for "+ToString(m_bunch_particles[beam])+" bunch."); if (pid==PDF::isr::hard_process) rpa->gen.SetPDF(beam,pdfbase); if (pdfbase==NULL) { isrbases[beam] = new Intact(flav); needs_resc = false; } else { pdfbase->SetBounds(); isrbases[beam] = new Structure_Function(pdfbase,flav); } ATOOLS::rpa->gen.SetBunch(m_bunch_particles[beam],beam); } else if (pid==PDF::isr::bunch_rescatter && needs_resc) { PDF_Base * pdfbase = PDF_Base::PDF_Getter_Function::GetObject(set,args); if (pdfbase==NULL) THROW(critical_error,"PDF '"+set+"' for rescattering does not exist in any of the loaded" +" libraries for "+ToString(m_bunch_particles[beam])+" bunch."); pdfbase->SetBounds(); isrbases[beam] = new Structure_Function(pdfbase,flav); } } if ((pid==PDF::isr::bunch_rescatter && needs_resc) || pid!=isr::bunch_rescatter) { ISR_Handler * isr = new ISR_Handler(isrbases,pid); for (size_t beam=0;beam<2;beam++) isr->SetBeam(p_beamspectra->GetBeam(beam),beam); isr->Init(); if (!(p_beamspectra->CheckConsistency(m_bunch_particles))) { msg_Error()<<"Error in Environment::InitializeThePDFs()"<Output(); m_isrhandlers[pid] = isr; /////////////////////////////////////////////////////////// // This is a bit of an ad-hoc fix for bunch rescattering in EPA configurations only. // The tags we fill here will by default be bunchtags = {0,0}, i.e. look at the // first pair of bunches coming out of the beam bases. // For EPA we will add {1,1} - effectively then also adding the hadrons. // These tags will be handed over to the remnants, making sure we have // one remnant per state that has a PDF. // TODO: Extend this for ions - this will need even more restructuring /////////////////////////////////////////////////////////// vector bunchtags; bunchtags.resize(2,pid==isr::bunch_rescatter ? 1 : 0); m_bunchtags[pid] = bunchtags; } } void Initialization_Handler::DefineBunchFlavours(Settings& settings) { std::vector bunches{ settings["BUNCHES"].GetVector() }; if (bunches.size() > 2) { THROW(fatal_error, "You can not specify more than two bunches."); } for (size_t beam=0;beam<2;beam++) { if (bunches.empty()) m_bunch_particles[beam] = p_beamspectra->GetBeam(beam)->Bunch(0); else { int flav = bunches[Min(beam,bunches.size()-1)]; m_bunch_particles[beam] = Flavour((kf_code)abs(flav)); if (flav<0) m_bunch_particles[beam] = m_bunch_particles[beam].Bar(); } } } bool Initialization_Handler::InitializeTheRemnants() { /////////////////////////////////////////////////////////// // define two sets of remnants, if necessary (i.e. if we have bunch rescattering). // we will have to make sure that the remnants for hard process and hard subprocess - // the MPI related to the hard process - are the same. // I have the feeling we will have to communicate the mode to the Remnant_Handler in question /////////////////////////////////////////////////////////// REMNANTS::Remnants_Parameters(); m_remnanthandlers[isr::hard_process] = new Remnant_Handler(m_isrhandlers[isr::hard_process],p_beamspectra, m_bunchtags[isr::hard_process]); m_remnanthandlers[isr::hard_subprocess] = m_remnanthandlers[isr::hard_process]; if (m_isrhandlers.find(isr::bunch_rescatter)!=m_isrhandlers.end()) { m_remnanthandlers[isr::bunch_rescatter] = new Remnant_Handler(m_isrhandlers[isr::bunch_rescatter],p_beamspectra, m_bunchtags[isr::bunch_rescatter]); } msg_Info()<<"Remnant_Handlers:\n" <<" hard process: " <GetRemnant(0)->GetBeam()->Bunch(0)<<": " <GetRemnant(0)->Type()<<" + " <GetRemnant(1)->GetBeam()->Bunch(0)<<": " <GetRemnant(1)->Type()<<"\n"; if (m_remnanthandlers.find(isr::bunch_rescatter)!=m_remnanthandlers.end()) msg_Info()<<" rescattering: " <GetRemnant(0)->GetBeam()->Bunch(1)<<": " <GetRemnant(0)->Type()<<" + " <GetRemnant(1)->GetBeam()->Bunch(1)<<": " <GetRemnant(1)->Type()<<"\n"; return true; } bool Initialization_Handler::InitializeTheHardDecays() { if (!Settings::GetMainSettings()["HARD_DECAYS"]["Enabled"].Get()) return true; if (p_harddecays) { delete p_harddecays; p_harddecays = NULL; } p_harddecays = new Hard_Decay_Handler(); return true; } bool Initialization_Handler::InitializeTheMatrixElements() { #ifdef USING__EWSud // in case that KFACTOR=EWSud is used we need to be ready when the ME handler // sets up the KFactor setters if (!s_loader->LoadLibrary("SherpaEWSud")) THROW(missing_module,"Cannot load EWsud library."); #endif if (p_mehandler) delete p_mehandler; p_mehandler = new Matrix_Element_Handler(p_model); p_mehandler->SetShowerHandler(m_showerhandlers[isr::hard_process]); p_mehandler->SetRemnantHandler(m_remnanthandlers[isr::hard_process]); auto ret = p_mehandler->InitializeProcesses(p_beamspectra, m_isrhandlers[isr::hard_process]); msg_Info()<<"Initialized the Matrix_Element_Handler for the hard processes." < isrtypes; isrtypes.push_back(isr::hard_process); isrtypes.push_back(isr::hard_subprocess); if (m_isrhandlers.find(isr::bunch_rescatter)!=m_isrhandlers.end()) isrtypes.push_back(isr::bunch_rescatter); for (size_t i=0; iSetActiveAs(id); Shower_Handler_Map::iterator it=m_showerhandlers.find(id); if (it!=m_showerhandlers.end()) delete it->second; m_showerhandlers[id] = new Shower_Handler(p_model, m_isrhandlers[id], i); m_showerhandlers[id]->SetRemnants(m_remnanthandlers[id]); for (size_t beam=0;beam<2;beam++) { m_isrhandlers[id]->SetRemnant(m_remnanthandlers[id]->GetRemnant(beam),beam); } } as->SetActiveAs(isr::hard_process); msg_Info()<<"Shower_Handler initialised.\n"< isrtypes; isrtypes.push_back(isr::hard_subprocess); if (m_isrhandlers.find(isr::bunch_rescatter)!=m_isrhandlers.end()) isrtypes.push_back(isr::bunch_rescatter); for (size_t i=0; iSetActiveAs(isr::hard_subprocess); MI_Handler * mih = new MI_Handler(p_model,m_isrhandlers[id],m_remnanthandlers[id]); mih->SetShowerHandler(m_showerhandlers[id]); as->SetActiveAs(isr::hard_process); m_mihandlers[id] = mih; } msg_Info()<<"Underlying event/multiple interactions handler:\n"; for (size_t i=0; iName()<<")\n"; } return true; } bool Initialization_Handler::InitializeTheSoftCollisions() { /////////////////////////////////////////////////////////// // define up to two soft collision handlers - // they will have to differ in how they fill blobs. // modify the beam remnants to take care of the rescatter /////////////////////////////////////////////////////////// std::vector isrtypes; isrtypes.push_back(isr::hard_subprocess); if (m_isrhandlers.find(isr::bunch_rescatter)!=m_isrhandlers.end()) isrtypes.push_back(isr::bunch_rescatter); for (size_t i=0; iOn() ? new Soft_Collision_Handler(mih->Amisic(),mih->Shrimps(), id==isr::bunch_rescatter) : NULL ); if (id==isr::bunch_rescatter) { p_beamremnants->AddBunchRescattering(m_remnanthandlers[isr::bunch_rescatter], m_schandlers[isr::bunch_rescatter]); } } msg_Info()<<"Soft-collision handlers:\n"; for (size_t i=0;iSoft_CollisionModel()<<"\n"; } return true; } bool Initialization_Handler::InitializeTheBeamRemnants() { if (p_beamremnants) delete p_beamremnants; p_beamremnants = new Beam_Remnant_Handler(p_beamspectra, m_remnanthandlers[isr::hard_process], m_schandlers[isr::hard_subprocess]); msg_Info()<<"Initialized the Beam_Remnant_Handler."<Output(); return 1; } bool Initialization_Handler::InitializeTheFragmentation() { if (p_fragmentation) { delete p_fragmentation; p_fragmentation = NULL; } as->SetActiveAs(isr::hard_subprocess); Settings& s = Settings::GetMainSettings(); string fragmentationmodel = s["FRAGMENTATION"].Get(); if (fragmentationmodel!="None") { if (!s["BEAM_REMNANTS"].Get()) msg_Error()<Tracking()); } p_fragmentation = Fragmentation_Getter::GetObject (fragmentationmodel, Fragmentation_Getter_Parameters(m_showerhandlers[isr::hard_process]->ShowerGenerator())); if (p_fragmentation==NULL) THROW(fatal_error, " Fragmentation model '"+fragmentationmodel+"' not found."); as->SetActiveAs(isr::hard_process); msg_Info()<<"Initialized the Fragmentation_Handler."<() == "None") return true; std::string decmodel{ s["HADRON_DECAYS"]["Model"].Get() }; msg_Tracking()<<"Decaymodel = "<SetActiveAs(isr::hard_subprocess); Hadron_Decay_Handler* hd=new Hadron_Decay_Handler(); as->SetActiveAs(isr::hard_process); p_hdhandler=hd; } else { THROW(fatal_error,"Hadron decay model '"+decmodel+"' not implemented."); } msg_Info()<<"Initialized the Hadron_Decay_Handler, Decay model = " <SetSoftPhotonHandler(p_softphotons); if (p_hdhandler) p_hdhandler->SetSoftPhotonHandler(p_softphotons); msg_Info()<<"Initialized the Soft_Photon_Handler."<(); const Analysis_Arguments args{ Analysis_Arguments(outpath) }; std::vector analyses=s["ANALYSIS"].GetVector(); for (size_t i=0; iLoadLibrary("SherpaAnalysis")) THROW(missing_module,"Cannot load Analysis library (-DSHERPA_ENABLE_ANALYSIS=ON)."); if (analyses[i]=="Rivet" || analyses[i]=="RivetME" || analyses[i]=="RivetShower") { bool hepmc_loaded {false}; #ifdef USING__HEPMC2 hepmc_loaded |= (s_loader->LoadLibrary("SherpaHepMCOutput") != nullptr); #endif #ifdef USING__HEPMC3 hepmc_loaded |= (s_loader->LoadLibrary("SherpaHepMC3Output") != nullptr); #endif if (!hepmc_loaded) { THROW(missing_module, "Cannot load HepMC library (-DSHERPA_ENABLE_HEPMC2=ON and/or -DSHERPA_ENABLE_HEPMC3=ON)."); } if (!s_loader->LoadLibrary("SherpaRivetAnalysis")) THROW(missing_module,"Cannot load RivetAnalysis library (-DSHERPA_ENABLE_RIVET=ON)."); } Analysis_Interface* ana = Analysis_Interface::Analysis_Getter_Function::GetObject(analyses[i], args); if (ana==NULL) { if (!s_loader->LoadLibrary("Sherpa"+analyses[i]+"Analysis")) THROW(missing_module,"Cannot load Analysis library '"+analyses[i]+"'."); ana = Analysis_Interface::Analysis_Getter_Function::GetObject( analyses[i], args); if (ana==NULL) THROW(fatal_error,"Cannot initialize Analysis "+analyses[i]); } m_analyses.push_back(ana); } return true; } bool Initialization_Handler::InitializeTheReweighting(Variations_Mode mode) { if (p_variations) { delete p_variations; } if (mode != Variations_Mode::nominal_only) Variations::CheckConsistencyWithBeamSpectra(p_beamspectra); p_variations = new Variations(mode); s_variations = p_variations; if (p_mehandler) p_mehandler->InitializeTheReweighting(mode); if (mode != Variations_Mode::nominal_only) msg_Info()<<"Initialized the Reweighting."<Init()) { delete p_filter; p_filter = NULL; } return true; } bool Initialization_Handler::CalculateTheHardProcesses() { if (m_mode!=eventtype::StandardPerturbative) return true; msg_Events()<<"===================================================================\n" <<"Start calculating the hard cross sections. This may take some time.\n"; as->SetActiveAs(isr::hard_process); p_variations->DisableVariations(); int ok = p_mehandler->CalculateTotalXSecs(); p_variations->EnableVariations(); if (ok) { msg_Events()<<"Calculating the hard cross sections has been successful.\n" <<"====================================================================\n"; } else { msg_Events()<<"Calculating the hard cross sections failed. Check this carefully.\n" <<"=======================================================================\n"; } return ok; } void Initialization_Handler::SetGlobalVariables() { Settings& s = Settings::GetMainSettings(); double sf(s["SCALE_FACTOR"].Get()); double fsf(sf*s["FACTORIZATION_SCALE_FACTOR"].Get()); double rsf(sf*s["RENORMALIZATION_SCALE_FACTOR"].Get()); double qsf(sf*s["RESUMMATION_SCALE_FACTOR"].Get()); rpa->gen.SetVariable("FACTORIZATION_SCALE_FACTOR", ToString(fsf)); rpa->gen.SetVariable("RENORMALIZATION_SCALE_FACTOR", ToString(rsf)); rpa->gen.SetVariable("RESUMMATION_SCALE_FACTOR", ToString(qsf)); msg_Debugging()<gen.Variable("FACTORIZATION_SCALE_FACTOR")<<"\n" <<" ren scale: "<gen.Variable("RENORMALIZATION_SCALE_FACTOR")<<"\n" <<" res scale: "<gen.Variable("RESUMMATION_SCALE_FACTOR")<<"\n}\n"; // TODO: remove from rpa? double virtfrac = s["VIRTUAL_EVALUATION_FRACTION"].Get(); rpa->gen.SetVariable("VIRTUAL_EVALUATION_FRACTION", ToString(virtfrac)); std::string evtm{ s["EVENT_GENERATION_MODE"].Get() }; int eventmode{ 0 }; if (evtm=="Unweighted" || evtm=="U") eventmode=1; else if (evtm=="PartiallyUnweighted" || evtm=="P") eventmode=2; rpa->gen.SetVariable("EVENT_GENERATION_MODE",ToString(eventmode)); }