Herwig6Hadronizer Class Reference

Inheritance diagram for Herwig6Hadronizer:

Public Member Functions
const char *	classname () const
bool	decay ()
bool	declareSpecialSettings (const std::vector< std::string > &)
bool	declareStableParticles (const std::vector< int > &pdgIds)
void	finalizeEvent ()
bool	generatePartonsAndHadronize ()
bool	hadronize ()
	Herwig6Hadronizer (const edm::ParameterSet &params)
bool	initialize (const lhef::HEPRUP *heprup)
bool	initializeForExternalPartons ()
bool	initializeForInternalPartons ()
bool	readSettings (int)
bool	residualDecay ()
void	setSLHAFromHeader (const std::vector< std::string > &lines)
void	statistics ()
	~Herwig6Hadronizer () override
Public Member Functions inherited from gen::BaseHadronizer
	BaseHadronizer (edm::ParameterSet const &ps)
void	cleanLHE ()
void	generateLHE (edm::LuminosityBlock const &lumi, CLHEP::HepRandomEngine *rengine, unsigned int ncpu)
edm::Event &	getEDMEvent () const
std::unique_ptr< HepMC::GenEvent >	getGenEvent ()
std::unique_ptr < GenEventInfoProduct >	getGenEventInfo ()
virtual std::unique_ptr < GenLumiInfoHeader >	getGenLumiInfoHeader () const
GenRunInfoProduct &	getGenRunInfo ()
std::unique_ptr< lhef::LHEEvent >	getLHEEvent ()
const std::shared_ptr < lhef::LHERunInfo > &	getLHERunInfo () const
const std::string &	gridpackPath () const
int	randomIndex () const
const std::string &	randomInitConfigDescription () const
void	randomizeIndex (edm::LuminosityBlock const &lumi, CLHEP::HepRandomEngine *rengine)
void	resetEvent (std::unique_ptr< HepMC::GenEvent > event)
void	resetEventInfo (std::unique_ptr< GenEventInfoProduct > eventInfo)
virtual bool	select (HepMC::GenEvent *) const
void	setEDMEvent (edm::Event &event)
void	setLHEEvent (std::unique_ptr< lhef::LHEEvent > event)
void	setLHERunInfo (std::unique_ptr< lhef::LHERunInfo > runInfo)
void	setRandomEngine (CLHEP::HepRandomEngine *v)
std::vector< std::string > const &	sharedResources () const
virtual	~BaseHadronizer () noexcept(false)
Public Member Functions inherited from gen::Herwig6Instance
bool	callWithTimeout (unsigned int secs, void(*fn)())
bool	give (const std::string &line)
	Herwig6Instance ()
void	openParticleSpecFile (const std::string fileName)
void	setHerwigRandomEngine (CLHEP::HepRandomEngine *v)
	~Herwig6Instance () override
Public Member Functions inherited from gen::FortranInstance
void	call (void(&fn)())
template<typename T >
T	call (T(&fn)())
template<typename A >
void	call (void(&fn)(A), A a)
template<typename T , typename A >
T	call (T(&fn)(A), A a)
template<typename A1 , typename A2 >
void	call (void(&fn)(A1, A2), A1 a1, A2 a2)
template<typename T , typename A1 , typename A2 >
T	call (T(&fn)(A1, A2), A1 a1, A2 a2)
virtual void	enter ()
	FortranInstance ()
virtual void	leave ()
virtual bool	upVeto ()
virtual	~FortranInstance () noexcept(false)

Private Member Functions
void	clear ()
void	doSetRandomEngine (CLHEP::HepRandomEngine *v) override
std::vector< std::string > const &	doSharedResources () const override
int	pythiaStatusCode (const HepMC::GenParticle *p) const
void	pythiaStatusCodes ()
void	upEvnt () override
void	upInit () override

Private Attributes
double	comEnergy
HepMC::IO_HERWIG	conv
bool	doMatching
bool	doMPInteraction
bool	emulatePythiaStatusCodes
bool	externalPartons
bool	fConvertToPDG
int	hepmcVerbosity
int	herwigVerbosity
bool	inclusiveMatching
double	matchingScale
int	maxEventsToPrint
bool	needClear
int	nMatch
int	numTrials
gen::ParameterCollector	parameters
std::string	particleSpecFileName
bool	printCards
bool	readMCatNLOfile
bool	readParticleSpecFile
bool	useJimmy

Static Private Attributes
static const std::vector < std::string >	theSharedResources

Additional Inherited Members
Static Public Member Functions inherited from gen::FortranInstance
template<typename T >
static T *	getInstance ()
Static Public Attributes inherited from gen::FortranInstance
static const std::string	kFortranInstance = "FortranInstance"
Protected Member Functions inherited from gen::BaseHadronizer
std::unique_ptr < HepMC::GenEvent > &	event ()
std::unique_ptr < GenEventInfoProduct > &	eventInfo ()
lhef::LHEEvent *	lheEvent ()
lhef::LHERunInfo *	lheRunInfo ()
GenRunInfoProduct &	runInfo ()
Protected Member Functions inherited from gen::Herwig6Instance
virtual bool	hwwarn (const std::string &fn, int code)
Protected Attributes inherited from gen::BaseHadronizer
std::string	lheFile_
int	randomIndex_

Detailed Description

Definition at line 83 of file Herwig6Hadronizer.cc.

Constructor & Destructor Documentation

Herwig6Hadronizer::Herwig6Hadronizer

(

const edm::ParameterSet &

params

)

Definition at line 163 of file Herwig6Hadronizer.cc.

