Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:09:21 

0001 #include "DQMOffline/EGamma/plugins/ElectronTagProbeAnalyzer.h"
0002 
0003 #include "DQMServices/Core/interface/DQMStore.h"
0004 
0005 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0006 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0007 #include "DataFormats/EgammaReco/interface/BasicClusterFwd.h"
0008 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0009 #include "DataFormats/VertexReco/interface/Vertex.h"
0010 #include "DataFormats/Common/interface/Handle.h"
0011 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0012 #include "DataFormats/Common/interface/TriggerResults.h"
0013 
0014 #include "FWCore/Common/interface/TriggerNames.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/Framework/interface/Event.h"
0017 #include "FWCore/Framework/interface/MakerMacros.h"
0018 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0019 #include "FWCore/ServiceRegistry/interface/Service.h"
0020 
0021 #include "CLHEP/Units/GlobalPhysicalConstants.h"
0022 #include "TMath.h"
0023 
0024 #include <iostream>
0025 #include <vector>
0026 
0027 using namespace reco;
0028 
0029 ElectronTagProbeAnalyzer::ElectronTagProbeAnalyzer(const edm::ParameterSet& conf) : ElectronDqmAnalyzerBase(conf) {
0030   // general, collections
0031   Selection_ = conf.getParameter<int>("Selection");
0032   electronCollection_ = consumes<GsfElectronCollection>(conf.getParameter<edm::InputTag>("ElectronCollection"));
0033   matchingObjectCollection_ =
0034       consumes<SuperClusterCollection>(conf.getParameter<edm::InputTag>("MatchingObjectCollection"));
0035   trackCollection_ = consumes<TrackCollection>(conf.getParameter<edm::InputTag>("TrackCollection"));
0036   vertexCollection_ = consumes<VertexCollection>(conf.getParameter<edm::InputTag>("VertexCollection"));
0037   gsftrackCollection_ = consumes<GsfTrackCollection>(conf.getParameter<edm::InputTag>("GsfTrackCollection"));
0038   beamSpotTag_ = consumes<BeamSpot>(conf.getParameter<edm::InputTag>("BeamSpot"));
0039   readAOD_ = conf.getParameter<bool>("ReadAOD");
0040 
0041   // tag and probe
0042   massLow_ = conf.getParameter<double>("MassLow");
0043   massHigh_ = conf.getParameter<double>("MassHigh");
0044   TPchecksign_ = conf.getParameter<bool>("TpCheckSign");
0045   TAGcheckclass_ = conf.getParameter<bool>("TagCheckClass");
0046   PROBEetcut_ = conf.getParameter<bool>("ProbeEtCut");
0047   PROBEcheckclass_ = conf.getParameter<bool>("ProbeCheckClass");
0048 
0049   // electron selection
0050   minEt_ = conf.getParameter<double>("MinEt");
0051   minPt_ = conf.getParameter<double>("MinPt");
0052   maxAbsEta_ = conf.getParameter<double>("MaxAbsEta");
0053   isEB_ = conf.getParameter<bool>("SelectEb");
0054   isEE_ = conf.getParameter<bool>("SelectEe");
0055   isNotEBEEGap_ = conf.getParameter<bool>("SelectNotEbEeGap");
0056   isEcalDriven_ = conf.getParameter<bool>("SelectEcalDriven");
0057   isTrackerDriven_ = conf.getParameter<bool>("SelectTrackerDriven");
0058   eOverPMinBarrel_ = conf.getParameter<double>("MinEopBarrel");
0059   eOverPMaxBarrel_ = conf.getParameter<double>("MaxEopBarrel");
0060   eOverPMinEndcaps_ = conf.getParameter<double>("MinEopEndcaps");
0061   eOverPMaxEndcaps_ = conf.getParameter<double>("MaxEopEndcaps");
0062   dEtaMinBarrel_ = conf.getParameter<double>("MinDetaBarrel");
0063   dEtaMaxBarrel_ = conf.getParameter<double>("MaxDetaBarrel");
0064   dEtaMinEndcaps_ = conf.getParameter<double>("MinDetaEndcaps");
0065   dEtaMaxEndcaps_ = conf.getParameter<double>("MaxDetaEndcaps");
0066   dPhiMinBarrel_ = conf.getParameter<double>("MinDphiBarrel");
0067   dPhiMaxBarrel_ = conf.getParameter<double>("MaxDphiBarrel");
0068   dPhiMinEndcaps_ = conf.getParameter<double>("MinDphiEndcaps");
0069   dPhiMaxEndcaps_ = conf.getParameter<double>("MaxDphiEndcaps");
0070   sigIetaIetaMinBarrel_ = conf.getParameter<double>("MinSigIetaIetaBarrel");
0071   sigIetaIetaMaxBarrel_ = conf.getParameter<double>("MaxSigIetaIetaBarrel");
0072   sigIetaIetaMinEndcaps_ = conf.getParameter<double>("MinSigIetaIetaEndcaps");
0073   sigIetaIetaMaxEndcaps_ = conf.getParameter<double>("MaxSigIetaIetaEndcaps");
0074   hadronicOverEmMaxBarrel_ = conf.getParameter<double>("MaxHoeBarrel");
0075   hadronicOverEmMaxEndcaps_ = conf.getParameter<double>("MaxHoeEndcaps");
0076   mvaMin_ = conf.getParameter<double>("MinMva");
0077   tipMaxBarrel_ = conf.getParameter<double>("MaxTipBarrel");
0078   tipMaxEndcaps_ = conf.getParameter<double>("MaxTipEndcaps");
0079   tkIso03Max_ = conf.getParameter<double>("MaxTkIso03");
0080   hcalIso03Depth1MaxBarrel_ = conf.getParameter<double>("MaxHcalIso03Depth1Barrel");
0081   hcalIso03Depth1MaxEndcaps_ = conf.getParameter<double>("MaxHcalIso03Depth1Endcaps");
0082   hcalIso03Depth2MaxEndcaps_ = conf.getParameter<double>("MaxHcalIso03Depth2Endcaps");
0083   ecalIso03MaxBarrel_ = conf.getParameter<double>("MaxEcalIso03Barrel");
0084   ecalIso03MaxEndcaps_ = conf.getParameter<double>("MaxEcalIso03Endcaps");
0085 
0086   // for trigger
0087   triggerResults_ = conf.getParameter<edm::InputTag>("TriggerResults");
0088   //  HLTPathsByName_= conf.getParameter<std::vector<std::string > >("HltPaths");
0089   //  HLTPathsByIndex_.resize(HLTPathsByName_.size());
0090 
0091   // histos limits and binning
0092   nbineta = conf.getParameter<int>("NbinEta");
0093   nbineta2D = conf.getParameter<int>("NbinEta2D");
0094   etamin = conf.getParameter<double>("EtaMin");
0095   etamax = conf.getParameter<double>("EtaMax");
0096   //
0097   nbinphi = conf.getParameter<int>("NbinPhi");
0098   nbinphi2D = conf.getParameter<int>("NbinPhi2D");
0099   phimin = conf.getParameter<double>("PhiMin");
0100   phimax = conf.getParameter<double>("PhiMax");
0101   //
0102   nbinpt = conf.getParameter<int>("NbinPt");
0103   nbinpteff = conf.getParameter<int>("NbinPtEff");
0104   nbinpt2D = conf.getParameter<int>("NbinPt2D");
0105   ptmax = conf.getParameter<double>("PtMax");
0106   //
0107   nbinp = conf.getParameter<int>("NbinP");
0108   nbinp2D = conf.getParameter<int>("NbinP2D");
0109   pmax = conf.getParameter<double>("PMax");
0110   //
0111   nbineop = conf.getParameter<int>("NbinEop");
0112   nbineop2D = conf.getParameter<int>("NbinEop2D");
0113   eopmax = conf.getParameter<double>("EopMax");
0114   eopmaxsht = conf.getParameter<double>("EopMaxSht");
0115   //
0116   nbindeta = conf.getParameter<int>("NbinDeta");
0117   detamin = conf.getParameter<double>("DetaMin");
0118   detamax = conf.getParameter<double>("DetaMax");
0119   //
0120   nbindphi = conf.getParameter<int>("NbinDphi");
0121   dphimin = conf.getParameter<double>("DphiMin");
0122   dphimax = conf.getParameter<double>("DphiMax");
0123   //
0124   nbindetamatch = conf.getParameter<int>("NbinDetaMatch");
0125   nbindetamatch2D = conf.getParameter<int>("NbinDetaMatch2D");
0126   detamatchmin = conf.getParameter<double>("DetaMatchMin");
0127   detamatchmax = conf.getParameter<double>("DetaMatchMax");
0128   //
0129   nbindphimatch = conf.getParameter<int>("NbinDphiMatch");
0130   nbindphimatch2D = conf.getParameter<int>("NbinDphiMatch2D");
0131   dphimatchmin = conf.getParameter<double>("DphiMatchMin");
0132   dphimatchmax = conf.getParameter<double>("DphiMatchMax");
0133   //
0134   nbinfhits = conf.getParameter<int>("NbinFhits");
0135   fhitsmax = conf.getParameter<double>("FhitsMax");
0136   //
0137   nbinlhits = conf.getParameter<int>("NbinLhits");
0138   lhitsmax = conf.getParameter<double>("LhitsMax");
0139   //
0140   nbinxyz = conf.getParameter<int>("NbinXyz");
0141   nbinxyz2D = conf.getParameter<int>("NbinXyz2D");
0142   //
0143   nbinpoptrue = conf.getParameter<int>("NbinPopTrue");
0144   poptruemin = conf.getParameter<double>("PopTrueMin");
0145   poptruemax = conf.getParameter<double>("PopTrueMax");
0146   //
0147   nbinmee = conf.getParameter<int>("NbinMee");
0148   meemin = conf.getParameter<double>("MeeMin");
0149   meemax = conf.getParameter<double>("MeeMax");
0150   //
0151   nbinhoe = conf.getParameter<int>("NbinHoe");
0152   hoemin = conf.getParameter<double>("HoeMin");
0153   hoemax = conf.getParameter<double>("HoeMax");
0154 }
0155 
0156 ElectronTagProbeAnalyzer::~ElectronTagProbeAnalyzer() {}
0157 
0158 void ElectronTagProbeAnalyzer::bookHistograms(DQMStore::IBooker& iBooker, edm::Run const&, edm::EventSetup const&) {
0159   iBooker.setCurrentFolder(outputInternalPath_);
0160   nEvents_ = 0;
0161   //nAfterTrigger_ = 0 ;
0162 
0163   // basic quantities
0164   h1_vertexPt_barrel = bookH1(
0165       iBooker, "vertexPt_barrel", "ele transverse momentum in barrel", nbinpt, 0., ptmax, "p_{T vertex} (GeV/c)");
0166   h1_vertexPt_endcaps = bookH1(
0167       iBooker, "vertexPt_endcaps", "ele transverse momentum in endcaps", nbinpt, 0., ptmax, "p_{T vertex} (GeV/c)");
0168   h1_vertexEta = bookH1(iBooker, "vertexEta", "ele momentum #eta", nbineta, etamin, etamax, "#eta");
0169   h2_vertexEtaVsPhi = bookH2(iBooker,
0170                              "vertexEtaVsPhi",
0171                              "ele momentum #eta vs #phi",
0172                              nbineta2D,
0173                              etamin,
0174                              etamax,
0175                              nbinphi2D,
0176                              phimin,
0177                              phimax,
0178                              "#eta",
0179                              "#phi (rad)");
0180   h2_vertexXvsY = bookH2(
0181       iBooker, "vertexXvsY", "ele vertex x vs y", nbinxyz2D, -0.1, 0.1, nbinxyz2D, -0.1, 0.1, "x (cm)", "y (cm)");
0182   h1_vertexZ = bookH1(iBooker, "vertexZ", "ele vertex z", nbinxyz, -25, 25, "z (cm)");
0183 
0184   // super-clusters
0185   //  h1_sclPhi = bookH1(iBooker, "sclPhi","ele supercluster phi",nbinphi,phimin,phimax);
0186   h1_sclEt = bookH1(iBooker, "sclEt", "ele supercluster transverse energy", nbinpt, 0., ptmax);
0187 
0188   // electron track
0189   h1_chi2 = bookH1(iBooker, "chi2", "ele track #chi^{2}", 100, 0., 15., "#Chi^{2}");
0190   h1_foundHits = bookH1(iBooker, "foundHits", "ele track # found hits", nbinfhits, 0., fhitsmax, "N_{hits}");
0191   h1_lostHits = bookH1(iBooker, "lostHits", "ele track # lost hits", 5, 0., 5., "N_{lost hits}");
0192 
0193   // electron matching and ID
0194   h1_Eop_barrel = bookH1(iBooker, "Eop_barrel", "ele E/P_{vertex} in barrel", nbineop, 0., eopmax, "E/P_{vertex}");
0195   h1_Eop_endcaps = bookH1(iBooker, "Eop_endcaps", "ele E/P_{vertex} in endcaps", nbineop, 0., eopmax, "E/P_{vertex}");
0196   h1_EeleOPout_barrel =
0197       bookH1(iBooker, "EeleOPout_barrel", "ele E_{ele}/P_{out} in barrel", nbineop, 0., eopmax, "E_{ele}/P_{out}");
0198   h1_EeleOPout_endcaps =
0199       bookH1(iBooker, "EeleOPout_endcaps", "ele E_{ele}/P_{out} in endcaps", nbineop, 0., eopmax, "E_{ele}/P_{out}");
0200   h1_dEtaSc_propVtx_barrel = bookH1(iBooker,
0201                                     "dEtaSc_propVtx_barrel",
0202                                     "ele #eta_{sc} - #eta_{tr}, prop from vertex, in barrel",
0203                                     nbindetamatch,
0204                                     detamatchmin,
0205                                     detamatchmax,
0206                                     "#eta_{sc} - #eta_{tr}");
0207   h1_dEtaSc_propVtx_endcaps = bookH1(iBooker,
0208                                      "dEtaSc_propVtx_endcaps",
0209                                      "ele #eta_{sc} - #eta_{tr}, prop from vertex, in endcaps",
0210                                      nbindetamatch,
0211                                      detamatchmin,
0212                                      detamatchmax,
0213                                      "#eta_{sc} - #eta_{tr}");
0214   h1_dEtaEleCl_propOut_barrel = bookH1(iBooker,
0215                                        "dEtaEleCl_propOut_barrel",
0216                                        "ele #eta_{EleCl} - #eta_{tr}, prop from outermost, in barrel",
0217                                        nbindetamatch,
0218                                        detamatchmin,
0219                                        detamatchmax,
0220                                        "#eta_{elecl} - #eta_{tr}");
0221   h1_dEtaEleCl_propOut_endcaps = bookH1(iBooker,
0222                                         "dEtaEleCl_propOut_endcaps",
0223                                         "ele #eta_{EleCl} - #eta_{tr}, prop from outermost, in endcaps",
0224                                         nbindetamatch,
0225                                         detamatchmin,
0226                                         detamatchmax,
0227                                         "#eta_{elecl} - #eta_{tr}");
0228   h1_dPhiSc_propVtx_barrel = bookH1(iBooker,
0229                                     "dPhiSc_propVtx_barrel",
0230                                     "ele #phi_{sc} - #phi_{tr}, prop from vertex, in barrel",
0231                                     nbindphimatch,
0232                                     dphimatchmin,
0233                                     dphimatchmax,
0234                                     "#phi_{sc} - #phi_{tr} (rad)");
0235   h1_dPhiSc_propVtx_endcaps = bookH1(iBooker,
0236                                      "dPhiSc_propVtx_endcaps",
0237                                      "ele #phi_{sc} - #phi_{tr}, prop from vertex, in endcaps",
0238                                      nbindphimatch,
0239                                      dphimatchmin,
0240                                      dphimatchmax,
0241                                      "#phi_{sc} - #phi_{tr} (rad)");
0242   h1_dPhiEleCl_propOut_barrel = bookH1(iBooker,
0243                                        "dPhiEleCl_propOut_barrel",
0244                                        "ele #phi_{EleCl} - #phi_{tr}, prop from outermost, in barrel",
0245                                        nbindphimatch,
0246                                        dphimatchmin,
0247                                        dphimatchmax,
0248                                        "#phi_{elecl} - #phi_{tr} (rad)");
0249   h1_dPhiEleCl_propOut_endcaps = bookH1(iBooker,
0250                                         "dPhiEleCl_propOut_endcaps",
0251                                         "ele #phi_{EleCl} - #phi_{tr}, prop from outermost, in endcaps",
0252                                         nbindphimatch,
0253                                         dphimatchmin,
0254                                         dphimatchmax,
0255                                         "#phi_{elecl} - #phi_{tr} (rad)");
0256   h1_Hoe_barrel = bookH1(iBooker,
0257                          "Hoe_barrel",
0258                          "ele hadronic energy / em energy, in barrel",
0259                          nbinhoe,
0260                          hoemin,
0261                          hoemax,
0262                          "H/E",
0263                          "Events",
0264                          "ELE_LOGY E1 P");
0265   h1_Hoe_endcaps = bookH1(iBooker,
0266                           "Hoe_endcaps",
0267                           "ele hadronic energy / em energy, in endcaps",
0268                           nbinhoe,
0269                           hoemin,
0270                           hoemax,
0271                           "H/E",
0272                           "Events",
0273                           "ELE_LOGY E1 P");
0274   h1_sclSigEtaEta_barrel =
0275       bookH1(iBooker, "sclSigEtaEta_barrel", "ele supercluster sigma eta eta in barrel", 100, 0., 0.05);
0276   h1_sclSigEtaEta_endcaps =
0277       bookH1(iBooker, "sclSigEtaEta_endcaps", "ele supercluster sigma eta eta in endcaps", 100, 0., 0.05);
0278 
0279   // fbrem
0280   h1_fbrem = bookH1(iBooker, "fbrem", "ele brem fraction", 100, 0., 1., "P_{in} - P_{out} / P_{in}");
0281   h1_classes = bookH1(iBooker, "classes", "ele electron classes", 10, 0.0, 10.);
0282 
0283   // pflow
0284   h1_mva = bookH1(iBooker, "mva", "ele identification mva", 100, -1., 1.);
0285   h1_provenance = bookH1(iBooker, "provenance", "ele provenance", 5, -2., 3.);
0286 
0287   // isolation
0288   h1_tkSumPt_dr03 = bookH1(iBooker,
0289                            "tkSumPt_dr03",
0290                            "tk isolation sum, dR=0.3",
0291                            100,
0292                            0.0,
0293                            20.,
0294                            "TkIsoSum, cone 0.3 (GeV/c)",
0295                            "Events",
0296                            "ELE_LOGY E1 P");
0297   h1_ecalRecHitSumEt_dr03 = bookH1(iBooker,
0298                                    "ecalRecHitSumEt_dr03",
0299                                    "ecal isolation sum, dR=0.3",
0300                                    100,
0301                                    0.0,
0302                                    20.,
0303                                    "EcalIsoSum, cone 0.3 (GeV)",
0304                                    "Events",
0305                                    "ELE_LOGY E1 P");
0306   h1_hcalTowerSumEt_dr03 = bookH1(iBooker,
0307                                   "hcalTowerSumEt_dr03",
0308                                   "hcal isolation sum, dR=0.3",
0309                                   100,
0310                                   0.0,
0311                                   20.,
0312                                   "HcalIsoSum, cone 0.3 (GeV)",
0313                                   "Events",
0314                                   "ELE_LOGY E1 P");
0315 
0316   // di-electron mass
0317   setBookIndex(200);
0318   h1_mee = bookH1(iBooker, "mesc", "Tag ele Probe SC invariant mass", nbinmee, meemin, meemax, "m_{eSC} (GeV/c^{2})");
0319   h1_mee_os = bookH1(iBooker,
0320                      "mee_os",
0321                      "ele pairs invariant mass, opposite sign",
0322                      nbinmee,
0323                      meemin,
0324                      meemax,
0325                      "m_{e^{+}e^{-}} (GeV/c^{2})");
0326 
0327   //===========================
0328   // histos for matching and matched matched objects
0329   //===========================
0330 
0331   // matching object
0332   std::string matchingObjectType;
0333   Labels l;
0334   labelsForToken(matchingObjectCollection_, l);
0335   if (std::string::npos != std::string(l.module).find("SuperCluster", 0)) {
0336     matchingObjectType = "SC";
0337   }
0338   if (matchingObjectType.empty()) {
0339     edm::LogError("ElectronMcFakeValidator::beginJob") << "Unknown matching object type !";
0340   } else {
0341     edm::LogInfo("ElectronMcFakeValidator::beginJob") << "Matching object type: " << matchingObjectType;
0342   }
0343 
0344   // matching object distributions
0345   h1_matchingObject_Eta = bookH1withSumw2(
0346       iBooker, "matchingObject_Eta", matchingObjectType + " #eta", nbineta, etamin, etamax, "#eta_{SC}");
0347   h1_matchingObject_Pt =
0348       bookH1withSumw2(iBooker, "matchingObject_Pt", matchingObjectType + " pt", nbinpteff, 5., ptmax);
0349   h1_matchingObject_Phi =
0350       bookH1withSumw2(iBooker, "matchingObject_Phi", matchingObjectType + " phi", nbinphi, phimin, phimax);
0351   //h1_matchingObject_Z = bookH1withSumw2(iBooker, "matchingObject_Z",matchingObjectType+" z",nbinxyz,-25,25);
0352 
0353   h1_matchedObject_Eta =
0354       bookH1withSumw2(iBooker, "matchedObject_Eta", "Efficiency vs matching SC #eta", nbineta, etamin, etamax);
0355   h1_matchedObject_Pt =
0356       bookH1withSumw2(iBooker, "matchedObject_Pt", "Efficiency vs matching SC E_{T}", nbinpteff, 5., ptmax);
0357   h1_matchedObject_Phi =
0358       bookH1withSumw2(iBooker, "matchedObject_Phi", "Efficiency vs matching SC phi", nbinphi, phimin, phimax);
0359   //h1_matchedObject_Z = bookH1withSumw2(iBooker, "matchedObject_Z","Efficiency vs matching SC z",nbinxyz,-25,25);
0360 }
0361 
0362 void ElectronTagProbeAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0363   nEvents_++;
0364 
0365   edm::Handle<GsfElectronCollection> gsfElectrons;
0366   iEvent.getByToken(electronCollection_, gsfElectrons);
0367   edm::Handle<reco::SuperClusterCollection> recoClusters;
0368   iEvent.getByToken(matchingObjectCollection_, recoClusters);
0369   edm::Handle<reco::TrackCollection> tracks;
0370   iEvent.getByToken(trackCollection_, tracks);
0371   edm::Handle<reco::GsfTrackCollection> gsfTracks;
0372   iEvent.getByToken(gsftrackCollection_, gsfTracks);
0373   edm::Handle<reco::VertexCollection> vertices;
0374   iEvent.getByToken(vertexCollection_, vertices);
0375   edm::Handle<reco::BeamSpot> recoBeamSpotHandle;
0376   iEvent.getByToken(beamSpotTag_, recoBeamSpotHandle);
0377   const BeamSpot bs = *recoBeamSpotHandle;
0378 
0379   edm::EventNumber_t ievt = iEvent.id().event();
0380   edm::RunNumber_t irun = iEvent.id().run();
0381   edm::LuminosityBlockNumber_t ils = iEvent.luminosityBlock();
0382 
0383   edm::LogInfo("ElectronMcSignalValidator::analyze")
0384       << "Treating " << gsfElectrons.product()->size() << " electrons"
0385       << " from event " << ievt << " in run " << irun << " and lumiblock " << ils;
0386   //h1_num_->Fill((*gsfElectrons).size()) ;
0387 
0388   std::vector<std::pair<double, double> > TTCheck;
0389   std::vector<std::pair<double, double> > TTscCheck;
0390 
0391   // selected rec electrons
0392   reco::GsfElectronCollection::const_iterator gsfIter;
0393   for (gsfIter = gsfElectrons->begin(); gsfIter != gsfElectrons->end(); gsfIter++) {
0394     // vertex TIP
0395     double vertexTIP = (gsfIter->vertex().x() - bs.position().x()) * (gsfIter->vertex().x() - bs.position().x()) +
0396                        (gsfIter->vertex().y() - bs.position().y()) * (gsfIter->vertex().y() - bs.position().y());
0397     vertexTIP = sqrt(vertexTIP);
0398 
0399     // select electrons
0400     if (!selected(gsfIter, vertexTIP))
0401       continue;
0402 
0403     reco::SuperClusterRef sclTagRef = gsfIter->superCluster();
0404     reco::SuperClusterCollection::const_iterator moIter;
0405     for (moIter = recoClusters->begin(); moIter != recoClusters->end(); moIter++) {
0406       if (moIter->eta() == sclTagRef->eta())
0407         continue;
0408 
0409       /*
0410       if
0411         ( moIter->energy()/cosh(moIter->eta())>maxPtMatchingObject_ ||
0412           std::abs(moIter->eta())> maxAbsEtaMatchingObject_ )
0413         { continue ; }
0414       */
0415 
0416       // Additional cuts on Tag
0417 
0418       // Additional cuts on Probe
0419       if (PROBEetcut_ && (moIter->energy() / cosh(moIter->eta()) < minEt_))
0420         continue;
0421 
0422       float SCenergy = moIter->energy();
0423       math::XYZPoint caloposition = moIter->position();
0424       float theta = caloposition.Theta();
0425       math::XYZVector momentum;
0426       float SCmomentumX = SCenergy * sin(theta) * cos(moIter->phi());
0427       float SCmomentumY = SCenergy * sin(theta) * sin(moIter->phi());
0428       float SCmomentumZ = SCenergy * cos(theta);
0429       const reco::Particle::LorentzVector pSCprobecandidate(SCmomentumX, SCmomentumY, SCmomentumZ, SCenergy);
0430 
0431       math::XYZTLorentzVector p12 = (*gsfIter).p4() + pSCprobecandidate;
0432       float mee2 = p12.Dot(p12);
0433       float invMass = mee2 > 0. ? sqrt(mee2) : 0.;
0434 
0435       if (invMass < massLow_ || invMass > massHigh_)
0436         continue;
0437 
0438       h1_matchingObject_Eta->Fill(moIter->eta());
0439       h1_matchingObject_Pt->Fill(moIter->energy() / cosh(moIter->eta()));
0440       h1_matchingObject_Phi->Fill(moIter->phi());
0441       //h1_matchingObject_Z->Fill(  moIter->z() );
0442 
0443       reco::GsfElectron bestGsfElectron;
0444       reco::SuperClusterRef sclProbeRef;
0445       bool okGsfFound = false;
0446       //double gsfOkRatio = 999999.;
0447       reco::GsfElectronCollection::const_iterator gsfIter2;
0448 
0449       for (gsfIter2 = gsfElectrons->begin(); gsfIter2 != gsfElectrons->end(); gsfIter2++) {
0450         // matching with ref
0451         sclProbeRef = gsfIter2->superCluster();
0452 
0453         if (sclProbeRef->eta() == moIter->eta()) {
0454           //std::cout << "les deux ref SC sont egales : " << std::endl;
0455           bestGsfElectron = *gsfIter2;
0456 
0457           // opposite sign checking
0458           bool opsign = (((gsfIter->charge()) * (bestGsfElectron.charge())) < 0.);
0459           if (TPchecksign_ && !opsign) {
0460             okGsfFound = false;
0461             h1_mee->Fill(invMass);
0462             break;
0463           } else {
0464             okGsfFound = true;
0465           }
0466         }  //fi on gsfEleSC.eta == probeSC.eta
0467 
0468         //        // matching with cone
0469         //        double dphi = gsfIter2->phi()-moIter->phi();
0470         //        if (std::abs(dphi)>CLHEP::pi)
0471         //         dphi = dphi < 0? (CLHEP::twopi) + dphi : dphi - CLHEP::twopi;
0472         //        double deltaR = sqrt(pow((moIter->eta()-gsfIter2->eta()),2) + pow(dphi,2));
0473         //        if ( deltaR < deltaR_ )
0474         //        {
0475         //         double tmpGsfRatio = gsfIter2->p()/moIter->energy();
0476         //         if ( std::abs(tmpGsfRatio-1) < std::abs(gsfOkRatio-1))
0477         //          {
0478         //           gsfOkRatio = tmpGsfRatio;
0479         //           bestGsfElectron=*gsfIter2;
0480         //           okGsfFound = true;
0481         //          }
0482         //        } // fi on deltaR
0483 
0484       }  // end of loop on gsfEle to find the best one which matches with probe SC
0485 
0486       if (okGsfFound) {
0487         // fill matched histos for eff
0488         fillMatchedHistos(moIter, bestGsfElectron);
0489 
0490         // basic quantities
0491         if (bestGsfElectron.isEB())
0492           h1_vertexPt_barrel->Fill(bestGsfElectron.pt());
0493         if (bestGsfElectron.isEE())
0494           h1_vertexPt_endcaps->Fill(bestGsfElectron.pt());
0495         h1_vertexEta->Fill(bestGsfElectron.eta());
0496         h2_vertexEtaVsPhi->Fill(bestGsfElectron.eta(), bestGsfElectron.phi());
0497         h2_vertexXvsY->Fill(bestGsfElectron.vertex().x(), bestGsfElectron.vertex().y());
0498         h1_vertexZ->Fill(bestGsfElectron.vertex().z());
0499 
0500         // supercluster related distributions
0501         reco::SuperClusterRef sclRef = bestGsfElectron.superCluster();
0502         double R = TMath::Sqrt(sclRef->x() * sclRef->x() + sclRef->y() * sclRef->y() + sclRef->z() * sclRef->z());
0503         double Rt = TMath::Sqrt(sclRef->x() * sclRef->x() + sclRef->y() * sclRef->y());
0504         h1_sclEt->Fill(sclRef->energy() * (Rt / R));
0505 
0506         if (!readAOD_) {  // track extra does not exist in AOD
0507           h1_foundHits->Fill(bestGsfElectron.gsfTrack()->numberOfValidHits());
0508           h1_lostHits->Fill(bestGsfElectron.gsfTrack()->numberOfLostHits());
0509           h1_chi2->Fill(bestGsfElectron.gsfTrack()->normalizedChi2());
0510         }
0511 
0512         // match distributions
0513         if (bestGsfElectron.isEB()) {
0514           h1_Eop_barrel->Fill(bestGsfElectron.eSuperClusterOverP());
0515           h1_EeleOPout_barrel->Fill(bestGsfElectron.eEleClusterOverPout());
0516           h1_dEtaSc_propVtx_barrel->Fill(bestGsfElectron.deltaEtaSuperClusterTrackAtVtx());
0517           h1_dEtaEleCl_propOut_barrel->Fill(bestGsfElectron.deltaEtaEleClusterTrackAtCalo());
0518           h1_dPhiSc_propVtx_barrel->Fill(bestGsfElectron.deltaPhiSuperClusterTrackAtVtx());
0519           h1_dPhiEleCl_propOut_barrel->Fill(bestGsfElectron.deltaPhiEleClusterTrackAtCalo());
0520           h1_Hoe_barrel->Fill(bestGsfElectron.hadronicOverEm());
0521           h1_sclSigEtaEta_barrel->Fill(bestGsfElectron.scSigmaEtaEta());
0522         }
0523         if (bestGsfElectron.isEE()) {
0524           h1_Eop_endcaps->Fill(bestGsfElectron.eSuperClusterOverP());
0525           h1_EeleOPout_endcaps->Fill(bestGsfElectron.eEleClusterOverPout());
0526           h1_dEtaSc_propVtx_endcaps->Fill(bestGsfElectron.deltaEtaSuperClusterTrackAtVtx());
0527           h1_dEtaEleCl_propOut_endcaps->Fill(bestGsfElectron.deltaEtaEleClusterTrackAtCalo());
0528           h1_dPhiSc_propVtx_endcaps->Fill(bestGsfElectron.deltaPhiSuperClusterTrackAtVtx());
0529           h1_dPhiEleCl_propOut_endcaps->Fill(bestGsfElectron.deltaPhiEleClusterTrackAtCalo());
0530           h1_Hoe_endcaps->Fill(bestGsfElectron.hadronicOverEm());
0531           h1_sclSigEtaEta_endcaps->Fill(bestGsfElectron.scSigmaEtaEta());
0532         }
0533 
0534         // fbrem
0535         h1_fbrem->Fill(bestGsfElectron.fbrem());
0536         int eleClass = bestGsfElectron.classification();
0537         if (bestGsfElectron.isEE())
0538           eleClass += 5;
0539         h1_classes->Fill(eleClass);
0540 
0541         // pflow
0542         h1_mva->Fill(bestGsfElectron.mva_e_pi());
0543         if (bestGsfElectron.ecalDrivenSeed())
0544           h1_provenance->Fill(1.);
0545         if (bestGsfElectron.trackerDrivenSeed())
0546           h1_provenance->Fill(-1.);
0547         if (bestGsfElectron.trackerDrivenSeed() || bestGsfElectron.ecalDrivenSeed())
0548           h1_provenance->Fill(0.);
0549         if (bestGsfElectron.trackerDrivenSeed() && !bestGsfElectron.ecalDrivenSeed())
0550           h1_provenance->Fill(-2.);
0551         if (!bestGsfElectron.trackerDrivenSeed() && bestGsfElectron.ecalDrivenSeed())
0552           h1_provenance->Fill(2.);
0553 
0554         // isolation
0555         h1_tkSumPt_dr03->Fill(bestGsfElectron.dr03TkSumPt());
0556         h1_ecalRecHitSumEt_dr03->Fill(bestGsfElectron.dr03EcalRecHitSumEt());
0557         h1_hcalTowerSumEt_dr03->Fill(bestGsfElectron.dr03HcalTowerSumEt());
0558 
0559         // inv Mass with opposite sign
0560         bool invMassTTAlreadyFilled = false;
0561         if (!TTCheck.empty()) {
0562           int TTCheckDim = TTCheck.size();
0563           for (int i = 0; i < TTCheckDim; i++) {
0564             if ((bestGsfElectron.eta() == TTCheck.at(i).first) && (gsfIter->eta() == TTCheck.at(i).second)) {
0565               invMassTTAlreadyFilled = true;
0566             }
0567           }
0568         }
0569 
0570         if (!invMassTTAlreadyFilled && (((gsfIter->charge()) * (bestGsfElectron.charge())) < 0.)) {
0571           h1_mee->Fill(invMass);
0572           math::XYZTLorentzVector p12bis = (*gsfIter).p4() + bestGsfElectron.p4();
0573           float mee2bis = p12.Dot(p12bis);
0574           float invMassEE = mee2bis > 0. ? sqrt(mee2bis) : 0.;
0575           if (invMassEE >= massLow_ && invMassEE <= massHigh_) {
0576             h1_mee_os->Fill(invMassEE);
0577           }
0578           std::pair<double, double> p(gsfIter->eta(), bestGsfElectron.eta());
0579           TTCheck.push_back(p);
0580         }
0581 
0582       }  // fi on OkGsfFound
0583       else {
0584         h1_mee->Fill(invMass);
0585       }
0586 
0587     }  // end of loop on SC to find probe SC
0588 
0589   }  // end of loop on Tag gsfEle
0590 }
0591 
0592 float ElectronTagProbeAnalyzer::computeInvMass(const reco::GsfElectron& e1, const reco::GsfElectron& e2) {
0593   math::XYZTLorentzVector p12 = e1.p4() + e2.p4();
0594   float mee2 = p12.Dot(p12);
0595   float invMass = mee2 > 0. ? sqrt(mee2) : 0.;
0596   return invMass;
0597 }
0598 
0599 void ElectronTagProbeAnalyzer::fillMatchedHistos(const reco::SuperClusterCollection::const_iterator& moIter,
0600                                                  const reco::GsfElectron& electron) {
0601   // generated distributions for matched electrons
0602   h1_matchedObject_Eta->Fill(moIter->eta());
0603   //  h1_matchedObject_AbsEta->Fill( std::abs(moIter->eta()) );
0604   h1_matchedObject_Pt->Fill(moIter->energy() / cosh(moIter->eta()));
0605   h1_matchedObject_Phi->Fill(moIter->phi());
0606   //h1_matchedObject_Z->Fill( moIter->z() );
0607 
0608   //classes
0609   //  int eleClass = electron.classification() ;
0610   //  h_classes->Fill(eleClass) ;
0611   //  h_matchedEle_eta->Fill(std::abs(electron.eta()));
0612   //  if (electron.classification() == GsfElectron::GOLDEN) h_matchedEle_eta_golden->Fill(std::abs(electron.eta()));
0613   //  if (electron.classification() == GsfElectron::SHOWERING) h_matchedEle_eta_shower->Fill(std::abs(electron.eta()));
0614   //  //if (electron.classification() == GsfElectron::BIGBREM) h_matchedEle_eta_bbrem->Fill(std::abs(electron.eta()));
0615   //  //if (electron.classification() == GsfElectron::OLDNARROW) h_matchedEle_eta_narrow->Fill(std::abs(electron.eta()));
0616 }
0617 
0618 bool ElectronTagProbeAnalyzer::selected(const reco::GsfElectronCollection::const_iterator& gsfIter, double vertexTIP) {
0619   if ((Selection_ > 0) && generalCut(gsfIter))
0620     return false;
0621   if ((Selection_ >= 1) && etCut(gsfIter))
0622     return false;
0623   if ((Selection_ >= 2) && isolationCut(gsfIter, vertexTIP))
0624     return false;
0625   if ((Selection_ >= 3) && idCut(gsfIter))
0626     return false;
0627   return true;
0628 }
0629 
0630 bool ElectronTagProbeAnalyzer::generalCut(const reco::GsfElectronCollection::const_iterator& gsfIter) {
0631   if (std::abs(gsfIter->eta()) > maxAbsEta_)
0632     return true;
0633   if (gsfIter->pt() < minPt_)
0634     return true;
0635 
0636   if (gsfIter->isEB() && isEE_)
0637     return true;
0638   if (gsfIter->isEE() && isEB_)
0639     return true;
0640   if (gsfIter->isEBEEGap() && isNotEBEEGap_)
0641     return true;
0642 
0643   if (gsfIter->ecalDrivenSeed() && isTrackerDriven_)
0644     return true;
0645   if (gsfIter->trackerDrivenSeed() && isEcalDriven_)
0646     return true;
0647 
0648   return false;
0649 }
0650 
0651 bool ElectronTagProbeAnalyzer::etCut(const reco::GsfElectronCollection::const_iterator& gsfIter) {
0652   if (gsfIter->superCluster()->energy() / cosh(gsfIter->superCluster()->eta()) < minEt_)
0653     return true;
0654 
0655   return false;
0656 }
0657 
0658 bool ElectronTagProbeAnalyzer::isolationCut(const reco::GsfElectronCollection::const_iterator& gsfIter,
0659                                             double vertexTIP) {
0660   if (gsfIter->isEB() && vertexTIP > tipMaxBarrel_)
0661     return true;
0662   if (gsfIter->isEE() && vertexTIP > tipMaxEndcaps_)
0663     return true;
0664 
0665   if (gsfIter->dr03TkSumPt() > tkIso03Max_)
0666     return true;
0667   if (gsfIter->isEB() && gsfIter->dr03HcalTowerSumEt(1) > hcalIso03Depth1MaxBarrel_)
0668     return true;
0669   if (gsfIter->isEE() && gsfIter->dr03HcalTowerSumEt(1) > hcalIso03Depth1MaxEndcaps_)
0670     return true;
0671   if (gsfIter->isEE() && gsfIter->dr03HcalTowerSumEt(2) > hcalIso03Depth2MaxEndcaps_)
0672     return true;
0673   if (gsfIter->isEB() && gsfIter->dr03EcalRecHitSumEt() > ecalIso03MaxBarrel_)
0674     return true;
0675   if (gsfIter->isEE() && gsfIter->dr03EcalRecHitSumEt() > ecalIso03MaxEndcaps_)
0676     return true;
0677 
0678   return false;
0679 }
0680 
0681 bool ElectronTagProbeAnalyzer::idCut(const reco::GsfElectronCollection::const_iterator& gsfIter) {
0682   if (gsfIter->isEB() && gsfIter->eSuperClusterOverP() < eOverPMinBarrel_)
0683     return true;
0684   if (gsfIter->isEB() && gsfIter->eSuperClusterOverP() > eOverPMaxBarrel_)
0685     return true;
0686   if (gsfIter->isEE() && gsfIter->eSuperClusterOverP() < eOverPMinEndcaps_)
0687     return true;
0688   if (gsfIter->isEE() && gsfIter->eSuperClusterOverP() > eOverPMaxEndcaps_)
0689     return true;
0690   if (gsfIter->isEB() && std::abs(gsfIter->deltaEtaSuperClusterTrackAtVtx()) < dEtaMinBarrel_)
0691     return true;
0692   if (gsfIter->isEB() && std::abs(gsfIter->deltaEtaSuperClusterTrackAtVtx()) > dEtaMaxBarrel_)
0693     return true;
0694   if (gsfIter->isEE() && std::abs(gsfIter->deltaEtaSuperClusterTrackAtVtx()) < dEtaMinEndcaps_)
0695     return true;
0696   if (gsfIter->isEE() && std::abs(gsfIter->deltaEtaSuperClusterTrackAtVtx()) > dEtaMaxEndcaps_)
0697     return true;
0698   if (gsfIter->isEB() && std::abs(gsfIter->deltaPhiSuperClusterTrackAtVtx()) < dPhiMinBarrel_)
0699     return true;
0700   if (gsfIter->isEB() && std::abs(gsfIter->deltaPhiSuperClusterTrackAtVtx()) > dPhiMaxBarrel_)
0701     return true;
0702   if (gsfIter->isEE() && std::abs(gsfIter->deltaPhiSuperClusterTrackAtVtx()) < dPhiMinEndcaps_)
0703     return true;
0704   if (gsfIter->isEE() && std::abs(gsfIter->deltaPhiSuperClusterTrackAtVtx()) > dPhiMaxEndcaps_)
0705     return true;
0706   if (gsfIter->isEB() && gsfIter->scSigmaIEtaIEta() < sigIetaIetaMinBarrel_)
0707     return true;
0708   if (gsfIter->isEB() && gsfIter->scSigmaIEtaIEta() > sigIetaIetaMaxBarrel_)
0709     return true;
0710   if (gsfIter->isEE() && gsfIter->scSigmaIEtaIEta() < sigIetaIetaMinEndcaps_)
0711     return true;
0712   if (gsfIter->isEE() && gsfIter->scSigmaIEtaIEta() > sigIetaIetaMaxEndcaps_)
0713     return true;
0714   if (gsfIter->isEB() && gsfIter->hadronicOverEm() > hadronicOverEmMaxBarrel_)
0715     return true;
0716   if (gsfIter->isEE() && gsfIter->hadronicOverEm() > hadronicOverEmMaxEndcaps_)
0717     return true;
0718 
0719   return false;
0720 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!