File indexing completed on 2024-04-06 12:13:44
0001 #include <iostream>
0002 #include "GeneratorInterface/HiGenCommon/interface/PartonHadronDecayGenEvtSelector.h"
0003 #include "FWCore/Utilities/interface/EDMException.h"
0004 using namespace std;
0005
0006 PartonHadronDecayGenEvtSelector::PartonHadronDecayGenEvtSelector(const edm::ParameterSet& pset)
0007 : BaseHiGenEvtSelector(pset) {
0008 hadronId_ = pset.getParameter<vector<int> >("hadrons");
0009 hadronStatus_ = pset.getParameter<vector<int> >("hadronStatus");
0010 hadronEtaMax_ = pset.getParameter<vector<double> >("hadronEtaMax");
0011 hadronEtaMin_ = pset.getParameter<vector<double> >("hadronEtaMin");
0012 hadronPMin_ = pset.getParameter<vector<double> >("hadronPMin");
0013 hadronPtMax_ = pset.getParameter<vector<double> >("hadronPtMax");
0014 hadronPtMin_ = pset.getParameter<vector<double> >("hadronPtMin");
0015
0016 decayId_ = pset.getParameter<int>("decays");
0017 decayStatus_ = pset.getParameter<int>("decayStatus");
0018 decayEtaMax_ = pset.getParameter<double>("decayEtaMax");
0019 decayEtaMin_ = pset.getParameter<double>("decayEtaMin");
0020 decayPMin_ = pset.getParameter<double>("decayPMin");
0021 decayPtMax_ = pset.getParameter<double>("decayPtMax");
0022 decayPtMin_ = pset.getParameter<double>("decayPtMin");
0023 decayNtrig_ = pset.getParameter<int>("decayNtrig");
0024
0025 partonId_ = pset.getParameter<vector<int> >("partons");
0026 partonStatus_ = pset.getParameter<vector<int> >("partonStatus");
0027 partonEtaMax_ = pset.getParameter<vector<double> >("partonEtaMax");
0028 partonPtMin_ = pset.getParameter<vector<double> >("partonPtMin");
0029
0030 int id = hadronId_.size();
0031 int st = hadronStatus_.size();
0032 int etamax = hadronEtaMax_.size();
0033 int etamin = hadronEtaMin_.size();
0034 int pmin = hadronPMin_.size();
0035 int ptmax = hadronPtMax_.size();
0036 int ptmin = hadronPtMin_.size();
0037
0038 if (id != st || id != etamax || id != etamin || id != ptmax || id != ptmin || id != pmin) {
0039 throw edm::Exception(edm::errors::LogicError)
0040 << "Hadron selection parameters: " << id << st << etamax << etamin << pmin << ptmax << ptmin << endl;
0041 }
0042
0043 id = partonId_.size();
0044 st = partonStatus_.size();
0045 etamax = partonEtaMax_.size();
0046 ptmin = partonPtMin_.size();
0047
0048 if (id != st || id != etamax || id != ptmin) {
0049 throw edm::Exception(edm::errors::LogicError)
0050 << "Parton selection parameters: " << id << st << etamax << ptmin << endl;
0051 }
0052 }
0053
0054
0055 bool PartonHadronDecayGenEvtSelector::filter(HepMC::GenEvent* evt) {
0056
0057
0058 HepMC::GenEvent::particle_const_iterator begin = evt->particles_begin();
0059 HepMC::GenEvent::particle_const_iterator end = evt->particles_end();
0060
0061 bool foundHadron = false;
0062 bool foundDecay = false;
0063 bool foundParton = false;
0064
0065 HepMC::GenEvent::particle_const_iterator it = begin;
0066 while (!foundParton && it != end) {
0067 for (unsigned i = 0; i < partonId_.size(); ++i) {
0068 if (selectParticle(*it, partonStatus_[i], partonId_[i], partonPtMin_[i], partonEtaMax_[i]))
0069 foundParton = true;
0070 }
0071 ++it;
0072 }
0073
0074 int foundtrig = 0;
0075 HepMC::GenEvent::particle_const_iterator it2 = begin;
0076
0077 if (foundParton) {
0078 while ((!foundHadron || !foundDecay) && it2 != end) {
0079 for (unsigned i = 0; i < hadronId_.size(); ++i) {
0080 if (selectParticle(*it2,
0081 hadronStatus_[i],
0082 hadronId_[i],
0083 hadronEtaMax_[i],
0084 hadronEtaMin_[i],
0085 hadronPMin_[i],
0086 hadronPtMax_[i],
0087 hadronPtMin_[i]))
0088 foundHadron = true;
0089 }
0090
0091 if (selectParticle(
0092 *it2, decayStatus_, decayId_, decayEtaMax_, decayEtaMin_, decayPMin_, decayPtMax_, decayPtMin_))
0093 foundtrig++;
0094 if (decayNtrig_ == foundtrig)
0095 foundDecay = true;
0096
0097 ++it2;
0098 }
0099 }
0100
0101 return (foundHadron && foundDecay && foundParton);
0102 }
0103
0104