CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
GenParticles2HepMCConverter Class Reference
Inheritance diagram for GenParticles2HepMCConverter:
edm::stream::EDProducer<>

Public Member Functions

void beginRun (edm::Run const &iRun, edm::EventSetup const &) override
 
 GenParticles2HepMCConverter (const edm::ParameterSet &pset)
 
void produce (edm::Event &event, const edm::EventSetup &eventSetup) override
 
 ~GenParticles2HepMCConverter () override
 
- Public Member Functions inherited from edm::stream::EDProducer<>
 EDProducer ()=default
 
 EDProducer (const EDProducer &)=delete
 
bool hasAbilityToProduceInBeginLumis () const final
 
bool hasAbilityToProduceInBeginProcessBlocks () const final
 
bool hasAbilityToProduceInBeginRuns () const final
 
bool hasAbilityToProduceInEndLumis () const final
 
bool hasAbilityToProduceInEndProcessBlocks () const final
 
bool hasAbilityToProduceInEndRuns () const final
 
const EDProduceroperator= (const EDProducer &)=delete
 

Private Member Functions

HepMC3::FourVector FourVector (const reco::Candidate::Point &point)
 
HepMC3::FourVector FourVector (const reco::Candidate::LorentzVector &lvec)
 

Private Attributes

const double cmEnergy_
 
edm::EDGetTokenT< GenEventInfoProductgenEventInfoToken_
 
edm::EDGetTokenT< reco::CandidateViewgenParticlesToken_
 
edm::EDGetTokenT< GenRunInfoProductgenRunInfoToken_
 
edm::ESGetToken< HepPDT::ParticleDataTable, PDTRecordpTable_
 
HepMC3::GenCrossSectionPtr xsec_
 

Additional Inherited Members

- Public Types inherited from edm::stream::EDProducer<>
using CacheTypes = CacheContexts< T... >
 
using GlobalCache = typename CacheTypes::GlobalCache
 
using HasAbility = AbilityChecker< T... >
 
using InputProcessBlockCache = typename CacheTypes::InputProcessBlockCache
 
using LuminosityBlockCache = typename CacheTypes::LuminosityBlockCache
 
using LuminosityBlockContext = LuminosityBlockContextT< LuminosityBlockCache, RunCache, GlobalCache >
 
using LuminosityBlockSummaryCache = typename CacheTypes::LuminosityBlockSummaryCache
 
using RunCache = typename CacheTypes::RunCache
 
using RunContext = RunContextT< RunCache, GlobalCache >
 
using RunSummaryCache = typename CacheTypes::RunSummaryCache
 

Detailed Description

Definition at line 31 of file GenParticles2HepMCConverter.cc.

Constructor & Destructor Documentation

◆ GenParticles2HepMCConverter()

GenParticles2HepMCConverter::GenParticles2HepMCConverter ( const edm::ParameterSet pset)
explicit

Definition at line 59 of file GenParticles2HepMCConverter.cc.

References genEventInfoToken_, genParticlesToken_, genRunInfoToken_, muonDTDigis_cfi::pset, and pTable_.

61  : cmEnergy_(pset.getUntrackedParameter<double>("cmEnergy", 13000)) {
62  genParticlesToken_ = consumes<reco::CandidateView>(pset.getParameter<edm::InputTag>("genParticles"));
63  genEventInfoToken_ = consumes<GenEventInfoProduct>(pset.getParameter<edm::InputTag>("genEventInfo"));
64  genRunInfoToken_ = consumes<GenRunInfoProduct, edm::InRun>(pset.getParameter<edm::InputTag>("genEventInfo"));
65  pTable_ = esConsumes<HepPDT::ParticleDataTable, PDTRecord>();
66 
67  produces<edm::HepMC3Product>("unsmeared");
68 }
edm::EDGetTokenT< GenRunInfoProduct > genRunInfoToken_
edm::EDGetTokenT< GenEventInfoProduct > genEventInfoToken_
edm::EDGetTokenT< reco::CandidateView > genParticlesToken_
edm::ESGetToken< HepPDT::ParticleDataTable, PDTRecord > pTable_

◆ ~GenParticles2HepMCConverter()

GenParticles2HepMCConverter::~GenParticles2HepMCConverter ( )
inlineoverride

Definition at line 34 of file GenParticles2HepMCConverter.cc.

34 {};

Member Function Documentation

◆ beginRun()

void GenParticles2HepMCConverter::beginRun ( edm::Run const &  iRun,
edm::EventSetup const &   
)
override

Definition at line 70 of file GenParticles2HepMCConverter.cc.

