---

edm::PythiaSource Class Reference

#include <GeneratorInterface/Pythia6Interface/interface/PythiaSource.h>

Inheritance diagram for edm::PythiaSource:

List of all members.

Public Member Functions
void	endRun (Run &r)
	PythiaSource (const ParameterSet &, const InputSourceDescription &)
	Constructor.
virtual	~PythiaSource ()
	Destructor.
Private Member Functions
bool	call_pygive (const std::string &iParm)
	Interface to the PYGIVE pythia routine, with add'l protections.
bool	call_slha_init ()
bool	call_slhagive (const std::string &iParm)
bool	call_txgive (const std::string &iParm)
bool	call_txgive_init ()
void	clear ()
virtual bool	produce (Event &e)
Private Attributes
double	comenergy
bool	doubleParticle
double	emax
double	emin
double	etamax
double	etamin
HepMC::GenEvent *	evt
double	extCrossSect
double	extFilterEff
bool	flatEnergy
PtYDistributor *	fPtYGenerator
CLHEP::HepRandomEngine &	fRandomEngine
CLHEP::RandFlat *	fRandomGenerator
bool	gluinoHadronsEnabled
bool	imposeProperTimes_
	Impose proper times for pions/kaons at generator level.
std::string	kinedata
unsigned int	maxEventsToPrint_
	Events to print if verbosity.
int	particleID
std::vector< int >	particleIDs
double	phimax
double	phimin
double	ptmax
double	ptmin
bool	pythiaHepMCVerbosity_
	HepMC verbosity flag.
unsigned int	pythiaPylistVerbosity_
	Pythia PYLIST Verbosity flag.
bool	stopHadronsEnabled
TauolaInterface	tauola_
bool	useExternalGenerators_
bool	useTauola_
bool	useTauolaPolarization_
double	ymax
double	ymin

---

Detailed Description

Definition at line 34 of file PythiaSource.h.

---

Constructor & Destructor Documentation

PythiaSource::PythiaSource	(	const ParameterSet &	pset,
		const InputSourceDescription &	desc
	)

Constructor.

Definition at line 140 of file PythiaSource.cc.