References edm::ParameterSet::exists(), fConvertToPDG, and edm::ParameterSet::getParameter().

    : BaseHadronizer(params),
      needClear(false),
      parameters(params.getParameter<edm::ParameterSet>("HerwigParameters")),
      herwigVerbosity(params.getUntrackedParameter<int>("herwigVerbosity", 0)),
      hepmcVerbosity(params.getUntrackedParameter<int>("hepmcVerbosity", 0)),
      maxEventsToPrint(params.getUntrackedParameter<int>("maxEventsToPrint", 0)),
      printCards(params.getUntrackedParameter<bool>("printCards", false)),
      emulatePythiaStatusCodes(params.getUntrackedParameter<bool>("emulatePythiaStatusCodes", false)),
      comEnergy(params.getParameter<double>("comEnergy")),
      useJimmy(params.getParameter<bool>("useJimmy")),
      doMPInteraction(params.getParameter<bool>("doMPInteraction")),
      numTrials(params.getUntrackedParameter<int>("numTrialsMPI", 100)),
      doMatching(params.getUntrackedParameter<bool>("doMatching", false)),
      inclusiveMatching(params.getUntrackedParameter<bool>("inclusiveMatching", true)),
      nMatch(params.getUntrackedParameter<int>("nMatch", 0)),
      matchingScale(params.getUntrackedParameter<double>("matchingScale", 0.0)),
      readMCatNLOfile(false),

      // added to be able to read external particle spectrum file
      particleSpecFileName(params.getUntrackedParameter<std::string>("ParticleSpectrumFileName", "")),
      readParticleSpecFile(params.getUntrackedParameter<bool>("readParticleSpecFile", false))

{
  fConvertToPDG = false;
  if (params.exists("doPDGConvert"))
    fConvertToPDG = params.getParameter<bool>("doPDGConvert");
}

Herwig6Hadronizer::~Herwig6Hadronizer ( )

override

Definition at line 192 of file Herwig6Hadronizer.cc.

References clear().

{ clear(); }

Member Function Documentation

const char* Herwig6Hadronizer::classname ( ) const

inline

Definition at line 110 of file Herwig6Hadronizer.cc.

{ return "Herwig6Hadronizer"; }

void Herwig6Hadronizer::clear ( void  )

private

Definition at line 196 of file Herwig6Hadronizer.cc.

References gen::FortranInstance::call(), jmefin, needClear, and useJimmy.

Referenced by initialize(), readSettings(), and ~Herwig6Hadronizer().

                              {
  if (!needClear)
    return;

  // teminate elementary process
  call(hwefin);
  if (useJimmy)
    call(jmefin);

  needClear = false;
}

bool Herwig6Hadronizer::decay

(

)

Definition at line 824 of file Herwig6Hadronizer.cc.

                              {
  // hadron decays

  //    InstanceWrapper wrapper(this);  // safe guard
  //
  //    //int iproc = hwproc.IPROC;
  //    //hwproc.IPROC = 312;
  //    hwdhad();   // unstable particle decays
  //    //hwproc.IPROC = iproc;
  //
  //    hwdhvy();   // heavy flavour decays
  //    hwmevt();   // soft underlying event
  //
  //    if (doMatching) {
  //      bool pass = call_hwmatch();
  //      if (!pass) {
  //        printf("Event failed MLM matching\n");
  //        hwufne();
  //        if (lheEvent()) lheEvent()->count(lhef::LHERunInfo::kKilled);
  //        return false;
  //      }
  //    }
  //
  //    hwufne();   // finalize event
  //
  //    if (hwevnt.IERROR)
  //        return false;
  //
  //    if (lheEvent()) lheEvent()->count(lhef::LHERunInfo::kAccepted);
  //
  //    event().reset(new HepMC::GenEvent);
  //    if (!conv.fill_next_event(event().get()))
  //        throw cms::Exception("Herwig6Error")
  //            << "HepMC Conversion problems in event." << std::endl;
  //
  //         // do particle ID conversion Herwig->PDG, if requested
  //         if ( fConvertToPDG ) {
  //            for ( HepMC::GenEvent::particle_iterator part = event()->particles_begin(); part != event()->particles_end(); ++part) {
  //              if ((*part)->pdg_id() != HepPID::translateHerwigtoPDT((*part)->pdg_id()))
  //                (*part)->set_pdg_id(HepPID::translateHerwigtoPDT((*part)->pdg_id()));
  //            }
  //         }

  return true;
}

bool Herwig6Hadronizer::declareSpecialSettings

(

const std::vector< std::string > &

)

Definition at line 612 of file Herwig6Hadronizer.cc.

{ return true; }

bool Herwig6Hadronizer::declareStableParticles

(

const std::vector< int > &

pdgIds

)

Definition at line 584 of file Herwig6Hadronizer.cc.

                                                                           {
  markStable(13);     // MU+
  markStable(-13);    // MU-
  markStable(3112);   // SIGMA+
  markStable(-3112);  // SIGMABAR+
  markStable(3222);   // SIGMA-
  markStable(-3222);  // SIGMABAR-
  markStable(3122);   // LAMBDA0
  markStable(-3122);  // LAMBDABAR0
  markStable(3312);   // XI-
  markStable(-3312);  // XIBAR+
  markStable(3322);   // XI0
  markStable(-3322);  // XI0BAR
  markStable(3334);   // OMEGA-
  markStable(-3334);  // OMEGABAR+
  markStable(211);    // PI+
  markStable(-211);   // PI-
  markStable(321);    // K+
  markStable(-321);   // K-
  markStable(310);    // K_S0
  markStable(130);    // K_L0

  for (std::vector<int>::const_iterator iter = pdgIds.begin(); iter != pdgIds.end(); ++iter)
    if (!markStable(*iter))
      return false;
  return true;
}

void Herwig6Hadronizer::doSetRandomEngine ( CLHEP::HepRandomEngine *  v )

overrideprivatevirtual

Reimplemented from gen::BaseHadronizer.

Definition at line 194 of file Herwig6Hadronizer.cc.

References gen::Herwig6Instance::setHerwigRandomEngine().

{ setHerwigRandomEngine(v); }

std::vector<std::string> const& Herwig6Hadronizer::doSharedResources ( ) const

inlineoverrideprivatevirtual

Reimplemented from gen::BaseHadronizer.

Definition at line 114 of file Herwig6Hadronizer.cc.

References theSharedResources.

{ return theSharedResources; }

void Herwig6Hadronizer::finalizeEvent

(

)

Definition at line 710 of file Herwig6Hadronizer.cc.

