Utils/src/CheckEventContent.cc

// -*- C++ -*-
//
// Package:    CheckEventContent
// Class:      CheckEventContent
// 
/**\class CheckEventContent CheckEventContent.h LJMet/Utils/interface/CheckEventContent.h

 Description: times various aspects of data processing

 Implementation:
     Keep this class as simple and fast as possible
*/
//
// Original Author:  "Gennadiy Kukartsev"
//         Created:  Tue Oct  6 13:39:12 CDT 2009
// $Id: CheckEventContent.cc,v 1.1 2009/10/27 23:24:12 kukartse Exp $
//
//

#include "LJMet/Utils/interface/CheckEventContent.h"
#include "FWCore/Framework/interface/TriggerNames.h"
#include "DataFormats/Common/interface/TriggerResults.h"
#include "DataFormats/PatCandidates/interface/Jet.h"
#include "DataFormats/PatCandidates/interface/Electron.h"
#include "DataFormats/PatCandidates/interface/Muon.h"
#include "DataFormats/PatCandidates/interface/MET.h"
#include "DataFormats/PatCandidates/interface/TriggerPath.h"
#include "DataFormats/BeamSpot/interface/BeamSpot.h"
#include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
#include "DataFormats/Math/interface/deltaR.h"

using namespace edm;

CheckEventContent::CheckEventContent(const edm::ParameterSet& iConfig){
  std::cout << "=======>: CheckEventContent::CheckEventContent()" << std::endl;
}


CheckEventContent::~CheckEventContent()
{
  std::cout << "=======>: CheckEventContent::~CheckEventContent()" << std::endl;
}


void
CheckEventContent::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
{
  std::cout << "=======>: CheckEventContent::analyze()" << std::endl;
  //
  //_____ event book keeping
  //
  unsigned int irun   = (unsigned int)iEvent.id().run();
  // this will work starting with CMSSW_340:
  //unsigned int ilumi  = (unsigned int)iEvent.id().luminosityBlock();
  unsigned int ilumi  = (unsigned int)iEvent.getLuminosityBlock().luminosityBlock();
  unsigned int ievent = (unsigned int)iEvent.id().event();
  cout << std::endl << "===> Provenance: " << std::endl;
  cout << "Run, lumi, event: " << irun << ", " << ilumi << ", "<< ievent << std::endl;
  //
  //_____ check trigger
  //
  cout << std::endl << "===> Check trigger: " << std::endl;
  Handle<std::vector<pat::TriggerPath> > trigs;
  iEvent.getByLabel("patTrigger", trigs);
  std::string triggerName = "HLT_Mu9";
  for (std::vector<pat::TriggerPath>::const_iterator trig = trigs->begin();
       trig!=trigs->end();
       ++trig){
    //cout << "Trigger name and number: " << trig->name() << ",  " << trig->index() << std::endl;
    if( trig->name().find(triggerName)!=std::string::npos ){
      cout << "Trigger " << triggerName << " is found in the event..." << std::endl;
      //trig->wasAccept();
    }
  }
  //
  //_____ check PAT muons _______________________________________________
  //
  cout << std::endl << "===> Check PAT muons: " << std::endl;
  Handle< vector< pat::Muon > > muons;
  iEvent . getByLabel( "selectedLayer1Muons", muons );
  int muon_coll_size = muons->size();
  cout << "muon collection size: "<< muon_coll_size << endl;
  if (muon_coll_size>0){
    cout << "Global muon selector: " << (*muons)[0].isGlobalMuon() << std::endl;
    double _pt = (*muons)[0].pt();
    double _eta = (*muons)[0].eta();
    double _d0 = (*muons)[0].dB(); // global track
    double _chi2ndof = (*muons)[0].normChi2();
    int _nhits = (*muons)[0].numberOfValidHits();
    cout << "pT= " << _pt << ", eta=" << _eta << std::endl;
    cout << "PAT d0= " << _d0 << ", PAT norm chi2=" << _chi2ndof << std::endl;
    cout << "PAT n of hits on track= " << _nhits << std::endl;
    //
    //_____ check that there is a track
    reco::TrackRef _track = (*muons)[0].track();
    reco::TrackRef _innertrack = (*muons)[0].innerTrack();
    if(_track.isNull()){
      cout << "there is no pat::muon::track()" << endl;
    }
    else{
      cout << "there is a pat::muon::track()" << endl;
    }
    if(_innertrack.isNull()){
      cout << "there is no pat::muon::innerTrack()" << endl;
    }
    else{
      cout << "there is a pat::muon::innerTrack()" << endl;
    }
    //
    //_____ check isolation
    if ( !((*muons)[0].isIsolationValid()) ) {
      cout << "Isolation info is invalid" << endl;
    }
    else{
      cout << "Isolation info is valid..." << endl;
      double _emveto = (*muons)[0].isolationR03().emVetoEt;
      double _hadveto = (*muons)[0].isolationR03().hadVetoEt;
      double _trackiso = (*muons)[0].isolationR03().sumPt;
      double _emcaloiso = (*muons)[0].isolationR03().emEt;
      double _hadcaloiso = (*muons)[0].isolationR03().hadEt;
      cout << "EM veto= " << _emveto << ", Had veto=" << _hadveto << std::endl;
      cout << "Track iso= " << _trackiso << std::endl;
      cout << "EM calo iso= " << _emcaloiso << ", Had calo iso=" << _hadcaloiso << std::endl;
    }
  }
  //
  //_____ check Beam spot _____________________________________________________
  //
  cout << std::endl << "===> Check Beam spot: " << std::endl;
  Handle< reco::BeamSpot >          beamSpotHandle;
  iEvent . getByLabel( "offlineBeamSpot", beamSpotHandle);
  reco::BeamSpot beamSpot;
  if ( beamSpotHandle.isValid() ){
    beamSpot = *beamSpotHandle;
    cout << "beam spot present in the Event..." << std::endl;
  }
  else{
    edm::LogInfo("TtLJetsAnalyzer")
      << "No beam spot available in the event \n";
  }
  cout << beamSpot.x0() << "   " << beamSpot.y0() << "   " << beamSpot.z0() << std::endl;
  //
  //_____ check PAT electrons _______________________________________________
  //
  std::cout << std::endl << "===> Check PAT electrons: " << std::endl;
  Handle< vector< pat::Electron > > electrons;
  iEvent . getByLabel( "selectedLayer1Electrons", electrons );
  int electron_coll_size = electrons->size();
  cout << "electron collection size: "<< electron_coll_size << endl;
  if (electron_coll_size>0){
    cout << "Robust tight electron selector: " << (*electrons)[0].electronID("eidRobustTight") << std::endl;
    double _pt = (*electrons)[0].pt();
    double _eta = (*electrons)[0].eta();
    double _d0 = (*electrons)[0].dB();
    cout << "pT= " << _pt << ", eta=" << _eta << std::endl;
    cout << "PAT d0= " << _d0 << std::endl;
    //
    //_____ check that there is a track
    /*
    reco::TrackRef _track = (*electrons)[0].track();
    if(_track.isNull()){
      cout << "there is no pat::electron::track()" << endl;
    }
    else{
      cout << "there is a pat::electron::track()" << endl;
    }
    */
    reco::GsfTrackRef _gsfTrack = (*electrons)[0].gsfTrack();
    if(_gsfTrack.isNull()){
      cout << "there is no pat::electron::gsfTrack()" << endl;
    }
    else{
      cout << "there is a pat::electron::gsfTrack()" << endl;
    }
    //
    //_____ check isolation
      cout << "Isolation info is valid..." << endl;
      double _trackiso = (*electrons)[0].trackIso();
      double _emcaloiso = (*electrons)[0].ecalIso();
      double _hadcaloiso = (*electrons)[0].hcalIso();
      double _caloiso = (*electrons)[0].caloIso();
      cout << "Track iso= " << _trackiso << std::endl;
      cout << "EM calo iso= " << _emcaloiso << ", Had calo iso=" << _hadcaloiso << std::endl;
      cout << "Calo iso= " << _caloiso << std::endl;
  }
  //
  //_____ check PAT jets _______________________________________________
  //
  std::cout << std::endl << "===> Check PAT jets: " << std::endl;
  Handle< vector< pat::Jet > > jets;
  iEvent . getByLabel( "selectedLayer1Jets", jets );
  int jet_coll_size = jets->size();
  cout << "jet collection size: "<< jet_coll_size << endl;
  if (jet_coll_size>0){
    double _pt = (*jets)[0].pt();
    double _eta = (*jets)[0].eta();
    //float _bDiscriminator = (*jets)[0].bDiscriminator("");
    const std::vector<std::pair<std::string, float> > & btag = (*jets)[0].getPairDiscri();
    cout << "B-tagging discriminators: " << endl;
    for (std::vector<std::pair<std::string, float> >::const_iterator _d = btag.begin();
         _d != btag.end();
         ++_d){
      cout << _d->first << " = " << _d->second << endl;
    }
    cout << "pT= " << _pt << ", eta=" << _eta << std::endl;
  }
}