References call_pygive(), call_slha_init(), call_slhagive(), call_txgive(), call_txgive_init(), cards, comenergy, edm::errors::Configuration, GenMuonPlsPt100GeV_cfg::cout, edm::TauolaInterface::disablePolarizationEffects(), doubleParticle, emax, emin, edm::TauolaInterface::enablePolarizationEffects(), lat::endl(), etamax, etamin, fPtYGenerator, fRandomEngine, fRandomGenerator, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), gluinoHadronsEnabled, i, edm::TauolaInterface::initialize(), kinedata, maxEventsToPrint_, pars, particleID, particleIDs, phimax, phimin, ptmax, ptmin, PYGLRHAD, PYSTRHAD, pythiaHepMCVerbosity_, pythiaPylistVerbosity_, randomEngine, stopHadronsEnabled, tauola_, useExternalGenerators_, useTauola_, useTauolaPolarization_, ymax, and ymin.

                                                                 :
  GeneratedInputSource(pset, desc), evt(0), 
  pythiaPylistVerbosity_ (pset.getUntrackedParameter<int>("pythiaPylistVerbosity",0)),
  pythiaHepMCVerbosity_ (pset.getUntrackedParameter<bool>("pythiaHepMCVerbosity",false)),
  imposeProperTimes_ (pset.getUntrackedParameter<bool>("imposeProperTimes",false)),
  maxEventsToPrint_ (pset.getUntrackedParameter<int>("maxEventsToPrint",1)),
  extCrossSect(pset.getUntrackedParameter<double>("crossSection", -1.)),
  extFilterEff(pset.getUntrackedParameter<double>("filterEfficiency", -1.)),
  comenergy(pset.getUntrackedParameter<double>("comEnergy",14000.)),
  useExternalGenerators_(false),
  useTauola_(false),
  useTauolaPolarization_(false),
  stopHadronsEnabled(false), gluinoHadronsEnabled(false),
  fRandomEngine(getEngineReference())
{
  
  // PYLIST Verbosity Level
  // Valid PYLIST arguments are: 1, 2, 3, 5, 7, 11, 12, 13
  pythiaPylistVerbosity_ = pset.getUntrackedParameter<int>("pythiaPylistVerbosity",0);
  
  // HepMC event verbosity Level
  pythiaHepMCVerbosity_ = pset.getUntrackedParameter<bool>("pythiaHepMCVerbosity",false);

  //Max number of events printed on verbosity level 
  maxEventsToPrint_ = pset.getUntrackedParameter<int>("maxEventsToPrint",0);
  
  particleID = pset.getUntrackedParameter<int>("ParticleID", 0);

  particleIDs = 
    pset.getUntrackedParameter<std::vector<int> >("ParticleIDs", 
                                                 std::vector<int>(0));

// Initialize the random engine unconditionally!
  randomEngine = &fRandomEngine;
  fRandomGenerator = new CLHEP::RandFlat(fRandomEngine) ;

  if(particleID) {

    cout <<" Particle ID = " << particleID << endl; 

    doubleParticle = pset.getUntrackedParameter<bool>("DoubleParticle",true);
    cout <<" double back-to-back " << doubleParticle << endl; 

    kinedata = pset.getUntrackedParameter<string>("kinematicsFile","");

    ptmin = pset.getUntrackedParameter<double>("Ptmin",20.);
    ptmax = pset.getUntrackedParameter<double>("Ptmax",420.);
    cout <<" ptmin = " << ptmin <<" ptmax = " << ptmax << endl;
  
    emin = pset.getUntrackedParameter<double>("Emin",-1);
    emax = pset.getUntrackedParameter<double>("Emax",-1);
    if ( emin > 0 && emax > 0 ) {
      cout <<" emin = " << emin <<" emax = " << emax << endl;
    }
  
    if(kinedata.size() < 1){                                                                                 
       etamin = pset.getUntrackedParameter<double>("Etamin",0.);                                             
       etamax = pset.getUntrackedParameter<double>("Etamax",2.2);                                            
       cout <<" etamin = " << etamin <<" etamax = " << etamax << endl;
    }else{                                                                                                   
       ymin = pset.getUntrackedParameter<double>("ymin",0.);                                                 
       ymax = pset.getUntrackedParameter<double>("ymax",10.);                                                
       cout <<" ymin = " << ymin <<" ymax = " << ymax << endl;                                               
    }
      
    
    phimin = pset.getUntrackedParameter<double>("Phimin",0.);
    phimax = pset.getUntrackedParameter<double>("Phimax",360.);
    cout <<" phimin = " << phimin <<" phimax = " << phimax << endl;

    if(kinedata.size() > 0)
    {
       int ptbins = pset.getUntrackedParameter<int>("ptBinning",1000);
       int ybins = pset.getUntrackedParameter<int>("yBinning",50);
       fPtYGenerator = new PtYDistributor(kinedata, fRandomEngine, 
                                          ptmax, ptmin, ymax, ymin, ptbins, ybins);
    }
  } else if ( particleIDs.size() > 1 ) {

    // Here a simple jet with a number of particles of your choice
    // is to be generated. The particles are given an energy between 
    // 500 MeV and 1 GeV (flat/random), and are then boosted with 
    // a |beta|*E between Emin and Emax and with eta and phi in the
    // ranges given.

    // Boost absolute value
    emin = pset.getUntrackedParameter<double>("Emin",0.5);
    emax = pset.getUntrackedParameter<double>("Emax",2.0);
    
    // Boost absolute value
    ptmin = pset.getUntrackedParameter<double>("Pmin",20.);
    ptmax = pset.getUntrackedParameter<double>("Pmax",20.);

    // Boost pseudo-rapidity range
    etamin = pset.getUntrackedParameter<double>("Etamin",0.); 
    etamax = pset.getUntrackedParameter<double>("Etamax",2.2);
    
    // Boost phi range
    phimin = pset.getUntrackedParameter<double>("Phimin",-3.1415926536);
    phimax = pset.getUntrackedParameter<double>("Phimax",+3.1415926536);

  }

  // Set PYTHIA parameters in a single ParameterSet
  ParameterSet pythia_params = 
    pset.getParameter<ParameterSet>("PythiaParameters") ;
  
  // The parameter sets to be read (default, min bias, user ...) in the
  // proper order.
  vector<string> setNames = 
    pythia_params.getParameter<vector<string> >("parameterSets");
  
  // Loop over the sets
  for ( unsigned i=0; i<setNames.size(); ++i ) {
    
    string mySet = setNames[i];
    
    // Read the PYTHIA parameters for each set of parameters
    vector<string> pars = 
      pythia_params.getParameter<vector<string> >(mySet);
    
    if (mySet != "SLHAParameters" && mySet != "CSAParameters"){
    
    // Loop over all parameters and stop in case of mistake
    for( vector<string>::const_iterator  
           itPar = pars.begin(); itPar != pars.end(); ++itPar ) {
      static string sRandomValueSetting("MRPY(1)");
      if( 0 == itPar->compare(0,sRandomValueSetting.size(),sRandomValueSetting) ) {
        throw edm::Exception(edm::errors::Configuration,"PythiaError")
          <<" attempted to set random number using pythia command 'MRPY(1)' this is not allowed.\n  Please use the RandomNumberGeneratorService to set the random number seed.";
      }
      if( ! call_pygive(*itPar) ) {
        throw edm::Exception(edm::errors::Configuration,"PythiaError") 
          <<" pythia did not accept the following \""<<*itPar<<"\"";
      }
    }
    } else if(mySet == "CSAParameters"){   

   // Read CSA parameter
  
   pars = pythia_params.getParameter<vector<string> >("CSAParameters");

   call_txgive_init();
  
  
   // Loop over all parameters and stop in case of a mistake
    for (vector<string>::const_iterator 
            itPar = pars.begin(); itPar != pars.end(); ++itPar) {
      call_txgive(*itPar); 
     
         } 
     
   } else if(mySet == "SLHAParameters"){   

   // Read SLHA parameter
  
   pars = pythia_params.getParameter<vector<string> >("SLHAParameters");

   // Loop over all parameters and stop in case of a mistake
    for (vector<string>::const_iterator 
            itPar = pars.begin(); itPar != pars.end(); ++itPar) {
      call_slhagive(*itPar); 
     
         } 
 
    call_slha_init(); 
  
  }
  }

   stopHadronsEnabled = pset.getUntrackedParameter<bool>("stopHadrons",false);
   gluinoHadronsEnabled = pset.getUntrackedParameter<bool>("gluinoHadrons",false);

  //Init names and pdg code of r-hadrons
   if(stopHadronsEnabled)  PYSTRHAD();
   if(gluinoHadronsEnabled)  PYGLRHAD();

  //In the future, we will get the random number seed on each event and tell 
  // pythia to use that new seed
  // The random engine has already been initialized.  DO NOT do it again!
#ifdef NOTYET
  edm::Service<RandomNumberGenerator> rng;
  uint32_t seed = rng->mySeed();
  ostringstream sRandomSet;
  sRandomSet <<"MRPY(1)="<<seed;
  call_pygive(sRandomSet.str());
#endif
  
  if(particleID || particleIDs.size() > 1 ) 
    {
      call_pyinit( "NONE", "p", "p", comenergy );
    } else {
      call_pyinit( "CMS", "p", "p", comenergy );
    }

  // TAUOLA, etc.
  //
  useExternalGenerators_ = pset.getUntrackedParameter<bool>("UseExternalGenerators",false);
//  useTauola_ = pset.getUntrackedParameter<bool>("UseTauola", false);
//  useTauolaPolarization_ = pset.getUntrackedParameter<bool>("UseTauolaPolarization", false);
  
  if ( useExternalGenerators_ ) {
 // read External Generator parameters
    ParameterSet ext_gen_params =
       pset.getParameter<ParameterSet>("ExternalGenerators") ;
    vector<string> extGenNames =
       ext_gen_params.getParameter< vector<string> >("parameterSets");
    for (unsigned int ip=0; ip<extGenNames.size(); ++ip )
    {
      string curSet = extGenNames[ip];
      ParameterSet gen_par_set =
         ext_gen_params.getUntrackedParameter< ParameterSet >(curSet);
/*
     cout << "----------------------------------------------" << endl;
     cout << "Read External Generator parameter set "  << endl;
     cout << "----------------------------------------------" << endl;
*/
     if ( curSet == "Tauola" )
     {
        useTauola_ = true;
        if ( useTauola_ ) {
           cout << "--> use TAUOLA" << endl;
        } 
        useTauolaPolarization_ = gen_par_set.getParameter<bool>("UseTauolaPolarization");
        if ( useTauolaPolarization_ ) 
        {
           cout << "(Polarization effects enabled)" << endl;
           tauola_.enablePolarizationEffects();
        } 
        else 
        {
           cout << "(Polarization effects disabled)" << endl;
           tauola_.disablePolarizationEffects();
        }
        vector<string> cards = gen_par_set.getParameter< vector<string> >("InputCards");
        cout << "----------------------------------------------" << endl;
        cout << "Initializing Tauola" << endl;
        for( vector<string>::const_iterator
                itPar = cards.begin(); itPar != cards.end(); ++itPar )
        {
           call_txgive(*itPar);
        }
        tauola_.initialize();
        //call_pretauola(-1); // initialize TAUOLA package for tau decays
     }
    }
    // cout << "----------------------------------------------" << endl;
  }


  cout << endl; // Stetically add for the output
  //********                                      
  
  produces<HepMCProduct>();
  produces<GenInfoProduct, edm::InRun>();
}