References funct::abs(), emulatePythiaStatusCodes, gen::BaseHadronizer::event(), gen::BaseHadronizer::eventInfo(), externalPartons, lhef::LHEEvent::fillEventInfo(), lhef::LHEEvent::fillPdfInfo(), lhef::LHEEvent::fixHepMCEventTimeOrdering(), GenParticle::GenParticle, hwuaem_(), hwualf_(), gen::BaseHadronizer::lheEvent(), gen::BaseHadronizer::lheRunInfo(), SiPixelPI::one, and pythiaStatusCodes().

                                      {
  lhef::LHEEvent::fixHepMCEventTimeOrdering(event().get());

  HepMC::PdfInfo pdfInfo;
  if (externalPartons) {
    lheEvent()->fillEventInfo(event().get());
    lheEvent()->fillPdfInfo(&pdfInfo);

    // for MC@NLO: IDWRUP is not filled...
    if (event()->signal_process_id() == 0)
      event()->set_signal_process_id(abs(hwproc.IPROC));
  }

  HepMC::GenParticle *incomingParton = nullptr;
  HepMC::GenParticle *targetParton = nullptr;

  HepMC::GenParticle *incomingProton = nullptr;
  HepMC::GenParticle *targetProton = nullptr;

  // find incoming parton (first entry with IST=121)
  for (HepMC::GenEvent::particle_const_iterator it = event()->particles_begin();
       (it != event()->particles_end() && incomingParton == nullptr);
       it++)
    if ((*it)->status() == 121)
      incomingParton = (*it);

  // find target parton (first entry with IST=122)
  for (HepMC::GenEvent::particle_const_iterator it = event()->particles_begin();
       (it != event()->particles_end() && targetParton == nullptr);
       it++)
    if ((*it)->status() == 122)
      targetParton = (*it);

  // find incoming Proton (first entry ID=2212, IST=101)
  for (HepMC::GenEvent::particle_const_iterator it = event()->particles_begin();
       (it != event()->particles_end() && incomingProton == nullptr);
       it++)
    if ((*it)->status() == 101 && (*it)->pdg_id() == 2212)
      incomingProton = (*it);

  // find target Proton (first entry ID=2212, IST=102)
  for (HepMC::GenEvent::particle_const_iterator it = event()->particles_begin();
       (it != event()->particles_end() && targetProton == nullptr);
       it++)
    if ((*it)->status() == 102 && (*it)->pdg_id() == 2212)
      targetProton = (*it);

  // find hard scale Q (computed from colliding partons)
  if (incomingParton && targetParton) {
    math::XYZTLorentzVector totMomentum(0, 0, 0, 0);
    totMomentum += incomingParton->momentum();
    totMomentum += targetParton->momentum();
    double evScale = totMomentum.mass();
    double evScale2 = evScale * evScale;

    // find alpha_QED & alpha_QCD
    int one = 1;
    double alphaQCD = hwualf_(&one, &evScale);
    double alphaQED = hwuaem_(&evScale2);

    if (!externalPartons || event()->event_scale() < 0)
      event()->set_event_scale(evScale);
    if (!externalPartons || event()->alphaQCD() < 0)
      event()->set_alphaQCD(alphaQCD);
    if (!externalPartons || event()->alphaQED() < 0)
      event()->set_alphaQED(alphaQED);

    if (!externalPartons || pdfInfo.x1() < 0) {
      // get the PDF information
      pdfInfo.set_id1(incomingParton->pdg_id() == 21 ? 0 : incomingParton->pdg_id());
      pdfInfo.set_id2(targetParton->pdg_id() == 21 ? 0 : targetParton->pdg_id());
      if (incomingProton && targetProton) {
        double x1 = incomingParton->momentum().pz() / incomingProton->momentum().pz();
        double x2 = targetParton->momentum().pz() / targetProton->momentum().pz();
        pdfInfo.set_x1(x1);
        pdfInfo.set_x2(x2);
      }
      // we do not fill pdf1 & pdf2, since they are not easily available (what are they needed for anyways???)
      pdfInfo.set_scalePDF(evScale);  // the same as Q above... does this make sense?
    }

    if (!externalPartons || event()->signal_process_id() < 0)
      event()->set_signal_process_id(abs(hwproc.IPROC));
    event()->set_pdf_info(pdfInfo);
  }

  // add event weight & PDF information
  if (lheRunInfo() != nullptr && std::abs(lheRunInfo()->getHEPRUP()->IDWTUP) == 4)
    // in LHE weighting mode 4 the weight is an xsec, so convert form HERWIG
    // to standard CMS unit "picobarn"
    event()->weights().push_back(1.0e3 * hwevnt.EVWGT);
  else
    event()->weights().push_back(hwevnt.EVWGT);

  // find final parton (first entry with IST=123)
  HepMC::GenParticle *finalParton = nullptr;
  for (HepMC::GenEvent::particle_const_iterator it = event()->particles_begin();
       (it != event()->particles_end() && finalParton == nullptr);
       it++)
    if ((*it)->status() == 123)
      finalParton = (*it);

  // add GenEventInfo & binning Values
  eventInfo() = std::make_unique<GenEventInfoProduct>(event().get());
  if (finalParton) {
    double thisPt = finalParton->momentum().perp();
    eventInfo()->setBinningValues(std::vector<double>(1, thisPt));
  }

  // emulate PY6 status codes, if switched on...
  if (emulatePythiaStatusCodes)
    pythiaStatusCodes();
}

bool Herwig6Hadronizer::generatePartonsAndHadronize ( )

inline

Definition at line 104 of file Herwig6Hadronizer.cc.

References hadronize().

{ return hadronize(); }

bool Herwig6Hadronizer::hadronize

(

)

Definition at line 629 of file Herwig6Hadronizer.cc.

References gen::Herwig6Instance::callWithTimeout(), conv, lhef::LHEEvent::count(), doMatching, doMPInteraction, relativeConstraints::error, gen::BaseHadronizer::event(), fConvertToPDG, hwwarn_(), lhef::LHERunInfo::kAccepted, lhef::LHERunInfo::kKilled, lhef::LHERunInfo::kSelected, gen::BaseHadronizer::lheEvent(), gpuVertexFinder::printf(), useJimmy, and wrapper.

