ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Utils/src/MVAMet.cc
Revision: 1.14
Committed: Tue May 1 16:53:38 2012 UTC (13 years ago) by pharris
Content type: text/plain
Branch: MAIN
Changes since 1.13: +42 -6 lines
Log Message:
Fixed Castor Jet Bug

File Contents

# User Rev Content
1 pharris 1.1 #include "MitPhysics/Utils/interface/MVAMet.h"
2     #include "MitPhysics/Utils/interface/JetTools.h"
3     #include "MitPhysics/Utils/interface/RecoilTools.h"
4     #include "MitAna/DataTree/interface/StableData.h"
5     #include <TFile.h>
6     #include <TRandom3.h>
7 bendavid 1.6 #include "CondFormats/EgammaObjects/interface/GBRForest.h"
8     #include "Cintex/Cintex.h"
9 pharris 1.1
10     ClassImp(mithep::MVAMet)
11    
12     using namespace mithep;
13    
14     //--------------------------------------------------------------------------------------------------
15     MVAMet::MVAMet() :
16     fRecoilTools(0),
17 pharris 1.4 fPhiMethodName ("PhiCorrection"),
18     fU1MethodName ("U1Correction"),
19     fCovU1MethodName ("CovU1"),
20     fCovU2MethodName ("CovU2"),
21 pharris 1.1 fIsInitialized(kFALSE),
22     fU (0),
23     fUPhi (0),
24     fTKSumEt(0),
25     fTKU (0),
26     fTKUPhi (0),
27     fNPSumEt(0),
28     fNPU (0),
29     fNPUPhi (0),
30     fPUSumEt(0),
31     fPUMet (0),
32     fPUMetPhi(0),
33     fPCSumEt(0),
34     fPCU (0),
35     fPCUPhi (0),
36     fJSPt1 (0),
37     fJSEta1 (0),
38     fJSPhi1 (0),
39     fJSPt2 (0),
40     fJSEta2 (0),
41     fJSPhi2 (0),
42     fNJet (0),
43     fNAllJet(0),
44     fNPV (0),
45     fPhiReader(0),
46 pharris 1.4 fU1Reader(0),
47     fCovU1Reader(0),
48     fCovU2Reader(0) { }
49 pharris 1.1 //--------------------------------------------------------------------------------------------------
50     MVAMet::~MVAMet()
51     {
52 pharris 1.4 delete fPhiReader;
53     delete fU1Reader;
54     delete fCovU1Reader;
55     delete fCovU2Reader;
56    
57 pharris 1.1 }
58     //--------------------------------------------------------------------------------------------------
59 pharris 1.2 /*
60 pharris 1.1 void MVAMet::setVariables(TMVA::Reader *iReader,bool iScale) {
61     iReader->AddVariable( "npv" , &fNPV );
62     iReader->AddVariable( "u" , &fU );
63     iReader->AddVariable( "uphi" , &fUPhi );
64     iReader->AddVariable( "chsumet/sumet" , &fTKSumEt );
65     iReader->AddVariable( "tku" , &fTKU );
66     iReader->AddVariable( "tkuphi" , &fTKUPhi );
67     iReader->AddVariable( "nopusumet/sumet" , &fNPSumEt );
68     iReader->AddVariable( "nopuu" , &fNPU );
69     iReader->AddVariable( "nopuuphi" , &fNPUPhi );
70     iReader->AddVariable( "pusumet/sumet" , &fPUSumEt );
71     iReader->AddVariable( "pumet" , &fPUMet );
72     iReader->AddVariable( "pumetphi" , &fPUMetPhi);
73     iReader->AddVariable( "pucsumet/sumet" , &fPCSumEt );
74     iReader->AddVariable( "pucu" , &fPCU );
75     iReader->AddVariable( "pucuphi" , &fPCUPhi );
76     iReader->AddVariable( "jspt_1" , &fJSPt1 );
77     iReader->AddVariable( "jseta_1" , &fJSEta1 );
78     iReader->AddVariable( "jsphi_1" , &fJSPhi1 );
79     iReader->AddVariable( "jspt_2" , &fJSPt2 );
80     iReader->AddVariable( "jseta_2" , &fJSEta2 );
81     iReader->AddVariable( "jsphi_2" , &fJSPhi2 );
82     iReader->AddVariable( "nalljet" , &fNAllJet );
83     iReader->AddVariable( "njet" , &fNJet );
84     if(iScale) iReader->AddVariable( "uphi_mva" , &fUPhiMVA );
85     }
86     //--------------------------------------------------------------------------------------------------
87     void MVAMet::Initialize( TString iU1MethodName,
88     TString iPhiMethodName,
89     TString iJetMVAFile,
90     TString iU1Weights,
91     TString iPhiWeights,
92     MVAMet::MVAType iType) {
93    
94     fIsInitialized = kTRUE;
95     fRecoilTools = new RecoilTools(iJetMVAFile);
96 pharris 1.2
97 pharris 1.1 fU1MethodName = iU1MethodName;
98     fPhiMethodName = iPhiMethodName;
99     fType = iType;
100     fPhiReader = new TMVA::Reader( "!Color:!Silent:Error" );
101     fU1Reader = new TMVA::Reader( "!Color:!Silent:Error" );
102     if (fType == kBaseline) {
103     setVariables(fPhiReader,false);
104     setVariables(fU1Reader ,true);
105     }
106     fPhiReader->BookMVA(fPhiMethodName , iPhiWeights );
107     fU1Reader ->BookMVA(fU1MethodName , iU1Weights );
108 pharris 1.2
109 pharris 1.1 std::cout << "Jet ID MVA Initialization\n";
110     std::cout << "Phi Method Name : " << fPhiMethodName << " , type == " << iType << std::endl;
111     std::cout << "U1 Method Name : " << fU1MethodName << " , type == " << iType << std::endl;
112 pharris 1.2
113    
114 pharris 1.1 }
115 pharris 1.2 */
116     //--------------------------------------------------------------------------------------------------
117 pharris 1.3 void MVAMet::Initialize(TString iJetLowPtFile,
118     TString iJetHighPtFile,
119     TString iJetCutFile,
120     TString iU1Weights,
121     TString iPhiWeights,
122 pharris 1.4 TString iCovU1Weights,
123     TString iCovU2Weights,
124 pharris 1.3 MVAMet::MVAType iType) {
125 pharris 1.2
126     fIsInitialized = kTRUE;
127 pharris 1.3 fRecoilTools = new RecoilTools(iJetLowPtFile,iJetHighPtFile,iJetCutFile);
128    
129 pharris 1.2 fType = iType;
130    
131 bendavid 1.6 ROOT::Cintex::Cintex::Enable();
132    
133 pharris 1.2 TFile *lPhiForest = new TFile(iPhiWeights,"READ");
134     fPhiReader = (GBRForest*)lPhiForest->Get(fPhiMethodName);
135     lPhiForest->Close();
136    
137     TFile *lU1Forest = new TFile(iU1Weights,"READ");
138     fU1Reader = (GBRForest*)lU1Forest->Get(fU1MethodName);
139     lU1Forest->Close();
140 pharris 1.1
141 pharris 1.4 TFile *lCovU1Forest = new TFile(iCovU1Weights,"READ");
142     fCovU1Reader = (GBRForest*)lCovU1Forest->Get(fCovU1MethodName);
143     lCovU1Forest->Close();
144    
145     TFile *lCovU2Forest = new TFile(iCovU2Weights,"READ");
146     fCovU2Reader = (GBRForest*)lCovU2Forest->Get(fCovU2MethodName);
147     lCovU2Forest->Close();
148 pharris 1.5
149     fCov = new TMatrixD(2,2);
150 pharris 1.7 fPhiVals = new Float_t[23];
151     fU1Vals = new Float_t[25];
152     fCovVals = new Float_t[26];
153 pharris 1.2 }
154     //--------------------------------------------------------------------------------------------------
155     Double_t MVAMet::evaluatePhi() {
156     fPhiVals[0] = fNPV ;
157     fPhiVals[1] = fU ;
158     fPhiVals[2] = fUPhi ;
159     fPhiVals[3] = fTKSumEt ;
160     fPhiVals[4] = fTKU ;
161     fPhiVals[5] = fTKUPhi ;
162     fPhiVals[6] = fNPSumEt ;
163     fPhiVals[7] = fNPU ;
164     fPhiVals[8] = fNPUPhi ;
165     fPhiVals[9] = fPUSumEt ;
166     fPhiVals[10] = fPUMet ;
167     fPhiVals[11] = fPUMetPhi;
168     fPhiVals[12] = fPCSumEt ;
169     fPhiVals[13] = fPCU ;
170     fPhiVals[14] = fPCUPhi ;
171     fPhiVals[15] = fJSPt1 ;
172     fPhiVals[16] = fJSEta1 ;
173     fPhiVals[17] = fJSPhi1 ;
174     fPhiVals[18] = fJSPt2 ;
175     fPhiVals[19] = fJSEta2 ;
176     fPhiVals[20] = fJSPhi2 ;
177     fPhiVals[21] = fNAllJet ;
178     fPhiVals[22] = fNJet ;
179     return fPhiReader->GetResponse(fPhiVals);
180     }
181     //--------------------------------------------------------------------------------------------------
182     Double_t MVAMet::evaluateU1() {
183 pharris 1.5 fU1Vals[0] = fSumEt ;
184     fU1Vals[1] = fNPV ;
185     fU1Vals[2] = fU ;
186     fU1Vals[3] = fUPhi ;
187     fU1Vals[4] = fTKSumEt ;
188     fU1Vals[5] = fTKU ;
189     fU1Vals[6] = fTKUPhi ;
190     fU1Vals[7] = fNPSumEt ;
191     fU1Vals[8] = fNPU ;
192     fU1Vals[9] = fNPUPhi ;
193     fU1Vals[10] = fPUSumEt ;
194     fU1Vals[11] = fPUMet ;
195     fU1Vals[12] = fPUMetPhi;
196     fU1Vals[13] = fPCSumEt ;
197     fU1Vals[14] = fPCU ;
198     fU1Vals[15] = fPCUPhi ;
199     fU1Vals[16] = fJSPt1 ;
200     fU1Vals[17] = fJSEta1 ;
201     fU1Vals[18] = fJSPhi1 ;
202     fU1Vals[19] = fJSPt2 ;
203     fU1Vals[20] = fJSEta2 ;
204     fU1Vals[21] = fJSPhi2 ;
205     fU1Vals[22] = fNAllJet ;
206     fU1Vals[23] = fNJet ;
207     fU1Vals[24] = fUPhiMVA ;
208 pharris 1.2 return fU1Reader->GetResponse(fU1Vals);
209     }
210 pharris 1.1 //--------------------------------------------------------------------------------------------------
211 pharris 1.4 Double_t MVAMet::evaluateCovU1() {
212 pharris 1.5 fCovVals[0] = fSumEt ;
213     fCovVals[1] = fNPV ;
214     fCovVals[2] = fU ;
215     fCovVals[3] = fUPhi ;
216     fCovVals[4] = fTKSumEt ;
217     fCovVals[5] = fTKU ;
218     fCovVals[6] = fTKUPhi ;
219     fCovVals[7] = fNPSumEt ;
220     fCovVals[8] = fNPU ;
221     fCovVals[9] = fNPUPhi ;
222     fCovVals[10] = fPUSumEt ;
223     fCovVals[11] = fPUMet ;
224     fCovVals[12] = fPUMetPhi;
225     fCovVals[13] = fPCSumEt ;
226     fCovVals[14] = fPCU ;
227     fCovVals[15] = fPCUPhi ;
228     fCovVals[16] = fJSPt1 ;
229     fCovVals[17] = fJSEta1 ;
230     fCovVals[18] = fJSPhi1 ;
231     fCovVals[19] = fJSPt2 ;
232     fCovVals[20] = fJSEta2 ;
233     fCovVals[21] = fJSPhi2 ;
234     fCovVals[22] = fNAllJet ;
235     fCovVals[23] = fNJet ;
236     fCovVals[24] = fUPhiMVA ;
237     fCovVals[25] = fUMVA ;
238 pharris 1.4 return fCovU1Reader->GetResponse(fCovVals);
239     }
240     //--------------------------------------------------------------------------------------------------
241     Double_t MVAMet::evaluateCovU2() {
242 pharris 1.5 fCovVals[0] = fSumEt ;
243     fCovVals[1] = fNPV ;
244     fCovVals[2] = fU ;
245     fCovVals[3] = fUPhi ;
246     fCovVals[4] = fTKSumEt ;
247     fCovVals[5] = fTKU ;
248     fCovVals[6] = fTKUPhi ;
249     fCovVals[7] = fNPSumEt ;
250     fCovVals[8] = fNPU ;
251     fCovVals[9] = fNPUPhi ;
252     fCovVals[10] = fPUSumEt ;
253     fCovVals[11] = fPUMet ;
254     fCovVals[12] = fPUMetPhi;
255     fCovVals[13] = fPCSumEt ;
256     fCovVals[14] = fPCU ;
257     fCovVals[15] = fPCUPhi ;
258     fCovVals[16] = fJSPt1 ;
259     fCovVals[17] = fJSEta1 ;
260     fCovVals[18] = fJSPhi1 ;
261     fCovVals[19] = fJSPt2 ;
262     fCovVals[20] = fJSEta2 ;
263     fCovVals[21] = fJSPhi2 ;
264     fCovVals[22] = fNAllJet ;
265     fCovVals[23] = fNJet ;
266     fCovVals[24] = fUPhiMVA ;
267     fCovVals[25] = fUMVA ;
268 pharris 1.4 return fCovU2Reader->GetResponse(fCovVals);
269     }
270     //--------------------------------------------------------------------------------------------------
271 pharris 1.1 Double_t MVAMet::MVAValue( Bool_t iPhi,
272     Float_t iPFSumEt,
273     Float_t iU ,
274     Float_t iUPhi ,
275     Float_t iTKSumEt,
276     Float_t iTKU ,
277     Float_t iTKUPhi ,
278     Float_t iNPSumEt,
279     Float_t iNPU ,
280     Float_t iNPUPhi ,
281     Float_t iPUSumEt,
282     Float_t iPUMet ,
283     Float_t iPUMetPhi,
284     Float_t iPCSumEt,
285     Float_t iPCU ,
286     Float_t iPCUPhi ,
287     Float_t iJSPt1 ,
288     Float_t iJSEta1 ,
289     Float_t iJSPhi1 ,
290     Float_t iJSPt2 ,
291     Float_t iJSEta2 ,
292     Float_t iJSPhi2 ,
293     Float_t iNAllJet,
294     Float_t iNJet ,
295     Float_t iNPV
296     ){
297    
298     if (!fIsInitialized) {
299     std::cout << "Error: MVA Met not properly initialized.\n";
300     return -9999;
301     }
302 pharris 1.5
303     fSumEt = iPFSumEt;
304 pharris 1.1 fU = iU ;
305     fUPhi = iUPhi ;
306     fTKSumEt = iTKSumEt/iPFSumEt;
307     fTKU = iTKU ;
308     fTKUPhi = iTKUPhi ;
309     fNPSumEt = iNPSumEt/iPFSumEt;
310     fNPU = iNPU ;
311     fNPUPhi = iNPUPhi ;
312     fPUSumEt = iPUSumEt/iPFSumEt;
313     fPUMet = iPUMet ;
314     fPUMetPhi = iPUMetPhi;
315     fPCSumEt = iPCSumEt/iPFSumEt;
316     fPCU = iPCU ;
317     fPCUPhi = iPCUPhi ;
318     fJSPt1 = iJSPt1 ;
319     fJSEta1 = iJSEta1 ;
320     fJSPhi1 = iJSPhi1 ;
321     fJSPt2 = iJSPt2 ;
322     fJSEta2 = iJSEta2 ;
323     fJSPhi2 = iJSPhi2 ;
324     fNAllJet = iNAllJet;
325     fNJet = iNJet ;
326     fNPV = iNPV ;
327 pharris 1.5
328 pharris 1.1 Double_t lMVA = -9999;
329 pharris 1.2 lMVA = evaluatePhi();
330 pharris 1.1 if(!iPhi) fUPhiMVA = iUPhi+lMVA;
331 pharris 1.2 //Not no nice feature of the training
332 pharris 1.4 //fTKSumEt /= iPFSumEt;
333     //fNPSumEt /= iPFSumEt;
334     //fPUSumEt /= iPFSumEt;
335     //fPCSumEt /= iPFSumEt;
336 pharris 1.2 if(!iPhi) lMVA = evaluateU1();
337 pharris 1.1 return lMVA;
338     }
339     //--------------------------------------------------------------------------------------------------
340     //====> Please not that the jet collection must be cleaned => all jets near leptons must be removed
341     Met MVAMet::GetMet( Bool_t iPhi,Float_t iPtVis,Float_t iPhiVis,Float_t iSumEtVis,
342     const PFMet *iMet ,
343 pharris 1.3 const PFCandidateCol *iCands,
344     const Vertex *iVertex,const VertexCol *iVertices,
345 pharris 1.1 const PFJetCol *iJets ,
346     FactorizedJetCorrector *iJetCorrector,
347     const PileupEnergyDensityCol *iPUEnergyDensity,
348     int iNPV,
349     Bool_t printDebug) {
350 pharris 1.14
351     int lNPV = 0;
352     for(unsigned int i0 = 0; i0 < iVertices->GetEntries(); i0++) {
353     const Vertex *lPV = iVertices->At(i0);
354     if(lPV->Ndof() < 4.0) continue;
355     if(fabs(lPV->Z()) > 24.) continue;
356     if(lPV->Position().Rho() > 2.) continue;
357     lNPV++;
358     }
359 pharris 1.1
360 pharris 1.7 Met lPFRec = fRecoilTools->pfRecoil (iPtVis,iPhiVis,iSumEtVis,iCands);
361 pharris 1.1 Met lTKRec = fRecoilTools->trackRecoil(iPtVis,iPhiVis,iSumEtVis,iCands,iVertex);
362 pharris 1.3 Met lNPRec = fRecoilTools->NoPURecoil (iPtVis,iPhiVis,iSumEtVis,iJets,iJetCorrector,iPUEnergyDensity,iCands,iVertex,iVertices);
363     Met lPCRec = fRecoilTools->PUCRecoil (iPtVis,iPhiVis,iSumEtVis,iJets,iJetCorrector,iPUEnergyDensity,iCands,iVertex,iVertices);
364     Met lPUMet = fRecoilTools->PUMet (iJets,iJetCorrector,iPUEnergyDensity,iCands,iVertex,iVertices);
365 pharris 1.1
366     Double_t lPt0 = 0; const PFJet *lLead = 0;
367     Double_t lPt1 = 0; const PFJet *l2nd = 0;
368     int lNAllJet = 0;
369     int lNJet = 0;
370     for(unsigned int i0 = 0; i0 < iJets->GetEntries(); i0++) {
371     const PFJet *pJet = iJets->At(i0);
372     Double_t pPt = fRecoilTools->fJetIDMVA->correctedPt(pJet,iJetCorrector,iPUEnergyDensity);
373     if(!JetTools::passPFLooseId(pJet)) continue;
374     lNAllJet++;
375     if(pPt > 30) lNJet++;
376     if(lPt0 < pPt) {lPt0 = pPt; lLead = pJet; continue;}
377     if(lPt1 < pPt) {lPt1 = pPt; l2nd = pJet; continue;}
378     }
379 pharris 1.5
380     fSumEt = lPFRec.SumEt();
381 pharris 1.1 fU = lPFRec.Pt();
382     fUPhi = lPFRec.Phi();
383     fTKSumEt = lTKRec.SumEt()/lPFRec.SumEt();
384     fTKU = lTKRec.Pt();
385     fTKUPhi = lTKRec.Phi();
386     fNPSumEt = lNPRec.SumEt()/lPFRec.SumEt();
387     fNPU = lNPRec.Pt();
388     fNPUPhi = lNPRec.Phi();
389     fPUSumEt = lPUMet.SumEt()/lPFRec.SumEt();
390     fPUMet = lPUMet.Pt();
391     fPUMetPhi = lPUMet.Phi();
392     fPCSumEt = lPCRec.SumEt()/lPFRec.SumEt();
393     fPCU = lPCRec.Pt() ;
394     fPCUPhi = lPCRec.Phi() ;
395     fJSPt1 = lPt0;
396     fJSEta1 = 0; if(lLead != 0) fJSEta1 = lLead->Eta();
397     fJSPhi1 = 0; if(lLead != 0) fJSPhi1 = lLead->Phi();
398     fJSPt2 = lPt1;
399     fJSEta2 = 0; if(l2nd != 0) fJSEta2 = l2nd ->Eta();
400     fJSPhi2 = 0; if(l2nd != 0) fJSPhi2 = l2nd ->Phi();
401     fNJet = lNJet ;
402     fNAllJet = lNAllJet;
403 pharris 1.14 fNPV = lNPV ;
404 pharris 1.5 Float_t lMVA = evaluatePhi();
405     if(!iPhi) fUPhiMVA = fUPhi + lMVA;
406 pharris 1.1 //Not no nice feature of teh training
407 pharris 1.4 //fTKSumEt /= lPFRec.SumEt();
408     //fNPSumEt /= lPFRec.SumEt();
409     //fPUSumEt /= lPFRec.SumEt();
410     //fPCSumEt /= lPFRec.SumEt();
411 pharris 1.7 if(!iPhi) lMVA = evaluateU1();//fU1Reader ->EvaluateMVA( fU1MethodName );
412     fUMVA = fU*lMVA;
413    
414     TLorentzVector lUVec(0,0,0,0); lUVec .SetPtEtaPhiM(fUMVA,0,fUPhiMVA,0);
415 pharris 1.1 TLorentzVector lVVec(0,0,0,0); lVVec .SetPtEtaPhiM(iPtVis ,0,iPhiVis ,0);
416     if(lMVA < 0) lUVec .RotateZ(TMath::Pi());
417     lUVec -= lVVec;
418     Met lMet(lUVec.Px(),lUVec.Py());
419 pharris 1.5
420 pharris 1.4 //TMatrixD lCov(2,2);
421     //Covariance matrix perpendicular and parallel to the recoil (sort of correct)
422     double lCovU1 = evaluateCovU1();
423     double lCovU2 = evaluateCovU2();
424    
425     double lCos2 = cos(fUPhiMVA)*cos(fUPhiMVA);
426     double lSin2 = sin(fUPhiMVA)*sin(fUPhiMVA);
427 pharris 1.5
428 pharris 1.4 //Now Compute teh covariance matrix in X and Y
429 pharris 1.5 (*fCov)(0,0) = lCovU1*lCos2+lCovU2*lSin2;
430     (*fCov)(1,0) = -lCovU1*sin(fUPhiMVA)*cos(fUPhiMVA)+lCovU2*sin(fUPhiMVA)*cos(fUPhiMVA);
431     (*fCov)(0,1) = (*fCov)(1,0);
432     (*fCov)(1,1) = lCovU1*lSin2+lCovU2*lCos2;
433 pharris 1.10 TMatrixD lInv(2,2);
434     lInv(0,0) = (*fCov)(0,0); lInv(1,1) = (*fCov)(1,1); lInv(1,0) = (*fCov)(1,0); lInv(0,1) = (*fCov)(0,1);
435     if(lInv.Determinant() != 0) lInv.Invert();
436     fSignificance = TMath::Sqrt(lUVec.Px()*lUVec.Px()*(lInv)(0,0) + 2.*lUVec.Px()*lUVec.Py()*(lInv)(1,0) + lUVec.Py()*lUVec.Py()*(lInv)(1,1));
437     fUncertainty = sqrt(lCovU1+lCovU2);
438 pharris 1.7
439 pharris 1.1 if (printDebug == kTRUE) {
440 ceballos 1.11 std::cout << "Debug Met MVA: "
441 pharris 1.1 << fU << " : "
442     << fUPhi << " : "
443     << fTKSumEt << " : "
444     << fTKU << " : "
445     << fTKUPhi << " : "
446     << fNPSumEt << " : "
447     << fNPU << " : "
448     << fNPUPhi << " : "
449     << fPUSumEt << " : "
450     << fPUMet << " : "
451     << fPUMetPhi << " : "
452     << fPCUPhi << " : "
453     << fPCSumEt << " : "
454     << fPCU << " : "
455     << fPCUPhi << " : "
456     << fJSPt1 << " : "
457     << fJSEta1 << " : "
458     << fJSPhi1 << " : "
459     << fJSPt2 << " : "
460     << fJSEta2 << " : "
461     << fJSPhi2 << " : "
462     << fNJet << " : "
463     << fNAllJet << " : "
464     << fNPV << " : "
465     << " === : === "
466 ceballos 1.12 << lMet.Pt() << " : "
467     << lMet.Phi() << " : "
468 pharris 1.1 << std::endl;
469     }
470    
471     return lMet;
472     }
473 pharris 1.2 //--------------------------------------------------------------------------------------------------
474     //====> Please not that the jet collection must be cleaned => all jets near leptons must be removed
475     //====> Corrected Jets
476     Met MVAMet::GetMet( Bool_t iPhi,Float_t iPtVis,Float_t iPhiVis,Float_t iSumEtVis,
477     const PFMet *iMet ,
478 pharris 1.3 const PFCandidateCol *iCands,
479 pharris 1.7 const Vertex *iVertex ,const VertexCol *iVertices,Double_t iRho,
480 pharris 1.2 const PFJetCol *iJets ,
481     int iNPV,
482     Bool_t printDebug) {
483 pharris 1.14
484     int lNPV = 0;
485     for(unsigned int i0 = 0; i0 < iVertices->GetEntries(); i0++) {
486     const Vertex *lPV = iVertices->At(i0);
487     if(lPV->Ndof() < 4.0) continue;
488     if(fabs(lPV->Z()) > 24.) continue;
489     if(lPV->Position().Rho() > 2.) continue;
490     lNPV++;
491     }
492    
493 pharris 1.7 Met lPFRec = fRecoilTools->pfRecoil (iPtVis,iPhiVis,iSumEtVis,iCands);
494 pharris 1.3 Met lTKRec = fRecoilTools->trackRecoil(iPtVis,iPhiVis,iSumEtVis, iCands,iVertex);
495 pharris 1.7 Met lNPRec = fRecoilTools->NoPURecoil (iPtVis,iPhiVis,iSumEtVis,iJets,iCands,iVertex,iVertices,iRho);
496     Met lPCRec = fRecoilTools->PUCRecoil (iPtVis,iPhiVis,iSumEtVis,iJets,iCands,iVertex,iVertices,iRho);
497     Met lPUMet = fRecoilTools->PUMet ( iJets,iCands,iVertex,iVertices,iRho);
498 pharris 1.2
499     Double_t lPt0 = 0; const PFJet *lLead = 0;
500     Double_t lPt1 = 0; const PFJet *l2nd = 0;
501     int lNAllJet = 0;
502     int lNJet = 0;
503     for(unsigned int i0 = 0; i0 < iJets->GetEntries(); i0++) {
504     const PFJet *pJet = iJets->At(i0);
505     Double_t pPt = pJet->Pt();
506     if(!JetTools::passPFLooseId(pJet)) continue;
507     lNAllJet++;
508     if(pPt > 30) lNJet++;
509     if(lPt0 < pPt) {lPt0 = pPt; lLead = pJet; continue;}
510     if(lPt1 < pPt) {lPt1 = pPt; l2nd = pJet; continue;}
511     }
512    
513 pharris 1.5 fSumEt = lPFRec.SumEt();
514 pharris 1.2 fU = lPFRec.Pt();
515     fUPhi = lPFRec.Phi();
516     fTKSumEt = lTKRec.SumEt()/lPFRec.SumEt();
517     fTKU = lTKRec.Pt();
518     fTKUPhi = lTKRec.Phi();
519     fNPSumEt = lNPRec.SumEt()/lPFRec.SumEt();
520     fNPU = lNPRec.Pt();
521     fNPUPhi = lNPRec.Phi();
522     fPUSumEt = lPUMet.SumEt()/lPFRec.SumEt();
523     fPUMet = lPUMet.Pt();
524     fPUMetPhi = lPUMet.Phi();
525     fPCSumEt = lPCRec.SumEt()/lPFRec.SumEt();
526     fPCU = lPCRec.Pt() ;
527     fPCUPhi = lPCRec.Phi() ;
528     fJSPt1 = lPt0;
529     fJSEta1 = 0; if(lLead != 0) fJSEta1 = lLead->Eta();
530     fJSPhi1 = 0; if(lLead != 0) fJSPhi1 = lLead->Phi();
531     fJSPt2 = lPt1;
532     fJSEta2 = 0; if(l2nd != 0) fJSEta2 = l2nd ->Eta();
533     fJSPhi2 = 0; if(l2nd != 0) fJSPhi2 = l2nd ->Phi();
534     fNJet = lNJet ;
535     fNAllJet = lNAllJet;
536 pharris 1.14 fNPV = lNPV ;
537 pharris 1.7
538 pharris 1.2 Float_t lMVA = evaluatePhi();
539    
540     if(!iPhi) fUPhiMVA = fUPhi + lMVA;
541     //Not no nice feature of teh training
542 pharris 1.4 //fTKSumEt /= lPFRec.SumEt();
543     //fNPSumEt /= lPFRec.SumEt();
544     //fPUSumEt /= lPFRec.SumEt();
545     //fPCSumEt /= lPFRec.SumEt();
546 pharris 1.7 if(!iPhi) lMVA = evaluateU1();//fU1Reader ->EvaluateMVA( fU1MethodName );
547     fUMVA = lMVA*fU;
548    
549     TLorentzVector lUVec(0,0,0,0); lUVec .SetPtEtaPhiM(fUMVA ,0,fUPhiMVA,0);
550 pharris 1.2 TLorentzVector lVVec(0,0,0,0); lVVec .SetPtEtaPhiM(iPtVis ,0,iPhiVis ,0);
551     if(lMVA < 0) lUVec .RotateZ(TMath::Pi());
552     lUVec -= lVVec;
553     Met lMet(lUVec.Px(),lUVec.Py());
554 pharris 1.4 //Cov matrix
555     double lCovU1 = evaluateCovU1();
556     double lCovU2 = evaluateCovU2();
557    
558     double lCos2 = cos(fUPhiMVA)*cos(fUPhiMVA);
559     double lSin2 = sin(fUPhiMVA)*sin(fUPhiMVA);
560    
561     //Now Compute teh covariance matrix in X and Y
562 pharris 1.5 (*fCov)(0,0) = lCovU1*lCos2+lCovU2*lSin2;
563     (*fCov)(1,0) = -lCovU1*sin(fUPhiMVA)*cos(fUPhiMVA)+lCovU2*sin(fUPhiMVA)*cos(fUPhiMVA);
564     (*fCov)(0,1) = (*fCov)(1,0);
565     (*fCov)(1,1) = lCovU1*lSin2+lCovU2*lCos2;
566 pharris 1.10 TMatrixD lInv(2,2);
567     lInv(0,0) = (*fCov)(0,0); lInv(1,1) = (*fCov)(1,1); lInv(1,0) = (*fCov)(1,0); lInv(0,1) = (*fCov)(0,1);
568     if(lInv.Determinant() != 0) lInv.Invert();
569     fSignificance = TMath::Sqrt(lUVec.Px()*lUVec.Px()*(lInv)(0,0) + 2.*lUVec.Px()*lUVec.Py()*(lInv)(1,0) + lUVec.Py()*lUVec.Py()*(lInv)(1,1));
570     fUncertainty = sqrt(lCovU1+lCovU2);
571 pharris 1.2
572     if (printDebug == kTRUE) {
573 ceballos 1.11 std::cout << "Debug Met MVA: "
574 pharris 1.2 << fU << " : "
575     << fUPhi << " : "
576     << fTKSumEt << " : "
577     << fTKU << " : "
578     << fTKUPhi << " : "
579     << fNPSumEt << " : "
580     << fNPU << " : "
581     << fNPUPhi << " : "
582     << fPUSumEt << " : "
583     << fPUMet << " : "
584     << fPUMetPhi << " : "
585     << fPCUPhi << " : "
586     << fPCSumEt << " : "
587     << fPCU << " : "
588     << fPCUPhi << " : "
589     << fJSPt1 << " : "
590     << fJSEta1 << " : "
591     << fJSPhi1 << " : "
592     << fJSPt2 << " : "
593     << fJSEta2 << " : "
594     << fJSPhi2 << " : "
595     << fNJet << " : "
596     << fNAllJet << " : "
597     << fNPV << " : "
598     << " === : === "
599 ceballos 1.12 << lMet.Pt() << " : "
600     << lMet.Phi() << " : "
601 pharris 1.2 << std::endl;
602     }
603    
604     return lMet;
605     }
606     //--------------------------------------------------------------------------------------------------
607 pharris 1.1 //Interms of the two candidates
608     Met MVAMet::GetMet( Bool_t iPhi,
609     Float_t iPt1,Float_t iPhi1,Float_t iEta1,
610     Float_t iPt2,Float_t iPhi2,Float_t iEta2,
611     const PFMet *iMet ,
612 pharris 1.3 const PFCandidateCol *iCands,
613     const Vertex *iVertex,const VertexCol *iVertices,
614 pharris 1.1 const PFJetCol *iJets ,
615     FactorizedJetCorrector *iJetCorrector,
616     const PileupEnergyDensityCol *iPUEnergyDensity,
617     int iNPV,
618     Bool_t printDebug) {
619 pharris 1.14
620     int lNPV = 0;
621     for(unsigned int i0 = 0; i0 < iVertices->GetEntries(); i0++) {
622     const Vertex *lPV = iVertices->At(i0);
623     if(lPV->Ndof() < 4.0) continue;
624     if(fabs(lPV->Z()) > 24.) continue;
625     if(lPV->Position().Rho() > 2.) continue;
626     lNPV++;
627     }
628    
629 pharris 1.1 TLorentzVector lVVec1(0,0,0,0); lVVec1.SetPtEtaPhiM(iPt1,0,iPhi1 ,0);
630     TLorentzVector lVVec2(0,0,0,0); lVVec2.SetPtEtaPhiM(iPt2,0,iPhi2 ,0);
631     lVVec1+=lVVec2;
632     Float_t lPtVis = lVVec1.Pt();
633     Float_t lPhiVis = lVVec1.Phi();
634     Float_t lSumEtVis = iPt1 + iPt2;
635 pharris 1.7 Met lPFRec = fRecoilTools->pfRecoil (lPtVis,lPhiVis,lSumEtVis,iCands);
636 pharris 1.1 Met lTKRec = fRecoilTools->trackRecoil(lPtVis,lPhiVis,lSumEtVis,iCands,iVertex);
637 pharris 1.3 Met lNPRec = fRecoilTools->NoPURecoil (lPtVis,lPhiVis,lSumEtVis,iJets,iJetCorrector,iPUEnergyDensity,iCands,iVertex,iVertices,
638     iPhi1,iEta1,iPhi2,iEta2);
639     Met lPCRec = fRecoilTools->PUCRecoil (lPtVis,lPhiVis,lSumEtVis,iJets,iJetCorrector,iPUEnergyDensity,iCands,iVertex,iVertices,
640     iPhi1,iEta1,iPhi2,iEta2);
641     Met lPUMet = fRecoilTools->PUMet ( iJets,iJetCorrector,iPUEnergyDensity,iCands,iVertex,iVertices,
642     iPhi1,iEta1,iPhi2,iEta2);
643 pharris 1.1
644     Double_t lPt0 = 0; const PFJet *lLead = 0;
645     Double_t lPt1 = 0; const PFJet *l2nd = 0;
646     int lNAllJet = 0;
647     int lNJet = 0;
648     for(unsigned int i0 = 0; i0 < iJets->GetEntries(); i0++) {
649     const PFJet *pJet = iJets->At(i0);
650     Double_t pPt = fRecoilTools->fJetIDMVA->correctedPt(pJet,iJetCorrector,iPUEnergyDensity);
651     double pDEta1 = pJet->Eta() - iEta1;
652     double pDPhi1 = fabs(pJet->Phi() - iPhi1); if(pDPhi1 > 2.*TMath::Pi()-pDPhi1) pDPhi1 = 2.*TMath::Pi()-pDPhi1;
653     double pDR1 = sqrt(pDEta1*pDEta1 + pDPhi1*pDPhi1);
654     if(pDR1 < 0.5) continue;
655     double pDEta2 = pJet->Eta() - iEta2;
656     double pDPhi2 = fabs(pJet->Phi() - iPhi2); if(pDPhi2 > 2.*TMath::Pi()-pDPhi2) pDPhi2 = 2.*TMath::Pi()-pDPhi2;
657     double pDR2 = sqrt(pDEta2*pDEta2 + pDPhi2*pDPhi2);
658     if(pDR2 < 0.5) continue;
659     if(!JetTools::passPFLooseId(pJet)) continue;
660     lNAllJet++;
661     if(pPt > 30.) lNJet++;
662     if(lPt0 < pPt) {lPt0 = pPt; lLead = pJet; continue;}
663     if(lPt1 < pPt) {lPt1 = pPt; l2nd = pJet; continue;}
664     }
665    
666 pharris 1.5 fSumEt = lPFRec.SumEt();
667 pharris 1.1 fU = lPFRec.Pt();
668     fUPhi = lPFRec.Phi();
669     fTKSumEt = lTKRec.SumEt()/lPFRec.SumEt();
670     fTKU = lTKRec.Pt();
671     fTKUPhi = lTKRec.Phi();
672     fNPSumEt = lNPRec.SumEt()/lPFRec.SumEt();
673     fNPU = lNPRec.Pt();
674     fNPUPhi = lNPRec.Phi();
675     fPUSumEt = lPUMet.SumEt()/lPFRec.SumEt();
676     fPUMet = lPUMet.Pt();
677     fPUMetPhi= lPUMet.Phi();
678     fPCSumEt = lPCRec.SumEt()/lPFRec.SumEt();
679     fPCU = lPCRec.Pt() ;
680     fPCUPhi = lPCRec.Phi() ;
681     fJSPt1 = lPt0;
682     fJSEta1 = 0; if(lLead != 0) fJSEta1 = lLead->Eta();
683     fJSPhi1 = 0; if(lLead != 0) fJSPhi1 = lLead->Phi();
684     fJSPt2 = lPt1;
685     fJSEta2 = 0; if(l2nd != 0) fJSEta2 = l2nd ->Eta();
686     fJSPhi2 = 0; if(l2nd != 0) fJSPhi2 = l2nd ->Phi();
687     fNJet = lNJet ;
688     fNAllJet = lNAllJet;
689 pharris 1.14 fNPV = lNPV ;
690 pharris 1.5
691 pharris 1.2 Float_t lMVA = evaluatePhi();//fPhiReader ->EvaluateMVA( fPhiMethodName );
692     if(!iPhi) fUPhiMVA = fUPhi + lMVA;
693 pharris 1.4 //fTKSumEt /= lPFRec.SumEt();
694     //fNPSumEt /= lPFRec.SumEt();
695     //fPUSumEt /= lPFRec.SumEt();
696     //fPCSumEt /= lPFRec.SumEt();
697 pharris 1.7 if(!iPhi) lMVA = evaluateU1();//fU1Reader ->EvaluateMVA( fU1MethodName );
698     fUMVA = fU*lMVA;
699    
700     TLorentzVector lUVec (0,0,0,0); lUVec .SetPtEtaPhiM(fUMVA ,0,fUPhiMVA,0);
701 pharris 1.2 TLorentzVector lVVec (0,0,0,0); lVVec .SetPtEtaPhiM(lPtVis ,0,lPhiVis ,0);
702     if(lMVA < 0) lUVec .RotateZ(TMath::Pi());
703     lUVec -= lVVec;
704     Met lMet(lUVec.Px(),lUVec.Py());
705 pharris 1.4 //Cov matrix
706     double lCovU1 = evaluateCovU1();
707     double lCovU2 = evaluateCovU2();
708    
709     double lCos2 = cos(fUPhiMVA)*cos(fUPhiMVA);
710     double lSin2 = sin(fUPhiMVA)*sin(fUPhiMVA);
711    
712     //Now Compute teh covariance matrix in X and Y
713 pharris 1.5 (*fCov)(0,0) = lCovU1*lCos2+lCovU2*lSin2;
714     (*fCov)(1,0) = -lCovU1*sin(fUPhiMVA)*cos(fUPhiMVA)+lCovU2*sin(fUPhiMVA)*cos(fUPhiMVA);
715     (*fCov)(0,1) = (*fCov)(1,0);
716     (*fCov)(1,1) = lCovU1*lSin2+lCovU2*lCos2;
717 pharris 1.10 TMatrixD lInv(2,2);
718     lInv(0,0) = (*fCov)(0,0); lInv(1,1) = (*fCov)(1,1); lInv(1,0) = (*fCov)(1,0); lInv(0,1) = (*fCov)(0,1);
719     if(lInv.Determinant() != 0) lInv.Invert();
720     fSignificance = TMath::Sqrt(lUVec.Px()*lUVec.Px()*(lInv)(0,0) + 2.*lUVec.Px()*lUVec.Py()*(lInv)(1,0) + lUVec.Py()*lUVec.Py()*(lInv)(1,1));
721     fUncertainty = sqrt(lCovU1+lCovU2);
722    
723 pharris 1.2
724     if (printDebug == kTRUE) {
725 ceballos 1.11 std::cout << "Debug Met MVA: "
726 pharris 1.14 << fSumEt << " : "
727 pharris 1.2 << fU << " : "
728     << fUPhi << " : "
729     << fTKSumEt << " : "
730     << fTKU << " : "
731     << fTKUPhi << " : "
732     << fNPSumEt << " : "
733     << fNPU << " : "
734     << fNPUPhi << " : "
735     << fPUSumEt << " : "
736     << fPUMet << " : "
737     << fPUMetPhi << " : "
738     << fPCSumEt << " : "
739     << fPCU << " : "
740     << fPCUPhi << " : "
741     << fJSPt1 << " : "
742     << fJSEta1 << " : "
743     << fJSPhi1 << " : "
744     << fJSPt2 << " : "
745     << fJSEta2 << " : "
746     << fJSPhi2 << " : "
747     << fNJet << " : "
748     << fNAllJet << " : "
749     << fNPV << " : "
750     << " === : === "
751 ceballos 1.12 << lMet.Pt() << " : "
752     << lMet.Phi() << " : "
753 pharris 1.2 << std::endl;
754     }
755    
756     return lMet;
757     }
758     //--------------------------------------------------------------------------------------------------
759     //Interms of the two candidates => corrected Jets
760     Met MVAMet::GetMet( Bool_t iPhi,
761     Float_t iPt1,Float_t iPhi1,Float_t iEta1,
762     Float_t iPt2,Float_t iPhi2,Float_t iEta2,
763     const PFMet *iMet ,
764 pharris 1.3 const PFCandidateCol *iCands,
765 pharris 1.7 const Vertex *iVertex,const VertexCol *iVertices,Double_t iRho,
766 pharris 1.2 const PFJetCol *iJets ,
767     int iNPV,
768     Bool_t printDebug) {
769 pharris 1.14 int lNPV = 0;
770     for(unsigned int i0 = 0; i0 < iVertices->GetEntries(); i0++) {
771     const Vertex *lPV = iVertices->At(i0);
772     if(lPV->Ndof() < 4.0) continue;
773     if(fabs(lPV->Z()) > 24.) continue;
774     if(lPV->Position().Rho() > 2.) continue;
775     lNPV++;
776     }
777 pharris 1.2
778     TLorentzVector lVVec1(0,0,0,0); lVVec1.SetPtEtaPhiM(iPt1,0,iPhi1 ,0);
779     TLorentzVector lVVec2(0,0,0,0); lVVec2.SetPtEtaPhiM(iPt2,0,iPhi2 ,0);
780     lVVec1+=lVVec2;
781     Float_t lPtVis = lVVec1.Pt();
782     Float_t lPhiVis = lVVec1.Phi();
783     Float_t lSumEtVis = iPt1 + iPt2;
784 pharris 1.7 Met lPFRec = fRecoilTools->pfRecoil (lPtVis,lPhiVis,lSumEtVis,iCands);
785 pharris 1.2 Met lTKRec = fRecoilTools->trackRecoil(lPtVis,lPhiVis,lSumEtVis,iCands,iVertex);
786 pharris 1.3 Met lNPRec = fRecoilTools->NoPURecoil (lPtVis,lPhiVis,lSumEtVis,iJets,iCands,iVertex,iVertices,iPhi1,iEta1,iPhi2,iEta2);
787     Met lPCRec = fRecoilTools->PUCRecoil (lPtVis,lPhiVis,lSumEtVis,iJets,iCands,iVertex,iVertices,iPhi1,iEta1,iPhi2,iEta2);
788     Met lPUMet = fRecoilTools->PUMet ( iJets,iCands,iVertex,iVertices,iPhi1,iEta1,iPhi2,iEta2);
789 pharris 1.2
790     Double_t lPt0 = 0; const PFJet *lLead = 0;
791     Double_t lPt1 = 0; const PFJet *l2nd = 0;
792     int lNAllJet = 0;
793     int lNJet = 0;
794     for(unsigned int i0 = 0; i0 < iJets->GetEntries(); i0++) {
795     const PFJet *pJet = iJets->At(i0);
796     Double_t pPt = pJet->Pt();
797     double pDEta1 = pJet->Eta() - iEta1;
798     double pDPhi1 = fabs(pJet->Phi() - iPhi1); if(pDPhi1 > 2.*TMath::Pi()-pDPhi1) pDPhi1 = 2.*TMath::Pi()-pDPhi1;
799     double pDR1 = sqrt(pDEta1*pDEta1 + pDPhi1*pDPhi1);
800     if(pDR1 < 0.5) continue;
801     double pDEta2 = pJet->Eta() - iEta2;
802     double pDPhi2 = fabs(pJet->Phi() - iPhi2); if(pDPhi2 > 2.*TMath::Pi()-pDPhi2) pDPhi2 = 2.*TMath::Pi()-pDPhi2;
803     double pDR2 = sqrt(pDEta2*pDEta2 + pDPhi2*pDPhi2);
804     if(pDR2 < 0.5) continue;
805     if(!JetTools::passPFLooseId(pJet)) continue;
806     lNAllJet++;
807     if(pPt > 30.) lNJet++;
808     if(lPt0 < pPt) {lPt0 = pPt; lLead = pJet; continue;}
809     if(lPt1 < pPt) {lPt1 = pPt; l2nd = pJet; continue;}
810     }
811 pharris 1.5 fSumEt = lPFRec.SumEt();
812 pharris 1.2 fU = lPFRec.Pt();
813     fUPhi = lPFRec.Phi();
814     fTKSumEt = lTKRec.SumEt()/lPFRec.SumEt();
815     fTKU = lTKRec.Pt();
816     fTKUPhi = lTKRec.Phi();
817     fNPSumEt = lNPRec.SumEt()/lPFRec.SumEt();
818     fNPU = lNPRec.Pt();
819     fNPUPhi = lNPRec.Phi();
820     fPUSumEt = lPUMet.SumEt()/lPFRec.SumEt();
821     fPUMet = lPUMet.Pt();
822     fPUMetPhi= lPUMet.Phi();
823     fPCSumEt = lPCRec.SumEt()/lPFRec.SumEt();
824     fPCU = lPCRec.Pt() ;
825     fPCUPhi = lPCRec.Phi() ;
826     fJSPt1 = lPt0;
827     fJSEta1 = 0; if(lLead != 0) fJSEta1 = lLead->Eta();
828     fJSPhi1 = 0; if(lLead != 0) fJSPhi1 = lLead->Phi();
829     fJSPt2 = lPt1;
830     fJSEta2 = 0; if(l2nd != 0) fJSEta2 = l2nd ->Eta();
831     fJSPhi2 = 0; if(l2nd != 0) fJSPhi2 = l2nd ->Phi();
832     fNJet = lNJet ;
833     fNAllJet = lNAllJet;
834 pharris 1.14 fNPV = lNPV ;
835 pharris 1.2
836     Float_t lMVA = evaluatePhi();//fPhiReader ->EvaluateMVA( fPhiMethodName );
837 pharris 1.1
838     if(!iPhi) fUPhiMVA = fUPhi + lMVA;
839 pharris 1.4 //fTKSumEt /= lPFRec.SumEt();
840     //fNPSumEt /= lPFRec.SumEt();
841     //fPUSumEt /= lPFRec.SumEt();
842     //fPCSumEt /= lPFRec.SumEt();
843 pharris 1.7 if(!iPhi) lMVA = evaluateU1();//fU1Reader ->EvaluateMVA( fU1MethodName );
844     fUMVA = lMVA*fU;
845    
846     TLorentzVector lUVec (0,0,0,0); lUVec .SetPtEtaPhiM(fUMVA ,0,fUPhiMVA,0);
847 pharris 1.1 TLorentzVector lVVec (0,0,0,0); lVVec .SetPtEtaPhiM(lPtVis ,0,lPhiVis ,0);
848 pharris 1.7 if(lMVA < 0) lUVec .RotateZ(TMath::Pi());
849 pharris 1.1 lUVec -= lVVec;
850     Met lMet(lUVec.Px(),lUVec.Py());
851 pharris 1.4 //Cov matrix
852     double lCovU1 = evaluateCovU1();
853     double lCovU2 = evaluateCovU2();
854    
855     double lCos2 = cos(fUPhiMVA)*cos(fUPhiMVA);
856     double lSin2 = sin(fUPhiMVA)*sin(fUPhiMVA);
857     //Now Compute teh covariance matrix in X and Y
858 pharris 1.5 (*fCov)(0,0) = lCovU1*lCos2+lCovU2*lSin2;
859     (*fCov)(1,0) = -lCovU1*sin(fUPhiMVA)*cos(fUPhiMVA)+lCovU2*sin(fUPhiMVA)*cos(fUPhiMVA);
860     (*fCov)(0,1) = (*fCov)(1,0);
861     (*fCov)(1,1) = lCovU1*lSin2+lCovU2*lCos2;
862 pharris 1.10 TMatrixD lInv(2,2);
863     lInv(0,0) = (*fCov)(0,0); lInv(1,1) = (*fCov)(1,1); lInv(1,0) = (*fCov)(1,0); lInv(0,1) = (*fCov)(0,1);
864     if(lInv.Determinant() != 0) lInv.Invert();
865     fSignificance = TMath::Sqrt(lUVec.Px()*lUVec.Px()*(lInv)(0,0) + 2.*lUVec.Px()*lUVec.Py()*(lInv)(1,0) + lUVec.Py()*lUVec.Py()*(lInv)(1,1));
866     fUncertainty = sqrt(lCovU1+lCovU2);
867 pharris 1.7
868 pharris 1.1 if (printDebug == kTRUE) {
869 ceballos 1.11 std::cout << "Debug Met MVA: "
870 pharris 1.1 << fU << " : "
871     << fUPhi << " : "
872     << fTKSumEt << " : "
873     << fTKU << " : "
874     << fTKUPhi << " : "
875     << fNPSumEt << " : "
876     << fNPU << " : "
877     << fNPUPhi << " : "
878     << fPUSumEt << " : "
879     << fPUMet << " : "
880     << fPUMetPhi << " : "
881     << fPCSumEt << " : "
882     << fPCU << " : "
883     << fPCUPhi << " : "
884     << fJSPt1 << " : "
885     << fJSEta1 << " : "
886     << fJSPhi1 << " : "
887     << fJSPt2 << " : "
888     << fJSEta2 << " : "
889     << fJSPhi2 << " : "
890     << fNJet << " : "
891     << fNAllJet << " : "
892     << fNPV << " : "
893     << " === : === "
894 ceballos 1.12 << lMet.Pt() << " : "
895     << lMet.Phi() << " : "
896 pharris 1.1 << std::endl;
897     }
898    
899     return lMet;
900     }
901