CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
DQMSourcePi0.cc
Go to the documentation of this file.
1 #include "FWCore/Framework/interface/Event.h"
2 #include "FWCore/Framework/interface/EventSetup.h"
3 #include "FWCore/ParameterSet/interface/ParameterSet.h"
4 #include "FWCore/ServiceRegistry/interface/Service.h"
5 
6 // DQM include files
7 
8 #include "DQMServices/Core/interface/MonitorElement.h"
9 
10 // work on collections
11 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
12 #include "DataFormats/EcalDetId/interface/EBDetId.h"
13 #include "DataFormats/EcalDetId/interface/EEDetId.h"
14 
15 #include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
16 
17 #include "DQMServices/Core/interface/DQMStore.h"
18 #include "DQMOffline/CalibCalo/src/DQMSourcePi0.h"
19 
20 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
21 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
22 #include "DataFormats/EcalDetId/interface/EBDetId.h"
23 #include "DataFormats/DetId/interface/DetId.h"
24 #include "DataFormats/Math/interface/Point3D.h"
25 
26 #include "FWCore/Framework/interface/Event.h"
27 #include "FWCore/Framework/interface/EventSetup.h"
28 #include "DataFormats/Common/interface/Handle.h"
29 #include "FWCore/Framework/interface/ESHandle.h"
30 #include "FWCore/MessageLogger/interface/MessageLogger.h"
31 
32 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
33 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidateFwd.h"
35 //#include "RecoEcal/EgammaCoreTools/interface/PositionCalc.h"
36 #include "RecoEcal/EgammaCoreTools/interface/EcalEtaPhiRegion.h"
37 
38 #include "TVector3.h"
39 
40 #define TWOPI 6.283185308
41 
42 
43 
44 using namespace std;
45 using namespace edm;
46 
47 
48 // ******************************************
49 // constructors
50 // *****************************************
51 
52 DQMSourcePi0::DQMSourcePi0( const edm::ParameterSet& ps ) :
53 eventCounter_(0)
54 {
55  folderName_ = ps.getUntrackedParameter<string>("FolderName","HLT/AlCaEcalPi0");
56  prescaleFactor_ = ps.getUntrackedParameter<int>("prescaleFactor",1);
57  productMonitoredEBpi0_= consumes<EcalRecHitCollection>(ps.getUntrackedParameter<edm::InputTag>("AlCaStreamEBpi0Tag"));
58  productMonitoredEBeta_= consumes<EcalRecHitCollection>(ps.getUntrackedParameter<edm::InputTag>("AlCaStreamEBetaTag"));
59  productMonitoredEEpi0_= consumes<EcalRecHitCollection>(ps.getUntrackedParameter<edm::InputTag>("AlCaStreamEEpi0Tag"));
60  productMonitoredEEeta_= consumes<EcalRecHitCollection>(ps.getUntrackedParameter<edm::InputTag>("AlCaStreamEEetaTag"));
61 
62  isMonEBpi0_ = ps.getUntrackedParameter<bool>("isMonEBpi0",false);
63  isMonEBeta_ = ps.getUntrackedParameter<bool>("isMonEBeta",false);
64  isMonEEpi0_ = ps.getUntrackedParameter<bool>("isMonEEpi0",false);
65  isMonEEeta_ = ps.getUntrackedParameter<bool>("isMonEEeta",false);
66 
67  saveToFile_=ps.getUntrackedParameter<bool>("SaveToFile",false);
68  fileName_= ps.getUntrackedParameter<string>("FileName","MonitorAlCaEcalPi0.root");
69 
70  clusSeedThr_ = ps.getParameter<double> ("clusSeedThr");
71  clusSeedThrEndCap_ = ps.getParameter<double> ("clusSeedThrEndCap");
72  clusEtaSize_ = ps.getParameter<int> ("clusEtaSize");
73  clusPhiSize_ = ps.getParameter<int> ("clusPhiSize");
74  if ( clusPhiSize_ % 2 == 0 || clusEtaSize_ % 2 == 0)
75  edm::LogError("AlCaPi0RecHitsProducerError") << "Size of eta/phi for simple clustering should be odd numbers";
76 
77 
78  seleXtalMinEnergy_ = ps.getParameter<double>("seleXtalMinEnergy");
79  seleXtalMinEnergyEndCap_ = ps.getParameter<double>("seleXtalMinEnergyEndCap");
80 
82  selePtGamma_ = ps.getParameter<double> ("selePtGamma");
83  selePtPi0_ = ps.getParameter<double> ("selePtPi0");
84  seleMinvMaxPi0_ = ps.getParameter<double> ("seleMinvMaxPi0");
85  seleMinvMinPi0_ = ps.getParameter<double> ("seleMinvMinPi0");
86  seleS4S9Gamma_ = ps.getParameter<double> ("seleS4S9Gamma");
87  selePi0Iso_ = ps.getParameter<double> ("selePi0Iso");
88  ptMinForIsolation_ = ps.getParameter<double> ("ptMinForIsolation");
89  selePi0BeltDR_ = ps.getParameter<double> ("selePi0BeltDR");
90  selePi0BeltDeta_ = ps.getParameter<double> ("selePi0BeltDeta");
91 
92 
94  selePtGammaEndCap_ = ps.getParameter<double> ("selePtGammaEndCap");
95  selePtPi0EndCap_ = ps.getParameter<double> ("selePtPi0EndCap");
96  seleS4S9GammaEndCap_ = ps.getParameter<double> ("seleS4S9GammaEndCap");
97  seleMinvMaxPi0EndCap_ = ps.getParameter<double> ("seleMinvMaxPi0EndCap");
98  seleMinvMinPi0EndCap_ = ps.getParameter<double> ("seleMinvMinPi0EndCap");
99  ptMinForIsolationEndCap_ = ps.getParameter<double> ("ptMinForIsolationEndCap");
100  selePi0BeltDREndCap_ = ps.getParameter<double> ("selePi0BeltDREndCap");
101  selePi0BeltDetaEndCap_ = ps.getParameter<double> ("selePi0BeltDetaEndCap");
102  selePi0IsoEndCap_ = ps.getParameter<double> ("selePi0IsoEndCap");
103 
104 
106  selePtGammaEta_ = ps.getParameter<double> ("selePtGammaEta");
107  selePtEta_ = ps.getParameter<double> ("selePtEta");
108  seleS4S9GammaEta_ = ps.getParameter<double> ("seleS4S9GammaEta");
109  seleS9S25GammaEta_ = ps.getParameter<double> ("seleS9S25GammaEta");
110  seleMinvMaxEta_ = ps.getParameter<double> ("seleMinvMaxEta");
111  seleMinvMinEta_ = ps.getParameter<double> ("seleMinvMinEta");
112  ptMinForIsolationEta_ = ps.getParameter<double> ("ptMinForIsolationEta");
113  seleEtaIso_ = ps.getParameter<double> ("seleEtaIso");
114  seleEtaBeltDR_ = ps.getParameter<double> ("seleEtaBeltDR");
115  seleEtaBeltDeta_ = ps.getParameter<double> ("seleEtaBeltDeta");
116 
117 
119  selePtGammaEtaEndCap_ = ps.getParameter<double> ("selePtGammaEtaEndCap");
120  selePtEtaEndCap_ = ps.getParameter<double> ("selePtEtaEndCap");
121  seleS4S9GammaEtaEndCap_ = ps.getParameter<double> ("seleS4S9GammaEtaEndCap");
122  seleS9S25GammaEtaEndCap_ = ps.getParameter<double> ("seleS9S25GammaEtaEndCap");
123  seleMinvMaxEtaEndCap_ = ps.getParameter<double> ("seleMinvMaxEtaEndCap");
124  seleMinvMinEtaEndCap_ = ps.getParameter<double> ("seleMinvMinEtaEndCap");
125  ptMinForIsolationEtaEndCap_ = ps.getParameter<double> ("ptMinForIsolationEtaEndCap");
126  seleEtaIsoEndCap_ = ps.getParameter<double> ("seleEtaIsoEndCap");
127  seleEtaBeltDREndCap_ = ps.getParameter<double> ("seleEtaBeltDREndCap");
128  seleEtaBeltDetaEndCap_ = ps.getParameter<double> ("seleEtaBeltDetaEndCap");
129 
130 
131 
132  // Parameters for the position calculation:
133  edm::ParameterSet posCalcParameters =
134  ps.getParameter<edm::ParameterSet>("posCalcParameters");
135  posCalculator_ = PositionCalc(posCalcParameters);
136 
137 }
138 
139 
140 DQMSourcePi0::~DQMSourcePi0()
141 {}
142 
143 
144 //--------------------------------------------------------
145 void DQMSourcePi0::bookHistograms(DQMStore::IBooker & ibooker, edm::Run const & irun, edm::EventSetup const & isetup) {
146 
147 
148  // create and cd into new folder
149  ibooker.setCurrentFolder(folderName_);
150 
151  // book some histograms 1D
152 
153  hiPhiDistrEBpi0_ = ibooker.book1D("iphiDistributionEBpi0", "RechitEB pi0 iphi", 361, 1,361);
154  hiPhiDistrEBpi0_->setAxisTitle("i#phi ", 1);
155  hiPhiDistrEBpi0_->setAxisTitle("# rechits", 2);
156 
157  hiXDistrEEpi0_ = ibooker.book1D("iXDistributionEEpi0", "RechitEE pi0 ix", 100, 0,100);
158  hiXDistrEEpi0_->setAxisTitle("ix ", 1);
159  hiXDistrEEpi0_->setAxisTitle("# rechits", 2);
160 
161  hiPhiDistrEBeta_ = ibooker.book1D("iphiDistributionEBeta", "RechitEB eta iphi", 361, 1,361);
162  hiPhiDistrEBeta_->setAxisTitle("i#phi ", 1);
163  hiPhiDistrEBeta_->setAxisTitle("# rechits", 2);
164 
165  hiXDistrEEeta_ = ibooker.book1D("iXDistributionEEeta", "RechitEE eta ix", 100, 0,100);
166  hiXDistrEEeta_->setAxisTitle("ix ", 1);
167  hiXDistrEEeta_->setAxisTitle("# rechits", 2);
168 
169 
170  hiEtaDistrEBpi0_ = ibooker.book1D("iEtaDistributionEBpi0", "RechitEB pi0 ieta", 171, -85, 86);
171  hiEtaDistrEBpi0_->setAxisTitle("i#eta", 1);
172  hiEtaDistrEBpi0_->setAxisTitle("#rechits", 2);
173 
174  hiYDistrEEpi0_ = ibooker.book1D("iYDistributionEEpi0", "RechitEE pi0 iY", 100, 0, 100);
175  hiYDistrEEpi0_->setAxisTitle("iy", 1);
176  hiYDistrEEpi0_->setAxisTitle("#rechits", 2);
177 
178  hiEtaDistrEBeta_ = ibooker.book1D("iEtaDistributionEBeta", "RechitEB eta ieta", 171, -85, 86);
179  hiEtaDistrEBeta_->setAxisTitle("i#eta", 1);
180  hiEtaDistrEBeta_->setAxisTitle("#rechits", 2);
181 
182  hiYDistrEEeta_ = ibooker.book1D("iYDistributionEEeta", "RechitEE eta iY", 100, 0, 100);
183  hiYDistrEEeta_->setAxisTitle("iy", 1);
184  hiYDistrEEeta_->setAxisTitle("#rechits", 2);
185 
186 
187  hRechitEnergyEBpi0_ = ibooker.book1D("rhEnergyEBpi0","Pi0 rechits energy EB",160,0.,2.0);
188  hRechitEnergyEBpi0_->setAxisTitle("energy (GeV) ",1);
189  hRechitEnergyEBpi0_->setAxisTitle("#rechits",2);
190 
191  hRechitEnergyEEpi0_ = ibooker.book1D("rhEnergyEEpi0","Pi0 rechits energy EE",160,0.,3.0);
192  hRechitEnergyEEpi0_->setAxisTitle("energy (GeV) ",1);
193  hRechitEnergyEEpi0_->setAxisTitle("#rechits",2);
194 
195  hRechitEnergyEBeta_ = ibooker.book1D("rhEnergyEBeta","Eta rechits energy EB",160,0.,2.0);
196  hRechitEnergyEBeta_->setAxisTitle("energy (GeV) ",1);
197  hRechitEnergyEBeta_->setAxisTitle("#rechits",2);
198 
199  hRechitEnergyEEeta_ = ibooker.book1D("rhEnergyEEeta","Eta rechits energy EE",160,0.,3.0);
200  hRechitEnergyEEeta_->setAxisTitle("energy (GeV) ",1);
201  hRechitEnergyEEeta_->setAxisTitle("#rechits",2);
202 
203  hEventEnergyEBpi0_ = ibooker.book1D("eventEnergyEBpi0","Pi0 event energy EB",100,0.,20.0);
204  hEventEnergyEBpi0_->setAxisTitle("energy (GeV) ",1);
205 
206  hEventEnergyEEpi0_ = ibooker.book1D("eventEnergyEEpi0","Pi0 event energy EE",100,0.,50.0);
207  hEventEnergyEEpi0_->setAxisTitle("energy (GeV) ",1);
208 
209  hEventEnergyEBeta_ = ibooker.book1D("eventEnergyEBeta","Eta event energy EB",100,0.,20.0);
210  hEventEnergyEBeta_->setAxisTitle("energy (GeV) ",1);
211 
212  hEventEnergyEEeta_ = ibooker.book1D("eventEnergyEEeta","Eta event energy EE",100,0.,50.0);
213  hEventEnergyEEeta_->setAxisTitle("energy (GeV) ",1);
214 
215  hNRecHitsEBpi0_ = ibooker.book1D("nRechitsEBpi0","#rechits in pi0 collection EB",100,0.,250.);
216  hNRecHitsEBpi0_->setAxisTitle("rechits ",1);
217 
218  hNRecHitsEEpi0_ = ibooker.book1D("nRechitsEEpi0","#rechits in pi0 collection EE",100,0.,250.);
219  hNRecHitsEEpi0_->setAxisTitle("rechits ",1);
220 
221  hNRecHitsEBeta_ = ibooker.book1D("nRechitsEBeta","#rechits in eta collection EB",100,0.,250.);
222  hNRecHitsEBeta_->setAxisTitle("rechits ",1);
223 
224  hNRecHitsEEeta_ = ibooker.book1D("nRechitsEEeta","#rechits in eta collection EE",100,0.,250.);
225  hNRecHitsEEeta_->setAxisTitle("rechits ",1);
226 
227  hMeanRecHitEnergyEBpi0_ = ibooker.book1D("meanEnergyEBpi0","Mean rechit energy in pi0 collection EB",50,0.,2.);
228  hMeanRecHitEnergyEBpi0_->setAxisTitle("Mean Energy [GeV] ",1);
229 
230  hMeanRecHitEnergyEEpi0_ = ibooker.book1D("meanEnergyEEpi0","Mean rechit energy in pi0 collection EE",100,0.,5.);
231  hMeanRecHitEnergyEEpi0_->setAxisTitle("Mean Energy [GeV] ",1);
232 
233  hMeanRecHitEnergyEBeta_ = ibooker.book1D("meanEnergyEBeta","Mean rechit energy in eta collection EB",50,0.,2.);
234  hMeanRecHitEnergyEBeta_->setAxisTitle("Mean Energy [GeV] ",1);
235 
236  hMeanRecHitEnergyEEeta_ = ibooker.book1D("meanEnergyEEeta","Mean rechit energy in eta collection EE",100,0.,5.);
237  hMeanRecHitEnergyEEeta_->setAxisTitle("Mean Energy [GeV] ",1);
238 
239  hMinvPi0EB_ = ibooker.book1D("Pi0InvmassEB","Pi0 Invariant Mass in EB",100,0.,0.5);
240  hMinvPi0EB_->setAxisTitle("Inv Mass [GeV] ",1);
241 
242  hMinvPi0EE_ = ibooker.book1D("Pi0InvmassEE","Pi0 Invariant Mass in EE",100,0.,0.5);
243  hMinvPi0EE_->setAxisTitle("Inv Mass [GeV] ",1);
244 
245  hMinvEtaEB_ = ibooker.book1D("EtaInvmassEB","Eta Invariant Mass in EB",100,0.,0.85);
246  hMinvEtaEB_->setAxisTitle("Inv Mass [GeV] ",1);
247 
248  hMinvEtaEE_ = ibooker.book1D("EtaInvmassEE","Eta Invariant Mass in EE",100,0.,0.85);
249  hMinvEtaEE_->setAxisTitle("Inv Mass [GeV] ",1);
250 
251 
252  hPt1Pi0EB_ = ibooker.book1D("Pt1Pi0EB","Pt 1st most energetic Pi0 photon in EB",100,0.,20.);
253  hPt1Pi0EB_->setAxisTitle("1st photon Pt [GeV] ",1);
254 
255  hPt1Pi0EE_ = ibooker.book1D("Pt1Pi0EE","Pt 1st most energetic Pi0 photon in EE",100,0.,20.);
256  hPt1Pi0EE_->setAxisTitle("1st photon Pt [GeV] ",1);
257 
258  hPt1EtaEB_ = ibooker.book1D("Pt1EtaEB","Pt 1st most energetic Eta photon in EB",100,0.,20.);
259  hPt1EtaEB_->setAxisTitle("1st photon Pt [GeV] ",1);
260 
261  hPt1EtaEE_ = ibooker.book1D("Pt1EtaEE","Pt 1st most energetic Eta photon in EE",100,0.,20.);
262  hPt1EtaEE_->setAxisTitle("1st photon Pt [GeV] ",1);
263 
264  hPt2Pi0EB_ = ibooker.book1D("Pt2Pi0EB","Pt 2nd most energetic Pi0 photon in EB",100,0.,20.);
265  hPt2Pi0EB_->setAxisTitle("2nd photon Pt [GeV] ",1);
266 
267  hPt2Pi0EE_ = ibooker.book1D("Pt2Pi0EE","Pt 2nd most energetic Pi0 photon in EE",100,0.,20.);
268  hPt2Pi0EE_->setAxisTitle("2nd photon Pt [GeV] ",1);
269 
270  hPt2EtaEB_ = ibooker.book1D("Pt2EtaEB","Pt 2nd most energetic Eta photon in EB",100,0.,20.);
271  hPt2EtaEB_->setAxisTitle("2nd photon Pt [GeV] ",1);
272 
273  hPt2EtaEE_ = ibooker.book1D("Pt2EtaEE","Pt 2nd most energetic Eta photon in EE",100,0.,20.);
274  hPt2EtaEE_->setAxisTitle("2nd photon Pt [GeV] ",1);
275 
276 
277  hPtPi0EB_ = ibooker.book1D("PtPi0EB","Pi0 Pt in EB",100,0.,20.);
278  hPtPi0EB_->setAxisTitle("Pi0 Pt [GeV] ",1);
279 
280  hPtPi0EE_ = ibooker.book1D("PtPi0EE","Pi0 Pt in EE",100,0.,20.);
281  hPtPi0EE_->setAxisTitle("Pi0 Pt [GeV] ",1);
282 
283  hPtEtaEB_ = ibooker.book1D("PtEtaEB","Eta Pt in EB",100,0.,20.);
284  hPtEtaEB_->setAxisTitle("Eta Pt [GeV] ",1);
285 
286  hPtEtaEE_ = ibooker.book1D("PtEtaEE","Eta Pt in EE",100,0.,20.);
287  hPtEtaEE_->setAxisTitle("Eta Pt [GeV] ",1);
288 
289  hIsoPi0EB_ = ibooker.book1D("IsoPi0EB","Pi0 Iso in EB",50,0.,1.);
290  hIsoPi0EB_->setAxisTitle("Pi0 Iso",1);
291 
292  hIsoPi0EE_ = ibooker.book1D("IsoPi0EE","Pi0 Iso in EE",50,0.,1.);
293  hIsoPi0EE_->setAxisTitle("Pi0 Iso",1);
294 
295  hIsoEtaEB_ = ibooker.book1D("IsoEtaEB","Eta Iso in EB",50,0.,1.);
296  hIsoEtaEB_->setAxisTitle("Eta Iso",1);
297 
298  hIsoEtaEE_ = ibooker.book1D("IsoEtaEE","Eta Iso in EE",50,0.,1.);
299  hIsoEtaEE_->setAxisTitle("Eta Iso",1);
300 
301  hS4S91Pi0EB_ = ibooker.book1D("S4S91Pi0EB","S4S9 1st most energetic Pi0 photon in EB",50,0.,1.);
302  hS4S91Pi0EB_->setAxisTitle("S4S9 of the 1st Pi0 Photon ",1);
303 
304  hS4S91Pi0EE_ = ibooker.book1D("S4S91Pi0EE","S4S9 1st most energetic Pi0 photon in EE",50,0.,1.);
305  hS4S91Pi0EE_->setAxisTitle("S4S9 of the 1st Pi0 Photon ",1);
306 
307  hS4S91EtaEB_ = ibooker.book1D("S4S91EtaEB","S4S9 1st most energetic Eta photon in EB",50,0.,1.);
308  hS4S91EtaEB_->setAxisTitle("S4S9 of the 1st Eta Photon ",1);
309 
310  hS4S91EtaEE_ = ibooker.book1D("S4S91EtaEE","S4S9 1st most energetic Eta photon in EE",50,0.,1.);
311  hS4S91EtaEE_->setAxisTitle("S4S9 of the 1st Eta Photon ",1);
312 
313  hS4S92Pi0EB_ = ibooker.book1D("S4S92Pi0EB","S4S9 2nd most energetic Pi0 photon in EB",50,0.,1.);
314  hS4S92Pi0EB_->setAxisTitle("S4S9 of the 2nd Pi0 Photon",1);
315 
316  hS4S92Pi0EE_ = ibooker.book1D("S4S92Pi0EE","S4S9 2nd most energetic Pi0 photon in EE",50,0.,1.);
317  hS4S92Pi0EE_->setAxisTitle("S4S9 of the 2nd Pi0 Photon",1);
318 
319  hS4S92EtaEB_ = ibooker.book1D("S4S92EtaEB","S4S9 2nd most energetic Pi0 photon in EB",50,0.,1.);
320  hS4S92EtaEB_->setAxisTitle("S4S9 of the 2nd Eta Photon",1);
321 
322  hS4S92EtaEE_ = ibooker.book1D("S4S92EtaEE","S4S9 2nd most energetic Pi0 photon in EE",50,0.,1.);
323  hS4S92EtaEE_->setAxisTitle("S4S9 of the 2nd Eta Photon",1);
324 
325 
326 
327 
328 }
329 
330 //--------------------------------------------------------
331 //void DQMSourcePi0::beginRun(const edm::Run& r, const EventSetup& context) {
332 //
333 //}
334 
335 //--------------------------------------------------------
336 void DQMSourcePi0::beginLuminosityBlock(const LuminosityBlock& lumiSeg,
337  const EventSetup& context) {
338 
339 }
340 
341 //-------------------------------------------------------------
342 
343 void DQMSourcePi0::analyze(const Event& iEvent,
344  const EventSetup& iSetup ){
345 
346  if (eventCounter_% prescaleFactor_ ) return;
347  eventCounter_++;
348 
349  edm::ESHandle<CaloTopology> theCaloTopology;
350  iSetup.get<CaloTopologyRecord>().get(theCaloTopology);
351 
352 
353  std::vector<EcalRecHit> seeds;
354  seeds.clear();
355 
356  vector<EBDetId> usedXtals;
357  usedXtals.clear();
358 
359  detIdEBRecHits.clear();
360  EBRecHits.clear();
361 
362 
363 
364 
365 
366 
367  edm::Handle<EcalRecHitCollection> rhEBpi0;
368  edm::Handle<EcalRecHitCollection> rhEBeta;
369  edm::Handle<EcalRecHitCollection> rhEEpi0;
370  edm::Handle<EcalRecHitCollection> rhEEeta;
371 
372  if(isMonEBpi0_) iEvent.getByToken(productMonitoredEBpi0_, rhEBpi0);
373  if(isMonEBeta_) iEvent.getByToken(productMonitoredEBeta_, rhEBeta);
374  if(isMonEEpi0_) iEvent.getByToken(productMonitoredEEpi0_, rhEEpi0);
375  if(isMonEEeta_) iEvent.getByToken(productMonitoredEEeta_, rhEEeta);
376 
377  // Initialize the Position Calc
378  const CaloSubdetectorGeometry *geometry_p;
379  const CaloSubdetectorGeometry *geometryEE_p;
380  const CaloSubdetectorGeometry *geometryES_p;
381 
382  const CaloSubdetectorTopology *topology_p;
383  const CaloSubdetectorTopology *topology_ee;
384 
385 
386 
387  edm::ESHandle<CaloGeometry> geoHandle;
388  iSetup.get<CaloGeometryRecord>().get(geoHandle);
389  geometry_p = geoHandle->getSubdetectorGeometry(DetId::Ecal,EcalBarrel);
390  geometryEE_p = geoHandle->getSubdetectorGeometry(DetId::Ecal,EcalEndcap);
391  geometryES_p = geoHandle->getSubdetectorGeometry(DetId::Ecal, EcalPreshower);
392  topology_p = theCaloTopology->getSubdetectorTopology(DetId::Ecal,EcalBarrel);
393  topology_ee = theCaloTopology->getSubdetectorTopology(DetId::Ecal,EcalEndcap);
394 
395  EcalRecHitCollection::const_iterator itb;
396 
397  // fill EB pi0 histos
398  if(isMonEBpi0_ ){
399  if (rhEBpi0.isValid() && (rhEBpi0->size() > 0)){
400 
401 
402  const EcalRecHitCollection *hitCollection_p = rhEBpi0.product();
403  float etot =0;
404  for(itb=rhEBpi0->begin(); itb!=rhEBpi0->end(); ++itb){
405 
406  EBDetId id(itb->id());
407  double energy = itb->energy();
408  if( energy < seleXtalMinEnergy_) continue;
409 
410  EBDetId det = itb->id();
411 
412 
413  detIdEBRecHits.push_back(det);
414  EBRecHits.push_back(*itb);
415 
416  if (energy > clusSeedThr_) seeds.push_back(*itb);
417 
418  hiPhiDistrEBpi0_->Fill(id.iphi());
419  hiEtaDistrEBpi0_->Fill(id.ieta());
420  hRechitEnergyEBpi0_->Fill(itb->energy());
421 
422  etot+= itb->energy();
423  } // Eb rechits
424 
425  hNRecHitsEBpi0_->Fill(rhEBpi0->size());
426  hMeanRecHitEnergyEBpi0_->Fill(etot/rhEBpi0->size());
427  hEventEnergyEBpi0_->Fill(etot);
428 
429  // cout << " EB RH Pi0 collection: #, mean rh_e, event E "<<rhEBpi0->size()<<" "<<etot/rhEBpi0->size()<<" "<<etot<<endl;
430 
431 
432  // Pi0 maker
433 
434  //cout<< " RH coll size: "<<rhEBpi0->size()<<endl;
435  //cout<< " Pi0 seeds: "<<seeds.size()<<endl;
436 
437  int nClus;
438  vector<float> eClus;
439  vector<float> etClus;
440  vector<float> etaClus;
441  vector<float> thetaClus;
442  vector<float> phiClus;
443  vector<EBDetId> max_hit;
444 
445  vector< vector<EcalRecHit> > RecHitsCluster;
446  vector< vector<EcalRecHit> > RecHitsCluster5x5;
447  vector<float> s4s9Clus;
448  vector<float> s9s25Clus;
449 
450 
451  nClus=0;
452 
453 
454  // Make own simple clusters (3x3, 5x5 or clusPhiSize_ x clusEtaSize_)
455  sort(seeds.begin(), seeds.end(), ecalRecHitLess());
456 
457  for (std::vector<EcalRecHit>::iterator itseed=seeds.begin(); itseed!=seeds.end(); itseed++) {
458  EBDetId seed_id = itseed->id();
459  std::vector<EBDetId>::const_iterator usedIds;
460 
461  bool seedAlreadyUsed=false;
462  for(usedIds=usedXtals.begin(); usedIds!=usedXtals.end(); usedIds++){
463  if(*usedIds==seed_id){
464  seedAlreadyUsed=true;
465  //cout<< " Seed with energy "<<itseed->energy()<<" was used !"<<endl;
466  break;
467  }
468  }
469  if(seedAlreadyUsed)continue;
470  std::vector<DetId> clus_v = topology_p->getWindow(seed_id,clusEtaSize_,clusPhiSize_);
471  std::vector<std::pair<DetId,float> > clus_used;
472  //Reject the seed if not able to build the cluster around it correctly
473  //if(clus_v.size() < clusEtaSize_*clusPhiSize_){cout<<" Not enough RecHits "<<endl; continue;}
474  vector<EcalRecHit> RecHitsInWindow;
475  vector<EcalRecHit> RecHitsInWindow5x5;
476 
477  double simple_energy = 0;
478 
479  for (std::vector<DetId >::iterator det=clus_v.begin(); det!=clus_v.end(); det++) {
480  EBDetId EBdet = *det;
481  // cout<<" det "<< EBdet<<" ieta "<<EBdet.ieta()<<" iphi "<<EBdet.iphi()<<endl;
482  bool HitAlreadyUsed=false;
483  for(usedIds=usedXtals.begin(); usedIds!=usedXtals.end(); usedIds++){
484  if(*usedIds==*det){
485  HitAlreadyUsed=true;
486  break;
487  }
488  }
489  if(HitAlreadyUsed)continue;
490 
491 
492  std::vector<EBDetId>::iterator itdet = find( detIdEBRecHits.begin(),detIdEBRecHits.end(),EBdet);
493  if(itdet == detIdEBRecHits.end()) continue;
494 
495  int nn = int(itdet - detIdEBRecHits.begin());
496  usedXtals.push_back(*det);
497  RecHitsInWindow.push_back(EBRecHits[nn]);
498  clus_used.push_back(std::make_pair(*det,1));
499  simple_energy = simple_energy + EBRecHits[nn].energy();
500 
501 
502  }
503 
504  if(simple_energy <= 0) continue;
505 
506  math::XYZPoint clus_pos = posCalculator_.Calculate_Location(clus_used,hitCollection_p,geometry_p,geometryES_p);
507  //cout<< " Simple Clustering: Total energy for this simple cluster : "<<simple_energy<<endl;
508  //cout<< " Simple Clustering: eta phi : "<<clus_pos.eta()<<" "<<clus_pos.phi()<<endl;
509  //cout<< " Simple Clustering: x y z : "<<clus_pos.x()<<" "<<clus_pos.y()<<" "<<clus_pos.z()<<endl;
510 
511  float theta_s = 2. * atan(exp(-clus_pos.eta()));
512  // float p0x_s = simple_energy * sin(theta_s) * cos(clus_pos.phi());
513  //float p0y_s = simple_energy * sin(theta_s) * sin(clus_pos.phi());
514  // float p0z_s = simple_energy * cos(theta_s);
515  //float et_s = sqrt( p0x_s*p0x_s + p0y_s*p0y_s);
516 
517  float et_s = simple_energy * sin(theta_s);
518  //cout << " Simple Clustering: E,Et,px,py,pz: "<<simple_energy<<" "<<et_s<<" "<<p0x_s<<" "<<p0y_s<<" "<<endl;
519 
520  //Compute S4/S9 variable
521  //We are not sure to have 9 RecHits so need to check eta and phi:
523  float s4s9_tmp[4];
524  for(int i=0;i<4;i++)s4s9_tmp[i]= 0;
525 
526  int seed_ieta = seed_id.ieta();
527  int seed_iphi = seed_id.iphi();
528 
529  convxtalid( seed_iphi,seed_ieta);
530 
531  float e3x3 = 0;
532  float e5x5 = 0;
533 
534  for(unsigned int j=0; j<RecHitsInWindow.size();j++){
535  EBDetId det = (EBDetId)RecHitsInWindow[j].id();
536 
537  int ieta = det.ieta();
538  int iphi = det.iphi();
539 
540  convxtalid(iphi,ieta);
541 
542  float en = RecHitsInWindow[j].energy();
543 
544  int dx = diff_neta_s(seed_ieta,ieta);
545  int dy = diff_nphi_s(seed_iphi,iphi);
546 
547  if(dx <= 0 && dy <=0) s4s9_tmp[0] += en;
548  if(dx >= 0 && dy <=0) s4s9_tmp[1] += en;
549  if(dx <= 0 && dy >=0) s4s9_tmp[2] += en;
550  if(dx >= 0 && dy >=0) s4s9_tmp[3] += en;
551 
552  if(std::abs(dx)<=1 && std::abs(dy)<=1) e3x3 += en;
553  if(std::abs(dx)<=2 && std::abs(dy)<=2) e5x5 += en;
554 
555 
556  }
557 
558  if(e3x3 <= 0) continue;
559 
560  float s4s9_max = *max_element( s4s9_tmp,s4s9_tmp+4)/e3x3;
561 
562 
564  std::vector<DetId> clus_v5x5 = topology_p->getWindow(seed_id,5,5);
565  for( std::vector<DetId>::const_iterator idItr = clus_v5x5.begin(); idItr != clus_v5x5.end(); idItr++){
566  EBDetId det = *idItr;
567 
568 
569  std::vector<EBDetId>::iterator itdet0 = find(usedXtals.begin(),usedXtals.end(),det);
570 
572  if(itdet0 != usedXtals.end()) continue;
573 
574  //inside collections
575  std::vector<EBDetId>::iterator itdet = find( detIdEBRecHits.begin(),detIdEBRecHits.end(),det);
576  if(itdet == detIdEBRecHits.end()) continue;
577 
578  int nn = int(itdet - detIdEBRecHits.begin());
579 
580  RecHitsInWindow5x5.push_back(EBRecHits[nn]);
581  e5x5 += EBRecHits[nn].energy();
582 
583  }
584 
585 
586 
587  if(e5x5 <= 0) continue;
588 
589  eClus.push_back(simple_energy);
590  etClus.push_back(et_s);
591  etaClus.push_back(clus_pos.eta());
592  thetaClus.push_back(theta_s);
593  phiClus.push_back(clus_pos.phi());
594  s4s9Clus.push_back(s4s9_max);
595  s9s25Clus.push_back(e3x3/e5x5);
596  RecHitsCluster.push_back(RecHitsInWindow);
597  RecHitsCluster5x5.push_back(RecHitsInWindow5x5);
598 
599  // std::cout<<" EB pi0 cluster (n,nxt,e,et eta,phi,s4s9) "<<nClus<<" "<<int(RecHitsInWindow.size())<<" "<<eClus[nClus]<<" "<<" "<<etClus[nClus]<<" "<<etaClus[nClus]<<" "<<phiClus[nClus]<<" "<<s4s9Clus[nClus]<<std::endl;
600 
601  nClus++;
602 
603 
604  }
605 
606  // cout<< " Pi0 clusters: "<<nClus<<endl;
607 
608  // Selection, based on Simple clustering
609  //pi0 candidates
610  int npi0_s=0;
611 
612 
613  // if (nClus <= 1) return;
614  for(Int_t i=0 ; i<nClus ; i++){
615  for(Int_t j=i+1 ; j<nClus ; j++){
616  // cout<<" i "<<i<<" etClus[i] "<<etClus[i]<<" j "<<j<<" etClus[j] "<<etClus[j]<<endl;
617  if( etClus[i]>selePtGamma_ && etClus[j]>selePtGamma_ && s4s9Clus[i]>seleS4S9Gamma_ && s4s9Clus[j]>seleS4S9Gamma_){
618 
619 
620  float p0x = etClus[i] * cos(phiClus[i]);
621  float p1x = etClus[j] * cos(phiClus[j]);
622  float p0y = etClus[i] * sin(phiClus[i]);
623  float p1y = etClus[j] * sin(phiClus[j]);
624  float p0z = eClus[i] * cos(thetaClus[i]);
625  float p1z = eClus[j] * cos(thetaClus[j]);
626 
627 
628  float pt_pair = sqrt( (p0x+p1x)*(p0x+p1x) + (p0y+p1y)*(p0y+p1y));
629 
630  if (pt_pair < selePtPi0_)continue;
631 
632  float m_inv = sqrt ( (eClus[i] + eClus[j])*(eClus[i] + eClus[j]) - (p0x+p1x)*(p0x+p1x) - (p0y+p1y)*(p0y+p1y) - (p0z+p1z)*(p0z+p1z) );
633  if ( (m_inv<seleMinvMaxPi0_) && (m_inv>seleMinvMinPi0_) ){
634 
635  //New Loop on cluster to measure isolation:
636  vector<int> IsoClus;
637  IsoClus.clear();
638  float Iso = 0;
639  TVector3 pairVect = TVector3((p0x+p1x), (p0y+p1y), (p0z+p1z));
640  for(Int_t k=0 ; k<nClus ; k++){
641 
642 
643  if(etClus[k] < ptMinForIsolation_) continue;
644 
645  if(k==i || k==j)continue;
646  TVector3 ClusVect = TVector3(etClus[k] *cos(phiClus[k]), etClus[k] * sin(phiClus[k]) , eClus[k] * cos(thetaClus[k]));
647 
648  float dretacl = fabs(etaClus[k] - pairVect.Eta());
649  float drcl = ClusVect.DeltaR(pairVect);
650  // cout<< " Iso: k, E, drclpi0, detaclpi0, dphiclpi0 "<<k<<" "<<eClus[k]<<" "<<drclpi0<<" "<<dretaclpi0<<endl;
651  if((drcl<selePi0BeltDR_) && (dretacl<selePi0BeltDeta_) ){
652  // cout<< " ... good iso cluster #: "<<k<<" etClus[k] "<<etClus[k] <<endl;
653  Iso = Iso + etClus[k];
654  IsoClus.push_back(k);
655  }
656  }
657 
658  // cout<<" Iso/pt_pi0 "<<Iso/pt_pi0<<endl;
659  if(Iso/pt_pair<selePi0Iso_){
660  //for(unsigned int Rec=0;Rec<RecHitsCluster[i].size();Rec++)pi0EBRecHitCollection->push_back(RecHitsCluster[i][Rec]);
661  //for(unsigned int Rec2=0;Rec2<RecHitsCluster[j].size();Rec2++)pi0EBRecHitCollection->push_back(RecHitsCluster[j][Rec2]);
662 
663 
664  hMinvPi0EB_->Fill(m_inv);
665  hPt1Pi0EB_->Fill(etClus[i]);
666  hPt2Pi0EB_->Fill(etClus[j]);
667  hPtPi0EB_->Fill(pt_pair);
668  hIsoPi0EB_->Fill(Iso/pt_pair);
669  hS4S91Pi0EB_->Fill(s4s9Clus[i]);
670  hS4S92Pi0EB_->Fill(s4s9Clus[j]);
671 
672  // cout <<" EB Simple Clustering: pi0 Candidate pt, eta, phi, Iso, m_inv, i, j : "<<pt_pair<<" "<<pairVect.Eta()<<" "<<pairVect.Phi()<<" "<<Iso<<" "<<m_inv<<" "<<i<<" "<<j<<" "<<endl;
673 
674  npi0_s++;
675 
676  }
677 
678 
679 
680 
681  }
682  }
683  } // End of the "j" loop over Simple Clusters
684  } // End of the "i" loop over Simple Clusters
685 
686  // cout<<" (Simple Clustering) EB Pi0 candidates #: "<<npi0_s<<endl;
687 
688  } // rhEBpi0.valid() ends
689 
690  } // isMonEBpi0 ends
691 
692  //------------------ End of pi0 in EB --------------------------//
693 
694  // fill EB eta histos
695  if(isMonEBeta_ ){
696  if (rhEBeta.isValid() && (rhEBeta->size() > 0)){
697 
698 
699  const EcalRecHitCollection *hitCollection_p = rhEBeta.product();
700  float etot =0;
701  for(itb=rhEBeta->begin(); itb!=rhEBeta->end(); ++itb){
702 
703  EBDetId id(itb->id());
704  double energy = itb->energy();
705  if( energy < seleXtalMinEnergy_) continue;
706 
707  EBDetId det = itb->id();
708 
709 
710  detIdEBRecHits.push_back(det);
711  EBRecHits.push_back(*itb);
712 
713  if (energy > clusSeedThr_) seeds.push_back(*itb);
714 
715  hiPhiDistrEBeta_->Fill(id.iphi());
716  hiEtaDistrEBeta_->Fill(id.ieta());
717  hRechitEnergyEBeta_->Fill(itb->energy());
718 
719  etot+= itb->energy();
720  } // Eb rechits
721 
722  hNRecHitsEBeta_->Fill(rhEBeta->size());
723  hMeanRecHitEnergyEBeta_->Fill(etot/rhEBeta->size());
724  hEventEnergyEBeta_->Fill(etot);
725 
726  // cout << " EB RH Eta collection: #, mean rh_e, event E "<<rhEBeta->size()<<" "<<etot/rhEBeta->size()<<" "<<etot<<endl;
727 
728 
729  // Eta maker
730 
731  //cout<< " RH coll size: "<<rhEBeta->size()<<endl;
732  //cout<< " Eta seeds: "<<seeds.size()<<endl;
733 
734  int nClus;
735  vector<float> eClus;
736  vector<float> etClus;
737  vector<float> etaClus;
738  vector<float> thetaClus;
739  vector<float> phiClus;
740  vector<EBDetId> max_hit;
741 
742  vector< vector<EcalRecHit> > RecHitsCluster;
743  vector< vector<EcalRecHit> > RecHitsCluster5x5;
744  vector<float> s4s9Clus;
745  vector<float> s9s25Clus;
746 
747 
748  nClus=0;
749 
750  // Make own simple clusters (3x3, 5x5 or clusPhiSize_ x clusEtaSize_)
751  sort(seeds.begin(), seeds.end(), ecalRecHitLess());
752 
753  for (std::vector<EcalRecHit>::iterator itseed=seeds.begin(); itseed!=seeds.end(); itseed++) {
754  EBDetId seed_id = itseed->id();
755  std::vector<EBDetId>::const_iterator usedIds;
756 
757  bool seedAlreadyUsed=false;
758  for(usedIds=usedXtals.begin(); usedIds!=usedXtals.end(); usedIds++){
759  if(*usedIds==seed_id){
760  seedAlreadyUsed=true;
761  //cout<< " Seed with energy "<<itseed->energy()<<" was used !"<<endl;
762  break;
763  }
764  }
765  if(seedAlreadyUsed)continue;
766  std::vector<DetId> clus_v = topology_p->getWindow(seed_id,clusEtaSize_,clusPhiSize_);
767  std::vector<std::pair<DetId,float> > clus_used;
768  //Reject the seed if not able to build the cluster around it correctly
769  //if(clus_v.size() < clusEtaSize_*clusPhiSize_){cout<<" Not enough RecHits "<<endl; continue;}
770  vector<EcalRecHit> RecHitsInWindow;
771  vector<EcalRecHit> RecHitsInWindow5x5;
772 
773  double simple_energy = 0;
774 
775  for (std::vector<DetId>::iterator det=clus_v.begin(); det!=clus_v.end(); det++) {
776  EBDetId EBdet = *det;
777  // cout<<" det "<< EBdet<<" ieta "<<EBdet.ieta()<<" iphi "<<EBdet.iphi()<<endl;
778  bool HitAlreadyUsed=false;
779  for(usedIds=usedXtals.begin(); usedIds!=usedXtals.end(); usedIds++){
780  if(*usedIds==*det){
781  HitAlreadyUsed=true;
782  break;
783  }
784  }
785  if(HitAlreadyUsed)continue;
786 
787 
788  std::vector<EBDetId>::iterator itdet = find( detIdEBRecHits.begin(),detIdEBRecHits.end(),EBdet);
789  if(itdet == detIdEBRecHits.end()) continue;
790 
791  int nn = int(itdet - detIdEBRecHits.begin());
792  usedXtals.push_back(*det);
793  RecHitsInWindow.push_back(EBRecHits[nn]);
794  clus_used.push_back(std::make_pair(*det,1));
795  simple_energy = simple_energy + EBRecHits[nn].energy();
796 
797 
798  }
799 
800  if(simple_energy <= 0) continue;
801 
802  math::XYZPoint clus_pos = posCalculator_.Calculate_Location(clus_used,hitCollection_p,geometry_p,geometryES_p);
803  //cout<< " Simple Clustering: Total energy for this simple cluster : "<<simple_energy<<endl;
804  //cout<< " Simple Clustering: eta phi : "<<clus_pos.eta()<<" "<<clus_pos.phi()<<endl;
805  //cout<< " Simple Clustering: x y z : "<<clus_pos.x()<<" "<<clus_pos.y()<<" "<<clus_pos.z()<<endl;
806 
807  float theta_s = 2. * atan(exp(-clus_pos.eta()));
808  // float p0x_s = simple_energy * sin(theta_s) * cos(clus_pos.phi());
809  //float p0y_s = simple_energy * sin(theta_s) * sin(clus_pos.phi());
810  // float p0z_s = simple_energy * cos(theta_s);
811  //float et_s = sqrt( p0x_s*p0x_s + p0y_s*p0y_s);
812 
813  float et_s = simple_energy * sin(theta_s);
814  //cout << " Simple Clustering: E,Et,px,py,pz: "<<simple_energy<<" "<<et_s<<" "<<p0x_s<<" "<<p0y_s<<" "<<endl;
815 
816  //Compute S4/S9 variable
817  //We are not sure to have 9 RecHits so need to check eta and phi:
819  float s4s9_tmp[4];
820  for(int i=0;i<4;i++)s4s9_tmp[i]= 0;
821 
822  int seed_ieta = seed_id.ieta();
823  int seed_iphi = seed_id.iphi();
824 
825  convxtalid( seed_iphi,seed_ieta);
826 
827  float e3x3 = 0;
828  float e5x5 = 0;
829 
830  for(unsigned int j=0; j<RecHitsInWindow.size();j++){
831  EBDetId det = (EBDetId)RecHitsInWindow[j].id();
832 
833  int ieta = det.ieta();
834  int iphi = det.iphi();
835 
836  convxtalid(iphi,ieta);
837 
838  float en = RecHitsInWindow[j].energy();
839 
840  int dx = diff_neta_s(seed_ieta,ieta);
841  int dy = diff_nphi_s(seed_iphi,iphi);
842 
843  if(dx <= 0 && dy <=0) s4s9_tmp[0] += en;
844  if(dx >= 0 && dy <=0) s4s9_tmp[1] += en;
845  if(dx <= 0 && dy >=0) s4s9_tmp[2] += en;
846  if(dx >= 0 && dy >=0) s4s9_tmp[3] += en;
847 
848  if(std::abs(dx)<=1 && std::abs(dy)<=1) e3x3 += en;
849  if(std::abs(dx)<=2 && std::abs(dy)<=2) e5x5 += en;
850 
851 
852  }
853 
854  if(e3x3 <= 0) continue;
855 
856  float s4s9_max = *max_element( s4s9_tmp,s4s9_tmp+4)/e3x3;
857 
858 
860  std::vector<DetId> clus_v5x5 = topology_p->getWindow(seed_id,5,5);
861  for( std::vector<DetId>::const_iterator idItr = clus_v5x5.begin(); idItr != clus_v5x5.end(); idItr++){
862  EBDetId det = *idItr;
863 
864 
865  std::vector<EBDetId>::iterator itdet0 = find(usedXtals.begin(),usedXtals.end(),det);
866 
868  if(itdet0 != usedXtals.end()) continue;
869 
870  //inside collections
871  std::vector<EBDetId>::iterator itdet = find( detIdEBRecHits.begin(),detIdEBRecHits.end(),det);
872  if(itdet == detIdEBRecHits.end()) continue;
873 
874  int nn = int(itdet - detIdEBRecHits.begin());
875 
876  RecHitsInWindow5x5.push_back(EBRecHits[nn]);
877  e5x5 += EBRecHits[nn].energy();
878 
879  }
880 
881 
882 
883  if(e5x5 <= 0) continue;
884 
885  eClus.push_back(simple_energy);
886  etClus.push_back(et_s);
887  etaClus.push_back(clus_pos.eta());
888  thetaClus.push_back(theta_s);
889  phiClus.push_back(clus_pos.phi());
890  s4s9Clus.push_back(s4s9_max);
891  s9s25Clus.push_back(e3x3/e5x5);
892  RecHitsCluster.push_back(RecHitsInWindow);
893  RecHitsCluster5x5.push_back(RecHitsInWindow5x5);
894 
895  // std::cout<<" EB Eta cluster (n,nxt,e,et eta,phi,s4s9) "<<nClus<<" "<<int(RecHitsInWindow.size())<<" "<<eClus[nClus]<<" "<<" "<<etClus[nClus]<<" "<<etaClus[nClus]<<" "<<phiClus[nClus]<<" "<<s4s9Clus[nClus]<<std::endl;
896 
897  nClus++;
898 
899 
900  }
901 
902  // cout<< " Eta clusters: "<<nClus<<endl;
903 
904  // Selection, based on Simple clustering
905  //eta candidates
906  int npi0_s=0;
907 
908 
909  // if (nClus <= 1) return;
910  for(Int_t i=0 ; i<nClus ; i++){
911  for(Int_t j=i+1 ; j<nClus ; j++){
912  // cout<<" i "<<i<<" etClus[i] "<<etClus[i]<<" j "<<j<<" etClus[j] "<<etClus[j]<<endl;
913  if( etClus[i]>selePtGammaEta_ && etClus[j]>selePtGammaEta_ && s4s9Clus[i]>seleS4S9GammaEta_ && s4s9Clus[j]>seleS4S9GammaEta_){
914 
915 
916  float p0x = etClus[i] * cos(phiClus[i]);
917  float p1x = etClus[j] * cos(phiClus[j]);
918  float p0y = etClus[i] * sin(phiClus[i]);
919  float p1y = etClus[j] * sin(phiClus[j]);
920  float p0z = eClus[i] * cos(thetaClus[i]);
921  float p1z = eClus[j] * cos(thetaClus[j]);
922 
923 
924  float pt_pair = sqrt( (p0x+p1x)*(p0x+p1x) + (p0y+p1y)*(p0y+p1y));
925 
926  if (pt_pair < selePtEta_)continue;
927 
928  float m_inv = sqrt ( (eClus[i] + eClus[j])*(eClus[i] + eClus[j]) - (p0x+p1x)*(p0x+p1x) - (p0y+p1y)*(p0y+p1y) - (p0z+p1z)*(p0z+p1z) );
929  if ( (m_inv<seleMinvMaxEta_) && (m_inv>seleMinvMinEta_) ){
930 
931  //New Loop on cluster to measure isolation:
932  vector<int> IsoClus;
933  IsoClus.clear();
934  float Iso = 0;
935  TVector3 pairVect = TVector3((p0x+p1x), (p0y+p1y), (p0z+p1z));
936  for(Int_t k=0 ; k<nClus ; k++){
937 
938 
939  if(etClus[k] < ptMinForIsolationEta_) continue;
940 
941  if(k==i || k==j)continue;
942  TVector3 ClusVect = TVector3(etClus[k] *cos(phiClus[k]), etClus[k] * sin(phiClus[k]) , eClus[k] * cos(thetaClus[k]));
943 
944  float dretacl = fabs(etaClus[k] - pairVect.Eta());
945  float drcl = ClusVect.DeltaR(pairVect);
946  // cout<< " Iso: k, E, drclpi0, detaclpi0, dphiclpi0 "<<k<<" "<<eClus[k]<<" "<<drclpi0<<" "<<dretaclpi0<<endl;
947  if((drcl<seleEtaBeltDR_) && (dretacl<seleEtaBeltDeta_) ){
948  // cout<< " ... good iso cluster #: "<<k<<" etClus[k] "<<etClus[k] <<endl;
949  Iso = Iso + etClus[k];
950  IsoClus.push_back(k);
951  }
952  }
953 
954  // cout<<" Iso/pt_pi0 "<<Iso/pt_pi0<<endl;
955  if(Iso/pt_pair<seleEtaIso_){
956  //for(unsigned int Rec=0;Rec<RecHitsCluster[i].size();Rec++)pi0EBRecHitCollection->push_back(RecHitsCluster[i][Rec]);
957  //for(unsigned int Rec2=0;Rec2<RecHitsCluster[j].size();Rec2++)pi0EBRecHitCollection->push_back(RecHitsCluster[j][Rec2]);
958 
959 
960  hMinvEtaEB_->Fill(m_inv);
961  hPt1EtaEB_->Fill(etClus[i]);
962  hPt2EtaEB_->Fill(etClus[j]);
963  hPtEtaEB_->Fill(pt_pair);
964  hIsoEtaEB_->Fill(Iso/pt_pair);
965  hS4S91EtaEB_->Fill(s4s9Clus[i]);
966  hS4S92EtaEB_->Fill(s4s9Clus[j]);
967 
968  // cout <<" EB Simple Clustering: Eta Candidate pt, eta, phi, Iso, m_inv, i, j : "<<pt_pair<<" "<<pairVect.Eta()<<" "<<pairVect.Phi()<<" "<<Iso<<" "<<m_inv<<" "<<i<<" "<<j<<" "<<endl;
969 
970  npi0_s++;
971 
972  }
973 
974 
975 
976 
977  }
978  }
979  } // End of the "j" loop over Simple Clusters
980  } // End of the "i" loop over Simple Clusters
981 
982  // cout<<" (Simple Clustering) EB Eta candidates #: "<<npi0_s<<endl;
983 
984  } // rhEBeta.valid() ends
985 
986  } // isMonEBeta ends
987 
988  //------------------ End of Eta in EB --------------------------//
989 
990 
991 
992  //----------------- End of the EB --------------------------//
993 
994 
995 
996 
997  //----------------- Start of the EE --------------------//
998 
999  // fill pi0 EE histos
1000  if(isMonEEpi0_){
1001  if (rhEEpi0.isValid() && (rhEEpi0->size() > 0)){
1002 
1003  const EcalRecHitCollection *hitCollection_ee = rhEEpi0.product();
1004  float etot =0;
1005 
1006 
1007  detIdEERecHits.clear();
1008  EERecHits.clear();
1009 
1010 
1011  std::vector<EcalRecHit> seedsEndCap;
1012  seedsEndCap.clear();
1013 
1014  vector<EEDetId> usedXtalsEndCap;
1015  usedXtalsEndCap.clear();
1016 
1017 
1019  EERecHitCollection::const_iterator ite;
1020  for (ite=rhEEpi0->begin(); ite!=rhEEpi0->end(); ite++) {
1021  double energy = ite->energy();
1022  if( energy < seleXtalMinEnergyEndCap_) continue;
1023 
1024  EEDetId det = ite->id();
1025  EEDetId id(ite->id());
1026 
1027 
1028  detIdEERecHits.push_back(det);
1029  EERecHits.push_back(*ite);
1030 
1031  hiXDistrEEpi0_->Fill(id.ix());
1032  hiYDistrEEpi0_->Fill(id.iy());
1033  hRechitEnergyEEpi0_->Fill(ite->energy());
1034 
1035  if (energy > clusSeedThrEndCap_) seedsEndCap.push_back(*ite);
1036 
1037  etot+= ite->energy();
1038  } // EE rechits
1039 
1040  hNRecHitsEEpi0_->Fill(rhEEpi0->size());
1041  hMeanRecHitEnergyEEpi0_->Fill(etot/rhEEpi0->size());
1042  hEventEnergyEEpi0_->Fill(etot);
1043 
1044  // cout << " EE RH Pi0 collection: #, mean rh_e, event E "<<rhEEpi0->size()<<" "<<etot/rhEEpi0->size()<<" "<<etot<<endl;
1045 
1046  int nClusEndCap;
1047  vector<float> eClusEndCap;
1048  vector<float> etClusEndCap;
1049  vector<float> etaClusEndCap;
1050  vector<float> thetaClusEndCap;
1051  vector<float> phiClusEndCap;
1052  vector< vector<EcalRecHit> > RecHitsClusterEndCap;
1053  vector< vector<EcalRecHit> > RecHitsCluster5x5EndCap;
1054  vector<float> s4s9ClusEndCap;
1055  vector<float> s9s25ClusEndCap;
1056 
1057 
1058  nClusEndCap=0;
1059 
1060 
1061  // Make own simple clusters (3x3, 5x5 or clusPhiSize_ x clusEtaSize_)
1062  sort(seedsEndCap.begin(), seedsEndCap.end(), ecalRecHitLess());
1063 
1064  for (std::vector<EcalRecHit>::iterator itseed=seedsEndCap.begin(); itseed!=seedsEndCap.end(); itseed++) {
1065  EEDetId seed_id = itseed->id();
1066  std::vector<EEDetId>::const_iterator usedIds;
1067 
1068  bool seedAlreadyUsed=false;
1069  for(usedIds=usedXtalsEndCap.begin(); usedIds!=usedXtalsEndCap.end(); usedIds++){
1070  if(*usedIds==seed_id){
1071  seedAlreadyUsed=true;
1072  break;
1073  }
1074  }
1075 
1076  if(seedAlreadyUsed)continue;
1077  std::vector<DetId> clus_v = topology_ee->getWindow(seed_id,clusEtaSize_,clusPhiSize_);
1078  std::vector<std::pair<DetId,float> > clus_used;
1079 
1080 
1081  vector<EcalRecHit> RecHitsInWindow;
1082  vector<EcalRecHit> RecHitsInWindow5x5;
1083 
1084  float simple_energy = 0;
1085 
1086  for (std::vector<DetId>::iterator det=clus_v.begin(); det!=clus_v.end(); det++) {
1087  EEDetId EEdet = *det;
1088 
1089  bool HitAlreadyUsed=false;
1090  for(usedIds=usedXtalsEndCap.begin(); usedIds!=usedXtalsEndCap.end(); usedIds++){
1091  if(*usedIds==*det){
1092  HitAlreadyUsed=true;
1093  break;
1094  }
1095  }
1096 
1097  if(HitAlreadyUsed)continue;
1098 
1099  std::vector<EEDetId>::iterator itdet = find( detIdEERecHits.begin(),detIdEERecHits.end(),EEdet);
1100  if(itdet == detIdEERecHits.end()) continue;
1101 
1102 
1103  int nn = int(itdet - detIdEERecHits.begin());
1104  usedXtalsEndCap.push_back(*det);
1105  RecHitsInWindow.push_back(EERecHits[nn]);
1106  clus_used.push_back(std::make_pair(*det,1));
1107  simple_energy = simple_energy + EERecHits[nn].energy();
1108 
1109 
1110  }
1111 
1112  if( simple_energy <= 0) continue;
1113 
1114  math::XYZPoint clus_pos = posCalculator_.Calculate_Location(clus_used,hitCollection_ee,geometryEE_p,geometryES_p);
1115 
1116 
1117  float theta_s = 2. * atan(exp(-clus_pos.eta()));
1118  float et_s = simple_energy * sin(theta_s);
1119  // float p0x_s = simple_energy * sin(theta_s) * cos(clus_pos.phi());
1120  //float p0y_s = simple_energy * sin(theta_s) * sin(clus_pos.phi());
1121  //float et_s = sqrt( p0x_s*p0x_s + p0y_s*p0y_s);
1122 
1123 
1124  //Compute S4/S9 variable
1125  //We are not sure to have 9 RecHits so need to check eta and phi:
1126  float s4s9_tmp[4];
1127  for(int i=0;i<4;i++) s4s9_tmp[i]= 0;
1128 
1129  int ixSeed = seed_id.ix();
1130  int iySeed = seed_id.iy();
1131  float e3x3 = 0;
1132  float e5x5 = 0;
1133 
1134  for(unsigned int j=0; j<RecHitsInWindow.size();j++){
1135  EEDetId det_this = (EEDetId)RecHitsInWindow[j].id();
1136  int dx = ixSeed - det_this.ix();
1137  int dy = iySeed - det_this.iy();
1138 
1139  float en = RecHitsInWindow[j].energy();
1140 
1141  if(dx <= 0 && dy <=0) s4s9_tmp[0] += en;
1142  if(dx >= 0 && dy <=0) s4s9_tmp[1] += en;
1143  if(dx <= 0 && dy >=0) s4s9_tmp[2] += en;
1144  if(dx >= 0 && dy >=0) s4s9_tmp[3] += en;
1145 
1146  if( std::abs(dx)<=1 && std::abs(dy)<=1) e3x3 += en;
1147  if( std::abs(dx)<=2 && std::abs(dy)<=2) e5x5 += en;
1148 
1149  }
1150 
1151 
1152  if(e3x3 <= 0) continue;
1153 
1154  eClusEndCap.push_back(simple_energy);
1155  etClusEndCap.push_back(et_s);
1156  etaClusEndCap.push_back(clus_pos.eta());
1157  thetaClusEndCap.push_back(theta_s);
1158  phiClusEndCap.push_back(clus_pos.phi());
1159  s4s9ClusEndCap.push_back(*max_element( s4s9_tmp,s4s9_tmp+4)/e3x3);
1160  s9s25ClusEndCap.push_back(e3x3/e5x5);
1161  RecHitsClusterEndCap.push_back(RecHitsInWindow);
1162  RecHitsCluster5x5EndCap.push_back(RecHitsInWindow5x5);
1163 
1164 
1165  // std::cout<<" EE pi0 cluster (n,nxt,e,et eta,phi,s4s9) "<<nClusEndCap<<" "<<int(RecHitsInWindow.size())<<" "<<eClusEndCap[nClusEndCap]<<" "<<" "<<etClusEndCap[nClusEndCap]<<" "<<etaClusEndCap[nClusEndCap]<<" "<<phiClusEndCap[nClusEndCap]<<" "<<s4s9ClusEndCap[nClusEndCap]<<std::endl;
1166 
1167 
1168  nClusEndCap++;
1169 
1170 
1171  }
1172 
1173 
1174 
1175  // Selection, based on Simple clustering
1176  //pi0 candidates
1177  int npi0_se=0;
1178 
1179 
1180  for(Int_t i=0 ; i<nClusEndCap ; i++){
1181  for(Int_t j=i+1 ; j<nClusEndCap ; j++){
1182 
1183  if( etClusEndCap[i]>selePtGammaEndCap_ && etClusEndCap[j]>selePtGammaEndCap_ && s4s9ClusEndCap[i]>seleS4S9GammaEndCap_ && s4s9ClusEndCap[j]>seleS4S9GammaEndCap_){
1184 
1185  float p0x = etClusEndCap[i] * cos(phiClusEndCap[i]);
1186  float p1x = etClusEndCap[j] * cos(phiClusEndCap[j]);
1187  float p0y = etClusEndCap[i] * sin(phiClusEndCap[i]);
1188  float p1y = etClusEndCap[j] * sin(phiClusEndCap[j]);
1189  float p0z = eClusEndCap[i] * cos(thetaClusEndCap[i]);
1190  float p1z = eClusEndCap[j] * cos(thetaClusEndCap[j]);
1191 
1192 
1193  float pt_pair = sqrt( (p0x+p1x)*(p0x+p1x) + (p0y+p1y)*(p0y+p1y));
1194  if (pt_pair < selePtPi0EndCap_)continue;
1195  float m_inv = sqrt ( (eClusEndCap[i] + eClusEndCap[j])*(eClusEndCap[i] + eClusEndCap[j]) - (p0x+p1x)*(p0x+p1x) - (p0y+p1y)*(p0y+p1y) - (p0z+p1z)*(p0z+p1z) );
1196 
1197 
1198  if ( (m_inv<seleMinvMaxPi0EndCap_) && (m_inv>seleMinvMinPi0EndCap_) ){
1199 
1200 
1201  //New Loop on cluster to measure isolation:
1202  vector<int> IsoClus;
1203  IsoClus.clear();
1204  float Iso = 0;
1205  TVector3 pairVect = TVector3((p0x+p1x), (p0y+p1y), (p0z+p1z));
1206  for(Int_t k=0 ; k<nClusEndCap ; k++){
1207 
1208  if(etClusEndCap[k] < ptMinForIsolationEndCap_) continue;
1209 
1210 
1211  if(k==i || k==j)continue;
1212 
1213 
1214  TVector3 clusVect = TVector3(etClusEndCap[k] * cos(phiClusEndCap[k]), etClusEndCap[k] * sin(phiClusEndCap[k]) , eClusEndCap[k] * cos(thetaClusEndCap[k]) ) ;
1215  float dretacl = fabs(etaClusEndCap[k] - pairVect.Eta());
1216  float drcl = clusVect.DeltaR(pairVect);
1217 
1218  if(drcl < selePi0BeltDREndCap_ && dretacl < selePi0BeltDetaEndCap_ ){
1219  Iso = Iso + etClusEndCap[k];
1220  IsoClus.push_back(k);
1221  }
1222  }
1223 
1224 
1225  if(Iso/pt_pair<selePi0IsoEndCap_){
1226  //cout <<" EE Simple Clustering: pi0 Candidate pt, eta, phi, Iso, m_inv, i, j : "<<pt_pair<<" "<<pairVect.Eta()<<" "<<pairVect.Phi()<<" "<<Iso<<" "<<m_inv<<" "<<i<<" "<<j<<" "<<endl;
1227 
1228  hMinvPi0EE_->Fill(m_inv);
1229  hPt1Pi0EE_->Fill(etClusEndCap[i]);
1230  hPt2Pi0EE_->Fill(etClusEndCap[j]);
1231  hPtPi0EE_->Fill(pt_pair);
1232  hIsoPi0EE_->Fill(Iso/pt_pair);
1233  hS4S91Pi0EE_->Fill(s4s9ClusEndCap[i]);
1234  hS4S92Pi0EE_->Fill(s4s9ClusEndCap[j]);
1235 
1236  npi0_se++;
1237  }
1238 
1239  }
1240  }
1241  } // End of the "j" loop over Simple Clusters
1242  } // End of the "i" loop over Simple Clusters
1243 
1244  // cout<<" (Simple Clustering) EE Pi0 candidates #: "<<npi0_se<<endl;
1245 
1246  } //rhEEpi0
1247  } //isMonEEpi0
1248 
1249  //================End of Pi0 endcap=======================//
1250 
1251 
1252  //================ Eta in EE===============================//
1253 
1254  // fill pi0 EE histos
1255  if(isMonEEeta_){
1256  if (rhEEeta.isValid() && (rhEEeta->size() > 0)){
1257 
1258  const EcalRecHitCollection *hitCollection_ee = rhEEeta.product();
1259  float etot =0;
1260 
1261  detIdEERecHits.clear();
1262  EERecHits.clear();
1263 
1264 
1265  std::vector<EcalRecHit> seedsEndCap;
1266  seedsEndCap.clear();
1267 
1268  vector<EEDetId> usedXtalsEndCap;
1269  usedXtalsEndCap.clear();
1270 
1271 
1273  EERecHitCollection::const_iterator ite;
1274  for (ite=rhEEeta->begin(); ite!=rhEEeta->end(); ite++) {
1275  double energy = ite->energy();
1276  if( energy < seleXtalMinEnergyEndCap_) continue;
1277 
1278  EEDetId det = ite->id();
1279  EEDetId id(ite->id());
1280 
1281 
1282  detIdEERecHits.push_back(det);
1283  EERecHits.push_back(*ite);
1284 
1285  hiXDistrEEeta_->Fill(id.ix());
1286  hiYDistrEEeta_->Fill(id.iy());
1287  hRechitEnergyEEeta_->Fill(ite->energy());
1288 
1289  if (energy > clusSeedThrEndCap_) seedsEndCap.push_back(*ite);
1290 
1291  etot+= ite->energy();
1292  } // EE rechits
1293 
1294  hNRecHitsEEeta_->Fill(rhEEeta->size());
1295  hMeanRecHitEnergyEEeta_->Fill(etot/rhEEeta->size());
1296  hEventEnergyEEeta_->Fill(etot);
1297 
1298  // cout << " EE RH Eta collection: #, mean rh_e, event E "<<rhEEeta->size()<<" "<<etot/rhEEeta->size()<<" "<<etot<<endl;
1299 
1300  int nClusEndCap;
1301  vector<float> eClusEndCap;
1302  vector<float> etClusEndCap;
1303  vector<float> etaClusEndCap;
1304  vector<float> thetaClusEndCap;
1305  vector<float> phiClusEndCap;
1306  vector< vector<EcalRecHit> > RecHitsClusterEndCap;
1307  vector< vector<EcalRecHit> > RecHitsCluster5x5EndCap;
1308  vector<float> s4s9ClusEndCap;
1309  vector<float> s9s25ClusEndCap;
1310 
1311 
1312  nClusEndCap=0;
1313 
1314 
1315  // Make own simple clusters (3x3, 5x5 or clusPhiSize_ x clusEtaSize_)
1316  sort(seedsEndCap.begin(), seedsEndCap.end(), ecalRecHitLess());
1317 
1318  for (std::vector<EcalRecHit>::iterator itseed=seedsEndCap.begin(); itseed!=seedsEndCap.end(); itseed++) {
1319  EEDetId seed_id = itseed->id();
1320  std::vector<EEDetId>::const_iterator usedIds;
1321 
1322  bool seedAlreadyUsed=false;
1323  for(usedIds=usedXtalsEndCap.begin(); usedIds!=usedXtalsEndCap.end(); usedIds++){
1324  if(*usedIds==seed_id){
1325  seedAlreadyUsed=true;
1326  break;
1327  }
1328  }
1329 
1330  if(seedAlreadyUsed)continue;
1331  std::vector<DetId> clus_v = topology_ee->getWindow(seed_id,clusEtaSize_,clusPhiSize_);
1332  std::vector<std::pair<DetId,float> > clus_used;
1333 
1334 
1335  vector<EcalRecHit> RecHitsInWindow;
1336  vector<EcalRecHit> RecHitsInWindow5x5;
1337 
1338  float simple_energy = 0;
1339 
1340  for (std::vector<DetId>::iterator det=clus_v.begin(); det!=clus_v.end(); det++) {
1341  EEDetId EEdet = *det;
1342 
1343  bool HitAlreadyUsed=false;
1344  for(usedIds=usedXtalsEndCap.begin(); usedIds!=usedXtalsEndCap.end(); usedIds++){
1345  if(*usedIds==*det){
1346  HitAlreadyUsed=true;
1347  break;
1348  }
1349  }
1350 
1351  if(HitAlreadyUsed)continue;
1352 
1353  std::vector<EEDetId>::iterator itdet = find( detIdEERecHits.begin(),detIdEERecHits.end(),EEdet);
1354  if(itdet == detIdEERecHits.end()) continue;
1355 
1356 
1357  int nn = int(itdet - detIdEERecHits.begin());
1358  usedXtalsEndCap.push_back(*det);
1359  RecHitsInWindow.push_back(EERecHits[nn]);
1360  clus_used.push_back(std::make_pair(*det,1));
1361  simple_energy = simple_energy + EERecHits[nn].energy();
1362 
1363 
1364  }
1365 
1366  if( simple_energy <= 0) continue;
1367 
1368  math::XYZPoint clus_pos = posCalculator_.Calculate_Location(clus_used,hitCollection_ee,geometryEE_p,geometryES_p);
1369 
1370 
1371  float theta_s = 2. * atan(exp(-clus_pos.eta()));
1372  float et_s = simple_energy * sin(theta_s);
1373  // float p0x_s = simple_energy * sin(theta_s) * cos(clus_pos.phi());
1374  //float p0y_s = simple_energy * sin(theta_s) * sin(clus_pos.phi());
1375  //float et_s = sqrt( p0x_s*p0x_s + p0y_s*p0y_s);
1376 
1377 
1378  //Compute S4/S9 variable
1379  //We are not sure to have 9 RecHits so need to check eta and phi:
1380  float s4s9_tmp[4];
1381  for(int i=0;i<4;i++) s4s9_tmp[i]= 0;
1382 
1383  int ixSeed = seed_id.ix();
1384  int iySeed = seed_id.iy();
1385  float e3x3 = 0;
1386  float e5x5 = 0;
1387 
1388  for(unsigned int j=0; j<RecHitsInWindow.size();j++){
1389  EEDetId det_this = (EEDetId)RecHitsInWindow[j].id();
1390  int dx = ixSeed - det_this.ix();
1391  int dy = iySeed - det_this.iy();
1392 
1393  float en = RecHitsInWindow[j].energy();
1394 
1395  if(dx <= 0 && dy <=0) s4s9_tmp[0] += en;
1396  if(dx >= 0 && dy <=0) s4s9_tmp[1] += en;
1397  if(dx <= 0 && dy >=0) s4s9_tmp[2] += en;
1398  if(dx >= 0 && dy >=0) s4s9_tmp[3] += en;
1399 
1400  if( std::abs(dx)<=1 && std::abs(dy)<=1) e3x3 += en;
1401  if( std::abs(dx)<=2 && std::abs(dy)<=2) e5x5 += en;
1402 
1403  }
1404 
1405 
1406  if(e3x3 <= 0) continue;
1407 
1408  eClusEndCap.push_back(simple_energy);
1409  etClusEndCap.push_back(et_s);
1410  etaClusEndCap.push_back(clus_pos.eta());
1411  thetaClusEndCap.push_back(theta_s);
1412  phiClusEndCap.push_back(clus_pos.phi());
1413  s4s9ClusEndCap.push_back(*max_element( s4s9_tmp,s4s9_tmp+4)/e3x3);
1414  s9s25ClusEndCap.push_back(e3x3/e5x5);
1415  RecHitsClusterEndCap.push_back(RecHitsInWindow);
1416  RecHitsCluster5x5EndCap.push_back(RecHitsInWindow5x5);
1417 
1418  // std::cout<<" EE Eta cluster (n,nxt,e,et eta,phi,s4s9) "<<nClusEndCap<<" "<<int(RecHitsInWindow.size())<<" "<<eClusEndCap[nClusEndCap]<<" "<<" "<<etClusEndCap[nClusEndCap]<<" "<<etaClusEndCap[nClusEndCap]<<" "<<phiClusEndCap[nClusEndCap]<<" "<<s4s9ClusEndCap[nClusEndCap]<<std::endl;
1419 
1420  nClusEndCap++;
1421 
1422  }
1423 
1424 
1425 
1426  // Selection, based on Simple clustering
1427  //pi0 candidates
1428  int npi0_se=0;
1429 
1430 
1431  for(Int_t i=0 ; i<nClusEndCap ; i++){
1432  for(Int_t j=i+1 ; j<nClusEndCap ; j++){
1433 
1434  if( etClusEndCap[i]>selePtGammaEtaEndCap_ && etClusEndCap[j]>selePtGammaEtaEndCap_ && s4s9ClusEndCap[i]>seleS4S9GammaEtaEndCap_ && s4s9ClusEndCap[j]>seleS4S9GammaEtaEndCap_){
1435 
1436  float p0x = etClusEndCap[i] * cos(phiClusEndCap[i]);
1437  float p1x = etClusEndCap[j] * cos(phiClusEndCap[j]);
1438  float p0y = etClusEndCap[i] * sin(phiClusEndCap[i]);
1439  float p1y = etClusEndCap[j] * sin(phiClusEndCap[j]);
1440  float p0z = eClusEndCap[i] * cos(thetaClusEndCap[i]);
1441  float p1z = eClusEndCap[j] * cos(thetaClusEndCap[j]);
1442 
1443 
1444  float pt_pair = sqrt( (p0x+p1x)*(p0x+p1x) + (p0y+p1y)*(p0y+p1y));
1445  if (pt_pair < selePtEtaEndCap_)continue;
1446  float m_inv = sqrt ( (eClusEndCap[i] + eClusEndCap[j])*(eClusEndCap[i] + eClusEndCap[j]) - (p0x+p1x)*(p0x+p1x) - (p0y+p1y)*(p0y+p1y) - (p0z+p1z)*(p0z+p1z) );
1447 
1448 
1449  if ( (m_inv<seleMinvMaxEtaEndCap_) && (m_inv>seleMinvMinEtaEndCap_) ){
1450 
1451 
1452  //New Loop on cluster to measure isolation:
1453  vector<int> IsoClus;
1454  IsoClus.clear();
1455  float Iso = 0;
1456  TVector3 pairVect = TVector3((p0x+p1x), (p0y+p1y), (p0z+p1z));
1457  for(Int_t k=0 ; k<nClusEndCap ; k++){
1458 
1459  if(etClusEndCap[k] < ptMinForIsolationEtaEndCap_) continue;
1460 
1461 
1462  if(k==i || k==j)continue;
1463 
1464 
1465  TVector3 clusVect = TVector3(etClusEndCap[k] * cos(phiClusEndCap[k]), etClusEndCap[k] * sin(phiClusEndCap[k]) , eClusEndCap[k] * cos(thetaClusEndCap[k]) ) ;
1466  float dretacl = fabs(etaClusEndCap[k] - pairVect.Eta());
1467  float drcl = clusVect.DeltaR(pairVect);
1468 
1469  if(drcl < seleEtaBeltDREndCap_ && dretacl < seleEtaBeltDetaEndCap_ ){
1470  Iso = Iso + etClusEndCap[k];
1471  IsoClus.push_back(k);
1472  }
1473  }
1474 
1475 
1476  if(Iso/pt_pair<seleEtaIsoEndCap_){
1477  // cout <<" EE Simple Clustering: Eta Candidate pt, eta, phi, Iso, m_inv, i, j : "<<pt_pair<<" "<<pairVect.Eta()<<" "<<pairVect.Phi()<<" "<<Iso<<" "<<m_inv<<" "<<i<<" "<<j<<" "<<endl;
1478 
1479  hMinvEtaEE_->Fill(m_inv);
1480  hPt1EtaEE_->Fill(etClusEndCap[i]);
1481  hPt2EtaEE_->Fill(etClusEndCap[j]);
1482  hPtEtaEE_->Fill(pt_pair);
1483  hIsoEtaEE_->Fill(Iso/pt_pair);
1484  hS4S91EtaEE_->Fill(s4s9ClusEndCap[i]);
1485  hS4S92EtaEE_->Fill(s4s9ClusEndCap[j]);
1486 
1487  npi0_se++;
1488  }
1489 
1490  }
1491  }
1492  } // End of the "j" loop over Simple Clusters
1493  } // End of the "i" loop over Simple Clusters
1494 
1495  // cout<<" (Simple Clustering) EE Eta candidates #: "<<npi0_se<<endl;
1496 
1497  } //rhEEeta
1498  } //isMonEEeta
1499 
1500  //================End of Pi0 endcap=======================//
1501 
1502 
1503 
1504 
1506 
1507 
1508 }
1509 
1510 
1511 
1512 
1513 //--------------------------------------------------------
1514 void DQMSourcePi0::endLuminosityBlock(const LuminosityBlock& lumiSeg,
1515  const EventSetup& context) {
1516 }
1517 //--------------------------------------------------------
1518 void DQMSourcePi0::endRun(const Run& r, const EventSetup& context){
1519 
1520 }
1521 //--------------------------------------------------------
1522 void DQMSourcePi0::endJob(){
1523 
1524 // if(dbe_) {
1525 // if (saveToFile_) {
1526 // dbe_->save(fileName_);
1527 // }
1528 // }
1529 }
1530 
1531 
1532 void DQMSourcePi0::convxtalid(Int_t &nphi,Int_t &neta)
1533 {
1534  // Barrel only
1535  // Output nphi 0...359; neta 0...84; nside=+1 (for eta>0), or 0 (for eta<0).
1536  // neta will be [-85,-1] , or [0,84], the minus sign indicates the z<0 side.
1537 
1538  if(neta > 0) neta -= 1;
1539  if(nphi > 359) nphi=nphi-360;
1540 
1541 } //end of convxtalid
1542 
1543 int DQMSourcePi0::diff_neta_s(Int_t neta1, Int_t neta2){
1544  Int_t mdiff;
1545  mdiff=(neta1-neta2);
1546  return mdiff;
1547 }
1548 
1549 // Calculate the distance in xtals taking into account the periodicity of the Barrel
1550 int DQMSourcePi0::diff_nphi_s(Int_t nphi1,Int_t nphi2) {
1551  Int_t mdiff;
1552  if(std::abs(nphi1-nphi2) < (360-std::abs(nphi1-nphi2))) {
1553  mdiff=nphi1-nphi2;
1554  }
1555  else {
1556  mdiff=360-std::abs(nphi1-nphi2);
1557  if(nphi1>nphi2) mdiff=-mdiff;
1558  }
1559  return mdiff;
1560 }
1561 
1562 
DQMSourcePi0::productMonitoredEEeta_
edm::EDGetTokenT< EcalRecHitCollection > productMonitoredEEeta_
Definition: DQMSourcePi0.h:248
DQMSourcePi0::seleEtaBeltDREndCap_
double seleEtaBeltDREndCap_
Definition: DQMSourcePi0.h:305
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
DQMSourcePi0::hiYDistrEEeta_
MonitorElement * hiYDistrEEeta_
Distribution of rechits in iy EE (eta)
Definition: DQMSourcePi0.h:103
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
i
int i
Definition: DBlmapReader.cc:9
DQMSourcePi0::hMinvPi0EE_
MonitorElement * hMinvPi0EE_
Pi0 invariant mass in EE.
Definition: DQMSourcePi0.h:157
DQMSourcePi0::seleMinvMinPi0EndCap_
double seleMinvMinPi0EndCap_
Definition: DQMSourcePi0.h:277
DQMSourcePi0::hiXDistrEEeta_
MonitorElement * hiXDistrEEeta_
Distribution of rechits in ix EE (eta)
Definition: DQMSourcePi0.h:91
DQMSourcePi0::hiPhiDistrEBpi0_
MonitorElement * hiPhiDistrEBpi0_
Distribution of rechits in iPhi (pi0)
Definition: DQMSourcePi0.h:82
create_public_lumi_plots.exp
tuple exp
Definition: create_public_lumi_plots.py:1088
DQMSourcePi0::hPt1Pi0EB_
MonitorElement * hPt1Pi0EB_
Pt of the 1st most energetic Pi0 photon in EB.
Definition: DQMSourcePi0.h:166
EEDetId::ix
int ix() const
Definition: EEDetId.h:76
DQMSourcePi0::hPt2Pi0EE_
MonitorElement * hPt2Pi0EE_
Pt of the 2nd most energetic Pi0 photon in EE.
Definition: DQMSourcePi0.h:182
DQMSourcePi0::hiYDistrEEpi0_
MonitorElement * hiYDistrEEpi0_
Distribution of rechits in iy EE (pi0)
Definition: DQMSourcePi0.h:97
DQMSourcePi0::ptMinForIsolationEtaEndCap_
double ptMinForIsolationEtaEndCap_
Definition: DQMSourcePi0.h:303
DQMSourcePi0::hPtPi0EE_
MonitorElement * hPtPi0EE_
Pi0 Pt in EE.
Definition: DQMSourcePi0.h:195
DQMSourcePi0::hMinvEtaEB_
MonitorElement * hMinvEtaEB_
Eta invariant mass in EB.
Definition: DQMSourcePi0.h:160
DQMSourcePi0::seleEtaBeltDR_
double seleEtaBeltDR_
Definition: DQMSourcePi0.h:293
edm::Event::getByToken
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:462
DQMSourcePi0::detIdEERecHits
std::vector< EEDetId > detIdEERecHits
Definition: DQMSourcePi0.h:322
DQMSourcePi0::endRun
void endRun(const edm::Run &r, const edm::EventSetup &c)
Definition: DQMSourcePi0.cc:1518
DQMSourcePi0::bookHistograms
virtual void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
Definition: DQMSourcePi0.cc:145
DQMSourcePi0::posCalculator_
PositionCalc posCalculator_
Definition: DQMSourcePi0.h:79
DQMSourcePi0::seleMinvMinEtaEndCap_
double seleMinvMinEtaEndCap_
Definition: DQMSourcePi0.h:302
DQMSourcePi0::hS4S92EtaEB_
MonitorElement * hS4S92EtaEB_
S4S9 of the 2nd most energetic eta photon.
Definition: DQMSourcePi0.h:234
DQMSourcePi0::detIdEBRecHits
std::vector< EBDetId > detIdEBRecHits
Definition: DQMSourcePi0.h:318
relval_2017.k
list k
Definition: relval_2017.py:23
DQMSourcePi0::productMonitoredEBeta_
edm::EDGetTokenT< EcalRecHitCollection > productMonitoredEBeta_
Definition: DQMSourcePi0.h:244
DQMSourcePi0::selePtGammaEndCap_
double selePtGammaEndCap_
for pi0-&gt;gg endcap
Definition: DQMSourcePi0.h:274
DQMSourcePi0::hMeanRecHitEnergyEBeta_
MonitorElement * hMeanRecHitEnergyEBeta_
Distribution of Mean energy per rechit EB (eta)
Definition: DQMSourcePi0.h:148
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
DQMSourcePi0::seleMinvMaxPi0EndCap_
double seleMinvMaxPi0EndCap_
Definition: DQMSourcePi0.h:276
DQMSourcePi0::folderName_
std::string folderName_
DQM folder name.
Definition: DQMSourcePi0.h:331
edm::SortedCollection< EcalRecHit >::const_iterator
std::vector< EcalRecHit >::const_iterator const_iterator
Definition: SortedCollection.h:82
DQMSourcePi0::hMeanRecHitEnergyEEpi0_
MonitorElement * hMeanRecHitEnergyEEpi0_
Distribution of Mean energy per rechit EE (pi0)
Definition: DQMSourcePi0.h:145
DQMSourcePi0::hNRecHitsEBpi0_
MonitorElement * hNRecHitsEBpi0_
Distribution of number of RecHits EB (pi0)
Definition: DQMSourcePi0.h:130
DQMSourcePi0::hPt2EtaEE_
MonitorElement * hPt2EtaEE_
Pt of the 2nd most energetic Eta photon in EE.
Definition: DQMSourcePi0.h:188
DQMSourcePi0::hRechitEnergyEBpi0_
MonitorElement * hRechitEnergyEBpi0_
Energy Distribution of rechits EB (pi0)
Definition: DQMSourcePi0.h:106
DQMSourcePi0::fileName_
std::string fileName_
Output file name if required.
Definition: DQMSourcePi0.h:343
DQMSourcePi0::clusSeedThr_
double clusSeedThr_
Definition: DQMSourcePi0.h:256
DQMSourcePi0::hRechitEnergyEEpi0_
MonitorElement * hRechitEnergyEEpi0_
Energy Distribution of rechits EE (pi0)
Definition: DQMSourcePi0.h:109
DQMSourcePi0::seleMinvMaxEtaEndCap_
double seleMinvMaxEtaEndCap_
Definition: DQMSourcePi0.h:301
DQMSourcePi0::EERecHits
std::vector< EcalRecHit > EERecHits
Definition: DQMSourcePi0.h:323
DQMSourcePi0::hPtEtaEB_
MonitorElement * hPtEtaEB_
Eta Pt in EB.
Definition: DQMSourcePi0.h:198
DQMSourcePi0::selePi0Iso_
double selePi0Iso_
Definition: DQMSourcePi0.h:270
DQMSourcePi0::hIsoEtaEB_
MonitorElement * hIsoEtaEB_
Eta Iso EB.
Definition: DQMSourcePi0.h:210
DQMSourcePi0::selePtEtaEndCap_
double selePtEtaEndCap_
Definition: DQMSourcePi0.h:300
DQMSourcePi0::saveToFile_
bool saveToFile_
Write to file.
Definition: DQMSourcePi0.h:334
DQMSourcePi0::seleS4S9GammaEtaEndCap_
double seleS4S9GammaEtaEndCap_
Definition: DQMSourcePi0.h:298
spr::find
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
DQMSourcePi0::hIsoPi0EB_
MonitorElement * hIsoPi0EB_
Pi0 Iso EB.
Definition: DQMSourcePi0.h:204
DQMSourcePi0::hS4S91Pi0EE_
MonitorElement * hS4S91Pi0EE_
S4S9 of the 1st most energetic pi0 photon EE.
Definition: DQMSourcePi0.h:219
DQMSourcePi0::hS4S91EtaEB_
MonitorElement * hS4S91EtaEB_
S4S9 of the 1st most energetic eta photon.
Definition: DQMSourcePi0.h:222
DQMSourcePi0::selePi0BeltDeta_
double selePi0BeltDeta_
Definition: DQMSourcePi0.h:269
DQMSourcePi0::seleS4S9Gamma_
double seleS4S9Gamma_
Definition: DQMSourcePi0.h:267
DQMSourcePi0::ptMinForIsolationEta_
double ptMinForIsolationEta_
Definition: DQMSourcePi0.h:291
MonitorElement::Fill
void Fill(long long x)
Definition: MonitorElement.h:170
EBDetId::iphi
int iphi() const
get the crystal iphi
Definition: EBDetId.h:53
DQMSourcePi0::seleMinvMinPi0_
double seleMinvMinPi0_
Definition: DQMSourcePi0.h:266
DQMSourcePi0::hMinvPi0EB_
MonitorElement * hMinvPi0EB_
Pi0 invariant mass in EB.
Definition: DQMSourcePi0.h:154
DQMSourcePi0::seleMinvMaxEta_
double seleMinvMaxEta_
Definition: DQMSourcePi0.h:289
DQMSourcePi0::hS4S92Pi0EE_
MonitorElement * hS4S92Pi0EE_
S4S9 of the 2nd most energetic pi0 photon EE.
Definition: DQMSourcePi0.h:231
DQMSourcePi0::convxtalid
void convxtalid(int &, int &)
Definition: DQMSourcePi0.cc:1532
DQMSourcePi0::ptMinForIsolationEndCap_
double ptMinForIsolationEndCap_
Definition: DQMSourcePi0.h:282
iEvent
int iEvent
Definition: GenABIO.cc:230
DQMSourcePi0::seleS4S9GammaEta_
double seleS4S9GammaEta_
Definition: DQMSourcePi0.h:287
DQMSourcePi0::hiEtaDistrEBpi0_
MonitorElement * hiEtaDistrEBpi0_
Distribution of rechits in iEta (pi0)
Definition: DQMSourcePi0.h:94
DQMSourcePi0::selePi0BeltDetaEndCap_
double selePi0BeltDetaEndCap_
Definition: DQMSourcePi0.h:281
DQMSourcePi0::hS4S91Pi0EB_
MonitorElement * hS4S91Pi0EB_
S4S9 of the 1st most energetic pi0 photon.
Definition: DQMSourcePi0.h:216
DQMSourcePi0::hiXDistrEEpi0_
MonitorElement * hiXDistrEEpi0_
Distribution of rechits in ix EE (pi0)
Definition: DQMSourcePi0.h:85
DQMSourcePi0::hiEtaDistrEBeta_
MonitorElement * hiEtaDistrEBeta_
Distribution of rechits in iEta (eta)
Definition: DQMSourcePi0.h:100
DQMSourcePi0::productMonitoredEBpi0_
edm::EDGetTokenT< EcalRecHitCollection > productMonitoredEBpi0_
object to monitor
Definition: DQMSourcePi0.h:243
DQMSourcePi0::isMonEBpi0_
bool isMonEBpi0_
which subdet will be monitored
Definition: DQMSourcePi0.h:337
DQMSourcePi0::selePtPi0EndCap_
double selePtPi0EndCap_
Definition: DQMSourcePi0.h:275
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:18
DQMSourcePi0::seleMinvMinEta_
double seleMinvMinEta_
Definition: DQMSourcePi0.h:290
DQMSourcePi0::selePtGammaEtaEndCap_
double selePtGammaEtaEndCap_
for eta-&gt;gg endcap
Definition: DQMSourcePi0.h:297
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
DQMSourcePi0::diff_nphi_s
int diff_nphi_s(int, int)
Definition: DQMSourcePi0.cc:1550
DQMSourcePi0::ptMinForIsolation_
double ptMinForIsolation_
Definition: DQMSourcePi0.h:271
DQMSourcePi0::diff_neta_s
int diff_neta_s(int, int)
Definition: DQMSourcePi0.cc:1543
DQMStore::IBooker::book1D
MonitorElement * book1D(Args &&...args)
Definition: DQMStore.h:115
DQMSourcePi0::seleXtalMinEnergy_
double seleXtalMinEnergy_
Definition: DQMSourcePi0.h:253
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
j
int j
Definition: DBlmapReader.cc:9
DQMSourcePi0::DQMSourcePi0
DQMSourcePi0(const edm::ParameterSet &)
Definition: DQMSourcePi0.cc:52
EEDetId::iy
int iy() const
Definition: EEDetId.h:82
EBDetId::ieta
int ieta() const
get the crystal ieta
Definition: EBDetId.h:51
HLT_FULL_cff.posCalcParameters
tuple posCalcParameters
Definition: HLT_FULL_cff.py:60262
DQMSourcePi0::seleEtaBeltDeta_
double seleEtaBeltDeta_
Definition: DQMSourcePi0.h:294
edm::HandleBase::isValid
bool isValid() const
Definition: HandleBase.h:75
DQMSourcePi0::hPtEtaEE_
MonitorElement * hPtEtaEE_
Eta Pt in EE.
Definition: DQMSourcePi0.h:201
DQMSourcePi0::seleXtalMinEnergyEndCap_
double seleXtalMinEnergyEndCap_
Definition: DQMSourcePi0.h:254
DQMSourcePi0::seleS9S25GammaEtaEndCap_
double seleS9S25GammaEtaEndCap_
Definition: DQMSourcePi0.h:299
DQMSourcePi0::selePi0BeltDR_
double selePi0BeltDR_
Definition: DQMSourcePi0.h:268
DQMSourcePi0::hEventEnergyEEpi0_
MonitorElement * hEventEnergyEEpi0_
Distribution of total event energy EE (pi0)
Definition: DQMSourcePi0.h:121
DQMSourcePi0::endLuminosityBlock
void endLuminosityBlock(const edm::LuminosityBlock &lumiSeg, const edm::EventSetup &c)
Definition: DQMSourcePi0.cc:1514
DQMSourcePi0::hNRecHitsEEeta_
MonitorElement * hNRecHitsEEeta_
Distribution of number of RecHits EE (eta)
Definition: DQMSourcePi0.h:139
DQMSourcePi0::selePtGammaEta_
double selePtGammaEta_
for eta-&gt;gg barrel
Definition: DQMSourcePi0.h:285
DQMSourcePi0::productMonitoredEEpi0_
edm::EDGetTokenT< EcalRecHitCollection > productMonitoredEEpi0_
object to monitor
Definition: DQMSourcePi0.h:247
DQMSourcePi0::analyze
void analyze(const edm::Event &e, const edm::EventSetup &c)
Definition: DQMSourcePi0.cc:343
DQMStore::IBooker::setCurrentFolder
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:276
CaloSubdetectorTopology::getWindow
virtual std::vector< DetId > getWindow(const DetId &id, const int &northSouthSize, const int &eastWestSize) const
Definition: CaloSubdetectorTopology.cc:4
DQMSourcePi0::seleEtaBeltDetaEndCap_
double seleEtaBeltDetaEndCap_
Definition: DQMSourcePi0.h:306
edm::Handle::product
T const * product() const
Definition: Handle.h:81
math::XYZPoint
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
DQMSourcePi0::hiPhiDistrEBeta_
MonitorElement * hiPhiDistrEBeta_
Distribution of rechits in iPhi (eta)
Definition: DQMSourcePi0.h:88
edm::EventSetup::get
const T & get() const
Definition: EventSetup.h:56
DQMSourcePi0::selePtEta_
double selePtEta_
Definition: DQMSourcePi0.h:286
DQMSourcePi0::hRechitEnergyEBeta_
MonitorElement * hRechitEnergyEBeta_
Energy Distribution of rechits EB (eta)
Definition: DQMSourcePi0.h:112
DQMSourcePi0::hPt1EtaEE_
MonitorElement * hPt1EtaEE_
Pt of the 1st most energetic Eta photon in EE.
Definition: DQMSourcePi0.h:175
DQMSourcePi0::hEventEnergyEBpi0_
MonitorElement * hEventEnergyEBpi0_
Distribution of total event energy EB (pi0)
Definition: DQMSourcePi0.h:118
DQMSourcePi0::seleEtaIso_
double seleEtaIso_
Definition: DQMSourcePi0.h:292
DQMSourcePi0::beginLuminosityBlock
void beginLuminosityBlock(const edm::LuminosityBlock &lumiSeg, const edm::EventSetup &context)
Definition: DQMSourcePi0.cc:336
DQMSourcePi0::hPt1EtaEB_
MonitorElement * hPt1EtaEB_
Pt of the 1st most energetic Eta photon in EB.
Definition: DQMSourcePi0.h:172
DQMSourcePi0::hEventEnergyEBeta_
MonitorElement * hEventEnergyEBeta_
Distribution of total event energy EB (eta)
Definition: DQMSourcePi0.h:124
DQMSourcePi0::hIsoEtaEE_
MonitorElement * hIsoEtaEE_
Eta Iso EE.
Definition: DQMSourcePi0.h:213
DQMSourcePi0::hS4S92Pi0EB_
MonitorElement * hS4S92Pi0EB_
S4S9 of the 2nd most energetic pi0 photon.
Definition: DQMSourcePi0.h:228
DQMSourcePi0::selePi0IsoEndCap_
double selePi0IsoEndCap_
Definition: DQMSourcePi0.h:279
DQMSourcePi0::hPt2EtaEB_
MonitorElement * hPt2EtaEB_
Pt of the 2nd most energetic Eta photon in EB.
Definition: DQMSourcePi0.h:185
DQMSourcePi0::seleMinvMaxPi0_
double seleMinvMaxPi0_
Definition: DQMSourcePi0.h:265
DQMSourcePi0::EBRecHits
std::vector< EcalRecHit > EBRecHits
Definition: DQMSourcePi0.h:319
DQMSourcePi0::hMinvEtaEE_
MonitorElement * hMinvEtaEE_
Eta invariant mass in EE.
Definition: DQMSourcePi0.h:163
DQMSourcePi0::hPtPi0EB_
MonitorElement * hPtPi0EB_
Pi0 Pt in EB.
Definition: DQMSourcePi0.h:192
PositionCalc::Calculate_Location
math::XYZPoint Calculate_Location(const HitsAndFractions &iDetIds, const edm::SortedCollection< HitType > *iRecHits, const CaloSubdetectorGeometry *iSubGeom, const CaloSubdetectorGeometry *iESGeom=0)
Definition: PositionCalc.h:68
DQMSourcePi0::seleEtaIsoEndCap_
double seleEtaIsoEndCap_
Definition: DQMSourcePi0.h:304
DQMSourcePi0::hS4S92EtaEE_
MonitorElement * hS4S92EtaEE_
S4S9 of the 2nd most energetic eta photon EE.
Definition: DQMSourcePi0.h:237
DQMSourcePi0::prescaleFactor_
unsigned int prescaleFactor_
Monitor every prescaleFactor_ events.
Definition: DQMSourcePi0.h:328
DQMSourcePi0::hMeanRecHitEnergyEEeta_
MonitorElement * hMeanRecHitEnergyEEeta_
Distribution of Mean energy per rechit EE (eta)
Definition: DQMSourcePi0.h:151
DQMSourcePi0::hIsoPi0EE_
MonitorElement * hIsoPi0EE_
Pi0 Iso EE.
Definition: DQMSourcePi0.h:207
DQMSourcePi0::selePtPi0_
double selePtPi0_
Definition: DQMSourcePi0.h:264
DQMSourcePi0::seleS4S9GammaEndCap_
double seleS4S9GammaEndCap_
Definition: DQMSourcePi0.h:278
DQMSourcePi0::hNRecHitsEBeta_
MonitorElement * hNRecHitsEBeta_
Distribution of number of RecHits EB (eta)
Definition: DQMSourcePi0.h:136
relval_parameters_module.energy
string energy
Definition: relval_parameters_module.py:29
DQMSourcePi0::hPt1Pi0EE_
MonitorElement * hPt1Pi0EE_
Pt of the 1st most energetic Pi0 photon in EE.
Definition: DQMSourcePi0.h:169
DQMSourcePi0::hS4S91EtaEE_
MonitorElement * hS4S91EtaEE_
S4S9 of the 1st most energetic eta photon EE.
Definition: DQMSourcePi0.h:225
DQMSourcePi0::hRechitEnergyEEeta_
MonitorElement * hRechitEnergyEEeta_
Energy Distribution of rechits EE (eta)
Definition: DQMSourcePi0.h:115
DQMSourcePi0::hMeanRecHitEnergyEBpi0_
MonitorElement * hMeanRecHitEnergyEBpi0_
Distribution of Mean energy per rechit EB (pi0)
Definition: DQMSourcePi0.h:142
DQMSourcePi0::hEventEnergyEEeta_
MonitorElement * hEventEnergyEEeta_
Distribution of total event energy EE (eta)
Definition: DQMSourcePi0.h:127
MonitorElement::setAxisTitle
void setAxisTitle(const std::string &title, int axis=1)
set x-, y- or z-axis title (axis=1, 2, 3 respectively)
Definition: MonitorElement.cc:947
DQMSourcePi0::hNRecHitsEEpi0_
MonitorElement * hNRecHitsEEpi0_
Distribution of number of RecHits EE (pi0)
Definition: DQMSourcePi0.h:133
DQMSourcePi0::selePi0BeltDREndCap_
double selePi0BeltDREndCap_
Definition: DQMSourcePi0.h:280
edm::Run
Definition: Run.h:43
DQMSourcePi0::seleS9S25GammaEta_
double seleS9S25GammaEta_
Definition: DQMSourcePi0.h:288
DQMSourcePi0::hPt2Pi0EB_
MonitorElement * hPt2Pi0EB_
Pt of the 2nd most energetic Pi0 photon in EB.
Definition: DQMSourcePi0.h:179
DQMSourcePi0::selePtGamma_
double selePtGamma_
Definition: DQMSourcePi0.h:263
DQMSourcePi0::clusSeedThrEndCap_
double clusSeedThrEndCap_
Definition: DQMSourcePi0.h:260
Generated for CMSSW Reference Manual by doxygen 1.8.5