Referenced by generatePartonsAndHadronize().

                                  {
  // hard process generation, parton shower, hadron formation

  InstanceWrapper wrapper(this);  // safe guard

  event().reset();

  // call herwig routines to create HEPEVT

  hwuine();  // initialize event

  if (callWithTimeout(10, hwepro)) {  // process event and PS
    // We hung for more than 10 seconds
    int error = 199;
    hwwarn_("HWHGUP", &error);
  }

  hwbgen();  // parton cascades

  // call jimmy ... only if event is not killed yet by HERWIG
  if (useJimmy && doMPInteraction && !hwevnt.IERROR && call_hwmsct()) {
    if (lheEvent())
      lheEvent()->count(lhef::LHERunInfo::kKilled);
    return false;
  }

  hwdhob();  // heavy quark decays
  hwcfor();  // cluster formation
  hwcdec();  // cluster decays

  //    // if event *not* killed by HERWIG, return true
  //    if (hwevnt.IERROR) {
  //      hwufne(); // finalize event, to keep system clean
  //      if (lheEvent()) lheEvent()->count(lhef::LHERunInfo::kKilled);
  //      return false;
  //    }

  //if (lheEvent()) lheEvent()->count(lhef::LHERunInfo::kAccepted);

  hwdhad();  // unstable particle decays
  hwdhvy();  // heavy flavour decays
  hwmevt();  // soft underlying event

  hwufne();  // finalize event

  // if event *not* killed by HERWIG, return true
  if (hwevnt.IERROR) {
    if (lheEvent())
      lheEvent()->count(lhef::LHERunInfo::kKilled);
    return false;
  }

  if (doMatching) {
    bool pass = call_hwmatch();
    if (!pass) {
      printf("Event failed MLM matching\n");
      if (lheEvent())
        lheEvent()->count(lhef::LHERunInfo::kSelected);
      return false;
    }
  }

  if (lheEvent())
    lheEvent()->count(lhef::LHERunInfo::kAccepted);

  event() = std::make_unique<HepMC::GenEvent>();
  if (!conv.fill_next_event(event().get()))
    throw cms::Exception("Herwig6Error") << "HepMC Conversion problems in event." << std::endl;

  // do particle ID conversion Herwig->PDG, if requested
  if (fConvertToPDG) {
    for (HepMC::GenEvent::particle_iterator part = event()->particles_begin(); part != event()->particles_end();
         ++part) {
      if ((*part)->pdg_id() != HepPID::translateHerwigtoPDT((*part)->pdg_id()))
        (*part)->set_pdg_id(HepPID::translateHerwigtoPDT((*part)->pdg_id()));
    }
  }

  return true;
}

bool Herwig6Hadronizer::initialize

(

const lhef::HEPRUP *

heprup

)

Definition at line 413 of file Herwig6Hadronizer.cc.

