#include <PFTauElecRejectionBenchmark.h>

Public Member Functions
	PFTauElecRejectionBenchmark ()
void	process (edm::Handle< edm::HepMCProduct > mcevt, edm::Handle< reco::PFTauCollection > pfTaus, edm::Handle< reco::PFTauDiscriminator > pfTauIsoDiscr, edm::Handle< reco::PFTauDiscriminator > pfTauElecDiscr)
void	setup (std::string Filename, std::string benchmarkLabel, double maxDeltaR, double minRecoPt, double maxRecoAbsEta, double minMCPt, double maxMCAbsEta, std::string sGenMatchObjectLabel, bool applyEcalCrackCut, DQMStore *db_store)
void	write ()
virtual	~PFTauElecRejectionBenchmark ()
Protected Attributes
DQMStore *	db_
Private Member Functions
bool	isInEcalCrack (double eta) const
Private Attributes
std::vector< TLorentzVector >	_GenObjects
bool	applyEcalCrackCut_
std::string	benchmarkLabel_
TFile *	file_
TH1F *	hElecMVA
TH1F *	hElecPreID
TH1F *	hEmfrac
TH1F *	hEmfrac_barrel
TH1F *	hEmfrac_endcap
TH1F *	hEmfrac_preid0
TH1F *	hEmfrac_preid1
TH2F *	hEmfracvsEoP
TH2F *	hEmfracvsEoP_preid0
TH2F *	hEmfracvsEoP_preid1
TH1F *	hEoverP
TH1F *	hEoverP_barrel
TH1F *	hEoverP_endcap
TH1F *	hEoverP_preid0
TH1F *	hEoverP_preid1
TH2F *	hHoPvsEoP
TH2F *	hHoPvsEoP_preid0
TH2F *	hHoPvsEoP_preid1
TH1F *	hHoverP
TH1F *	hHoverP_barrel
TH1F *	hHoverP_endcap
TH1F *	hHoverP_preid0
TH1F *	hHoverP_preid1
TH1F *	hleadGsfTk_eta
TH1F *	hleadGsfTk_phi
TH1F *	hleadGsfTk_pt
TH1F *	hleadTk_eta
TH1F *	hleadTk_phi
TH1F *	hleadTk_pt
TH1F *	hpfcand_deltaEta
TH1F *	hpfcand_deltaEta_weightE
TH1F *	hpfcand_deltaPhiOverQ
TH1F *	hpfcand_deltaPhiOverQ_weightE
TH1F *	hTauElecDiscriminant
double	maxDeltaR_
double	maxMCAbsEta_
double	maxRecoAbsEta_
double	minMCPt_
double	minRecoPt_
std::string	outputFile_
std::string	sGenMatchObjectLabel_

Detailed Description

Definition at line 36 of file PFTauElecRejectionBenchmark.h.

Constructor & Destructor Documentation

PFTauElecRejectionBenchmark::PFTauElecRejectionBenchmark ( )

Definition at line 24 of file PFTauElecRejectionBenchmark.cc.

: file_(0) {}

PFTauElecRejectionBenchmark::~PFTauElecRejectionBenchmark ( ) [virtual]

Definition at line 26 of file PFTauElecRejectionBenchmark.cc.

References file_.

                                                          {
  if(file_) file_->Close();
}

Member Function Documentation

bool PFTauElecRejectionBenchmark::isInEcalCrack ( double eta ) const [private]

Definition at line 357 of file PFTauElecRejectionBenchmark.cc.

Referenced by process().

                                                               {  
  return (eta < 0.018 || 
          (eta>0.423 && eta<0.461) ||
          (eta>0.770 && eta<0.806) ||
          (eta>1.127 && eta<1.163) ||
          (eta>1.460 && eta<1.558));
}

void PFTauElecRejectionBenchmark::process	(	edm::Handle< edm::HepMCProduct >	mcevt,
		edm::Handle< reco::PFTauCollection >	pfTaus,
		edm::Handle< reco::PFTauDiscriminator >	pfTauIsoDiscr,
		edm::Handle< reco::PFTauDiscriminator >	pfTauElecDiscr
	)

Definition at line 216 of file PFTauElecRejectionBenchmark.cc.