PythiaSource::~PythiaSource

(

)

[virtual]

Destructor.

Definition at line 399 of file PythiaSource.cc.

References clear().

                           {
  clear(); 
}

---

Member Function Documentation

bool PythiaSource::call_pygive

(

const std::string &

iParm

)

[private]

Interface to the PYGIVE pythia routine, with add'l protections.

Definition at line 678 of file PythiaSource.cc.

References pydat1, and PYGIVE.

Referenced by produce(), and PythiaSource().

                                                 {

  int numWarn = pydat1.mstu[26]; //# warnings
  int numErr = pydat1.mstu[22];// # errors
  
//call the fortran routine pygive with a fortran string
  PYGIVE( iParm.c_str(), iParm.length() );  
  //  PYGIVE( iParm );  
//if an error or warning happens it is problem
  return pydat1.mstu[26] == numWarn && pydat1.mstu[22] == numErr;   
 
}

bool PythiaSource::call_slha_init

(

)

[private]

Definition at line 745 of file PythiaSource.cc.

References SLHA_INIT.

Referenced by PythiaSource().

                             {
  
   SLHA_INIT();
   return 1;  
}

bool PythiaSource::call_slhagive

(

const std::string &

iParm

)

[private]

Definition at line 707 of file PythiaSource.cc.