References gen::ParameterCollector::begin(), gen::FortranInstance::call(), clear(), comEnergy, edm::errors::Configuration, doMPInteraction, alignCSCRings::e, lhef::HEPRUP::EBMUP, gen::ParameterCollector::end(), Exception, externalPartons, lhef::LHERunInfo::findHeader(), gen::Herwig6Instance::give(), hepmcVerbosity, herwigVerbosity, hwdspn, hwprch, mps_fire::i, lhef::HEPRUP::IDBMUP, info(), jimmin, jminit, jmparm, gen::BaseHadronizer::lheRunInfo(), geometryCSVtoXML::line, maxEventsToPrint, needClear, gen::Herwig6Instance::openParticleSpecFile(), parameters, particleSpecFileName, lhef::LHERunInfo::pdfSetTranslation(), printCards, gpuVertexFinder::printf(), readParticleSpecFile, setSLHAFromHeader(), and useJimmy.

                                                           {
  clear();

  externalPartons = (heprup != nullptr);

  std::ostringstream info;
  info << "---------------------------------------------------\n";
  info << "Initializing Herwig6Hadronizer for " << (externalPartons ? "external" : "internal") << " partons\n";
  info << "---------------------------------------------------\n";

  info << "   Herwig verbosity level         = " << herwigVerbosity << "\n";
  info << "   HepMC verbosity                = " << hepmcVerbosity << "\n";
  info << "   Number of events to be printed = " << maxEventsToPrint << "\n";

  // Call hwudat to set up HERWIG block data
  hwudat();

  // Setting basic parameters
  if (externalPartons) {
    hwproc.PBEAM1 = heprup->EBMUP.first;
    hwproc.PBEAM2 = heprup->EBMUP.second;
    pdgToHerwig(heprup->IDBMUP.first, hwbmch.PART1);
    pdgToHerwig(heprup->IDBMUP.second, hwbmch.PART2);
  } else {
    hwproc.PBEAM1 = 0.5 * comEnergy;
    hwproc.PBEAM2 = 0.5 * comEnergy;
    pdgToHerwig(2212, hwbmch.PART1);
    pdgToHerwig(2212, hwbmch.PART2);
  }

  if (useJimmy) {
    info << "   HERWIG will be using JIMMY for UE/MI.\n";
    jmparm.MSFLAG = 1;
    if (doMPInteraction)
      info << "   JIMMY trying to generate multiple interactions.\n";
  }

  // set the IPROC already here... needed for VB pairs

  bool iprocFound = false;

  for (gen::ParameterCollector::const_iterator line = parameters.begin(); line != parameters.end(); ++line) {
    if (!strcmp((line->substr(0, 5)).c_str(), "IPROC")) {
      if (!give(*line))
        throw edm::Exception(edm::errors::Configuration)
            << "Herwig 6 did not accept the following: \"" << *line << "\"." << std::endl;
      else
        iprocFound = true;
    }
  }

  if (!iprocFound && !externalPartons)
    throw edm::Exception(edm::errors::Configuration) << "You have to define the process with IPROC." << std::endl;

  // initialize other common blocks ...
  call(hwigin);
  hwevnt.MAXER = 100000000;  // O(inf)
  hwpram.LWSUD = 0;          // don't write Sudakov form factors
  hwdspn.LWDEC = 0;          // don't write three/four body decays
                             // (no fort.77 and fort.88 ...)

  // init LHAPDF glue

  std::memset(hwprch.AUTPDF, ' ', sizeof hwprch.AUTPDF);
  for (unsigned int i = 0; i < 2; i++) {
    hwpram.MODPDF[i] = -111;
    std::memcpy(hwprch.AUTPDF[i], "HWLHAPDF", 8);
  }

  if (useJimmy)
    call(jimmin);

  hwevnt.MAXPR = maxEventsToPrint;
  hwpram.IPRINT = herwigVerbosity;
  //    hwprop.RMASS[6] = 175.0;    //FIXME

  if (printCards) {
    info << "\n";
    info << "------------------------------------\n";
    info << "Reading HERWIG parameters\n";
    info << "------------------------------------\n";
  }
  for (gen::ParameterCollector::const_iterator line = parameters.begin(); line != parameters.end(); ++line) {
    if (printCards)
      info << "   " << *line << "\n";
    if (!give(*line))
      throw edm::Exception(edm::errors::Configuration)
          << "Herwig 6 did not accept the following: \"" << *line << "\"." << std::endl;
  }

  if (printCards)
    info << "\n";

  if (externalPartons) {
    std::vector<std::string> slha = lheRunInfo()->findHeader("slha");
    if (!slha.empty())
      setSLHAFromHeader(slha);
  }

  needClear = true;

  std::pair<int, int> pdfs(-1, -1);
  if (externalPartons)
    pdfs = lheRunInfo()->pdfSetTranslation();

  if (hwpram.MODPDF[0] != -111 || hwpram.MODPDF[1] != -111) {
    for (unsigned int i = 0; i < 2; i++)
      if (hwpram.MODPDF[i] == -111)
        hwpram.MODPDF[i] = -1;

    if (pdfs.first != -1 || pdfs.second != -1)
      edm::LogError("Generator|Herwig6Hadronzier") << "Both external Les Houches event and "
                                                      "config file specify a PDF set.  "
                                                      "User PDF will override external one."
                                                   << std::endl;

    pdfs.first = hwpram.MODPDF[0] != -111 ? hwpram.MODPDF[0] : -1;
    pdfs.second = hwpram.MODPDF[1] != -111 ? hwpram.MODPDF[1] : -1;
  }

  printf("pdfs.first = %i, pdfs.second = %i\n", pdfs.first, pdfs.second);

  hwpram.MODPDF[0] = pdfs.first;
  hwpram.MODPDF[1] = pdfs.second;

  if (externalPartons)
    hwproc.IPROC = -1;

  //Lars: lower EFFMIN threshold, to continue execution of IPROC=4000, lambda'_211=0.01 at LM7,10
  if (readParticleSpecFile) {
    openParticleSpecFile(particleSpecFileName);
    hwpram.EFFMIN = 1e-5;
  }

  // HERWIG preparations ...
  call(hwuinc);
  markStable(13);     // MU+
  markStable(-13);    // MU-
  markStable(3112);   // SIGMA+
  markStable(-3112);  // SIGMABAR+
  markStable(3222);   // SIGMA-
  markStable(-3222);  // SIGMABAR-
  markStable(3122);   // LAMBDA0
  markStable(-3122);  // LAMBDABAR0
  markStable(3312);   // XI-
  markStable(-3312);  // XIBAR+
  markStable(3322);   // XI0
  markStable(-3322);  // XI0BAR
  markStable(3334);   // OMEGA-
  markStable(-3334);  // OMEGABAR+
  markStable(211);    // PI+
  markStable(-211);   // PI-
  markStable(321);    // K+
  markStable(-321);   // K-
  markStable(310);    // K_S0
  markStable(130);    // K_L0

  // better: merge with declareStableParticles
  // and get the list from configuration / Geant4 / Core somewhere

  // initialize HERWIG event generation
  call(hweini);

  if (useJimmy)
    call(jminit);

  edm::LogInfo(info.str());

  return true;
}

bool Herwig6Hadronizer::initializeForExternalPartons ( )

inline

Definition at line 97 of file Herwig6Hadronizer.cc.

References initializeForInternalPartons().

{ return initializeForInternalPartons(); }

bool Herwig6Hadronizer::initializeForInternalPartons

(

)

Definition at line 397 of file Herwig6Hadronizer.cc.

References gen::FortranInstance::call(), jimmin, jminit, and useJimmy.

Referenced by initializeForExternalPartons().

                                                     {
  if (useJimmy)
    call(jimmin);

  call(hwuinc);

  // initialize HERWIG event generation
  call(hweini);

  if (useJimmy) {
    call(jminit);
  }

  return true;
}

int Herwig6Hadronizer::pythiaStatusCode ( const HepMC::GenParticle *  p ) const

private

Definition at line 911 of file Herwig6Hadronizer.cc.

References externalPartons, GenParticle::GenParticle, submitPVResolutionJobs::q, and mps_update::status.

