ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/LeptonSelection/src/FSR.cc
Revision: 1.2
Committed: Mon May 28 19:04:28 2012 UTC (12 years, 11 months ago) by khahn
Content type: text/plain
Branch: MAIN
Changes since 1.1: +76 -3 lines
Log Message:
make FSR switchable, add betaCorrection

File Contents

# User Rev Content
1 khahn 1.1 #include <vector>
2     #include "MathUtils.h"
3     #include "FSR.h"
4     #include "CommonDefs.h"
5     #include "TLorentzVector.h"
6    
7     using namespace std;
8    
9     extern vector<bool> PFnoPUflag;
10    
11     //--------------------------------------------------------------------------------------------------
12 khahn 1.2 double betaCorrectedIso(ControlFlags & ctrl,
13     const mithep::PFCandidate * photon,
14     const mithep::Array<mithep::PFCandidate> * fPFCandidates )
15     //--------------------------------------------------------------------------------------------------
16     {
17    
18     //
19     // final iso
20     //
21     Double_t fChargedIso = 0.0;
22     Double_t fGammaIso = 0.0;
23     Double_t fNeutralHadronIso = 0.0;
24    
25     //
26     // Loop over PF Candidates
27     //
28     for(int k=0; k<fPFCandidates->GetEntries(); ++k) {
29    
30     if( !(PFnoPUflag[k]) ) continue; // my PF no PU hack
31     const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*fPFCandidates)[k]);
32    
33     Double_t deta = (photon->Eta() - pf->Eta());
34     Double_t dphi = mithep::MathUtils::DeltaPhi(Double_t(photon->Phi()),Double_t(pf->Phi()));
35     Double_t dr = mithep::MathUtils::DeltaR(photon->Phi(),photon->Eta(), pf->Phi(), pf->Eta());
36     if (dr > 0.4) continue;
37    
38     // skip this photon
39     if( abs(pf->PFType()) == mithep::PFCandidate::eGamma &&
40     pf->Et() == photon->Et() ) continue;
41    
42     if (dr < 1.0) {
43    
44     //
45     // Charged Iso
46     //
47     if (pf->Charge() != 0 && (pf->HasTrackerTrk()||pf->HasGsfTrk()) ) {
48     if( dr < 0.01 ) continue;
49     // if (abs(pf->PFType()) == mithep::PFCandidate::eElectron ||
50     // abs(pf->PFType()) == mithep::PFCandidate::eMuon) continue;
51     fChargedIso += pf->Pt();
52     }
53    
54     //
55     // Gamma Iso
56     //
57     else if (abs(pf->PFType()) == mithep::PFCandidate::eGamma) {
58     if( pf->Pt() > 0.5 && dr > 0.08) // need the inner veto?
59     fGammaIso += pf->Pt();
60     }
61    
62     //
63     // Other Neutrals
64     //
65     else {
66     // KH, add to sync
67     if( pf->Pt() > 0.5 )
68     fNeutralHadronIso += pf->Pt();
69     }
70    
71     }
72    
73     }
74    
75     double pfIso = fChargedIso + fGammaIso + fNeutralHadronIso;
76     return pfIso;
77     }
78    
79    
80    
81     //--------------------------------------------------------------------------------------------------
82 khahn 1.1 // typeI = PF IDed photons. NB : repurpose PFnoPUflag, flip for recovered photons
83     // so that they are skipped in the isolation calculation
84     //--------------------------------------------------------------------------------------------------
85 khahn 1.2 bool recover_typeI_Photon( ControlFlags & ctrl,
86     mithep::Muon * mu,
87 khahn 1.1 const mithep::Array<mithep::PFCandidate> * pfArr )
88     //--------------------------------------------------------------------------------------------------
89     {
90     vector<int> photonIndices;
91     for( int i=0; i<pfArr->GetEntries(); i++ ) {
92     if( !(PFnoPUflag[i])) continue; // my PF no PU hack
93     const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*pfArr)[i]);
94     if( abs(pf->PFType()) == mithep::PFCandidate::eGamma &&
95     pf->Pt() > 2.0 && fabs(pf->Eta()) < 2.4 ) {
96    
97     float dR = mithep::MathUtils::DeltaR(pf->Phi(),pf->Eta(), mu->Phi(), mu->Eta());
98    
99     // "keep all photons close to one of the 4L muons ..."
100     if( dR < 0.07 ) photonIndices.push_back(i);
101    
102     // "need tighter cuts for other photons ..."
103 khahn 1.2 if( dR < 0.5 && pf->Pt() > 4. && betaCorrectedIso(ctrl, pf, pfArr)/pf->Pt() < 1.0)
104     photonIndices.push_back(i);
105 khahn 1.1 }
106    
107     int index=-1;
108     float highest_pt = -1;
109     const mithep::PFCandidate * thepf;
110     for( int i=0; i<photonIndices.size(); i++ ) {
111     const mithep::PFCandidate *pf = (mithep::PFCandidate*)(pfArr->At(photonIndices[i]));
112     if( pf->Pt() > highest_pt ) {
113     thepf = pf;
114     highest_pt = pf->Pt();
115     index = photonIndices[i];
116     }
117     }
118    
119     if( thepf != NULL ) {
120     // "... remove it from lepton isolation ..."
121     PFnoPUflag[index] = 1;
122     // add to the muon
123     TLorentzVector muvec,phvec,newmuvec;
124     muvec.SetPtEtaPhiM( mu->Pt(), mu->Eta(), mu->Phi(), MUON_MASS);
125     phvec.SetPtEtaPhiM( thepf->Pt(), thepf->Eta(), thepf->Phi(), 0.);
126     newmuvec = muvec+phvec;
127     mu->SetPtEtaPhi (newmuvec.Pt(),
128     newmuvec.Eta(),
129     newmuvec.Phi());
130     return true;
131     }
132     }
133     return false;
134     }
135    
136     //--------------------------------------------------------------------------------------------------
137     // typeII = "PFClusters linked to muons"
138     //--------------------------------------------------------------------------------------------------
139 khahn 1.2 bool recover_typeII_Photon( ControlFlags & ctrl,
140     mithep::Muon * mu,
141 khahn 1.1 const mithep::Array<mithep::PFCandidate> * pfArr )
142     //--------------------------------------------------------------------------------------------------
143     {
144    
145     bool foundPF=false;
146     const mithep::PFCandidate * thepf;
147     for( int i=0; i<pfArr->GetEntries(); i++ ) {
148     if( !(PFnoPUflag[i]) ) continue; // my PF no PU hack
149     const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*pfArr)[i]);
150     if( (pf->TrackerTrk() == mu->TrackerTrk()) && abs(pf->PFType()) == mithep::PFCandidate::eMuon ) {
151     foundPF = true;
152     thepf = pf;
153     break;
154     }
155     }
156    
157     if( foundPF ) {
158     double sintet = thepf->Pt()/thepf->E();
159     double phpt = thepf->EECal() * sintet;
160     if ( thepf->EECal() >= 2.0 && phpt >= 2.0 ) {
161     mu->SetPtEtaPhi (mu->Pt()+phpt,
162     mu->Eta(),
163     mu->Phi());
164     return true;
165     }
166     }
167    
168    
169     return false;
170     }