diff --git a/Analysis/NLORivetAnalysis.cc b/Analysis/NLORivetAnalysis.cc --- a/Analysis/NLORivetAnalysis.cc +++ b/Analysis/NLORivetAnalysis.cc @@ -1,309 +1,310 @@ // -*- C++ -*- // // This is the implementation of the non-inlined, non-templated member // functions of the NLORivetAnalysis class. // +#include #include "NLORivetAnalysis.h" #include "ThePEG/Interface/ClassDocumentation.h" #include "ThePEG/Interface/ParVector.h" #include "ThePEG/Interface/Parameter.h" #include "ThePEG/Interface/Switch.h" #include "ThePEG/Persistency/PersistentOStream.h" #include "ThePEG/Persistency/PersistentIStream.h" #include "ThePEG/Vectors/HepMCConverter.h" #include "ThePEG/Config/HepMCHelper.h" #include "ThePEG/EventRecord/Event.h" #include "ThePEG/EventRecord/SubProcess.h" #include "ThePEG/EventRecord/SubProcessGroup.h" #include "ThePEG/Repository/EventGenerator.h" #include "ThePEG/Repository/CurrentGenerator.h" -#include "HepMC/GenEvent.h" #include "Rivet/AnalysisHandler.hh" #include "Rivet/Tools/Logging.hh" -#include "HepMC/IO_AsciiParticles.h" -#include "HepMC/IO_GenEvent.h" -#include using namespace ThePEG; NLORivetAnalysis::NLORivetAnalysis() : _remnantId(82), _format(1),_unitchoice(), _geneventPrecision(16), debug(false), _rivet(), _nevent(0) {} -HepMC::GenEvent * NLORivetAnalysis::makeEvent(tEventPtr event, tSubProPtr sub, long no, +namespace { + +// Special anonymous function for creating a genEvent. +HepMC::GenEvent * makeEvent(tEventPtr event, tSubProPtr sub, long no, long remnantId, Energy eUnit, Length lUnit, - CrossSection xsec, CrossSection xsecErr) const { + CrossSection xsec, CrossSection xsecErr) { // generate beam particles const PPair& beam = event->incoming(); - HepMC::GenParticle * b1 = + HepMC::GenParticlePtr b1 = HepMCTraits::newParticle(beam.first->momentum(),beam.first->id(), 1,eUnit); - HepMC::GenParticle * b2 = + HepMC::GenParticlePtr b2 = HepMCTraits::newParticle(beam.second->momentum(),beam.second->id(), 1,eUnit); // generate remnants - HepMC::GenParticle * r1 = + HepMC::GenParticlePtr r1 = HepMCTraits::newParticle(beam.first->momentum() - sub->incoming().first->momentum(), - _remnantId,1,eUnit); - HepMC::GenParticle * r2 = + remnantId,1,eUnit); + HepMC::GenParticlePtr r2 = HepMCTraits::newParticle(beam.second->momentum() - sub->incoming().second->momentum(), - _remnantId,1,eUnit); + remnantId,1,eUnit); // generate outgoing particles - vector outgoing; + vector outgoing; for ( ParticleVector::const_iterator p = sub->outgoing().begin(); p != sub->outgoing().end(); ++p ) { outgoing.push_back(HepMCTraits::newParticle((**p).momentum(),(**p).id(), 1,eUnit)); } // generate one blob vertex - HepMC::GenVertex * vertex = HepMCTraits::newVertex(); + HepMC::GenVertexPtr vertex = HepMCTraits::newVertex(); HepMCTraits::addIncoming(*vertex,b1); HepMCTraits::addIncoming(*vertex,b2); HepMCTraits::addOutgoing(*vertex,r1); HepMCTraits::addOutgoing(*vertex,r2); - for ( vector::const_iterator p = outgoing.begin(); + for ( vector::const_iterator p = outgoing.begin(); p != outgoing.end(); ++p ) HepMCTraits::addOutgoing(*vertex,*p); HepMC::GenEvent * ev = HepMCTraits::newEvent(no,event->weight()*sub->groupWeight(), event->optionalWeights()); HepMCTraits::setUnits(*ev,eUnit,lUnit); HepMCTraits::setBeamParticles(*ev,b1,b2); HepMCTraits::addVertex(*ev,vertex); HepMCTraits::setCrossSection(*ev,xsec/picobarn, xsecErr/picobarn); return ev; } +} void NLORivetAnalysis::analyze(ThePEG::tEventPtr event, long ieve, int loop, int state) { Energy eUnit; Length lUnit; switch (_unitchoice) { default: eUnit = GeV; lUnit = millimeter; break; case 1: eUnit = MeV; lUnit = millimeter; break; case 2: eUnit = GeV; lUnit = centimeter; break; case 3: eUnit = MeV; lUnit = centimeter; break; } tcEHPtr eh = dynamic_ptr_cast(event->primaryCollision()->handler()); assert(eh); CrossSection xsec = eh->integratedXSec(); CrossSection xsecErr = eh->integratedXSecErr(); tSubProPtr sub = event->primarySubProcess(); Ptr::tptr grp = dynamic_ptr_cast::tptr>(sub); AnalysisHandler::analyze(event, ieve, loop, state); // Rotate to CMS, extract final state particles and call analyze(particles). // convert to hepmc HepMC::GenEvent * hepmc = - makeEvent(event,sub,_nevent,eUnit,lUnit,xsec,xsecErr); + makeEvent(event,sub,_nevent,_remnantId,eUnit,lUnit,xsec,xsecErr); CurrentGenerator::Redirect stdout(cout); if(_rivet){ #if ThePEG_RIVET_VERSION == 1 _rivet->analyze(*hepmc); #elif ThePEG_RIVET_VERSION > 1 try { _rivet->analyze(*hepmc); } catch (const YODA::Exception & e) { Throw() << "Warning: Rivet/Yoda got the exception: "<< e.what()<<"\n" << Exception::warning; } #else #error "Unknown ThePEG_RIVET_VERSION" #endif } // delete hepmc event delete hepmc; if ( grp ) { for ( SubProcessVector::const_iterator s = grp->dependent().begin(); s != grp->dependent().end(); ++s ) { - hepmc = makeEvent(event,*s,_nevent,eUnit,lUnit,xsec,xsecErr); + hepmc = makeEvent(event,*s,_nevent,_remnantId,eUnit,lUnit,xsec,xsecErr); if ( _rivet ){ #if ThePEG_RIVET_VERSION == 1 _rivet->analyze(*hepmc); #elif ThePEG_RIVET_VERSION > 1 try { _rivet->analyze(*hepmc); } catch (const YODA::Exception & e) { Throw() << "Warning: Rivet/Yoda got the exception: "<< e.what()<<"\n" << Exception::warning; } #else #error "Unknown ThePEG_RIVET_VERSION" #endif } // delete hepmc event delete hepmc; } } ++_nevent; } ThePEG::IBPtr NLORivetAnalysis::clone() const { return new_ptr(*this); } ThePEG::IBPtr NLORivetAnalysis::fullclone() const { return new_ptr(*this); } void NLORivetAnalysis::persistentOutput(ThePEG::PersistentOStream & os) const { os << _analyses << filename << debug; } void NLORivetAnalysis::persistentInput(ThePEG::PersistentIStream & is, int) { is >> _analyses >> filename >> debug; } ThePEG::ClassDescription NLORivetAnalysis::initNLORivetAnalysis; // Definition of the static class description member. void NLORivetAnalysis::Init() { static ThePEG::ClassDocumentation documentation ("The NLORivetAnalysis class is a simple class to allow analyses" " from the Rivet library to be called from ThePEG"); static ThePEG::ParVector interfaceAnalyses ("Analyses", "The names of the Rivet analyses to use", &NLORivetAnalysis::_analyses, -1, "", "","" "", false, false, ThePEG::Interface::nolimits); static Parameter interfaceRemnantId ("RemnantId", "Set the PDG id to be used for remnants.", &NLORivetAnalysis::_remnantId, 82, 0, 0, false, false, Interface::nolimits); static Parameter interfaceFilename ("Filename", #if ThePEG_RIVET_VERSION == 1 "The name of the file where the AIDA histograms are put. If empty, " "the run name will be used instead. '.aida' will in any case be " "appended to the file name.", #elif ThePEG_RIVET_VERSION > 1 "The name of the file where the YODA histograms are put. If empty, " "the run name will be used instead. '.yoda' will in any case be " "appended to the file name.", #else #error "Unknown ThePEG_RIVET_VERSION" #endif &NLORivetAnalysis::filename, "", true, false); static Switch interfaceDebug ("Debug", "Enable debug information from Rivet", &NLORivetAnalysis::debug, false, true, false); static SwitchOption interfaceDebugNo (interfaceDebug, "No", "Disable debug information.", false); static SwitchOption interfaceDebugYes (interfaceDebug, "Yes", "Enable debug information from Rivet.", true); interfaceAnalyses.rank(10); } void NLORivetAnalysis::dofinish() { AnalysisHandler::dofinish(); if( _nevent > 0 && _rivet ) { CurrentGenerator::Redirect stdout(cout); #if ThePEG_RIVET_VERSION > 2 - _rivet->setCrossSection(generator()->integratedXSec()/picobarn, - generator()->integratedXSecErr()/picobarn); + _rivet->setCrossSection(make_pair(generator()->integratedXSec()/picobarn, + generator()->integratedXSecErr()/picobarn)); #else _rivet->setCrossSection(generator()->integratedXSec()/picobarn); #endif _rivet->finalize(); string fname = filename; #if ThePEG_RIVET_VERSION == 1 if ( fname.empty() ) fname = generator()->path() + "/" + generator()->runName() + ".aida"; #elif ThePEG_RIVET_VERSION > 1 if ( fname.empty() ) fname = generator()->path() + "/" + generator()->runName() + ".yoda"; #else #error "Unknown ThePEG_RIVET_VERSION" #endif _rivet->writeData(fname); } delete _rivet; _rivet = nullptr; } void NLORivetAnalysis::doinit() { AnalysisHandler::doinit(); if(_analyses.empty()) throw ThePEG::Exception() << "Must have at least one analysis loaded in " << "NLORivetAnalysis::doinitrun()" << ThePEG::Exception::runerror; // check that analysis list is available _rivet = new Rivet::AnalysisHandler; //(fname); _rivet->addAnalyses(_analyses); if ( _rivet->analysisNames().size() != _analyses.size() ) { throw ThePEG::Exception() << "Rivet could not find all requested analyses.\n" << "Use 'rivet --list-analyses' to check availability.\n" << ThePEG::Exception::runerror; } delete _rivet; _rivet = 0; } void NLORivetAnalysis::doinitrun() { AnalysisHandler::doinitrun(); // create NLORivet analysis handler CurrentGenerator::Redirect stdout(cout); _rivet = new Rivet::AnalysisHandler; //(fname); _rivet->addAnalyses(_analyses); // check that analysis list is still available if ( _rivet->analysisNames().size() != _analyses.size() ) { throw ThePEG::Exception() << "Rivet could not find all requested analyses.\n" << "Use 'rivet --list-analyses' to check availability.\n" << ThePEG::Exception::runerror; } if ( debug ) Rivet::Log::setLevel("Rivet",Rivet::Log::DEBUG); } diff --git a/Analysis/NLORivetAnalysis.h b/Analysis/NLORivetAnalysis.h --- a/Analysis/NLORivetAnalysis.h +++ b/Analysis/NLORivetAnalysis.h @@ -1,223 +1,216 @@ // -*- C++ -*- #ifndef THEPEG_NLORivetAnalysis_H #define THEPEG_NLORivetAnalysis_H // // This is the declaration of the NLORivetAnalysis class. // #include "ThePEG/Handlers/AnalysisHandler.h" #include "Rivet/AnalysisHandler.hh" namespace ThePEG { /** * Here is the documentation of the NLORivetAnalysis class. * * @see \ref NLORivetAnalysisInterfaces "The interfaces" * defined for NLORivetAnalysis. */ class NLORivetAnalysis: public ThePEG::AnalysisHandler { public: /** * The default constructor. */ NLORivetAnalysis(); public: /** @name Virtual functions required by the AnalysisHandler class. */ //@{ /** * Analyze a given Event. Note that a fully generated event * may be presented several times, if it has been manipulated in * between. The default version of this function will call transform * to make a lorentz transformation of the whole event, then extract * all final state particles and call analyze(tPVector) of this * analysis object and those of all associated analysis objects. The * default version will not, however, do anything on events which * have not been fully generated, or have been manipulated in any * way. * @param event pointer to the Event to be analyzed. * @param ieve the event number. * @param loop the number of times this event has been presented. * If negative the event is now fully generated. * @param state a number different from zero if the event has been * manipulated in some way since it was last presented. */ virtual void analyze(ThePEG::tEventPtr event, long ieve, int loop, int state); - /** - * Produca a HepMC event for the given subprocess - */ - HepMC::GenEvent * makeEvent(tEventPtr event, tSubProPtr sub, long no, - Energy eUnit, Length lUnit, - CrossSection xsec, CrossSection xsecErr) const; - //@} public: /** @name Functions used by the persistent I/O system. */ //@{ /** * Function used to write out object persistently. * @param os the persistent output stream written to. */ void persistentOutput(ThePEG::PersistentOStream & os) const; /** * Function used to read in object persistently. * @param is the persistent input stream read from. * @param version the version number of the object when written. */ void persistentInput(ThePEG::PersistentIStream & is, int version); //@} /** * The standard Init function used to initialize the interfaces. * Called exactly once for each class by the class description system * before the main function starts or * when this class is dynamically loaded. */ static void Init(); protected: /** @name Clone Methods. */ //@{ /** * Make a simple clone of this object. * @return a pointer to the new object. */ virtual ThePEG::IBPtr clone() const; /** Make a clone of this object, possibly modifying the cloned object * to make it sane. * @return a pointer to the new object. */ virtual ThePEG::IBPtr fullclone() const; //@} protected: /** @name Standard Interfaced functions. */ //@{ /** * Initialize this object. Called in the read phase. */ virtual void doinit(); /** * Initialize this object. Called in the run phase just before * a run begins. */ virtual void doinitrun(); /** * Finalize this object. Called in the run phase just after a * run has ended. Used eg. to write out statistics. */ virtual void dofinish(); //@} private: /** * The static object used to initialize the description of this class. * Indicates that this is a concrete class with persistent data. */ static ThePEG::ClassDescription initNLORivetAnalysis; /** * The assignment operator is private and must never be called. * In fact, it should not even be implemented. */ NLORivetAnalysis & operator=(const NLORivetAnalysis &) = delete; private: /** * The PDG ID to be used for remnants */ long _remnantId; /** * The HepMC format */ int _format; /** * Selector for the choice of units */ int _unitchoice; /** * Choice of output precision in GenEvent format */ unsigned int _geneventPrecision; /** * The Analyses to use */ vector _analyses; /** * The base name of the output file. */ string filename; /** * Enable debugging information from NLORivet */ bool debug; /** * The NLORivetAnalysisHandler */ Rivet::AnalysisHandler * _rivet; /** * Event count */ unsigned long _nevent; }; } #include "ThePEG/Utilities/ClassTraits.h" namespace ThePEG { /** @cond TRAITSPECIALIZATIONS */ /** This template specialization informs ThePEG about the * base classes of NLORivetAnalysis. */ template <> struct BaseClassTrait { /** Typedef of the first base class of NLORivetAnalysis. */ typedef AnalysisHandler NthBase; }; /** This template specialization informs ThePEG about the name of * the NLORivetAnalysis class and the shared object where it is defined. */ template <> struct ClassTraits : public ClassTraitsBase { /** Return a platform-independent class name */ static string className() { return "ThePEG::NLORivetAnalysis"; } /** * The name of a file containing the dynamic library where the class * NLORivetAnalysis is implemented. It may also include several, space-separated, * libraries if the class NLORivetAnalysis depends on other classes (base classes * excepted). In this case the listed libraries will be dynamically * linked in the order they are specified. */ static string library() { return "RivetAnalysis.so"; } }; /** @endcond */ } #endif /* THEPEG_NLORivetAnalysis_H */ diff --git a/Analysis/RivetAnalysis.cc b/Analysis/RivetAnalysis.cc --- a/Analysis/RivetAnalysis.cc +++ b/Analysis/RivetAnalysis.cc @@ -1,203 +1,201 @@ // -*- C++ -*- // // This is the implementation of the non-inlined, non-templated member // functions of the RivetAnalysis class. // - +#include #include "ThePEG/Interface/Interfaced.h" #include "ThePEG/Interface/ClassDocumentation.h" #include "ThePEG/Interface/ParVector.h" #include "ThePEG/Interface/Parameter.h" #include "ThePEG/Interface/Switch.h" #include "ThePEG/Persistency/PersistentOStream.h" #include "ThePEG/Persistency/PersistentIStream.h" #include "ThePEG/Vectors/HepMCConverter.h" #include "ThePEG/Config/HepMCHelper.h" #include "ThePEG/Repository/EventGenerator.h" #include "ThePEG/Repository/CurrentGenerator.h" -#include "HepMC/GenEvent.h" #include "RivetAnalysis.h" #include "Rivet/AnalysisHandler.hh" #include "Rivet/Tools/RivetPaths.hh" #include "Rivet/Tools/Logging.hh" -#include using namespace ThePEG; RivetAnalysis::RivetAnalysis() : debug(false), _rivet(), _nevent(0) {} void RivetAnalysis::analyze(ThePEG::tEventPtr event, long ieve, int loop, int state) { ++_nevent; AnalysisHandler::analyze(event, ieve, loop, state); // Rotate to CMS, extract final state particles and call analyze(particles). // convert to hepmc HepMC::GenEvent * hepmc = ThePEG::HepMCConverter::convert(*event); // analyse the event if(_nevent>1) CurrentGenerator::Redirect stdout(cout); if ( _rivet ){ #if ThePEG_RIVET_VERSION == 1 _rivet->analyze(*hepmc); #elif ThePEG_RIVET_VERSION > 1 try { _rivet->analyze(*hepmc); } catch (const YODA::Exception & e) { Throw() << "Warning: Rivet/Yoda got the exception: "<< e.what()<<"\n" << Exception::warning; } #else #error "Unknown ThePEG_RIVET_VERSION" #endif } if(_nevent<=1) { // check that analysis list is still available if ( _rivet->analysisNames().size() != _analyses.size() ) { throw ThePEG::Exception() << "Rivet could not find all requested analyses.\n" << "Use 'rivet --list-analyses' to check availability.\n" << ThePEG::Exception::runerror; } } // delete hepmc event delete hepmc; } ThePEG::IBPtr RivetAnalysis::clone() const { return new_ptr(*this); } ThePEG::IBPtr RivetAnalysis::fullclone() const { return new_ptr(*this); } void RivetAnalysis::persistentOutput(ThePEG::PersistentOStream & os) const { os << _analyses << _paths << filename << debug; } void RivetAnalysis::persistentInput(ThePEG::PersistentIStream & is, int) { is >> _analyses >> _paths >> filename >> debug; } ThePEG::ClassDescription RivetAnalysis::initRivetAnalysis; // Definition of the static class description member. void RivetAnalysis::Init() { static ThePEG::ClassDocumentation documentation ("The RivetAnalysis class is a simple class to allow analyses" " from the Rivet library to be called from ThePEG"); static ThePEG::ParVector interfaceAnalyses ("Analyses", "The names of the Rivet analyses to use", &RivetAnalysis::_analyses, -1, "", "","" "", false, false, ThePEG::Interface::nolimits); static ThePEG::ParVector interfacePaths ("Paths", "The directory paths where Rivet should look for analyses.", &RivetAnalysis::_paths, -1, "", "","" "", false, false, ThePEG::Interface::nolimits); static Parameter interfaceFilename ("Filename", #if ThePEG_RIVET_VERSION == 1 "The name of the file where the AIDA histograms are put. If empty, " "the run name will be used instead. '.aida' will in any case be " "appended to the file name.", #elif ThePEG_RIVET_VERSION > 1 "The name of the file where the YODA histograms are put. If empty, " "the run name will be used instead. '.yoda' will in any case be " "appended to the file name.", #else #error "Unknown ThePEG_RIVET_VERSION" #endif &RivetAnalysis::filename, "", true, false); static Switch interfaceDebug ("Debug", "Enable debug information from Rivet", &RivetAnalysis::debug, false, true, false); static SwitchOption interfaceDebugNo (interfaceDebug, "No", "Disable debug information.", false); static SwitchOption interfaceDebugYes (interfaceDebug, "Yes", "Enable debug information from Rivet.", true); interfaceAnalyses.rank(10); } void RivetAnalysis::dofinish() { AnalysisHandler::dofinish(); if( _nevent > 0 && _rivet ) { CurrentGenerator::Redirect stdout(cout); #if ThePEG_RIVET_VERSION > 2 - _rivet->setCrossSection(generator()->integratedXSec()/picobarn, - generator()->integratedXSecErr()/picobarn); + _rivet->setCrossSection(make_pair(generator()->integratedXSec()/picobarn, + generator()->integratedXSecErr()/picobarn)); #else _rivet->setCrossSection(generator()->integratedXSec()/picobarn); #endif _rivet->finalize(); string fname = filename; #if ThePEG_RIVET_VERSION == 1 if ( fname.empty() ) fname = generator()->path() + "/" + generator()->runName() + ".aida"; #elif ThePEG_RIVET_VERSION > 1 if ( fname.empty() ) fname = generator()->path() + "/" + generator()->runName() + ".yoda"; #else #error "Unknown ThePEG_RIVET_VERSION" #endif _rivet->writeData(fname); } delete _rivet; _rivet = nullptr; } void RivetAnalysis::doinit() { AnalysisHandler::doinit(); if(_analyses.empty()) throw ThePEG::Exception() << "Must have at least one analysis loaded in " << "RivetAnalysis::doinitrun()" << ThePEG::Exception::runerror; // check that analysis list is available _rivet = new Rivet::AnalysisHandler; //(fname); for ( int i = 0, N = _paths.size(); i < N; ++i ) Rivet::addAnalysisLibPath(_paths[i]); _rivet->addAnalyses(_analyses); if ( _rivet->analysisNames().size() != _analyses.size() ) { throw ThePEG::Exception() << "Rivet could not find all requested analyses.\n" << "Use 'rivet --list-analyses' to check availability.\n" << ThePEG::Exception::runerror; } delete _rivet; _rivet = 0; } void RivetAnalysis::doinitrun() { AnalysisHandler::doinitrun(); // create Rivet analysis handler CurrentGenerator::Redirect stdout(cout); _rivet = new Rivet::AnalysisHandler; //(fname); for ( int i = 0, N = _paths.size(); i < N; ++i ) Rivet::addAnalysisLibPath(_paths[i]); _rivet->addAnalyses(_analyses); // check that analysis list is still available if ( _rivet->analysisNames().size() != _analyses.size() ) { throw ThePEG::Exception() << "Rivet could not find all requested analyses.\n" << "Use 'rivet --list-analyses' to check availability.\n" << ThePEG::Exception::runerror; } if ( debug ) Rivet::Log::setLevel("Rivet",Rivet::Log::DEBUG); } diff --git a/Config/HepMCHelper.h b/Config/HepMCHelper.h --- a/Config/HepMCHelper.h +++ b/Config/HepMCHelper.h @@ -1,174 +1,180 @@ // -*- C++ -*- // // HepMCHelper_HepMC.h is a part of ThePEG - A multi-purpose Monte Carlo event generator // Copyright (C) 2002-2017 The Herwig Collaboration // // ThePEG is licenced under version 3 of the GPL, see COPYING for details. // Please respect the MCnet academic guidelines, see GUIDELINES for details. // // // This is a helper header to implement HepMC conversions // #include "ThePEG/Vectors/HepMCTraits.h" #ifdef HAVE_HEPMC3 #include "HepMC3/GenEvent.h" #include "HepMC3/GenVertex.h" #include "HepMC3/GenParticle.h" #include "HepMC3/Version.h" #include "HepMC3/WriterAscii.h" #include "HepMC3/WriterHEPEVT.h" #include "HepMC3/WriterAsciiHepMC2.h" #ifdef HAVE_HEPMC3_WRITERROOT_H #include "HepMC3/WriterRoot.h" #endif #ifdef HAVE_HEPMC3_WRITERROOTTREE_H #include "HepMC3/WriterRootTree.h" #endif namespace HepMC3 { using PdfInfo=GenPdfInfo; using Polarization=std::pair; } namespace HepMC=HepMC3; #else #include "HepMC/GenEvent.h" #include "HepMC/GenVertex.h" #include "HepMC/GenParticle.h" #include "HepMC/Version.h" #include "HepMC/IO_BaseClass.h" #include "HepMC/IO_GenEvent.h" #include "HepMC/IO_AsciiParticles.h" +namespace HepMC { +#undef GenParticlePtr +#undef GenVertexPtr +typedef GenParticle * GenParticlePtr; +typedef GenVertex * GenVertexPtr; +} #endif namespace ThePEG { /** * Struct for HepMC conversion */ #ifndef HAVE_HEPMC3 template<> struct HepMCTraits : public HepMCTraitsBase { #else // This is version 3! template<> struct HepMCTraits : public HepMCTraitsBase { /** Create an event object with number \a evno and \a weight. */ static EventT * newEvent(long evno, double weight, const map& optionalWeights) { EventT * e = new EventT(HepMC::Units::GEV, HepMC::Units::MM); e->set_event_number(evno); e->set_event_number(evno); //std::vector wnames; std::vector wvalues; //wnames.push_back("Default"); wvalues.push_back(weight); for ( map::const_iterator w = optionalWeights.begin(); w != optionalWeights.end(); ++w ) { //wnames.push_back(w->first); wvalues.push_back(w->second); } //e->run_info()->set_weight_names(wnames); e->weights()=wvalues; return e; } /** Create a new vertex. */ static VertexPtrT newVertex() { return std::make_shared(VertexT()); } /** Set the \a scale, \f$\alpha_S\f$ (\a aS) and \f$\alpha_{EM}\f$ (\a aEM) for the event \a e. The scale will be scaled with \a unit before given to the GenEvent. */ static void setScaleAndAlphas(EventT & e, Energy2 scale, double aS, double aEM, Energy unit) { e.add_attribute("event_scale",std::make_shared(sqrt(scale)/unit)); e.add_attribute("mpi",std::make_shared(-1));//Please fix it later, once ThePEG authors respond e.add_attribute("signal_process_id",std::make_shared(0));//Please fix it later, once ThePEG authors respond e.add_attribute("alphaQCD",std::make_shared(aS)); e.add_attribute("alphaQED",std::make_shared(aEM)); } /** Set the colour line (with index \a indx) to \a coline for particle \a p. */ static void setColourLine(ParticleT & p, int indx, int coline) { p.add_attribute("flow"+std::to_string(indx),std::make_shared(coline)); } /** Add an incoming particle, \a p, to the vertex, \a v. */ static void addIncoming(VertexT & v, ParticlePtrT p) { v.add_particle_in(p); } /** Add an outgoing particle, \a p, to the vertex, \a v. */ static void addOutgoing(VertexT & v, ParticlePtrT p) { v.add_particle_out(p); } /** Set the primary vertex, \a v, for the event \a e. */ static void setSignalProcessVertex(EventT & e, VertexPtrT v) { e.add_vertex(v); e.add_attribute("signal_process_vertex", std::make_shared(v->id())); } /** Set a vertex, \a v, for the event \a e. */ static void addVertex(EventT & e, VertexPtrT v) { e.add_vertex(v); } /** Set the beam particles for the event.*/ static void setBeamParticles(EventT & e, ParticlePtrT p1, ParticlePtrT p2) { // e.set_beam_particles(p1,p2); p1->set_status(4); p2->set_status(4); e.set_beam_particles(p1, p2); } /** Create a new particle object with momentum \a p, PDG number \a id and status code \a status. The momentum will be scaled with \a unit which according to the HepMC documentation should be GeV. */ static ParticlePtrT newParticle(const Lorentz5Momentum & p, long id, int status, Energy unit) { // Note that according to the documentation the momentum is stored in a // HepLorentzVector in GeV (event though the CLHEP standard is MeV). HepMC::FourVector p_scalar(p.x()/unit, p.y()/unit, p.z()/unit, p.e()/unit); ParticlePtrT genp = std::make_shared(ParticleT(p_scalar, id, status)); genp->set_generated_mass(p.mass()/unit); return genp; } /** Set the polarization directions, \a the and \a phi, for particle \a p. */ static void setPolarization(ParticleT & genp, double the, double phi) { genp.add_attribute("theta",std::make_shared(the)); genp.add_attribute("phi",std::make_shared(phi)); } /** Set the position \a p for the vertex, \a v. The length will be scaled with \a unit which normally should be millimeters. */ static void setPosition(VertexT & v, const LorentzPoint & p, Length unit) { HepMC::FourVector v_scaled(p.x()/unit, p.y()/unit, p.z()/unit, p.e()/unit); v.set_position(v_scaled); } #endif }; }