Referenced by pythiaStatusCodes().

                                                                       {
  int status = p->status();

  // weird 9922212 particles...
  if (status == 3 && !p->end_vertex())
    return 2;

  if (status >= 1 && status <= 3)
    return status;

  if (!p->end_vertex())
    return 1;

  // let's prevent particles having status 3, if the identical
  // particle downstream is a better status 3 candidate
  int currentId = p->pdg_id();
  int orig = status;
  if (status == 123 || status == 124 || status == 155 || status == 156 || status == 160 ||
      (status >= 195 && status <= 197)) {
    for (const HepMC::GenParticle *q = p;;) {
      const HepMC::GenVertex *vtx = q->end_vertex();
      if (!vtx)
        break;

      HepMC::GenVertex::particles_out_const_iterator iter;
      for (iter = vtx->particles_out_const_begin(); iter != vtx->particles_out_const_end(); ++iter)
        if ((*iter)->pdg_id() == currentId)
          break;

      if (iter == vtx->particles_out_const_end())
        break;

      q = *iter;
      if (q->status() == 3 || ((status == 120 || status == 123 || status == 124) && orig > 124))
        return 4;
    }
  }

  int nesting = 0;
  for (;;) {
    if ((status >= 120 && status <= 122) || status == 3) {
      // avoid flagging status 3 if there is a
      // better status 3 candidate upstream
      if (externalPartons)
        return ((orig >= 121 && orig <= 124) || orig == 3) ? 3 : 4;
      else
        return (nesting || (status != 3 && orig <= 124)) ? 3 : 4;
    }

    // check whether we are leaving the hard process
    // including heavy resonance decays
    if (!(status == 4 || status == 123 || status == 124 || status == 155 || status == 156 || status == 160 ||
          (status >= 195 && status <= 197)))
      break;

    const HepMC::GenVertex *vtx = p->production_vertex();
    if (!vtx || !vtx->particles_in_size())
      break;

    p = *vtx->particles_in_const_begin();
    status = p->status();

    int newId = p->pdg_id();

    if (!newId)
      break;

    // nesting increases if we move to the next-best mother
    if (newId != currentId) {
      if (++nesting > 1 && externalPartons)
        break;
      currentId = newId;
    }
  }

  return 2;
}

void Herwig6Hadronizer::pythiaStatusCodes ( )

private

Definition at line 989 of file Herwig6Hadronizer.cc.

References gen::BaseHadronizer::event(), and pythiaStatusCode().

Referenced by finalizeEvent().

                                          {
  for (HepMC::GenEvent::particle_iterator iter = event()->particles_begin(); iter != event()->particles_end(); iter++)
    (*iter)->set_status(pythiaStatusCode(*iter));

  for (HepMC::GenEvent::particle_iterator iter = event()->particles_begin(); iter != event()->particles_end(); iter++)
    if ((*iter)->status() == 4)
      (*iter)->set_status(2);
}

bool Herwig6Hadronizer::readSettings

(

int

key

)

Definition at line 246 of file Herwig6Hadronizer.cc.

References gen::ParameterCollector::begin(), gen::FortranInstance::call(), clear(), comEnergy, edm::errors::Configuration, doMatching, doMPInteraction, alignCSCRings::e, lhef::HEPRUP::EBMUP, gen::ParameterCollector::end(), Exception, externalPartons, lhef::LHERunInfo::findHeader(), lhef::LHERunInfo::getHEPRUP(), gen::Herwig6Instance::give(), hepmcVerbosity, herwigVerbosity, hwdspn, hwmatchpram, hwprch, mps_fire::i, lhef::HEPRUP::IDBMUP, inclusiveMatching, info(), jmparm, gen::BaseHadronizer::lheRunInfo(), geometryCSVtoXML::line, matchingScale, maxEventsToPrint, needClear, nMatch, gen::Herwig6Instance::openParticleSpecFile(), parameters, particleSpecFileName, lhef::LHERunInfo::pdfSetTranslation(), printCards, gpuVertexFinder::printf(), readParticleSpecFile, setSLHAFromHeader(), and useJimmy.

                                            {
  clear();
  const lhef::HEPRUP *heprup = lheRunInfo()->getHEPRUP();
  externalPartons = (heprup != nullptr);

  if (key == 0 && externalPartons)
    return false;
  if (key > 0 && !externalPartons)
    return false;

  std::ostringstream info;
  info << "---------------------------------------------------\n";
  info << "Taking in settinsg for Herwig6Hadronizer for " << (externalPartons ? "external" : "internal")
       << " partons\n";
  info << "---------------------------------------------------\n";

  info << "   Herwig verbosity level         = " << herwigVerbosity << "\n";
  info << "   HepMC verbosity                = " << hepmcVerbosity << "\n";
  info << "   Number of events to be printed = " << maxEventsToPrint << "\n";

  // Call hwudat to set up HERWIG block data
  hwudat();

  // Setting basic parameters
  if (externalPartons) {
    hwproc.PBEAM1 = heprup->EBMUP.first;
    hwproc.PBEAM2 = heprup->EBMUP.second;
    pdgToHerwig(heprup->IDBMUP.first, hwbmch.PART1);
    pdgToHerwig(heprup->IDBMUP.second, hwbmch.PART2);
  } else {
    hwproc.PBEAM1 = 0.5 * comEnergy;
    hwproc.PBEAM2 = 0.5 * comEnergy;
    pdgToHerwig(2212, hwbmch.PART1);
    pdgToHerwig(2212, hwbmch.PART2);
  }

  if (doMatching) {
    hwmatchpram.n_match = nMatch;
    hwmatchpram.matching_scale = matchingScale;

    if (inclusiveMatching)
      hwmatchpram.max_multiplicity_flag = 1;
    else
      hwmatchpram.max_multiplicity_flag = 0;
  }

  if (useJimmy) {
    info << "   HERWIG will be using JIMMY for UE/MI.\n";
    jmparm.MSFLAG = 1;
    if (doMPInteraction)
      info << "   JIMMY trying to generate multiple interactions.\n";
  }

  // set the IPROC already here... needed for VB pairs

  bool iprocFound = false;

  for (gen::ParameterCollector::const_iterator line = parameters.begin(); line != parameters.end(); ++line) {
    if (!strcmp((line->substr(0, 5)).c_str(), "IPROC")) {
      if (!give(*line))
        throw edm::Exception(edm::errors::Configuration)
            << "Herwig 6 did not accept the following: \"" << *line << "\"." << std::endl;
      else
        iprocFound = true;
    }
  }

  if (!iprocFound && !externalPartons)
    throw edm::Exception(edm::errors::Configuration) << "You have to define the process with IPROC." << std::endl;

  // initialize other common blocks ...
  call(hwigin);  // default init

  hwevnt.MAXER = 100000000;  // O(inf)
  hwpram.LWSUD = 0;          // don't write Sudakov form factors
  hwdspn.LWDEC = 0;          // don't write three/four body decays
                             // (no fort.77 and fort.88 ...)
  // init LHAPDF glue
  std::memset(hwprch.AUTPDF, ' ', sizeof hwprch.AUTPDF);
  for (unsigned int i = 0; i < 2; i++) {
    hwpram.MODPDF[i] = -111;
    std::memcpy(hwprch.AUTPDF[i], "HWLHAPDF", 8);
  }

  hwevnt.MAXPR = maxEventsToPrint;
  hwpram.IPRINT = herwigVerbosity;
  //    hwprop.RMASS[6] = 175.0;    //FIXME

  if (printCards) {
    info << "\n";
    info << "------------------------------------\n";
    info << "Reading HERWIG parameters\n";
    info << "------------------------------------\n";
  }
  for (gen::ParameterCollector::const_iterator line = parameters.begin(); line != parameters.end(); ++line) {
    if (printCards)
      info << "   " << *line << "\n";
    if (!give(*line))
      throw edm::Exception(edm::errors::Configuration)
          << "Herwig 6 did not accept the following: \"" << *line << "\"." << std::endl;
  }

  if (printCards)
    info << "\n";

  if (externalPartons) {
    std::vector<std::string> slha = lheRunInfo()->findHeader("slha");
    if (!slha.empty())
      setSLHAFromHeader(slha);
  }

  needClear = true;

  std::pair<int, int> pdfs(-1, -1);
  if (externalPartons)
    pdfs = lheRunInfo()->pdfSetTranslation();

  if (hwpram.MODPDF[0] != -111 || hwpram.MODPDF[1] != -111) {
    for (unsigned int i = 0; i < 2; i++)
      if (hwpram.MODPDF[i] == -111)
        hwpram.MODPDF[i] = -1;

    if (pdfs.first != -1 || pdfs.second != -1)
      edm::LogError("Generator|Herwig6Hadronzier") << "Both external Les Houches event and "
                                                      "config file specify a PDF set.  "
                                                      "User PDF will override external one."
                                                   << std::endl;

    pdfs.first = hwpram.MODPDF[0] != -111 ? hwpram.MODPDF[0] : -1;
    pdfs.second = hwpram.MODPDF[1] != -111 ? hwpram.MODPDF[1] : -1;
  }

  printf("pdfs.first = %i, pdfs.second = %i\n", pdfs.first, pdfs.second);

  hwpram.MODPDF[0] = pdfs.first;
  hwpram.MODPDF[1] = pdfs.second;

  if (externalPartons && hwproc.IPROC >= 0)
    hwproc.IPROC = -1;

  //Lars: lower EFFMIN threshold, to continue execution of IPROC=4000, lambda'_211=0.01 at LM7,10
  if (readParticleSpecFile) {
    openParticleSpecFile(particleSpecFileName);
    hwpram.EFFMIN = 1e-5;
  }

  edm::LogInfo(info.str());

  return true;
}

