00001
00002 #include <memory>
00003 #include <string>
00004 #include <iostream>
00005 #include <algorithm>
00006
00007
00008 #include "FWCore/Framework/interface/Frameworkfwd.h"
00009 #include "FWCore/Framework/interface/EDAnalyzer.h"
00010 #include "FWCore/Framework/interface/Event.h"
00011 #include "FWCore/Framework/interface/EventSetup.h"
00012 #include "FWCore/Framework/interface/MakerMacros.h"
00013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00014 #include "Calibration/EcalAlCaRecoProducers/interface/AlCaElectronsProducer.h"
00015 #include "DataFormats/EgammaCandidates/interface/SiStripElectron.h"
00016 #include "DataFormats/EgammaCandidates/interface/ElectronFwd.h"
00017 #include "DataFormats/EgammaCandidates/interface/Electron.h"
00018 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00019 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
00020 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
00021 #include "DataFormats/EcalDetId/interface/EBDetId.h"
00022 #include "DataFormats/EcalDetId/interface/EEDetId.h"
00023 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
00024
00025
00026 #include "TFile.h"
00027 #include "TH1.h"
00028 #include "TH2.h"
00029
00030 #include <Math/VectorUtil.h>
00031
00032 using namespace reco ;
00033 using namespace edm ;
00034
00035
00036 class AlCaElectronsTest : public edm::EDAnalyzer {
00037 public:
00038 explicit AlCaElectronsTest (const edm::ParameterSet&) ;
00039 ~AlCaElectronsTest () {}
00040 virtual void analyze (const edm::Event& iEvent,
00041 const edm::EventSetup& iSetup) ;
00042 virtual void beginJob (const edm::EventSetup& iSetup) ;
00043 virtual void endJob () ;
00044
00045 private:
00046
00047 EcalRecHit getMaximum (const EcalRecHitCollection * recHits) ;
00048 void fillAroundBarrel (const EcalRecHitCollection * recHits, int eta, int phi) ;
00049 void fillAroundEndcap (const EcalRecHitCollection * recHits, int ics, int ips) ;
00050
00051
00052 private:
00053
00054 edm::InputTag m_barrelAlCa ;
00055 edm::InputTag m_endcapAlCa ;
00056 std::string m_outputFileName ;
00057
00059 TH2F * m_barrelGlobalCrystalsMap ;
00061 TH2F * m_barrelLocalCrystalsMap ;
00063 TH2F * m_endcapGlobalCrystalsMap ;
00065 TH2F * m_endcapLocalCrystalsMap ;
00067 TH2F * m_barrelGlobalCrystalsEnergy ;
00069 TH2F * m_barrelLocalCrystalsEnergy ;
00071 TH2F * m_endcapGlobalCrystalsEnergy ;
00073 TH2F * m_endcapLocalCrystalsEnergy ;
00075 TH2F * m_barrelGlobalCrystalsEnergyMap ;
00077 TH2F * m_endcapGlobalCrystalsEnergyMap ;
00078 } ;
00079
00080
00081
00082
00083
00084 AlCaElectronsTest::AlCaElectronsTest (const edm::ParameterSet& iConfig) :
00085 m_barrelAlCa (iConfig.getParameter<edm::InputTag> ("alcaBarrelHitCollection")) ,
00086 m_endcapAlCa (iConfig.getParameter<edm::InputTag> ("alcaEndcapHitCollection")) ,
00087 m_outputFileName (iConfig.getUntrackedParameter<std::string>
00088 ("HistOutFile",std::string ("AlCaElectronsTest.root")))
00089 {}
00090
00091
00092
00093
00094
00095 void
00096 AlCaElectronsTest::beginJob ( const edm::EventSetup& iSetup)
00097 {
00098 m_barrelGlobalCrystalsMap = new TH2F ("m_barrelGlobalCrystalsMap","m_barrelGlobalCrystalsMap",171,-85,86,360,0,360) ;
00099 m_barrelLocalCrystalsMap = new TH2F ("m_barrelLocalCrystalsMap","m_barrelLocalCrystalsMap",20,-10,10,20,-10,10) ;
00100 m_endcapGlobalCrystalsMap = new TH2F ("m_endcapGlobalCrystalsMap","m_endcapGlobalCrystalsMap",100,0,100,100,0,100) ;
00101 m_endcapLocalCrystalsMap = new TH2F ("m_endcapLocalCrystalsMap","m_endcapLocalCrystalsMap",20,-10,10,20,-10,10) ;
00102 m_barrelGlobalCrystalsEnergy = new TH2F ("m_barrelGlobalCrystalsEnergy","m_barrelGlobalCrystalsEnergy",171,-85,86,360,0,360) ;
00103 m_barrelLocalCrystalsEnergy = new TH2F ("m_barrelLocalCrystalsEnergy","m_barrelLocalCrystalsEnergy",20,-10,10,20,-10,10) ;
00104 m_endcapGlobalCrystalsEnergy = new TH2F ("m_endcapGlobalCrystalsEnergy","m_endcapGlobalCrystalsEnergy",100,0,100,100,0,100) ;
00105 m_endcapLocalCrystalsEnergy = new TH2F ("m_endcapLocalCrystalsEnergy","m_endcapLocalCrystalsEnergy",20,-10,10,20,-10,10) ;
00106 m_barrelGlobalCrystalsEnergyMap = new TH2F ("m_barrelGlobalCrystalsEnergyMap","m_barrelGlobalCrystalsEnergyMap",171,-85,86,360,0,360) ;
00107 m_endcapGlobalCrystalsEnergyMap = new TH2F ("m_endcapGlobalCrystalsEnergyMap","m_endcapGlobalCrystalsEnergyMap",100,0,100,100,0,100) ;
00108 return ;
00109 }
00110
00111
00112
00113
00114
00115 void
00116 AlCaElectronsTest::endJob ()
00117 {
00118 TFile output (m_outputFileName.c_str (),"recreate") ;
00119 m_barrelGlobalCrystalsMap->Write () ;
00120 m_barrelLocalCrystalsMap->Write () ;
00121 m_endcapGlobalCrystalsMap->Write () ;
00122 m_endcapLocalCrystalsMap->Write () ;
00123 m_barrelGlobalCrystalsEnergy->Write () ;
00124 m_barrelLocalCrystalsEnergy->Write () ;
00125 m_endcapGlobalCrystalsEnergy->Write () ;
00126 m_endcapLocalCrystalsEnergy->Write () ;
00127 m_barrelGlobalCrystalsEnergyMap->Write () ;
00128 m_endcapGlobalCrystalsEnergyMap->Write () ;
00129 output.Close () ;
00130
00131 return ;
00132 }
00133
00134
00135
00136
00137
00138 void
00139 AlCaElectronsTest::analyze (const edm::Event& iEvent,
00140 const edm::EventSetup& iSetup)
00141 {
00142
00143 std::cout << "[AlCaElectronsTest] analysing event "
00144 << iEvent.id () << std::endl ;
00145
00146
00147
00148 Handle<EBRecHitCollection> barrelRecHitsHandle ;
00149 iEvent.getByLabel (m_barrelAlCa, barrelRecHitsHandle) ;
00150 if (!barrelRecHitsHandle.isValid()) {
00151 std::cerr << "[AlCaElectronsTest] caught std::exception "
00152 << " in rertieving " << m_barrelAlCa
00153 << std::endl ;
00154 return ;
00155 } else {
00156 const EBRecHitCollection* barrelHitsCollection = barrelRecHitsHandle.product () ;
00157
00158 EcalRecHit barrelMax = getMaximum (barrelHitsCollection) ;
00159 EBDetId barrelMaxId (barrelMax.id ()) ;
00160 m_barrelGlobalCrystalsMap->Fill (
00161 barrelMaxId.ieta () ,
00162 barrelMaxId.iphi ()
00163 ) ;
00164 m_barrelGlobalCrystalsEnergy->Fill (
00165 barrelMaxId.ieta () ,
00166 barrelMaxId.iphi () ,
00167 barrelMax.energy ()
00168 ) ;
00169 fillAroundBarrel (
00170 barrelHitsCollection,
00171 barrelMaxId.ieta (),
00172 barrelMaxId.iphi ()
00173 ) ;
00174 }
00175
00176
00177 Handle<EERecHitCollection> endcapRecHitsHandle ;
00178 iEvent.getByLabel (m_endcapAlCa,endcapRecHitsHandle) ;
00179 if (!endcapRecHitsHandle.isValid()) {
00180 std::cerr << "[AlCaElectronsTest] caught std::exception "
00181 << " in rertieving " << m_endcapAlCa
00182 << std::endl ;
00183 return ;
00184 } else {
00185 const EERecHitCollection* endcapHitsCollection = endcapRecHitsHandle.product () ;
00186
00187 EcalRecHit endcapMax = getMaximum (endcapHitsCollection) ;
00188 EEDetId endcapMaxId (endcapMax.id ()) ;
00189 m_endcapGlobalCrystalsMap->Fill (
00190 endcapMaxId.ix () ,
00191 endcapMaxId.iy ()
00192 ) ;
00193 m_endcapGlobalCrystalsEnergy->Fill (
00194 endcapMaxId.ix () ,
00195 endcapMaxId.iy () ,
00196 endcapMax.energy ()
00197 ) ;
00198 fillAroundEndcap (
00199 endcapHitsCollection,
00200 endcapMaxId.ix (),
00201 endcapMaxId.iy ()
00202 ) ;
00203 }
00204 }
00205
00206
00207
00208
00209
00210 EcalRecHit
00211 AlCaElectronsTest::getMaximum (const EcalRecHitCollection * recHits)
00212 {
00213 double energy = 0. ;
00214 EcalRecHit max ;
00215 for (EcalRecHitCollection::const_iterator elem = recHits->begin () ;
00216 elem != recHits->end () ;
00217 ++elem)
00218 {
00219 if (elem->energy () > energy)
00220 {
00221 energy = elem->energy () ;
00222 max = *elem ;
00223 }
00224 }
00225 return max ;
00226 }
00227
00228
00229
00230
00231
00232 void
00233 AlCaElectronsTest::fillAroundBarrel (const EcalRecHitCollection * recHits, int eta, int phi)
00234 {
00235 for (EcalRecHitCollection::const_iterator elem = recHits->begin () ;
00236 elem != recHits->end () ;
00237 ++elem)
00238 {
00239 EBDetId elementId = elem->id () ;
00240 m_barrelLocalCrystalsMap->Fill (
00241 elementId.ieta () - eta ,
00242 elementId.iphi () - phi
00243 ) ;
00244 m_barrelLocalCrystalsEnergy->Fill (
00245 elementId.ieta () - eta ,
00246 elementId.iphi () - phi ,
00247 elem->energy ()
00248 ) ;
00249 m_barrelGlobalCrystalsEnergyMap->Fill (
00250 elementId.ieta () ,
00251 elementId.iphi () ,
00252 elem->energy ()
00253 ) ;
00254
00255 }
00256 return ;
00257 }
00258
00259
00260
00261
00262
00263 void
00264 AlCaElectronsTest::fillAroundEndcap (const EcalRecHitCollection * recHits, int ics, int ips)
00265 {
00266 for (EcalRecHitCollection::const_iterator elem = recHits->begin () ;
00267 elem != recHits->end () ;
00268 ++elem)
00269 {
00270 EEDetId elementId = elem->id () ;
00271 m_endcapLocalCrystalsMap->Fill (
00272 elementId.ix () - ics ,
00273 elementId.iy () - ips
00274 ) ;
00275 m_endcapLocalCrystalsEnergy->Fill (
00276 elementId.ix () - ics ,
00277 elementId.iy () - ips ,
00278 elem->energy ()
00279 ) ;
00280 m_endcapGlobalCrystalsEnergyMap->Fill (
00281 elementId.ix () ,
00282 elementId.iy () ,
00283 elem->energy ()
00284 ) ;
00285 }
00286 return ;
00287 }