File indexing completed on 2024-04-06 12:06:48
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include <memory>
0010 #include <vector>
0011 #include <map>
0012 #include <set>
0013
0014
0015
0016 #include "FWCore/Utilities/interface/InputTag.h"
0017 #include "FWCore/Framework/interface/Frameworkfwd.h"
0018 #include "FWCore/Framework/interface/Event.h"
0019 #include "FWCore/Framework/interface/MakerMacros.h"
0020 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0021 #include "FWCore/Framework/interface/ESHandle.h"
0022 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0023
0024 #include "DPGAnalysis/Skims/interface/FilterOutScraping.h"
0025
0026 using namespace edm;
0027 using namespace std;
0028
0029 FilterOutScraping::FilterOutScraping(const edm::ParameterSet& iConfig) {
0030 applyfilter = iConfig.getUntrackedParameter<bool>("applyfilter", true);
0031 debugOn = iConfig.getUntrackedParameter<bool>("debugOn", false);
0032 thresh = iConfig.getUntrackedParameter<double>("thresh", 0.2);
0033 numtrack = iConfig.getUntrackedParameter<unsigned int>("numtrack", 10);
0034 tracks_ = consumes<reco::TrackCollection>(
0035 iConfig.getUntrackedParameter<edm::InputTag>("src", edm::InputTag("generalTracks")));
0036 }
0037
0038 FilterOutScraping::~FilterOutScraping() = default;
0039
0040 bool FilterOutScraping::filter(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0041 bool accepted = false;
0042 float fraction = 0;
0043
0044
0045 edm::Handle<reco::TrackCollection> tkRef;
0046 iEvent.getByToken(tracks_, tkRef);
0047 const reco::TrackCollection* tkColl = tkRef.product();
0048
0049
0050
0051 int numhighpurity = 0;
0052 _trackQuality = reco::TrackBase::qualityByName("highPurity");
0053
0054 if (tkColl->size() > numtrack) {
0055 reco::TrackCollection::const_iterator itk = tkColl->begin();
0056 reco::TrackCollection::const_iterator itk_e = tkColl->end();
0057 for (; itk != itk_e; ++itk) {
0058
0059 if (itk->quality(_trackQuality))
0060 numhighpurity++;
0061 }
0062 fraction = (float)numhighpurity / (float)tkColl->size();
0063 if (fraction > thresh)
0064 accepted = true;
0065 } else {
0066
0067 accepted = true;
0068 }
0069
0070 if (debugOn) {
0071 int ievt = iEvent.id().event();
0072 int irun = iEvent.id().run();
0073 int ils = iEvent.luminosityBlock();
0074 int bx = iEvent.bunchCrossing();
0075
0076 std::cout << "FilterOutScraping_debug: Run " << irun << " Event " << ievt << " Lumi Block " << ils
0077 << " Bunch Crossing " << bx << " Fraction " << fraction << " NTracks " << tkColl->size() << " Accepted "
0078 << accepted << std::endl;
0079 }
0080
0081 if (applyfilter)
0082 return accepted;
0083 else
0084 return true;
0085 }
0086
0087
0088 DEFINE_FWK_MODULE(FilterOutScraping);