File indexing completed on 2023-03-17 10:59:16
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include <memory>
0021
0022
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/stream/EDAnalyzer.h"
0025
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/MakerMacros.h"
0028
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030
0031 #include "DataFormats/Common/interface/ValueMap.h"
0032 #include "DataFormats/Common/interface/View.h"
0033
0034 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0035 #include "DataFormats/VertexReco/interface/Vertex.h"
0036 #include "DataFormats/PatCandidates/interface/Electron.h"
0037
0038 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0039 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0040 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
0041
0042 #include "DataFormats/Candidate/interface/Candidate.h"
0043 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0044
0045 #include "DataFormats/EgammaCandidates/interface/ConversionFwd.h"
0046 #include "DataFormats/EgammaCandidates/interface/Conversion.h"
0047 #include "CommonTools/Egamma/interface/ConversionTools.h"
0048
0049 #include "FWCore/ServiceRegistry/interface/Service.h"
0050 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0051
0052 #include "TTree.h"
0053 #include "Math/VectorUtil.h"
0054
0055
0056
0057
0058
0059 class ElectronIDValidationAnalyzer : public edm::stream::EDAnalyzer<> {
0060 public:
0061 explicit ElectronIDValidationAnalyzer(const edm::ParameterSet &);
0062 ~ElectronIDValidationAnalyzer();
0063
0064 static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0065
0066 enum ElectronMatchType {
0067 UNMATCHED = 0,
0068 TRUE_PROMPT_ELECTRON,
0069 TRUE_ELECTRON_FROM_TAU,
0070 TRUE_NON_PROMPT_ELECTRON
0071 };
0072
0073 private:
0074 virtual void analyze(const edm::Event &, const edm::EventSetup &) override;
0075
0076 int matchToTruth(const reco::GsfElectron &el, const edm::Handle<edm::View<reco::GenParticle>> &genParticles);
0077 void findFirstNonElectronMother(const reco::Candidate *particle, int &ancestorPID, int &ancestorStatus);
0078
0079
0080 edm::EDGetTokenT<reco::VertexCollection> vtxToken_;
0081 edm::EDGetTokenT<edm::View<reco::GenParticle>> genToken_;
0082 edm::EDGetTokenT<reco::ConversionCollection> convToken_;
0083 edm::EDGetTokenT<reco::BeamSpot> beamToken_;
0084 edm::EDGetTokenT<edm::ValueMap<float>> full5x5SigmaIEtaIEtaMapToken_;
0085 edm::EDGetTokenT<edm::View<reco::GsfElectron>> electronCollectionToken_;
0086 edm::EDGetTokenT<edm::ValueMap<bool>> electronIdToken_;
0087
0088 TTree *electronTree_;
0089 Float_t pt_;
0090 Float_t etaSC_;
0091
0092 Float_t dEtaIn_;
0093 Float_t dPhiIn_;
0094 Float_t hOverE_;
0095 Float_t sigmaIetaIeta_;
0096 Float_t full5x5_sigmaIetaIeta_;
0097 Float_t relIsoWithDBeta_;
0098 Float_t ooEmooP_;
0099 Float_t d0_;
0100 Float_t dz_;
0101 Int_t expectedMissingInnerHits_;
0102 Int_t passConversionVeto_;
0103
0104 Int_t isTrue_;
0105 Int_t isPass_;
0106 };
0107
0108
0109
0110
0111
0112
0113
0114
0115
0116
0117
0118
0119 ElectronIDValidationAnalyzer::ElectronIDValidationAnalyzer(const edm::ParameterSet &iConfig)
0120 : vtxToken_(consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("vertices"))),
0121 genToken_(consumes<edm::View<reco::GenParticle>>(iConfig.getParameter<edm::InputTag>("genparticles"))),
0122 convToken_(consumes<reco::ConversionCollection>(iConfig.getParameter<edm::InputTag>("convcollection"))),
0123 beamToken_(consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("beamspot"))),
0124 full5x5SigmaIEtaIEtaMapToken_(
0125 consumes<edm::ValueMap<float>>(iConfig.getParameter<edm::InputTag>("full5x5SigmaIEtaIEtaMap"))),
0126 electronCollectionToken_(
0127 consumes<edm::View<reco::GsfElectron>>(iConfig.getParameter<edm::InputTag>("electrons"))),
0128 electronIdToken_(consumes<edm::ValueMap<bool>>(iConfig.getParameter<edm::InputTag>("electronIDs"))) {
0129
0130 edm::Service<TFileService> fs;
0131 electronTree_ = fs->make<TTree>("ElectronTree", "Electron data");
0132
0133 electronTree_->Branch("pt", &pt_, "pt/F");
0134 electronTree_->Branch("etaSC", &etaSC_, "etaSC/F");
0135
0136 electronTree_->Branch("dEtaIn", &dEtaIn_, "dEtaIn/F");
0137 electronTree_->Branch("dPhiIn", &dPhiIn_, "dPhiIn/F");
0138 electronTree_->Branch("hOverE", &hOverE_, "hOverE/F");
0139 electronTree_->Branch("sigmaIetaIeta", &sigmaIetaIeta_, "sigmaIetaIeta/F");
0140 electronTree_->Branch("full5x5_sigmaIetaIeta", &full5x5_sigmaIetaIeta_, "full5x5_sigmaIetaIeta/F");
0141 electronTree_->Branch("relIsoWithDBeta", &relIsoWithDBeta_, "relIsoWithDBeta/F");
0142 electronTree_->Branch("ooEmooP", &ooEmooP_, "ooEmooP/F");
0143 electronTree_->Branch("d0", &d0_, "d0/F");
0144 electronTree_->Branch("dz", &dz_, "dz/F");
0145 electronTree_->Branch("expectedMissingInnerHits", &expectedMissingInnerHits_, "expectedMissingInnerHits/I");
0146 electronTree_->Branch("passConversionVeto", &passConversionVeto_, "passConversionVeto/I");
0147
0148 electronTree_->Branch("isTrue", &isTrue_, "isTrue/I");
0149 electronTree_->Branch("isPass", &isPass_, "isPass/I");
0150 }
0151
0152 ElectronIDValidationAnalyzer::~ElectronIDValidationAnalyzer() {
0153
0154
0155 }
0156
0157
0158
0159
0160
0161
0162 void ElectronIDValidationAnalyzer::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0163 edm::Handle<edm::ValueMap<float>> full5x5sieie;
0164 edm::Handle<edm::View<reco::GsfElectron>> collection;
0165 edm::Handle<edm::ValueMap<bool>> id_decisions;
0166 iEvent.getByToken(full5x5SigmaIEtaIEtaMapToken_, full5x5sieie);
0167 iEvent.getByToken(electronCollectionToken_, collection);
0168 iEvent.getByToken(electronIdToken_, id_decisions);
0169
0170
0171 edm::Handle<reco::VertexCollection> vertices;
0172 iEvent.getByToken(vtxToken_, vertices);
0173 if (vertices->empty())
0174 return;
0175 const reco::Vertex &pv = vertices->front();
0176
0177
0178 edm::Handle<edm::View<reco::GenParticle>> genParticles;
0179 iEvent.getByToken(genToken_, genParticles);
0180
0181
0182 edm::Handle<reco::ConversionCollection> convs;
0183 edm::Handle<reco::BeamSpot> thebs;
0184 iEvent.getByToken(convToken_, convs);
0185 iEvent.getByToken(beamToken_, thebs);
0186
0187 for (size_t i = 0; i < collection->size(); ++i) {
0188 auto el = collection->refAt(i);
0189 pt_ = el->pt();
0190 etaSC_ = el->superCluster()->eta();
0191
0192
0193 dEtaIn_ = el->deltaEtaSuperClusterTrackAtVtx();
0194 dPhiIn_ = el->deltaPhiSuperClusterTrackAtVtx();
0195 hOverE_ = el->hcalOverEcal();
0196 sigmaIetaIeta_ = el->sigmaIetaIeta();
0197 full5x5_sigmaIetaIeta_ = (*full5x5sieie)[el];
0198
0199
0200
0201 if (el->ecalEnergy() == 0) {
0202 std::cout << "Electron energy is zero!\n";
0203 ooEmooP_ = 1e30;
0204 } else if (!std::isfinite(el->ecalEnergy())) {
0205 std::cout << "Electron energy is not finite!\n";
0206 ooEmooP_ = 1e30;
0207 } else {
0208 ooEmooP_ = fabs(1.0 / el->ecalEnergy() - el->eSuperClusterOverP() / el->ecalEnergy());
0209 }
0210
0211
0212 reco::GsfElectron::PflowIsolationVariables pfIso = el->pfIsolationVariables();
0213
0214 float absiso =
0215 pfIso.sumChargedHadronPt + std::max(0.0, pfIso.sumNeutralHadronEt + pfIso.sumPhotonEt - 0.5 * pfIso.sumPUPt);
0216 relIsoWithDBeta_ = absiso / pt_;
0217
0218
0219 d0_ = (-1) * el->gsfTrack()->dxy(pv.position());
0220 dz_ = el->gsfTrack()->dz(pv.position());
0221
0222
0223 constexpr reco::HitPattern::HitCategory missingHitType = reco::HitPattern::MISSING_INNER_HITS;
0224 expectedMissingInnerHits_ = el->gsfTrack()->hitPattern().numberOfLostHits(missingHitType);
0225 passConversionVeto_ = false;
0226 if (thebs.isValid() && convs.isValid()) {
0227 passConversionVeto_ = !ConversionTools::hasMatchedConversion(*el, *convs, thebs->position());
0228 } else {
0229 std::cout << "\n\nERROR!!! conversions not found!!!\n";
0230 }
0231
0232 isTrue_ = matchToTruth(*el, genParticles);
0233 isPass_ = (*id_decisions)[el];
0234
0235 electronTree_->Fill();
0236 }
0237 }
0238
0239
0240 void ElectronIDValidationAnalyzer::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0241
0242
0243 edm::ParameterSetDescription desc;
0244 desc.setUnknown();
0245 descriptions.addDefault(desc);
0246 }
0247
0248
0249
0250 int ElectronIDValidationAnalyzer::matchToTruth(const reco::GsfElectron &el,
0251 const edm::Handle<edm::View<reco::GenParticle>> &genParticles) {
0252
0253
0254
0255
0256
0257 double dR = 999;
0258 const reco::Candidate *closestElectron = 0;
0259 for (size_t i = 0; i < genParticles->size(); i++) {
0260 const reco::Candidate *particle = &(*genParticles)[i];
0261
0262 if (abs(particle->pdgId()) != 11 || particle->status() != 1)
0263 continue;
0264
0265 double dRtmp = ROOT::Math::VectorUtil::DeltaR(el.p4(), particle->p4());
0266 if (dRtmp < dR) {
0267 dR = dRtmp;
0268 closestElectron = particle;
0269 }
0270 }
0271
0272
0273 if (!(closestElectron != 0 && dR < 0.1)) {
0274 return UNMATCHED;
0275 }
0276
0277
0278 int ancestorPID = -999;
0279 int ancestorStatus = -999;
0280 findFirstNonElectronMother(closestElectron, ancestorPID, ancestorStatus);
0281
0282 if (ancestorPID == -999 && ancestorStatus == -999) {
0283
0284
0285 std::cout << "ElectronNtupler: ERROR! Electron does not apper to have a non-electron parent\n";
0286 return UNMATCHED;
0287 }
0288
0289 if (abs(ancestorPID) > 50 && ancestorStatus == 2)
0290 return TRUE_NON_PROMPT_ELECTRON;
0291
0292 if (abs(ancestorPID) == 15 && ancestorStatus == 2)
0293 return TRUE_ELECTRON_FROM_TAU;
0294
0295
0296 return TRUE_PROMPT_ELECTRON;
0297 }
0298
0299 void ElectronIDValidationAnalyzer::findFirstNonElectronMother(const reco::Candidate *particle,
0300 int &ancestorPID,
0301 int &ancestorStatus) {
0302 if (particle == 0) {
0303 std::cout << "ElectronNtupler: ERROR! null candidate pointer, this should never happen\n";
0304 return;
0305 }
0306
0307
0308
0309 if (abs(particle->pdgId()) == 11) {
0310 findFirstNonElectronMother(particle->mother(0), ancestorPID, ancestorStatus);
0311 } else {
0312 ancestorPID = particle->pdgId();
0313 ancestorStatus = particle->status();
0314 }
0315
0316 return;
0317 }
0318
0319
0320 DEFINE_FWK_MODULE(ElectronIDValidationAnalyzer);