bool Herwig6Hadronizer::residualDecay

(

)

Definition at line 870 of file Herwig6Hadronizer.cc.

{ return true; }

void Herwig6Hadronizer::setSLHAFromHeader

(

const std::vector< std::string > &

lines

)

Definition at line 208 of file Herwig6Hadronizer.cc.

References createPayload::block, gather_cfg::blocks, geometryCSVtoXML::line, submitPVValidationJobs::split(), AlCaHLTBitMon_QueryRunRegistry::string, and HcalDetIdTransform::transform().

Referenced by initialize(), and readSettings().

                                                                           {
  std::set<std::string> blocks;
  std::string block;
  for (std::vector<std::string>::const_iterator iter = lines.begin(); iter != lines.end(); ++iter) {
    std::string line = *iter;
    std::transform(line.begin(), line.end(), line.begin(), (int (*)(int))std::toupper);
    std::string::size_type pos = line.find('#');
    if (pos != std::string::npos)
      line.resize(pos);

    if (line.empty())
      continue;

    if (!boost::algorithm::is_space()(line[0])) {
      std::vector<std::string> tokens;
      boost::split(tokens, line, boost::algorithm::is_space(), boost::token_compress_on);
      if (tokens.empty())
        continue;
      block.clear();
      if (tokens.size() < 2)
        continue;
      if (tokens[0] == "BLOCK")
        block = tokens[1];
      else if (tokens[0] == "DECAY")
        block = "DECAY";

      if (block.empty())
        continue;

      if (!blocks.count(block)) {
        blocks.insert(block);
        edm::LogWarning("Generator|Herwig6Hadronzier")
            << "Unsupported SLHA block \"" << block << "\".  It will be ignored." << std::endl;
      }
    }
  }
}

void Herwig6Hadronizer::statistics

(

)

Definition at line 614 of file Herwig6Hadronizer.cc.

References gen::BaseHadronizer::runInfo(), and GenRunInfoProduct::setInternalXSec().

                                   {
  if (!runInfo().internalXSec()) {
    // not set via LHE, so get it from HERWIG
    // the reason is that HERWIG doesn't compute the xsec
    // in all LHE modes

    double RNWGT = 1. / hwevnt.NWGTS;
    double AVWGT = hwevnt.WGTSUM * RNWGT;

    double xsec = 1.0e3 * AVWGT;

    runInfo().setInternalXSec(xsec);
  }
}

void Herwig6Hadronizer::upEvnt ( )

overrideprivatevirtual

Reimplemented from gen::FortranInstance.

Definition at line 895 of file Herwig6Hadronizer.cc.

References SplitLinear::begin, edm::errors::Configuration, lhef::CommonBlocks::fillHEPEUP(), lhef::LHEEvent::getComments(), gen::Herwig6Instance::give(), gen::BaseHadronizer::lheEvent(), readMCatNLOfile, and AlCaHLTBitMon_QueryRunRegistry::string.

                               {
  lhef::CommonBlocks::fillHEPEUP(lheEvent()->getHEPEUP());

  // if MCatNLO external file is read, read comment & pass IHPRO to HERWIG
  if (readMCatNLOfile) {
    for (std::vector<std::string>::const_iterator iter = lheEvent()->getComments().begin();
         iter != lheEvent()->getComments().end();
         ++iter) {
      std::string toParse(iter->substr(1));
      if (!give(toParse))
        throw edm::Exception(edm::errors::Configuration)
            << "Herwig 6 did not accept the following: \"" << toParse << "\"." << std::endl;
    }
  }
}