void 
CheckEventContent::beginJob()
{
  std::cout << "=======>: CheckEventContent::beginJob()" << std::endl;
}


void 
CheckEventContent::endJob() {
  std::cout << "=======>: CheckEventContent::endJob()" << std::endl;
}


/*
bool CheckEventContent::is_muon(std::string mLabel,
                      double mPt,
                      double mEta, 
                      int mNHits,
                      double mD0,
                      double mChi2Ndof,
                      double mEmVeto,
                      double mHadVeto,
                      double mRelIso,
                      const edm::Event& iEvent){
  bool result = false;
  //cout << "Cut Flow: passed stage 1 (mu+jets) or stage 2 (e+jets): ";
  Handle< vector< pat::Muon > >     muons;
  iEvent . getByLabel( "selectedLayer1Muons", muons );
  //
  //_____ loop over muons _______________________________________________
  //
  double _pt = 0.0;
  double _abseta = 0.0;
  int _nhits = 0;
  double _d0 = 0.0;
  double _chi2ndof = 0.0;
  double _emveto = 10000.0;
  double _hadveto = 10000.0;
  double _trackiso = 10000.0;
  double _caloiso = 10000.0;
  double _reliso = 10000.0;
  MeanCounter passed_muons("");
  for ( vector<pat::Muon>::const_iterator mu = muons -> begin(); mu != muons -> end(); mu++){
    _reliso = (_trackiso + _caloiso)/_pt;
    if (_pt > mPt &&
        _abseta < mEta &&
        _nhits >= mNHits &&
        fabs(_d0) < mD0 &&
        _chi2ndof < mChi2Ndof &&
        fabs(_emveto) < mEmVeto &&
        fabs(_hadveto) < mHadVeto &&
        fabs(_reliso) < mRelIso
        ){
      passed_muons.count();
      muon_index = (int)(mu - muons->begin());
      if (passed_muons.getCount()==1){
        cout << _pt << ", ";
        cout << _abseta << ", ";
        cout << _nhits << ", ";
        cout << _d0 << ", ";
        cout << _chi2ndof << ", ";
        cout << _emveto << ", ";
        cout << _hadveto << ", ";
        cout << _trackiso << ", ";
        cout << _caloiso << ", ";
        cout << _reliso;
        //cout << std::endl;
      }
    }
  }
  cout << "Number of good muons: " << passed_muons.getCount();
  result = passed_muons.getCount()==1;
  return result;
}
  

bool CheckEventContent::is_jets(std::string mLabel,
                      double mPt,
                      double mEta,
                      const edm::Event& iEvent){
  bool result = false;

  //cout << "Cut Flow: passed stage 2: ";

  Handle< vector< pat::Jet > >      jets;
  iEvent . getByLabel( mLabel, jets );
  MeanCounter c_jets("");
  MeanCounter passed_jets("");
  for ( vector<pat::Jet>::const_iterator jet = jets -> begin(); jet != jets -> end(); jet++ ){
    c_jets.count();
    double _pt = jet->pt();
    double _eta = jet->eta();
    if (
        ( fabs( _eta ) < mEta ) &&
        ( _pt > mPt )
        )
      {
        passed_jets.count();
      }
    if (c_jets.getCount()<=4){
      cout << _pt << ", ";
      cout << _eta;
    }
  }
  //cout << std::endl;
  result = passed_jets.getCount() >=4;
  return result;
}


bool CheckEventContent::no_loose_muon(std::string mLabel,
                            double mEta,
                            double mPt,
                            double mRelIso,
                            const edm::Event& iEvent){
  bool result = false;
  //cout << "Cut Flow: passed stage 3: ";
  Handle< vector< pat::Muon > >     muons;
  iEvent . getByLabel( "selectedLayer1Muons", muons );
  //
  //_____ loop over muons _______________________________________________
  //
  double _pt = 0.0;
  double _abseta = 0.0;
  double _trackiso = 10000.0;
  double _caloiso = 10000.0;
  double _reliso = 10000.0;
  MeanCounter passed_muons("");
  for ( vector<pat::Muon>::const_iterator mu = muons -> begin(); mu != muons -> end(); mu++){
    //_____ check that the muon is global
    if( !(mu->isGlobalMuon()) ) continue;
    _pt = mu->pt();
    _abseta = fabs(mu->eta());
    //_____ check that isolation is valid
    if ( !(mu->isIsolationValid()) ) continue;
    _trackiso = mu->isolationR03().sumPt;
    _caloiso = mu->isolationR03().emEt + mu->isolationR03().hadEt;
    _reliso = (_trackiso + _caloiso)/_pt;
    int _index = (int)(mu - muons->begin());
    if (_pt > mPt &&
        _abseta < mEta &&
        _reliso < mRelIso &&
        _index != muon_index
        ){
      passed_muons.count();
    }
  }
  result = passed_muons.getCount()==0;
  //cout << std::endl;
  return result;
}


bool CheckEventContent::no_loose_electron(std::string mLabel,
                                double mEta,
                                double mPt,
                                double mRelIso,
                                const edm::Event& iEvent){
  bool result = false;
  //cout << "Cut Flow: passed stage 4: ";
  Handle< vector< pat::Electron > >     electrons;
  iEvent . getByLabel( "selectedLayer1Electrons", electrons );
  //
  //_____ loop over electrons _______________________________________________
  //
  double _pt = 0.0;
  double _abseta = 0.0;
  double _trackiso = 10000.0;
  double _caloiso = 10000.0;
  double _reliso = 10000.0;
  MeanCounter passed_electrons("");
  for ( vector<pat::Electron>::const_iterator e = electrons -> begin(); e != electrons -> end(); e++){
    _pt = e->pt();
    _abseta = fabs(e->eta());
    _trackiso = e->dr03IsolationVariables().tkSumPt;
    _caloiso = e->dr03IsolationVariables().ecalRecHitSumEt + e->dr03IsolationVariables().hcalDepth1TowerSumEt + e->dr03IsolationVariables().hcalDepth2TowerSumEt;
    _reliso = (_trackiso + _caloiso)/_pt;
    int _index = (int)(e - electrons->begin());
    if (_pt > mPt &&
        _abseta < mEta &&
        _reliso < mRelIso &&
        _index != electron_index
        ){
      passed_electrons.count();
    }
  }
  result = passed_electrons.getCount()==0;
  //cout << std::endl;
  return result;
}


bool CheckEventContent::has_electron(std::string mLabel,
                           const edm::Event& iEvent){
  bool result = false;
  Handle< vector< pat::Electron > >     electrons;
  iEvent . getByLabel( mLabel, electrons );
  result = electrons->size()>0;
  return result;
}


bool CheckEventContent::is_electron(std::string mLabel,
                      double mPt,
                      double mEta, 
                      double mD0,
                      double mRelIso,
                      const edm::Event& iEvent){
  bool result = false;
  //cout << "Cut Flow: passed stage 1: ";
  Handle< vector< pat::Electron > >     electrons;
  iEvent . getByLabel( "selectedLayer1Electrons", electrons );
  Handle< reco::BeamSpot >          beamSpotHandle;
  iEvent . getByLabel( "offlineBeamSpot", beamSpotHandle);
  //
  //_____ Beam spot _____________________________________________________
  //
  reco::BeamSpot beamSpot;
  if ( beamSpotHandle.isValid() ){
    beamSpot = *beamSpotHandle;
  }
  else{
    edm::LogInfo("TtLJetsAnalyzer")
      << "No beam spot available from EventSetup \n";
  }
  const reco::BeamSpot::Point _bs(beamSpot.x0(), beamSpot.y0(),beamSpot.z0());
  //
  //_____ loop over electrons _______________________________________________
  //
  double _pt = 0.0;
  double _abseta = 0.0;
  double _d0 = 0.0;
  double _d0_pat = 0.0;
  double _trackiso = 10000.0;
  double _ecaloiso = 10000.0;
  double _ecaloiso2 = 10000.0;
  double _pat_ecaloiso = 10000.0;
  double _hcaloiso1 = 10000.0;
  double _hcaloiso2 = 10000.0;
  double _caloiso = 10000.0;
  double _reliso = 10000.0;
  MeanCounter passed_electrons("");
  for ( vector<pat::Electron>::const_iterator e = electrons -> begin(); e != electrons -> end(); e++){
    //_____ check that the electron is robust-tight
    if( !(e->electronID("eidRobustTight")) ) continue;
    _pt = e->pt();
    _abseta = fabs(e->eta());
    _d0_pat = e->dB();
    _d0 = -(e->gsfTrack()->dxy(_bs));
    _trackiso = e->dr03IsolationVariables().tkSumPt;
    _ecaloiso = e->dr04IsolationVariables().ecalRecHitSumEt;
    _ecaloiso2 = e->dr03IsolationVariables().ecalRecHitSumEt;
    _pat_ecaloiso = e->isolation(pat::ECalIso);
    _hcaloiso1 = e->dr04IsolationVariables().hcalDepth1TowerSumEt;
    _hcaloiso2 = e->dr04IsolationVariables().hcalDepth2TowerSumEt;
    _caloiso = _ecaloiso + _hcaloiso1 + _hcaloiso2;
    //_caloiso = _pat_ecaloiso +_hcaloiso1 + _hcaloiso2;
    _reliso = (_trackiso + _caloiso)/_pt;
    if (_pt > mPt &&
        _abseta < mEta &&
        fabs(_d0) < mD0 &&
        fabs(_reliso) < mRelIso
        ){
      passed_electrons.count();
      electron_index = (int)(e - electrons->begin());
      if (passed_electrons.getCount()==1){
        cout << "pt=" << _pt << ", abseta=";
        cout << _abseta << ", d0=";
        cout << _d0 << ", trackiso=";
        cout << _trackiso << ", ecaloiso04=";
        cout << _ecaloiso << ", ecaloiso03=";
        cout << _ecaloiso2 << ", pat_ecaloiso=";
        cout << _pat_ecaloiso << ", hcaloiso1=";
        cout << _hcaloiso1 << ", hcaloiso2=";
        cout << _hcaloiso2 << ", caloiso04=";
        cout << _caloiso << ", reliso=";
        cout << _reliso;
        //cout << std::endl;
      }
    }
  }
  //cout << "Number of good electrons: " << passed_electrons.getCount();
  result = passed_electrons.getCount()==1;
  return result;
}
  

int CheckEventContent::count_jets(std::string jLabel,
                        double jPt,
                        double jEta,
                        std::string lLabel,
                        double lPt,
                        double lEta,
                        double lD0,
                        double lDR,
                        const edm::Event & iEvent){
  MeanCounter _njets("");
  Handle< vector< pat::Jet > >     jets;
  iEvent . getByLabel( jLabel, jets );
  Handle< vector< pat::Electron > >     electrons;
  iEvent . getByLabel( lLabel, electrons );
  Handle< reco::BeamSpot >          beamSpotHandle;
  iEvent . getByLabel( "offlineBeamSpot", beamSpotHandle);
  //
  //_____ Beam spot _____________________________________________________
  //
  reco::BeamSpot beamSpot;
  if ( beamSpotHandle.isValid() ){
    beamSpot = *beamSpotHandle;
  }
  else{
    edm::LogInfo("TtLJetsAnalyzer")
      << "No beam spot available from EventSetup \n";
  }
  const reco::BeamSpot::Point _bs(beamSpot.x0(), beamSpot.y0(),beamSpot.z0());
  //
  //_____ loop over jets _______________________________________________
  //
  double j_pt = 0.0;
  double j_eta = 0.0;
  for ( vector<pat::Jet>::const_iterator jet = jets -> begin(); jet != jets -> end(); ++jet){
    j_pt = jet->pt();
    j_eta = jet->eta();
    //
    //_____ loop over electrons _______________________________________________
    //
    double e_pt = 0.0;
    double e_abseta = 0.0;
    double e_d0 = 0.0;
    double e_dr = 0.0;
    bool jet_removed = false;
    for ( vector<pat::Electron>::const_iterator e = electrons -> begin(); e != electrons -> end(); e++){
      //_____ check that the electron is robust-tight
      if( !(e->electronID("eidRobustTight")) ) continue;
      e_pt = e->pt();
      e_abseta = fabs(e->eta());
      e_d0 = -(e->gsfTrack()->dxy(_bs));
      e_dr = reco::deltaR(*jet, *e);
      if ( e_pt > lPt && e_abseta < lEta && fabs(e_d0)<lD0 && e_dr<lDR){
        jet_removed = true;
      }
    }
    if (j_pt>jPt && fabs(j_eta)<jEta && !jet_removed) _njets.count();
  }
  cout << "Number of jets after cleaning: " << _njets.getCount();
  return _njets.getCount();
}

*/

