CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_6_1_2_SLHC4_patch1/src/HLTriggerOffline/Tau/src/L25TauAnalyzer.cc

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 //
00003 // Package:    L25TauAnalyzer
00004 // Class:      L25TauAnalyzer
00005 // 
00013 //
00014 // Original Author:  Eduardo Luiggi
00015 //         Created:  Mon Apr 19 16:09:13 CDT 2010
00016 // $Id: L25TauAnalyzer.cc,v 1.11 2011/03/01 22:54:26 eluiggi Exp $
00017 //
00018 //
00019 
00020 
00021 #include "HLTriggerOffline/Tau/interface/L25TauAnalyzer.h"
00022 #include "Math/GenVector/VectorUtil.h"
00023 #include "FWCore/ServiceRegistry/interface/Service.h"
00024 #include "CommonTools/UtilAlgos/interface/TFileService.h"
00025 #include "DataFormats/Math/interface/deltaR.h"
00026 
00027 using namespace edm;
00028 using namespace reco;
00029 using namespace std;
00030 
00031 L25TauAnalyzer::L25TauAnalyzer(const edm::ParameterSet& iConfig){
00032    //now do what ever initialization is needed
00033 
00034   _l25JetSource = iConfig.getParameter<InputTag>("L25JetSource");
00035   _l2TauInfoAssoc = iConfig.getParameter<InputTag>("L2TauSource");
00036   _pfTauSource = iConfig.getParameter<InputTag>("PFTauSource");
00037   _pfTauIsoSource = iConfig.getParameter<InputTag>("PFTauIsoSource");
00038   _pfTauMuonDiscSource = iConfig.getParameter<InputTag>("PFTauMuonDiscSource");
00039   _pVtxSource = iConfig.getParameter<InputTag>("PrimaryVtxSource");
00040   _l2l25MatchingCone = iConfig.getParameter<double>("L2L25MatchingCone");
00041   
00042   _l25JetLeadTkMacthingCone =  iConfig.getParameter<double>("L25JetLeadTkMatchingCone");
00043   _minTrackPt = iConfig.getParameter<double>("MinTrackPt");
00044   _signalCone = iConfig.getParameter<double>("SignalCone");
00045   _isolationCone = iConfig.getParameter<double>("IsolationCone");
00046   _l25Dz = iConfig.getParameter<double>("L25DzCut");
00047   _nTrkIso = iConfig.getParameter<int>("NTrkIso");
00048   _l25LeadTkPtMin = iConfig.getParameter<double>("L25LeadTkMinPt");
00049   
00050      
00051   l25JetEt=0.;
00052   l25JetEta=10.;
00053   l25JetPhi=10.;
00054   hasLeadTrk=0;
00055   leadSignalTrackPt = 0.;
00056   leadTrkJetDeltaR = 0.;
00057   numPixTrkInJet = 0;
00058   numQPixTrkInJet = 0;
00059   numQPixTrkInSignalCone = 0;
00060   numQPixTrkInAnnulus = 10;
00061   
00062   pfTauPt = 0.;
00063   pfTauEt = 0.;
00064   pfTauEta = 10.;
00065   pfTauPhi = 10.;
00066   pfTauLTPt = 0.;
00067   pfTauTrkIso = 100.;
00068   pfTauGammaIso = 100.;
00069   
00070   l2JetEt = 0;
00071   l2JetEta = 10;
00072   l2JetPhi = 10;
00073   
00074   l25MatchedToPFTau = 0;
00075   L2MatchedToPFtau = 0;
00076   L25MatchedToL2 = 0;
00077   
00078   initializeVectors();
00079   edm::Service<TFileService> fs;
00080   l25tree = fs->make<TTree>("l25tree","Level 2.5 Tau Tree");
00081   
00082   l25tree->Branch("pfTauHasLeadTrk", &pfTauHasLeadTrk, "pfTauHasLeadTrk/B");
00083   l25tree->Branch("pfTauInBounds", &pfTauInBounds, "pfTauInBounds/B");
00084   l25tree->Branch("pfTauPt", &pfTauPt, "pfTauPt/F");
00085   l25tree->Branch("pfTauEt", &pfTauEt, "pfTauEt/F");
00086   l25tree->Branch("pfTauEta", &pfTauEta, "pfTauEta/F");
00087   l25tree->Branch("pfTauphi", &pfTauPhi, "pfTauPhi/F");
00088   l25tree->Branch("pfTauLTPt", &pfTauLTPt, "pfTauLTPt/F");
00089   l25tree->Branch("pfTauIsoDisc", &pfTauIsoDisc, "pfTauIsoDisc/F");
00090   l25tree->Branch("pfTauMuonDisc", &pfTauMuonDisc, "pfTauMuonDisc/F");
00091   l25tree->Branch("pfTauNProngs", &pfTauNProngs, "pfTauNProngs/I");
00092   l25tree->Branch("pfTauTrkIso", &pfTauTrkIso, "pfTauTrkIso/F");
00093   l25tree->Branch("pfTauNTrkIso", &pfTauNTrkIso, "pfTauNTrkIso/I");
00094   l25tree->Branch("pfTauIsoTrkPt", &pfTauIsoTrkPt, "pfTauIsoTrkPt/F");
00095   l25tree->Branch("pfTauGammaIso", &pfTauGammaIso, "pfTauGammaIso/F");
00096   
00097   l25tree->Branch("pftauL25DeltaR", &pftauL25DeltaR, "pftauL25DeltaR/F");
00098   l25tree->Branch("pftauSignalTrkDeltaR", &pftauSignalTrkDeltaR, "pftauSignalTrkDeltaR/F");
00099   l25tree->Branch("pftauIsoTrkDeltaR", &pftauIsoTrkDeltaR, "pftauIsoTrkDeltaR/F");
00100   l25tree->Branch("leadTrkPtRes", &leadTrkPtRes, "leadTrkPtRes/F");
00101   l25tree->Branch("leadTrkDeltaR", &leadTrkDeltaR, "leadTrkDeltaR/F");
00102   l25tree->Branch("leadIsoTrkPtRes", &leadIsoTrkPtRes, "leadIsoTrkPtRes/F");
00103   l25tree->Branch("leadIsoTrkDeltaR", &leadIsoTrkDeltaR, "leadIsoTrkDeltaR/F");
00104   l25tree->Branch("l25Disc_LeadTkDir", &l25Disc_LeadTkDir, "l25Disc_LeadTkDir/B");
00105   l25tree->Branch("l25Disc_JetDir", &l25Disc_JetDir, "l25Disc_JetDir/B");
00106  
00107   l25tree->Branch("l2JetEt", &l2JetEt, "l2JetEt/F");
00108   l25tree->Branch("l2JetEta", &l2JetEta, "l2JetEta/F");
00109   l25tree->Branch("l2JetPhi", &l2JetPhi, "l2JetPhi/F");
00110  
00111   l25tree->Branch("l25MatchedToPFTau", &l25MatchedToPFTau, "l25MatchedToPFTau/B");
00112   l25tree->Branch("L2MatchedToPFtau", &L2MatchedToPFtau, "L2MatchedToPFtau/B");
00113   l25tree->Branch("L25MatchedToL2", &L25MatchedToL2, "L25MatchedToL2/B");
00114   
00115   l25tree->Branch("l25JetEt", &l25JetEt, "l25JetEt/F");
00116   l25tree->Branch("l25JetEta", &l25JetEta, "l25JetEta/F");
00117   l25tree->Branch("l25JetPhi", &l25JetPhi, "l25JetPhi/F");
00118   l25tree->Branch("hasLeadTrack", &hasLeadTrk, "hasLeadTrack/B");
00119   l25tree->Branch("leadTrackPt", &leadSignalTrackPt, "leadTrackPt/F");
00120   l25tree->Branch("nTracks", &numPixTrkInJet, "nTracks/I");
00121   l25tree->Branch("nQTracks", &numQPixTrkInJet, "nQTracks/I");
00122   l25tree->Branch("nQTracksInSignal", &numQPixTrkInSignalCone, "nQTracksInSignal/I");
00123   l25tree->Branch("nQTracksInAnnulus", &numQPixTrkInAnnulus, "nQTracksInAnnulus/I");
00124   
00125   l25tree->Branch("l25TrkPt", &l25TrkPt);
00126   l25tree->Branch("l25TrkEta", &l25TrkEta);
00127   l25tree->Branch("l25TrkPhi", &l25TrkPhi);
00128   l25tree->Branch("l25TrkDz", &l25TrkDz);
00129   l25tree->Branch("l25TrkDxy", &l25TrkDxy);
00130   l25tree->Branch("l25TrkChi2", &l25TrkChi2);
00131   l25tree->Branch("l25TrkChi2NdF", &l25TrkChi2NdF);
00132   l25tree->Branch("l25TrkNRecHits", &l25TrkNRecHits);
00133   l25tree->Branch("l25TrkNValidPixelHits", &l25TrkNValidPixelHits);
00134   
00135   l25tree->Branch("l25SignalTrkPt", &l25SignalTrkPt);
00136   l25tree->Branch("l25SignalTrkEta", &l25SignalTrkEta);
00137   l25tree->Branch("l25SignalTrkPhi", &l25SignalTrkPhi);
00138   l25tree->Branch("l25SignalTrkDz", &l25SignalTrkDz);
00139   l25tree->Branch("l25SignalTrkDxy", &l25SignalTrkDxy);
00140   l25tree->Branch("l25SignalTrkChi2", &l25SignalTrkChi2);
00141   l25tree->Branch("l25SignalTrkChi2NdF", &l25SignalTrkChi2NdF);
00142   l25tree->Branch("l25SignalTrkNRecHits", &l25SignalTrkNRecHits);
00143   l25tree->Branch("l25SignalTrkNValidHits", &l25SignalTrkNValidHits);
00144   l25tree->Branch("l25SignalTrkNValidPixelHits", &l25SignalTrkNValidPixelHits);
00145   l25tree->Branch("l25SignalTrkNLostHits", &l25SignalTrkNLostHits);
00146   
00147   l25tree->Branch("l25IsoTrkPt", &l25IsoTrkPt);
00148   l25tree->Branch("l25IsoTrkEta", &l25IsoTrkEta);
00149   l25tree->Branch("l25IsoTrkPhi", &l25IsoTrkPhi);
00150   l25tree->Branch("l25IsoTrkDz", &l25IsoTrkDz);
00151   l25tree->Branch("l25IsoTrkDxy", &l25IsoTrkDxy);
00152   l25tree->Branch("l25IsoTrkChi2", &l25IsoTrkChi2);
00153   l25tree->Branch("l25IsoTrkChi2NdF", &l25IsoTrkChi2NdF);
00154   l25tree->Branch("l25IsoTrkNRecHits", &l25IsoTrkNRecHits);
00155   l25tree->Branch("l25IsoTrkNValidHits", &l25IsoTrkNValidHits);
00156   l25tree->Branch("l25IsoTrkNValidPixelHits", &l25IsoTrkNValidPixelHits);
00157   l25tree->Branch("l25IsoTrkNLostHits", &l25IsoTrkNLostHits);
00158   
00159   l25tree->Branch("myNtrkIso", &myNtrkIso, "myNtrkIso/I");
00160 }
00161 
00162 
00163 L25TauAnalyzer::~L25TauAnalyzer(){
00164  
00165    // do anything here that needs to be done at desctruction time
00166    // (e.g. close files, deallocate resources etc.)
00167 
00168 }
00169 
00170 
00171 //
00172 // member functions
00173 //
00174 
00175 // ------------ method called to for each event  ------------
00176 void L25TauAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup){
00177    using namespace edm;
00178 
00179   //get l2 jet; match to pftau
00180   Handle<PFTauCollection> thePFTaus;
00181   iEvent.getByLabel(_pfTauSource,thePFTaus);
00182   
00183   Handle<PFTauDiscriminator> thePFTauDiscriminatorByIsolation;
00184   iEvent.getByLabel(_pfTauIsoSource,thePFTauDiscriminatorByIsolation);  
00185   
00186   Handle<PFTauDiscriminator> thePFTauDiscriminatorAgainstMuon;
00187   iEvent.getByLabel(_pfTauMuonDiscSource, thePFTauDiscriminatorAgainstMuon); 
00188   
00189   Handle<VertexCollection> thePrimaryVertices;
00190   iEvent.getByLabel(_pVtxSource, thePrimaryVertices);
00191   const reco::Vertex& theHltPrimaryVertex = (*thePrimaryVertices->begin());
00192   theVertexPosition = math::XYZPoint(0.,0.,0.);
00193   if(thePrimaryVertices->size() > 0) theVertexPosition = math::XYZPoint(theHltPrimaryVertex.position().x(), 
00194                                                      theHltPrimaryVertex.position().y(),
00195                                                      theHltPrimaryVertex.position().z());  
00196   
00197   //Loop over PFtaus
00198   if(thePFTaus.isValid()){
00199     edm::LogInfo("L25Analysis")  << "Inside PFTau\n";
00200     for(unsigned int pfTauIt = 0; pfTauIt < thePFTaus->size(); ++pfTauIt){
00201       const PFTauRef thisTauRef(thePFTaus,pfTauIt);
00202       pfTauIsoDisc = 0;
00203       pfTauMuonDisc = 0;
00204       if(thePFTauDiscriminatorByIsolation.isValid()){
00205         const PFTauDiscriminator& disc = *(thePFTauDiscriminatorByIsolation.product());
00206         pfTauIsoDisc = disc[thisTauRef];
00207       }
00208       if(thePFTauDiscriminatorAgainstMuon.isValid()){
00209         const PFTauDiscriminator& disc = *(thePFTauDiscriminatorAgainstMuon.product());
00210         pfTauMuonDisc = disc[thisTauRef];
00211       }
00212       
00213       L2MatchedToPFtau = false;      
00214       pfTauHasLeadTrk = false;
00215       pfTauInBounds = true;
00216       
00217       pfTauPt = thePFTaus->at(pfTauIt).pt();                                                      
00218       pfTauEt = thePFTaus->at(pfTauIt).et();                                                      
00219       pfTauEta = thePFTaus->at(pfTauIt).eta();                                            
00220       pfTauPhi = thePFTaus->at(pfTauIt).eta(); 
00221       
00222       const PFCandidateRef& thePFTauLeadTrack = thePFTaus->at(pfTauIt).leadPFChargedHadrCand(); 
00223       if(thePFTauLeadTrack.isNonnull()){
00224         pfTauHasLeadTrk = true;
00225         pfTauNProngs = thePFTaus->at(pfTauIt).signalPFChargedHadrCands().size();                              
00226         pfTauLTPt = thePFTaus->at(pfTauIt).leadPFChargedHadrCand()->pt();
00227         pfTauTrkIso = thePFTaus->at(pfTauIt).isolationPFChargedHadrCandsPtSum();
00228         pfTauGammaIso = thePFTaus->at(pfTauIt).isolationPFGammaCandsEtSum();
00229       }
00230       
00231       const PFCandidateRefVector& theSignalCands = thePFTaus->at(pfTauIt).signalPFChargedHadrCands();
00232       for(PFCandidateRefVector::const_iterator vIt = theSignalCands.begin(); vIt != theSignalCands.end(); ++vIt){
00233         pftauSignalTrkDeltaR = ROOT::Math::VectorUtil::DeltaR((*vIt)->trackRef()->momentum(), thePFTauLeadTrack->momentum());
00234         if(pftauSignalTrkDeltaR > _signalCone) pfTauInBounds = false;
00235       }
00236     
00237       const PFCandidateRefVector& theIsoCands = thePFTaus->at(pfTauIt).isolationPFChargedHadrCands();
00238       pfTauNTrkIso = theIsoCands.size();
00239       float PFTauLeadIsoPt = 0.;
00240       Track thePFTauLeadIsoTrk;      
00241       for(PFCandidateRefVector::const_iterator vIt = theIsoCands.begin(); vIt != theIsoCands.end(); ++vIt){
00242         pftauIsoTrkDeltaR = ROOT::Math::VectorUtil::DeltaR((*vIt)->trackRef()->momentum(), thePFTauLeadTrack->momentum());
00243         if(pftauIsoTrkDeltaR < _isolationCone) pfTauInBounds = false;
00244         pfTauIsoTrkPt = (*vIt)->trackRef()->pt();                     
00245         if((*vIt)->trackRef()->pt() > PFTauLeadIsoPt){
00246           thePFTauLeadIsoTrk = *((*vIt)->trackRef());
00247           PFTauLeadIsoPt = (*vIt)->trackRef()->pt();
00248         }
00249       }
00250       
00251       // matched PFtau to l2 Jet with Et > 5 ... originially 15 but to strong for min bias
00252       Handle<L2TauInfoAssociation> l2TauInfoAssoc; //Handle to the input (L2 Tau Info Association) 
00253       iEvent.getByLabel(_l2TauInfoAssoc,l2TauInfoAssoc);                                           
00254       reco::CaloJet theMatchedL2Jet;
00255       IsolatedTauTagInfo theMatchedL25TauTagInfo;
00256       if(iEvent.getByLabel(_l2TauInfoAssoc,l2TauInfoAssoc) && l2TauInfoAssoc->size()>0){
00257         double matchDr = _l2l25MatchingCone;
00258         for(L2TauInfoAssociation::const_iterator it = l2TauInfoAssoc->begin(); it != l2TauInfoAssoc->end(); ++it){
00259           const reco::CaloJet& l2Jet =*(it->key);
00260           if(l2Jet.et() < 15.) continue;
00261           double delta = ROOT::Math::VectorUtil::DeltaR(thePFTaus->at(pfTauIt).p4(),l2Jet.p4());
00262           if(delta < matchDr){
00263             matchDr = delta;
00264             L2MatchedToPFtau = true;
00265             theMatchedL2Jet = l2Jet;
00266           }
00267         }
00268       }
00269       
00270       //if(L2MatchedToPFtau) match to l25 and fill l25 variables
00271       Handle<IsolatedTauTagInfoCollection> tauTagInfo;
00272       if(L2MatchedToPFtau){
00273         l2JetEt = theMatchedL2Jet.et();
00274         l2JetEta = theMatchedL2Jet.eta();
00275         l2JetPhi = theMatchedL2Jet.phi();
00276         if(iEvent.getByLabel(_l25JetSource, tauTagInfo) && tauTagInfo->size()>0){
00277           L25MatchedToL2 = false;                                                                              
00278           double minDeltaR = _l2l25MatchingCone;                                                                       
00279           //declare l25 tauTagInfo ...                                                                         
00280           //IsolatedTauTagInfo theMatchedL25TauTagInfo;                                                        
00281           for(IsolatedTauTagInfoCollection::const_iterator i = tauTagInfo->begin();i!=tauTagInfo->end();++i){ 
00282             double delta = ROOT::Math::VectorUtil::DeltaR(theMatchedL2Jet.p4(), i->jet()->p4());               
00283             if(delta < minDeltaR){                                                                             
00284               minDeltaR = delta;                                                                               
00285               L25MatchedToL2 = true;                                                                           
00286               theMatchedL25TauTagInfo = *i;
00287             }                                                                                                  
00288           }                                                                                                    
00289         }
00290         
00291         pftauL25DeltaR = ROOT::Math::VectorUtil::DeltaR(thePFTaus->at(pfTauIt).p4(), theMatchedL25TauTagInfo.jet()->p4());
00292         if(pftauL25DeltaR < _l2l25MatchingCone)l25MatchedToPFTau = true;
00293         
00294         const TrackRefVector theTauTagTracks = theMatchedL25TauTagInfo.allTracks();
00295         
00296         l25JetEt = theMatchedL25TauTagInfo.jet()->et();                                                                                                              
00297         l25JetEta = theMatchedL25TauTagInfo.jet()->eta();
00298         l25JetPhi = theMatchedL25TauTagInfo.jet()->phi();
00299         numPixTrkInJet = theTauTagTracks.size();                                                                                             
00300         numQPixTrkInJet = theMatchedL25TauTagInfo.selectedTracks().size();
00301         //get info for all tracks in tauTagInfo
00302         for(TrackRefVector::const_iterator trkIt = theTauTagTracks.begin(); trkIt != theTauTagTracks.end(); ++trkIt){
00303           l25TrkPt->push_back((*trkIt)->pt());
00304           l25TrkEta->push_back((*trkIt)->eta());
00305           l25TrkPhi->push_back((*trkIt)->phi());
00306           l25TrkDz->push_back((*trkIt)->dz(theVertexPosition));
00307           l25TrkDxy->push_back((*trkIt)->dxy(theVertexPosition));
00308           l25TrkChi2->push_back((*trkIt)->chi2());
00309           l25TrkChi2NdF->push_back((*trkIt)->normalizedChi2());
00310           l25TrkNRecHits->push_back((*trkIt)->recHitsSize());
00311           l25TrkNValidPixelHits->push_back((*trkIt)->hitPattern().numberOfValidPixelHits());
00312         }
00313         const TrackRef leadTk = theMatchedL25TauTagInfo.leadingSignalTrack(_l25JetLeadTkMacthingCone, _minTrackPt);                                                                      
00314         if(!leadTk){                                                                                                                                
00315           hasLeadTrk=false;                                                                                                                 
00316           numQPixTrkInSignalCone=0;                                                                                                         
00317           numQPixTrkInAnnulus=0;                                                                                                            
00318           leadSignalTrackPt=0;                                                                                                              
00319           leadTrkJetDeltaR=10;
00320           leadTrkDeltaR=10;                                                                                                         
00321         }                                                                                                                                   
00322         else{                                                                                                                               
00323           hasLeadTrk=true;                                                                                                                
00324           leadTrkJetDeltaR = ROOT::Math::VectorUtil::DeltaR(theMatchedL25TauTagInfo.jet()->p4().Vect(), leadTk->momentum());
00325           leadSignalTrackPt=leadTk->pt();
00326           if(pfTauLTPt != 0)leadTrkPtRes = (leadTk->pt() - pfTauLTPt) / pfTauLTPt;
00327           leadTrkDeltaR = ROOT::Math::VectorUtil::DeltaR(thePFTauLeadTrack->momentum(), leadTk->momentum());
00328           
00329           //print info in the case where the pftau is isolated but the tauTag is not
00330           bool l25IsoDisc = theMatchedL25TauTagInfo.discriminator(0.2,0.15,0.5,5.,1.,0,0.2);
00331           if(pfTauTrkIso <= 1. && l25MatchedToPFTau && !l25IsoDisc) printInfo(thePFTaus->at(pfTauIt), theMatchedL25TauTagInfo);
00332           
00333           const TrackRefVector theSignalTracks = theMatchedL25TauTagInfo.tracksInCone(leadTk->momentum(), _signalCone, _minTrackPt);
00334           const TrackRefVector theIsoTracks = theMatchedL25TauTagInfo.tracksInCone(leadTk->momentum(), _isolationCone, _minTrackPt);
00335           numQPixTrkInSignalCone = theSignalTracks.size();                                
00336           if(numQPixTrkInSignalCone > 0) numQPixTrkInAnnulus = theIsoTracks.size() - theSignalTracks.size();
00337           //get the lead track in isolation 
00338           float l25LeadIsoPt = 0.;
00339           Track theL25LeadIsoTrk;
00340           
00341           myNtrkIso = 0;
00342           for(TrackRefVector::const_iterator isoIt = theIsoTracks.begin(); isoIt != theIsoTracks.end(); ++isoIt){
00343             if(ROOT::Math::VectorUtil::DeltaR(leadTk->momentum(), (*isoIt)->momentum()) <= _signalCone){
00344           
00345               l25SignalTrkPt->push_back((*isoIt)->pt());
00346               l25SignalTrkEta->push_back((*isoIt)->eta());
00347               l25SignalTrkPhi->push_back((*isoIt)->phi());
00348               l25SignalTrkDz->push_back((*isoIt)->dz(theVertexPosition));
00349               l25SignalTrkDxy->push_back((*isoIt)->dxy(theVertexPosition));
00350               l25SignalTrkChi2->push_back((*isoIt)->chi2());
00351               l25SignalTrkChi2NdF->push_back((*isoIt)->normalizedChi2());
00352               l25SignalTrkNRecHits->push_back((*isoIt)->recHitsSize());
00353               l25SignalTrkNValidHits->push_back((*isoIt)->numberOfValidHits());
00354               l25SignalTrkNValidPixelHits->push_back((*isoIt)->hitPattern().numberOfValidPixelHits());
00355               l25SignalTrkNLostHits->push_back((*isoIt)->lost());
00356 
00357             }
00358             else{
00359               myNtrkIso++;
00360 
00361               l25IsoTrkPt->push_back((*isoIt)->pt());
00362               l25IsoTrkChi2NdF->push_back((*isoIt)->normalizedChi2());
00363               l25IsoTrkChi2->push_back((*isoIt)->chi2());
00364               l25IsoTrkNValidHits->push_back((*isoIt)->numberOfValidHits());
00365               l25IsoTrkNRecHits->push_back((*isoIt)->recHitsSize());
00366               l25IsoTrkNValidPixelHits->push_back((*isoIt)->hitPattern().numberOfValidPixelHits());
00367               l25IsoTrkNLostHits->push_back((*isoIt)->lost());
00368               l25IsoTrkDxy->push_back((*isoIt)->dxy(theVertexPosition));
00369               l25IsoTrkDz->push_back((*isoIt)->dz(theVertexPosition));
00370               l25IsoTrkEta->push_back((*isoIt)->eta());
00371               l25IsoTrkPhi->push_back((*isoIt)->phi());
00372               if((*isoIt)->pt() > l25LeadIsoPt){
00373                 theL25LeadIsoTrk = **isoIt;
00374                 l25LeadIsoPt = (*isoIt)->pt();
00375               }
00376             }
00377           }
00378           //if((nIsoTraks > numQPixTrkInAnnulus) ? std::cout << "FUCK \n" : std::cout << "GREAT\n" );
00379           if(thePFTauLeadIsoTrk.pt() != 0) leadIsoTrkPtRes = (theL25LeadIsoTrk.pt() - thePFTauLeadIsoTrk.pt()) /thePFTauLeadIsoTrk.pt();
00380           leadIsoTrkDeltaR = ROOT::Math::VectorUtil::DeltaR(theL25LeadIsoTrk.momentum(), thePFTauLeadIsoTrk.momentum());
00381           l25Disc_LeadTkDir = theMatchedL25TauTagInfo.discriminator(leadTk->momentum(), _l25JetLeadTkMacthingCone, _signalCone, _isolationCone, _l25LeadTkPtMin, _minTrackPt, _nTrkIso, _l25Dz);
00382           l25Disc_JetDir = theMatchedL25TauTagInfo.discriminator(theMatchedL25TauTagInfo.jet()->momentum(), _l25JetLeadTkMacthingCone, _signalCone, _isolationCone, _l25LeadTkPtMin, _minTrackPt, _nTrkIso, _l25Dz);
00383         }
00384       }
00385       //Fill here per pfTau ...
00386      l25tree->Fill();
00387      clearVectors();
00388     }
00389   }
00390   else std::cout << "Invalid PFTau Collection\n"; 
00391 }
00392 
00393 
00394 
00395 reco::PFTau L25TauAnalyzer::match(const reco::Jet& jet, const reco::PFTauCollection& thePFTauColl){
00396 
00397   //Loop On the Collection and see if your tau jet is matched to one there
00398   //Also find the nearest Matched MC Particle to your Jet (to be complete)
00399   // switch reference collection to pftau;  match to this.
00400  
00401  reco::PFTau theMatchedPFTau;
00402  double matchDr=0.3;
00403  
00404  if(thePFTauColl.size()>0)
00405   for(reco::PFTauCollection::const_iterator it = thePFTauColl.begin(); it != thePFTauColl.end(); ++it){
00406     double delta = ROOT::Math::VectorUtil::DeltaR(jet.p4(),(*it).p4());          
00407     if(delta<matchDr){
00408       matchDr = delta;
00409       theMatchedPFTau = *it;                                                          
00410     }                                                                    
00411   }
00412   return theMatchedPFTau;
00413 }
00414 
00415 reco::CaloJet L25TauAnalyzer::matchedToPFTau(const reco::PFTau& thePFTau, const L2TauInfoAssociation& theL2Info){
00416   double matchDr = 0.5;
00417   reco::CaloJet theMatchedJet;
00418   for(L2TauInfoAssociation::const_iterator it = theL2Info.begin(); it != theL2Info.end(); ++it){
00419     const reco::CaloJet& l2Jet =*(it->key);
00420     if(l2Jet.et() < 15.) continue;
00421     double delta = ROOT::Math::VectorUtil::DeltaR(thePFTau.p4(),l2Jet.p4());
00422     if(delta < matchDr){
00423       matchDr = delta;
00424       theMatchedJet = l2Jet;
00425     }
00426   }
00427   return theMatchedJet;
00428 }
00429 
00430 void L25TauAnalyzer::printInfo(const reco::PFTau& thePFTau, const reco::IsolatedTauTagInfo& theTauTagInfo){
00431   const TrackRef theLeadTrack = theTauTagInfo.leadingSignalTrack(_l25JetLeadTkMacthingCone, _minTrackPt);
00432   const TrackRefVector theSignalTracks = theTauTagInfo.tracksInCone(theLeadTrack->momentum(), _signalCone, _minTrackPt);
00433   const TrackRefVector theIsoTracks = theTauTagInfo.tracksInCone(theLeadTrack->momentum(), _isolationCone, _minTrackPt);
00434   
00435   std::cout << " Isolated PFTau Matched to Non-Isolated L25 Object( PFTau:L25)"
00436   << "\n Et\t" << thePFTau.et() << "\t" << theTauTagInfo.jet()->et()
00437   << "\n Eta\t" << thePFTau.eta() << "\t" << theTauTagInfo.jet()->eta()
00438   << "\n Phi\t" << thePFTau.phi() << "\t" << theTauTagInfo.jet()->phi()
00439   << "\n LTPt\t" << thePFTau.leadPFChargedHadrCand()->pt() << "\t" << theLeadTrack->pt()
00440   << "\n LTEta\t" << thePFTau.leadPFChargedHadrCand()->eta() << "\t" << theLeadTrack->eta()
00441   << "\n LTPhi\t" << thePFTau.leadPFChargedHadrCand()->phi() << "\t" << theLeadTrack->phi()
00442   << "\n NIso\t" << thePFTau.isolationPFChargedHadrCands().size() << "\t" << theIsoTracks.size() - theSignalTracks.size()
00443   <<"\n";
00444 
00445   unsigned int nIsoTrk = 0;
00446   std::cout << " Tracks in L2.5 Iso: (Pt:Eta:Phi:DR:Chi2:Chi2/NdF:PxlHits:dxy:dz)\n";
00447   for(TrackRefVector::const_iterator isoIt = theIsoTracks.begin(); isoIt != theIsoTracks.end(); ++isoIt){
00448     double isoTrackLeadTrkDeltaR = deltaR(theLeadTrack->eta(), theLeadTrack->phi(), (*isoIt)->eta(), (*isoIt)->phi());
00449     if(isoTrackLeadTrkDeltaR > _signalCone){
00450       nIsoTrk++;
00451       cout << nIsoTrk 
00452            << "\t" << (*isoIt)->pt() 
00453            << "\t" << (*isoIt)->eta() 
00454            << "\t" << (*isoIt)->phi()
00455            << "\t" << isoTrackLeadTrkDeltaR 
00456            << "\t" << (*isoIt)->chi2()
00457            << "\t" << (*isoIt)->normalizedChi2()
00458            << "\t" << (*isoIt)->hitPattern().numberOfValidPixelHits()
00459            << "\t" << (*isoIt)->dxy(theVertexPosition)
00460            << "\t" << (*isoIt)->dz(theVertexPosition)
00461            << "\t" << theVertexPosition
00462            << "\n";
00463     }
00464   }
00465   nIsoTrk = 0;
00466   std::cout << "Tracks in PFTau Iso: (Pt:Eta:Phi)\n";
00467   for(PFCandidateRefVector::const_iterator isoIt = thePFTau.isolationPFChargedHadrCands().begin(); 
00468      isoIt != thePFTau.isolationPFChargedHadrCands().end(); ++isoIt){
00469     nIsoTrk++;
00470     cout << nIsoTrk << "\t" 
00471          << (*isoIt)->pt() << "\t" 
00472          << (*isoIt)->eta() << "\t" 
00473          << (*isoIt)->phi() << "\n";
00474   }
00475 }
00476 
00477 void L25TauAnalyzer::clearVectors(){
00478  l25TrkPt->clear();
00479  l25TrkEta->clear();
00480  l25TrkPhi->clear();
00481  l25TrkDz->clear();
00482  l25TrkDxy->clear();
00483  l25TrkChi2->clear();
00484  l25TrkChi2NdF->clear();
00485  l25TrkNRecHits->clear();
00486  l25TrkNValidPixelHits->clear();
00487  
00488  l25SignalTrkPt->clear();
00489  l25SignalTrkEta->clear();
00490  l25SignalTrkPhi->clear();
00491  l25SignalTrkDz->clear();
00492  l25SignalTrkDxy->clear();
00493  l25SignalTrkChi2NdF->clear();
00494  l25SignalTrkChi2->clear();
00495  l25SignalTrkNRecHits->clear();
00496  l25SignalTrkNValidHits->clear();
00497  l25SignalTrkNValidPixelHits->clear();
00498  l25SignalTrkNLostHits->clear();
00499  
00500  l25IsoTrkPt->clear();
00501  l25IsoTrkEta->clear();
00502  l25IsoTrkPhi->clear();
00503  l25IsoTrkDz->clear();
00504  l25IsoTrkDxy->clear();
00505  l25IsoTrkChi2->clear();
00506  l25IsoTrkChi2NdF->clear();
00507  l25IsoTrkNRecHits->clear();
00508  l25IsoTrkNValidHits->clear();
00509  l25IsoTrkNValidPixelHits->clear();
00510  l25IsoTrkNLostHits->clear();
00511 }
00512 
00513 void L25TauAnalyzer::initializeVectors(){
00514  l25TrkPt = NULL;
00515  l25TrkEta= NULL;
00516  l25TrkPhi= NULL;
00517  l25TrkDz= NULL;
00518  l25TrkDxy= NULL;
00519  l25TrkChi2= NULL;
00520  l25TrkChi2NdF= NULL;
00521  l25TrkNRecHits= NULL;
00522  l25TrkNValidPixelHits= NULL;
00523  
00524  l25SignalTrkPt= NULL;
00525  l25SignalTrkEta= NULL;
00526  l25SignalTrkPhi= NULL;
00527  l25SignalTrkDz= NULL;
00528  l25SignalTrkDxy= NULL;
00529  l25SignalTrkChi2= NULL;
00530  l25SignalTrkChi2NdF= NULL;
00531  l25SignalTrkNRecHits= NULL;
00532  l25SignalTrkNValidHits= NULL;
00533  l25SignalTrkNValidPixelHits= NULL;
00534  l25SignalTrkNLostHits= NULL;
00535  
00536  l25IsoTrkPt= NULL;
00537  l25IsoTrkChi2NdF= NULL;
00538  l25IsoTrkChi2= NULL;
00539  l25IsoTrkNValidHits= NULL;
00540  l25IsoTrkNRecHits= NULL;
00541  l25IsoTrkNValidPixelHits= NULL;
00542  l25IsoTrkNLostHits= NULL;
00543  l25IsoTrkDxy= NULL;
00544  l25IsoTrkDz= NULL;
00545  l25IsoTrkEta= NULL;
00546  l25IsoTrkPhi= NULL;
00547 }
00548 
00549 // ------------ method called once each job just before starting event loop  ------------
00550 void L25TauAnalyzer::beginJob(){
00551 }
00552 
00553 // ------------ method called once each job just after ending the event loop  ------------
00554 void L25TauAnalyzer::endJob() {
00555 }
00556 
00557 //define this as a plug-in
00558 //DEFINE_FWK_MODULE(L25TauAnalyzer);