File indexing completed on 2023-03-17 11:15:56
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026 #include <memory>
0027
0028
0029 #include "FWCore/Framework/interface/Frameworkfwd.h"
0030 #include "FWCore/Framework/interface/stream/EDFilter.h"
0031
0032 #include "FWCore/Framework/interface/Event.h"
0033 #include "FWCore/Framework/interface/MakerMacros.h"
0034 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0035
0036 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0037 #include "FWCore/Utilities/interface/StreamID.h"
0038 #include "FWCore/Utilities/interface/InputTag.h"
0039
0040
0041
0042
0043
0044 class ttHFGenFilter : public edm::stream::EDFilter<> {
0045 public:
0046 explicit ttHFGenFilter(const edm::ParameterSet&);
0047 ~ttHFGenFilter() override;
0048
0049 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0050
0051 private:
0052 void beginStream(edm::StreamID) override;
0053 bool filter(edm::Event&, const edm::EventSetup&) override;
0054 void endStream() override;
0055
0056 virtual bool HasAdditionalBHadron(const std::vector<int>&,
0057 const std::vector<int>&,
0058 const std::vector<reco::GenParticle>&,
0059 std::vector<const reco::Candidate*>&);
0060 virtual bool analyzeMothersRecursive(const reco::Candidate*, std::vector<const reco::Candidate*>& AllTopMothers);
0061 virtual std::vector<const reco::Candidate*> GetTops(const std::vector<reco::GenParticle>&,
0062 std::vector<const reco::Candidate*>& AllTopMothers);
0063 virtual void FindAllTopMothers(const reco::Candidate* particle, std::vector<const reco::Candidate*>&);
0064
0065 const edm::EDGetTokenT<reco::GenParticleCollection> genParticlesToken_;
0066 const edm::EDGetTokenT<std::vector<int> > genBHadFlavourToken_;
0067 const edm::EDGetTokenT<std::vector<int> > genBHadFromTopWeakDecayToken_;
0068 const edm::EDGetTokenT<std::vector<reco::GenParticle> > genBHadPlusMothersToken_;
0069 const edm::EDGetTokenT<std::vector<std::vector<int> > > genBHadPlusMothersIndicesToken_;
0070 const edm::EDGetTokenT<std::vector<int> > genBHadIndexToken_;
0071 bool OnlyHardProcessBHadrons_;
0072 bool taggingMode_;
0073 };
0074
0075
0076
0077
0078
0079
0080
0081
0082
0083
0084
0085
0086 ttHFGenFilter::ttHFGenFilter(const edm::ParameterSet& iConfig)
0087 : genParticlesToken_(consumes<reco::GenParticleCollection>(iConfig.getParameter<edm::InputTag>("genParticles"))),
0088 genBHadFlavourToken_(consumes<std::vector<int> >(iConfig.getParameter<edm::InputTag>("genBHadFlavour"))),
0089 genBHadFromTopWeakDecayToken_(
0090 consumes<std::vector<int> >(iConfig.getParameter<edm::InputTag>("genBHadFromTopWeakDecay"))),
0091 genBHadPlusMothersToken_(
0092 consumes<std::vector<reco::GenParticle> >(iConfig.getParameter<edm::InputTag>("genBHadPlusMothers"))),
0093 genBHadPlusMothersIndicesToken_(
0094 consumes<std::vector<std::vector<int> > >(iConfig.getParameter<edm::InputTag>("genBHadPlusMothersIndices"))),
0095 genBHadIndexToken_(consumes<std::vector<int> >(iConfig.getParameter<edm::InputTag>("genBHadIndex"))) {
0096
0097 OnlyHardProcessBHadrons_ = iConfig.getParameter<bool>("OnlyHardProcessBHadrons");
0098 taggingMode_ = iConfig.getParameter<bool>("taggingMode");
0099
0100 produces<bool>();
0101 }
0102
0103 ttHFGenFilter::~ttHFGenFilter() {
0104
0105
0106 }
0107
0108
0109
0110
0111
0112
0113 bool ttHFGenFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0114 using namespace std;
0115
0116
0117 edm::Handle<reco::GenParticleCollection> genParticles;
0118 iEvent.getByToken(genParticlesToken_, genParticles);
0119
0120
0121
0122
0123 edm::Handle<std::vector<int> > genBHadFlavour;
0124 iEvent.getByToken(genBHadFlavourToken_, genBHadFlavour);
0125
0126 edm::Handle<std::vector<int> > genBHadFromTopWeakDecay;
0127 iEvent.getByToken(genBHadFromTopWeakDecayToken_, genBHadFromTopWeakDecay);
0128
0129 edm::Handle<std::vector<reco::GenParticle> > genBHadPlusMothers;
0130 iEvent.getByToken(genBHadPlusMothersToken_, genBHadPlusMothers);
0131
0132 edm::Handle<std::vector<std::vector<int> > > genBHadPlusMothersIndices;
0133 iEvent.getByToken(genBHadPlusMothersIndicesToken_, genBHadPlusMothersIndices);
0134
0135 edm::Handle<std::vector<int> > genBHadIndex;
0136 iEvent.getByToken(genBHadIndexToken_, genBHadIndex);
0137
0138
0139
0140 std::vector<const reco::Candidate*> AllTopMothers;
0141 std::vector<const reco::Candidate*> Tops = GetTops(*genParticles, AllTopMothers);
0142
0143 bool pass = HasAdditionalBHadron(*genBHadIndex, *genBHadFlavour, *genBHadPlusMothers, AllTopMothers);
0144
0145 iEvent.put(std::make_unique<bool>(pass));
0146
0147 return taggingMode_ || pass;
0148 }
0149
0150 bool ttHFGenFilter::HasAdditionalBHadron(const std::vector<int>& genBHadIndex,
0151 const std::vector<int>& genBHadFlavour,
0152 const std::vector<reco::GenParticle>& genBHadPlusMothers,
0153 std::vector<const reco::Candidate*>& AllTopMothers) {
0154 for (unsigned int i = 0; i < genBHadIndex.size(); i++) {
0155 const reco::GenParticle* bhadron = genBHadIndex[i] >= 0 && genBHadIndex[i] < int(genBHadPlusMothers.size())
0156 ? &(genBHadPlusMothers[genBHadIndex[i]])
0157 : nullptr;
0158 int motherflav = genBHadFlavour[i];
0159 bool from_tth = (abs(motherflav) == 6 || abs(motherflav) == 25);
0160 bool fromhp = false;
0161
0162 if (!OnlyHardProcessBHadrons_) {
0163 if (bhadron != nullptr && !from_tth) {
0164 return true;
0165 }
0166 }
0167
0168 if (OnlyHardProcessBHadrons_) {
0169 if (bhadron != nullptr && !from_tth) {
0170
0171
0172 fromhp = analyzeMothersRecursive(bhadron, AllTopMothers);
0173 if (fromhp) {
0174 return true;
0175 }
0176 }
0177 if (i == genBHadIndex.size() - 1) {
0178 return false;
0179 }
0180 }
0181 }
0182
0183 return false;
0184 }
0185
0186
0187 bool ttHFGenFilter::analyzeMothersRecursive(const reco::Candidate* particle,
0188 std::vector<const reco::Candidate*>& AllTopMothers) {
0189 if (particle->status() > 20 && particle->status() < 30) {
0190 return true;
0191 }
0192 for (unsigned int k = 0; k < AllTopMothers.size(); k++) {
0193 if (particle == AllTopMothers[k]) {
0194 return true;
0195 }
0196 }
0197 bool isFromHardProcess = false;
0198 for (unsigned int i = 0; i < particle->numberOfMothers(); i++) {
0199 const reco::Candidate* mother = particle->mother(i);
0200 isFromHardProcess = analyzeMothersRecursive(mother, AllTopMothers);
0201 if (isFromHardProcess) {
0202 return true;
0203 }
0204 }
0205 return isFromHardProcess;
0206 }
0207
0208 std::vector<const reco::Candidate*> ttHFGenFilter::GetTops(const std::vector<reco::GenParticle>& genParticles,
0209 std::vector<const reco::Candidate*>& AllTopMothers) {
0210
0211 const reco::GenParticle* firstTop = nullptr;
0212 const reco::GenParticle* firstTopBar = nullptr;
0213 bool foundTop = false;
0214 bool foundTopBar = false;
0215 std::vector<const reco::GenParticle*> Tops;
0216
0217
0218
0219
0220
0221 for (reco::GenParticleCollection::const_iterator i_particle = genParticles.begin(); i_particle != genParticles.end();
0222 ++i_particle) {
0223 const reco::GenParticle* thisParticle = &*i_particle;
0224 if (!foundTop && thisParticle->pdgId() == 6) {
0225 firstTop = thisParticle;
0226 for (unsigned int i = 0; i < firstTop->numberOfMothers(); i++) {
0227 FindAllTopMothers(thisParticle->mother(i), AllTopMothers);
0228 }
0229 foundTop = true;
0230 } else if (!foundTopBar && thisParticle->pdgId() == -6) {
0231 firstTopBar = thisParticle;
0232 for (unsigned int i = 0; i < firstTopBar->numberOfMothers(); i++) {
0233 FindAllTopMothers(firstTopBar->mother(i), AllTopMothers);
0234 }
0235 foundTopBar = true;
0236 }
0237 if (foundTopBar && foundTop) {
0238
0239
0240
0241 break;
0242 }
0243 }
0244 return AllTopMothers;
0245 }
0246
0247
0248
0249 void ttHFGenFilter::FindAllTopMothers(const reco::Candidate* particle,
0250 std::vector<const reco::Candidate*>& AllTopMothers) {
0251
0252 for (unsigned int i = 0; i < particle->numberOfMothers(); i++) {
0253 if (abs(particle->mother(i)->pdgId()) != 6 && particle->mother(i)->pdgId() != 2212) {
0254 AllTopMothers.push_back(particle->mother(i));
0255 if (particle->mother(i)->pdgId() != 2212 || particle->mother(i)->numberOfMothers() > 1) {
0256
0257 FindAllTopMothers(particle->mother(i), AllTopMothers);
0258 }
0259 }
0260 }
0261 }
0262
0263
0264 void ttHFGenFilter::beginStream(edm::StreamID) {}
0265
0266
0267 void ttHFGenFilter::endStream() {}
0268
0269
0270
0271
0272
0273
0274
0275
0276
0277
0278
0279
0280
0281
0282
0283
0284
0285
0286
0287
0288
0289
0290
0291
0292
0293
0294
0295
0296
0297
0298
0299
0300
0301
0302 void ttHFGenFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0303
0304
0305 edm::ParameterSetDescription desc;
0306 desc.setUnknown();
0307 descriptions.addDefault(desc);
0308 }
0309
0310 DEFINE_FWK_MODULE(ttHFGenFilter);