CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_1/src/GeneratorInterface/Pythia6Interface/plugins/Pythia6PtGun.cc

Go to the documentation of this file.
00001 
00002 #include <iostream>
00003 
00004 #include "Pythia6PtGun.h"
00005 
00006 #include "FWCore/Utilities/interface/Exception.h"
00007 
00008 #include "FWCore/Framework/interface/EDProducer.h"
00009 #include "FWCore/Framework/interface/EventSetup.h"
00010 
00011 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
00012 
00013 #include "FWCore/Framework/interface/MakerMacros.h"
00014 
00015 using namespace edm;
00016 using namespace gen;
00017 
00018 Pythia6PtGun::Pythia6PtGun( const ParameterSet& pset ) :
00019    Pythia6ParticleGun(pset)
00020 {
00021    
00022    // ParameterSet defpset ;
00023    ParameterSet pgun_params = 
00024       pset.getParameter<ParameterSet>("PGunParameters"); //, defpset ) ;
00025    fMinEta     = pgun_params.getParameter<double>("MinEta"); // ,-2.2);
00026    fMaxEta     = pgun_params.getParameter<double>("MaxEta"); // , 2.2);
00027    fMinPt       = pgun_params.getParameter<double>("MinPt"); // ,  20.);
00028    fMaxPt       = pgun_params.getParameter<double>("MaxPt"); // , 420.);
00029    fAddAntiParticle = pgun_params.getParameter<bool>("AddAntiParticle"); //, false) ;  
00030 
00031 }
00032 
00033 Pythia6PtGun::~Pythia6PtGun()
00034 {
00035 }
00036 
00037 void Pythia6PtGun::generateEvent()
00038 {
00039    
00040    Pythia6Service::InstanceWrapper guard(fPy6Service);  // grab Py6 instance
00041 
00042    // now actualy, start cooking up the event gun 
00043    //
00044 
00045    // 1st, primary vertex
00046    //
00047    HepMC::GenVertex* Vtx = new HepMC::GenVertex( HepMC::FourVector(0.,0.,0.));
00048 
00049    // here re-create fEvt (memory)
00050    //
00051    fEvt = new HepMC::GenEvent() ;
00052      
00053    int ip=1;
00054    for ( size_t i=0; i<fPartIDs.size(); i++ )
00055    {
00056          int particleID = fPartIDs[i]; // this is PDG - need to convert to Py6 !!!
00057          int py6PID = HepPID::translatePDTtoPythia( particleID );
00058          int dum = 0;
00059          double pt=0, mom=0, ee=0, the=0, eta=0;
00060          double mass = pymass_(py6PID);
00061          
00062          // fill p(ip,5) (in PYJETS) with mass value right now,
00063          // because the (hardcoded) mstu(10)=1 will make py1ent
00064          // pick the mass from there
00065          pyjets.p[4][ip-1]=mass; 
00066                  
00067          double phi = (fMaxPhi-fMinPhi)*pyr_(&dum)+fMinPhi;
00068 
00069          eta  = (fMaxEta-fMinEta)*pyr_(&dum)+fMinEta;                                                      
00070 
00071          the  = 2.*atan(exp(-eta));                                                                          
00072          
00073          pt = (fMaxPt-fMinPt)*pyr_(&dum)+fMinPt;
00074          
00075          mom = pt/sin(the);
00076          ee = sqrt(mom*mom+mass*mass);
00077 
00078          py1ent_(ip, py6PID, ee, the, phi);
00079          
00080          double px     = pyjets.p[0][ip-1]; // pt*cos(phi) ;
00081          double py     = pyjets.p[1][ip-1]; // pt*sin(phi) ;
00082          double pz     = pyjets.p[2][ip-1]; // mom*cos(the) ;
00083          
00084          HepMC::FourVector p(px,py,pz,ee) ;
00085          HepMC::GenParticle* Part = 
00086              new HepMC::GenParticle(p,particleID,1);
00087          Part->suggest_barcode( ip ) ;
00088          Vtx->add_particle_out(Part);
00089          
00090          if(fAddAntiParticle)
00091          {
00092             ip = ip + 1;
00093             HepMC::GenParticle* APart = addAntiParticle( ip, particleID, ee, eta, phi );
00094             if ( APart ) Vtx->add_particle_out(APart) ;     
00095          }
00096          ip++;
00097    }
00098    
00099    fEvt->add_vertex(Vtx);
00100      
00101    // run pythia
00102    pyexec_();
00103    
00104    return;
00105 }
00106 
00107 DEFINE_FWK_MODULE(Pythia6PtGun);