CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ZeeCalibration.cc
Go to the documentation of this file.
1 #include "FWCore/ParameterSet/interface/ParameterSet.h"
2 
3 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
4 #include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
5 #include "DataFormats/EcalDetId/interface/EBDetId.h"
6 #include "DataFormats/EcalDetId/interface/EEDetId.h"
7 #include "DataFormats/DetId/interface/DetId.h"
8 #include "CondFormats/EcalObjects/interface/EcalIntercalibConstants.h"
9 #include "CondFormats/DataRecord/interface/EcalIntercalibConstantsRcd.h"
10 #include "Calibration/Tools/interface/calibXMLwriter.h"
11 
12 #include "Calibration/Tools/interface/CalibrationCluster.h"
13 #include "Calibration/Tools/interface/HouseholderDecomposition.h"
14 #include "Calibration/Tools/interface/MinL3Algorithm.h"
15 #include "Calibration/Tools/interface/EcalRingCalibrationTools.h"
16 #include "Calibration/Tools/interface/EcalIndexingTools.h"
17 
18 #include "CLHEP/Vector/LorentzVector.h"
19 
20 #include "FWCore/Framework/interface/EventSetup.h"
21 #include "FWCore/Framework/interface/ESHandle.h"
22 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
23 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
24 #include "Geometry/Records/interface/CaloGeometryRecord.h"
25 
26 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
27 
28 #include "Calibration/EcalCalibAlgos/interface/ZeeCalibration.h"
29 #include "Calibration/EcalCalibAlgos/interface/ZeeKinematicTools.h"
30 
31 #include "CalibCalorimetry/CaloMiscalibTools/interface/MiscalibReaderFromXMLEcalBarrel.h"
32 #include "CalibCalorimetry/CaloMiscalibTools/interface/MiscalibReaderFromXMLEcalEndcap.h"
33 #include "CalibCalorimetry/CaloMiscalibTools/interface/CaloMiscalibMapEcal.h"
34 
35 #include "DataFormats/EgammaCandidates/interface/Electron.h"
36 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
37 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
38 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
39 #include "DataFormats/TrackReco/interface/Track.h"
40 #include "DataFormats/TrackReco/interface/TrackFwd.h"
41 #include "DataFormats/TrackReco/interface/TrackExtraFwd.h"
42 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
43 
44 #include "DataFormats/CaloRecHit/interface/CaloRecHit.h"
45 
47 #include "HLTrigger/HLTanalyzers/interface/HLTrigReport.h"
48 #include "DataFormats/Common/interface/Handle.h"
49 #include "DataFormats/Common/interface/TriggerResults.h"
50 #include "FWCore/MessageLogger/interface/MessageLogger.h"
51 #include "FWCore/ServiceRegistry/interface/Service.h"
52 #include "FWCore/Framework/interface/TriggerNamesService.h"
53 
54 #include "TTree.h"
55 #include "TBranch.h"
56 #include "TCanvas.h"
57 #include "TFile.h"
58 #include "TProfile.h"
59 #include "TH1.h"
60 #include "TH2.h"
61 #include "TF1.h"
62 #include "TGraph.h"
63 #include "TGraphErrors.h"
64 #include "TRandom.h"
65 
66 #include <iostream>
67 #include <string>
68 #include <stdexcept>
69 #include <vector>
70 #include <utility>
71 #include <map>
72 #include <fstream>
73 
74 #define MZ 91.1876
75 
76 #define DEBUG 1
77 
78 ZeeCalibration::ZeeCalibration(const edm::ParameterSet& iConfig)
79 {
80 
81 #ifdef DEBUG
82  std::cout<<"[ZeeCalibration] Starting the ctor"<<std::endl;
83 #endif
84 
85  theMaxLoops = iConfig.getUntrackedParameter<unsigned int>("maxLoops",0);
86 
87  wantEtaCorrection_ = iConfig.getUntrackedParameter<bool>("wantEtaCorrection",true);
88 
89  outputFileName_ = iConfig.getParameter<std::string>("outputFile");
90 
91  minInvMassCut_ = iConfig.getUntrackedParameter<double>("minInvMassCut", 70.);
92  maxInvMassCut_ = iConfig.getUntrackedParameter<double>("maxInvMassCut", 110.);
93 
94  rechitProducer_ = iConfig.getParameter<std::string>("rechitProducer");
95  rechitCollection_ = iConfig.getParameter<std::string>("rechitCollection");
96 
97  erechitProducer_ = iConfig.getParameter<std::string>("erechitProducer");
98  erechitCollection_ = iConfig.getParameter<std::string>("erechitCollection");
99 
100  scProducer_ = iConfig.getParameter<std::string>("scProducer");
101  scCollection_ = iConfig.getParameter<std::string>("scCollection");
102 
103  scIslandProducer_ = iConfig.getParameter<std::string>("scIslandProducer");
104  scIslandCollection_ = iConfig.getParameter<std::string>("scIslandCollection");
105 
106  calibMode_ = iConfig.getUntrackedParameter<std::string>("ZCalib_CalibType");
107 
108  mcProducer_ = iConfig.getUntrackedParameter<std::string>("mcProducer","");
109 
110 
111  electronProducer_ = iConfig.getParameter<std::string > ("electronProducer");
112  electronCollection_ = iConfig.getParameter<std::string > ("electronCollection");
113 
114  outputFile_ = TFile::Open(outputFileName_.c_str(),"RECREATE"); // open output file to store histograms
115 
116  myTree = new TTree("myTree","myTree");
117  // myTree->Branch("zMass","zMass", &mass);
118  myTree->Branch("zMass",&mass4tree,"mass/F");
119  myTree->Branch("zMassDiff",&massDiff4tree,"massDiff/F");
120 
121  barrelfile_=iConfig.getUntrackedParameter<std::string> ("initialMiscalibrationBarrel","");
122  endcapfile_=iConfig.getUntrackedParameter<std::string> ("initialMiscalibrationEndcap","");
123 
124  electronSelection_=iConfig.getUntrackedParameter<unsigned int> ("electronSelection",0);//option for electron selection
125 
126  etaBins_ = iConfig.getUntrackedParameter<unsigned int>("etaBins", 10);
127  etBins_ = iConfig.getUntrackedParameter<unsigned int>("etBins", 10);
128 
129  etaMin_ = iConfig.getUntrackedParameter<double>("etaMin", 0.);
130  etMin_ = iConfig.getUntrackedParameter<double>("etMin", 0.);
131  etaMax_ = iConfig.getUntrackedParameter<double>("etaMax", 3.);
132  etMax_ = iConfig.getUntrackedParameter<double>("etMax", 100.);
133 
134 
135  // new ZeePlots("zeePlots.root");
136  // ZeePlots->bookHistos();
137 
138  //ZeeCalibrationPLots("zeeCalibPlots");
139  //ZeecaPlots->bookHistos(maxsIter);
140 
141  hlTriggerResults_ = iConfig.getParameter<edm::InputTag> ("HLTriggerResults");
142 
143  theParameterSet=iConfig;
144  EcalIndexingTools* myIndexTool=0;
145 
146 
147  myIndexTool = EcalIndexingTools::getInstance();
148 
149  myIndexTool->setBinRange( etaBins_, etaMin_, etaMax_, etBins_, etMin_, etMax_ );
150 
151  //creating the algorithm
152  theAlgorithm_ = new ZIterativeAlgorithmWithFit(iConfig);
153 
154  // Tell the framework what data is being produced
155  //setWhatProduced(this);
156  setWhatProduced (this, &ZeeCalibration::produceEcalIntercalibConstants ) ;
157  findingRecord<EcalIntercalibConstantsRcd> () ;
158 
159  for(int i = 0; i<50; i++){
160 
161  coefficientDistanceAtIteration[i] = -1.;
162  loopArray[i] = -1.;
163  sigmaArray[i] = -1.;
164  sigmaErrorArray[i] = -1.;
165 
166  }
167 
168 #ifdef DEBUG
169  std::cout<<"[ZeeCalibration] Done with the ctor"<<std::endl;
170 #endif
171 
172 }
173 
174 
175 ZeeCalibration::~ZeeCalibration()
176 {
177 // if (theAlgorithm_)
178 // delete theAlgorithm_;
179 }
180 
181 //_____________________________________________________________________________
182 // Produce EcalIntercalibConstants
183 boost::shared_ptr<EcalIntercalibConstants>
184 ZeeCalibration::produceEcalIntercalibConstants( const EcalIntercalibConstantsRcd& iRecord )
185 {
186  std::cout << "@SUB=ZeeCalibration::produceEcalIntercalibConstants" << std::endl;
187  return ical;
188 }
189 
190 void ZeeCalibration::beginOfJob(){isfirstcall_=true;}
191 
192 
193 
194 
195 //========================================================================
196 void
197 ZeeCalibration::endOfJob() {
198 
199 
200  printStatistics();
201 
202  if(calibMode_ != "ETA_ET_MODE"){
203 
205 
206  //Writing out calibration coefficients
207  calibXMLwriter* barrelWriter = new calibXMLwriter(EcalBarrel);
208  for(int ieta=-EBDetId::MAX_IETA; ieta<=EBDetId::MAX_IETA ;++ieta) {
209  if(ieta==0) continue;
210  for(int iphi=EBDetId::MIN_IPHI; iphi<=EBDetId::MAX_IPHI; ++iphi) {
211  // make an EBDetId since we need EBDetId::rawId() to be used as the key for the pedestals
212  if (EBDetId::validDetId(ieta,iphi))
213  {
214  EBDetId ebid(ieta,iphi);
215  barrelWriter->writeLine(ebid,* (ical->getMap().find(ebid.rawId()) ));
216  }
217  }
218  }
219 
220 
221 
222  calibXMLwriter* endcapWriter = new calibXMLwriter(EcalEndcap);
223  for(int iX=EEDetId::IX_MIN; iX<=EEDetId::IX_MAX ;++iX) {
224  for(int iY=EEDetId::IY_MIN; iY<=EEDetId::IY_MAX; ++iY) {
225  // make an EEDetId since we need EEDetId::rawId() to be used as the key for the pedestals
226  if (EEDetId::validDetId(iX,iY,1))
227  {
228  EEDetId eeid(iX,iY,1);
229  endcapWriter->writeLine(eeid,*(ical->getMap().find(eeid.rawId()) ) );
230  }
231  if (EEDetId::validDetId(iX,iY,-1))
232  {
233  EEDetId eeid(iX,iY,-1);
234  endcapWriter->writeLine(eeid, *(ical->getMap().find(eeid.rawId())) );
235  }
236 
237  }
238  }
239 
240 
241  }
242 
243  std::cout<<"Writing histos..."<<std::endl;
244  outputFile_->cd();
245 
246  // zeeplts->Write();
247 
248  h1_eventsBeforeEWKSelection_ ->Write();
249  h1_eventsAfterEWKSelection_ ->Write();
250 
251  h1_eventsBeforeBorderSelection_ ->Write();
252  h1_eventsAfterBorderSelection_ ->Write();
253 
254  h1_borderElectronClassification_->Write();
255 
256  h2_xtalMiscalibCoeffBarrel_->Write();
257  h2_xtalMiscalibCoeffEndcapMinus_->Write();
258  h2_xtalMiscalibCoeffEndcapPlus_->Write();
259 
260  h1_electronCosTheta_SC_->Write();
261  h1_electronCosTheta_TK_->Write();
262  h1_electronCosTheta_SC_TK_->Write();
263 
264  h1_zMassResol_->Write();
265  h1_zEtaResol_->Write();
266  h1_zPhiResol_->Write();
267  h1_eleEtaResol_->Write();
268  h1_elePhiResol_->Write();
269  h1_seedOverSC_ ->Write();
270  h1_preshowerOverSC_ ->Write();
271 
272  for(unsigned int i =0; i<25; i++){
273  if( i < theMaxLoops ){
274 
275  h_ESCEtrueVsEta_[i]->Write();
276  h_ESCEtrue_[i]->Write();
277 
278  h_ESCcorrEtrueVsEta_[i]->Write();
279  h_ESCcorrEtrue_[i]->Write();
280 
281  h2_chi2_[i]->Write();
282  h2_iterations_[i]->Write();
283 
284  // h_DiffZMassDistr_[i]->Write();
285 
286  //h_ZMassDistr_[i]->Write();
287  }
288  }
289 
290  h2_fEtaBarrelGood_->Write();
291  h2_fEtaBarrelBad_->Write();
292  h2_fEtaEndcapGood_->Write();
293  h2_fEtaEndcapBad_->Write();
294  h1_eleClasses_->Write();
295 
296  h_eleEffEta_[0]->Write();
297  h_eleEffPhi_[0]->Write();
298  h_eleEffPt_[0]->Write();
299 
300  h_eleEffEta_[1]->Write();
301  h_eleEffPhi_[1]->Write();
302  h_eleEffPt_[1]->Write();
303 
304 
305  int j = 0;
306 
307  int flag=0;
308 
309  Double_t mean[25] = {0.};
310  Double_t num[25] = {0.};
311  Double_t meanErr[25] = {0.};
312  Float_t rms[25] = {0.};
313  Float_t tempRms[10][25];
314 
315  for(int ia = 0; ia<10; ia++){
316  for(int ib = 0; ib<25; ib++){
317 
318  tempRms[ia][ib] = 0.;
319 
320  }
321  }
322 
323  int aa = 0;
324 
325  for( int k = 0; k<theAlgorithm_->getNumberOfChannels(); k++ )
326  {
327 
328 
329 
331  bool isNearCrack = false;
332 
333  if( calibMode_ == "RING"){
334 
335  isNearCrack = ( abs( ringNumberCorrector(k) ) == 1 || abs( ringNumberCorrector(k) ) == 25 ||
336  abs( ringNumberCorrector(k) ) == 26 || abs( ringNumberCorrector(k) ) == 45 ||
337  abs( ringNumberCorrector(k) ) == 46 || abs( ringNumberCorrector(k) ) == 65 ||
338  abs( ringNumberCorrector(k) ) == 66 || abs( ringNumberCorrector(k) ) == 85 ||
339  abs( ringNumberCorrector(k) ) == 86 || abs( ringNumberCorrector(k) ) == 124 );
340  }
341 
342  if(k<85)
343  {
344 
345  if((k+1)%5!=0)
346  {
347 
348  if(!isNearCrack){
349  mean[j]+=calibCoeff[k];
350  mean[j]+=calibCoeff[169 - k];
351 
352  num[j] += 2.;
353 
354  //meanErr[j]+= calibCoeffError[k];
355  //meanErr[j]+= calibCoeffError[169 - k];
356 
357  meanErr[j]+= 1./ pow ( calibCoeffError[k], 2 );
358  meanErr[j]+= 1./ pow ( calibCoeffError[169 - k], 2);
359 
360 
361  tempRms[aa][j]+=calibCoeff[k];
362  aa++;
363  tempRms[aa][j]+=calibCoeff[169 - k];
364  aa++;
365  }
366  }
367 
368  else {
369  if(!isNearCrack){
370  mean[j]+=calibCoeff[k];
371  mean[j]+=calibCoeff[169 - k];
372 
373  num[j] += 2.;
374 
375  //meanErr[j]+= calibCoeffError[k];
376  //meanErr[j]+= calibCoeffError[169 - k];
377 
378  meanErr[j]+= 1./ pow ( calibCoeffError[k], 2 );
379  meanErr[j]+= 1./ pow ( calibCoeffError[169 - k], 2);
380 
381  tempRms[aa][j]+=calibCoeff[k];
382  aa++;
383  tempRms[aa][j]+=calibCoeff[169 - k];
384  aa++;
385 
386  }
387  j++;
388  aa = 0;
389 
390  }
391 
392  }
393  //EE begin
394 
395 
396  if(k>=170 && k<=204){
397 
398  if(flag<4){
399  //make groups of 5 Xtals in #eta
400  mean[j]+=calibCoeff[k]/10.;
401  mean[j]+=calibCoeff[k+39]/10.;
402 
403  meanErr[j]+= calibCoeffError[k]/30.;
404  meanErr[j]+= calibCoeffError[k + 39]/30.;
405 
406 
407  tempRms[aa][j]+=calibCoeff[k];
408  aa++;
409  tempRms[aa][j]+=calibCoeff[k + 39];
410  aa++;
411 
412  flag++;
413  }
414 
415  else if(flag==4){
416  //make groups of 5 Xtals in #eta
417  mean[j]+=calibCoeff[k]/10.;
418  mean[j]+=calibCoeff[k+39]/10.;
419 
420  meanErr[j]+= calibCoeffError[k]/30.;
421  meanErr[j]+= calibCoeffError[k + 39]/30.;
422 
423  tempRms[aa][j]+=calibCoeff[k];
424  aa++;
425  tempRms[aa][j]+=calibCoeff[k + 39];
426  aa++;
427 
428  flag=0;
429  // std::cout<<" index(>85) "<<k<<" j is "<<j<<" mean[j] is "<<mean[j]<<std::endl;
430  j++;
431  aa = 0;
432 
433  }
434 
435  }
436  if(k>=205 && k<=208){
437  mean[j]+=calibCoeff[k]/8.;
438  mean[j]+=calibCoeff[k+39]/8.;
439 
440  meanErr[j]+= calibCoeffError[k]/30.;
441  meanErr[j]+= calibCoeffError[k + 39]/30.;
442 
443 
444  tempRms[aa][j]+=calibCoeff[k];
445  aa++;
446  tempRms[aa][j]+=calibCoeff[k + 39];
447  aa++;
448  }
449  //EE end
450 
451  /*
452  for(int jj =0; jj< 25; jj++){
453  if(meanErr[jj] > 0.)
454  std::cout<<" meanErr[jj] before sqrt: "<<meanErr[jj]<<std::endl;
455 
456  meanErr[jj] = 1./sqrt( meanErr[jj] );
457 
458  std::cout<<" meanErr[jj] after sqrt: "<<meanErr[jj]<<std::endl;
459  }
460  */
461 
462 
463 
464  if(!isNearCrack){
465  h2_coeffVsEta_->Fill( ringNumberCorrector(k), calibCoeff[k] );
466  h2_miscalRecal_->Fill( initCalibCoeff[k], 1./calibCoeff[k] );
467  h1_mc_->Fill( initCalibCoeff[k]*calibCoeff[k] -1. );
468 
469 
470 
471  if(k<170){
472  h2_miscalRecalEB_->Fill( initCalibCoeff[k], 1./calibCoeff[k] );
473  h1_mcEB_->Fill( initCalibCoeff[k]*calibCoeff[k] -1. );
474  }
475 
476  if(k>=170){
477  h2_miscalRecalEE_->Fill( initCalibCoeff[k], 1./calibCoeff[k] );
478  h1_mcEE_->Fill( initCalibCoeff[k]*calibCoeff[k] -1. );
479  }
480 
481  }
482  }
483 
484  for(int ic = 0; ic< 17; ic++){
485 
486  mean[ic] = mean[ic] / num[ic]; //find mean of recalib coeff on group of rings
487  //meanErr[ic] = meanErr[ic] / ( sqrt( num[ic] ) * num[ic] ); //find mean of recalib coeff on group of rings
488  meanErr[ic] = 1. / sqrt(meanErr[ic]); //find mean of recalib coeff on group of rings
489 
490  }
491 
492 
493  //build array of RMS
494  for(int ic = 0; ic< 25; ic++){
495  for(int id = 0; id< 10; id++){
496 
497  if(tempRms[id][ic] > 0.){
498 
499  rms[ic] += (tempRms[id][ic] - mean[j])*(tempRms[id][ic] - mean[j]);
500 
501  }
502  }
503  rms[ic]/= 10.;//this is approximate
504  rms[ic] = sqrt(rms[ic]);
505  }
506 
507  //build array of RMS
508 
509 
510 
511  Double_t xtalEta[25] = {1.4425, 1.3567,1.2711,1.1855,
512  1.10,1.01,0.92,0.83,
513  0.7468,0.6612,0.5756,0.4897,0.3985,0.3117,0.2250,0.1384,0.0487,
514  1.546, 1.651, 1.771, 1.908, 2.071, 2.267, 2.516, 2.8};
515 
516  Double_t zero[25] = {0.026};//interval/sqrt(12)
517 
518  for(int j = 0; j <25; j++)
519  h2_coeffVsEtaGrouped_->Fill( xtalEta[j],mean[j]);
520 
521  // for(int sho = 0; sho <25; sho++)
522  //cout<<"xtalEta[j] "<< xtalEta[sho]<<" mean[j] "<<mean[sho]<<" err[j] "<<meanErr[sho]<<std::endl;
523 
524  TProfile *px = h2_coeffVsEta_->ProfileX("coeffVsEtaProfile");
525  px->SetXTitle("Eta channel");
526  px->SetYTitle("recalibCoeff");
527  px->Write();
528 
529  h2_coeffVsEta_->Write();
530  h2_coeffVsEtaGrouped_->Write();
531  h2_zMassVsLoop_->Write();
532  h2_zMassDiffVsLoop_->Write();
533  h2_zWidthVsLoop_->Write();
534  h2_coeffVsLoop_->Write();
535  h2_miscalRecal_->Write();
536  h1_mc_->Write();
537  h2_miscalRecalEB_->Write();
538  h1_mcEB_->Write();
539  h2_miscalRecalEE_->Write();
540  h1_mcEE_->Write();
541 
542  h2_residualSigma_->Write();
543 
544  const ZIterativeAlgorithmWithFit::ZIterativeAlgorithmWithFitPlots* algoHistos=theAlgorithm_->getHistos();
545 
546  double weightSumMeanBarrel = 0.;
547  double weightSumMeanEndcap = 0.;
548 
549  for (int iIteration=0;iIteration<theAlgorithm_->getNumberOfIterations();iIteration++)
550  for (int iChannel=0;iChannel<theAlgorithm_->getNumberOfChannels();iChannel++)
551  {
552 
553  if( iIteration==(theAlgorithm_->getNumberOfIterations()-1) ){
554 
555  if(iChannel < 170)
556  weightSumMeanBarrel += algoHistos->weightedRescaleFactor[iIteration][iChannel]->Integral()/170.;
557 
558  if(iChannel >= 170)
559  weightSumMeanEndcap += algoHistos->weightedRescaleFactor[iIteration][iChannel]->Integral()/78.;
560 
561  h1_occupancyVsEta_->Fill((Double_t)ringNumberCorrector(iChannel), algoHistos->weightedRescaleFactor[iIteration][iChannel]->Integral() );
562 
563 
564  h1_occupancy_->Fill( algoHistos->weightedRescaleFactor[iIteration][iChannel]->Integral() );
565 
566  if(iChannel < 170)
567  h1_occupancyBarrel_->Fill( algoHistos->weightedRescaleFactor[iIteration][iChannel]->Integral() );
568 
569  if(iChannel >= 170)
570  h1_occupancyEndcap_->Fill( algoHistos->weightedRescaleFactor[iIteration][iChannel]->Integral() );
571 
572 #ifdef DEBUG
573  std::cout<<"Writing weighted integral for channel "<<ringNumberCorrector(iChannel)<<" ,value "<<algoHistos->weightedRescaleFactor[iIteration][iChannel]->Integral()<<std::endl;
574 #endif
575 
576  }
577 
578  }
579 
580  // std::cout<<"Done! Closing output file... "<<std::endl;
581 
582  h1_weightSumMeanBarrel_ ->Fill(weightSumMeanBarrel);
583  h1_weightSumMeanEndcap_ ->Fill(weightSumMeanEndcap);
584 
585  std::cout<<"Weight sum mean on channels in Barrel is :"<<weightSumMeanBarrel<<std::endl;
586  std::cout<<"Weight sum mean on channels in Endcap is :"<<weightSumMeanEndcap<<std::endl;
587 
588  h1_weightSumMeanBarrel_ ->Write();
589  h1_weightSumMeanEndcap_ ->Write();
590 
591  h1_occupancyVsEta_->Write();
592  h1_occupancy_->Write();
593  h1_occupancyBarrel_->Write();
594  h1_occupancyEndcap_->Write();
595 
596  myTree->Write();
597 
598  TGraphErrors* graph = new TGraphErrors(25,xtalEta,mean,zero,meanErr);
599  graph->Draw("APL");
600  graph->Write();
601 
602  double zero50[50] = { 0. };
603 
604  TGraphErrors* residualSigmaGraph = new TGraphErrors(50,loopArray,sigmaArray,zero50,sigmaErrorArray);
605  residualSigmaGraph->SetName("residualSigmaGraph");
606  residualSigmaGraph->Draw("APL");
607  residualSigmaGraph->Write();
608 
609  TGraphErrors* coefficientDistanceAtIterationGraph = new TGraphErrors(50,loopArray,coefficientDistanceAtIteration,zero50,zero50);
610  coefficientDistanceAtIterationGraph->SetName("coefficientDistanceAtIterationGraph");
611  coefficientDistanceAtIterationGraph->Draw("APL");
612  coefficientDistanceAtIterationGraph->Write();
613 
614  Float_t noError[250] = {0.};
615 
616  Float_t ringInd[250];
617  for(int i =0; i<250; i++)
618  ringInd[i]=ringNumberCorrector(i);
619 
620  TGraphErrors* graphCoeff = new TGraphErrors(theAlgorithm_->getNumberOfChannels(),ringInd,calibCoeff,noError,calibCoeffError);
621  graphCoeff->SetName("graphCoeff");
622  graphCoeff->Draw("APL");
623  graphCoeff->Write();
624 
625  // outputFile_->Write();//this automatically writes all histos on file
626 
627 
628  h1_ZCandMult_->Write();
629  h1_reco_ZMass_->Write();
630 
631  h1_reco_ZMassCorr_->Write();
632  h1_reco_ZMassCorrBB_->Write();
633  h1_reco_ZMassCorrEE_->Write();
634 
635  outputFile_->Close();
636 
637  myZeePlots_ ->writeEleHistograms();
638  myZeePlots_ ->writeMCEleHistograms();
639  myZeePlots_ ->writeZHistograms();
640  myZeePlots_ ->writeMCZHistograms();
641 
642  // myZeeRescaleFactorPlots_ = new ZeeRescaleFactorPlots("zeeRescaleFactorPlots.root");
643  //myZeeRescaleFactorPlots_->writeHistograms( theAlgorithm_ );
644 
645  // delete myZeeRescaleFactorPlots_;
646 
647 
648 
649 }
650 
651 //_____________________________________________________________________________
652 // Called at each event
653 //________________________________________
654 
655 edm::EDLooper::Status
656 ZeeCalibration::duringLoop( const edm::Event& iEvent, const edm::EventSetup& iSetup )
657 {
658  using namespace edm;
659 
660 #ifdef DEBUG
661  std::cout<<"[ZeeCalibration] Entering duringLoop"<<std::endl;
662 #endif
663 
664 
665  // code that used to be in beginJob
666  if (isfirstcall_){
667 
668  //inizializzare la geometria di ecal
669  edm::ESHandle<CaloGeometry> pG;
670  iSetup.get<CaloGeometryRecord>().get(pG);
671  EcalRingCalibrationTools::setCaloGeometry(&(*pG));
672 
673  myZeePlots_ = new ZeePlots( "zeePlots.root" );
674  // myZeeRescaleFactorPlots_ = new ZeeRescaleFactorPlots("zeeRescaleFactorPlots.root");
675 
676  // go to *OUR* rootfile and book histograms
677  outputFile_->cd();
678  bookHistograms();
679 
680  std::cout<<"[ZeeCalibration::beginOfJob] Histograms booked "<<std::endl;
681 
682  loopFlag_ = 0;
683 
684  //Read miscalibration map if requested
685  CaloMiscalibMapEcal* miscalibMap=0;
686  if (!barrelfile_.empty() || !barrelfile_.empty())
687  {
688  miscalibMap=new CaloMiscalibMapEcal();
689  miscalibMap->prefillMap();
690  }
691 
692 
693  if(!barrelfile_.empty())
694  {
695  MiscalibReaderFromXMLEcalBarrel barrelreader_(*miscalibMap);
696  barrelreader_.parseXMLMiscalibFile(barrelfile_);
697 #ifdef DEBUG
698  std::cout<<"[ZeeCalibration::beginOfJob] Parsed EB miscal file"<<std::endl;
699 #endif
700  }
701 
702  if(!endcapfile_.empty())
703  {
704  MiscalibReaderFromXMLEcalEndcap endcapreader_(*miscalibMap);
705  endcapreader_.parseXMLMiscalibFile(endcapfile_);
706 #ifdef DEBUG
707  std::cout<<"[ZeeCalibration::beginOfJob] Parsed EE miscal file"<<std::endl;
708 #endif
709  }
710 
711  std::cout << " theAlgorithm_->getNumberOfChannels() "
712  << theAlgorithm_->getNumberOfChannels() << std::endl;
713 
714 
716  for(int k = 0; k < theAlgorithm_->getNumberOfChannels(); k++)
717  {
718  calibCoeff[k]=1.;
719  calibCoeffError[k]=0.;
720 
721  std::vector<DetId> ringIds;
722 
723  if(calibMode_ == "RING")
724  ringIds = EcalRingCalibrationTools::getDetIdsInRing(k);
725 
726  if(calibMode_ == "MODULE")
727  ringIds = EcalRingCalibrationTools::getDetIdsInModule(k);
728 
729  if(calibMode_ == "ABS_SCALE" || calibMode_ == "ETA_ET_MODE" )
730  ringIds = EcalRingCalibrationTools::getDetIdsInECAL();
731 
732  if (miscalibMap)
733  {
734  initCalibCoeff[k]=0.;
735  for (unsigned int iid=0; iid<ringIds.size();++iid)
736  {
737  float miscalib=* (miscalibMap->get().getMap().find(ringIds[iid]) );
738  // float miscalib=miscalibMap->get().getMap().find(ringIds[iid])->second; ////////AP
739  initCalibCoeff[k]+=miscalib;
740  }
741  initCalibCoeff[k]/=(float)ringIds.size();
742  std::cout << k << " " << initCalibCoeff[k] << " " << ringIds.size() << std::endl;
743  }
744  else
745  {
746  initCalibCoeff[k]=1.;
747  }
748  }
749 
750  ical = boost::shared_ptr<EcalIntercalibConstants>( new EcalIntercalibConstants() );
751 
752  for(int k = 0; k < theAlgorithm_->getNumberOfChannels(); k++)
753  {
754  // std::vector<DetId> ringIds = EcalRingCalibrationTools::getDetIdsInRing(k);
755 
756  std::vector<DetId> ringIds;
757 
758  if(calibMode_ == "RING")
759  ringIds = EcalRingCalibrationTools::getDetIdsInRing(k);
760 
761  if(calibMode_ == "MODULE")
762  ringIds = EcalRingCalibrationTools::getDetIdsInModule(k);
763 
764  if(calibMode_ == "ABS_SCALE" || calibMode_ == "ETA_ET_MODE")
765  ringIds = EcalRingCalibrationTools::getDetIdsInECAL();
766 
767 
768  for (unsigned int iid=0; iid<ringIds.size();++iid){
769  // ical->setValue( ringIds[iid], 1. * initCalibCoeff[k] );
770 
771  if(ringIds[iid].subdetId() == EcalBarrel){
772  EBDetId myEBDetId(ringIds[iid]);
773  h2_xtalMiscalibCoeffBarrel_->SetBinContent( myEBDetId.ieta() + 86, myEBDetId.iphi(), * (miscalibMap->get().getMap().find(ringIds[iid]) ) );//fill TH2 with miscalibCoeff
774 
775  }
776 
777  if(ringIds[iid].subdetId() == EcalEndcap){
778  EEDetId myEEDetId(ringIds[iid]);
779  if(myEEDetId.zside() < 0)
780  h2_xtalMiscalibCoeffEndcapMinus_->SetBinContent( myEEDetId.ix(), myEEDetId.iy(), * ( miscalibMap->get().getMap().find(ringIds[iid]) ) );//fill TH2 with miscalibCoeff
781 
782  if(myEEDetId.zside() > 0)
783  h2_xtalMiscalibCoeffEndcapPlus_->SetBinContent( myEEDetId.ix(), myEEDetId.iy(), * (miscalibMap->get().getMap().find(ringIds[iid]) ) );//fill TH2 with miscalibCoeff
784 
785  }
786 
787  ical->setValue( ringIds[iid], *(miscalibMap->get().getMap().find(ringIds[iid]) ) );
788 
789  }
790 
791 
792  read_events = 0;
793  init_ = false;
794 
795 
796  }
797  isfirstcall_=false;
798  }// if isfirstcall
799 
800 
801 
802 
803 
804 
805 
806 
807 
809 
810  for(unsigned int iHLT=0; iHLT<200; ++iHLT) {
811  aHLTResults[iHLT] = false;
812  }
813 
814 #ifdef DEBUG
815  std::cout<<"[ZeeCalibration::duringLoop] Done with initializing aHLTresults[] "<<std::endl;
816 #endif
817 
818  edm::Handle<edm::TriggerResults> hltTriggerResultHandle;
819  iEvent.getByLabel(hlTriggerResults_, hltTriggerResultHandle);
820 
821  if(!hltTriggerResultHandle.isValid()) {
822  //std::cout << "invalid handle for HLT TriggerResults" << std::endl;
823  } else {
824 
825  hltCount = hltTriggerResultHandle->size();
826 
827  if (loopFlag_ == 0)
828  myZeePlots_->fillHLTInfo(hltTriggerResultHandle);
829 
830 #ifdef DEBUG
831  std::cout<<"[ZeeCalibration::duringLoop] Done with myZeePlots_->fillHLTInfo(hltTriggerResultHandle); "<<std::endl;
832 #endif
833 
834  for(int i = 0 ; i < hltCount ; i++) {
835  aHLTResults[i] = hltTriggerResultHandle->accept(i);
836 
837  //HLT bit 32 = HLT1Electron
838  //HLT bit 34 = HLT2Electron
839  //HLT bit 35 = HLT2ElectronRelaxed
840 
841  }
842 
843  if(!aHLTResults[32] && !aHLTResults[34] && !aHLTResults[35])
844  return kContinue;
845 
846  }
847 
848 #ifdef DEBUG
849  std::cout<<"[ZeeCalibration::duringLoop] End HLT section"<<std::endl;
850 #endif
851 
853 
854 
855  std::vector<HepMC::GenParticle*> mcEle;
856 
857  float myGenZMass(-1);
858 
859  if (!mcProducer_.empty())
860  {
861 
862  //DUMP GENERATED Z MASS - BEGIN
863  Handle< HepMCProduct > hepProd ;
864  // iEvent.getByLabel( "source", hepProd ) ;
865  iEvent.getByLabel( mcProducer_.c_str(), hepProd ) ;
866 
867  const HepMC::GenEvent * myGenEvent = hepProd->GetEvent();
868 
869  if (loopFlag_ == 0)
870  myZeePlots_->fillZMCInfo( & (*myGenEvent) );
871 
872 #ifdef DEBUG
873  std::cout<<"[ZeeCalibration::duringLoop] Done with myZeePlots_->fillZMCInfo( & (*myGenEvent) ); "<<std::endl;
874 #endif
875 
876 
877  for ( HepMC::GenEvent::particle_const_iterator p = myGenEvent->particles_begin();
878  p != myGenEvent->particles_end(); ++p ) {
879  //return a pointer to MC Z in the event
880  if ( (*p)->pdg_id() == 23 && (*p)->status()==2){
881 
882  myGenZMass = (*p)->momentum().m();
883  }
884  }
885  //DUMP GENERATED Z MASS - END
886 
887 
888  if (loopFlag_ == 0)
889  myZeePlots_ ->fillEleMCInfo( & (*myGenEvent) );
890 
891 
892  //loop over MC positrons and find the closest MC positron in (eta,phi) phace space - begin
893  HepMC::GenParticle MCele;
894 
895  for ( HepMC::GenEvent::particle_const_iterator p = myGenEvent->particles_begin();
896  p != myGenEvent->particles_end(); ++p ) {
897 
898  if ( abs( (*p)->pdg_id() ) == 11 )
899  {
900  mcEle.push_back( (*p) );
901  MCele=*(*p);
902 
903  }
904  }
905 
906 
907  if(mcEle.size()==2 && fabs(mcEle[0]->momentum().eta())<2.4 && fabs(mcEle[1]->momentum().eta())<2.4 ){
908  NEVT++;
909 
910  if( fabs(mcEle[0]->momentum().eta())<1.479 && fabs(mcEle[1]->momentum().eta())<1.479 )MCZBB++;
911 
912  if( (fabs(mcEle[0]->momentum().eta())>1.479 && fabs(mcEle[1]->momentum().eta())<1.479) || (fabs(mcEle[0]->momentum().eta())<1.479 && fabs(mcEle[1]->momentum().eta())>1.479) )MCZEB++;
913 
914  if( fabs(mcEle[0]->momentum().eta())>1.479 && fabs(mcEle[1]->momentum().eta())>1.479 )MCZEE++;
915 
916 
917  }
918 
919  }
920 
921 
922 
923  // Get EBRecHits
924  Handle<EBRecHitCollection> phits;
925  try {
926  iEvent.getByLabel( rechitProducer_, rechitCollection_, phits);
927  } catch (std::exception& ex) {
928  std::cerr << "Error! can't get the product EBRecHitCollection " << std::endl;
929  }
930  const EBRecHitCollection* hits = phits.product(); // get a ptr to the product
931 
932  // Get EERecHits
933  Handle<EERecHitCollection> ephits;
934  try {
935  iEvent.getByLabel( erechitProducer_, erechitCollection_, ephits);
936  } catch (std::exception& ex) {
937  std::cerr << "Error! can't get the product EERecHitCollection " << std::endl;
938  }
939  const EERecHitCollection* ehits = ephits.product(); // get a ptr to the product
940 
941 
942  //Get Hybrid SuperClusters
943  Handle<reco::SuperClusterCollection> pSuperClusters;
944  try {
945  iEvent.getByLabel(scProducer_, scCollection_, pSuperClusters);
946  } catch (std::exception& ex ) {
947  std::cerr << "Error! can't get the product SuperClusterCollection "<< std::endl;
948  }
949  const reco::SuperClusterCollection* scCollection = pSuperClusters.product();
950 
951 #ifdef DEBUG
952  std::cout<<"scCollection->size()"<<scCollection->size()<<std::endl;
953  for(reco::SuperClusterCollection::const_iterator scIt = scCollection->begin(); scIt != scCollection->end(); scIt++)
954  {
955  std::cout<<scIt->energy()<<std::endl;
956  }
957 #endif
958 
959  //Get Island SuperClusters
960  Handle<reco::SuperClusterCollection> pIslandSuperClusters;
961  try {
962  iEvent.getByLabel(scIslandProducer_, scIslandCollection_, pIslandSuperClusters);
963  } catch (std::exception& ex ) {
964  std::cerr << "Error! can't get the product IslandSuperClusterCollection "<< std::endl;
965  }
966  const reco::SuperClusterCollection* scIslandCollection = pIslandSuperClusters.product();
967 
968 #ifdef DEBUG
969  std::cout<<"scCollection->size()"<<scIslandCollection->size()<<std::endl;
970 #endif
971 
972  if( ( scCollection->size()+scIslandCollection->size() ) < 2)
973  return kContinue;
974 
975  // Get Electrons
976  Handle<reco::GsfElectronCollection> pElectrons;
977  try {
978  iEvent.getByLabel(electronProducer_, electronCollection_, pElectrons);
979  } catch (std::exception& ex ) {
980  std::cerr << "Error! can't get the product ElectronCollection "<< std::endl;
981  }
982  const reco::GsfElectronCollection* electronCollection = pElectrons.product();
983 
984  /*
985  //reco-mc association map
986  std::map<HepMC::GenParticle*,const reco::PixelMatchGsfElectron*> myMCmap;
987 
988  fillMCmap(&(*electronCollection),mcEle,myMCmap);
989 
990  fillEleInfo(mcEle,myMCmap);
991  */
992 
993  if(electronCollection->size() < 2)
994  return kContinue;
995 
996  if ( !hits && !ehits){
997  std::cout << "!hits" << std::endl;
998  return kContinue;
999  }
1000 
1001  if (hits->size() == 0 && ehits->size() == 0){
1002  std::cout << "hits->size() == 0" << std::endl;
1003  return kContinue;
1004  }
1005 
1006  if (!electronCollection){
1007  std::cout << "!electronCollection" << std::endl;
1008  return kContinue;
1009  }
1010 
1011  if (electronCollection->size() == 0){
1012  std::cout << "electronCollection->size() == 0" << std::endl;
1013  return kContinue;
1014  }
1015 
1016 
1017 
1021 
1022  read_events++;
1023 
1024  // std::cout << "read_events = " << read_events << std::endl;
1025 
1027 
1028 #ifdef DEBUG
1029  std::cout <<" Starting with myZeePlots_->fillEleInfo(electronCollection); " << std::endl;
1030 #endif
1031 
1032  if (loopFlag_ == 0)
1033  myZeePlots_->fillEleInfo(electronCollection);
1034 
1035 #ifdef DEBUG
1036  std::cout <<" Done with myZeePlots_->fillEleInfo(electronCollection); " << std::endl;
1037 #endif
1038 
1039  //FILL an electron vector - end
1040  //###################################Electron-SC association: begin#####################################################
1041  //Filling new ElectronCollection with new SC ref and calibElectron container
1042  std::vector<calib::CalibElectron> calibElectrons;
1043  //std::map< const calib::CalibElectron* , const reco::SuperCluster* > eleScMap;
1044 
1045 
1046 
1047  //#####################################Electron-SC association map: end#####################################################
1048  for(unsigned int e_it = 0 ; e_it != electronCollection->size() ; e_it++)
1049  {
1050  calibElectrons.push_back(calib::CalibElectron(&((*electronCollection)[e_it]),hits,ehits));
1051 #ifdef DEBUG
1052  std::cout << calibElectrons.back().getRecoElectron()->superCluster()->energy() << " " << calibElectrons.back().getRecoElectron()->energy() << std::endl;
1053 #endif
1054  // h1_recoEleEnergy_->Fill(calibElectrons.back().getRecoElectron()->superCluster()->energy());
1055  }
1056  // if (iLoop == 0)
1057  //fillCalibElectrons(calibElectrons);
1058 
1059 #ifdef DEBUG
1060  std::cout << "Filled histos" << std::endl;
1061 #endif
1062 
1063  //COMBINATORY FOR Z MASS - begin
1064  std::vector<std::pair<calib::CalibElectron*,calib::CalibElectron*> > zeeCandidates;
1065  int myBestZ=-1;
1066 
1067  mass = -1.;
1068  double DeltaMinvMin(5000.);
1069 
1070  if (calibElectrons.size() < 2)
1071  return kContinue;
1072 
1073  for(unsigned int e_it = 0 ; e_it != calibElectrons.size() - 1 ; e_it++){
1074  for(unsigned int p_it = e_it + 1 ; p_it != calibElectrons.size() ; p_it++)
1075  {
1076 #ifdef DEBUG
1077  std::cout << e_it << " " << calibElectrons[e_it].getRecoElectron()->charge() << " " << p_it << " " << calibElectrons[p_it].getRecoElectron()->charge() << std::endl;
1078 #endif
1079  if (calibElectrons[e_it].getRecoElectron()->charge() * calibElectrons[p_it].getRecoElectron()->charge() != -1)
1080  continue;
1081 
1082  mass = ZeeKinematicTools::calculateZMass_withTK(std::pair<calib::CalibElectron*,calib::CalibElectron*>(&(calibElectrons[e_it]),&(calibElectrons[p_it])));
1083 
1084  if (mass<0)
1085  continue;
1086 
1087 #ifdef DEBUG
1088  std::cout << "#######################mass "<<mass << std::endl;
1089 #endif
1090 
1091  zeeCandidates.push_back(std::pair<calib::CalibElectron*,calib::CalibElectron*>(&(calibElectrons[e_it]),&(calibElectrons[p_it])));
1092  double DeltaMinv = fabs(mass - MZ);
1093 
1094  if( DeltaMinv < DeltaMinvMin)
1095  {
1096  DeltaMinvMin = DeltaMinv;
1097  myBestZ=zeeCandidates.size()-1;
1098  }
1099  }
1100  }
1101 
1102  // h_DeltaZMassDistr_[loopFlag_]->Fill( (mass-MZ) / MZ );
1103 
1104  // zeeCa->Fill(zeeCandidates);
1105  //
1106  h1_ZCandMult_->Fill(zeeCandidates.size());
1107 
1108  if(zeeCandidates.size()==0 || myBestZ==-1 )
1109  return kContinue;
1110 
1111  if (loopFlag_ == 0)
1112  myZeePlots_->fillZInfo( zeeCandidates[myBestZ] );
1113 
1114 #ifdef DEBUG
1115  std::cout << "Found ZCandidates " << myBestZ << std::endl;
1116 #endif
1117 
1118  // h1_zMassResol_ ->Fill(mass-myGenZMass);
1119 
1121 
1122 
1123  h1_eleClasses_->Fill(zeeCandidates[myBestZ].first->getRecoElectron()->classification());
1124  h1_eleClasses_->Fill(zeeCandidates[myBestZ].second->getRecoElectron()->classification());
1125 
1126  int class1 = zeeCandidates[myBestZ].first->getRecoElectron()->classification();
1127  int class2 = zeeCandidates[myBestZ].second->getRecoElectron()->classification();
1128 
1129  std::cout << "BEFORE "<<std::endl;
1130 
1131  // myZeePlots_->fillEleClassesPlots( zeeCandidates[myBestZ].first );
1132  //myZeePlots_->fillEleClassesPlots( zeeCandidates[myBestZ].second );
1133 
1134  std::cout << "AFTER "<<std::endl;
1135 
1137 
1138  if(class1 < 100)
1139  // h1_Elec_->Fill(1);
1140  TOTAL_ELECTRONS_IN_BARREL++;
1141  if(class1 >= 100)
1142  TOTAL_ELECTRONS_IN_ENDCAP++;
1143 
1144  if(class2 < 100)
1145  TOTAL_ELECTRONS_IN_BARREL++;
1146  if(class2 >= 100)
1147  TOTAL_ELECTRONS_IN_ENDCAP++;
1148 
1149 
1150  if( class1==0)
1151  GOLDEN_ELECTRONS_IN_BARREL++;
1152  if( class1==100)
1153  GOLDEN_ELECTRONS_IN_ENDCAP++;
1154  if( class1==10 || class1 ==20)
1155  SILVER_ELECTRONS_IN_BARREL++;
1156  if( class1==110 || class1 ==120)
1157  SILVER_ELECTRONS_IN_ENDCAP++;
1158  if( class1>=30 && class1 <=34)
1159  SHOWER_ELECTRONS_IN_BARREL++;
1160  if( class1>=130 && class1 <=134)
1161  SHOWER_ELECTRONS_IN_ENDCAP++;
1162  if( class1==40)
1163  CRACK_ELECTRONS_IN_BARREL++;
1164  if( class1==140)
1165  CRACK_ELECTRONS_IN_ENDCAP++;
1166 
1167  if( class2==0)
1168  GOLDEN_ELECTRONS_IN_BARREL++;
1169  if( class2==100)
1170  GOLDEN_ELECTRONS_IN_ENDCAP++;
1171  if( class2==10 || class2 ==20)
1172  SILVER_ELECTRONS_IN_BARREL++;
1173  if( class2==110 || class2 ==120)
1174  SILVER_ELECTRONS_IN_ENDCAP++;
1175  if( class2>=30 && class2 <=34)
1176  SHOWER_ELECTRONS_IN_BARREL++;
1177  if( class2>=130 && class2 <=134)
1178  SHOWER_ELECTRONS_IN_ENDCAP++;
1179  if( class2==40)
1180  CRACK_ELECTRONS_IN_BARREL++;
1181  if( class2==140)
1182  CRACK_ELECTRONS_IN_ENDCAP++;
1183 
1185 
1187 
1188 
1189  DetId firstElehottestDetId = getHottestDetId( zeeCandidates[myBestZ].first->getRecoElectron()->superCluster()->seed()->hitsAndFractions() , hits, ehits ).first;
1190  DetId secondElehottestDetId = getHottestDetId( zeeCandidates[myBestZ].second->getRecoElectron()->superCluster()->seed()->hitsAndFractions() , hits, ehits ).first;
1191 
1192  bool firstElectronIsOnModuleBorder(false);
1193  bool secondElectronIsOnModuleBorder(false);
1194 
1195  h1_eventsBeforeBorderSelection_->Fill(1);
1196 
1197  if(class1<100){
1198 
1199  if( firstElehottestDetId.subdetId() == EcalBarrel)
1200  firstElectronIsOnModuleBorder = xtalIsOnModuleBorder( firstElehottestDetId );
1201 
1202  BARREL_ELECTRONS_BEFORE_BORDER_CUT++;
1203 
1204  if( firstElehottestDetId.subdetId() == EcalBarrel && !firstElectronIsOnModuleBorder )
1205  BARREL_ELECTRONS_AFTER_BORDER_CUT++;
1206 
1207  }
1208 
1209  if(class2<100){
1210 
1211  if( secondElehottestDetId.subdetId() == EcalBarrel)
1212  secondElectronIsOnModuleBorder = xtalIsOnModuleBorder( secondElehottestDetId );
1213 
1214  BARREL_ELECTRONS_BEFORE_BORDER_CUT++;
1215 
1216  if( secondElehottestDetId.subdetId() == EcalBarrel && !secondElectronIsOnModuleBorder )
1217  BARREL_ELECTRONS_AFTER_BORDER_CUT++;
1218 
1219  }
1220 
1221 
1222  if(class1<100){
1223  if ( firstElehottestDetId.subdetId() == EcalBarrel && firstElectronIsOnModuleBorder ){
1224  h1_borderElectronClassification_ -> Fill( zeeCandidates[myBestZ].first->getRecoElectron()->classification() );
1225  return kContinue;
1226  }
1227  }
1228 
1229  if(class2<100){
1230  if ( secondElehottestDetId.subdetId() == EcalBarrel && secondElectronIsOnModuleBorder ){
1231  h1_borderElectronClassification_ -> Fill( zeeCandidates[myBestZ].second->getRecoElectron()->classification() );
1232  return kContinue;
1233  }
1234  }
1235 
1236  h1_eventsAfterBorderSelection_->Fill(1);
1238 
1239 
1240  if(class1<100 && class2<100){
1241  BBZN++;
1242  if(class1==0 && class2==0)BBZN_gg++;
1243  if(class1<21 && class2<21)BBZN_tt++;
1244  if(class1<21 || class2<21)BBZN_t0++;
1245 
1246  }
1247 
1248  if(class1>=100 && class2>=100){
1249  EEZN++;
1250  if(class1==100 && class2==100)EEZN_gg++;
1251  if(class1<121 && class2<121)EEZN_tt++;
1252  if(class1<121 || class2<121)EEZN_t0++;
1253 
1254  }
1255 
1256  if( (class1<100 && class2>=100) || (class2<100 && class1>=100)){
1257  EBZN++;
1258  if( (class1==0 && class2==100)||(class2==0 && class1==100) )EBZN_gg++;
1259  if( ( class1<21 && class2<121) ||(class2<21 && class1<121) )EBZN_tt++;
1260  if( class2<21 || class1<21 || class2<121 || class1<121 )EBZN_t0++;
1261  }
1262 
1263 
1265 
1266  if(myBestZ == -1)
1267  return kContinue;
1268 
1269 
1270  bool invMassBool = ( (mass > minInvMassCut_) && (mass < maxInvMassCut_) );
1271 
1272  bool selectionBool=false;
1273  //0 = all electrons (but no crack)
1274 
1276 
1277  float theta1 = 2. * atan( exp(- zeeCandidates[myBestZ].first->getRecoElectron()->superCluster()->eta()) );
1278  bool ET_1 = ( (zeeCandidates[myBestZ].first->getRecoElectron()->superCluster()->energy() * sin( theta1) ) > 20.);
1279 
1280  float theta2 = 2. * atan( exp(- zeeCandidates[myBestZ].second->getRecoElectron()->superCluster()->eta()) );
1281  bool ET_2 = ( (zeeCandidates[myBestZ].second->getRecoElectron()->superCluster()->energy() * sin( theta2) ) > 20.);
1282 
1283 
1284  bool HoE_1 = (zeeCandidates[myBestZ].first->getRecoElectron()->hadronicOverEm() < 0.115);
1285  bool HoE_2 = (zeeCandidates[myBestZ].second->getRecoElectron()->hadronicOverEm() < 0.115);
1286 
1287  bool DeltaPhiIn_1 = ( zeeCandidates[myBestZ].first->getRecoElectron()->deltaPhiSuperClusterTrackAtVtx() < 0.090);
1288  bool DeltaPhiIn_2 = ( zeeCandidates[myBestZ].second->getRecoElectron()->deltaPhiSuperClusterTrackAtVtx() < 0.090);
1289 
1290  bool DeltaEtaIn_1 = ( zeeCandidates[myBestZ].first->getRecoElectron()->deltaEtaSuperClusterTrackAtVtx() < 0.0090);
1291  bool DeltaEtaIn_2 = ( zeeCandidates[myBestZ].second->getRecoElectron()->deltaEtaSuperClusterTrackAtVtx() < 0.0090);
1292 
1293  h1_eventsBeforeEWKSelection_->Fill(1);
1294 
1295  if(! (invMassBool &&
1296  ET_1 && ET_2 &&
1297  HoE_1 && HoE_2 &&
1298  DeltaPhiIn_1 && DeltaPhiIn_2 &&
1299  DeltaEtaIn_1 && DeltaEtaIn_2
1300  ) )
1301  return kContinue;
1303 
1304  h1_eventsAfterEWKSelection_->Fill(1);
1305 
1307 
1308 
1309 
1310  if(electronSelection_==0)selectionBool=( myBestZ != -1 &&
1311  zeeCandidates[myBestZ].first->getRecoElectron()->classification()!= 40 &&
1312  zeeCandidates[myBestZ].first->getRecoElectron()->classification()!= 40 &&
1313  zeeCandidates[myBestZ].second->getRecoElectron()->classification()!= 40 &&
1314  zeeCandidates[myBestZ].second->getRecoElectron()->classification()!= 140);
1315 
1316  //1 = all electrons are Golden, BB or Narrow
1317 
1318  if(electronSelection_==1)selectionBool=( myBestZ != -1 &&
1319  (zeeCandidates[myBestZ].first->getRecoElectron()->classification() ==0 ||
1320  zeeCandidates[myBestZ].first->getRecoElectron()->classification() ==10 ||
1321  zeeCandidates[myBestZ].first->getRecoElectron()->classification() ==20 ||
1322  zeeCandidates[myBestZ].first->getRecoElectron()->classification() ==100 ||
1323  zeeCandidates[myBestZ].first->getRecoElectron()->classification() ==110 ||
1324  zeeCandidates[myBestZ].first->getRecoElectron()->classification() ==120
1325  ) &&
1326  (zeeCandidates[myBestZ].second->getRecoElectron()->classification() == 0 ||
1327  zeeCandidates[myBestZ].second->getRecoElectron()->classification() == 10 ||
1328  zeeCandidates[myBestZ].second->getRecoElectron()->classification() == 20 ||
1329  zeeCandidates[myBestZ].second->getRecoElectron()->classification() == 100 ||
1330  zeeCandidates[myBestZ].second->getRecoElectron()->classification() == 110 ||
1331  zeeCandidates[myBestZ].second->getRecoElectron()->classification() == 120
1332  )
1333  );
1334 
1335  //2 = all electrons are Golden
1336  if(electronSelection_==2)selectionBool=( myBestZ != -1 &&
1337  (zeeCandidates[myBestZ].first->getRecoElectron()->classification() == 0 ||
1338  zeeCandidates[myBestZ].first->getRecoElectron()->classification() == 100
1339  ) &&
1340  (zeeCandidates[myBestZ].second->getRecoElectron()->classification() == 0 ||
1341  zeeCandidates[myBestZ].second->getRecoElectron()->classification() == 100
1342  )
1343  );
1344  //3 = all electrons are showering
1345  if(electronSelection_==3)selectionBool=( myBestZ != -1 &&
1346  (
1347  (zeeCandidates[myBestZ].first->getRecoElectron()->classification() >=30 &&
1348  zeeCandidates[myBestZ].first->getRecoElectron()->classification() <=34)
1349  ||
1350  ((zeeCandidates[myBestZ].first->getRecoElectron()->classification() >=130 &&
1351  zeeCandidates[myBestZ].first->getRecoElectron()->classification() <=134))
1352  )
1353  &&
1354  ( (zeeCandidates[myBestZ].second->getRecoElectron()->classification() >=30 &&
1355  zeeCandidates[myBestZ].second->getRecoElectron()->classification() <=34)
1356  ||
1357  ((zeeCandidates[myBestZ].second->getRecoElectron()->classification() >=130 &&
1358  zeeCandidates[myBestZ].second->getRecoElectron()->classification() <=134))
1359  )
1360 
1361  );
1362 
1363  //4 = all Barrel electrons are Golden, BB or Narrow; take all Endcap electrons
1364 
1365  if(electronSelection_==4)selectionBool=( myBestZ != -1 &&
1366  (
1367 
1368  (
1369  (zeeCandidates[myBestZ].first->getRecoElectron()->classification() ==0 ||
1370  zeeCandidates[myBestZ].first->getRecoElectron()->classification() ==10 ||
1371  zeeCandidates[myBestZ].first->getRecoElectron()->classification() ==20
1372  ) && zeeCandidates[myBestZ].second->getRecoElectron()->classification()>=100
1373  && zeeCandidates[myBestZ].second->getRecoElectron()->classification()!=140
1374  )
1375 
1376  ||
1377 
1378  (
1379  (zeeCandidates[myBestZ].second->getRecoElectron()->classification() ==0 ||
1380  zeeCandidates[myBestZ].second->getRecoElectron()->classification() ==10 ||
1381  zeeCandidates[myBestZ].second->getRecoElectron()->classification() ==20
1382  ) && zeeCandidates[myBestZ].first->getRecoElectron()->classification()>=100
1383  && zeeCandidates[myBestZ].first->getRecoElectron()->classification()!=140
1384  )
1385 
1386 
1387  )
1388  );
1389 
1390  //5 = all Endcap electrons (but no crack)
1391 
1392  if(electronSelection_==5)selectionBool=( myBestZ != -1 &&
1393  zeeCandidates[myBestZ].first->getRecoElectron()->classification()>=100 &&
1394  zeeCandidates[myBestZ].second->getRecoElectron()->classification()>= 100 &&
1395  zeeCandidates[myBestZ].first->getRecoElectron()->classification()!= 140 &&
1396  zeeCandidates[myBestZ].second->getRecoElectron()->classification()!= 140);
1397 
1398  //6 = all Barrel electrons (but no crack)
1399 
1400  if(electronSelection_==6)selectionBool=( myBestZ != -1 &&
1401  zeeCandidates[myBestZ].first->getRecoElectron()->classification()<100 &&
1402  zeeCandidates[myBestZ].second->getRecoElectron()->classification()< 100 &&
1403  zeeCandidates[myBestZ].first->getRecoElectron()->classification()!= 40 &&
1404  zeeCandidates[myBestZ].second->getRecoElectron()->classification()!= 40);
1405 
1406  //7 = this eliminates the events which have 1 ele in the Barrel and 1 in the Endcap
1407 
1408  if(electronSelection_==7)selectionBool=( myBestZ != -1 &&
1409  !(zeeCandidates[myBestZ].first->getRecoElectron()->classification()<100 &&
1410  zeeCandidates[myBestZ].second->getRecoElectron()->classification()>=100) &&
1411  !(zeeCandidates[myBestZ].first->getRecoElectron()->classification()>=100 &&
1412  zeeCandidates[myBestZ].second->getRecoElectron()->classification()<100) );
1413 
1414 
1415  float ele1EnergyCorrection(1.);
1416  float ele2EnergyCorrection(1.);
1417 
1418  if(invMassBool && selectionBool && wantEtaCorrection_){
1419 
1420  ele1EnergyCorrection=getEtaCorrection(zeeCandidates[myBestZ].first->getRecoElectron());
1421  ele2EnergyCorrection=getEtaCorrection(zeeCandidates[myBestZ].second->getRecoElectron());
1422 
1423  }
1424 
1425  if (invMassBool && selectionBool)
1426  {
1427 
1428  h1_electronCosTheta_SC_ -> Fill( ZeeKinematicTools::cosThetaElectrons_SC(zeeCandidates[myBestZ],ele1EnergyCorrection,ele2EnergyCorrection) );
1429  h1_electronCosTheta_TK_ -> Fill( ZeeKinematicTools::cosThetaElectrons_TK(zeeCandidates[myBestZ],ele1EnergyCorrection,ele2EnergyCorrection) );
1430  h1_electronCosTheta_SC_TK_ -> Fill( ZeeKinematicTools::cosThetaElectrons_SC(zeeCandidates[myBestZ],ele1EnergyCorrection,ele2EnergyCorrection)/ZeeKinematicTools::cosThetaElectrons_TK(zeeCandidates[myBestZ],ele1EnergyCorrection,ele2EnergyCorrection) - 1. );
1431 
1432  if (!mcProducer_.empty())
1433  {
1434  h1_zMassResol_ ->Fill(mass-myGenZMass);
1435 
1436  //reco-mc association map - begin
1437 
1438  std::map<HepMC::GenParticle*,const reco::GsfElectron*> myMCmap;
1439 
1440  std::vector<const reco::GsfElectron*> dauElectronCollection;
1441 
1442  dauElectronCollection.push_back(zeeCandidates[myBestZ].first->getRecoElectron() );
1443  dauElectronCollection.push_back(zeeCandidates[myBestZ].second->getRecoElectron() );
1444 
1445  fillMCmap(&dauElectronCollection,mcEle,myMCmap);
1446  fillEleInfo(mcEle,myMCmap);
1447  //h_DiffZMassDistr_[loopFlag_]->Fill( (mass-myGenZMass) );
1448  }
1449 
1450  //PUT f(eta) IN OUR Zee ALGORITHM
1451  theAlgorithm_->addEvent(zeeCandidates[myBestZ].first, zeeCandidates[myBestZ].second,MZ*sqrt(ele1EnergyCorrection*ele2EnergyCorrection) );
1452 
1453  h1_reco_ZMass_->Fill(ZeeKinematicTools::calculateZMass_withTK(zeeCandidates[myBestZ]));
1454 
1455  h1_reco_ZMassCorr_->Fill(ZeeKinematicTools::calculateZMassWithCorrectedElectrons_withTK(zeeCandidates[myBestZ],ele1EnergyCorrection,ele2EnergyCorrection));
1456 
1457  if(zeeCandidates[myBestZ].first->getRecoElectron()->classification()<100 && zeeCandidates[myBestZ].second->getRecoElectron()->classification()<100 )
1458  h1_reco_ZMassCorrBB_->Fill(ZeeKinematicTools::calculateZMassWithCorrectedElectrons_withTK(zeeCandidates[myBestZ],ele1EnergyCorrection,ele2EnergyCorrection));
1459 
1460 
1461  if(zeeCandidates[myBestZ].first->getRecoElectron()->classification()>=100 && zeeCandidates[myBestZ].second->getRecoElectron()->classification()>=100 )
1462  h1_reco_ZMassCorrEE_->Fill(ZeeKinematicTools::calculateZMassWithCorrectedElectrons_withTK(zeeCandidates[myBestZ],ele1EnergyCorrection,ele2EnergyCorrection));
1463 
1464 
1465  mass4tree = ZeeKinematicTools::calculateZMassWithCorrectedElectrons_withTK(zeeCandidates[myBestZ],ele1EnergyCorrection,ele2EnergyCorrection);
1466 
1467  massDiff4tree = ZeeKinematicTools::calculateZMassWithCorrectedElectrons_withTK(zeeCandidates[myBestZ],ele1EnergyCorrection,ele2EnergyCorrection) - myGenZMass;
1468 
1469  // h_ZMassDistr_[loopFlag_]->Fill(ZeeKinematicTools::calculateZMassWithCorrectedElectrons_withTK(zeeCandidates[myBestZ],ele1EnergyCorrection,ele2EnergyCorrection));
1470 
1471  myTree->Fill();
1472 
1473  }
1474 
1475 
1476 
1477 #ifdef DEBUG
1478  std::cout << "Added event to algorithm" << std::endl;
1479 #endif
1480 
1481  return kContinue;
1482  }
1483 //end of ZeeCalibration::duringLoop
1484 
1485 
1486 // Called at beginning of loop
1487 void ZeeCalibration::startingNewLoop ( unsigned int iLoop )
1488 {
1489 
1490 std::cout<< "[ZeeCalibration] Starting loop number " << iLoop<<std::endl;
1491 
1492  theAlgorithm_->resetIteration();
1493 
1494  resetVariables();
1495 
1496  resetHistograms();
1497 
1498 #ifdef DEBUG
1499  std::cout<< "[ZeeCalibration] exiting from startingNewLoop" << std::endl;
1500 #endif
1501 
1502 }
1503 
1504 
1505 
1506 // Called at end of loop
1507 edm::EDLooper::Status
1508 ZeeCalibration::endOfLoop(const edm::EventSetup& iSetup, unsigned int iLoop)
1509 {
1510 
1511 
1512 
1513  double par[3];
1514  double errpar[3];
1515  double zChi2;
1516  int zIters;
1517 
1518  ZIterativeAlgorithmWithFit::gausfit(h1_reco_ZMass_,par,errpar,2.,2., &zChi2, &zIters );
1519 
1520  h2_zMassVsLoop_ -> Fill(loopFlag_, par[1] );
1521 
1522  h2_zMassDiffVsLoop_ -> Fill(loopFlag_, (par[1]-MZ)/MZ );
1523 
1524  h2_zWidthVsLoop_ -> Fill(loopFlag_, par[2] );
1525 
1526 
1528 
1529  std::cout<< "[ZeeCalibration] Ending loop " << iLoop<<std::endl;
1530  //RUN the algorithm
1531  theAlgorithm_->iterate();
1532 
1533  const std::vector<float>& optimizedCoefficients = theAlgorithm_->getOptimizedCoefficients();
1534  const std::vector<float>& optimizedCoefficientsError = theAlgorithm_->getOptimizedCoefficientsError();
1535  //const std::vector<float>& weightSum = theAlgorithm_->getWeightSum();
1536  const std::vector<float>& optimizedChi2 = theAlgorithm_->getOptimizedChiSquare();
1537  const std::vector<int>& optimizedIterations = theAlgorithm_->getOptimizedIterations();
1538 
1539 
1540  //#ifdef DEBUG
1541  std::cout<< "Optimized coefficients " << optimizedCoefficients.size() <<std::endl;
1542  //#endif
1543 
1544  // h2_coeffVsLoop_->Fill(loopFlag_, optimizedCoefficients[75]); //show the evolution of just 1 ring coefficient (well chosen...)
1545 
1547  for (unsigned int ieta=0;ieta<optimizedCoefficients.size();ieta++)
1548  {
1549  NewCalibCoeff[ieta] = calibCoeff[ieta] * optimizedCoefficients[ieta];
1550 
1551  h2_chi2_[loopFlag_]->Fill( ringNumberCorrector( ieta ), optimizedChi2[ieta] );
1552  h2_iterations_[loopFlag_]->Fill( ringNumberCorrector( ieta ), optimizedIterations[ieta] );
1553 
1554  }
1556 
1557 
1558  coefficientDistanceAtIteration[loopFlag_]= computeCoefficientDistanceAtIteration(calibCoeff, NewCalibCoeff, optimizedCoefficients.size() );
1559 
1560  std::cout<<"Iteration # : "<< loopFlag_ << " CoefficientDistanceAtIteration "<< coefficientDistanceAtIteration[loopFlag_] <<std::endl;
1561  std::cout<<"size "<<optimizedCoefficients.size()<<std::endl;
1562 
1563  for (unsigned int ieta=0;ieta<optimizedCoefficients.size();ieta++)
1564  {
1565  calibCoeff[ieta] *= optimizedCoefficients[ieta];
1566  calibCoeffError[ieta] = calibCoeff[ieta] * sqrt ( pow( optimizedCoefficientsError[ieta]/optimizedCoefficients[ieta], 2 ) + pow( calibCoeffError[ieta]/calibCoeff[ieta] , 2 ) );
1567  //calibCoeffError[ieta] = optimizedCoefficientsError[ieta];
1568 
1569 
1570 #ifdef DEBUG
1571  std::cout<< ieta << " " << optimizedCoefficients[ieta] <<std::endl;
1572 #endif
1573 
1574  std::vector<DetId> ringIds;
1575 
1576  if(calibMode_ == "RING")
1577  ringIds = EcalRingCalibrationTools::getDetIdsInRing(ieta);
1578 
1579  if(calibMode_ == "MODULE")
1580  ringIds = EcalRingCalibrationTools::getDetIdsInModule(ieta);
1581 
1582  if(calibMode_ == "ABS_SCALE" || calibMode_ == "ETA_ET_MODE" )
1583  ringIds = EcalRingCalibrationTools::getDetIdsInECAL();
1584 
1585 
1586  for (unsigned int iid=0; iid<ringIds.size();++iid){
1587 
1588  if(ringIds[iid].subdetId() == EcalBarrel){
1589  EBDetId myEBDetId(ringIds[iid]);
1590  h2_xtalRecalibCoeffBarrel_[loopFlag_]->SetBinContent( myEBDetId.ieta() + 86, myEBDetId.iphi(), 100 * (calibCoeff[ieta]*initCalibCoeff[ieta] - 1.) );//fill TH2 with recalibCoeff
1591 
1592  }
1593 
1594  if(ringIds[iid].subdetId() == EcalEndcap){
1595  EEDetId myEEDetId(ringIds[iid]);
1596  if(myEEDetId.zside() < 0)
1597  h2_xtalRecalibCoeffEndcapMinus_[loopFlag_]->SetBinContent( myEEDetId.ix(), myEEDetId.iy(), 100 * (calibCoeff[ieta]*initCalibCoeff[ieta] - 1.) );//fill TH2 with recalibCoeff
1598 
1599  if(myEEDetId.zside() > 0)
1600  h2_xtalRecalibCoeffEndcapPlus_[loopFlag_]->SetBinContent( myEEDetId.ix(), myEEDetId.iy(), 100 * (calibCoeff[ieta]*initCalibCoeff[ieta] - 1.) );//fill TH2 with recalibCoeff
1601 
1602  }
1603 
1604 
1605  ical->setValue( ringIds[iid], *(ical->getMap().find(ringIds[iid]) ) * optimizedCoefficients[ieta] );
1606  }
1607 
1608  }
1609 
1610 
1612 
1613  for( int k = 0; k<theAlgorithm_->getNumberOfChannels(); k++ )
1614  {
1615  bool isNearCrack = ( abs( ringNumberCorrector(k) ) == 1 || abs( ringNumberCorrector(k) ) == 25 ||
1616  abs( ringNumberCorrector(k) ) == 26 || abs( ringNumberCorrector(k) ) == 45 ||
1617  abs( ringNumberCorrector(k) ) == 46 || abs( ringNumberCorrector(k) ) == 65 ||
1618  abs( ringNumberCorrector(k) ) == 66 || abs( ringNumberCorrector(k) ) == 85 ||
1619  abs( ringNumberCorrector(k) ) == 86 || abs( ringNumberCorrector(k) ) == 124 );
1620 
1621  if(!isNearCrack){
1622 
1623  // h2_miscalRecalParz_[iLoop]->Fill( initCalibCoeff[k], 1./calibCoeff[k] );
1624  h1_mcParz_[iLoop]->Fill( initCalibCoeff[k]*calibCoeff[k] -1. );
1625 
1626  if(k<170){
1627  //h2_miscalRecalEBParz_[iLoop]->Fill( initCalibCoeff[k], 1./calibCoeff[k] );
1628  h1_mcEBParz_[iLoop]->Fill( initCalibCoeff[k]*calibCoeff[k] -1. );
1629 
1630  }
1631 
1632  if(k>=170){
1633  //h2_miscalRecalEEParz_[iLoop]->Fill( initCalibCoeff[k], 1./calibCoeff[k] );
1634  h1_mcEEParz_[iLoop]->Fill( initCalibCoeff[k]*calibCoeff[k] -1. );
1635  }
1636 
1637  }
1638  }
1639 
1640 
1642  double parResidual[3];
1643  double errparResidual[3];
1644  double zResChi2;
1645  int zResIters;
1646 
1647  ZIterativeAlgorithmWithFit::gausfit(h1_mcParz_[iLoop],parResidual,errparResidual,3.,3., &zResChi2, &zResIters);
1648  //h1_mcParz_[iLoop]->Fit("gaus");
1649 
1650  h2_residualSigma_ -> Fill(loopFlag_ + 1, parResidual[2]);
1651  loopArray[loopFlag_] = loopFlag_ + 1;
1652  sigmaArray[loopFlag_] = parResidual[2];
1653  sigmaErrorArray[loopFlag_] = errparResidual[2];
1654 
1655  std::cout<<"Fit on residuals, sigma is "<<parResidual[2]<<" +/- "<<errparResidual[2]<<std::endl;
1656 
1658  outputFile_->cd();
1659 
1660 
1661  // h2_miscalRecalParz_[iLoop]->Write();
1662  h1_mcParz_[iLoop]->Write();
1663 
1664  //h2_miscalRecalEBParz_[iLoop]->Write();
1665  h1_mcEBParz_[iLoop]->Write();
1666 
1667  //h2_miscalRecalEEParz_[iLoop]->Write();
1668  h1_mcEEParz_[iLoop]->Write();
1669  h2_xtalRecalibCoeffBarrel_[loopFlag_] -> Write();
1670  h2_xtalRecalibCoeffEndcapPlus_[loopFlag_] -> Write();
1671  h2_xtalRecalibCoeffEndcapMinus_[loopFlag_] -> Write();
1672 
1674 
1675  loopFlag_++;
1676 
1677 #ifdef DEBUG
1678  std::cout<<" loopFlag_ is "<<loopFlag_<<std::endl;
1679 #endif
1680 
1681  if ( iLoop == theMaxLoops-1 || iLoop >= theMaxLoops ) return kStop;
1682  else return kContinue;
1683 }
1684 
1685 void ZeeCalibration::bookHistograms()
1686 {
1687 
1688  h1_eventsBeforeEWKSelection_= new TH1F("h1_eventsBeforeEWKSelection", "h1_eventsBeforeEWKSelection", 5,0,5);
1689  h1_eventsAfterEWKSelection_ = new TH1F("h1_eventsAfterEWKSelection", "h1_eventsAfterEWKSelection", 5,0,5);
1690 
1691  h1_eventsBeforeBorderSelection_= new TH1F("h1_eventsBeforeBorderSelection", "h1_eventsBeforeBorderSelection", 5,0,5);
1692  h1_eventsAfterBorderSelection_ = new TH1F("h1_eventsAfterBorderSelection", "h1_eventsAfterBorderSelection", 5,0,5);
1693 
1694  h1_seedOverSC_= new TH1F("h1_seedOverSC", "h1_seedOverSC", 400, 0., 2.);
1695 
1696  myZeePlots_ -> bookHLTHistograms();
1697 
1698  h1_borderElectronClassification_ = new TH1F("h1_borderElectronClassification", "h1_borderElectronClassification", 55, -5 , 50);
1699  h1_preshowerOverSC_= new TH1F("h1_preshowerOverSC", "h1_preshowerOverSC", 400, 0., 1.);
1700 
1701  h2_fEtaBarrelGood_ = new TH2F("fEtaBarrelGood","fEtaBarrelGood",800,-4.,4.,800,0.8,1.2);
1702  h2_fEtaBarrelGood_->SetXTitle("Eta");
1703  h2_fEtaBarrelGood_->SetYTitle("1/fEtaBarrelGood");
1704 
1705  h2_fEtaBarrelBad_ = new TH2F("fEtaBarrelBad","fEtaBarrelBad",800,-4.,4.,800,0.8,1.2);
1706  h2_fEtaBarrelBad_->SetXTitle("Eta");
1707  h2_fEtaBarrelBad_->SetYTitle("1/fEtaBarrelBad");
1708 
1709  h2_fEtaEndcapGood_ = new TH2F("fEtaEndcapGood","fEtaEndcapGood",800,-4.,4.,800,0.8,1.2);
1710  h2_fEtaEndcapGood_->SetXTitle("Eta");
1711  h2_fEtaEndcapGood_->SetYTitle("1/fEtaEndcapGood");
1712 
1713  h2_fEtaEndcapBad_ = new TH2F("fEtaEndcapBad","fEtaEndcapBad",800,-4.,4.,800,0.8,1.2);
1714  h2_fEtaEndcapBad_->SetXTitle("Eta");
1715  h2_fEtaEndcapBad_->SetYTitle("1/fEtaEndcapBad");
1716 
1717  for (int i=0;i<2;i++)
1718  {
1719  char histoName[50];
1720 
1721  sprintf(histoName,"h_eleEffEta_%d",i);
1722  h_eleEffEta_[i] = new TH1F(histoName,histoName, 150, 0., 2.7);
1723  h_eleEffEta_[i]->SetXTitle("|#eta|");
1724 
1725  sprintf(histoName,"h_eleEffPhi_%d",i);
1726  h_eleEffPhi_[i] = new TH1F(histoName,histoName, 400, -4., 4.);
1727  h_eleEffPhi_[i]->SetXTitle("Phi");
1728 
1729  sprintf(histoName,"h_eleEffPt_%d",i);
1730  h_eleEffPt_[i] = new TH1F(histoName,histoName, 200, 0., 200.);
1731  h_eleEffPt_[i]->SetXTitle("p_{T}(GeV/c)");
1732  }
1733 
1734  h2_xtalMiscalibCoeffBarrel_ = new TH2F("h2_xtalMiscalibCoeffBarrel","h2_xtalMiscalibCoeffBarrel", 171, -85, 85, 360, 0, 360);
1735  h2_xtalMiscalibCoeffEndcapMinus_ = new TH2F("h2_xtalMiscalibCoeffEndcapMinus", "h2_xtalMiscalibCoeffEndcapMinus", 100, 0,100, 100, 0, 100);
1736  h2_xtalMiscalibCoeffEndcapPlus_ = new TH2F("h2_xtalMiscalibCoeffEndcapPlus", "h2_xtalMiscalibCoeffEndcapPlus", 100, 0,100, 100, 0, 100);
1737 
1738  h2_xtalMiscalibCoeffBarrel_ ->SetXTitle("ieta");
1739  h2_xtalMiscalibCoeffBarrel_ ->SetYTitle("iphi");
1740 
1741  h2_xtalMiscalibCoeffEndcapMinus_->SetXTitle("ix");
1742  h2_xtalMiscalibCoeffEndcapMinus_->SetYTitle("iy");
1743 
1744  for (int i=0;i<25;i++)
1745  {
1746 
1747  char histoName[50];
1748  sprintf(histoName,"h_ESCEtrueVsEta_%d",i);
1749 
1750  h_ESCEtrueVsEta_[i] = new TH2F(histoName,histoName, 150, 0., 2.7, 300,0.,1.5);
1751  h_ESCEtrueVsEta_[i]->SetXTitle("|#eta|");
1752  h_ESCEtrueVsEta_[i]->SetYTitle("E_{SC,raw}/E_{MC}");
1753 
1754  sprintf(histoName,"h_ESCEtrue_%d",i);
1755 
1756  h_ESCEtrue_[i] = new TH1F(histoName,histoName, 300,0.,1.5);
1757 
1758  sprintf(histoName,"h2_chi2_%d",i);
1759  h2_chi2_[i] = new TH2F(histoName,histoName, 1000,-150,150, 1000, -1, 5);
1760 
1761  sprintf(histoName,"h2_iterations_%d",i);
1762  h2_iterations_[i] = new TH2F(histoName,histoName, 1000,-150,150, 1000, -1, 15);
1763 
1764  sprintf(histoName,"h_ESCcorrEtrueVsEta_%d",i);
1765 
1766  h_ESCcorrEtrueVsEta_[i] = new TH2F(histoName,histoName, 150, 0., 2.7, 300,0.,1.5);
1767  h_ESCcorrEtrueVsEta_[i]->SetXTitle("|#eta|");
1768  h_ESCcorrEtrueVsEta_[i]->SetYTitle("E_{SC,#eta-corr}/E_{MC}");
1769 
1770  sprintf(histoName,"h_ESCcorrEtrue_%d",i);
1771 
1772  h_ESCcorrEtrue_[i] = new TH1F(histoName,histoName, 300,0.,1.5);
1773 
1774  sprintf(histoName,"h2_xtalRecalibCoeffBarrel_%d",i);
1775  h2_xtalRecalibCoeffBarrel_[i] = new TH2F(histoName,histoName, 171, -85, 85, 360, 0, 360);
1776 
1777  h2_xtalRecalibCoeffBarrel_[i]->SetXTitle("ieta");
1778  h2_xtalRecalibCoeffBarrel_[i]->SetYTitle("iphi");
1779 
1780  sprintf(histoName,"h2_xtalRecalibCoeffEndcapMinus_%d",i);
1781  h2_xtalRecalibCoeffEndcapMinus_[i] = new TH2F(histoName,histoName, 100, 0,100, 100, 0, 100);
1782  h2_xtalRecalibCoeffEndcapMinus_[i]->SetXTitle("ix");
1783  h2_xtalRecalibCoeffEndcapMinus_[i]->SetYTitle("iy");
1784 
1785  sprintf(histoName,"h2_xtalRecalibCoeffEndcapPlus_%d",i);
1786  h2_xtalRecalibCoeffEndcapPlus_[i] = new TH2F(histoName,histoName, 100, 0,100, 100, 0, 100);
1787  h2_xtalRecalibCoeffEndcapPlus_[i]->SetXTitle("ix");
1788  h2_xtalRecalibCoeffEndcapPlus_[i]->SetYTitle("iy");
1789 
1790  }
1791 
1792  /*
1793  for (int i=0;i<15;i++)
1794  {
1795 
1796  char histoName[50];
1797 
1798  sprintf(histoName,"h_DiffZMassDistr_%d",i);
1799  h_DiffZMassDistr_[i] = new TH1F(histoName,histoName, 400, -20., 20.);
1800  h_DiffZMassDistr_[i]->SetXTitle("M_{Z, reco} - M_{Z, MC}");
1801  h_DiffZMassDistr_[i]->SetYTitle("events");
1802 
1803  sprintf(histoName,"h_ZMassDistr_%d",i);
1804  h_ZMassDistr_[i] = new TH1F(histoName,histoName, 200, 0., 150.);
1805  h_ZMassDistr_[i]->SetXTitle("RecoZmass (GeV)");
1806  h_ZMassDistr_[i]->SetYTitle("events");
1807 
1808  }
1809  */
1810 
1811  h1_zMassResol_ = new TH1F("zMassResol", "zMassResol", 200, -50., 50.);
1812  h1_zMassResol_->SetXTitle("M_{Z, reco} - M_{Z, MC}");
1813  h1_zMassResol_->SetYTitle("events");
1814 
1815  h1_eleEtaResol_ = new TH1F("eleEtaResol", "eleEtaResol", 100, -0.01, 0.01);
1816  h1_eleEtaResol_->SetXTitle("#eta_{reco} - #eta_{MC}");
1817  h1_eleEtaResol_->SetYTitle("events");
1818 
1819  h1_electronCosTheta_TK_ = new TH1F("electronCosTheta_TK", "electronCosTheta_TK", 100, -1, 1);
1820  h1_electronCosTheta_TK_->SetXTitle("cos #theta_{12}");
1821  h1_electronCosTheta_TK_->SetYTitle("events");
1822 
1823  h1_electronCosTheta_SC_ = new TH1F("electronCosTheta_SC", "electronCosTheta_SC", 100, -1, 1);
1824  h1_electronCosTheta_SC_->SetXTitle("cos #theta_{12}");
1825  h1_electronCosTheta_SC_->SetYTitle("events");
1826 
1827  h1_electronCosTheta_SC_TK_ = new TH1F("electronCosTheta_SC_TK", "electronCosTheta_SC_TK", 200, -0.1, 0.1);
1828  h1_electronCosTheta_SC_TK_->SetXTitle("cos #theta_{12}^{SC}/ cos #theta_{12}^{TK} - 1");
1829  h1_electronCosTheta_SC_TK_->SetYTitle("events");
1830 
1831  h1_elePhiResol_ = new TH1F("elePhiResol", "elePhiResol", 100, -0.01, 0.01);
1832  h1_elePhiResol_->SetXTitle("#phi_{reco} - #phi_{MC}");
1833  h1_elePhiResol_->SetYTitle("events");
1834 
1835 
1836  h1_zEtaResol_ = new TH1F("zEtaResol", "zEtaResol", 200, -1., 1.);
1837  h1_zEtaResol_->SetXTitle("#eta_{Z, reco} - #eta_{Z, MC}");
1838  h1_zEtaResol_->SetYTitle("events");
1839 
1840 
1841  h1_zPhiResol_ = new TH1F("zPhiResol", "zPhiResol", 200, -1., 1.);
1842  h1_zPhiResol_->SetXTitle("#phi_{Z, reco} - #phi_{Z, MC}");
1843  h1_zPhiResol_->SetYTitle("events");
1844 
1845  h1_nEleReco_ = new TH1F("nEleReco","Number of reco electrons",10,-0.5,10.5);
1846  h1_nEleReco_->SetXTitle("nEleReco");
1847  h1_nEleReco_->SetYTitle("events");
1848 
1849  // h1_occupancyVsEta_ = new TH1F("occupancyVsEta","occupancyVsEta",EcalRingCalibrationTools::N_RING_TOTAL,0,(float)EcalRingCalibrationTools::N_RING_TOTAL);
1850 
1851  h1_occupancyVsEta_ = new TH1F("occupancyVsEta","occupancyVsEta",249, -124, 124);
1852  h1_occupancyVsEta_->SetYTitle("Weighted electron statistics");
1853  h1_occupancyVsEta_->SetXTitle("Eta channel");
1854 
1855  h1_weightSumMeanBarrel_= new TH1F("weightSumMeanBarrel","weightSumMeanBarrel",10000, 0, 10000);
1856  h1_weightSumMeanEndcap_= new TH1F("weightSumMeanEndcap","weightSumMeanEndcap",10000, 0, 10000);
1857 
1858  h1_occupancy_ = new TH1F("occupancy","occupancy",1000,0,10000);
1859  h1_occupancy_->SetXTitle("Weighted electron statistics");
1860 
1861  h1_occupancyBarrel_ = new TH1F("occupancyBarrel","occupancyBarrel",1000,0,10000);
1862  h1_occupancyBarrel_->SetXTitle("Weighted electron statistics");
1863 
1864  h1_occupancyEndcap_ = new TH1F("occupancyEndcap","occupancyEndcap",1000,0,10000);
1865  h1_occupancyEndcap_->SetXTitle("Weighted electron statistics");
1866 
1867 
1868  h1_eleClasses_= new TH1F("eleClasses","eleClasses",301,-1,300);
1869  h1_eleClasses_->SetXTitle("classCode");
1870  h1_eleClasses_->SetYTitle("#");
1871 
1872 
1873 
1874  myZeePlots_ ->bookZMCHistograms();
1875 
1876  myZeePlots_ ->bookZHistograms();
1877 
1878  myZeePlots_ ->bookEleMCHistograms();
1879 
1880  myZeePlots_ ->bookEleHistograms();
1881 
1882 
1883  h1_ZCandMult_ =new TH1F("ZCandMult","Multiplicity of Z candidates in one event",10,-0.5,10.5);
1884  h1_ZCandMult_ ->SetXTitle("ZCandMult");
1885 
1886  h1_reco_ZMass_ = new TH1F("reco_ZMass","Inv. mass of 2 reco Electrons",200,0.,150.);
1887  h1_reco_ZMass_->SetXTitle("reco_ZMass (GeV)");
1888  h1_reco_ZMass_->SetYTitle("events");
1889 
1890  h1_reco_ZMassCorr_ = new TH1F("reco_ZMassCorr","Inv. mass of 2 corrected reco Electrons",200,0.,150.);
1891  h1_reco_ZMassCorr_->SetXTitle("reco_ZMass (GeV)");
1892  h1_reco_ZMassCorr_->SetYTitle("events");
1893 
1894  h1_reco_ZMassCorrBB_ = new TH1F("reco_ZMassCorrBB","Inv. mass of 2 corrected reco Electrons",200,0.,150.);
1895  h1_reco_ZMassCorrBB_->SetXTitle("reco_ZMass (GeV)");
1896  h1_reco_ZMassCorrBB_->SetYTitle("events");
1897 
1898  h1_reco_ZMassCorrEE_ = new TH1F("reco_ZMassCorrEE","Inv. mass of 2 corrected reco Electrons",200,0.,150.);
1899  h1_reco_ZMassCorrEE_->SetXTitle("reco_ZMass (GeV)");
1900  h1_reco_ZMassCorrEE_->SetYTitle("events");
1901 
1902  // h2_coeffVsEta_= new TH2F("h2_calibCoeffVsEta","h2_calibCoeffVsEta",EcalRingCalibrationTools::N_RING_TOTAL,0, (double)EcalRingCalibrationTools::N_RING_TOTAL, 200, 0., 2.);
1903 
1904  h2_coeffVsEta_= new TH2F("h2_calibCoeffVsEta","h2_calibCoeffVsEta",249,-124,125, 200, 0., 2.);
1905  h2_coeffVsEta_->SetXTitle("Eta channel");
1906  h2_coeffVsEta_->SetYTitle("recalibCoeff");
1907 
1908  h2_coeffVsEtaGrouped_= new TH2F("h2_calibCoeffVsEtaGrouped","h2_calibCoeffVsEtaGrouped", 200, 0., 3., 200, 0.6, 1.4);
1909  h2_coeffVsEtaGrouped_->SetXTitle("|#eta|");
1910  h2_coeffVsEtaGrouped_->SetYTitle("recalibCoeff");
1911 
1912  h2_zMassVsLoop_= new TH2F("h2_zMassVsLoop","h2_zMassVsLoop",1000,0,40, 90, 80.,95.);
1913 
1914  h2_zMassDiffVsLoop_= new TH2F("h2_zMassDiffVsLoop","h2_zMassDiffVsLoop",1000,0,40, 100, -1., 1.);
1915  h2_zMassDiffVsLoop_->SetXTitle("Iteration");
1916  h2_zMassDiffVsLoop_->SetYTitle("M_{Z, reco peak} - M_{Z, true}");
1917 
1918  h2_zWidthVsLoop_= new TH2F("h2_zWidthVsLoop","h2_zWidthVsLoop",1000,0,40, 100, 0.,10.);
1919 
1920  h2_coeffVsLoop_= new TH2F("h2_coeffVsLoop","h2_coeffVsLoop",1000,0,40, 100, 0., 2.);
1921 
1922  h2_residualSigma_= new TH2F("h2_residualSigma","h2_residualSigma",1000, 0, 40, 100, 0., .5);
1923 
1924  h2_miscalRecal_ = new TH2F("h2_miscalRecal","h2_miscalRecal", 500, 0., 2., 500, 0., 2.);
1925  h2_miscalRecal_->SetXTitle("initCalibCoeff");
1926  h2_miscalRecal_->SetYTitle("1/RecalibCoeff");
1927 
1928  h2_miscalRecalEB_ = new TH2F("h2_miscalRecalEB","h2_miscalRecalEB", 500, 0., 2., 500, 0., 2.);
1929  h2_miscalRecalEB_->SetXTitle("initCalibCoeff");
1930  h2_miscalRecalEB_->SetYTitle("1/RecalibCoeff");
1931 
1932  h2_miscalRecalEE_ = new TH2F("h2_miscalRecalEE","h2_miscalRecalEE", 500, 0., 2., 500, 0., 2.);
1933  h2_miscalRecalEE_->SetXTitle("initCalibCoeff");
1934  h2_miscalRecalEE_->SetYTitle("1/RecalibCoeff");
1935 
1936  h1_mc_ = new TH1F("h1_residualMiscalib","h1_residualMiscalib", 200, -0.2, 0.2);
1937  h1_mcEB_ = new TH1F("h1_residualMiscalibEB","h1_residualMiscalibEB", 200, -0.2, 0.2);
1938  h1_mcEE_ = new TH1F("h1_residualMiscalibEE","h1_residualMiscalibEE", 200, -0.2, 0.2);
1939 
1940  for (int i=0;i<25;i++)
1941  {
1942  char histoName[50];
1943  /*
1944  sprintf(histoName,"h2_miscalRecalParz_%d",i);
1945  h2_miscalRecalParz_[i] = new TH2F(histoName,histoName,500, 0., 2., 500, 0., 2.);
1946  h2_miscalRecalParz_[i]->SetXTitle("initCalibCoeff");
1947  h2_miscalRecalParz_[i]->SetYTitle("1/recalibCoeff");
1948 
1949  sprintf(histoName,"h2_miscalRecalEBParz_%d",i);
1950  h2_miscalRecalEBParz_[i] = new TH2F(histoName,histoName,500, 0., 2., 500, 0., 2.);
1951  h2_miscalRecalEBParz_[i]->SetXTitle("initCalibCoeff");
1952  h2_miscalRecalEBParz_[i]->SetYTitle("1/recalibCoeff");
1953 
1954  sprintf(histoName,"h2_miscalRecalEEParz_%d",i);
1955  h2_miscalRecalEEParz_[i] = new TH2F(histoName,histoName,500, 0., 2., 500, 0., 2.);
1956  h2_miscalRecalEEParz_[i]->SetXTitle("initCalibCoeff");
1957  h2_miscalRecalEEParz_[i]->SetYTitle("1/recalibCoeff");
1958  */
1959 
1960  sprintf(histoName,"h1_residualMiscalibParz_%d",i);
1961  h1_mcParz_[i] = new TH1F(histoName,histoName, 200, -0.2, 0.2);
1962  sprintf(histoName,"h1_residualMiscalibEBParz_%d",i);
1963  h1_mcEBParz_[i] = new TH1F(histoName,histoName, 200, -0.2, 0.2);
1964  sprintf(histoName,"h1_residualMiscalibEEParz_%d",i);
1965  h1_mcEEParz_[i] = new TH1F(histoName,histoName, 200, -0.2, 0.2);
1966 
1967  }
1968 
1969 
1970 }
1971 
1972 
1973 double ZeeCalibration::fEtaBarrelBad(double scEta) const{
1974 
1975  float p0 = 1.00153e+00;
1976  float p1 = 3.29331e-02;
1977  float p2 = 1.21187e-03;
1978 
1979  double x = (double) fabs(scEta);
1980 
1981  return 1. / ( p0 + p1*x*x + p2*x*x*x*x );
1982 
1983 }
1984 
1985 double ZeeCalibration::fEtaEndcapGood(double scEta) const{
1986 
1987  // f(eta) for the first 3 classes (100, 110 and 120)
1988  // Ivica's new corrections 01/06
1989  float p0 = 1.06819e+00;
1990  float p1 = -1.53189e-02;
1991  float p2 = 4.01707e-04 ;
1992 
1993  double x = (double) fabs(scEta);
1994 
1995  return 1. / ( p0 + p1*x*x + p2*x*x*x*x );
1996 
1997 }
1998 
1999 double ZeeCalibration::fEtaEndcapBad(double scEta) const{
2000 
2001  float p0 = 1.17382e+00;
2002  float p1 = -6.52319e-02;
2003  float p2 = 6.26108e-03;
2004 
2005  double x = (double) fabs(scEta);
2006 
2007  return 1. / ( p0 + p1*x*x + p2*x*x*x*x );
2008 
2009 }
2010 
2011 double ZeeCalibration::fEtaBarrelGood(double scEta) const{
2012 
2013  float p0 = 9.99782e-01 ;
2014  float p1 = 1.26983e-02;
2015  float p2 = 2.16344e-03;
2016 
2017  double x = (double) fabs(scEta);
2018 
2019  return 1. / ( p0 + p1*x*x + p2*x*x*x*x );
2020 
2021 }
2022 
2023 
2025 
2026 void ZeeCalibration::fillMCmap(const std::vector<const reco::GsfElectron*>* electronCollection,const std::vector<HepMC::GenParticle*>& mcEle,std::map<HepMC::GenParticle*,const reco::GsfElectron*>& myMCmap)
2027 {
2028  for (unsigned int i=0;i<mcEle.size();i++)
2029  {
2030  float minDR=0.1;
2031  const reco::GsfElectron* myMatchEle=0;
2032  for (unsigned int j=0;j<electronCollection->size();j++)
2033  {
2034  float dr=EvalDR(mcEle[i]->momentum().pseudoRapidity(),(*(*electronCollection)[j]).eta(),mcEle[i]->momentum().phi(),(*(*electronCollection)[j]).phi());
2035  if (dr < minDR )
2036  {
2037  myMatchEle = (*electronCollection)[j];
2038  minDR = dr;
2039  }
2040  }
2041  myMCmap.insert(std::pair<HepMC::GenParticle*,const reco::GsfElectron*>(mcEle[i],myMatchEle));
2042 
2043  }
2044 }
2045 
2046 float ZeeCalibration::EvalDR(float Eta,float Eta_ref,float Phi,float Phi_ref)
2047 {
2048  if (Phi<0) Phi = 2*TMath::Pi() + Phi;
2049  if (Phi_ref<0) Phi_ref = 2*TMath::Pi() + Phi_ref;
2050  float DPhi = Phi - Phi_ref ;
2051  if (fabs(DPhi)>TMath::Pi()) DPhi = 2*TMath::Pi() - fabs(DPhi);
2052 
2053  float DEta = Eta - Eta_ref ;
2054 
2055  float DR = sqrt( DEta*DEta + DPhi*DPhi );
2056  return DR;
2057 }
2058 
2059 float ZeeCalibration::EvalDPhi(float Phi,float Phi_ref)
2060 {
2061  if (Phi<0) Phi = 2*TMath::Pi() + Phi;
2062  if (Phi_ref<0) Phi_ref = 2*TMath::Pi() + Phi_ref;
2063  return (Phi - Phi_ref);
2064 }
2065 
2066 void ZeeCalibration::fillEleInfo( std::vector<HepMC::GenParticle*>& mcEle, std::map<HepMC::GenParticle*,const reco::GsfElectron*>& associationMap)
2067 {
2068 
2069  for (unsigned int i=0;i<mcEle.size();i++)
2070  {
2071 
2072  h_eleEffEta_[0]->Fill(fabs(mcEle[i]->momentum().pseudoRapidity()));
2073  h_eleEffPhi_[0]->Fill(mcEle[i]->momentum().phi());
2074  h_eleEffPt_[0]->Fill(mcEle[i]->momentum().perp());
2075 
2076  std::map<HepMC::GenParticle*,const reco::GsfElectron*>::const_iterator mIter = associationMap.find(mcEle[i]);
2077  if (mIter == associationMap.end() )
2078  continue;
2079 
2080  if((*mIter).second)
2081  {
2082  const reco::GsfElectron* myEle=(*mIter).second;
2083 
2084  h_eleEffEta_[1]->Fill(fabs(mcEle[i]->momentum().pseudoRapidity()));
2085  h_eleEffPhi_[1]->Fill(mcEle[i]->momentum().phi());
2086  h_eleEffPt_[1]->Fill(mcEle[i]->momentum().perp());
2087  h1_eleEtaResol_->Fill( myEle->eta() - mcEle[i]->momentum().eta() );
2088  h1_elePhiResol_->Fill( myEle->phi() - mcEle[i]->momentum().phi() );
2089 
2090  const reco::SuperCluster* mySC=&(*(myEle->superCluster()));
2091  if (/*fabs(mySC->position().eta()) < 2.4*/1)
2092  {
2093  // if(myEle->classification()>=100)std::cout<<"mySC->preshowerEnergy()"<<mySC->preshowerEnergy()<<std::endl;
2094 
2095  h_ESCEtrue_[loopFlag_]->Fill(mySC->energy()/mcEle[i]->momentum().e());
2096  h_ESCEtrueVsEta_[loopFlag_]->Fill(fabs(mySC->position().eta()),mySC->energy()/mcEle[i]->momentum().e());
2097 
2098  double corrSCenergy = ( mySC->energy() )/getEtaCorrection(myEle);
2099  h_ESCcorrEtrue_[loopFlag_]->Fill(corrSCenergy/mcEle[i]->momentum().e());
2100  h_ESCcorrEtrueVsEta_[loopFlag_]->Fill(fabs(mySC->position().eta()),corrSCenergy/mcEle[i]->momentum().e());
2101 
2102 // std::vector<DetId> mySCRecHits = mySC->seed()->getHitsByDetId();
2103 
2104  h1_seedOverSC_->Fill( mySC->seed()->energy() / mySC->energy() );
2105  h1_preshowerOverSC_->Fill( mySC->preshowerEnergy() / mySC->energy() );
2106 
2107  }
2108 
2109  }
2110  }
2111 }
2112 
2113 int ZeeCalibration::ringNumberCorrector(int k)
2114 {
2115 
2116  int index=-999;
2117 
2118  if( calibMode_ == "RING"){
2119  if(k>=0 && k<=84)index = k - 85;
2120 
2121  if(k>=85 && k<=169)index = k - 84;
2122 
2123  if(k>=170 && k<=208)index = - k + 84;
2124 
2125  if(k>=209 && k<=247)index = k - 123;
2126 
2127  }
2128 
2129  else if( calibMode_ == "MODULE"){
2130 
2131  if(k>=0 && k<=71)index = k - 72;
2132 
2133  if(k>=72 && k<=143)index = k - 71;
2134 
2135  }
2136  return index;
2137 
2138 }
2139 
2140 
2141 double ZeeCalibration::getEtaCorrection(const reco::GsfElectron* ele){
2142 
2143  double correction(1.);
2144 
2145  if(ele->classification() ==0 ||
2146  ele->classification() ==10 ||
2147  ele->classification() ==20)
2148  correction = fEtaBarrelGood(ele->superCluster()->eta());
2149 
2150  if(ele->classification() ==100 ||
2151  ele->classification() ==110 ||
2152  ele->classification() ==120)
2153  correction = fEtaEndcapGood(ele->superCluster()->eta());
2154 
2155  if(ele->classification() ==30 ||
2156  ele->classification() ==31 ||
2157  ele->classification() ==32 ||
2158  ele->classification() ==33 ||
2159  ele->classification() ==34)
2160  correction = fEtaBarrelBad(ele->superCluster()->eta());
2161 
2162 
2163  if(ele->classification() ==130 ||
2164  ele->classification() ==131 ||
2165  ele->classification() ==132 ||
2166  ele->classification() ==133 ||
2167  ele->classification() ==134)
2168  correction = fEtaEndcapBad(ele->superCluster()->eta());
2169 
2170  return correction;
2171 }
2172 
2173 std::pair<DetId, double> ZeeCalibration::getHottestDetId(const std::vector<std::pair< DetId,float > >& mySCRecHits, const EBRecHitCollection* ebhits, const EERecHitCollection* eehits){
2174 
2175 
2176  double maxEnergy = -9999.;
2177  const EcalRecHit* hottestRecHit=0;
2178 
2179  std::pair<DetId, double> myPair (DetId(0), -9999.);
2180 
2181 
2182  for( std::vector<std::pair<DetId,float> >::const_iterator idIt=mySCRecHits.begin(); idIt != mySCRecHits.end(); idIt++){
2183 
2184  if (idIt->first.subdetId() == EcalBarrel )
2185  {
2186  hottestRecHit = & (* ( ebhits->find((*idIt).first) ) );
2187 
2188  if( hottestRecHit == & (*( ebhits->end())) )
2189  {
2190  std::cout<<"@@@@@@@@@@@@@@@@@@@@@@@@@@@ NO RECHIT FOUND SHOULD NEVER HAPPEN"<<std::endl;
2191  continue;
2192  }
2193  }
2194  else if (idIt->first.subdetId() == EcalEndcap )
2195  {
2196  hottestRecHit = & (* ( eehits->find((*idIt).first) ) );
2197  if( hottestRecHit == & (*( eehits->end())) )
2198  {
2199  std::cout<<"@@@@@@@@@@@@@@@@@@@@@@@@@@@ NO RECHIT FOUND SHOULD NEVER HAPPEN"<<std::endl;
2200  continue;
2201  }
2202  }
2203 
2204  //std::cout<<"[getHottestDetId] hottestRecHit->energy() "<<hottestRecHit->energy()<<std::endl;
2205 
2206  if(hottestRecHit && hottestRecHit->energy() > maxEnergy){
2207 
2208  maxEnergy = hottestRecHit->energy();
2209 
2210  myPair.first = hottestRecHit ->id();
2211  myPair.second = maxEnergy;
2212 
2213  }
2214 
2215  }//end loop to find hottest RecHit
2216 
2217  //std::cout<<"[ZeeCalibration::getHottestDetId] going to return..."<<std::endl;
2218 
2219  return myPair;
2220 
2221 }
2222 
2223 
2224 bool ZeeCalibration::xtalIsOnModuleBorder( EBDetId myEBDetId ){
2225 
2226  bool myBool(false);
2227 
2228  short ieta = myEBDetId.ieta();
2229  short iphi = myEBDetId.iphi();
2230 
2231  // std::cout<<"[xtalIsOnModuleBorder] ieta: "<<ieta<<" iphi "<<iphi<<std::endl;
2232 
2233  myBool = ( abs( ieta ) == 1 || abs( ieta ) == 25
2234  || abs( ieta ) ==26 || abs( ieta ) == 45
2235  || abs( ieta ) ==46 || abs( ieta ) == 65
2236  || abs( ieta ) ==66 || abs( ieta ) == 85 );
2237 
2238  for(int i = 0; i < 19; i++){
2239 
2240  if(iphi == ( 20*i + 1 ) || iphi == 20*i )
2241  myBool = true;
2242 
2243  }
2244 
2245  return myBool;
2246 }
2247 
2248 float ZeeCalibration::computeCoefficientDistanceAtIteration( float v1[250], float v2[250], int size ){
2249 
2250  float dist(0.);
2251 
2252  for(int i =0; i < size; i++){
2253 
2254  // std::cout<< "[ZeeCalibration::computeCoefficientDistanceAtIteration] Adding term "<<pow( v1[i]-v2[i], 2 )<<" from v1 "<<v1[i]<<" and v2 "<<v2[i]<<std::endl;
2255 
2256  bool isNearCrack = false;
2257 
2258  if( calibMode_ == "RING"){//exclude non-calibrated rings from computation
2259 
2260  isNearCrack = ( abs( ringNumberCorrector(i) ) == 1 || abs( ringNumberCorrector(i) ) == 25 ||
2261  abs( ringNumberCorrector(i) ) == 26 || abs( ringNumberCorrector(i) ) == 45 ||
2262  abs( ringNumberCorrector(i) ) == 46 || abs( ringNumberCorrector(i) ) == 65 ||
2263  abs( ringNumberCorrector(i) ) == 66 || abs( ringNumberCorrector(i) ) == 85 ||
2264  abs( ringNumberCorrector(i) ) == 86 || abs( ringNumberCorrector(i) ) == 124 );
2265  }
2266 
2267  if(!isNearCrack)
2268  dist += pow( v1[i]-v2[i], 2 );
2269  }
2270 
2271  dist = sqrt(dist) / size;
2272 
2273  return dist;
2274 
2275 }
2276 
2277 
2278 void ZeeCalibration::resetVariables(){
2279 
2280  BBZN=0;
2281  EBZN=0;
2282  EEZN=0;
2283  BBZN_gg=0;
2284  EBZN_gg=0;
2285  EEZN_gg=0;
2286 
2287  BBZN_tt=0;
2288  EBZN_tt=0;
2289  EEZN_tt=0;
2290 
2291  BBZN_t0=0;
2292  EBZN_t0=0;
2293  EEZN_t0=0;
2294 
2295  TOTAL_ELECTRONS_IN_BARREL=0;
2296  TOTAL_ELECTRONS_IN_ENDCAP=0;
2297 
2298  GOLDEN_ELECTRONS_IN_BARREL=0;
2299  GOLDEN_ELECTRONS_IN_ENDCAP=0;
2300  SILVER_ELECTRONS_IN_BARREL=0;
2301  SILVER_ELECTRONS_IN_ENDCAP=0;
2302  SHOWER_ELECTRONS_IN_BARREL=0;
2303  SHOWER_ELECTRONS_IN_ENDCAP=0;
2304  CRACK_ELECTRONS_IN_BARREL=0;
2305  CRACK_ELECTRONS_IN_ENDCAP=0;
2306 
2307 
2308  BARREL_ELECTRONS_BEFORE_BORDER_CUT = 0;
2309  BARREL_ELECTRONS_AFTER_BORDER_CUT = 0;
2310 
2311  return;
2312 
2313 }
2314 
2315 
2316 void ZeeCalibration::resetHistograms(){
2317 
2318  h1_eventsBeforeEWKSelection_ ->Reset();
2319  h1_eventsAfterEWKSelection_ ->Reset();
2320 
2321  h1_eventsBeforeBorderSelection_ ->Reset();
2322  h1_eventsAfterBorderSelection_ ->Reset();
2323 
2324  for (int i=0;i<2;i++)
2325  {
2326  h_eleEffEta_[i] ->Reset();
2327  h_eleEffPhi_[i] ->Reset();
2328  h_eleEffPt_[i] ->Reset();
2329  }
2330 
2331  h1_seedOverSC_ ->Reset();
2332  h1_preshowerOverSC_ ->Reset();
2333 
2334  h1_eleEtaResol_->Reset();
2335  h1_elePhiResol_->Reset();
2336 
2337  h1_zMassResol_->Reset();
2338 
2339  h1_electronCosTheta_TK_->Reset();
2340  h1_electronCosTheta_SC_->Reset();
2341  h1_electronCosTheta_SC_TK_->Reset();
2342 
2343  h2_fEtaBarrelGood_->Reset();
2344  h2_fEtaBarrelBad_->Reset();
2345  h2_fEtaEndcapGood_->Reset();
2346  h2_fEtaEndcapBad_->Reset();
2347  h1_eleClasses_->Reset();
2348 
2349  h1_ZCandMult_-> Reset();
2350  h1_reco_ZMass_-> Reset();
2351  h1_reco_ZMassCorr_-> Reset();
2352  h1_reco_ZMassCorrBB_-> Reset();
2353  h1_reco_ZMassCorrEE_-> Reset();
2354  h1_occupancyVsEta_-> Reset();
2355  h1_occupancy_-> Reset();
2356  h1_occupancyBarrel_-> Reset();
2357  h1_occupancyEndcap_-> Reset();
2358 
2359  return;
2360 
2361 }
2362 
2363 
2364 void ZeeCalibration::printStatistics(){
2365 
2366 
2367  std::cout<< "[ CHECK ON BARREL ELECTRON NUMBER ]"<<" first "<<BARREL_ELECTRONS_BEFORE_BORDER_CUT<<" second "<<TOTAL_ELECTRONS_IN_BARREL << std::endl;
2368 
2369 
2370  std::cout<< "[ EFFICIENCY OF THE BORDER SELECTION ]" << (float)BARREL_ELECTRONS_AFTER_BORDER_CUT / (float) BARREL_ELECTRONS_BEFORE_BORDER_CUT << std::endl;
2371 
2372  std::cout<< "[ EFFICIENCY OF THE GOLDEN SELECTION ] BARREL: " << (float)GOLDEN_ELECTRONS_IN_BARREL / (float) TOTAL_ELECTRONS_IN_BARREL << " ENDCAP: "<< (float)GOLDEN_ELECTRONS_IN_ENDCAP / (float) TOTAL_ELECTRONS_IN_ENDCAP << std::endl;
2373 
2374  std::cout<< "[ EFFICIENCY OF THE SILVER SELECTION ] BARREL: " << (float)SILVER_ELECTRONS_IN_BARREL / (float) TOTAL_ELECTRONS_IN_BARREL << " ENDCAP: "<< (float)SILVER_ELECTRONS_IN_ENDCAP / (float) TOTAL_ELECTRONS_IN_ENDCAP << std::endl;
2375 
2376  std::cout<< "[ EFFICIENCY OF THE SHOWER SELECTION ] BARREL: " << (float)SHOWER_ELECTRONS_IN_BARREL / (float) TOTAL_ELECTRONS_IN_BARREL << " ENDCAP: "<< (float)SHOWER_ELECTRONS_IN_ENDCAP / (float) TOTAL_ELECTRONS_IN_ENDCAP << std::endl;
2377 
2378  std::cout<< "[ EFFICIENCY OF THE CRACK SELECTION ] BARREL: " << (float)CRACK_ELECTRONS_IN_BARREL / (float) TOTAL_ELECTRONS_IN_BARREL << " ENDCAP: "<< (float)CRACK_ELECTRONS_IN_ENDCAP / (float) TOTAL_ELECTRONS_IN_ENDCAP << std::endl;
2379 
2380 
2381 
2382  std::ofstream fout("ZeeStatistics.txt");
2383 
2384  if(!fout) {
2385  std::cout << "Cannot open output file.\n";
2386  }
2387 
2388  fout<<"ZeeStatistics"<<std::endl;
2389 
2390  fout<<"##########################RECO#########################"<<std::endl;
2391  fout<<"##################Zee with Barrel-Barrel electrons: "<<BBZN<<std::endl;
2392 
2393  fout<<"Golden-Golden fraction: "<<(float)BBZN_gg/BBZN<<" 3-3 fraction is "<<(float)BBZN_tt/BBZN<<" 3-whatever fraction is "<<(float)BBZN_t0/BBZN<<std::endl;
2394  fout<<"##################Zee with Barrel-Endcap electrons: "<<EBZN<<std::endl;
2395  fout<<"Golden-Golden fraction: "<<(float)EBZN_gg/EBZN<<" 3-3 fraction is "<<(float)EBZN_tt/EBZN<<" 3-whatever fraction is "<<(float)EBZN_t0/EBZN<<std::endl;
2396  fout<<"##################Zee with Endcap-Endcap electrons: "<<EEZN<<std::endl;
2397  fout<<"Golden-Golden fraction: "<<(float)EEZN_gg/EEZN<<" 3-3 fraction is "<<(float)EEZN_tt/EEZN<<" 3-whatever fraction is "<<(float)EEZN_t0/EEZN<<std::endl;
2398 
2399  fout<<"\n"<<std::endl;
2400 
2401  fout<<"##########################GEN#########################"<<std::endl;
2402  fout<<"##################Zee with Barrel-Barrel electrons: "<<(float)MCZBB/NEVT<<std::endl;
2403  fout<<"##################Zee with Barrel-Endcap electrons: "<<(float)MCZEB/NEVT<<std::endl;
2404  fout<<"##################Zee with Endcap-Endcap electrons: "<<(float)MCZEE/NEVT<<std::endl;
2405 
2406  fout.close();
2407 
2408 
2409 
2410 }
ZeeCalibration::etBins_
unsigned int etBins_
Definition: ZeeCalibration.h:169
ZeeCalibration::endOfJob
virtual void endOfJob()
Called at end of job.
Definition: ZeeCalibration.cc:197
Pi
const double Pi
Definition: CosmicMuonParameters.h:18
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
reco::CaloCluster::position
const math::XYZPoint & position() const
cluster centroid position
Definition: CaloCluster.h:124
i
int i
Definition: DBlmapReader.cc:9
ZeeCalibration::h2_fEtaBarrelBad_
TH2F * h2_fEtaBarrelBad_
Definition: ZeeCalibration.h:216
ZeeCalibration::h1_mcEEParz_
TH1F * h1_mcEEParz_[25]
Definition: ZeeCalibration.h:279
ZIterativeAlgorithmWithFit::getHistos
const ZIterativeAlgorithmWithFitPlots * getHistos() const
Definition: ZIterativeAlgorithmWithFit.h:60
ZeeCalibration::calibMode_
std::string calibMode_
Definition: ZeeCalibration.h:161
ZeeCalibration::outputFileName_
std::string outputFileName_
Definition: ZeeCalibration.h:148
EBDetId::MIN_IPHI
static const int MIN_IPHI
Definition: EBDetId.h:142
ZeeCalibration::fEtaEndcapGood
double fEtaEndcapGood(double scEta) const
Definition: ZeeCalibration.cc:1985
create_public_lumi_plots.exp
tuple exp
Definition: create_public_lumi_plots.py:1087
EcalIndexingTools::setBinRange
void setBinRange(int, double, double, int, double, double)
Definition: EcalIndexingTools.cc:9
EEDetId::ix
int ix() const
Definition: EEDetId.h:76
cuy.ib
int ib
Definition: cuy.py:660
ZeeCalibration::scIslandCollection_
std::string scIslandCollection_
Definition: ZeeCalibration.h:158
ZeeCalibration::h_eleEffPhi_
TH1F * h_eleEffPhi_[2]
Definition: ZeeCalibration.h:248
ZeeCalibration::electronCollection_
std::string electronCollection_
Definition: ZeeCalibration.h:164
ZeeCalibration::BARREL_ELECTRONS_AFTER_BORDER_CUT
int BARREL_ELECTRONS_AFTER_BORDER_CUT
Definition: ZeeCalibration.h:328
ZeeCalibration::h2_miscalRecal_
TH2F * h2_miscalRecal_
Definition: ZeeCalibration.h:263
ZeeCalibration::h1_eventsBeforeEWKSelection_
TH1F * h1_eventsBeforeEWKSelection_
Definition: ZeeCalibration.h:208
ZeeCalibration::fEtaBarrelBad
double fEtaBarrelBad(double scEta) const
Definition: ZeeCalibration.cc:1973
ZeeCalibration::computeCoefficientDistanceAtIteration
float computeCoefficientDistanceAtIteration(float v1[250], float v2[250], int size)
Definition: ZeeCalibration.cc:2248
ZeeCalibration::endOfLoop
virtual Status endOfLoop(const edm::EventSetup &, unsigned int iLoop)
Called at end of loop.
Definition: ZeeCalibration.cc:1508
ZeeCalibration::beginOfJob
virtual void beginOfJob()
Called at beginning of job.
Definition: ZeeCalibration.cc:190
ZeeKinematicTools::calculateZMassWithCorrectedElectrons_withTK
static float calculateZMassWithCorrectedElectrons_withTK(const std::pair< calib::CalibElectron *, calib::CalibElectron * > &aZCandidate, float ele1EnergyCorrection, float ele2EnergyCorrection)
Definition: ZeeKinematicTools.cc:98
ZeeCalibration::calibCoeff
float calibCoeff[250]
Definition: ZeeCalibration.h:190
EcalCondObjectContainer::getMap
const self & getMap() const
Definition: EcalCondObjectContainer.h:100
ZeeKinematicTools::cosThetaElectrons_TK
static float cosThetaElectrons_TK(const std::pair< calib::CalibElectron *, calib::CalibElectron * > &aZCandidate, float ele1EnergyCorrection, float ele2EnergyCorrection)
Definition: ZeeKinematicTools.cc:67
ZeePlots::fillEleInfo
void fillEleInfo(const reco::GsfElectronCollection *)
Definition: ZeePlots.cc:261
ZeePlots::fillEleMCInfo
void fillEleMCInfo(const HepMC::GenEvent *)
Definition: ZeePlots.cc:210
ZIterativeAlgorithmWithFit::getOptimizedChiSquare
const std::vector< float > & getOptimizedChiSquare() const
Definition: ZIterativeAlgorithmWithFit.h:70
ZeeCalibration::theMaxLoops
unsigned int theMaxLoops
Definition: ZeeCalibration.h:316
ZeeCalibration::GOLDEN_ELECTRONS_IN_ENDCAP
int GOLDEN_ELECTRONS_IN_ENDCAP
Definition: ZeeCalibration.h:334
ZeePlots::fillZMCInfo
void fillZMCInfo(const HepMC::GenEvent *)
Definition: ZeePlots.cc:158
reco::LeafCandidate::phi
virtual float phi() const
momentum azimuthal angle
Definition: LeafCandidate.h:180
ZeeCalibration::h_ESCEtrue_
TH1F * h_ESCEtrue_[25]
Definition: ZeeCalibration.h:250
ZeeCalibration::h1_eventsAfterBorderSelection_
TH1F * h1_eventsAfterBorderSelection_
Definition: ZeeCalibration.h:212
ZeeCalibration::hlTriggerResults_
edm::InputTag hlTriggerResults_
Definition: ZeeCalibration.h:346
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
MiscalibReaderFromXML::parseXMLMiscalibFile
bool parseXMLMiscalibFile(std::string configFile)
Definition: MiscalibReaderXML.cc:80
ZeeCalibration::fillMCmap
void fillMCmap(const std::vector< const reco::GsfElectron * > *electronCollection, const std::vector< HepMC::GenParticle * > &mcEle, std::map< HepMC::GenParticle *, const reco::GsfElectron * > &myMCmap)
Definition: ZeeCalibration.cc:2026
ZeeKinematicTools::calculateZMass_withTK
static float calculateZMass_withTK(const std::pair< calib::CalibElectron *, calib::CalibElectron * > &aZCandidate)
Definition: ZeeKinematicTools.cc:115
ZIterativeAlgorithmWithFit::getOptimizedCoefficients
const std::vector< float > & getOptimizedCoefficients() const
Definition: ZIterativeAlgorithmWithFit.h:66
ZeeCalibration::h1_occupancyEndcap_
TH1F * h1_occupancyEndcap_
Definition: ZeeCalibration.h:302
MZ
#define MZ
Definition: ZeeCalibration.cc:74
ZIterativeAlgorithmWithFit::addEvent
bool addEvent(calib::CalibElectron *, calib::CalibElectron *, float)
Definition: ZIterativeAlgorithmWithFit.cc:253
ZeePlots::writeMCEleHistograms
void writeMCEleHistograms()
Definition: ZeePlots.cc:300
ZeeCalibration::duringLoop
virtual Status duringLoop(const edm::Event &, const edm::EventSetup &)
Called at each event.
Definition: ZeeCalibration.cc:656
ZeeCalibration::h1_eventsBeforeBorderSelection_
TH1F * h1_eventsBeforeBorderSelection_
Definition: ZeeCalibration.h:211
ZeeCalibration::h1_electronCosTheta_SC_
TH1F * h1_electronCosTheta_SC_
Definition: ZeeCalibration.h:305
ZeePlots::writeMCZHistograms
void writeMCZHistograms()
Definition: ZeePlots.cc:146
eta
T eta() const
Definition: Basic3DVectorLD.h:177
ZeeCalibration::fEtaEndcapBad
double fEtaEndcapBad(double scEta) const
Definition: ZeeCalibration.cc:1999
ZeeKinematicTools::cosThetaElectrons_SC
static float cosThetaElectrons_SC(const std::pair< calib::CalibElectron *, calib::CalibElectron * > &aZCandidate, float ele1EnergyCorrection, float ele2EnergyCorrection)
Definition: ZeeKinematicTools.cc:42
ZeeCalibration::h1_reco_ZMassCorr_
TH1F * h1_reco_ZMassCorr_
Definition: ZeeCalibration.h:234
ZIterativeAlgorithmWithFit::gausfit
static void gausfit(TH1F *histoou, double *par, double *errpar, float nsigmalow, float nsigmaup, double *mychi2, int *iterations)
Definition: ZIterativeAlgorithmWithFit.cc:367
ZeeCalibration::GOLDEN_ELECTRONS_IN_BARREL
int GOLDEN_ELECTRONS_IN_BARREL
Definition: ZeeCalibration.h:333
ZeeCalibration::h2_coeffVsEtaGrouped_
TH2F * h2_coeffVsEtaGrouped_
Definition: ZeeCalibration.h:257
ZeeCalibration::h1_occupancyVsEta_
TH1F * h1_occupancyVsEta_
Definition: ZeeCalibration.h:295
ZeeCalibration::endcapfile_
std::string endcapfile_
Definition: ZeeCalibration.h:177
DeDxDiscriminatorTools::charge
double charge(const std::vector< uint8_t > &Ampls)
Definition: DeDxDiscriminatorTools.cc:43
reco::GsfElectronCollection
std::vector< GsfElectron > GsfElectronCollection
collection of GsfElectron objects
Definition: GsfElectronFwd.h:14
EcalRingCalibrationTools::getDetIdsInRing
static std::vector< DetId > getDetIdsInRing(short aRingIndex)
Retrieve the DetIds in a phi-ring.
Definition: EcalRingCalibrationTools.cc:71
EEDetId::IX_MIN
static const int IX_MIN
Definition: EEDetId.h:294
ZeeCalibration::rechitCollection_
std::string rechitCollection_
Definition: ZeeCalibration.h:151
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:255
edm::ESProducer::setWhatProduced
void setWhatProduced(T *iThis, const es::Label &iLabel=es::Label())
Definition: ESProducer.h:115
ZeeCalibration::h1_occupancyBarrel_
TH1F * h1_occupancyBarrel_
Definition: ZeeCalibration.h:301
EEDetId::IY_MIN
static const int IY_MIN
Definition: EEDetId.h:298
ZeeCalibration::h2_miscalRecalEE_
TH2F * h2_miscalRecalEE_
Definition: ZeeCalibration.h:276
EBDetId::validDetId
static bool validDetId(int i, int j)
check if a valid index combination
Definition: EBDetId.h:124
ZeeCalibration::EvalDR
float EvalDR(float Eta, float Eta_ref, float Phi, float Phi_ref)
Definition: ZeeCalibration.cc:2046
EBDetId::iphi
int iphi() const
get the crystal iphi
Definition: EBDetId.h:53
DetId::rawId
uint32_t rawId() const
get the raw id
Definition: DetId.h:43
ZeePlots::bookEleHistograms
void bookEleHistograms()
Definition: ZeePlots.cc:231
edm::second
U second(std::pair< T, U > const &p)
Definition: ParameterSet.cc:250
ZeeCalibration::h1_electronCosTheta_TK_
TH1F * h1_electronCosTheta_TK_
Definition: ZeeCalibration.h:304
ZeeCalibration::~ZeeCalibration
~ZeeCalibration()
Destructor.
Definition: ZeeCalibration.cc:175
ZeeCalibration::h_ESCcorrEtrueVsEta_
TH2F * h_ESCcorrEtrueVsEta_[25]
Definition: ZeeCalibration.h:254
ZeeCalibration::ical
boost::shared_ptr< EcalIntercalibConstants > ical
Definition: ZeeCalibration.h:195
iEvent
int iEvent
Definition: GenABIO.cc:230
ZeeCalibration::h2_fEtaEndcapBad_
TH2F * h2_fEtaEndcapBad_
Definition: ZeeCalibration.h:218
ZeeCalibration::rechitProducer_
std::string rechitProducer_
Definition: ZeeCalibration.h:150
CaloMiscalibMapEcal::get
const EcalIntercalibConstants & get()
Definition: CaloMiscalibMapEcal.h:111
ZeeCalibration::ringNumberCorrector
int ringNumberCorrector(int k)
Definition: ZeeCalibration.cc:2113
ZeeCalibration::scIslandProducer_
std::string scIslandProducer_
Definition: ZeeCalibration.h:157
ZeeCalibration::SILVER_ELECTRONS_IN_BARREL
int SILVER_ELECTRONS_IN_BARREL
Definition: ZeeCalibration.h:336
reco::SuperClusterCollection
std::vector< SuperCluster > SuperClusterCollection
collection of SuperCluser objectr
Definition: SuperClusterFwd.h:9
HcalObjRepresent::Fill
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
Definition: HcalObjRepresent.h:611
ZeeCalibration::sigmaArray
double sigmaArray[50]
Definition: ZeeCalibration.h:323
ZeeCalibration::TOTAL_ELECTRONS_IN_ENDCAP
int TOTAL_ELECTRONS_IN_ENDCAP
Definition: ZeeCalibration.h:331
reco::LeafCandidate::eta
virtual float eta() const
momentum pseudorapidity
Definition: LeafCandidate.h:184
ZeeCalibration::BARREL_ELECTRONS_BEFORE_BORDER_CUT
int BARREL_ELECTRONS_BEFORE_BORDER_CUT
Definition: ZeeCalibration.h:327
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:48
ZeePlots::fillHLTInfo
void fillHLTInfo(edm::Handle< edm::TriggerResults >)
Definition: ZeePlots.cc:341
ZeeCalibration::SHOWER_ELECTRONS_IN_ENDCAP
int SHOWER_ELECTRONS_IN_ENDCAP
Definition: ZeeCalibration.h:340
reco::GsfElectron::superCluster
virtual SuperClusterRef superCluster() const
reference to a SuperCluster
Definition: GsfElectron.h:182
EcalIndexingTools::getInstance
static EcalIndexingTools * getInstance()
Definition: EcalIndexingTools.h:27
ZeeCalibration::scProducer_
std::string scProducer_
Definition: ZeeCalibration.h:154
ZeePlots::writeZHistograms
void writeZHistograms()
Definition: ZeePlots.cc:132
EEDetId::zside
int zside() const
Definition: EEDetId.h:70
ZeeCalibration::erechitProducer_
std::string erechitProducer_
Definition: ZeeCalibration.h:152
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
j
int j
Definition: DBlmapReader.cc:9
ZeeCalibration::h1_mcEBParz_
TH1F * h1_mcEBParz_[25]
Definition: ZeeCalibration.h:275
EcalRecHit::energy
float energy() const
Definition: EcalRecHit.h:68
reco::CaloCluster::energy
double energy() const
cluster energy
Definition: CaloCluster.h:120
ZeeCalibration::theAlgorithm_
ZIterativeAlgorithmWithFit * theAlgorithm_
Definition: ZeeCalibration.h:197
EEDetId::iy
int iy() const
Definition: EEDetId.h:82
EBDetId::ieta
int ieta() const
get the crystal ieta
Definition: EBDetId.h:51
ZeeCalibration::h1_elePhiResol_
TH1F * h1_elePhiResol_
Definition: ZeeCalibration.h:245
ZeeCalibration::electronSelection_
unsigned int electronSelection_
Definition: ZeeCalibration.h:320
first
bool first
Definition: L1TdeRCT.cc:75
ZeeCalibration::xtalIsOnModuleBorder
bool xtalIsOnModuleBorder(EBDetId myEBDetId)
Definition: ZeeCalibration.cc:2224
ZeeCalibration::h1_eleEtaResol_
TH1F * h1_eleEtaResol_
Definition: ZeeCalibration.h:244
ZeeCalibration::h2_xtalRecalibCoeffBarrel_
TH2F * h2_xtalRecalibCoeffBarrel_[25]
Definition: ZeeCalibration.h:284
p2
double p2[4]
Definition: TauolaWrapper.h:90
ZeeCalibration::EvalDPhi
float EvalDPhi(float Phi, float Phi_ref)
Definition: ZeeCalibration.cc:2059
ZeeCalibration::h1_eventsAfterEWKSelection_
TH1F * h1_eventsAfterEWKSelection_
Definition: ZeeCalibration.h:209
edm::Event::getByLabel
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:390
EcalRingCalibrationTools::getDetIdsInModule
static std::vector< DetId > getDetIdsInModule(short int)
Definition: EcalRingCalibrationTools.cc:139
ZeePlots::bookZMCHistograms
void bookZMCHistograms()
Definition: ZeePlots.cc:62
ZeeCalibration::h1_preshowerOverSC_
TH1F * h1_preshowerOverSC_
Definition: ZeeCalibration.h:227
ZeeCalibration::h1_electronCosTheta_SC_TK_
TH1F * h1_electronCosTheta_SC_TK_
Definition: ZeeCalibration.h:306
ZeeCalibration::coefficientDistanceAtIteration
double coefficientDistanceAtIteration[50]
Definition: ZeeCalibration.h:325
ZeeCalibration::h_eleEffPt_
TH1F * h_eleEffPt_[2]
Definition: ZeeCalibration.h:249
gen::k
int k[5][pyjets_maxn]
Definition: Cascade2Hadronizer.cc:79
EEDetId::IX_MAX
static const int IX_MAX
Definition: EEDetId.h:302
edm::SortedCollection::end
const_iterator end() const
Definition: SortedCollection.h:310
ZeeCalibration::h2_zMassDiffVsLoop_
TH2F * h2_zMassDiffVsLoop_
Definition: ZeeCalibration.h:259
ZeeCalibration::scCollection_
std::string scCollection_
Definition: ZeeCalibration.h:155
graph
Definition: adjgraph.h:12
DetId
Definition: DetId.h:18
ZeeCalibration::h1_mcParz_
TH1F * h1_mcParz_[25]
Definition: ZeeCalibration.h:266
ZIterativeAlgorithmWithFit::getOptimizedCoefficientsError
const std::vector< float > & getOptimizedCoefficientsError() const
Definition: ZIterativeAlgorithmWithFit.h:68
EBDetId::MAX_IPHI
static const int MAX_IPHI
Definition: EBDetId.h:144
EcalRecHit::id
DetId id() const
get the id
Definition: EcalRecHit.h:76
ZeeCalibration::myZeePlots_
ZeePlots * myZeePlots_
Definition: ZeeCalibration.h:199
ZeePlots::writeEleHistograms
void writeEleHistograms()
Definition: ZeePlots.cc:283
ZeeCalibration::etaBins_
unsigned int etaBins_
Definition: ZeeCalibration.h:168
ZeeCalibration::theParameterSet
edm::ParameterSet theParameterSet
Definition: ZeeCalibration.h:204
EEDetId::validDetId
static bool validDetId(int crystal_ix, int crystal_iy, int iz)
Definition: EEDetId.h:248
ZeeCalibration::h2_miscalRecalEB_
TH2F * h2_miscalRecalEB_
Definition: ZeeCalibration.h:272
reco::GsfElectron::classification
Classification classification() const
Definition: GsfElectron.h:680
ZeeCalibration::h_ESCcorrEtrue_
TH1F * h_ESCcorrEtrue_[25]
Definition: ZeeCalibration.h:253
ZeeCalibration::getHottestDetId
std::pair< DetId, double > getHottestDetId(const std::vector< std::pair< DetId, float > > &mySCRecHits, const EBRecHitCollection *ebhits, const EERecHitCollection *eehits)
Definition: ZeeCalibration.cc:2173
ZeeCalibration::loopArray
double loopArray[50]
Definition: ZeeCalibration.h:322
edm::EventSetup::get
const T & get() const
Definition: EventSetup.h:55
EcalRingCalibrationTools::setCaloGeometry
static void setCaloGeometry(const CaloGeometry *geometry)
Definition: EcalRingCalibrationTools.h:42
ZeeCalibration::h2_xtalMiscalibCoeffEndcapPlus_
TH2F * h2_xtalMiscalibCoeffEndcapPlus_
Definition: ZeeCalibration.h:290
EBDetId::MAX_IETA
static const int MAX_IETA
Definition: EBDetId.h:143
ZeeCalibration::produceEcalIntercalibConstants
virtual boost::shared_ptr< EcalIntercalibConstants > produceEcalIntercalibConstants(const EcalIntercalibConstantsRcd &iRecord)
Produce Ecal interCalibrations.
Definition: ZeeCalibration.cc:184
dtNoiseDBValidation_cfg.cerr
tuple cerr
Definition: dtNoiseDBValidation_cfg.py:22
EcalRingCalibrationTools::getDetIdsInECAL
static std::vector< DetId > getDetIdsInECAL()
Definition: EcalRingCalibrationTools.cc:113
ZeeCalibration::initCalibCoeff
float initCalibCoeff[250]
Definition: ZeeCalibration.h:193
ZeeCalibration::SILVER_ELECTRONS_IN_ENDCAP
int SILVER_ELECTRONS_IN_ENDCAP
Definition: ZeeCalibration.h:337
ZeeCalibration::mcProducer_
std::string mcProducer_
Definition: ZeeCalibration.h:160
AlCaHLTBitMon_ParallelJobs.p
tuple p
Definition: AlCaHLTBitMon_ParallelJobs.py:152
EcalIntercalibConstants
EcalIntercalibConstantMap EcalIntercalibConstants
Definition: EcalIntercalibConstants.h:12
ZeeCalibration::h1_reco_ZMassCorrBB_
TH1F * h1_reco_ZMassCorrBB_
Definition: ZeeCalibration.h:235
calibXMLwriter::writeLine
void writeLine(EBDetId const &, float)
Definition: calibXMLwriter.cc:36
edm::SortedCollection::find
iterator find(key_type k)
Definition: SortedCollection.h:271
perp
T perp() const
Magnitude of transverse component.
Definition: Basic3DVectorLD.h:153
ZeeCalibration::startingNewLoop
virtual void startingNewLoop(unsigned int iLoop)
Called at beginning of loop.
Definition: ZeeCalibration.cc:1487
p1
double p1[4]
Definition: TauolaWrapper.h:89
ZeeCalibration::h1_reco_ZMassCorrEE_
TH1F * h1_reco_ZMassCorrEE_
Definition: ZeeCalibration.h:236
ZeeCalibration::h2_xtalMiscalibCoeffBarrel_
TH2F * h2_xtalMiscalibCoeffBarrel_
Definition: ZeeCalibration.h:288
ZeeCalibration::h2_chi2_
TH2F * h2_chi2_[25]
Definition: ZeeCalibration.h:281
ZeeCalibration::h2_coeffVsLoop_
TH2F * h2_coeffVsLoop_
Definition: ZeeCalibration.h:261
ZeePlots::fillZInfo
void fillZInfo(std::pair< calib::CalibElectron *, calib::CalibElectron * > myZeeCandidate)
Definition: ZeePlots.cc:118
ZeeCalibration::fillEleInfo
void fillEleInfo(std::vector< HepMC::GenParticle * > &a, std::map< HepMC::GenParticle *, const reco::GsfElectron * > &b)
Definition: ZeeCalibration.cc:2066
EcalCondObjectContainer::find
const_iterator find(uint32_t rawId) const
Definition: EcalCondObjectContainer.h:60
ZeeCalibration::barrelfile_
std::string barrelfile_
Definition: ZeeCalibration.h:176
ZeeCalibration::h1_weightSumMeanBarrel_
TH1F * h1_weightSumMeanBarrel_
Definition: ZeeCalibration.h:292
ZeeCalibration::fEtaBarrelGood
double fEtaBarrelGood(double scEta) const
Definition: ZeeCalibration.cc:2011
ZeeCalibration::h_ESCEtrueVsEta_
TH2F * h_ESCEtrueVsEta_[25]
Definition: ZeeCalibration.h:251
ZeeCalibration::h2_zMassVsLoop_
TH2F * h2_zMassVsLoop_
Definition: ZeeCalibration.h:258
ZeePlots::bookZHistograms
void bookZHistograms()
Definition: ZeePlots.cc:89
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
EEDetId::IY_MAX
static const int IY_MAX
Definition: EEDetId.h:306
ZeeCalibration::h2_iterations_
TH2F * h2_iterations_[25]
Definition: ZeeCalibration.h:282
reco::SuperCluster::seed
const CaloClusterPtr & seed() const
seed BasicCluster
Definition: SuperCluster.h:66
ZeeCalibration::ZeeCalibration
ZeeCalibration(const edm::ParameterSet &iConfig)
Constructor.
Definition: ZeeCalibration.cc:78
HcalObjRepresent::Reset
void Reset(std::vector< TH2F > &depth)
Definition: HcalObjRepresent.h:620
ZeeCalibration::NewCalibCoeff
float NewCalibCoeff[250]
Definition: ZeeCalibration.h:191
reco::SuperCluster::preshowerEnergy
double preshowerEnergy() const
energy deposited in preshower
Definition: SuperCluster.h:50
ZeeCalibration::calibCoeffError
float calibCoeffError[250]
Definition: ZeeCalibration.h:192
x
Definition: DDAxes.h:10
ZeeCalibration::h2_xtalRecalibCoeffEndcapPlus_
TH2F * h2_xtalRecalibCoeffEndcapPlus_[25]
Definition: ZeeCalibration.h:286
ZeeCalibration::CRACK_ELECTRONS_IN_ENDCAP
int CRACK_ELECTRONS_IN_ENDCAP
Definition: ZeeCalibration.h:343
ZeeCalibration::CRACK_ELECTRONS_IN_BARREL
int CRACK_ELECTRONS_IN_BARREL
Definition: ZeeCalibration.h:342
ZeeCalibration::h2_fEtaEndcapGood_
TH2F * h2_fEtaEndcapGood_
Definition: ZeeCalibration.h:217
ZeeCalibration::h2_xtalRecalibCoeffEndcapMinus_
TH2F * h2_xtalRecalibCoeffEndcapMinus_[25]
Definition: ZeeCalibration.h:285
ZeePlots::bookEleMCHistograms
void bookEleMCHistograms()
Definition: ZeePlots.cc:186
ZeeCalibration::electronProducer_
std::string electronProducer_
Definition: ZeeCalibration.h:163
ZeeCalibration::getEtaCorrection
double getEtaCorrection(const reco::GsfElectron *)
Definition: ZeeCalibration.cc:2141
ZeeCalibration::aHLTResults
bool aHLTResults[200]
Definition: ZeeCalibration.h:369
ZeeCalibration::h2_fEtaBarrelGood_
TH2F * h2_fEtaBarrelGood_
Definition: ZeeCalibration.h:215
ZeeCalibration::erechitCollection_
std::string erechitCollection_
Definition: ZeeCalibration.h:153
findQualityFiles.size
tuple size
Write out results.
Definition: findQualityFiles.py:442
ZeeCalibration::h_eleEffEta_
TH1F * h_eleEffEta_[2]
Definition: ZeeCalibration.h:247
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
ZeeCalibration::h1_borderElectronClassification_
TH1F * h1_borderElectronClassification_
Definition: ZeeCalibration.h:308
ZeeCalibration::h2_xtalMiscalibCoeffEndcapMinus_
TH2F * h2_xtalMiscalibCoeffEndcapMinus_
Definition: ZeeCalibration.h:289
ZeeCalibration::h2_zWidthVsLoop_
TH2F * h2_zWidthVsLoop_
Definition: ZeeCalibration.h:260
ZeeCalibration::sigmaErrorArray
double sigmaErrorArray[50]
Definition: ZeeCalibration.h:324
ZeeCalibration::h2_residualSigma_
TH2F * h2_residualSigma_
Definition: ZeeCalibration.h:271
ZeeCalibration::SHOWER_ELECTRONS_IN_BARREL
int SHOWER_ELECTRONS_IN_BARREL
Definition: ZeeCalibration.h:339
ZIterativeAlgorithmWithFit::getOptimizedIterations
const std::vector< int > & getOptimizedIterations() const
Definition: ZIterativeAlgorithmWithFit.h:72
ZeeCalibration::TOTAL_ELECTRONS_IN_BARREL
int TOTAL_ELECTRONS_IN_BARREL
Definition: ZeeCalibration.h:330
phi
Definition: DDAxes.h:10
ZeeCalibration::h1_weightSumMeanEndcap_
TH1F * h1_weightSumMeanEndcap_
Definition: ZeeCalibration.h:293
Generated for CMSSW Reference Manual by doxygen 1.8.5