References GenMuonPlsPt100GeV_cfg::cout, end, lat::endl(), f1, file, edm::FileInPath::fullPath(), SLHAGIVE, and pyDBSRunClass::temp.

Referenced by PythiaSource().

                                                   {
        if( iParm.find( "SLHAFILE", 0 ) != string::npos ) {
                string::size_type start = iParm.find_first_of( "=" ) + 1;
                string::size_type end = iParm.length() - 1;
                string::size_type temp = iParm.find_first_of( "'", start );
                if( temp != string::npos ) {
                        start = temp + 1;
                        end = iParm.find_last_of( "'" ) - 1;
                } 
                start = iParm.find_first_not_of( " ", start );
                end = iParm.find_last_not_of( " ", end );               
                string shortfile = iParm.substr( start, end - start + 1 );
                string file;
                if( shortfile[0] == '/' ) {
                        cout << "SLHA file given with absolut path." << endl;
                        file = shortfile;
                } else {
                //      try {
                                FileInPath f1( shortfile );
                                file = f1.fullPath();
                //      } catch(...) {
                //              cout << "SLHA file not in path. Trying anyway." << endl;
                //              file = shortfile;
                //      }
                }
                file = "SLHAFILE = '" + file + "'";
                SLHAGIVE( file.c_str(), file.length() );
                cout << "     " <<  file.c_str() << endl;
                
        } else {
                SLHAGIVE( iParm.c_str(), iParm.length() );
                cout << "     " <<  iParm.c_str() << endl; 
        }
        return 1;
}

bool PythiaSource::call_txgive

(

const std::string &

iParm

)

[private]

Definition at line 692 of file PythiaSource.cc.

References GenMuonPlsPt100GeV_cfg::cout, lat::endl(), and TXGIVE.

