File indexing completed on 2024-04-06 12:13:54
0001
0002 #include <iostream>
0003
0004 #include "Pythia6PtGun.h"
0005
0006 #include "FWCore/Utilities/interface/Exception.h"
0007
0008 #include "FWCore/Framework/interface/EventSetup.h"
0009
0010 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0011
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013
0014 using namespace edm;
0015 using namespace gen;
0016
0017 Pythia6PtGun::Pythia6PtGun(const ParameterSet& pset) : Pythia6ParticleGun(pset) {
0018
0019 ParameterSet pgun_params = pset.getParameter<ParameterSet>("PGunParameters");
0020 fMinEta = pgun_params.getParameter<double>("MinEta");
0021 fMaxEta = pgun_params.getParameter<double>("MaxEta");
0022 fMinPt = pgun_params.getParameter<double>("MinPt");
0023 fMaxPt = pgun_params.getParameter<double>("MaxPt");
0024 fAddAntiParticle = pgun_params.getParameter<bool>("AddAntiParticle");
0025 }
0026
0027 Pythia6PtGun::~Pythia6PtGun() {}
0028
0029 void Pythia6PtGun::generateEvent(CLHEP::HepRandomEngine*) {
0030 Pythia6Service::InstanceWrapper guard(fPy6Service);
0031
0032
0033
0034
0035
0036
0037 HepMC::GenVertex* Vtx = new HepMC::GenVertex(HepMC::FourVector(0., 0., 0.));
0038
0039
0040
0041 fEvt = new HepMC::GenEvent();
0042
0043 int ip = 1;
0044 for (size_t i = 0; i < fPartIDs.size(); i++) {
0045 int particleID = fPartIDs[i];
0046 int py6PID = HepPID::translatePDTtoPythia(particleID);
0047 int dum = 0;
0048 double pt = 0, mom = 0, ee = 0, the = 0, eta = 0;
0049 double mass = pymass_(py6PID);
0050
0051
0052
0053
0054 pyjets.p[4][ip - 1] = mass;
0055
0056 double phi = (fMaxPhi - fMinPhi) * pyr_(&dum) + fMinPhi;
0057
0058 eta = (fMaxEta - fMinEta) * pyr_(&dum) + fMinEta;
0059
0060 the = 2. * atan(exp(-eta));
0061
0062 pt = (fMaxPt - fMinPt) * pyr_(&dum) + fMinPt;
0063
0064 mom = pt / sin(the);
0065 ee = sqrt(mom * mom + mass * mass);
0066
0067 py1ent_(ip, py6PID, ee, the, phi);
0068
0069 double px = pyjets.p[0][ip - 1];
0070 double py = pyjets.p[1][ip - 1];
0071 double pz = pyjets.p[2][ip - 1];
0072
0073 HepMC::FourVector p(px, py, pz, ee);
0074 HepMC::GenParticle* Part = new HepMC::GenParticle(p, particleID, 1);
0075 Part->suggest_barcode(ip);
0076 Vtx->add_particle_out(Part);
0077
0078 if (fAddAntiParticle) {
0079 ip = ip + 1;
0080 HepMC::GenParticle* APart = addAntiParticle(ip, particleID, ee, eta, phi);
0081 if (APart)
0082 Vtx->add_particle_out(APart);
0083 }
0084 ip++;
0085 }
0086
0087 fEvt->add_vertex(Vtx);
0088
0089
0090 pyexec_();
0091
0092 return;
0093 }
0094
0095 DEFINE_FWK_MODULE(Pythia6PtGun);