EcalAlCaRecoProducers/src/AlCaElectronsTest.h

0001 // system include files
0002 #include <memory>
0003 #include <string>
0004 #include <iostream>
0005 #include <algorithm>
0006
0007 // user include files
0008 #include "FWCore/Framework/interface/Frameworkfwd.h"
0009 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/ServiceRegistry/interface/Service.h"
0016 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0017 #include "FWCore/Utilities/interface/EDGetToken.h"
0018 //#include "Calibration/EcalAlCaRecoProducers/interface/AlCaElectronsProducer.h"
0019 #include "DataFormats/EgammaCandidates/interface/SiStripElectron.h"
0020 #include "DataFormats/EgammaCandidates/interface/ElectronFwd.h"
0021 #include "DataFormats/EgammaCandidates/interface/Electron.h"
0022 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0023 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0024 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0025 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0026 #include "DataFormats/EcalDetId/interface/EEDetId.h"
0027 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
0028
0029 #include "TFile.h"
0030 #include "TH1.h"
0031 #include "TH2.h"
0032
0033 #include <Math/VectorUtil.h>
0034
0035 class AlCaElectronsTest : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0036 public:
0037   explicit AlCaElectronsTest(const edm::ParameterSet&);
0038   ~AlCaElectronsTest() = default;
0039   void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0040   void beginJob() override;
0041
0042 private:
0043   EcalRecHit getMaximum(const EcalRecHitCollection* recHits);
0044   void fillAroundBarrel(const EcalRecHitCollection* recHits, int eta, int phi);
0045   void fillAroundEndcap(const EcalRecHitCollection* recHits, int ics, int ips);
0046
0047   edm::EDGetTokenT<EBRecHitCollection> m_barrelAlCa;
0048   edm::EDGetTokenT<EERecHitCollection> m_endcapAlCa;
0049   std::string m_outputFileName;
0050
0051   //! ECAL map
0052   TH2F* m_barrelGlobalCrystalsMap;
0053   //! local map
0054   TH2F* m_barrelLocalCrystalsMap;
0055   //! ECAL map
0056   TH2F* m_endcapGlobalCrystalsMap;
0057   //! local map
0058   TH2F* m_endcapLocalCrystalsMap;
0059   //! ECAL Energy
0060   TH2F* m_barrelGlobalCrystalsEnergy;
0061   //! local Energy
0062   TH2F* m_barrelLocalCrystalsEnergy;
0063   //! ECAL Energy
0064   TH2F* m_endcapGlobalCrystalsEnergy;
0065   //! local Energy
0066   TH2F* m_endcapLocalCrystalsEnergy;
0067   //! ECAL EnergyMap
0068   TH2F* m_barrelGlobalCrystalsEnergyMap;
0069   //! ECAL EnergyMap
0070   TH2F* m_endcapGlobalCrystalsEnergyMap;
0071 };
0072
0073 // ----------------------------------------------------------------
0074
0075 AlCaElectronsTest::AlCaElectronsTest(const edm::ParameterSet& iConfig)
0076     : m_barrelAlCa(consumes<EBRecHitCollection>(iConfig.getParameter<edm::InputTag>("alcaBarrelHitCollection"))),
0077       m_endcapAlCa(consumes<EERecHitCollection>(iConfig.getParameter<edm::InputTag>("alcaEndcapHitCollection"))),
0078       m_outputFileName(
0079           iConfig.getUntrackedParameter<std::string>("HistOutFile", std::string("AlCaElectronsTest.root"))) {
0080   usesResource(TFileService::kSharedResource);
0081 }
0082
0083 // ----------------------------------------------------------------
0084
0085 void AlCaElectronsTest::beginJob() {
0086   edm::Service<TFileService> fs;
0087   m_barrelGlobalCrystalsMap =
0088       fs->make<TH2F>("m_barrelGlobalCrystalsMap", "m_barrelGlobalCrystalsMap", 171, -85, 86, 360, 0, 360);
0089   m_barrelLocalCrystalsMap =
0090       fs->make<TH2F>("m_barrelLocalCrystalsMap", "m_barrelLocalCrystalsMap", 20, -10, 10, 20, -10, 10);
0091   m_endcapGlobalCrystalsMap =
0092       fs->make<TH2F>("m_endcapGlobalCrystalsMap", "m_endcapGlobalCrystalsMap", 100, 0, 100, 100, 0, 100);
0093   m_endcapLocalCrystalsMap =
0094       fs->make<TH2F>("m_endcapLocalCrystalsMap", "m_endcapLocalCrystalsMap", 20, -10, 10, 20, -10, 10);
0095   m_barrelGlobalCrystalsEnergy =
0096       fs->make<TH2F>("m_barrelGlobalCrystalsEnergy", "m_barrelGlobalCrystalsEnergy", 171, -85, 86, 360, 0, 360);
0097   m_barrelLocalCrystalsEnergy =
0098       fs->make<TH2F>("m_barrelLocalCrystalsEnergy", "m_barrelLocalCrystalsEnergy", 20, -10, 10, 20, -10, 10);
0099   m_endcapGlobalCrystalsEnergy =
0100       fs->make<TH2F>("m_endcapGlobalCrystalsEnergy", "m_endcapGlobalCrystalsEnergy", 100, 0, 100, 100, 0, 100);
0101   m_endcapLocalCrystalsEnergy =
0102       fs->make<TH2F>("m_endcapLocalCrystalsEnergy", "m_endcapLocalCrystalsEnergy", 20, -10, 10, 20, -10, 10);
0103   m_barrelGlobalCrystalsEnergyMap =
0104       fs->make<TH2F>("m_barrelGlobalCrystalsEnergyMap", "m_barrelGlobalCrystalsEnergyMap", 171, -85, 86, 360, 0, 360);
0105   m_endcapGlobalCrystalsEnergyMap =
0106       fs->make<TH2F>("m_endcapGlobalCrystalsEnergyMap", "m_endcapGlobalCrystalsEnergyMap", 100, 0, 100, 100, 0, 100);
0107 }
0108
0109 // ----------------------------------------------------------------
0110
0111 void AlCaElectronsTest::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0112   //FIXME replace with msg logger
0113   edm::LogVerbatim("ElectronsTest") << "[AlCaElectronsTest] analysing event " << iEvent.id();
0114
0115   //PG get the collections
0116   // get Barrel RecHits
0117   edm::Handle<EBRecHitCollection> barrelRecHitsHandle;
0118   iEvent.getByToken(m_barrelAlCa, barrelRecHitsHandle);
0119   if (!barrelRecHitsHandle.isValid()) {
0120     edm::EDConsumerBase::Labels labels;
0121     labelsForToken(m_barrelAlCa, labels);
0122     edm::LogError("ElectronsTest") << "[AlCaElectronsTest] caught std::exception in rertieving " << labels.module;
0123     return;
0124   } else {
0125     const EBRecHitCollection* barrelHitsCollection = barrelRecHitsHandle.product();
0126     //PG fill the histo with the maximum
0127     EcalRecHit barrelMax = getMaximum(barrelHitsCollection);
0128     EBDetId barrelMaxId(barrelMax.id());
0129     m_barrelGlobalCrystalsMap->Fill(barrelMaxId.ieta(), barrelMaxId.iphi());
0130     m_barrelGlobalCrystalsEnergy->Fill(barrelMaxId.ieta(), barrelMaxId.iphi(), barrelMax.energy());
0131     fillAroundBarrel(barrelHitsCollection, barrelMaxId.ieta(), barrelMaxId.iphi());
0132   }
0133
0134   // get Endcap RecHits
0135   edm::Handle<EERecHitCollection> endcapRecHitsHandle;
0136   iEvent.getByToken(m_endcapAlCa, endcapRecHitsHandle);
0137   if (!endcapRecHitsHandle.isValid()) {
0138     edm::EDConsumerBase::Labels labels;
0139     labelsForToken(m_endcapAlCa, labels);
0140     edm::LogError("ElectronsTest") << "[AlCaElectronsTest] caught std::exception in rertieving " << labels.module;
0141     return;
0142   } else {
0143     const EERecHitCollection* endcapHitsCollection = endcapRecHitsHandle.product();
0144     //PG fill the histo with the maximum
0145     EcalRecHit endcapMax = getMaximum(endcapHitsCollection);
0146     EEDetId endcapMaxId(endcapMax.id());
0147     m_endcapGlobalCrystalsMap->Fill(endcapMaxId.ix(), endcapMaxId.iy());
0148     m_endcapGlobalCrystalsEnergy->Fill(endcapMaxId.ix(), endcapMaxId.iy(), endcapMax.energy());
0149     fillAroundEndcap(endcapHitsCollection, endcapMaxId.ix(), endcapMaxId.iy());
0150   }
0151 }
0152
0153 // ----------------------------------------------------------------
0154
0155 EcalRecHit AlCaElectronsTest::getMaximum(const EcalRecHitCollection* recHits) {
0156   double energy = 0.;
0157   EcalRecHit max;
0158   for (EcalRecHitCollection::const_iterator elem = recHits->begin(); elem != recHits->end(); ++elem) {
0159     if (elem->energy() > energy) {
0160       energy = elem->energy();
0161       max = *elem;
0162     }
0163   }
0164   return max;
0165 }
0166
0167 // ----------------------------------------------------------------
0168
0169 void AlCaElectronsTest::fillAroundBarrel(const EcalRecHitCollection* recHits, int eta, int phi) {
0170   for (EcalRecHitCollection::const_iterator elem = recHits->begin(); elem != recHits->end(); ++elem) {
0171     EBDetId elementId = elem->id();
0172     m_barrelLocalCrystalsMap->Fill(elementId.ieta() - eta, elementId.iphi() - phi);
0173     m_barrelLocalCrystalsEnergy->Fill(elementId.ieta() - eta, elementId.iphi() - phi, elem->energy());
0174     m_barrelGlobalCrystalsEnergyMap->Fill(elementId.ieta(), elementId.iphi(), elem->energy());
0175   }
0176   return;
0177 }
0178
0179 // ----------------------------------------------------------------
0180
0181 void AlCaElectronsTest::fillAroundEndcap(const EcalRecHitCollection* recHits, int ics, int ips) {
0182   for (EcalRecHitCollection::const_iterator elem = recHits->begin(); elem != recHits->end(); ++elem) {
0183     EEDetId elementId = elem->id();
0184     m_endcapLocalCrystalsMap->Fill(elementId.ix() - ics, elementId.iy() - ips);
0185     m_endcapLocalCrystalsEnergy->Fill(elementId.ix() - ics, elementId.iy() - ips, elem->energy());
0186     m_endcapGlobalCrystalsEnergyMap->Fill(elementId.ix(), elementId.iy(), elem->energy());
0187   }
0188   return;
0189 }