Referenced by PythiaSource().

                                                 {
  
   TXGIVE( iParm.c_str(), iParm.length() );
   cout << "     " <<  iParm.c_str() << endl; 
   return 1;  
}

bool PythiaSource::call_txgive_init

(

)

[private]

Definition at line 700 of file PythiaSource.cc.

References TXGIVE_INIT.

Referenced by PythiaSource().

                               {
  
   TXGIVE_INIT();
   return 1;  
}

void PythiaSource::clear

(

void

)

[private]

Definition at line 403 of file PythiaSource.cc.

Referenced by ~PythiaSource().

                         {
 
}

void PythiaSource::endRun

(

Run &

r

)

[virtual]

Reimplemented from edm::ConfigurableInputSource.

Definition at line 407 of file PythiaSource.cc.

References extCrossSect, extFilterEff, edm::TauolaInterface::print(), edm::Run::put(), pypars, tauola_, and useTauola_.

                                 {
 
 double cs = pypars.pari[0]; // cross section in mb
 auto_ptr<GenInfoProduct> giprod (new GenInfoProduct());
 giprod->set_cross_section(cs);
 giprod->set_external_cross_section(extCrossSect);
 giprod->set_filter_efficiency(extFilterEff);
 r.put(giprod);

  call_pystat(1);
  if ( useTauola_ ) {
    tauola_.print();
    //call_pretauola(1); // print TAUOLA decay statistics output
  }

}

bool PythiaSource::produce

(

Event &

e

)

[private, virtual]

Implements edm::ConfigurableInputSource.

Definition at line 424 of file PythiaSource.cc.