References GenRunInfoProduct::XSec::error(), genRunInfoToken_, edm::Run::getByToken(), GenRunInfoProduct::internalXSec(), edm::HandleBase::isValid(), GenRunInfoProduct::XSec::value(), and xsec_.

70  {
71  edm::Handle<GenRunInfoProduct> genRunInfoHandle;
72  iRun.getByToken(genRunInfoToken_, genRunInfoHandle);
73 
74  xsec_ = make_shared<HepMC3::GenCrossSection>();
75  if (genRunInfoHandle.isValid()) {
76  xsec_->set_cross_section(genRunInfoHandle->internalXSec().value(), genRunInfoHandle->internalXSec().error());
77  } else {
78  // dummy cross section
79  xsec_->set_cross_section(1., 0.);
80  }
81 }
edm::EDGetTokenT< GenRunInfoProduct > genRunInfoToken_
const XSec & internalXSec() const
bool isValid() const
Definition: HandleBase.h:70

◆ FourVector() [1/2]

HepMC3::FourVector GenParticles2HepMCConverter::FourVector ( const reco::Candidate::Point point)
inlineprivate

Definition at line 49 of file GenParticles2HepMCConverter.cc.

References point.

Referenced by produce().

49  {
50  return HepMC3::FourVector(10 * point.x(), 10 * point.y(), 10 * point.z(), 0);
51  };
*vegas h *****************************************************used in the default bin number in original ***version of VEGAS is ***a higher bin number might help to derive a more precise ***grade subtle point
Definition: invegas.h:5

◆ FourVector() [2/2]

HepMC3::FourVector GenParticles2HepMCConverter::FourVector ( const reco::Candidate::LorentzVector lvec)
inlineprivate

Definition at line 53 of file GenParticles2HepMCConverter.cc.

References Matriplex::hypot(), and SiStripPI::max.

53  {
54  // Avoid negative mass, set minimum m^2 = 0
55  return HepMC3::FourVector(lvec.px(), lvec.py(), lvec.pz(), std::hypot(lvec.P(), std::max(0., lvec.mass())));
56  };
MPlex< T, D1, D2, N > hypot(const MPlex< T, D1, D2, N > &a, const MPlex< T, D1, D2, N > &b)
Definition: Matriplex.h:436

◆ produce()

void GenParticles2HepMCConverter::produce ( edm::Event event,
const edm::EventSetup eventSetup 
)
override

Definition at line 83 of file GenParticles2HepMCConverter.cc.

References funct::abs(), GenEventInfoProduct::alphaQCD(), GenEventInfoProduct::alphaQED(), edm::View< T >::at(), cmEnergy_, options_cfi::eventSetup, FourVector(), genEventInfoToken_, genParticlesToken_, mps_fire::i, gen::PdfInfo::id, dqmiolumiharvest::j, EgHLTOffHistBins_cfi::mass, eostools::move(), dqmiodumpmetadata::n, or, AlCaHLTBitMon_ParallelJobs::p, GenEventInfoProduct::pdf(), pTable_, GenEventInfoProduct::qScale(), gen::PdfInfo::scalePDF, GenEventInfoProduct::signalProcessID(), edm::View< T >::size(), bphysicsOniaDQM_cfi::vertex, reco::Candidate::vertex(), GenEventInfoProduct::weights(), gen::PdfInfo::x, gen::PdfInfo::xPDF, and xsec_.

