ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Utils/src/RecoilTools.cc
Revision: 1.24
Committed: Sat Jan 12 11:49:50 2013 UTC (12 years, 3 months ago) by pharris
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_029c, Mit_029b, Mit_029a, HEAD
Changes since 1.23: +9 -3 lines
Log Message:
Updated MVA Met and Jet ID

File Contents

# Content
1 #include "MitPhysics/Utils/interface/RecoilTools.h"
2 #include "MitPhysics/Utils/interface/JetTools.h"
3 #include <TFile.h>
4
5 ClassImp(mithep::RecoilTools)
6
7 using namespace mithep;
8
9 RecoilTools::RecoilTools(TString iJetLowPtMVAFile,TString iJetHighPtMVAFile,TString iCutFile,bool i42,JetIDMVA::MVAType iType) {
10 fJetIDMVA = new JetIDMVA();
11 fJetIDMVA->Initialize( JetIDMVA::kMET,iJetLowPtMVAFile,iJetHighPtMVAFile,iType,iCutFile);
12 f42 = i42;
13 }
14 //--------------------------------------------------------------------------------------------------
15 RecoilTools::~RecoilTools() {
16 delete fJetIDMVA;
17 }
18 //--------------------------------------------------------------------------------------------------
19 bool RecoilTools::filter(const PFJet *iJet,Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2) {
20 double pDEta1 = iJet->Eta() - iEta1;
21 double pDPhi1 = fabs(iJet->Phi() - iPhi1); if(pDPhi1 > 2.*TMath::Pi()-pDPhi1) pDPhi1 = 2.*TMath::Pi()-pDPhi1;
22 double pDR1 = sqrt(pDEta1*pDEta1 + pDPhi1*pDPhi1);
23 if(pDR1 < 0.5) return false;
24 double pDEta2 = iJet->Eta() - iEta2;
25 double pDPhi2 = fabs(iJet->Phi() - iPhi2); if(pDPhi2 > 2.*TMath::Pi()-pDPhi2) pDPhi2 = 2.*TMath::Pi()-pDPhi2;
26 double pDR2 = sqrt(pDEta2*pDEta2 + pDPhi2*pDPhi2);
27 if(pDR2 < 0.5) return false;
28 return true;
29 }
30 //--------------------------------------------------------------------------------------------------
31 bool RecoilTools::filter(const PFCandidate *iCand,Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2) {
32 double pDEta1 = iCand->Eta() - iEta1;
33 double pDPhi1 = fabs(iCand->Phi() - iPhi1); if(pDPhi1 > 2.*TMath::Pi()-pDPhi1) pDPhi1 = 2.*TMath::Pi()-pDPhi1;
34 double pDR1 = sqrt(pDEta1*pDEta1 + pDPhi1*pDPhi1);
35 if(pDR1 < 0.3) return false;
36 double pDEta2 = iCand->Eta() - iEta2;
37 double pDPhi2 = fabs(iCand->Phi() - iPhi2); if(pDPhi2 > 2.*TMath::Pi()-pDPhi2) pDPhi2 = 2.*TMath::Pi()-pDPhi2;
38 double pDR2 = sqrt(pDEta2*pDEta2 + pDPhi2*pDPhi2);
39 if(pDR2 < 0.3) return false;
40 return true;
41 }
42 //--------------------------------------------------------------------------------------------------
43 Met RecoilTools::pfRecoil(Double_t iVisPt,Double_t iVisPhi,Double_t iVisSumEt,
44 const PFCandidateCol *iCands) {
45 double lSumEt = 0;
46 FourVectorM lVec(0,0,0,0);
47 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
48 lVec -= iCands->At(i0)->Mom();
49 lSumEt += iCands->At(i0)->Pt();
50 }
51 Met lPFMet(lVec.Px(),lVec.Py());
52 lPFMet.SetSumEt(lSumEt);
53 lPFMet.SetMex (lPFMet.Mex()+iVisPt*cos(iVisPhi));
54 lPFMet.SetMey (lPFMet.Mey()+iVisPt*sin(iVisPhi));
55 lPFMet.SetSumEt(lPFMet.SumEt()-iVisSumEt);
56 return lPFMet;
57 }
58 //--------------------------------------------------------------------------------------------------
59 Met RecoilTools::pfRecoil(double iPhi1,double iEta1,double iPhi2,double iEta2,
60 const PFCandidateCol *iCands) {
61 double lSumEt = 0;
62 FourVectorM lVec(0,0,0,0);
63 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
64 const PFCandidate *pfcand = iCands->At(i0);
65 if(!filter(pfcand,iPhi1,iEta1,iPhi2,iEta2)) continue;
66 lVec -= pfcand->Mom();
67 lSumEt += pfcand->Pt();
68 }
69 Met lPFMet(lVec.Px(),lVec.Py());
70 lPFMet.SetSumEt(lSumEt);
71 return lPFMet;
72
73 }
74 //--------------------------------------------------------------------------------------------------
75 void RecoilTools::addType1(FourVectorM &iVec,Double_t &iSumEt,
76 const PFJetCol *iJets,FactorizedJetCorrector *iJetCorrector,const PileupEnergyDensityCol *iPUEnergyDensity,
77 double iPhi1,double iEta1,double iPhi2,double iEta2) {
78 for(UInt_t i0 = 0; i0 < iJets->GetEntries(); i0++) {
79 const PFJet *pJet = iJets->At(i0);
80 if(!JetTools::passPFLooseId(pJet)) continue;
81 if(!filter(pJet,iPhi1,iEta1,iPhi2,iEta2)) continue;
82 if(fJetIDMVA->correctedPt(pJet,iJetCorrector,iPUEnergyDensity) < 10) continue;
83 double lPt = fJetIDMVA->correctedPt(pJet,iJetCorrector,iPUEnergyDensity,RhoUtilities::DEFAULT,100);
84 FourVectorM pVec(0,0,0,0);
85 pVec.SetPt(lPt); pVec.SetEta(pJet->Eta()); pVec.SetPhi(pJet->Phi()); pVec.SetM(pJet->Mass());
86 iVec -= pVec;
87 iSumEt += pVec.Pt();
88 }
89 }
90 //--------------------------------------------------------------------------------------------------
91 Met RecoilTools::pfRecoilType1(Double_t iVisPt,Double_t iVisPhi,Double_t iVisSumEt,
92 const PFCandidateCol *iCands,const PFJetCol *iJets,
93 FactorizedJetCorrector *iJetCorrector,const PileupEnergyDensityCol *iPUEnergyDensity,
94 double iPhi1,double iEta1,double iPhi2,double iEta2) {
95 double lSumEt = 0;
96 FourVectorM lVec(0,0,0,0);
97 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
98 lVec -= iCands->At(i0)->Mom();
99 lSumEt += iCands->At(i0)->Pt();
100 }
101 addType1(lVec,lSumEt,iJets,iJetCorrector,iPUEnergyDensity,iPhi1,iEta1,iPhi2,iEta2);
102 Met lPFMet(lVec.Px(),lVec.Py());
103 lPFMet.SetSumEt(lSumEt);
104 lPFMet.SetMex (lPFMet.Mex()+iVisPt*cos(iVisPhi));
105 lPFMet.SetMey (lPFMet.Mey()+iVisPt*sin(iVisPhi));
106 lPFMet.SetSumEt(lPFMet.SumEt()-iVisSumEt);
107 return lPFMet;
108 }
109 //--------------------------------------------------------------------------------------------------
110 Met RecoilTools::pfCone(double iPhi1,double iEta1,
111 const PFCandidateCol *iCands,const Vertex *iVertex,bool iCharge,Double_t iDZCut) {
112 double lSumEt = 0;
113 FourVectorM lVec(0,0,0,0);
114 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
115 const PFCandidate *pfcand = iCands->At(i0);
116 if(filter(pfcand,iPhi1,iEta1,100,100)) continue;
117 if(iCharge) {
118 double lDZ = -999;
119 if(pfcand->HasTrackerTrk()) lDZ = fabs(pfcand->TrackerTrk()->DzCorrected(*iVertex));//
120 else if(pfcand->HasGsfTrk()) lDZ = fabs(pfcand->GsfTrk() ->DzCorrected(*iVertex));
121 if( fabs(lDZ) > iDZCut) continue;
122 }
123 lVec += pfcand->Mom();
124 lSumEt += pfcand->Pt();
125 }
126 Met lPFMet(lVec.Px(),lVec.Py());
127 lPFMet.SetSumEt(lSumEt);
128 return lPFMet;
129 }
130 //--------------------------------------------------------------------------------------------------
131 Met RecoilTools::trackMet(const PFCandidateCol *iCands,const Vertex *iVertex,Double_t iDZCut) {
132 double trkMetx = 0;
133 double trkMety = 0;
134 double trkSumEt = 0;
135 for(UInt_t i0=0; i0<iCands->GetEntries(); i0++) {
136 const PFCandidate *pfcand = iCands->At(i0);
137 double lDZ = -999;
138 if(pfcand->HasTrackerTrk()) lDZ = fabs(pfcand->TrackerTrk()->DzCorrected(*iVertex));//
139 else if(pfcand->HasGsfTrk()) lDZ = fabs(pfcand->GsfTrk() ->DzCorrected(*iVertex));
140 if( fabs(lDZ) > iDZCut) continue;
141 trkMetx -= pfcand->Px();
142 trkMety -= pfcand->Py();
143 trkSumEt += pfcand->Pt();
144 }
145 Met lMet(trkMetx,trkMety);
146 lMet.SetSumEt(trkSumEt);
147 return lMet;
148 }
149 //--------------------------------------------------------------------------------------------------
150 Met RecoilTools::trackRecoil(double iPhi1,double iEta1,double iPhi2,double iEta2,const PFCandidateCol *iCands,const Vertex *iVertex,Double_t iDZCut) {
151 double trkMetx = 0;
152 double trkMety = 0;
153 double trkSumEt = 0;
154 for(UInt_t i0=0; i0<iCands->GetEntries(); i0++) {
155 const PFCandidate *pfcand = iCands->At(i0);
156 if(!filter(pfcand,iPhi1,iEta1,iPhi2,iEta2)) continue;
157 if( (pfcand->HasTrackerTrk() && (fabs(pfcand->TrackerTrk()->DzCorrected(*iVertex))< iDZCut)) ||
158 (pfcand->HasGsfTrk() && (fabs(pfcand->GsfTrk()->DzCorrected(*iVertex)) < iDZCut)) ) {
159 trkMetx -= pfcand->Px();
160 trkMety -= pfcand->Py();
161 trkSumEt += pfcand->Pt();
162 }
163 }
164 Met lMet(trkMetx,trkMety);
165 lMet.SetSumEt(trkSumEt);
166 return lMet;
167 }
168 //--------------------------------------------------------------------------------------------------
169 //Compute the recoil => here this requires the vector sum of the visible components
170 //VisPt => Vector sum pT of the visible non-recoiling components
171 //VisPhi => Vector sum Phi of the visible non-recoiling components
172 //visSumEt => Vector sum Et =f the visible non-recoiling components
173 Met RecoilTools::trackRecoil(Double_t iVisPt,Double_t iVisPhi,Double_t iVisSumEt,
174 const PFCandidateCol *iCands,const Vertex *iVertex,double iDZCut) {
175 Met lTrkMet = trackMet(iCands,iVertex,iDZCut);
176 lTrkMet.SetMex (lTrkMet.Mex()+iVisPt*cos(iVisPhi));
177 lTrkMet.SetMey (lTrkMet.Mey()+iVisPt*sin(iVisPhi));
178 lTrkMet.SetSumEt(lTrkMet.SumEt()-iVisSumEt);
179 return lTrkMet;
180 }
181 //--------------------------------------------------------------------------------------------------
182 void RecoilTools::addNeut(const PFJet *iJet,FourVectorM &iVec,Double_t &iSumEt,
183 FactorizedJetCorrector *iJetCorrector,const PileupEnergyDensityCol *iPUEnergyDensity,
184 int iSign) {
185 FourVectorM lVec(0,0,0,0);
186 double lPt = fJetIDMVA->correctedPt(iJet,iJetCorrector,iPUEnergyDensity);
187 double lFrac = (iJet->NeutralEmEnergy()/iJet->RawMom().E() + iJet->NeutralHadronEnergy()/iJet->RawMom().E());
188 if(fabs(iJet->Eta()) > 2.5 && !f42) lFrac = 1.;
189 lPt *= lFrac;
190 lVec.SetPt(lPt); lVec.SetEta(iJet->Eta()); lVec.SetPhi(iJet->Phi()); lVec.SetM(iJet->Mass());
191 if(iSign > 0) iVec -= lVec;
192 if(iSign < 0) iVec += lVec;
193 //iSumEt += lPt;
194 //=== Above was a bug in the training
195 if(iSign > 0) iSumEt += lPt;
196 if(iSign < 0) iSumEt -= lPt;
197 }
198
199 //--------------------------------------------------------------------------------------------------
200 //Corrected Jets
201 void RecoilTools::addNeut(const PFJet *iJet,FourVectorM &iVec,Double_t &iSumEt,double iRho,int iSign) {
202 FourVectorM lVec(0,0,0,0);
203 double lPt = iJet->Pt();
204 //if(iJet->RawMom().Pt() < 10) lPt = TMath::Max(iJet->RawMom().Pt()-iJet->JetArea()*iRho,0.);//to be fixed
205 //if(iJet->RawMom().Pt() < 10) lPt = iJet->RawMom().Pt()*iJet->L1OffsetCorrectionScale();
206 lPt *= (iJet->NeutralEmEnergy()/iJet->RawMom().E() + iJet->NeutralHadronEnergy()/iJet->RawMom().E());
207 lVec.SetPt(lPt); lVec.SetEta(iJet->Eta()); lVec.SetPhi(iJet->Phi()); lVec.SetM(iJet->Mass());
208 if(iSign > 0) iVec -= lVec;
209 if(iSign < 0) iVec += lVec;
210 //iSumEt += lPt;
211 //=== Above was a bug in the training
212 if(iSign > 0) iSumEt += lPt;
213 if(iSign < 0) iSumEt -= lPt;
214 }
215
216 //--------------------------------------------------------------------------------------------------
217 Met RecoilTools::NoPUMet( const PFJetCol *iJets,FactorizedJetCorrector *iJetCorrector,
218 const PileupEnergyDensityCol *iPileupEnergyDensity,
219 const PFCandidateCol *iCands,const Vertex *iVertex,const VertexCol *iVertices,
220 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,Double_t iDZCut) {
221
222 FourVectorM lVec (0,0,0,0); double lSumEt = 0;
223 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
224 const PFCandidate *pPF = iCands->At(i0);
225 double lDZ = 999;
226 if(pPF->HasTrackerTrk()) lDZ = fabs(pPF->TrackerTrk()->DzCorrected(*iVertex));//
227 else if(pPF->HasGsfTrk()) lDZ = fabs(pPF->GsfTrk() ->DzCorrected(*iVertex));
228 if( fabs(lDZ) > iDZCut) continue;
229 lSumEt += pPF->Pt();
230 lVec -= pPF->Mom();
231 }
232 int lNPass = 0;
233 for(UInt_t i0 = 0; i0 < iJets->GetEntries(); i0++) {
234 const PFJet *pJet = iJets->At(i0);
235 if(!JetTools::passPFLooseId(pJet)) continue;
236 if(!filter(pJet,iPhi1,iEta1,iPhi2,iEta2)) continue;
237 if(!fJetIDMVA->passPt(pJet,iJetCorrector,iPileupEnergyDensity)) continue;
238 if(!fJetIDMVA->pass(pJet,iVertex,iVertices,iJetCorrector,iPileupEnergyDensity)) continue;
239 double pDPhi1 = fabs(pJet->Phi() - iPhi1); if(pDPhi1 > 2.*TMath::Pi()-pDPhi1) pDPhi1 = 2.*TMath::Pi()-pDPhi1;
240 double pDPhi2 = fabs(pJet->Phi() - iPhi2); if(pDPhi2 > 2.*TMath::Pi()-pDPhi2) pDPhi2 = 2.*TMath::Pi()-pDPhi2;
241 addNeut(pJet,lVec,lSumEt,iJetCorrector,iPileupEnergyDensity);
242 lNPass++;
243 }
244 Met lMet(lVec.Px(),lVec.Py());
245 lMet.SetSumEt( lSumEt);
246 return lMet;
247 }
248 //-------------------------------------------------------------------------------------------------
249 //Corrected Jets
250 Met RecoilTools::NoPUMet( const PFJetCol *iJets,
251 const PFCandidateCol *iCands,const Vertex *iVertex,const VertexCol *iVertices,Double_t iRho,
252 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,Double_t iDZCut) {
253
254 FourVectorM lVec (0,0,0,0); double lSumEt = 0;
255 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
256 const PFCandidate *pPF = iCands->At(i0);
257 double lDZ = 999;
258 if(pPF->HasTrackerTrk()) lDZ = fabs(pPF->TrackerTrk()->DzCorrected(*iVertex));//
259 else if(pPF->HasGsfTrk()) lDZ = fabs(pPF->GsfTrk() ->DzCorrected(*iVertex));
260 if( fabs(lDZ) > iDZCut) continue;
261 lVec -= pPF->Mom();
262 lSumEt += pPF->Pt();
263 }
264 for(UInt_t i0 = 0; i0 < iJets->GetEntries(); i0++) {
265 const PFJet *pJet = iJets->At(i0);
266 if(!filter(pJet,iPhi1,iEta1,iPhi2,iEta2)) continue;
267 if(!fJetIDMVA->passPt(pJet)) continue;
268 if(!fJetIDMVA->pass(pJet,iVertex,iVertices)) continue;
269 addNeut(pJet,lVec,lSumEt,iRho);
270 }
271 Met lMet(lVec.Px(),lVec.Py());
272 lMet.SetSumEt(lSumEt);
273 return lMet;
274 }
275 //--------------------------------------------------------------------------------------------------
276 Met RecoilTools::NoPURecoil(Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
277 const PFJetCol *iJets,
278 const PFCandidateCol *iCands ,
279 const Vertex *iVertex,const VertexCol *iVertices,
280 FactorizedJetCorrector *iJetCorrector,
281 const PileupEnergyDensityCol *iPileupEnergyDensity,Double_t iDZCut) {
282
283 FourVectorM lVec (0,0,0,0); double lSumEt = 0;
284 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
285 const PFCandidate *pPF = iCands->At(i0);
286 if(!filter(pPF,iPhi1,iEta1,iPhi2,iEta2)) continue;
287 double lDZ = 999;
288 if(pPF->HasTrackerTrk()) lDZ = fabs(pPF->TrackerTrk()->DzCorrected(*iVertex));//
289 else if(pPF->HasGsfTrk()) lDZ = fabs(pPF->GsfTrk() ->DzCorrected(*iVertex));
290 if( fabs(lDZ) > iDZCut) continue;
291 lSumEt += pPF->Pt();
292 lVec -= pPF->Mom();
293 }
294 int lNPass = 0;
295 for(UInt_t i0 = 0; i0 < iJets->GetEntries(); i0++) {
296 const PFJet *pJet = iJets->At(i0);
297 if(!filter(pJet,iPhi1,iEta1,iPhi2,iEta2)) continue;
298 if(!fJetIDMVA->passPt(pJet,iJetCorrector,iPileupEnergyDensity)) continue;
299 if(iJetCorrector != 0 && iPileupEnergyDensity != 0) if(!fJetIDMVA->pass(pJet,iVertex,iVertices,iJetCorrector,iPileupEnergyDensity)) continue;
300 if(iJetCorrector == 0 || iPileupEnergyDensity == 0) if(!fJetIDMVA->pass(pJet,iVertex,iVertices)) continue;
301 if(iJetCorrector != 0 && iPileupEnergyDensity != 0) addNeut(pJet,lVec,lSumEt,iJetCorrector,iPileupEnergyDensity);
302 if(iJetCorrector == 0 || iPileupEnergyDensity == 0) addNeut(pJet,lVec,lSumEt,1.);
303 lNPass++;
304 }
305 Met lMet(lVec.Px(),lVec.Py());
306 lMet.SetSumEt( lSumEt);
307 return lMet;
308 }
309
310 //--------------------------------------------------------------------------------------------------
311 Met RecoilTools::NoPURecoil(Double_t iVisPt,Double_t iVisPhi,Double_t iVisSumEt,
312 const PFJetCol *iJets,FactorizedJetCorrector *iJetCorrector,
313 const PileupEnergyDensityCol *iPileupEnergyDensity,
314 const PFCandidateCol *iCands ,const Vertex *iVertex,
315 const VertexCol *iVertices,
316 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
317 Double_t iDZCut) {
318
319 Met lNoPUMet = NoPUMet(iJets,iJetCorrector,iPileupEnergyDensity,iCands,iVertex,iVertices,iPhi1,iEta1,iPhi2,iEta2,iDZCut);
320 lNoPUMet.SetMex (lNoPUMet.Mex()+iVisPt*cos(iVisPhi));
321 lNoPUMet.SetMey (lNoPUMet.Mey()+iVisPt*sin(iVisPhi));
322 lNoPUMet.SetSumEt(lNoPUMet.SumEt()-iVisSumEt);
323 return lNoPUMet;
324 }
325 //--------------------------------------------------------------------------------------------------
326 //Corrected Jets
327 Met RecoilTools::NoPURecoil(Double_t iVisPt,Double_t iVisPhi,Double_t iVisSumEt,
328 const PFJetCol *iJets,const PFCandidateCol *iCands,
329 const Vertex *iVertex,const VertexCol *iVertices,Double_t iRho,
330 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
331 Double_t iDZCut) {
332
333 Met lNoPUMet = NoPUMet(iJets,iCands,iVertex,iVertices,iRho,iPhi1,iEta1,iPhi2,iEta2,iDZCut);
334 lNoPUMet.SetMex (lNoPUMet.Mex()+iVisPt*cos(iVisPhi));
335 lNoPUMet.SetMey (lNoPUMet.Mey()+iVisPt*sin(iVisPhi));
336 lNoPUMet.SetSumEt(lNoPUMet.SumEt()-iVisSumEt);
337 return lNoPUMet;
338 }
339 //--------------------------------------------------------------------------------------------------
340 Met RecoilTools::PUCMet( const PFJetCol *iJets,FactorizedJetCorrector *iJetCorrector,
341 const PileupEnergyDensityCol *iPileupEnergyDensity,
342 const PFCandidateCol *iCands,
343 const Vertex *iVertex,const VertexCol *iVertices,
344 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
345 Double_t iDZCut) {
346
347 FourVectorM lVec (0,0,0,0); double lSumEt = 0;
348 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
349 const PFCandidate *pPF = iCands->At(i0);
350 const Track* pTrack = pPF->TrackerTrk();
351 if(pPF->GsfTrk()) pTrack = pPF->GsfTrk();
352 if(pTrack == 0 &&
353 (pPF->PFType() == PFCandidate::eGamma ||
354 pPF->PFType() == PFCandidate::eEGammaHF ||
355 pPF->PFType() == PFCandidate::eNeutralHadron ||
356 pPF->PFType() == PFCandidate::eHadronHF ))
357 {lVec -= pPF->Mom(); lSumEt += pPF->Pt();}
358 if(pTrack == 0 ) continue;
359 double lDZ = 999;
360 if(pPF->HasTrackerTrk()) lDZ = fabs(pPF->TrackerTrk()->DzCorrected(*iVertex));//
361 else if(pPF->HasGsfTrk()) lDZ = fabs(pPF->GsfTrk() ->DzCorrected(*iVertex));
362 if( fabs(lDZ) > iDZCut) continue;
363 lVec -= pPF->Mom();
364 lSumEt += pPF->Pt();
365 }
366 for(UInt_t i0 = 0; i0 < iJets->GetEntries(); i0++) {
367 const PFJet *pJet = iJets->At(i0);
368 if(!JetTools::passPFLooseId(pJet)) continue;
369 if(!fJetIDMVA->passPt(pJet,iJetCorrector,iPileupEnergyDensity)) continue;
370 if(fJetIDMVA->pass(pJet,iVertex,iVertices,iJetCorrector,iPileupEnergyDensity)) continue;
371 if(!filter(pJet,iPhi1,iEta1,iPhi2,iEta2)) continue;
372 addNeut(pJet,lVec,lSumEt,iJetCorrector,iPileupEnergyDensity,-1);
373 }
374 Met lMet(lVec.Px(),lVec.Py());
375 lMet.SetSumEt(lSumEt);
376 return lMet;
377 }
378 //--------------------------------------------------------------------------------------------------
379 //Corrected jets
380 Met RecoilTools::PUCMet( const PFJetCol *iJets,const PFCandidateCol *iCands,
381 const Vertex *iVertex,const VertexCol *iVertices,Double_t iRho,
382 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
383 Double_t iDZCut) {
384
385 FourVectorM lVec (0,0,0,0); double lSumEt = 0;
386 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
387 const PFCandidate *pPF = iCands->At(i0);
388 const Track* pTrack = pPF->TrackerTrk();
389 if(pPF->GsfTrk()) pTrack = pPF->GsfTrk();
390 if(pTrack == 0 &&
391 (pPF->PFType() == PFCandidate::eGamma ||
392 pPF->PFType() == PFCandidate::eEGammaHF ||
393 pPF->PFType() == PFCandidate::eNeutralHadron ||
394 pPF->PFType() == PFCandidate::eHadronHF ))
395 {lVec -= pPF->Mom(); lSumEt += pPF->Pt();}
396 if(pTrack == 0 ) continue;
397 double lDZ = 999;
398 if(pPF->HasTrackerTrk()) lDZ = fabs(pPF->TrackerTrk()->DzCorrected(*iVertex));//
399 else if(pPF->HasGsfTrk()) lDZ = fabs(pPF->GsfTrk() ->DzCorrected(*iVertex));
400 if( fabs(lDZ) > iDZCut) continue;
401 lVec -= pPF->Mom();
402 lSumEt += pPF->Pt();
403 }
404 for(UInt_t i0 = 0; i0 < iJets->GetEntries(); i0++) {
405 const PFJet *pJet = iJets->At(i0);
406 if(!JetTools::passPFLooseId(pJet)) continue;
407 if(!fJetIDMVA->passPt(pJet)) continue;
408 if(fJetIDMVA->pass(pJet,iVertex,iVertices)) continue;
409 if(!filter(pJet,iPhi1,iEta1,iPhi2,iEta2)) continue; //Quick cleaning==> if not done already
410 addNeut(pJet,lVec,lSumEt,iRho,-1);
411 }
412 Met lMet(lVec.Px(),lVec.Py());
413 lMet.SetSumEt(lSumEt);
414 return lMet;
415 }
416 //--------------------------------------------------------------------------------------------------
417 Met RecoilTools::PUCRecoil( Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
418 const PFJetCol *iJets,
419 const PFCandidateCol *iCands ,
420 const Vertex *iVertex,const VertexCol *iVertices,
421 FactorizedJetCorrector *iJetCorrector,
422 const PileupEnergyDensityCol *iPileupEnergyDensity,Double_t iDZCut) {
423
424
425 FourVectorM lVec (0,0,0,0); double lSumEt = 0;
426 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
427 const PFCandidate *pPF = iCands->At(i0);
428 if(!filter(pPF,iPhi1,iEta1,iPhi2,iEta2)) continue;
429 const Track* pTrack = pPF->TrackerTrk();
430 if(pPF->GsfTrk()) pTrack = pPF->GsfTrk();
431 if(pTrack == 0 &&
432 (pPF->PFType() == PFCandidate::eGamma ||
433 pPF->PFType() == PFCandidate::eEGammaHF ||
434 pPF->PFType() == PFCandidate::eNeutralHadron ||
435 pPF->PFType() == PFCandidate::eHadronHF ))
436 {lVec -= pPF->Mom(); lSumEt += pPF->Pt();}
437 if(pTrack == 0 ) continue;
438 double lDZ = 999;
439 if(pPF->HasTrackerTrk()) lDZ = fabs(pPF->TrackerTrk()->DzCorrected(*iVertex));//
440 else if(pPF->HasGsfTrk()) lDZ = fabs(pPF->GsfTrk() ->DzCorrected(*iVertex));
441 if( fabs(lDZ) > iDZCut) continue;
442 lVec -= pPF->Mom();
443 lSumEt += pPF->Pt();
444 }
445 for(UInt_t i0 = 0; i0 < iJets->GetEntries(); i0++) {
446 const PFJet *pJet = iJets->At(i0);
447 if(!filter(pJet,iPhi1,iEta1,iPhi2,iEta2)) continue; //Quick cleaning==> if not done already
448 if(!fJetIDMVA->passPt(pJet,iJetCorrector,iPileupEnergyDensity)) continue;
449 if(!JetTools::passPFLooseId(pJet)) continue;
450 if(iJetCorrector != 0 && iPileupEnergyDensity != 0) if(fJetIDMVA->pass(pJet,iVertex,iVertices,iJetCorrector,iPileupEnergyDensity)) continue;
451 if(iJetCorrector == 0 || iPileupEnergyDensity == 0) if(fJetIDMVA->pass(pJet,iVertex,iVertices)) continue;
452 if(iJetCorrector != 0 && iPileupEnergyDensity != 0) addNeut(pJet,lVec,lSumEt,iJetCorrector,iPileupEnergyDensity,-1);
453 if(iJetCorrector == 0 || iPileupEnergyDensity == 0) addNeut(pJet,lVec,lSumEt,1,-1);
454 }
455 Met lMet(lVec.Px(),lVec.Py());
456 lMet.SetSumEt(lSumEt);
457 return lMet;
458 }
459 //----> This MET is a bug need to fix it
460 //--------------------------------------------------------------------------------------------------
461 Met RecoilTools::PUCRecoil(Double_t iVisPt,Double_t iVisPhi,Double_t iVisSumEt,
462 const PFJetCol *iJets,FactorizedJetCorrector *iJetCorrector,
463 const PileupEnergyDensityCol *iPileupEnergyDensity,
464 const PFCandidateCol *iCands,const Vertex *iVertex,
465 const VertexCol *iVertices,
466 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
467 Double_t iDZCut) {
468 Met lPUCMet = PUCMet(iJets,iJetCorrector,iPileupEnergyDensity,iCands,iVertex,iVertices,iPhi1,iEta1,iPhi2,iEta2,iDZCut);
469 lPUCMet.SetMex (lPUCMet.Mex()+iVisPt*cos(iVisPhi));
470 lPUCMet.SetMey (lPUCMet.Mey()+iVisPt*sin(iVisPhi));
471 lPUCMet.SetSumEt(lPUCMet.SumEt()-iVisSumEt);
472 return lPUCMet;
473 }
474 //--------------------------------------------------------------------------------------------------
475 //Corrected Jets
476 Met RecoilTools::PUCRecoil(Double_t iVisPt,Double_t iVisPhi,Double_t iVisSumEt,
477 const PFJetCol *iJets,
478 const PFCandidateCol *iCands,const Vertex *iVertex,
479 const VertexCol *iVertices,Double_t iRho,
480 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
481 Double_t iDZCut) {
482 Met lPUCMet = PUCMet(iJets,iCands,iVertex,iVertices,iRho,iPhi1,iEta1,iPhi2,iEta2,iDZCut);
483 lPUCMet.SetMex (lPUCMet.Mex()+iVisPt*cos(iVisPhi));
484 lPUCMet.SetMey (lPUCMet.Mey()+iVisPt*sin(iVisPhi));
485 lPUCMet.SetSumEt(lPUCMet.SumEt()-iVisSumEt);
486 return lPUCMet;
487 }
488 //--------------------------------------------------------------------------------------------------
489 Met RecoilTools::PUMet( const PFJetCol *iJets,FactorizedJetCorrector *iJetCorrector,
490 const PileupEnergyDensityCol *iPileupEnergyDensity,
491 const PFCandidateCol *iCands,const Vertex *iVertex,
492 const VertexCol *iVertices,
493 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
494 Double_t iDZCut) {
495
496 FourVectorM lVec (0,0,0,0); double lSumEt = 0;
497 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
498 const PFCandidate *pPF = iCands->At(i0);
499 const Track* pTrack = pPF->TrackerTrk();
500 if(pPF->GsfTrk()) pTrack = pPF->GsfTrk();
501 if(pTrack == 0 ) continue;
502 double lDZ = 999;
503 if(pPF->HasTrackerTrk()) lDZ = fabs(pPF->TrackerTrk()->DzCorrected(*iVertex));//
504 else if(pPF->HasGsfTrk()) lDZ = fabs(pPF->GsfTrk() ->DzCorrected(*iVertex));
505 if( fabs(lDZ) < iDZCut) continue;
506 lVec -= pPF->Mom();
507 lSumEt += pPF->Pt();
508 }
509 for(UInt_t i0 = 0; i0 < iJets->GetEntries(); i0++) {
510 const PFJet *pJet = iJets->At(i0);
511 if(!JetTools::passPFLooseId(pJet)) continue;
512 if(!fJetIDMVA->passPt(pJet,iJetCorrector,iPileupEnergyDensity)) continue;
513 if(!filter(pJet,iPhi1,iEta1,iPhi2,iEta2)) continue;
514 if(fJetIDMVA->pass(pJet,iVertex,iVertices,iJetCorrector,iPileupEnergyDensity)) continue;
515 addNeut(pJet,lVec,lSumEt,iJetCorrector,iPileupEnergyDensity);
516 }
517 Met lMet(lVec.Px(),lVec.Py());
518 lMet.SetSumEt(lSumEt);
519 return lMet;
520 }
521 //--------------------------------------------------------------------------------------------------
522 //Corrected Jets
523 Met RecoilTools::PUMet( const PFJetCol *iJets,
524 const PFCandidateCol *iCands,const Vertex *iVertex,const VertexCol *iVertices,Double_t iRho,
525 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
526 Double_t iDZCut) {
527
528 FourVectorM lVec (0,0,0,0); double lSumEt = 0;
529 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
530 const PFCandidate *pPF = iCands->At(i0);
531 const Track* pTrack = pPF->TrackerTrk();
532 if(pPF->GsfTrk()) pTrack = pPF->GsfTrk();
533 if(pTrack == 0 ) continue;
534 double lDZ = 999;
535 if(pPF->HasTrackerTrk()) lDZ = fabs(pPF->TrackerTrk()->DzCorrected(*iVertex));//
536 else if(pPF->HasGsfTrk()) lDZ = fabs(pPF->GsfTrk() ->DzCorrected(*iVertex));
537 if( fabs(lDZ) < iDZCut) continue;
538 lVec -= pPF->Mom();
539 lSumEt += pPF->Pt();
540 }
541 for(UInt_t i0 = 0; i0 < iJets->GetEntries(); i0++) {
542 const PFJet *pJet = iJets->At(i0);
543 if(!JetTools::passPFLooseId(pJet)) continue;
544 if(!fJetIDMVA->passPt(pJet)) continue;
545 if(!filter(pJet,iPhi1,iEta1,iPhi2,iEta2)) continue; //Quick cleaning
546 if(fJetIDMVA->pass(pJet,iVertex,iVertices)) continue;
547 addNeut(pJet,lVec,lSumEt,iRho);
548 }
549 Met lMet(lVec.Px(),lVec.Py());
550 lMet.SetSumEt(lSumEt);
551 return lMet;
552 }
553 //--------------------------------------------------------------------------------------------------
554 Met RecoilTools::PUMet( Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
555 const PFJetCol *iJets,
556 const PFCandidateCol *iCands ,
557 const Vertex *iVertex,const VertexCol *iVertices,
558 FactorizedJetCorrector *iJetCorrector,
559 const PileupEnergyDensityCol *iPileupEnergyDensity,Double_t iDZCut) {
560
561 FourVectorM lVec (0,0,0,0); double lSumEt = 0;
562 for(UInt_t i0 = 0; i0 < iCands->GetEntries(); i0++) {
563 const PFCandidate *pPF = iCands->At(i0);
564 if(!filter(pPF,iPhi1,iEta1,iPhi2,iEta2)) continue;
565 const Track* pTrack = pPF->TrackerTrk();
566 if(pPF->GsfTrk()) pTrack = pPF->GsfTrk();
567 if(pTrack == 0 ) continue;
568 double lDZ = 999;
569 if(pPF->HasTrackerTrk()) lDZ = fabs(pPF->TrackerTrk()->DzCorrected(*iVertex));//
570 else if(pPF->HasGsfTrk()) lDZ = fabs(pPF->GsfTrk() ->DzCorrected(*iVertex));
571 if( fabs(lDZ) < iDZCut) continue;
572 lVec -= pPF->Mom();
573 lSumEt += pPF->Pt();
574 }
575 for(UInt_t i0 = 0; i0 < iJets->GetEntries(); i0++) {
576 const PFJet *pJet = iJets->At(i0);
577 if(!JetTools::passPFLooseId(pJet)) continue;
578 if(!fJetIDMVA->passPt(pJet,iJetCorrector,iPileupEnergyDensity)) continue;
579 if(!filter(pJet,iPhi1,iEta1,iPhi2,iEta2)) continue; //Quick cleaning
580 if(fJetIDMVA->pass(pJet,iVertex,iVertices,iJetCorrector,iPileupEnergyDensity)) continue;
581 addNeut(pJet,lVec,lSumEt,iJetCorrector,iPileupEnergyDensity);
582 }
583 Met lMet(lVec.Px(),lVec.Py());
584 lMet.SetSumEt(lSumEt);
585 return lMet;
586 }
587 //--------------------------------------------------------------------------------------------------
588 Met RecoilTools::PURecoil(Double_t iVisPt,Double_t iVisPhi,Double_t iVisSumEt,
589 const PFJetCol *iJets,FactorizedJetCorrector *iJetCorrector,
590 const PileupEnergyDensityCol *iPileupEnergyDensity,
591 const PFCandidateCol *iCands,const Vertex *iVertex,
592 const VertexCol *iVertices,
593 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
594 Double_t iDZCut) {
595 Met lPUMet = PUMet(iJets,iJetCorrector,iPileupEnergyDensity,iCands,iVertex,iVertices,iPhi1,iEta1,iPhi2,iEta2,iDZCut);
596 lPUMet.SetMex (lPUMet.Mex()+iVisPt*cos(iVisPhi));
597 lPUMet.SetMey (lPUMet.Mey()+iVisPt*sin(iVisPhi));
598 lPUMet.SetSumEt(lPUMet.SumEt()-iVisSumEt);
599 return lPUMet;
600 }
601 //--------------------------------------------------------------------------------------------------
602 //Corrected Jets
603 Met RecoilTools::PURecoil(Double_t iVisPt,Double_t iVisPhi,Double_t iVisSumEt,
604 const PFJetCol *iJets,
605 const PFCandidateCol *iCands,const Vertex *iVertex,const VertexCol *iVertices,Double_t iRho,
606 Double_t iPhi1,Double_t iEta1,Double_t iPhi2,Double_t iEta2,
607 Double_t iDZCut) {
608 Met lPUMet = PUMet(iJets,iCands,iVertex,iVertices,iRho,iPhi1,iEta1,iPhi2,iEta2,iDZCut);
609 lPUMet.SetMex (lPUMet.Mex()+iVisPt*cos(iVisPhi));
610 lPUMet.SetMey (lPUMet.Mey()+iVisPt*sin(iVisPhi));
611 lPUMet.SetSumEt(lPUMet.SumEt()-iVisSumEt);
612 return lPUMet;
613 }