References funct::abs(), call_pygive(), call_pylist(), funct::cos(), GenMuonPlsPt100GeV_cfg::cout, doubleParticle, emax, emin, lat::endl(), eta, etamax, etamin, edm::ConfigurableInputSource::event(), evt, funct::exp(), edm::PtYDistributor::firePt(), edm::PtYDistributor::fireY(), edm::first(), fPtYGenerator, gluinoHadronsEnabled, id, imposeProperTimes_, kinedata, prof2calltree::last, funct::log(), maxEventsToPrint_, edm::ConfigurableInputSource::numberEventsInRun(), P, particleID, particleIDs, phi, phimax, phimin, pi, edm::TauolaInterface::processEvent(), ptmax, ptmin, edm::Event::put(), PY1ENT, PYCOMP, pydat1, pydat2, pydat3, PYEXEC, PYGLFR, pyint1, PYMASS, pypars, pyr_(), PYROBO, PYSTFR, pythiaHepMCVerbosity_, pythiaPylistVerbosity_, edm::InputSource::remainingEvents(), funct::sin(), funct::sqrt(), stopHadronsEnabled, t, funct::tan(), tauola_, useTauola_, vbegin, vend, x, y, and z.

                                    {

    auto_ptr<HepMCProduct> bare_product(new HepMCProduct());  

    //********                                         
    //  
   if(particleID) 
      {    
         double pi = 3.1415927;
         int ip = 1;  
         int dum;
         double ee=0,the=0,eta=0;
         double pmass = PYMASS(particleID);
         double phi = (phimax-phimin)*pyr_(&dum)+phimin; 
         
         if(kinedata.size() < 1){  // no kinematics input specified, use flat distribution, pt and eta         
            double pt  = (ptmax-ptmin)*pyr_(&dum)+ptmin;                                                 
            double e   = (emax-emin)*pyr_(&dum)+emin;
            eta = (etamax-etamin)*pyr_(&dum)+etamin;                                                      
            the = 2.*atan(exp(-eta));                                                                          
            if ( emin > pmass && emax > pmass ) { // generate single particle distribution flat in energy      
               ee = e;                                                                                         
            } else { // generate single particle distribution flat in pt                                       
               double pe = pt/sin(the);                                                                        
               ee = sqrt(pe*pe+pmass*pmass);                                                                   
            }                                                                                                  
         }else{ // kinematics from input file, pt and y                                                        
            double pt  = fPtYGenerator->firePt();                                                
            double y = fPtYGenerator->fireY();                                                       
            double u = exp(y);                                                                                 
            ee = 0.5*sqrt(pmass*pmass+pt*pt)*(u*u+1)/u;                                                        
            double pz = sqrt(ee*ee-pt*pt-pmass*pmass);                   
            if(y<0) pz = -pz;
            the = atan(pt/pz);                                                                                 
            if(pz < 0) the = pi + the;
            eta = -log(tan(the/2));                                                                            
         }                                 
         /*
           cout <<" pt = " << pt
           <<" eta = " << eta                                                                                  
           <<" the = " << the                                                                                  
           <<" pe = " << pe                                                                                    
           <<" phi = " << phi                                                                                  
           <<" pmass = " << pmass                                                                              
           <<" ee = " << ee << endl; 
        */

        phi = phi * (3.1415927/180.);

        PY1ENT(ip, particleID, ee, the, phi);

        if(doubleParticle)
          {
            ip = ip + 1;
// Check if particle is its own anti-particle.
            // int particleID2 = -1 * particleID;
            int pythiaCode = PYCOMP(particleID);
            int has_antipart = pydat2.kchg[3-1][pythiaCode-1];
            int particleID2 = has_antipart ? -1 * particleID : particleID;          
            the = 2.*atan(exp(eta));
            phi  = phi + 3.1415927;
            if (phi > 2.* 3.1415927) {phi = phi - 2.* 3.1415927;}         
            PY1ENT(ip, particleID2, ee, the, phi);
          }
        PYEXEC();

      } else if ( particleIDs.size() > 1 ) { 

        // Initialize some variables
        int dum;
        int ip = 0;
        double px = 0;  
        double py = 0;
        double pz = 0;
        double pe = 0;
        // The particles in the "c.o.m. frame"
        for ( unsigned id=0; id<particleIDs.size(); ++id ) {
          ++ip;
          int pythiaCode = particleIDs[id];
          // Generate the particles with emin -> emax Fermi motion
          double e = emin + (emax-emin) * pyr_(&pythiaCode);
          // Any direction
          double phi = 2. * 3.1415927 * pyr_(&pythiaCode);
          double ctt = -1. + 2.*pyr_(&pythiaCode);
          double the = std::acos(ctt);
          // Add the entry in PYJETS
          PY1ENT(ip,pythiaCode,e,the,phi);
          // To compute the overall mass
          px += P(ip,1);
          py += P(ip,2);
          pz += P(ip,3);
          pe += P(ip,4);
        }

        // The boost
        // The overall mass
        double mm = std::sqrt(pe*pe-px*px-py*py-pz*pz);
        // The boost absolute value
        double pp = ptmin + (ptmax-ptmin)*pyr_(&dum);
        double ee = sqrt(pp*pp+mm*mm);
        // The boost direction
        double fhi = phimin + (phimax-phimin)*pyr_(&dum);
        double eta = etamin + (etamax-etamin)*pyr_(&dum);
        double tet = 2.*atan(exp(-eta));
        // betax, betay, betaz
        double betax = pp/ee * std::sin(tet) * std::cos(fhi);
        double betay = pp/ee * std::sin(tet) * std::sin(fhi);
        double betaz = pp/ee * std::cos(tet);
        // No rotation
        double rothe = 0.;
        double rophi = 0.;
        // Boost all particles
        int first = -1;
        int last = -1;
        PYROBO(first, last, rothe, rophi, betax, betay, betaz);
        // And now decay boosted particles
        PYEXEC();

      } else {
          if(!gluinoHadronsEnabled && !stopHadronsEnabled)
          {
             call_pyevnt();      // generate one event with Pythia
          }
          else
          {
             call_pygive("MSTJ(14)=-1");
             call_pyevnt();      // generate one event with Pythia
             call_pygive("MSTJ(14)=1");
             if(gluinoHadronsEnabled)  PYGLFR();
             if(stopHadronsEnabled)  PYSTFR();
          }
          
      }

    if ( useTauola_ ) {
      tauola_.processEvent();
      //call_pretauola(0); // generate tau decays with TAUOLA
    }

    // convert to stdhep format
    //
    call_pyhepc( 1 );
    
    // convert stdhep (hepevt) to hepmc
    //
    //HepMC::GenEvent* evt = conv.getGenEventfromHEPEVT();
    HepMC::GenEvent* evt = conv.read_next_event();
    
    // fix for 1-part events
    if ( particleID ) evt->set_beam_particles(0,0);
    
    evt->set_signal_process_id(pypars.msti[0]);
    evt->set_event_scale(pypars.pari[16]);
    evt->set_event_number(numberEventsInRun() - remainingEvents() - 1);

    // int id1 = pypars.msti[14];
    // int id2 = pypars.msti[15];
    int id1 = pyint1.mint[14];
    int id2 = pyint1.mint[15];
    if ( id1 == 21 ) id1 = 0;
    if ( id2 == 21 ) id2 = 0; 
    double x1 = pyint1.vint[40];
    double x2 = pyint1.vint[41];  
    double Q  = pyint1.vint[50];
    double pdf1 = pyint1.vint[38];
    pdf1 /= x1 ;
    double pdf2 = pyint1.vint[39];
    pdf2 /= x2 ;
    evt->set_pdf_info( HepMC::PdfInfo(id1,id2,x1,x2,Q,pdf1,pdf2) ) ;
    
    evt->weights().push_back( pyint1.vint[96] );

    if (imposeProperTimes_ || pydat1.mstj[21]==3 || pydat1.mstj[21]==4 ) {
      int dumm;
      HepMC::GenEvent::vertex_const_iterator vbegin = evt->vertices_begin();
      HepMC::GenEvent::vertex_const_iterator vend = evt->vertices_end();
      HepMC::GenEvent::vertex_const_iterator vitr = vbegin;
      for (; vitr != vend; ++vitr ) {
            HepMC::GenVertex::particle_iterator pbegin = (*vitr)->particles_begin(HepMC::children);
            HepMC::GenVertex::particle_iterator pend = (*vitr)->particles_end(HepMC::children);
            HepMC::GenVertex::particle_iterator pitr = pbegin;
            for (; pitr != pend; ++pitr) {
                  if ((*pitr)->end_vertex()) continue;
                  if ((*pitr)->status()!=1) continue;
                  int pdgcode= abs((*pitr)->pdg_id());
                  // Do nothing if the particle is not expected to decay
                  if (pydat3.mdcy[0][PYCOMP(pdgcode)-1]!=1) continue;

                  double ctau = pydat2.pmas[3][PYCOMP(pdgcode)-1];
                  HepMC::FourVector mom = (*pitr)->momentum();
                  HepMC::FourVector vin = (*vitr)->position();
                  double x = 0.;
                  double y = 0.;
                  double z = 0.;
                  double t = 0.;
                  bool decayInRange = false;
                  while (!decayInRange) {
                        double unif_rand = pyr_(&dumm);
                        // Value of 0 is excluded, so following line is OK
                        double proper_length = - ctau * log(unif_rand);
                        double factor = proper_length/mom.m();
                        x = vin.x() + factor * mom.px();
                        y = vin.y() + factor * mom.py();
                        z = vin.z() + factor * mom.pz();
                        t = vin.t() + factor * mom.e();
                        // Decay must be happen outside a cylindrical region
                        if (pydat1.mstj[21]==4) {
                              if (sqrt(x*x+y*y)>pydat1.parj[72] || abs(z)>pydat1.parj[73]) decayInRange = true;
                        // Decay must be happen outside a given sphere
                        } else if (pydat1.mstj[21]==3) {
                              if (sqrt(x*x+y*y+z*z)>pydat1.parj[71]) decayInRange = true;
                        // Decay is always OK otherwise
                        } else {
                              decayInRange = true;
                        }
                  }
                  
                  HepMC::GenVertex* vdec = new HepMC::GenVertex(HepMC::FourVector(x,y,z,t));
                  evt->add_vertex(vdec);
                  vdec->add_particle_in((*pitr));
            }
      }
    }
    
    //******** Verbosity ********
    
    if(event() <= maxEventsToPrint_ &&
       (pythiaPylistVerbosity_ || pythiaHepMCVerbosity_)) {

      // Prints PYLIST info
      if(pythiaPylistVerbosity_) {
        call_pylist(pythiaPylistVerbosity_);
      }
      
      // Prints HepMC event
      if(pythiaHepMCVerbosity_) {
        cout << "Event process = " << pypars.msti[0] << endl 
        << "----------------------" << endl;
        evt->print();
      }
    }
    
    
    //evt = reader_->fillCurrentEventData(); 
    //********                                      

    if(evt)  bare_product->addHepMCData(evt );

    e.put(bare_product);

    return true;
}