83  {
84  edm::Handle<reco::CandidateView> genParticlesHandle;
85  event.getByToken(genParticlesToken_, genParticlesHandle);
86 
87  edm::Handle<GenEventInfoProduct> genEventInfoHandle;
88  event.getByToken(genEventInfoToken_, genEventInfoHandle);
89 
90  auto const& pTableData = eventSetup.getData(pTable_);
91 
92  HepMC3::GenEvent* hepmc_event = new HepMC3::GenEvent();
93  hepmc_event->set_event_number(event.id().event());
94  hepmc_event->add_attribute("signal_process_id",
95  std::make_shared<HepMC3::IntAttribute>(genEventInfoHandle->signalProcessID()));
96  hepmc_event->add_attribute("event_scale", std::make_shared<HepMC3::DoubleAttribute>(genEventInfoHandle->qScale()));
97  hepmc_event->add_attribute("alphaQCD", std::make_shared<HepMC3::DoubleAttribute>(genEventInfoHandle->alphaQCD()));
98  hepmc_event->add_attribute("alphaQED", std::make_shared<HepMC3::DoubleAttribute>(genEventInfoHandle->alphaQED()));
99 
100  hepmc_event->weights() = genEventInfoHandle->weights();
101  // add dummy weight if necessary
102  if (hepmc_event->weights().size() == 0) {
103  hepmc_event->weights().push_back(1.);
104  }
105 
106  // resize cross section to number of weights
107  if (xsec_->xsecs().size() < hepmc_event->weights().size()) {
108  xsec_->set_cross_section(std::vector<double>(hepmc_event->weights().size(), xsec_->xsec(0)),
109  std::vector<double>(hepmc_event->weights().size(), xsec_->xsec_err(0)));
110  }
111  hepmc_event->set_cross_section(xsec_);
112 
113  // Set PDF
114  const gen::PdfInfo* pdf = genEventInfoHandle->pdf();
115  if (pdf != nullptr) {
116  const int pdf_id1 = pdf->id.first, pdf_id2 = pdf->id.second;
117  const double pdf_x1 = pdf->x.first, pdf_x2 = pdf->x.second;
118  const double pdf_scalePDF = pdf->scalePDF;
119  const double pdf_xPDF1 = pdf->xPDF.first, pdf_xPDF2 = pdf->xPDF.second;
120  HepMC3::GenPdfInfoPtr hepmc_pdfInfo = make_shared<HepMC3::GenPdfInfo>();
121  hepmc_pdfInfo->set(pdf_id1, pdf_id2, pdf_x1, pdf_x2, pdf_scalePDF, pdf_xPDF1, pdf_xPDF2);
122  hepmc_event->set_pdf_info(hepmc_pdfInfo);
123  }
124 
125  // Prepare list of HepMC3::GenParticles
126  std::map<const reco::Candidate*, HepMC3::GenParticlePtr> genCandToHepMCMap;
127  HepMC3::GenParticlePtr hepmc_parton1, hepmc_parton2;
128  std::vector<HepMC3::GenParticlePtr> hepmc_particles;
129  const reco::Candidate *parton1 = nullptr, *parton2 = nullptr;
130  for (unsigned int i = 0, n = genParticlesHandle->size(); i < n; ++i) {
131  const reco::Candidate* p = &genParticlesHandle->at(i);
132  HepMC3::GenParticlePtr hepmc_particle =
133  std::make_shared<HepMC3::GenParticle>(FourVector(p->p4()), p->pdgId(), p->status());
134 
135  // Assign particle's generated mass from the standard particle data table
136  double particleMass;
137  if (pTableData.particle(p->pdgId()))
138  particleMass = pTableData.particle(p->pdgId())->mass();
139  else
140  particleMass = p->mass();
141 
142  hepmc_particle->set_generated_mass(particleMass);
143 
144  hepmc_particles.push_back(hepmc_particle);
145  genCandToHepMCMap[p] = hepmc_particle;
146 
147  // Find incident proton pair
148  if (p->mother() == nullptr and std::abs(p->eta()) > 5 and std::abs(p->pz()) > 1000) {
149  if (!parton1 and p->pz() > 0) {
150  parton1 = p;
151  hepmc_parton1 = hepmc_particle;
152  } else if (!parton2 and p->pz() < 0) {
153  parton2 = p;
154  hepmc_parton2 = hepmc_particle;
155  }
156  }
157  }
158 
159  HepMC3::GenVertexPtr vertex1;
160  HepMC3::GenVertexPtr vertex2;
161  if (parton1 == nullptr || parton2 == nullptr) {
162  // Particle gun samples do not have incident partons. Put dummy incident particle and prod vertex
163  // Note: leave parton1 and parton2 as nullptr since it is not used anymore after creating hepmc_parton1 and 2
164  const reco::Candidate::LorentzVector nullP4(0, 0, 0, 0);
165  const reco::Candidate::LorentzVector beamP4(0, 0, cmEnergy_ / 2, cmEnergy_ / 2);
166  vertex1 = make_shared<HepMC3::GenVertex>(FourVector(nullP4));
167  vertex2 = make_shared<HepMC3::GenVertex>(FourVector(nullP4));
168  hepmc_parton1 = make_shared<HepMC3::GenParticle>(FourVector(+beamP4), 2212, 4);
169  hepmc_parton2 = make_shared<HepMC3::GenParticle>(FourVector(-beamP4), 2212, 4);
170  } else {
171  // Put incident beam particles : proton -> parton vertex
172  vertex1 = make_shared<HepMC3::GenVertex>(FourVector(parton1->vertex()));
173  vertex2 = make_shared<HepMC3::GenVertex>(FourVector(parton2->vertex()));
174  }
175  hepmc_event->add_vertex(vertex1);
176  hepmc_event->add_vertex(vertex2);
177  vertex1->add_particle_in(hepmc_parton1);
178  vertex2->add_particle_in(hepmc_parton2);
179  //hepmc_event->set_beam_particles(hepmc_parton1, hepmc_parton2);
180 
181  // Prepare vertex list
182  typedef std::map<const reco::Candidate*, HepMC3::GenVertexPtr> ParticleToVertexMap;
183  ParticleToVertexMap particleToVertexMap;
184  particleToVertexMap[parton1] = vertex1;
185  particleToVertexMap[parton2] = vertex2;
186  for (unsigned int i = 0, n = genParticlesHandle->size(); i < n; ++i) {
187  const reco::Candidate* p = &genParticlesHandle->at(i);
188  if (p == parton1 or p == parton2)
189  continue;
190 
191  // Connect mother-daughters for the other cases
192  for (unsigned int j = 0, nMothers = p->numberOfMothers(); j < nMothers; ++j) {
193  // Mother-daughter hierarchy defines vertex
194  const reco::Candidate* elder = p->mother(j)->daughter(0);
195  HepMC3::GenVertexPtr vertex;
196  if (particleToVertexMap.find(elder) == particleToVertexMap.end()) {
197  vertex = make_shared<HepMC3::GenVertex>(FourVector(elder->vertex()));
198  hepmc_event->add_vertex(vertex);
199  particleToVertexMap[elder] = vertex;
200  } else {
201  vertex = particleToVertexMap[elder];
202  }
203 
204  // Vertex is found. Now connect each other
205  const reco::Candidate* mother = p->mother(j);
206  vertex->add_particle_in(genCandToHepMCMap[mother]);
207  vertex->add_particle_out(hepmc_particles[i]);
208  }
209  }
210 
211  // Finalize HepMC event record
212  std::unique_ptr<edm::HepMC3Product> hepmc_product(new edm::HepMC3Product());
213  hepmc_product->addHepMCData(hepmc_event);
214  event.put(std::move(hepmc_product), "unsmeared");
215 }
edm::EDGetTokenT< GenEventInfoProduct > genEventInfoToken_
std::pair< double, double > x
Definition: PdfInfo.h:13
edm::EDGetTokenT< reco::CandidateView > genParticlesToken_
size_type size() const
edm::ESGetToken< HepPDT::ParticleDataTable, PDTRecord > pTable_
std::pair< double, double > xPDF
Definition: PdfInfo.h:14
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
virtual const Point & vertex() const =0
vertex position
unsigned int signalProcessID() const
const_reference at(size_type pos) const
edm::Ptr< GenParticle > GenParticlePtr
persistent reference to a GenParticle
std::pair< int, int > id
Definition: PdfInfo.h:12
HepMC3::FourVector FourVector(const reco::Candidate::Point &point)
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Candidate.h:36
const PDF * pdf() const
std::vector< double > & weights()
double scalePDF
Definition: PdfInfo.h:15
def move(src, dest)
Definition: eostools.py:511
Definition: event.py:1

