![]() |
![]() |
#include <Validation/RecoMET/interface/DumpEvent.h>
Definition at line 83 of file DumpEvent.h.
DumpEvent::DumpEvent | ( | const edm::ParameterSet & | iConfig | ) |
Definition at line 7 of file DumpEvent.cc.
References debug_, EnergyThreshold, FirstEvent, edm::ParameterSet::getParameter(), LastEvent, and theEvent.
00008 { 00009 EnergyThreshold = iConfig.getParameter<double>("EnergyThreshold"); 00010 theEvent = iConfig.getParameter<int>("theEvent"); 00011 FirstEvent = iConfig.getParameter<int>("FirstEvent"); 00012 LastEvent = iConfig.getParameter<int>("LastEvent"); 00013 debug_ = iConfig.getParameter<bool>("Debug"); 00014 }
void DumpEvent::analyze | ( | const edm::Event & | iEvent, | |
const edm::EventSetup & | iSetup | |||
) | [virtual] |
Implements edm::EDAnalyzer.
Definition at line 61 of file DumpEvent.cc.
References CurrentEvent, DEBUG, theEvent, WriteElectrons(), WriteJets(), WriteMET(), WriteMuons(), WritePhotons(), and WriteSCs().
00062 { 00063 CurrentEvent++; 00064 DEBUG( "Event: " << CurrentEvent); 00065 if (CurrentEvent == theEvent) 00066 { 00067 WritePhotons(iEvent, iSetup); 00068 WriteElectrons(iEvent, iSetup); 00069 WriteJets(iEvent, iSetup); 00070 WriteMuons(iEvent, iSetup); 00071 WriteMET(iEvent, iSetup); 00072 WriteSCs(iEvent, iSetup); 00073 } 00074 }
void DumpEvent::beginJob | ( | const edm::EventSetup & | iSetup | ) | [virtual] |
Reimplemented from edm::EDAnalyzer.
Definition at line 22 of file DumpEvent.cc.
References BookHistos(), CurrentEvent, and m_DataFile.
00022 { 00023 CurrentEvent = -1; 00024 // Make the output files 00025 m_DataFile=new TFile("DumpEvent_data.root" ,"RECREATE"); 00026 // Book the Histograms 00027 BookHistos(); 00028 }
void DumpEvent::BookHistos | ( | ) |
Definition at line 30 of file DumpEvent.cc.
References hCaloTowerToJetMap_ieta_iphi, hEB_SC_ieta_iphi, hEEmZ_SC_ix_iy, hEEpZ_SC_ix_iy, hElectron_energy, hElectron_eta, hElectron_phi, hGenMET_phi, hJet_energy, hJet_eta, hJet_phi, hMuon_eta, hMuon_phi, hMuon_pt, hPhoton_energy, hPhoton_eta, hPhoton_phi, hRecoMET_MET, hRecoMET_phi, and m_DataFile.
Referenced by beginJob().
00031 { 00032 // Book Data Histograms 00033 m_DataFile->cd(); 00034 hElectron_eta = new TH1F("hElectron_eta", "", 100,0,100); 00035 hElectron_phi = new TH1F("hElectron_phi", "", 100,0,100); 00036 hElectron_energy = new TH1F("hElectron_energy", "", 100,0,100); 00037 00038 hPhoton_eta = new TH1F("hPhoton_eta", "", 100,0,100); 00039 hPhoton_phi = new TH1F("hPhoton_phi", "", 100,0,100); 00040 hPhoton_energy = new TH1F("hPhoton_energy", "", 100,0,100); 00041 00042 hJet_eta = new TH1F("hJet_eta", "", 100,0,100); 00043 hJet_phi = new TH1F("hJet_phi", "", 100,0,100); 00044 hJet_energy = new TH1F("hJet_energy", "", 100,0,100); 00045 hCaloTowerToJetMap_ieta_iphi = new TH2F("hCaloTowerToJetMap_ieta_iphi","",83,-41,42, 73,0,73); 00046 00047 hMuon_eta = new TH1F("hMuon_eta", "", 100,0,100); 00048 hMuon_phi = new TH1F("hMuon_phi", "", 100,0,100); 00049 hMuon_pt = new TH1F("hMuon_pt", "", 100,0,100); 00050 00051 hRecoMET_phi = new TH1F("hRecoMET_phi", "", 1,0,1); 00052 hRecoMET_MET = new TH1F("hRecoMET_MET", "", 1,0,1); 00053 hGenMET_phi = new TH1F("hGenMET_phi", "", 1,0,1); 00054 00055 // SuperCluster Histos 00056 hEEpZ_SC_ix_iy = new TH2F("hEEpZ_SC_ix_iy","", 100,1,101, 100,1,101); 00057 hEEmZ_SC_ix_iy = new TH2F("hEEmZ_SC_ix_iy","", 100,1,101, 100,1,101); 00058 hEB_SC_ieta_iphi = new TH2F("hEB_SC_ieta_iphi","",171, -85, 86, 360, 1, 361); 00059 }
Reimplemented from edm::EDAnalyzer.
Definition at line 16 of file DumpEvent.cc.
References m_DataFile.
00017 { 00018 //Write out the histogram files. 00019 m_DataFile->Write(); 00020 }
void DumpEvent::WriteElectrons | ( | const edm::Event & | iEvent, | |
const edm::EventSetup & | iSetup | |||
) |
Definition at line 76 of file DumpEvent.cc.
References DEBUG, relval_parameters_module::energy, eta, edm::Event::getByLabel(), hElectron_energy, hElectron_eta, hElectron_phi, phi, and edm::Handle< T >::product().
Referenced by analyze().
00077 { 00078 // Get Reco Electrons 00079 edm::Handle<ElectronCollection> ElectronHandle; 00080 iEvent.getByLabel("pixelMatchElectrons",ElectronHandle); 00081 const ElectronCollection *Electrons = ElectronHandle.product(); 00082 ElectronCollection::const_iterator elec; 00083 DEBUG("***ELECTRONS***" ); 00084 DEBUG("Event has " << Electrons->size() << " reconstructed electrons"); 00085 int ele=0; 00086 for (elec = Electrons->begin(); elec != Electrons->end(); elec++) 00087 { 00088 Float_t eta = elec->eta(); 00089 Float_t phi = elec->phi(); 00090 Float_t energy = elec->energy(); 00091 hElectron_eta->SetBinContent(ele+1, eta); 00092 hElectron_phi->SetBinContent(ele+1, phi); 00093 hElectron_energy->SetBinContent(ele+1, energy); 00094 ele++; 00095 DEBUG( "Electron # " << ele << " : eta =" << eta << ", phi = " << phi << ", energy = " << energy ); 00096 } 00097 }
void DumpEvent::WriteJets | ( | const edm::Event & | iEvent, | |
const edm::EventSetup & | iSetup | |||
) |
Definition at line 122 of file DumpEvent.cc.
References IterativeConePu5Jets_PbPb_cff::caloTowers, DEBUG, relval_parameters_module::energy, eta, edm::Event::getByLabel(), hCaloTowerToJetMap_ieta_iphi, hJet_energy, hJet_eta, hJet_phi, i, CaloTower::id(), CaloTowerDetId::ieta(), CaloTowerDetId::iphi(), metsig::jet, njet, phi, and edm::Handle< T >::product().
Referenced by analyze().
00123 { 00124 // Jets 00125 edm::Handle<CaloJetCollection> JetHandle; 00126 iEvent.getByLabel( "iterativeCone5CaloJets", JetHandle ); 00127 const CaloJetCollection *Jets = JetHandle.product(); 00128 // Towers 00129 edm::Handle<CaloTowerCollection> caloTowers; 00130 iEvent.getByLabel( "towerMaker", caloTowers ); 00131 00132 CaloJetCollection::const_iterator jet; 00133 DEBUG("***JETS***"); 00134 DEBUG("Event has " << Jets->size() << " reconstructed jets"); 00135 int njet = 0; 00136 for (jet = Jets->begin(); jet != Jets->end(); jet++) 00137 { 00138 Float_t eta = jet->eta(); 00139 Float_t phi = jet->phi(); 00140 Float_t energy = jet->energy(); 00141 hJet_eta->SetBinContent(njet+1, eta); 00142 hJet_phi->SetBinContent(njet+1, phi); 00143 hJet_energy->SetBinContent(njet+1, energy); 00144 njet++; 00145 DEBUG( "Jet # " << njet << " : eta =" << eta << ", phi = " << phi << ", energy = " << energy ); 00146 00147 int nConstituents= jet->nConstituents(); 00148 for (int i = 0; i <nConstituents ; i++) 00149 { 00150 const CaloTower& tower = *(jet->getCaloConstituent (i)); 00151 int ietaTower = tower.id().ieta(); 00152 int iphiTower = tower.id().iphi(); 00153 hCaloTowerToJetMap_ieta_iphi->Fill(ietaTower, iphiTower, njet); 00154 } 00155 00156 } 00157 }
void DumpEvent::WriteMET | ( | const edm::Event & | iEvent, | |
const edm::EventSetup & | iSetup | |||
) |
Definition at line 184 of file DumpEvent.cc.
References DEBUG, ET, edm::Event::getByLabel(), hGenMET_phi, hRecoMET_MET, hRecoMET_phi, phi, reco::Particle::phi(), edm::Handle< T >::product(), and reco::Particle::pt().
Referenced by analyze().
00185 { 00186 { 00187 // Reco MET 00188 edm::Handle<CaloMETCollection> metHandle; 00189 iEvent.getByLabel("met", metHandle); 00190 const CaloMETCollection *MET = metHandle.product(); 00191 const CaloMET caloMET = MET->front(); 00192 Float_t phi = caloMET.phi(); 00193 Float_t ET = caloMET.pt(); 00194 hRecoMET_phi->SetBinContent(1, phi); 00195 hRecoMET_MET->SetBinContent(1, ET); 00196 DEBUG("***RecoMET***"); 00197 DEBUG("RecoMET: ET =" << ET << ", phi = " << phi); 00198 } 00199 00200 { 00201 // Get GenMET objects from event 00202 edm::Handle<GenMETCollection> genmetHandle; 00203 iEvent.getByLabel("genMet", genmetHandle); 00204 const GenMETCollection *gmet = genmetHandle.product(); 00205 const GenMET genMET = gmet->front(); 00206 Float_t phi = genMET.phi(); 00207 Float_t ET = genMET.pt(); 00208 hGenMET_phi->SetBinContent(1, phi); 00209 DEBUG("***GenMET***"); 00210 DEBUG("GenMET: MET = " << ET << ", phi = " << phi); 00211 } 00212 }
void DumpEvent::WriteMuons | ( | const edm::Event & | iEvent, | |
const edm::EventSetup & | iSetup | |||
) |
Definition at line 160 of file DumpEvent.cc.
References DEBUG, eta, edm::Event::getByLabel(), hMuon_eta, hMuon_phi, hMuon_pt, metsig::muon, phi, and edm::Handle< T >::product().
Referenced by analyze().
00161 { 00162 // Get Reco Muons 00163 edm::Handle<MuonCollection> MuonHandle; 00164 iEvent.getByLabel("globalMuons",MuonHandle); 00165 const MuonCollection *Muons = MuonHandle.product(); 00166 MuonCollection::const_iterator muon; 00167 DEBUG("***MUONS***"); 00168 DEBUG("Event has " << Muons->size() << " reconstructed muons"); 00169 int mu=0; 00170 for (muon = Muons->begin(); muon != Muons->end(); muon++) 00171 { 00172 Float_t eta = muon->eta(); 00173 Float_t phi = muon->phi(); 00174 Float_t pt = muon->pt(); 00175 hMuon_eta->SetBinContent(mu+1, eta); 00176 hMuon_phi->SetBinContent(mu+1, phi); 00177 hMuon_pt->SetBinContent(mu+1, pt); 00178 mu++; 00179 DEBUG( "Muon # " << mu << " : eta =" << eta << ", phi = " << phi << ", Pt = " << pt ); 00180 } 00181 }
void DumpEvent::WritePhotons | ( | const edm::Event & | iEvent, | |
const edm::EventSetup & | iSetup | |||
) |
Definition at line 99 of file DumpEvent.cc.
References DEBUG, relval_parameters_module::energy, eta, edm::Event::getByLabel(), hPhoton_energy, hPhoton_eta, hPhoton_phi, phi, and edm::Handle< T >::product().
Referenced by analyze().
00100 { 00101 // Get Reco Photons 00102 edm::Handle<PhotonCollection> PhotonHandle; 00103 iEvent.getByLabel("photons",PhotonHandle); 00104 const PhotonCollection *Photons = PhotonHandle.product(); 00105 PhotonCollection::const_iterator gamma; 00106 DEBUG("***PHOTONS***"); 00107 DEBUG("Event has " << Photons->size() << " reconstructed photons"); 00108 int photon=0; 00109 for (gamma = Photons->begin(); gamma != Photons->end(); gamma++) 00110 { 00111 Float_t eta = gamma->eta(); 00112 Float_t phi = gamma->phi(); 00113 Float_t energy = gamma->energy(); 00114 hPhoton_eta->SetBinContent(photon+1, eta); 00115 hPhoton_phi->SetBinContent(photon+1, phi); 00116 hPhoton_energy->SetBinContent(photon+1, energy); 00117 photon++; 00118 DEBUG( "Photon # " << photon << " : eta =" << eta << ", phi = " << phi << ", energy = " << energy ); 00119 } 00120 }
void DumpEvent::WriteSCs | ( | const edm::Event & | iEvent, | |
const edm::EventSetup & | iSetup | |||
) |
Definition at line 215 of file DumpEvent.cc.
References DEBUG, relval_parameters_module::energy, edm::Event::getByLabel(), hEB_SC_ieta_iphi, hEEmZ_SC_ix_iy, hEEpZ_SC_ix_iy, EBDetId::ieta(), EBDetId::iphi(), EEDetId::ix(), EEDetId::iy(), edm::Handle< T >::product(), SubDet, and EEDetId::zside().
Referenced by analyze().
00216 { 00217 // Get ECAL Barrel identified SuperClusters 00218 edm::Handle<SuperClusterCollection> BSCHandle; 00219 iEvent.getByLabel("correctedHybridSuperClusters",BSCHandle); 00220 const SuperClusterCollection *BSuperClusters = BSCHandle.product(); 00221 DEBUG("***Barrel Super Clusters***"); 00222 Int_t whichSC = 0; 00223 SuperClusterCollection::const_iterator BSCC;// int num = 0; 00224 DEBUG("Event has " << BSuperClusters->size() << " Barrel SuperClusters" ); 00225 int Bscs=0; 00226 for (BSCC = BSuperClusters->begin(); BSCC != BSuperClusters->end(); BSCC++) { 00227 whichSC++; 00228 Bscs++; 00229 Float_t energy = BSCC->energy(); 00230 DEBUG("SC# " << Bscs << " : eta = " << BSCC->eta() << ", phi = " << BSCC->phi() << ", Energy = " << energy ); 00231 00232 vector<DetId> crystals = BSCC->getHitsByDetId(); 00233 vector<DetId>::const_iterator crystal; 00234 //int cry = 0; 00235 for (crystal = crystals.begin(); crystal != crystals.end(); crystal++) 00236 { 00237 // get the subdetector 00238 int SubDet = (*crystal).subdetId(); 00239 //DetId::Detector DetNum=(*crystal).det(); 00240 if (SubDet == 1) 00241 { 00242 EBDetId EcalID = *crystal; 00243 float ieta = EcalID.ieta(); 00244 float iphi = EcalID.iphi(); 00245 hEB_SC_ieta_iphi->Fill(ieta, iphi, whichSC); 00246 } 00247 } 00248 } 00249 00250 // Get ECAL EndCap identified SuperClusters 00251 edm::Handle<SuperClusterCollection> ESCHandle; 00252 iEvent.getByLabel("correctedIslandEndcapSuperClusters",ESCHandle); 00253 const SuperClusterCollection *ESuperClusters = ESCHandle.product(); 00254 DEBUG("***EndCap Super Clusters***"); 00255 SuperClusterCollection::const_iterator ESCC;// int num = 0; 00256 DEBUG("Event has " << ESuperClusters->size() << " EndCap SuperClusters" ); 00257 int Escs=0; 00258 for (ESCC = ESuperClusters->begin(); ESCC != ESuperClusters->end(); ESCC++) { 00259 whichSC++; 00260 Escs++; 00261 Float_t energy = ESCC->energy(); 00262 DEBUG("SC# " << Bscs << " : eta = " << BSCC->eta() << ", phi = " << BSCC->phi() << ", Energy = " << energy ); 00263 vector<DetId> crystals = ESCC->getHitsByDetId(); 00264 vector<DetId>::const_iterator crystal; 00265 //int cry = 0; 00266 for (crystal = crystals.begin(); crystal != crystals.end(); crystal++) 00267 { 00268 // get the subdetector 00269 int SubDet = (*crystal).subdetId(); 00270 //DetId::Detector DetNum=(*crystal).det(); 00271 if (SubDet == 2) 00272 { 00273 EEDetId EcalID = *crystal; 00274 float ix = EcalID.ix(); 00275 float iy = EcalID.iy(); 00276 float Crystal_zside = EcalID.zside(); 00277 if (Crystal_zside == -1) 00278 { 00279 hEEmZ_SC_ix_iy->Fill(ix, iy, whichSC); 00280 } 00281 if (Crystal_zside == 1) 00282 { 00283 hEEpZ_SC_ix_iy->Fill(ix, iy, whichSC); 00284 } 00285 } 00286 } // end loop over crystals 00287 } 00288 }
int DumpEvent::CurrentEvent [private] |
bool DumpEvent::debug_ [private] |
float DumpEvent::EnergyThreshold [private] |
int DumpEvent::FirstEvent [private] |
TH2F* DumpEvent::hCaloTowerToJetMap_ieta_iphi [private] |
TH2F* DumpEvent::hEB_SC_ieta_iphi [private] |
TH2F* DumpEvent::hEEmZ_SC_ix_iy [private] |
TH2F* DumpEvent::hEEpZ_SC_ix_iy [private] |
TH1F * DumpEvent::hElectron_energy [private] |
TH1F* DumpEvent::hElectron_eta [private] |
TH1F * DumpEvent::hElectron_phi [private] |
TH1F* DumpEvent::hGenMET_phi [private] |
TH1F * DumpEvent::hJet_energy [private] |
TH1F* DumpEvent::hJet_eta [private] |
TH1F * DumpEvent::hJet_phi [private] |
TH1F* DumpEvent::hMuon_eta [private] |
TH1F * DumpEvent::hMuon_phi [private] |
TH1F * DumpEvent::hMuon_pt [private] |
TH1F * DumpEvent::hPhoton_energy [private] |
TH1F* DumpEvent::hPhoton_eta [private] |
TH1F * DumpEvent::hPhoton_phi [private] |
TH1F * DumpEvent::hRecoMET_MET [private] |
TH1F* DumpEvent::hRecoMET_phi [private] |
int DumpEvent::LastEvent [private] |
TFile* DumpEvent::m_DataFile [private] |
int DumpEvent::theEvent [private] |