---

Member Data Documentation

double edm::PythiaSource::comenergy [private]

Definition at line 82 of file PythiaSource.h.

Referenced by PythiaSource().

bool edm::PythiaSource::doubleParticle [private]

Definition at line 77 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

double edm::PythiaSource::emax [private]

Definition at line 83 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

double edm::PythiaSource::emin [private]

Definition at line 83 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

double edm::PythiaSource::etamax [private]

Definition at line 80 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

double edm::PythiaSource::etamin [private]

Definition at line 80 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

HepMC::GenEvent* edm::PythiaSource::evt [private]

Definition at line 59 of file PythiaSource.h.

Referenced by produce().

double edm::PythiaSource::extCrossSect [private]

Definition at line 71 of file PythiaSource.h.

Referenced by endRun().

double edm::PythiaSource::extFilterEff [private]

Definition at line 72 of file PythiaSource.h.

Referenced by endRun().

bool edm::PythiaSource::flatEnergy [private]

Definition at line 85 of file PythiaSource.h.

PtYDistributor* edm::PythiaSource::fPtYGenerator [private]

Definition at line 98 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

CLHEP::HepRandomEngine& edm::PythiaSource::fRandomEngine [private]

Definition at line 96 of file PythiaSource.h.

Referenced by PythiaSource().

