EcalCalibAlgos/interface/ZeeCalibration.h

0001 #ifndef CALIBRATION_ECALCALIBALGOS_ZEECALIBRATION
0002 #define CALIBRATION_ECALCALIBALGOS_ZEECALIBRATION
0003
0004 // -*- C++ -*-
0005 //
0006 // Package:    ZeeCalibration
0007 // Class:      ZeeCalibration
0008 //
0009 /**\class ZeeCalibration ZeeCalibration.cc Calibration/EcalCalibAlgos/src/ZeeCalibration.cc
0010
0011  Description: Perform single electron calibration (tested on TB data only).
0012
0013  Implementation:
0014      <Notes on implementation>
0015 */
0016 //
0017 //
0018 //
0019
0020 // system include files
0021 #include <memory>
0022 #include <vector>
0023 #include <map>
0024 #include <string>
0025
0026 // user include files
0027 #include "FWCore/Framework/interface/LooperFactory.h"
0028 #include "FWCore/Framework/interface/ESProducerLooper.h"
0029 #include "FWCore/Framework/interface/Frameworkfwd.h"
0030 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0031 #include "FWCore/Framework/interface/Event.h"
0032 #include "FWCore/Framework/interface/EventSetup.h"
0033 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0034
0035 #include "Calibration/Tools/interface/ZIterativeAlgorithmWithFit.h"
0036 #include "Calibration/Tools/interface/CalibElectron.h"
0037
0038 #include "Calibration/EcalCalibAlgos/interface/ZeePlots.h"
0039 #include "Calibration/EcalCalibAlgos/interface/ZeeRescaleFactorPlots.h"
0040
0041 #include "CondFormats/EcalObjects/interface/EcalIntercalibConstants.h"
0042 #include "CondFormats/DataRecord/interface/EcalIntercalibConstantsRcd.h"
0043
0044 #include "DataFormats/Common/interface/TriggerResults.h"
0045 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0046 #include "DataFormats/DetId/interface/DetId.h"
0047
0048 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0049 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0050
0051 #include "TTree.h"
0052 #include "TFile.h"
0053 #include "TGraph.h"
0054 #include "TGraphErrors.h"
0055 #include "TH1.h"
0056 #include "TH2.h"
0057
0058 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0059 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0060 #include "DataFormats/EgammaCandidates/interface/ElectronFwd.h"
0061 #include "DataFormats/EgammaCandidates/interface/Electron.h"
0062 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0063 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
0064
0065 // class declaration
0066 //
0067
0068 class ZeeCalibration : public edm::ESProducerLooper {
0069 public:
0070   /// Constructor
0071   ZeeCalibration(const edm::ParameterSet& iConfig);
0072
0073   /// Destructor
0074   ~ZeeCalibration() override;
0075
0076   /// Dummy implementation (job done in duringLoop)
0077   virtual void produce(edm::Event&, const edm::EventSetup&) {}
0078
0079   /// Called at beginning of job
0080   void beginOfJob() override;
0081
0082   /// Called at end of job
0083   void endOfJob() override;
0084
0085   /// Called at beginning of loop
0086   void startingNewLoop(unsigned int iLoop) override;
0087
0088   /// Called at end of loop
0089   Status endOfLoop(const edm::EventSetup&, unsigned int iLoop) override;
0090
0091   /// Called at each event
0092   Status duringLoop(const edm::Event&, const edm::EventSetup&) override;
0093
0094   /// Produce Ecal interCalibrations
0095   virtual std::shared_ptr<EcalIntercalibConstants> produceEcalIntercalibConstants(
0096       const EcalIntercalibConstantsRcd& iRecord);
0097
0098 private:
0099   /*   ElectronEnergyCorrector myCorrector; */
0100   /*   ElectronClassification myClassificator; */
0101
0102   double fEtaBarrelBad(double scEta) const;
0103   double fEtaBarrelGood(double scEta) const;
0104   double fEtaEndcapBad(double scEta) const;
0105   double fEtaEndcapGood(double scEta) const;
0106
0107   int ringNumberCorrector(int k);
0108   double getEtaCorrection(const reco::GsfElectron*);
0109
0110   void fillEleInfo(std::vector<HepMC::GenParticle*>& a, std::map<HepMC::GenParticle*, const reco::GsfElectron*>& b);
0111   void fillMCInfo(HepMC::GenParticle* mcele);
0112
0113   void fillMCmap(const std::vector<const reco::GsfElectron*>* electronCollection,
0114                  const std::vector<HepMC::GenParticle*>& mcEle,
0115                  std::map<HepMC::GenParticle*, const reco::GsfElectron*>& myMCmap);
0116   //  void fillMCmap(const reco::ElectronCollection* electronCollection, const std::vector<HepMC::GenParticle*>& mcEle,std::map<HepMC::GenParticle*,const reco::Electron*>& myMCmap);
0117
0118   float EvalDPhi(float Phi, float Phi_ref);
0119   float EvalDR(float Eta, float Eta_ref, float Phi, float Phi_ref);
0120
0121   void bookHistograms();
0122
0123   void resetVariables();
0124
0125   void resetHistograms();
0126
0127   void printStatistics();
0128
0129   std::pair<DetId, double> getHottestDetId(const std::vector<std::pair<DetId, float> >& mySCRecHits,
0130                                            const EBRecHitCollection* ebhits,
0131                                            const EERecHitCollection* eehits);
0132
0133   bool xtalIsOnModuleBorder(EBDetId myEBDetId);
0134
0135   float computeCoefficientDistanceAtIteration(float v1[250], float v2[250], int size);
0136
0137   //  float Calculate_SigmaEtaEta(const reco::SuperCluster &passedCluster);
0138
0139   // ----------member data ---------------------------
0140
0141   TTree* myTree;
0142
0143   std::string outputFileName_;
0144
0145   const edm::InputTag hlTriggerResults_;
0146
0147   const std::string mcProducer_;
0148   const std::string rechitProducer_;
0149   const std::string rechitCollection_;
0150   const std::string erechitProducer_;
0151   const std::string erechitCollection_;
0152   const std::string scProducer_;
0153   const std::string scCollection_;
0154
0155   const std::string scIslandProducer_;
0156   const std::string scIslandCollection_;
0157
0158   const std::string electronProducer_;
0159   const std::string electronCollection_;
0160
0161   std::string calibMode_;
0162
0163   const edm::EDGetTokenT<edm::TriggerResults> trigResultsToken_;
0164   const edm::EDGetTokenT<edm::HepMCProduct> hepMCToken_;
0165   const edm::EDGetTokenT<EBRecHitCollection> ebRecHitToken_;
0166   const edm::EDGetTokenT<EERecHitCollection> eeRecHitToken_;
0167   const edm::EDGetTokenT<reco::SuperClusterCollection> scToken_;
0168   const edm::EDGetTokenT<reco::SuperClusterCollection> islandSCToken_;
0169   const edm::EDGetTokenT<reco::GsfElectronCollection> gsfElectronToken_;
0170
0171   const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> geometryToken_;
0172
0173   unsigned int etaBins_;
0174   unsigned int etBins_;
0175
0176   double etaMin_;
0177   double etMin_;
0178   double etaMax_;
0179   double etMax_;
0180
0181   std::string barrelfile_;
0182   std::string endcapfile_;
0183
0184   double minInvMassCut_;
0185   double maxInvMassCut_;
0186   double mass;
0187
0188   float mass4tree;
0189   float massDiff4tree;
0190
0191   int read_events;
0192
0193   int loopFlag_;
0194
0195   float calibCoeff[nMaxChannels];
0196   float NewCalibCoeff[nMaxChannels];
0197   float calibCoeffError[nMaxChannels];
0198   float initCalibCoeff[nMaxChannels];
0199
0200   std::shared_ptr<EcalIntercalibConstants> ical;
0201
0202   ZIterativeAlgorithmWithFit* theAlgorithm_;
0203
0204   ZeePlots* myZeePlots_;
0205   ZeeRescaleFactorPlots* myZeeRescaleFactorPlots_;
0206
0207   // steering parameters
0208
0209   edm::ParameterSet theParameterSet;
0210
0211   //  TGraph* graph;
0212
0213   TH1F* h1_eventsBeforeEWKSelection_;
0214   TH1F* h1_eventsAfterEWKSelection_;
0215
0216   TH1F* h1_eventsBeforeBorderSelection_;
0217   TH1F* h1_eventsAfterBorderSelection_;
0218
0219   TH2F* h2_fEtaBarrelGood_;
0220   TH2F* h2_fEtaBarrelBad_;
0221   TH2F* h2_fEtaEndcapGood_;
0222   TH2F* h2_fEtaEndcapBad_;
0223   TH1F* h1_nEleReco_;
0224   TH1F* h1_eleClasses_;
0225
0226   TH1F* h_eleEffEta[2];
0227   TH1F* h_eleEffPhi[2];
0228   TH1F* h_eleEffPt[2];
0229
0230   TH1F* h1_seedOverSC_;
0231   TH1F* h1_preshowerOverSC_;
0232
0233   TH1F* h1_zMassResol_;
0234   TH1F* h1_zEtaResol_;
0235   TH1F* h1_zPhiResol_;
0236   TH1F* h1_reco_ZMass_;
0237
0238   TH1F* h1_reco_ZMassCorr_;
0239   TH1F* h1_reco_ZMassCorrBB_;
0240   TH1F* h1_reco_ZMassCorrEE_;
0241   TH1F* h1_reco_ZMassGood_;
0242   TH1F* h1_reco_ZMassBad_;
0243   TH1F* h1_ZCandMult_;
0244   TH1F* h1_RMin_;
0245   TH1F* h1_RMinZ_;
0246   TH1F* h1_eleERecoOverEtrue_;
0247
0248   TH1F* h1_eleEtaResol_;
0249   TH1F* h1_elePhiResol_;
0250
0251   TH1F* h_eleEffEta_[2];
0252   TH1F* h_eleEffPhi_[2];
0253   TH1F* h_eleEffPt_[2];
0254   TH1F* h_ESCEtrue_[25];
0255   TH2F* h_ESCEtrueVsEta_[25];
0256
0257   TH1F* h_ESCcorrEtrue_[25];
0258   TH2F* h_ESCcorrEtrueVsEta_[25];
0259
0260   TH2F* h2_coeffVsEta_;
0261   TH2F* h2_coeffVsEtaGrouped_;
0262   TH2F* h2_zMassVsLoop_;
0263   TH2F* h2_zMassDiffVsLoop_;
0264   TH2F* h2_zWidthVsLoop_;
0265   TH2F* h2_coeffVsLoop_;
0266
0267   TH2F* h2_miscalRecal_;
0268   //  TH2F* h2_miscalRecalParz_[25];
0269   TH1F* h1_mc_;
0270   TH1F* h1_mcParz_[25];
0271   /*
0272   TH1F* h_DiffZMassDistr_[25];
0273   TH1F* h_ZMassDistr_[25];
0274   */
0275   TH2F* h2_residualSigma_;
0276   TH2F* h2_miscalRecalEB_;
0277   //TH2F* h2_miscalRecalEBParz_[25];
0278   TH1F* h1_mcEB_;
0279   TH1F* h1_mcEBParz_[25];
0280   TH2F* h2_miscalRecalEE_;
0281   //TH2F* h2_miscalRecalEEParz_[25];
0282   TH1F* h1_mcEE_;
0283   TH1F* h1_mcEEParz_[25];
0284
0285   TH2F* h2_chi2_[25];
0286   TH2F* h2_iterations_[25];
0287
0288   TH2F* h2_xtalRecalibCoeffBarrel_[25];
0289   TH2F* h2_xtalRecalibCoeffEndcapMinus_[25];
0290   TH2F* h2_xtalRecalibCoeffEndcapPlus_[25];
0291
0292   TH2F* h2_xtalMiscalibCoeffBarrel_;
0293   TH2F* h2_xtalMiscalibCoeffEndcapMinus_;
0294   TH2F* h2_xtalMiscalibCoeffEndcapPlus_;
0295
0296   TH1F* h1_weightSumMeanBarrel_;
0297   TH1F* h1_weightSumMeanEndcap_;
0298
0299   TH1F* h1_occupancyVsEta_;
0300   TH1F* h1_occupancyVsEtaGold_;
0301   TH1F* h1_occupancyVsEtaSilver_;
0302   TH1F* h1_occupancyVsEtaCrack_;
0303   TH1F* h1_occupancyVsEtaShower_;
0304   TH1F* h1_occupancy_;
0305   TH1F* h1_occupancyBarrel_;
0306   TH1F* h1_occupancyEndcap_;
0307
0308   TH1F* h1_electronCosTheta_TK_;
0309   TH1F* h1_electronCosTheta_SC_;
0310   TH1F* h1_electronCosTheta_SC_TK_;
0311
0312   TH1F* h1_borderElectronClassification_;
0313
0314   Int_t BBZN, EBZN, EEZN, BBZN_gg, EBZN_gg, EEZN_gg, BBZN_tt, EBZN_tt, EEZN_tt, BBZN_t0, EBZN_t0, EEZN_t0;
0315   Int_t NEVT, MCZBB, MCZEB, MCZEE;
0316
0317   TFile* outputFile_;
0318
0319   unsigned int theMaxLoops;  // Number of loops to loop
0320
0321   bool wantEtaCorrection_;
0322
0323   unsigned int electronSelection_;
0324
0325   double loopArray[50];
0326   double sigmaArray[50];
0327   double sigmaErrorArray[50];
0328   double coefficientDistanceAtIteration[50];
0329
0330   int BARREL_ELECTRONS_BEFORE_BORDER_CUT;
0331   int BARREL_ELECTRONS_AFTER_BORDER_CUT;
0332
0333   int TOTAL_ELECTRONS_IN_BARREL;
0334   int TOTAL_ELECTRONS_IN_ENDCAP;
0335
0336   int GOLDEN_ELECTRONS_IN_BARREL;
0337   int GOLDEN_ELECTRONS_IN_ENDCAP;
0338
0339   int SILVER_ELECTRONS_IN_BARREL;
0340   int SILVER_ELECTRONS_IN_ENDCAP;
0341
0342   int SHOWER_ELECTRONS_IN_BARREL;
0343   int SHOWER_ELECTRONS_IN_ENDCAP;
0344
0345   int CRACK_ELECTRONS_IN_BARREL;
0346   int CRACK_ELECTRONS_IN_ENDCAP;
0347
0348   unsigned int nEvents_;  // number of events processed
0349
0350   unsigned int nWasRun_;  // # where at least one HLT was run
0351   unsigned int nAccept_;  // # of accepted events
0352   unsigned int nErrors_;  // # where at least one HLT had error
0353
0354   std::vector<unsigned int> hlWasRun_;  // # where HLT[i] was run
0355   std::vector<unsigned int> hlAccept_;  // # of events accepted by HLT[i]
0356   std::vector<unsigned int> hlErrors_;  // # of events with error in HLT[i]
0357
0358   std::vector<std::string> hlNames_;  // name of each HLT algorithm
0359   bool init_;                         // vectors initialised or not
0360
0361   Int_t triggerCount;
0362   char aTriggerNames[200][30];
0363   bool aTriggerResults[200];
0364
0365   Int_t hltCount;
0366   char aHLTNames[6000];
0367   Int_t hltNamesLen;
0368   TString aNames[200];
0369   bool aHLTResults[200];
0370
0371   bool isfirstcall_;
0372 };
0373 #endif