Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:26:36

0001 // -*- C++ -*-
0002 //
0003 // Package:    GsfGEDElectronAnalyzer
0004 // Class:      GsfGEDElectronAnalyzer
0005 //
0006 /**\class GsfGEDElectronAnalyzer
0007 
0008  Description: <one line class summary>
0009 
0010  Implementation:
0011      <Notes on implementation>
0012 */
0013 //
0014 // Original Author:  Daniele Benedetti
0015 
0016 // system include files
0017 #include <memory>
0018 
0019 // user include files
0020 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0021 #include "FWCore/Framework/interface/Event.h"
0022 #include "FWCore/Framework/interface/MakerMacros.h"
0023 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0024 #include "FWCore/ServiceRegistry/interface/Service.h"
0025 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0026 #include "DataFormats/ParticleFlowReco/interface/GsfPFRecTrack.h"
0027 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0028 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0029 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0030 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0031 #include "DataFormats/EgammaReco/interface/ElectronSeed.h"
0032 #include "DataFormats/METReco/interface/PFMET.h"
0033 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0034 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0035 #include "DataFormats/Math/interface/normalizedPhi.h"
0036 
0037 #include <vector>
0038 #include <TROOT.h>
0039 #include <TFile.h>
0040 #include <TTree.h>
0041 #include <TH1F.h>
0042 #include <TH2F.h>
0043 #include <TLorentzVector.h>
0044 //
0045 // class decleration
0046 //
0047 
0048 using namespace edm;
0049 using namespace reco;
0050 using namespace std;
0051 class GsfGEDElectronAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0052 public:
0053   explicit GsfGEDElectronAnalyzer(const edm::ParameterSet &);
0054   ~GsfGEDElectronAnalyzer() override;
0055 
0056 private:
0057   void beginJob() override;
0058   void analyze(const edm::Event &, const edm::EventSetup &) override;
0059   void endJob() override;
0060 
0061   ParameterSet conf_;
0062 
0063   unsigned int ev;
0064   // ----------member data ---------------------------
0065 
0066   TH1F *h_etamc_ele, *h_etaec_ele, *h_etaecTrDr_ele, *h_etaecEcDr_ele, *h_etaeg_ele, *h_etaegTrDr_ele, *h_etaegEcDr_ele,
0067       *h_etaged_ele, *h_etagedTrDr_ele, *h_etagedEcDr_ele;
0068 
0069   TH1F *h_ptmc_ele, *h_ptec_ele, *h_ptecEcDr_ele, *h_ptecTrDr_ele, *h_pteg_ele, *h_ptegTrDr_ele, *h_ptegEcDr_ele,
0070       *h_ptged_ele, *h_ptgedTrDr_ele, *h_ptgedEcDr_ele;
0071 
0072   TH1F *h_mva_ele;
0073   TH2F *etaphi_nonreco;
0074 
0075   TH1F *eg_EEcalEtrue_1, *eg_EEcalEtrue_2, *eg_EEcalEtrue_3, *eg_EEcalEtrue_4, *eg_EEcalEtrue_5;
0076   TH1F *pf_EEcalEtrue_1, *pf_EEcalEtrue_2, *pf_EEcalEtrue_3, *pf_EEcalEtrue_4, *pf_EEcalEtrue_5;
0077   TH1F *ged_EEcalEtrue_1, *ged_EEcalEtrue_2, *ged_EEcalEtrue_3, *ged_EEcalEtrue_4, *ged_EEcalEtrue_5;
0078 
0079   TH1F *eg_e1x5_all, *eg_e1x5_eb, *eg_e1x5_ee, *ged_e1x5_all, *ged_e1x5_eb, *ged_e1x5_ee, *pf_e1x5_all, *pf_e1x5_eb,
0080       *pf_e1x5_ee;
0081 
0082   TH1F *eg_sihih_all, *eg_sihih_eb, *eg_sihih_ee, *ged_sihih_all, *ged_sihih_eb, *ged_sihih_ee, *pf_sihih_all,
0083       *pf_sihih_eb, *pf_sihih_ee;
0084 
0085   TH1F *eg_r9_all, *eg_r9_eb, *eg_r9_ee, *ged_r9_all, *ged_r9_eb, *ged_r9_ee, *pf_r9_all, *pf_r9_eb, *pf_r9_ee;
0086 
0087   TH1F *eg_scshh_all, *eg_scshh_eb, *eg_scshh_ee, *ged_scshh_all, *ged_scshh_eb, *ged_scshh_ee, *pf_scshh_all,
0088       *pf_scshh_eb, *pf_scshh_ee;
0089 
0090   TH1F *eg_scsff_all, *eg_scsff_eb, *eg_scsff_ee, *ged_scsff_all, *ged_scsff_eb, *ged_scsff_ee, *pf_scsff_all,
0091       *pf_scsff_eb, *pf_scsff_ee;
0092 
0093   TH1F *pf_MET_reco, *pf_MET_rereco;
0094 };
0095 
0096 //
0097 // constants, enums and typedefs
0098 //
0099 
0100 //
0101 // static data member definitions
0102 //
0103 
0104 //
0105 // constructors and destructor
0106 //
0107 GsfGEDElectronAnalyzer::GsfGEDElectronAnalyzer(const edm::ParameterSet &iConfig)
0108     : conf_(iConfig)
0109 
0110 {
0111   usesResource(TFileService::kSharedResource);
0112 
0113   edm::Service<TFileService> fs;
0114 
0115   // Efficiency plots
0116 
0117   h_etamc_ele = fs->make<TH1F>("h_etamc_ele", " ", 50, -2.5, 2.5);
0118   h_etaec_ele = fs->make<TH1F>("h_etaec_ele", " ", 50, -2.5, 2.5);
0119   h_etaecTrDr_ele = fs->make<TH1F>("h_etaecTrDr_ele", " ", 50, -2.5, 2.5);
0120   h_etaecEcDr_ele = fs->make<TH1F>("h_etaecEcDr_ele", " ", 50, -2.5, 2.5);
0121   h_etaeg_ele = fs->make<TH1F>("h_etaeg_ele", " ", 50, -2.5, 2.5);
0122   h_etaegTrDr_ele = fs->make<TH1F>("h_etaegTrDr_ele", " ", 50, -2.5, 2.5);
0123   h_etaegEcDr_ele = fs->make<TH1F>("h_etaegEcDr_ele", " ", 50, -2.5, 2.5);
0124   h_etaged_ele = fs->make<TH1F>("h_etaged_ele", " ", 50, -2.5, 2.5);
0125   h_etagedTrDr_ele = fs->make<TH1F>("h_etagedTrDr_ele", " ", 50, -2.5, 2.5);
0126   h_etagedEcDr_ele = fs->make<TH1F>("h_etagedEcDr_ele", " ", 50, -2.5, 2.5);
0127 
0128   h_ptmc_ele = fs->make<TH1F>("h_ptmc_ele", " ", 50, 0, 100.);
0129   h_ptec_ele = fs->make<TH1F>("h_ptec_ele", " ", 50, 0, 100.);
0130   h_ptecTrDr_ele = fs->make<TH1F>("h_ptecTrDr_ele", " ", 50, 0, 100.);
0131   h_ptecEcDr_ele = fs->make<TH1F>("h_ptecEcDr_ele", " ", 50, 0, 100.);
0132   h_pteg_ele = fs->make<TH1F>("h_pteg_ele", " ", 50, 0, 100.);
0133   h_ptegTrDr_ele = fs->make<TH1F>("h_ptegTrDr_ele", " ", 50, 0, 100.);
0134   h_ptegEcDr_ele = fs->make<TH1F>("h_ptegEcDr_ele", " ", 50, 0, 100.);
0135   h_ptged_ele = fs->make<TH1F>("h_ptged_ele", " ", 50, 0, 100.);
0136   h_ptgedTrDr_ele = fs->make<TH1F>("h_ptgedTrDr_ele", " ", 50, 0, 100.);
0137   h_ptgedEcDr_ele = fs->make<TH1F>("h_ptgedEcDr_ele", " ", 50, 0, 100.);
0138 
0139   h_mva_ele = fs->make<TH1F>("h_mva_ele", " ", 50, -1.1, 1.1);
0140 
0141   // Ecal map
0142   etaphi_nonreco = fs->make<TH2F>("etaphi_nonreco", " ", 50, -2.5, 2.5, 50, -3.2, 3.2);
0143 
0144   // Energies
0145   eg_EEcalEtrue_1 = fs->make<TH1F>("eg_EEcalEtrue_1", "  ", 50, 0., 2.);
0146   eg_EEcalEtrue_2 = fs->make<TH1F>("eg_EEcalEtrue_2", "  ", 50, 0., 2.);
0147   eg_EEcalEtrue_3 = fs->make<TH1F>("eg_EEcalEtrue_3", "  ", 50, 0., 2.);
0148   eg_EEcalEtrue_4 = fs->make<TH1F>("eg_EEcalEtrue_4", "  ", 50, 0., 2.);
0149   eg_EEcalEtrue_5 = fs->make<TH1F>("eg_EEcalEtrue_5", "  ", 50, 0., 2.);
0150 
0151   pf_EEcalEtrue_1 = fs->make<TH1F>("pf_EEcalEtrue_1", "  ", 50, 0., 2.);
0152   pf_EEcalEtrue_2 = fs->make<TH1F>("pf_EEcalEtrue_2", "  ", 50, 0., 2.);
0153   pf_EEcalEtrue_3 = fs->make<TH1F>("pf_EEcalEtrue_3", "  ", 50, 0., 2.);
0154   pf_EEcalEtrue_4 = fs->make<TH1F>("pf_EEcalEtrue_4", "  ", 50, 0., 2.);
0155   pf_EEcalEtrue_5 = fs->make<TH1F>("pf_EEcalEtrue_5", "  ", 50, 0., 2.);
0156 
0157   ged_EEcalEtrue_1 = fs->make<TH1F>("ged_EEcalEtrue_1", "  ", 50, 0., 2.);
0158   ged_EEcalEtrue_2 = fs->make<TH1F>("ged_EEcalEtrue_2", "  ", 50, 0., 2.);
0159   ged_EEcalEtrue_3 = fs->make<TH1F>("ged_EEcalEtrue_3", "  ", 50, 0., 2.);
0160   ged_EEcalEtrue_4 = fs->make<TH1F>("ged_EEcalEtrue_4", "  ", 50, 0., 2.);
0161   ged_EEcalEtrue_5 = fs->make<TH1F>("ged_EEcalEtrue_5", "  ", 50, 0., 2.);
0162 
0163   //shower shapes
0164   eg_e1x5_all = fs->make<TH1F>("eg_e1x5_all", "  ", 60, 0., 300.);
0165   eg_e1x5_eb = fs->make<TH1F>("eg_e1x5_eb", "  ", 60, 0., 300.);
0166   eg_e1x5_ee = fs->make<TH1F>("eg_e1x5_ee", "  ", 60, 0., 300.);
0167 
0168   ged_e1x5_all = fs->make<TH1F>("ged_e1x5_all", "  ", 60, 0., 300.);
0169   ged_e1x5_eb = fs->make<TH1F>("ged_e1x5_eb", "  ", 60, 0., 300.);
0170   ged_e1x5_ee = fs->make<TH1F>("ged_e1x5_ee", "  ", 60, 0., 300.);
0171 
0172   pf_e1x5_all = fs->make<TH1F>("pf_e1x5_all", "  ", 60, 0., 300.);
0173   pf_e1x5_eb = fs->make<TH1F>("pf_e1x5_eb", "  ", 60, 0., 300.);
0174   pf_e1x5_ee = fs->make<TH1F>("pf_e1x5_ee", "  ", 60, 0., 300.);
0175 
0176   eg_sihih_all = fs->make<TH1F>("eg_sihih_all", "  ", 100, 0.0, 0.05);
0177   eg_sihih_eb = fs->make<TH1F>("eg_sihih_eb", "  ", 100, 0.0, 0.05);
0178   eg_sihih_ee = fs->make<TH1F>("eg_sihih_ee", "  ", 100, 0.0, 0.05);
0179 
0180   ged_sihih_all = fs->make<TH1F>("ged_sihih_all", "  ", 100, 0.0, 0.05);
0181   ged_sihih_eb = fs->make<TH1F>("ged_sihih_eb", "  ", 100, 0.0, 0.05);
0182   ged_sihih_ee = fs->make<TH1F>("ged_sihih_ee", "  ", 100, 0.0, 0.05);
0183 
0184   pf_sihih_all = fs->make<TH1F>("pf_sihih_all", "  ", 100, 0.0, 0.05);
0185   pf_sihih_eb = fs->make<TH1F>("pf_sihih_eb", "  ", 100, 0.0, 0.05);
0186   pf_sihih_ee = fs->make<TH1F>("pf_sihih_ee", "  ", 100, 0.0, 0.05);
0187 
0188   eg_r9_all = fs->make<TH1F>("eg_r9_all", "  ", 200, 0.0, 2.0);
0189   eg_r9_eb = fs->make<TH1F>("eg_r9_eb", "  ", 200, 0.0, 2.0);
0190   eg_r9_ee = fs->make<TH1F>("eg_r9_ee", "  ", 200, 0.0, 2.0);
0191 
0192   ged_r9_all = fs->make<TH1F>("ged_r9_all", "  ", 200, 0.0, 2.0);
0193   ged_r9_eb = fs->make<TH1F>("ged_r9_eb", "  ", 200, 0.0, 2.0);
0194   ged_r9_ee = fs->make<TH1F>("ged_r9_ee", "  ", 200, 0.0, 2.0);
0195 
0196   pf_r9_all = fs->make<TH1F>("pf_r9_all", "  ", 200, 0.0, 2.0);
0197   pf_r9_eb = fs->make<TH1F>("pf_r9_eb", "  ", 200, 0.0, 2.0);
0198   pf_r9_ee = fs->make<TH1F>("pf_r9_ee", "  ", 200, 0.0, 2.0);
0199 
0200   eg_scshh_all = fs->make<TH1F>("eg_scshh_all", "  ", 100, 0.0, 0.06);
0201   eg_scshh_eb = fs->make<TH1F>("eg_scshh_eb", "  ", 100, 0.0, 0.06);
0202   eg_scshh_ee = fs->make<TH1F>("eg_scshh_ee", "  ", 100, 0.0, 0.06);
0203 
0204   ged_scshh_all = fs->make<TH1F>("ged_scshh_all", "  ", 100, 0.0, 0.06);
0205   ged_scshh_eb = fs->make<TH1F>("ged_scshh_eb", "  ", 100, 0.0, 0.06);
0206   ged_scshh_ee = fs->make<TH1F>("ged_scshh_ee", "  ", 100, 0.0, 0.06);
0207 
0208   pf_scshh_all = fs->make<TH1F>("pf_scshh_all", "  ", 100, 0.0, 0.06);
0209   pf_scshh_eb = fs->make<TH1F>("pf_scshh_eb", "  ", 100, 0.0, 0.06);
0210   pf_scshh_ee = fs->make<TH1F>("pf_scshh_ee", "  ", 100, 0.0, 0.06);
0211 
0212   eg_scsff_all = fs->make<TH1F>("eg_scsff_all", "  ", 150, 0.0, 0.15);
0213   eg_scsff_eb = fs->make<TH1F>("eg_scsff_eb", "  ", 150, 0.0, 0.15);
0214   eg_scsff_ee = fs->make<TH1F>("eg_scsff_ee", "  ", 150, 0.0, 0.15);
0215 
0216   ged_scsff_all = fs->make<TH1F>("ged_scsff_all", "  ", 150, 0.0, 0.15);
0217   ged_scsff_eb = fs->make<TH1F>("ged_scsff_eb", "  ", 150, 0.0, 0.15);
0218   ged_scsff_ee = fs->make<TH1F>("ged_scsff_ee", "  ", 150, 0.0, 0.15);
0219 
0220   pf_scsff_all = fs->make<TH1F>("pf_scsff_all", "  ", 150, 0.0, 0.15);
0221   pf_scsff_eb = fs->make<TH1F>("pf_scsff_eb", "  ", 150, 0.0, 0.15);
0222   pf_scsff_ee = fs->make<TH1F>("pf_scsff_ee", "  ", 150, 0.0, 0.15);
0223 
0224   // MET
0225   pf_MET_reco = fs->make<TH1F>("pf_MET_reco", "  ", 10, 0, 1000);
0226   pf_MET_rereco = fs->make<TH1F>("pf_MET_rereco", "  ", 10, 0, 1000);
0227 }
0228 
0229 GsfGEDElectronAnalyzer::~GsfGEDElectronAnalyzer() {
0230   // do anything here that needs to be done at desctruction time
0231   // (e.g. close files, deallocate resources etc.)
0232 }
0233 
0234 //
0235 // member functions
0236 //
0237 
0238 // ------------ method called to for each event  ------------
0239 void GsfGEDElectronAnalyzer::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0240   // Candidate info
0241   Handle<reco::PFCandidateCollection> collection;
0242   InputTag label("particleFlow::reRECO");  // <- Special electron coll.
0243   iEvent.getByLabel(label, collection);
0244   std::vector<reco::PFCandidate> candidates = (*collection.product());
0245 
0246   InputTag egammaLabel(string("gsfElectrons"));
0247   Handle<GsfElectronCollection> theGsfEleCollection;
0248   iEvent.getByLabel(egammaLabel, theGsfEleCollection);
0249   const GsfElectronCollection theGsfEle = *(theGsfEleCollection.product());
0250 
0251   InputTag MCTruthCollection(string("genParticles"));
0252   edm::Handle<std::vector<reco::GenParticle> > pMCTruth;
0253   iEvent.getByLabel(MCTruthCollection, pMCTruth);
0254 
0255   InputTag gedEleLabel(string("gedGsfElectrons::reRECO"));
0256   Handle<GsfElectronCollection> theGedEleCollection;
0257   iEvent.getByLabel(gedEleLabel, theGedEleCollection);
0258   const GsfElectronCollection theGedEle = *(theGedEleCollection.product());
0259 
0260   //get and plot the reco met
0261   InputTag pfMETRecoLabel("pfMet::RECO");
0262   Handle<std::vector<reco::PFMET> > recoMet;
0263   iEvent.getByLabel(pfMETRecoLabel, recoMet);
0264   pf_MET_reco->Fill(recoMet->at(0).et());
0265 
0266   //get and plot the rereco met
0267   InputTag pfMETReRecoLabel("pfMet::reRECO");
0268   Handle<std::vector<reco::PFMET> > rerecoMet;
0269   iEvent.getByLabel(pfMETReRecoLabel, rerecoMet);
0270   pf_MET_rereco->Fill(rerecoMet->at(0).et());
0271 
0272   bool debug = true;
0273 
0274   ev++;
0275 
0276   if (debug)
0277     cout << "************************* New Event:: " << ev << " *************************" << endl;
0278 
0279   // Validation from generator events
0280 
0281   for (std::vector<reco::GenParticle>::const_iterator cP = pMCTruth->begin(); cP != pMCTruth->end(); cP++) {
0282     float etamc = cP->eta();
0283     float phimc = cP->phi();
0284     float ptmc = cP->pt();
0285     float Emc = cP->energy();
0286 
0287     if (abs(cP->pdgId()) == 11 && cP->status() == 1 && cP->pt() > 2. && fabs(cP->eta()) < 2.5) {
0288       h_etamc_ele->Fill(etamc);
0289       h_ptmc_ele->Fill(ptmc);
0290 
0291       float MindR = 1000.;
0292       float mvaCutEle = -1.;
0293       bool isPfTrDr = false;
0294       bool isPfEcDr = false;
0295 
0296       if (debug)
0297         cout << " MC particle:  pt " << ptmc << " eta,phi " << etamc << ", " << phimc << endl;
0298 
0299       std::vector<reco::PFCandidate>::iterator it;
0300       for (it = candidates.begin(); it != candidates.end(); ++it) {
0301         reco::PFCandidate::ParticleType type = (*it).particleId();
0302 
0303         if (type == reco::PFCandidate::e) {
0304           float eta = (*it).eta();
0305           float phi = (*it).phi();
0306           float pfmva = (*it).mva_e_pi();
0307 
0308           reco::GsfTrackRef refGsf = (*it).gsfTrackRef();
0309           //ElectronSeedRef seedRef= refGsf->extra()->seedRef().castTo<ElectronSeedRef>();
0310 
0311           float deta = etamc - eta;
0312           float dphi = normalizedPhi(phimc - phi);
0313           float dR = sqrt(deta * deta + dphi * dphi);
0314 
0315           if (dR < 0.05) {
0316             MindR = dR;
0317             mvaCutEle = (*it).mva_e_pi();
0318             /*
0319         if(seedRef->isEcalDriven())
0320           isPfEcDr = true;
0321         if(seedRef->isTrackerDriven())
0322           isPfTrDr = true;
0323         */
0324 
0325             if (debug)
0326               cout << " PF ele matched:  pt " << (*it).pt() << " (" << (*it).ecalEnergy() / std::cosh(eta) << ") "
0327                    << " eta,phi " << eta << ", " << phi << " pfmva " << pfmva << endl;
0328             // all for the moment
0329           }
0330         }  // End PFCandidates Electron Selection
0331       }    // End Loop PFCandidates
0332 
0333       if (MindR < 0.05) {
0334         h_mva_ele->Fill(mvaCutEle);
0335         h_etaec_ele->Fill(etamc);
0336         h_ptec_ele->Fill(ptmc);
0337 
0338         if (isPfEcDr) {
0339           h_etaecEcDr_ele->Fill(etamc);
0340           h_ptecEcDr_ele->Fill(ptmc);
0341         }
0342         if (isPfTrDr) {
0343           h_etaecTrDr_ele->Fill(etamc);
0344           h_ptecTrDr_ele->Fill(ptmc);
0345         }
0346       } else {
0347         etaphi_nonreco->Fill(etamc, phimc);
0348       }
0349 
0350       float MindREG = 1000;
0351       bool isEgTrDr = false;
0352       bool isEgEcDr = false;
0353 
0354       for (uint j = 0; j < theGsfEle.size(); j++) {
0355         reco::GsfTrackRef egGsfTrackRef = (theGsfEle[j]).gsfTrack();
0356         float etareco = theGsfEle[j].eta();
0357         float phireco = theGsfEle[j].phi();
0358 
0359         float pfmva = theGsfEle[j].mva_e_pi();
0360 
0361         reco::GsfTrackRef refGsf = theGsfEle[j].gsfTrack();
0362         //ElectronSeedRef seedRef= refGsf->extra()->seedRef().castTo<ElectronSeedRef>();
0363 
0364         float deta = etamc - etareco;
0365         float dphi = normalizedPhi(phimc - phireco);
0366         float dR = sqrt(deta * deta + dphi * dphi);
0367 
0368         float SCEnergy = (theGsfEle[j]).superCluster()->energy();
0369         float ErecoEtrue = SCEnergy / Emc;
0370 
0371         if (dR < 0.05) {
0372           if (debug)
0373             cout << " EG ele matched: pt " << theGsfEle[j].pt() << " (" << SCEnergy / std::cosh(etareco) << ") "
0374                  << " eta,phi " << etareco << ", " << phireco << " pfmva " << pfmva << endl;
0375 
0376           if (fabs(etamc) < 0.5) {
0377             eg_EEcalEtrue_1->Fill(ErecoEtrue);
0378           }
0379           if (fabs(etamc) >= 0.5 && fabs(etamc) < 1.0) {
0380             eg_EEcalEtrue_2->Fill(ErecoEtrue);
0381           }
0382           if (fabs(etamc) >= 1.0 && fabs(etamc) < 1.5) {
0383             eg_EEcalEtrue_3->Fill(ErecoEtrue);
0384           }
0385           if (fabs(etamc) >= 1.5 && fabs(etamc) < 2.0) {
0386             eg_EEcalEtrue_4->Fill(ErecoEtrue);
0387           }
0388           if (fabs(etamc) >= 2.0 && fabs(etamc) < 2.5) {
0389             eg_EEcalEtrue_5->Fill(ErecoEtrue);
0390           }
0391 
0392           if ((theGsfEle[j].parentSuperCluster()).isNonnull()) {
0393             float SCPF = (theGsfEle[j]).parentSuperCluster()->rawEnergy();
0394             float EpfEtrue = SCPF / Emc;
0395             if (fabs(etamc) < 0.5) {
0396               pf_EEcalEtrue_1->Fill(EpfEtrue);
0397             }
0398             if (fabs(etamc) >= 0.5 && fabs(etamc) < 1.0) {
0399               pf_EEcalEtrue_2->Fill(EpfEtrue);
0400             }
0401             if (fabs(etamc) >= 1.0 && fabs(etamc) < 1.5) {
0402               pf_EEcalEtrue_3->Fill(EpfEtrue);
0403             }
0404             if (fabs(etamc) >= 1.5 && fabs(etamc) < 2.0) {
0405               pf_EEcalEtrue_4->Fill(EpfEtrue);
0406             }
0407             if (fabs(etamc) >= 2.0 && fabs(etamc) < 2.5) {
0408               pf_EEcalEtrue_5->Fill(EpfEtrue);
0409             }
0410             const reco::GsfElectron::ShowerShape &pfshapes = theGsfEle[j].showerShape();
0411 
0412             pf_e1x5_all->Fill(pfshapes.e1x5);
0413             pf_sihih_all->Fill(pfshapes.sigmaIetaIeta);
0414             pf_r9_all->Fill(pfshapes.r9);
0415             pf_scshh_all->Fill((theGsfEle[j]).parentSuperCluster()->etaWidth());
0416             pf_scsff_all->Fill((theGsfEle[j]).parentSuperCluster()->phiWidth());
0417             if (std::abs(etareco) < 1.479) {
0418               pf_e1x5_eb->Fill(pfshapes.e1x5);
0419               pf_sihih_eb->Fill(pfshapes.sigmaIetaIeta);
0420               pf_r9_eb->Fill(pfshapes.r9);
0421               pf_scshh_eb->Fill((theGsfEle[j]).parentSuperCluster()->etaWidth());
0422               pf_scsff_eb->Fill((theGsfEle[j]).parentSuperCluster()->phiWidth());
0423             }
0424             if (std::abs(etareco) >= 1.479) {
0425               pf_e1x5_ee->Fill(pfshapes.e1x5);
0426               pf_sihih_ee->Fill(pfshapes.sigmaIetaIeta);
0427               pf_r9_ee->Fill(pfshapes.r9);
0428               pf_scshh_ee->Fill((theGsfEle[j]).parentSuperCluster()->etaWidth());
0429               pf_scsff_ee->Fill((theGsfEle[j]).parentSuperCluster()->phiWidth());
0430             }
0431           }
0432 
0433           eg_e1x5_all->Fill(theGsfEle[j].e1x5());
0434           eg_sihih_all->Fill(theGsfEle[j].sigmaIetaIeta());
0435           eg_r9_all->Fill(theGsfEle[j].r9());
0436           eg_scshh_all->Fill((theGsfEle[j]).superCluster()->etaWidth());
0437           eg_scsff_all->Fill((theGsfEle[j]).superCluster()->phiWidth());
0438           if (std::abs(etareco) < 1.479) {
0439             eg_e1x5_eb->Fill(theGsfEle[j].e1x5());
0440             eg_sihih_eb->Fill(theGsfEle[j].sigmaIetaIeta());
0441             eg_r9_eb->Fill(theGsfEle[j].r9());
0442             eg_scshh_eb->Fill((theGsfEle[j]).superCluster()->etaWidth());
0443             eg_scsff_eb->Fill((theGsfEle[j]).superCluster()->phiWidth());
0444           }
0445           if (std::abs(etareco) >= 1.479) {
0446             eg_e1x5_ee->Fill(theGsfEle[j].e1x5());
0447             eg_sihih_ee->Fill(theGsfEle[j].sigmaIetaIeta());
0448             eg_r9_ee->Fill(theGsfEle[j].r9());
0449             eg_scshh_ee->Fill((theGsfEle[j]).superCluster()->etaWidth());
0450             eg_scsff_ee->Fill((theGsfEle[j]).superCluster()->phiWidth());
0451           }
0452 
0453           MindREG = dR;
0454           /*
0455       if(seedRef->isEcalDriven())
0456         isEgEcDr = true;
0457       if(seedRef->isTrackerDriven())
0458         isEgTrDr = true;
0459       */
0460         }
0461       }
0462       if (MindREG < 0.05) {
0463         h_etaeg_ele->Fill(etamc);
0464         h_pteg_ele->Fill(ptmc);
0465 
0466         // This is to understand the contribution of each seed type.
0467         if (isEgEcDr) {
0468           h_etaegEcDr_ele->Fill(etamc);
0469           h_ptegEcDr_ele->Fill(ptmc);
0470         }
0471         if (isEgTrDr) {
0472           h_etaegTrDr_ele->Fill(etamc);
0473           h_ptegTrDr_ele->Fill(ptmc);
0474         }
0475       }
0476 
0477       //Add loop on the GED electrons
0478 
0479       float MindRGedEg = 1000;
0480       bool isGedEgTrDr = false;
0481       bool isGedEgEcDr = false;
0482 
0483       for (uint j = 0; j < theGedEle.size(); j++) {
0484         reco::GsfTrackRef egGsfTrackRef = (theGedEle[j]).gsfTrack();
0485         float etareco = theGedEle[j].eta();
0486         float phireco = theGedEle[j].phi();
0487         float pfmva = theGedEle[j].mva_e_pi();
0488 
0489         reco::GsfTrackRef refGsf = theGedEle[j].gsfTrack();
0490         //ElectronSeedRef seedRef= refGsf->extra()->seedRef().castTo<ElectronSeedRef>();
0491 
0492         float deta = etamc - etareco;
0493         float dphi = normalizedPhi(phimc - phireco);
0494         float dR = sqrt(deta * deta + dphi * dphi);
0495 
0496         float SCEnergy = (theGedEle[j]).superCluster()->energy();
0497         float ErecoEtrue = SCEnergy / Emc;
0498 
0499         if (dR < 0.05) {
0500           const reco::PFCandidate *matchPF = NULL;
0501           if (debug)
0502             cout << " GED ele matched: pt " << theGedEle[j].pt() << " (" << SCEnergy / std::cosh(etareco) << ") "
0503                  << " eta,phi " << etareco << ", " << phireco << " pfmva " << pfmva << endl;
0504 
0505           for (unsigned k = 0; k < candidates.size(); ++k) {
0506             if (std::abs(candidates[k].pdgId()) == 11) {
0507               reco::GsfTrackRef gsfEleGsfTrackRef = candidates[k].gsfTrackRef();
0508               if (gsfEleGsfTrackRef->ptMode() == egGsfTrackRef->ptMode()) {
0509                 matchPF = &candidates[k];
0510               }
0511             }
0512           }
0513 
0514           if (debug && matchPF)
0515             std::cout << "GED-PF match!" << std::endl;
0516 
0517           if (ErecoEtrue < 0.6)
0518             std::cout << "++bad Ereco/Etrue: " << ErecoEtrue << std::endl;
0519 
0520           if (fabs(etamc) < 0.5) {
0521             ged_EEcalEtrue_1->Fill(ErecoEtrue);
0522           }
0523           if (fabs(etamc) >= 0.5 && fabs(etamc) < 1.0) {
0524             ged_EEcalEtrue_2->Fill(ErecoEtrue);
0525           }
0526           if (fabs(etamc) >= 1.0 && fabs(etamc) < 1.5) {
0527             ged_EEcalEtrue_3->Fill(ErecoEtrue);
0528           }
0529           if (fabs(etamc) >= 1.5 && fabs(etamc) < 2.0) {
0530             ged_EEcalEtrue_4->Fill(ErecoEtrue);
0531           }
0532           if (fabs(etamc) >= 2.0 && fabs(etamc) < 2.5) {
0533             ged_EEcalEtrue_5->Fill(ErecoEtrue);
0534             /*
0535         if( debug && matchPF && ErecoEtrue > 1.30 ) {
0536           cout << " -PF ele matched: pt " 
0537            << matchPF->pt()
0538          << " eta,phi " <<  matchPF->eta() 
0539          << ", " << matchPF->phi() << " pfmva " 
0540          <<  matchPF->mva_e_pi() << endl;
0541           std::cout << "GED Clusters: " << std::endl;
0542           for( const auto& c : theGedEle[j].superCluster()->clusters() ) {
0543         std::cout << " E/eta/phi : " << c->energy() 
0544               << '/' << c->eta() 
0545               << '/' << c->phi() << std::endl;      
0546           }
0547           std::cout << "-PF Clusters: " << std::endl;
0548           for( const auto& c : matchPF->superClusterRef()->clusters() ) {
0549         std::cout << " E/eta/phi : " << c->energy() 
0550               << '/' << c->eta() 
0551               << '/' << c->phi() << std::endl;
0552           }
0553           
0554         }
0555         */
0556           }
0557 
0558           ged_e1x5_all->Fill(theGedEle[j].e1x5());
0559           ged_sihih_all->Fill(theGedEle[j].sigmaIetaIeta());
0560           ged_r9_all->Fill(theGedEle[j].r9());
0561           ged_scshh_all->Fill((theGedEle[j]).superCluster()->etaWidth());
0562           ged_scsff_all->Fill((theGedEle[j]).superCluster()->phiWidth());
0563           if (std::abs(etareco) < 1.479) {
0564             ged_e1x5_eb->Fill(theGedEle[j].e1x5());
0565             ged_sihih_eb->Fill(theGedEle[j].sigmaIetaIeta());
0566             ged_r9_eb->Fill(theGedEle[j].r9());
0567             ged_scshh_eb->Fill((theGedEle[j]).superCluster()->etaWidth());
0568             ged_scsff_eb->Fill((theGedEle[j]).superCluster()->phiWidth());
0569           }
0570           if (std::abs(etareco) >= 1.479) {
0571             ged_e1x5_ee->Fill(theGedEle[j].e1x5());
0572             ged_sihih_ee->Fill(theGedEle[j].sigmaIetaIeta());
0573             ged_r9_ee->Fill(theGedEle[j].r9());
0574             ged_scshh_ee->Fill((theGedEle[j]).superCluster()->etaWidth());
0575             ged_scsff_eb->Fill((theGedEle[j]).superCluster()->phiWidth());
0576           }
0577 
0578           MindRGedEg = dR;
0579           /*
0580       if(seedRef->isEcalDriven())
0581         isGedEgEcDr = true;
0582       if(seedRef->isTrackerDriven())
0583         isGedEgTrDr = true;
0584       */
0585         }
0586       }
0587       if (MindRGedEg < 0.05) {
0588         h_etaged_ele->Fill(etamc);
0589         h_ptged_ele->Fill(ptmc);
0590 
0591         // This is to understand the contribution of each seed type.
0592         if (isGedEgEcDr) {
0593           h_etagedEcDr_ele->Fill(etamc);
0594           h_ptgedEcDr_ele->Fill(ptmc);
0595         }
0596         if (isGedEgTrDr) {
0597           h_etagedTrDr_ele->Fill(etamc);
0598           h_ptgedTrDr_ele->Fill(ptmc);
0599         }
0600       }  //End Loop Generator Particles
0601 
0602     }  //End IF Generator Particles
0603 
0604   }  //End Loop Generator Particles
0605 }
0606 // ------------ method called once each job just before starting event loop  ------------
0607 void GsfGEDElectronAnalyzer::beginJob() { ev = 0; }
0608 
0609 // ------------ method called once each job just after ending the event loop  ------------
0610 void GsfGEDElectronAnalyzer::endJob() { cout << " endJob:: #events " << ev << endl; }
0611 //define this as a plug-in
0612 DEFINE_FWK_MODULE(GsfGEDElectronAnalyzer);