//define this as a plug-in
DEFINE_FWK_MODULE(CheckEventContent);
Revision:	1.2
Committed:	Wed Nov 4 14:12:22 2009 UTC (15 years, 6 months ago) by kukartse
Content type:	text/plain
Branch:	MAIN
Changes since 1.1:	+125 -24 lines
Log Message:	* empty log message *
#	User	Rev	Content
1	kukartse	1.1	// -- C++ --
2			//
3			// Package: CheckEventContent
4			// Class: CheckEventContent
5			//
6			/**\class CheckEventContent CheckEventContent.h LJMet/Utils/interface/CheckEventContent.h
7
8			Description: times various aspects of data processing
9
10			Implementation:
11			Keep this class as simple and fast as possible
12			*/
13			//
14			// Original Author: "Gennadiy Kukartsev"
15			// Created: Tue Oct 6 13:39:12 CDT 2009
16	kukartse	1.2	// $Id: CheckEventContent.cc,v 1.1 2009/10/27 23:24:12 kukartse Exp $
17	kukartse	1.1	//
18			//
19
20			#include "LJMet/Utils/interface/CheckEventContent.h"
21			#include "FWCore/Framework/interface/TriggerNames.h"
22			#include "DataFormats/Common/interface/TriggerResults.h"
23			#include "DataFormats/PatCandidates/interface/Jet.h"
24			#include "DataFormats/PatCandidates/interface/Electron.h"
25			#include "DataFormats/PatCandidates/interface/Muon.h"
26			#include "DataFormats/PatCandidates/interface/MET.h"
27			#include "DataFormats/PatCandidates/interface/TriggerPath.h"
28			#include "DataFormats/BeamSpot/interface/BeamSpot.h"
29			#include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
30			#include "DataFormats/Math/interface/deltaR.h"
31
32			using namespace edm;
33
34			CheckEventContent::CheckEventContent(const edm::ParameterSet& iConfig){
35	kukartse	1.2	std::cout << "=======>: CheckEventContent::CheckEventContent()" << std::endl;
36	kukartse	1.1	}
37
38
39			CheckEventContent::~CheckEventContent()
40			{
41	kukartse	1.2	std::cout << "=======>: CheckEventContent::~CheckEventContent()" << std::endl;
42	kukartse	1.1	}
43
44
45			void
46			CheckEventContent::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
47			{
48	kukartse	1.2	std::cout << "=======>: CheckEventContent::analyze()" << std::endl;
49			//
50			//_____ event book keeping
51			//
52			unsigned int irun = (unsigned int)iEvent.id().run();
53			// this will work starting with CMSSW_340:
54			//unsigned int ilumi = (unsigned int)iEvent.id().luminosityBlock();
55			unsigned int ilumi = (unsigned int)iEvent.getLuminosityBlock().luminosityBlock();
56			unsigned int ievent = (unsigned int)iEvent.id().event();
57			cout << std::endl << "===> Provenance: " << std::endl;
58			cout << "Run, lumi, event: " << irun << ", " << ilumi << ", "<< ievent << std::endl;
59	kukartse	1.1	//
60			//_____ check trigger
61			//
62			cout << std::endl << "===> Check trigger: " << std::endl;
63			Handle<std::vector<pat::TriggerPath> > trigs;
64			iEvent.getByLabel("patTrigger", trigs);
65			std::string triggerName = "HLT_Mu9";
66			for (std::vector<pat::TriggerPath>::const_iterator trig = trigs->begin();
67			trig!=trigs->end();
68			++trig){
69			//cout << "Trigger name and number: " << trig->name() << ", " << trig->index() << std::endl;
70			if( trig->name().find(triggerName)!=std::string::npos ){
71			cout << "Trigger " << triggerName << " is found in the event..." << std::endl;
72			//trig->wasAccept();
73			}
74			}
75			//
76			//_____ check PAT muons _______________________________________________
77			//
78			cout << std::endl << "===> Check PAT muons: " << std::endl;
79			Handle< vector< pat::Muon > > muons;
80			iEvent . getByLabel( "selectedLayer1Muons", muons );
81	kukartse	1.2	int muon_coll_size = muons->size();
82			cout << "muon collection size: "<< muon_coll_size << endl;
83			if (muon_coll_size>0){
84			cout << "Global muon selector: " << (*muons)[0].isGlobalMuon() << std::endl;
85			double _pt = (*muons)[0].pt();
86			double _eta = (*muons)[0].eta();
87			double _d0 = (*muons)[0].dB(); // global track
88			double _chi2ndof = (*muons)[0].normChi2();
89			int _nhits = (*muons)[0].numberOfValidHits();
90			cout << "pT= " << _pt << ", eta=" << _eta << std::endl;
91			cout << "PAT d0= " << _d0 << ", PAT norm chi2=" << _chi2ndof << std::endl;
92			cout << "PAT n of hits on track= " << _nhits << std::endl;
93			//
94			//_____ check that there is a track
95			reco::TrackRef _track = (*muons)[0].track();
96			reco::TrackRef _innertrack = (*muons)[0].innerTrack();
97			if(_track.isNull()){
98			cout << "there is no pat::muon::track()" << endl;
99			}
100			else{
101			cout << "there is a pat::muon::track()" << endl;
102			}
103			if(_innertrack.isNull()){
104			cout << "there is no pat::muon::innerTrack()" << endl;
105			}
106			else{
107			cout << "there is a pat::muon::innerTrack()" << endl;
108			}
109			//
110			//_____ check isolation
111			if ( !((*muons)[0].isIsolationValid()) ) {
112			cout << "Isolation info is invalid" << endl;
113			}
114			else{
115			cout << "Isolation info is valid..." << endl;
116			double _emveto = (*muons)[0].isolationR03().emVetoEt;
117			double _hadveto = (*muons)[0].isolationR03().hadVetoEt;
118			double _trackiso = (*muons)[0].isolationR03().sumPt;
119			double _emcaloiso = (*muons)[0].isolationR03().emEt;
120			double _hadcaloiso = (*muons)[0].isolationR03().hadEt;
121			cout << "EM veto= " << _emveto << ", Had veto=" << _hadveto << std::endl;
122			cout << "Track iso= " << _trackiso << std::endl;
123			cout << "EM calo iso= " << _emcaloiso << ", Had calo iso=" << _hadcaloiso << std::endl;
124			}
125	kukartse	1.1	}
126			//
127			//_____ check Beam spot _____________________________________________________
128			//
129			cout << std::endl << "===> Check Beam spot: " << std::endl;
130			Handle< reco::BeamSpot > beamSpotHandle;
131			iEvent . getByLabel( "offlineBeamSpot", beamSpotHandle);
132			reco::BeamSpot beamSpot;
133			if ( beamSpotHandle.isValid() ){
134			beamSpot = *beamSpotHandle;
135			cout << "beam spot present in the Event..." << std::endl;
136			}
137			else{
138			edm::LogInfo("TtLJetsAnalyzer")
139			<< "No beam spot available in the event \n";
140			}
141			cout << beamSpot.x0() << " " << beamSpot.y0() << " " << beamSpot.z0() << std::endl;
142	kukartse	1.2	//
143			//_____ check PAT electrons _______________________________________________
144			//
145			std::cout << std::endl << "===> Check PAT electrons: " << std::endl;
146			Handle< vector< pat::Electron > > electrons;
147			iEvent . getByLabel( "selectedLayer1Electrons", electrons );
148			int electron_coll_size = electrons->size();
149			cout << "electron collection size: "<< electron_coll_size << endl;
150			if (electron_coll_size>0){
151			cout << "Robust tight electron selector: " << (*electrons)[0].electronID("eidRobustTight") << std::endl;
152			double _pt = (*electrons)[0].pt();
153			double _eta = (*electrons)[0].eta();
154			double _d0 = (*electrons)[0].dB();
155			cout << "pT= " << _pt << ", eta=" << _eta << std::endl;
156			cout << "PAT d0= " << _d0 << std::endl;
157			//
158			//_____ check that there is a track
159			/*
160			reco::TrackRef _track = (*electrons)[0].track();
161			if(_track.isNull()){
162			cout << "there is no pat::electron::track()" << endl;
163			}
164			else{
165			cout << "there is a pat::electron::track()" << endl;
166			}
167			*/
168			reco::GsfTrackRef _gsfTrack = (*electrons)[0].gsfTrack();
169			if(_gsfTrack.isNull()){
170			cout << "there is no pat::electron::gsfTrack()" << endl;
171			}
172			else{
173			cout << "there is a pat::electron::gsfTrack()" << endl;
174			}
175			//
176			//_____ check isolation
177			cout << "Isolation info is valid..." << endl;
178			double _trackiso = (*electrons)[0].trackIso();
179			double _emcaloiso = (*electrons)[0].ecalIso();
180			double _hadcaloiso = (*electrons)[0].hcalIso();
181			double _caloiso = (*electrons)[0].caloIso();
182			cout << "Track iso= " << _trackiso << std::endl;
183			cout << "EM calo iso= " << _emcaloiso << ", Had calo iso=" << _hadcaloiso << std::endl;
184			cout << "Calo iso= " << _caloiso << std::endl;
185			}
186			//
187			//_____ check PAT jets _______________________________________________
188			//
189			std::cout << std::endl << "===> Check PAT jets: " << std::endl;
190			Handle< vector< pat::Jet > > jets;
191			iEvent . getByLabel( "selectedLayer1Jets", jets );
192			int jet_coll_size = jets->size();
193			cout << "jet collection size: "<< jet_coll_size << endl;
194			if (jet_coll_size>0){
195			double _pt = (*jets)[0].pt();
196			double _eta = (*jets)[0].eta();
197			//float _bDiscriminator = (*jets)[0].bDiscriminator("");
198			const std::vector<std::pair<std::string, float> > & btag = (*jets)[0].getPairDiscri();
199			cout << "B-tagging discriminators: " << endl;
200			for (std::vector<std::pair<std::string, float> >::const_iterator _d = btag.begin();
201			_d != btag.end();
202			++_d){
203			cout << _d->first << " = " << _d->second << endl;
204			}
205			cout << "pT= " << _pt << ", eta=" << _eta << std::endl;
206			}
207	kukartse	1.1	}
208
209
210			void
211			CheckEventContent::beginJob()
212			{
213	kukartse	1.2	std::cout << "=======>: CheckEventContent::beginJob()" << std::endl;
214	kukartse	1.1	}
215
216
217			void
218			CheckEventContent::endJob() {
219	kukartse	1.2	std::cout << "=======>: CheckEventContent::endJob()" << std::endl;
220	kukartse	1.1	}
221
222
223			/*
224			bool CheckEventContent::is_muon(std::string mLabel,
225			double mPt,
226			double mEta,
227			int mNHits,
228			double mD0,
229			double mChi2Ndof,
230			double mEmVeto,
231			double mHadVeto,
232			double mRelIso,
233			const edm::Event& iEvent){
234			bool result = false;
235			//cout << "Cut Flow: passed stage 1 (mu+jets) or stage 2 (e+jets): ";
236			Handle< vector< pat::Muon > > muons;
237			iEvent . getByLabel( "selectedLayer1Muons", muons );
238			//
239			//_____ loop over muons _______________________________________________
240			//
241			double _pt = 0.0;
242			double _abseta = 0.0;
243			int _nhits = 0;
244			double _d0 = 0.0;
245			double _chi2ndof = 0.0;
246			double _emveto = 10000.0;
247			double _hadveto = 10000.0;
248			double _trackiso = 10000.0;
249			double _caloiso = 10000.0;
250			double _reliso = 10000.0;
251			MeanCounter passed_muons("");
252			for ( vector<pat::Muon>::const_iterator mu = muons -> begin(); mu != muons -> end(); mu++){
253			_reliso = (_trackiso + _caloiso)/_pt;
254			if (_pt > mPt &&
255			_abseta < mEta &&
256			_nhits >= mNHits &&
257			fabs(_d0) < mD0 &&
258			_chi2ndof < mChi2Ndof &&
259			fabs(_emveto) < mEmVeto &&
260			fabs(_hadveto) < mHadVeto &&
261			fabs(_reliso) < mRelIso
262			){
263			passed_muons.count();
264			muon_index = (int)(mu - muons->begin());
265			if (passed_muons.getCount()==1){
266			cout << _pt << ", ";
267			cout << _abseta << ", ";
268			cout << _nhits << ", ";
269			cout << _d0 << ", ";
270			cout << _chi2ndof << ", ";
271			cout << _emveto << ", ";
272			cout << _hadveto << ", ";
273			cout << _trackiso << ", ";
274			cout << _caloiso << ", ";
275			cout << _reliso;
276			//cout << std::endl;
277			}
278			}
279			}
280			cout << "Number of good muons: " << passed_muons.getCount();
281			result = passed_muons.getCount()==1;
282			return result;
283			}
284
285
286			bool CheckEventContent::is_jets(std::string mLabel,
287			double mPt,
288			double mEta,
289			const edm::Event& iEvent){
290			bool result = false;
291
292			//cout << "Cut Flow: passed stage 2: ";
293
294			Handle< vector< pat::Jet > > jets;
295			iEvent . getByLabel( mLabel, jets );
296			MeanCounter c_jets("");
297			MeanCounter passed_jets("");
298			for ( vector<pat::Jet>::const_iterator jet = jets -> begin(); jet != jets -> end(); jet++ ){
299			c_jets.count();
300			double _pt = jet->pt();
301			double _eta = jet->eta();
302			if (
303			( fabs( _eta ) < mEta ) &&
304			( _pt > mPt )
305			)
306			{
307			passed_jets.count();
308			}
309			if (c_jets.getCount()<=4){
310			cout << _pt << ", ";
311			cout << _eta;
312			}
313			}
314			//cout << std::endl;
315			result = passed_jets.getCount() >=4;
316			return result;
317			}
318
319
320			bool CheckEventContent::no_loose_muon(std::string mLabel,
321			double mEta,
322			double mPt,
323			double mRelIso,
324			const edm::Event& iEvent){
325			bool result = false;
326			//cout << "Cut Flow: passed stage 3: ";
327			Handle< vector< pat::Muon > > muons;
328			iEvent . getByLabel( "selectedLayer1Muons", muons );
329			//
330			//_____ loop over muons _______________________________________________
331			//
332			double _pt = 0.0;
333			double _abseta = 0.0;
334			double _trackiso = 10000.0;
335			double _caloiso = 10000.0;
336			double _reliso = 10000.0;
337			MeanCounter passed_muons("");
338			for ( vector<pat::Muon>::const_iterator mu = muons -> begin(); mu != muons -> end(); mu++){
339			//_____ check that the muon is global
340			if( !(mu->isGlobalMuon()) ) continue;
341			_pt = mu->pt();
342			_abseta = fabs(mu->eta());
343			//_____ check that isolation is valid
344			if ( !(mu->isIsolationValid()) ) continue;
345			_trackiso = mu->isolationR03().sumPt;
346			_caloiso = mu->isolationR03().emEt + mu->isolationR03().hadEt;
347			_reliso = (_trackiso + _caloiso)/_pt;
348			int _index = (int)(mu - muons->begin());
349			if (_pt > mPt &&
350			_abseta < mEta &&
351			_reliso < mRelIso &&
352			_index != muon_index
353			){
354			passed_muons.count();
355			}
356			}
357			result = passed_muons.getCount()==0;
358			//cout << std::endl;
359			return result;
360			}
361
362
363			bool CheckEventContent::no_loose_electron(std::string mLabel,
364			double mEta,
365			double mPt,
366			double mRelIso,
367			const edm::Event& iEvent){
368			bool result = false;
369			//cout << "Cut Flow: passed stage 4: ";
370			Handle< vector< pat::Electron > > electrons;
371			iEvent . getByLabel( "selectedLayer1Electrons", electrons );
372			//
373			//_____ loop over electrons _______________________________________________
374			//
375			double _pt = 0.0;
376			double _abseta = 0.0;
377			double _trackiso = 10000.0;
378			double _caloiso = 10000.0;
379			double _reliso = 10000.0;
380			MeanCounter passed_electrons("");
381			for ( vector<pat::Electron>::const_iterator e = electrons -> begin(); e != electrons -> end(); e++){
382			_pt = e->pt();
383			_abseta = fabs(e->eta());
384			_trackiso = e->dr03IsolationVariables().tkSumPt;
385			_caloiso = e->dr03IsolationVariables().ecalRecHitSumEt + e->dr03IsolationVariables().hcalDepth1TowerSumEt + e->dr03IsolationVariables().hcalDepth2TowerSumEt;
386			_reliso = (_trackiso + _caloiso)/_pt;
387			int _index = (int)(e - electrons->begin());
388			if (_pt > mPt &&
389			_abseta < mEta &&
390			_reliso < mRelIso &&
391			_index != electron_index
392			){
393			passed_electrons.count();
394			}
395			}
396			result = passed_electrons.getCount()==0;
397			//cout << std::endl;
398			return result;
399			}
400
401
402			bool CheckEventContent::has_electron(std::string mLabel,
403			const edm::Event& iEvent){
404			bool result = false;
405			Handle< vector< pat::Electron > > electrons;
406			iEvent . getByLabel( mLabel, electrons );
407			result = electrons->size()>0;
408			return result;
409			}
410
411
412			bool CheckEventContent::is_electron(std::string mLabel,
413			double mPt,
414			double mEta,
415			double mD0,
416			double mRelIso,
417			const edm::Event& iEvent){
418			bool result = false;
419			//cout << "Cut Flow: passed stage 1: ";
420			Handle< vector< pat::Electron > > electrons;
421			iEvent . getByLabel( "selectedLayer1Electrons", electrons );
422			Handle< reco::BeamSpot > beamSpotHandle;
423			iEvent . getByLabel( "offlineBeamSpot", beamSpotHandle);
424			//
425			//_____ Beam spot _____________________________________________________
426			//
427			reco::BeamSpot beamSpot;
428			if ( beamSpotHandle.isValid() ){
429			beamSpot = *beamSpotHandle;
430			}
431			else{
432			edm::LogInfo("TtLJetsAnalyzer")
433			<< "No beam spot available from EventSetup \n";
434			}
435			const reco::BeamSpot::Point _bs(beamSpot.x0(), beamSpot.y0(),beamSpot.z0());
436			//
437			//_____ loop over electrons _______________________________________________
438			//
439			double _pt = 0.0;
440			double _abseta = 0.0;
441			double _d0 = 0.0;
442			double _d0_pat = 0.0;
443			double _trackiso = 10000.0;
444			double _ecaloiso = 10000.0;
445			double _ecaloiso2 = 10000.0;
446			double _pat_ecaloiso = 10000.0;
447			double _hcaloiso1 = 10000.0;
448			double _hcaloiso2 = 10000.0;
449			double _caloiso = 10000.0;
450			double _reliso = 10000.0;
451			MeanCounter passed_electrons("");
452			for ( vector<pat::Electron>::const_iterator e = electrons -> begin(); e != electrons -> end(); e++){
453			//_____ check that the electron is robust-tight
454			if( !(e->electronID("eidRobustTight")) ) continue;
455			_pt = e->pt();
456			_abseta = fabs(e->eta());
457			_d0_pat = e->dB();
458			_d0 = -(e->gsfTrack()->dxy(_bs));
459			_trackiso = e->dr03IsolationVariables().tkSumPt;
460			_ecaloiso = e->dr04IsolationVariables().ecalRecHitSumEt;
461			_ecaloiso2 = e->dr03IsolationVariables().ecalRecHitSumEt;
462			_pat_ecaloiso = e->isolation(pat::ECalIso);
463			_hcaloiso1 = e->dr04IsolationVariables().hcalDepth1TowerSumEt;
464			_hcaloiso2 = e->dr04IsolationVariables().hcalDepth2TowerSumEt;
465			_caloiso = _ecaloiso + _hcaloiso1 + _hcaloiso2;
466			//_caloiso = _pat_ecaloiso +_hcaloiso1 + _hcaloiso2;
467			_reliso = (_trackiso + _caloiso)/_pt;
468			if (_pt > mPt &&
469			_abseta < mEta &&
470			fabs(_d0) < mD0 &&
471			fabs(_reliso) < mRelIso
472			){
473			passed_electrons.count();
474			electron_index = (int)(e - electrons->begin());
475			if (passed_electrons.getCount()==1){
476			cout << "pt=" << _pt << ", abseta=";
477			cout << _abseta << ", d0=";
478			cout << _d0 << ", trackiso=";
479			cout << _trackiso << ", ecaloiso04=";
480			cout << _ecaloiso << ", ecaloiso03=";
481			cout << _ecaloiso2 << ", pat_ecaloiso=";
482			cout << _pat_ecaloiso << ", hcaloiso1=";
483			cout << _hcaloiso1 << ", hcaloiso2=";
484			cout << _hcaloiso2 << ", caloiso04=";
485			cout << _caloiso << ", reliso=";
486			cout << _reliso;
487			//cout << std::endl;
488			}
489			}
490			}
491			//cout << "Number of good electrons: " << passed_electrons.getCount();
492			result = passed_electrons.getCount()==1;
493			return result;
494			}
495
496
497			int CheckEventContent::count_jets(std::string jLabel,
498			double jPt,
499			double jEta,
500			std::string lLabel,
501			double lPt,
502			double lEta,
503			double lD0,
504			double lDR,
505			const edm::Event & iEvent){
506			MeanCounter _njets("");
507			Handle< vector< pat::Jet > > jets;
508			iEvent . getByLabel( jLabel, jets );
509			Handle< vector< pat::Electron > > electrons;
510			iEvent . getByLabel( lLabel, electrons );
511			Handle< reco::BeamSpot > beamSpotHandle;
512			iEvent . getByLabel( "offlineBeamSpot", beamSpotHandle);
513			//
514			//_____ Beam spot _____________________________________________________
515			//
516			reco::BeamSpot beamSpot;
517			if ( beamSpotHandle.isValid() ){
518			beamSpot = *beamSpotHandle;
519			}
520			else{
521			edm::LogInfo("TtLJetsAnalyzer")
522			<< "No beam spot available from EventSetup \n";
523			}
524			const reco::BeamSpot::Point _bs(beamSpot.x0(), beamSpot.y0(),beamSpot.z0());
525			//
526			//_____ loop over jets _______________________________________________
527			//
528			double j_pt = 0.0;
529			double j_eta = 0.0;
530			for ( vector<pat::Jet>::const_iterator jet = jets -> begin(); jet != jets -> end(); ++jet){
531			j_pt = jet->pt();
532			j_eta = jet->eta();
533			//
534			//_____ loop over electrons _______________________________________________
535			//
536			double e_pt = 0.0;
537			double e_abseta = 0.0;
538			double e_d0 = 0.0;
539			double e_dr = 0.0;
540			bool jet_removed = false;
541			for ( vector<pat::Electron>::const_iterator e = electrons -> begin(); e != electrons -> end(); e++){
542			//_____ check that the electron is robust-tight
543			if( !(e->electronID("eidRobustTight")) ) continue;
544			e_pt = e->pt();
545			e_abseta = fabs(e->eta());
546			e_d0 = -(e->gsfTrack()->dxy(_bs));
547			e_dr = reco::deltaR(jet, e);
548			if ( e_pt > lPt && e_abseta < lEta && fabs(e_d0)<lD0 && e_dr<lDR){
549			jet_removed = true;
550			}
551			}
552			if (j_pt>jPt && fabs(j_eta)<jEta && !jet_removed) _njets.count();
553			}
554			cout << "Number of jets after cleaning: " << _njets.getCount();
555			return _njets.getCount();
556			}
557
558			*/
559
560			//define this as a plug-in
561			DEFINE_FWK_MODULE(CheckEventContent);