References _GenObjects, abs, applyEcalCrackCut_, HLTFastRecoForTau_cff::deltaEta, SiPixelRawToDigiRegional_cfi::deltaPhi, reco::tau::disc::Eta(), hElecMVA, hEmfrac, hEmfrac_barrel, hEmfrac_endcap, hEmfrac_preid0, hEmfrac_preid1, hEmfracvsEoP, hEmfracvsEoP_preid0, hEmfracvsEoP_preid1, hEoverP, hEoverP_barrel, hEoverP_endcap, hEoverP_preid0, hEoverP_preid1, hHoPvsEoP, hHoPvsEoP_preid0, hHoPvsEoP_preid1, hHoverP, hHoverP_barrel, hHoverP_endcap, hHoverP_preid0, hHoverP_preid1, hleadTk_eta, hleadTk_phi, hleadTk_pt, hpfcand_deltaEta, hpfcand_deltaEta_weightE, hpfcand_deltaPhiOverQ, hpfcand_deltaPhiOverQ_weightE, hTauElecDiscriminant, i, isInEcalCrack(), edm::Ref< C, T, F >::isNonnull(), maxDeltaR_, maxMCAbsEta_, maxRecoAbsEta_, minMCPt_, minRecoPt_, AlCaHLTBitMon_ParallelJobs::p, sGenMatchObjectLabel_, edm::RefVector< C, T, F >::size(), metsig::tau, and z.

                                                                                            {


  // Find Gen Objects to be matched with
  HepMC::GenEvent * generated_event = new HepMC::GenEvent(*(mcevt->GetEvent()));
  _GenObjects.clear();
  
  TLorentzVector taunet;
  HepMC::GenEvent::particle_iterator p;
  for (p = generated_event->particles_begin(); p != generated_event->particles_end(); p++) {
    if(std::abs((*p)->pdg_id()) == 15&&(*p)->status()==2) { 
      bool lept_decay = false;     
      TLorentzVector tau((*p)->momentum().px(),(*p)->momentum().py(),(*p)->momentum().pz(),(*p)->momentum().e());
      HepMC::GenVertex::particle_iterator z = (*p)->end_vertex()->particles_begin(HepMC::descendants);
      for(; z != (*p)->end_vertex()->particles_end(HepMC::descendants); z++) {
        if(std::abs((*z)->pdg_id()) == 11 || std::abs((*z)->pdg_id()) == 13) lept_decay=true;
        if(std::abs((*z)->pdg_id()) == 16)
          taunet.SetPxPyPzE((*z)->momentum().px(),(*z)->momentum().py(),(*z)->momentum().pz(),(*z)->momentum().e());
        
      }
      if(lept_decay==false) {
        TLorentzVector jetMom=tau-taunet;
        if (sGenMatchObjectLabel_=="tau") _GenObjects.push_back(jetMom);
      }
    } else if(std::abs((*p)->pdg_id()) == 11&&(*p)->status()==1) { 
      TLorentzVector elec((*p)->momentum().px(),(*p)->momentum().py(),(*p)->momentum().pz(),(*p)->momentum().e());
      if (sGenMatchObjectLabel_=="e") _GenObjects.push_back(elec);
    } 
  }
 

  
  // Loop over all PFTaus
  math::XYZPointF myleadTkEcalPos;
  for (PFTauCollection::size_type iPFTau=0;iPFTau<pfTaus->size();iPFTau++) { 
    PFTauRef thePFTau(pfTaus,iPFTau); 
    if ((*pfTauIsoDiscr)[thePFTau] == 1) {
      if ((*thePFTau).et() > minRecoPt_ && std::abs((*thePFTau).eta()) < maxRecoAbsEta_) {

        // Check if track goes to Ecal crack
        TrackRef myleadTk;
        if(thePFTau->leadPFChargedHadrCand().isNonnull()){
          myleadTk=thePFTau->leadPFChargedHadrCand()->trackRef();
          myleadTkEcalPos = thePFTau->leadPFChargedHadrCand()->positionAtECALEntrance();
          
          if(myleadTk.isNonnull()){ 
            if (applyEcalCrackCut_ && isInEcalCrack(std::abs((double)myleadTkEcalPos.eta()))) {
              continue; // do nothing
            } else {

              // Match with gen object
              for (unsigned int i = 0; i<_GenObjects.size();i++) {
                if (_GenObjects[i].Et() >= minMCPt_ && std::abs(_GenObjects[i].Eta()) < maxMCAbsEta_ ) {
                  TLorentzVector pftau((*thePFTau).px(),(*thePFTau).py(),(*thePFTau).pz(),(*thePFTau).energy());
                  double GenDeltaR = pftau.DeltaR(_GenObjects[i]);
                  if (GenDeltaR<maxDeltaR_) {
                    
                    hleadTk_pt->Fill((float)myleadTk->pt());
                    hleadTk_eta->Fill((float)myleadTk->eta());
                    hleadTk_phi->Fill((float)myleadTk->phi());

                    hEoverP->Fill((*thePFTau).ecalStripSumEOverPLead());
                    hHoverP->Fill((*thePFTau).hcal3x3OverPLead());
                    hEmfrac->Fill((*thePFTau).emFraction());

                    if (std::abs(myleadTk->eta())<1.5) {
                      hEoverP_barrel->Fill((*thePFTau).ecalStripSumEOverPLead());
                      hHoverP_barrel->Fill((*thePFTau).hcal3x3OverPLead());
                      hEmfrac_barrel->Fill((*thePFTau).emFraction());
                    } else if (std::abs(myleadTk->eta())>1.5 && std::abs(myleadTk->eta())<2.5) {
                      hEoverP_endcap->Fill((*thePFTau).ecalStripSumEOverPLead());
                      hHoverP_endcap->Fill((*thePFTau).hcal3x3OverPLead());
                      hEmfrac_endcap->Fill((*thePFTau).emFraction());
                    }

                    // if -999 fill in -1 bin!
                    if ((*thePFTau).electronPreIDOutput()<-1)
                      hElecMVA->Fill(-1);
                    else
                      hElecMVA->Fill((*thePFTau).electronPreIDOutput());

                    hTauElecDiscriminant->Fill((*pfTauElecDiscr)[thePFTau]);

                    hHoPvsEoP->Fill((*thePFTau).ecalStripSumEOverPLead(),(*thePFTau).hcal3x3OverPLead());
                    hEmfracvsEoP->Fill((*thePFTau).emFraction(),(*thePFTau).hcal3x3OverPLead());

                    if ((*thePFTau).electronPreIDDecision()==1) {
                      hEoverP_preid1->Fill((*thePFTau).ecalStripSumEOverPLead());
                      hHoverP_preid1->Fill((*thePFTau).hcal3x3OverPLead());
                      hEmfrac_preid1->Fill((*thePFTau).emFraction());
                      hHoPvsEoP_preid1->Fill((*thePFTau).ecalStripSumEOverPLead(),(*thePFTau).hcal3x3OverPLead());
                      hEmfracvsEoP_preid1->Fill((*thePFTau).emFraction(),(*thePFTau).hcal3x3OverPLead());
                    } else {
                      hEoverP_preid0->Fill((*thePFTau).ecalStripSumEOverPLead());
                      hHoverP_preid0->Fill((*thePFTau).hcal3x3OverPLead());
                      hEmfrac_preid0->Fill((*thePFTau).emFraction());
                      hHoPvsEoP_preid0->Fill((*thePFTau).ecalStripSumEOverPLead(),(*thePFTau).hcal3x3OverPLead());
                      hEmfracvsEoP_preid0->Fill((*thePFTau).emFraction(),(*thePFTau).hcal3x3OverPLead());
                    }


                  }
                }

                // Loop over all PFCands for cluster plots  
                PFCandidateRefVector myPFCands=(*thePFTau).pfTauTagInfoRef()->PFCands();
                for(int i=0;i<(int)myPFCands.size();i++){

                  math::XYZPointF candPos;
                  if (myPFCands[i]->particleId()==1 || myPFCands[i]->particleId()==2) // if charged hadron or electron
                    candPos = myPFCands[i]->positionAtECALEntrance();
                  else
                    candPos = math::XYZPointF(myPFCands[i]->px(),myPFCands[i]->py(),myPFCands[i]->pz());

                  //double deltaR   = ROOT::Math::VectorUtil::DeltaR(myleadTkEcalPos,candPos);
                  double deltaPhi = ROOT::Math::VectorUtil::DeltaPhi(myleadTkEcalPos,candPos);
                  double deltaEta = std::abs(myleadTkEcalPos.eta()-myPFCands[i]->eta());
                  double deltaPhiOverQ = deltaPhi/(double)myleadTk->charge();
                  
                  hpfcand_deltaEta->Fill(deltaEta);
                  hpfcand_deltaEta_weightE->Fill(deltaEta*myPFCands[i]->ecalEnergy());
                  hpfcand_deltaPhiOverQ->Fill(deltaPhiOverQ);
                  hpfcand_deltaPhiOverQ_weightE->Fill(deltaPhiOverQ*myPFCands[i]->ecalEnergy());        
                  
                }

              }

            }
          }
        }

      }
    }
  }
}

void PFTauElecRejectionBenchmark::setup	(	std::string	Filename,
		std::string	benchmarkLabel,
		double	maxDeltaR,
		double	minRecoPt,
		double	maxRecoAbsEta,
		double	minMCPt,
		double	maxMCAbsEta,
		std::string	sGenMatchObjectLabel,
		bool	applyEcalCrackCut,
		DQMStore *	db_store
	)

void PFTauElecRejectionBenchmark::write ( void )

Definition at line 30 of file PFTauElecRejectionBenchmark.cc.

References gather_cfg::cout, db_, file_, outputFile_, and DQMStore::save().

                                        {
   // Store the DAQ Histograms 
  if (outputFile_.size() != 0) {
    if (db_)
          db_->save(outputFile_.c_str());
    // use bare Root if no DQM (FWLite applications)
    else if (file_) {
       file_->Write(outputFile_.c_str());
       cout << "Benchmark output written to file " << outputFile_.c_str() << endl;
       file_->Close();
       }
  }
  else 
    cout << "No output file specified ("<<outputFile_<<"). Results will not be saved!" << endl;
    
}