CLHEP::RandFlat* edm::PythiaSource::fRandomGenerator [private]

Definition at line 97 of file PythiaSource.h.

Referenced by PythiaSource().

bool edm::PythiaSource::gluinoHadronsEnabled [private]

Definition at line 94 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

bool edm::PythiaSource::imposeProperTimes_ [private]

Impose proper times for pions/kaons at generator level.

Definition at line 66 of file PythiaSource.h.

Referenced by produce().

std::string edm::PythiaSource::kinedata [private]

Definition at line 78 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

unsigned int edm::PythiaSource::maxEventsToPrint_ [private]

Events to print if verbosity.

Definition at line 68 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

int edm::PythiaSource::particleID [private]

Definition at line 75 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

std::vector<int> edm::PythiaSource::particleIDs [private]

Definition at line 76 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

double edm::PythiaSource::phimax [private]

Definition at line 81 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

double edm::PythiaSource::phimin [private]

Definition at line 81 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

double edm::PythiaSource::ptmax [private]

Definition at line 79 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

double edm::PythiaSource::ptmin [private]

Definition at line 79 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

bool edm::PythiaSource::pythiaHepMCVerbosity_ [private]

HepMC verbosity flag.

Definition at line 64 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

unsigned int edm::PythiaSource::pythiaPylistVerbosity_ [private]

Pythia PYLIST Verbosity flag.

Definition at line 62 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

bool edm::PythiaSource::stopHadronsEnabled [private]

Definition at line 93 of file PythiaSource.h.

Referenced by produce(), and PythiaSource().

TauolaInterface edm::PythiaSource::tauola_ [private]

Definition at line 91 of file PythiaSource.h.

Referenced by endRun(), produce(), and PythiaSource().

bool edm::PythiaSource::useExternalGenerators_ [private]

Definition at line 88 of file PythiaSource.h.

Referenced by PythiaSource().

bool edm::PythiaSource::useTauola_ [private]

Definition at line 89 of file PythiaSource.h.

Referenced by endRun(), produce(), and PythiaSource().

bool edm::PythiaSource::useTauolaPolarization_ [private]

Definition at line 90 of file PythiaSource.h.

Referenced by PythiaSource().

double edm::PythiaSource::ymax [private]

Definition at line 84 of file PythiaSource.h.

Referenced by PythiaSource().

double edm::PythiaSource::ymin [private]

Definition at line 84 of file PythiaSource.h.

Referenced by PythiaSource().

---

The documentation for this class was generated from the following files:

• GeneratorInterface/Pythia6Interface/interface/PythiaSource.h
• GeneratorInterface/Pythia6Interface/src/PythiaSource.cc

---