Member Data Documentation

◆ cmEnergy_

const double GenParticles2HepMCConverter::cmEnergy_
private

Definition at line 45 of file GenParticles2HepMCConverter.cc.

Referenced by produce().

◆ genEventInfoToken_

edm::EDGetTokenT<GenEventInfoProduct> GenParticles2HepMCConverter::genEventInfoToken_
private

Definition at line 41 of file GenParticles2HepMCConverter.cc.

Referenced by GenParticles2HepMCConverter(), and produce().

◆ genParticlesToken_

edm::EDGetTokenT<reco::CandidateView> GenParticles2HepMCConverter::genParticlesToken_
private

Definition at line 40 of file GenParticles2HepMCConverter.cc.

Referenced by GenParticles2HepMCConverter(), and produce().

◆ genRunInfoToken_

edm::EDGetTokenT<GenRunInfoProduct> GenParticles2HepMCConverter::genRunInfoToken_
private

Definition at line 42 of file GenParticles2HepMCConverter.cc.

Referenced by beginRun(), and GenParticles2HepMCConverter().

◆ pTable_

edm::ESGetToken<HepPDT::ParticleDataTable, PDTRecord> GenParticles2HepMCConverter::pTable_
private

Definition at line 43 of file GenParticles2HepMCConverter.cc.

Referenced by GenParticles2HepMCConverter(), and produce().

◆ xsec_

HepMC3::GenCrossSectionPtr GenParticles2HepMCConverter::xsec_
private

Definition at line 46 of file GenParticles2HepMCConverter.cc.

Referenced by beginRun(), and produce().