void Herwig6Hadronizer::upInit ( )

overrideprivatevirtual

Reimplemented from gen::FortranInstance.

Definition at line 872 of file Herwig6Hadronizer.cc.

References edm::errors::Configuration, lhef::CommonBlocks::fillHEPRUP(), lhef::LHERunInfo::getHeaders(), gen::Herwig6Instance::give(), getRunAppsInfo::headers, heprup_, gen::BaseHadronizer::lheRunInfo(), HEPRUP_::pdfgup, HEPRUP_::pdfsup, readMCatNLOfile, and AlCaHLTBitMon_QueryRunRegistry::string.

                               {
  lhef::CommonBlocks::fillHEPRUP(lheRunInfo()->getHEPRUP());
  heprup_.pdfgup[0] = heprup_.pdfgup[1] = -1;
  heprup_.pdfsup[0] = heprup_.pdfsup[1] = -1;
  // we set up the PDFs ourselves

  // pass HERWIG paramaters fomr header (if present)
  std::string mcnloHeader = "herwig6header";
  std::vector<lhef::LHERunInfo::Header> headers = lheRunInfo()->getHeaders();
  for (std::vector<lhef::LHERunInfo::Header>::const_iterator hIter = headers.begin(); hIter != headers.end(); ++hIter) {
    if (hIter->tag() == mcnloHeader) {
      readMCatNLOfile = true;
      for (lhef::LHERunInfo::Header::const_iterator lIter = hIter->begin(); lIter != hIter->end(); ++lIter) {
        if ((lIter->c_str())[0] != '#' && (lIter->c_str())[0] != '\n') {  // it's not a comment)
          if (!give(*lIter))
            throw edm::Exception(edm::errors::Configuration)
                << "Herwig 6 did not accept the following: \"" << *lIter << "\"." << std::endl;
        }
      }
    }
  }
}

Member Data Documentation

double Herwig6Hadronizer::comEnergy

private

Definition at line 136 of file Herwig6Hadronizer.cc.

Referenced by initialize(), and readSettings().

HepMC::IO_HERWIG Herwig6Hadronizer::conv

private

Definition at line 126 of file Herwig6Hadronizer.cc.

Referenced by hadronize().

bool Herwig6Hadronizer::doMatching

private

Definition at line 141 of file Herwig6Hadronizer.cc.

Referenced by hadronize(), and readSettings().

bool Herwig6Hadronizer::doMPInteraction

private

Definition at line 138 of file Herwig6Hadronizer.cc.

Referenced by hadronize(), initialize(), and readSettings().

bool Herwig6Hadronizer::emulatePythiaStatusCodes

private

Definition at line 135 of file Herwig6Hadronizer.cc.

Referenced by finalizeEvent().

bool Herwig6Hadronizer::externalPartons

private

Definition at line 128 of file Herwig6Hadronizer.cc.

Referenced by finalizeEvent(), initialize(), pythiaStatusCode(), and readSettings().

bool Herwig6Hadronizer::fConvertToPDG

private

Definition at line 140 of file Herwig6Hadronizer.cc.

Referenced by hadronize(), and Herwig6Hadronizer().

int Herwig6Hadronizer::hepmcVerbosity

private

Definition at line 132 of file Herwig6Hadronizer.cc.

Referenced by initialize(), and readSettings().

int Herwig6Hadronizer::herwigVerbosity

private

Definition at line 131 of file Herwig6Hadronizer.cc.

Referenced by initialize(), and readSettings().

bool Herwig6Hadronizer::inclusiveMatching

private

Definition at line 142 of file Herwig6Hadronizer.cc.

Referenced by readSettings().

double Herwig6Hadronizer::matchingScale

private

Definition at line 144 of file Herwig6Hadronizer.cc.

Referenced by readSettings().

int Herwig6Hadronizer::maxEventsToPrint

private

Definition at line 133 of file Herwig6Hadronizer.cc.

Referenced by initialize(), and readSettings().

bool Herwig6Hadronizer::needClear

private

Definition at line 127 of file Herwig6Hadronizer.cc.

Referenced by clear(), initialize(), and readSettings().

int Herwig6Hadronizer::nMatch

private

Definition at line 143 of file Herwig6Hadronizer.cc.

Referenced by readSettings().

int Herwig6Hadronizer::numTrials

private

Definition at line 139 of file Herwig6Hadronizer.cc.

gen::ParameterCollector Herwig6Hadronizer::parameters

private

Definition at line 130 of file Herwig6Hadronizer.cc.

Referenced by initialize(), and readSettings().

std::string Herwig6Hadronizer::particleSpecFileName

private

Definition at line 149 of file Herwig6Hadronizer.cc.

Referenced by initialize(), and readSettings().

bool Herwig6Hadronizer::printCards

private

Definition at line 134 of file Herwig6Hadronizer.cc.

Referenced by initialize(), and readSettings().

bool Herwig6Hadronizer::readMCatNLOfile

private

Definition at line 146 of file Herwig6Hadronizer.cc.

Referenced by upEvnt(), and upInit().

bool Herwig6Hadronizer::readParticleSpecFile

private

Definition at line 150 of file Herwig6Hadronizer.cc.

Referenced by initialize(), and readSettings().

const std::vector< std::string > Herwig6Hadronizer::theSharedResources

staticprivate

Initial value:

= {edm::SharedResourceNames::kHerwig6,
                                                                        gen::FortranInstance::kFortranInstance}

Definition at line 124 of file Herwig6Hadronizer.cc.

Referenced by doSharedResources().

bool Herwig6Hadronizer::useJimmy

private

Definition at line 137 of file Herwig6Hadronizer.cc.

Referenced by clear(), hadronize(), initialize(), initializeForInternalPartons(), and readSettings().