ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Utils/src/JetIDMVA.cc
(Generate patch)

Comparing UserCode/MitPhysics/Utils/src/JetIDMVA.cc (file contents):
Revision 1.7 by pharris, Tue May 1 16:53:38 2012 UTC vs.
Revision 1.11 by pharris, Mon May 14 09:42:27 2012 UTC

# Line 15 | Line 15 | using namespace mithep;
15  
16   //--------------------------------------------------------------------------------------------------
17   JetIDMVA::JetIDMVA() :
18 <  fJetPtMin(0)  , //We need to lower this
18 >  fJetPtMin(0.)  , //We need to lower this
19    fDZCut   (0.2),
20    fLowPtMethodName ("JetIDMVALowPt" ),
21    fHighPtMethodName("JetIDMVAHighPt"),
# Line 35 | Line 35 | JetIDMVA::JetIDMVA() :
35    fFrac02   (0),
36    fFrac03   (0),
37    fFrac04   (0),
38 <  fFrac05   (0)
38 >  fFrac05   (0),
39 >  fDR2Mean  (0)
40   {    
41    fReader = 0;
42   }
43   //--------------------------------------------------------------------------------------------------
44   JetIDMVA::~JetIDMVA() {
45  
46 <  fReader = 0;
46 >  delete fReader;
47 >  delete fLowPtReader;
48   }
49  
50   //--------------------------------------------------------------------------------------------------
# Line 56 | Line 58 | void JetIDMVA::Initialize( JetIDMVA::Cut
58    fIsInitialized = kTRUE;
59    fType          = iType;
60    fCutType       = iCutType;
61 +  fLowPtReader   = 0;
62 +  fLowPtReader   = new TMVA::Reader( "!Color:!Silent:Error" );  
63 +  fLowPtReader->AddVariable( "nvtx"     , &fNVtx      );
64 +  fLowPtReader->AddVariable( "jetPt"    , &fJPt1      );  
65 +  fLowPtReader->AddVariable( "jetEta"   , &fJEta1     );
66 +  fLowPtReader->AddVariable( "jetPhi"   , &fJPhi1     );            
67 +  fLowPtReader->AddVariable( "dZ"       , &fJDZ1      );
68 +  fLowPtReader->AddVariable( "d0"       , &fJD01      );
69 +  fLowPtReader->AddVariable( "beta"     , &fBeta      );
70 +  fLowPtReader->AddVariable( "betaStar" , &fBetaStar  );
71 +  fLowPtReader->AddVariable( "nCharged" , &fNCharged  );
72 +  fLowPtReader->AddVariable( "nNeutrals", &fNNeutrals );
73 +  fLowPtReader->AddVariable( "dRMean"   , &fDRMean    );
74 +  fLowPtReader->AddVariable( "frac01"   , &fFrac01    );
75 +  fLowPtReader->AddVariable( "frac02"   , &fFrac02    );
76 +  fLowPtReader->AddVariable( "frac03"   , &fFrac03    );
77 +  fLowPtReader->AddVariable( "frac04"   , &fFrac04    );
78 +  fLowPtReader->AddVariable( "frac05"   , &fFrac05    );
79 +  
80    fReader        = 0;
81    fReader        = new TMVA::Reader( "!Color:!Silent:Error" );  
82    if (fType == kBaseline) {
# Line 76 | Line 97 | void JetIDMVA::Initialize( JetIDMVA::Cut
97      fReader->AddVariable( "frac04"   , &fFrac04    );
98      fReader->AddVariable( "frac05"   , &fFrac05    );
99    }
100 <  fReader->BookMVA(fLowPtMethodName  , iLowPtWeights );
100 >  if (fType == k42) {
101 >    fReader->AddVariable( "frac01"   , &fFrac01    );
102 >    fReader->AddVariable( "frac02"   , &fFrac02    );
103 >    fReader->AddVariable( "frac03"   , &fFrac03    );
104 >    fReader->AddVariable( "frac04"   , &fFrac04    );
105 >    fReader->AddVariable( "frac05"   , &fFrac05    );
106 >    fReader->AddVariable( "nvtx"     , &fNVtx      );
107 >    fReader->AddVariable( "nNeutrals", &fNNeutrals );
108 >    fReader->AddVariable( "beta"     , &fBeta      );
109 >    fReader->AddVariable( "betaStar" , &fBetaStar  );
110 >    fReader->AddVariable( "dZ"       , &fJDZ1      );
111 >    fReader->AddVariable( "nCharged" , &fNCharged  );
112 >    fReader->AddSpectator( "jetPt"    , &fJPt1      );  
113 >    fReader->AddSpectator( "jetEta"   , &fJEta1     );
114 >  }
115 >  if (fType == k52) {
116 >    fReader->AddVariable( "frac01"   , &fFrac01    );
117 >    fReader->AddVariable( "frac02"   , &fFrac02    );
118 >    fReader->AddVariable( "frac03"   , &fFrac03    );
119 >    fReader->AddVariable( "frac04"   , &fFrac04    );
120 >    fReader->AddVariable( "frac05"   , &fFrac05    );
121 >    fReader->AddVariable( "dR2Mean"  , &fDR2Mean   );
122 >    fReader->AddVariable( "nvtx"     , &fNVtx      );
123 >    fReader->AddVariable( "nNeutrals", &fNNeutrals );
124 >    fReader->AddVariable( "beta"     , &fBeta      );
125 >    fReader->AddVariable( "betaStar" , &fBetaStar  );
126 >    fReader->AddVariable( "dZ"       , &fJDZ1      );
127 >    fReader->AddVariable( "nCharged" , &fNCharged  );
128 >    fReader->AddSpectator( "jetPt"    , &fJPt1      );  
129 >    fReader->AddSpectator( "jetEta"   , &fJEta1     );
130 >  }
131 >
132 >  fLowPtReader->BookMVA(fLowPtMethodName  , iLowPtWeights );
133    fReader->BookMVA(fHighPtMethodName , iHighPtWeights );
134    std::cout << "Jet ID MVA Initialization\n";
135    std::cout << "MethodName : " << fLowPtMethodName << " , type == " << fType << std::endl;
136  
137 +  std::string lCutId = "JetIdParams";
138 +  if(fType == k42) lCutId = "PuJetIdOptMVA_wp";
139 +  if(fType == k52) lCutId = "full_5x_wp";
140    //Load Cut Matrix
141 <  edm::ParameterSet lConfig = edm::readPSetsFrom(iCutFileName.Data())->getParameter<edm::ParameterSet>("JetIdParams");
141 >  edm::ParameterSet lDefConfig = edm::readPSetsFrom(iCutFileName.Data())->getParameter<edm::ParameterSet>("JetIdParams");
142 >  edm::ParameterSet lConfig    = edm::readPSetsFrom(iCutFileName.Data())->getParameter<edm::ParameterSet>(lCutId);
143    std::string lCutType = "Tight";
144    if(fCutType == kMedium) lCutType = "Medium";
145    if(fCutType == kLoose ) lCutType = "Loose";
146    if(fCutType == kMET   ) lCutType = "MET";
147 <  std::vector<double> lPt010  = lConfig.getParameter<std::vector<double> >(("Pt010_" +lCutType).c_str());
147 >  std::string lLowPtCut = "MET";
148 >  std::vector<double> lPt010  = lDefConfig.getParameter<std::vector<double> >(("Pt010_" +lLowPtCut).c_str());
149    std::vector<double> lPt1020 = lConfig.getParameter<std::vector<double> >(("Pt1020_"+lCutType).c_str());
150    std::vector<double> lPt2030 = lConfig.getParameter<std::vector<double> >(("Pt2030_"+lCutType).c_str());
151    std::vector<double> lPt3050 = lConfig.getParameter<std::vector<double> >(("Pt3050_"+lCutType).c_str());
# Line 117 | Line 175 | Double_t JetIDMVA::MVAValue(
175                              Float_t iFrac02  ,
176                              Float_t iFrac03  ,
177                              Float_t iFrac04  ,
178 <                            Float_t iFrac05  
178 >                            Float_t iFrac05  ,
179 >                            Float_t iDR2Mean  
180                              ){
181    
182    if(!fIsInitialized) {
# Line 141 | Line 200 | Double_t JetIDMVA::MVAValue(
200    fFrac03    = iFrac03;
201    fFrac04    = iFrac04;
202    fFrac05    = iFrac05;
203 +  fDR2Mean   = iDR2Mean;
204  
205    Double_t lMVA = -9999;  
206 <  if(iJPt1 < 10) lMVA = fReader->EvaluateMVA( fLowPtMethodName  );
206 >  if(iJPt1 < 10) lMVA = fLowPtReader->EvaluateMVA( fLowPtMethodName  );
207    if(iJPt1 > 10) lMVA = fReader->EvaluateMVA( fHighPtMethodName );
208  
209    return lMVA;
# Line 151 | Line 211 | Double_t JetIDMVA::MVAValue(
211   //--------------------------------------------------------------------------------------------------
212   Bool_t JetIDMVA::pass(const PFJet *iJet,const Vertex *iVertex,const VertexCol *iVertices,
213                        FactorizedJetCorrector *iJetCorrector,
214 <                      const PileupEnergyDensityCol *iPileupEnergyDensity) {
214 >                      const PileupEnergyDensityCol *iPileupEnergyDensity,
215 >                      RhoUtilities::RhoType type) {
216    
217    if(!JetTools::passPFLooseId(iJet))                 return false;
157 –  if(iJet->Pt() < fJetPtMin)                         return false;
218    if(fabs(iJet->Eta()) > 4.99)                       return false;
219    
220    double lMVA = MVAValue   (iJet,iVertex,iVertices,iJetCorrector,iPileupEnergyDensity);
221 <  double lPt  = correctedPt(iJet,                  iJetCorrector,iPileupEnergyDensity);
222 <  if(lPt > 50)                                return true; //==> we can raise this
221 >  double lPt  = correctedPt(iJet,                  iJetCorrector,iPileupEnergyDensity,type);
222 >  if(lPt < fJetPtMin)                         return false;
223    
224    int lPtId = 0;
225    if(lPt > 10 && lPt < 20) lPtId = 1;
# Line 174 | Line 234 | Bool_t JetIDMVA::pass(const PFJet *iJet,
234    double lMVACut = fMVACut[lPtId][lEtaId];
235    if(lMVA < lMVACut) return false;
236    return true;
177 –   //if( fabs(JetTools::impactParameter(iJet,iVertex,true)) < 0.2) return true;
178 –  //if(correctedPt(iJet,iJetCorrector,iPileupEnergyDensity) < fJetPtMin && iJet->TrackCountingHighEffBJetTagsDisc() == -100) return false;
179 –  //This line is a bug in the Met training
180 –  //if(lMVA < -0.8)                            return false;
181 –  //if(lMVA < -0.5 && fabs(iJet->Eta()) > 3.0) return false;
237   }
238   //--------------------------------------------------------------------------------------------------
239   Bool_t JetIDMVA::pass(const PFJet *iJet,const Vertex *iVertex,const VertexCol *iVertices) {
240    if(!JetTools::passPFLooseId(iJet))                 return false;
241    if(iJet->Pt()        < fJetPtMin) return false;
187 –  if(iJet->Pt()        > 50)        return true; //==> we can raise this
242    if(fabs(iJet->Eta()) > 4.99)      return false;
243    double lMVA = MVAValue(iJet,iVertex,iVertices);
244    
# Line 231 | Line 285 | Double_t JetIDMVA::MVAValue(const PFJet
285    fNNeutrals = iJet->NeutralMultiplicity();
286  
287    fDRMean    = JetTools::dRMean(iJet,-1);
288 +  fDR2Mean   = JetTools::dR2Mean(iJet,-1);
289    fFrac01    = JetTools::frac  (iJet,0.1,0. ,-1);
290    fFrac02    = JetTools::frac  (iJet,0.2,0.1,-1);
291    fFrac03    = JetTools::frac  (iJet,0.3,0.2,-1);
# Line 238 | Line 293 | Double_t JetIDMVA::MVAValue(const PFJet
293    fFrac05    = JetTools::frac  (iJet,0.5,0.4,-1);
294  
295    double lMVA = 0;
296 <  if(fJPt1 < 10) lMVA = fReader->EvaluateMVA( fLowPtMethodName  );
296 >  if(fJPt1 < 10) lMVA = fLowPtReader->EvaluateMVA( fLowPtMethodName  );
297    if(fJPt1 > 10) lMVA = fReader->EvaluateMVA( fHighPtMethodName );
298    if (printDebug == kTRUE) {
299      std::cout << "Debug Jet MVA: "
# Line 257 | Line 312 | Double_t JetIDMVA::MVAValue(const PFJet
312                << fFrac02    << " "
313                << fFrac03    << " "
314                << fFrac04    << " "
315 <              << fFrac05    
315 >              << fFrac05    << " "
316 >              << fDRMean    
317                << " === : === "
318                << lMVA << " "    
319                << std::endl;
# Line 286 | Line 342 | Double_t JetIDMVA::MVAValue(const PFJet
342    fNCharged  = iJet->ChargedMultiplicity();
343    fNNeutrals = iJet->NeutralMultiplicity();
344  
345 <  fDRMean    = JetTools::dRMean(iJet,-1);
346 <  fFrac01    = JetTools::frac  (iJet,0.1,0. ,-1);
347 <  fFrac02    = JetTools::frac  (iJet,0.2,0.1,-1);
348 <  fFrac03    = JetTools::frac  (iJet,0.3,0.2,-1);
349 <  fFrac04    = JetTools::frac  (iJet,0.4,0.3,-1);
350 <  fFrac05    = JetTools::frac  (iJet,0.5,0.4,-1);
345 >  fDRMean    = JetTools::dRMean (iJet,-1);
346 >  fDR2Mean   = JetTools::dR2Mean(iJet,-1);
347 >  fFrac01    = JetTools::frac   (iJet,0.1,0. ,-1);
348 >  fFrac02    = JetTools::frac   (iJet,0.2,0.1,-1);
349 >  fFrac03    = JetTools::frac   (iJet,0.3,0.2,-1);
350 >  fFrac04    = JetTools::frac   (iJet,0.4,0.3,-1);
351 >  fFrac05    = JetTools::frac   (iJet,0.5,0.4,-1);
352  
353    double lMVA = 0;
354 <  if(fJPt1 < 10) lMVA = fReader->EvaluateMVA( fLowPtMethodName  );
355 <  if(fJPt1 > 10) lMVA = fReader->EvaluateMVA( fHighPtMethodName );
354 >  if(fJPt1 < 10) lMVA = fLowPtReader->EvaluateMVA( fLowPtMethodName  );
355 >  if(fJPt1 > 10) lMVA = fReader     ->EvaluateMVA( fHighPtMethodName );
356    
357    if (printDebug == kTRUE) {
358      std::cout << "Debug Jet MVA: "
# Line 314 | Line 371 | Double_t JetIDMVA::MVAValue(const PFJet
371                << fFrac02    << " "
372                << fFrac03    << " "
373                << fFrac04    << " "
374 <              << fFrac05    
374 >              << fFrac05    << " "
375 >              << fDRMean    
376                << " === : === "
377                << lMVA << " "    
378                << std::endl;
# Line 322 | Line 380 | Double_t JetIDMVA::MVAValue(const PFJet
380  
381    return lMVA;
382   }
383 < Double_t JetIDMVA::correctedPt(const PFJet *iJet, FactorizedJetCorrector *iJetCorrector,const PileupEnergyDensityCol *iPUEnergyDensity) {
383 > Double_t JetIDMVA::correctedPt(const PFJet *iJet, FactorizedJetCorrector *iJetCorrector,
384 >                               const PileupEnergyDensityCol *iPUEnergyDensity,
385 >                               RhoUtilities::RhoType type) {
386 >  Double_t Rho = 0.0;
387 >  switch(type) {
388 >  case RhoUtilities::MIT_RHO_VORONOI_HIGH_ETA:
389 >    Rho = iPUEnergyDensity->At(0)->Rho();
390 >    break;
391 >  case RhoUtilities::MIT_RHO_VORONOI_LOW_ETA:
392 >    Rho = iPUEnergyDensity->At(0)->RhoLowEta();
393 >    break;
394 >  case RhoUtilities::MIT_RHO_RANDOM_HIGH_ETA:
395 >    Rho = iPUEnergyDensity->At(0)->RhoRandom();
396 >    break;
397 >  case RhoUtilities::MIT_RHO_RANDOM_LOW_ETA:
398 >    Rho = iPUEnergyDensity->At(0)->RhoRandomLowEta();
399 >    break;
400 >  case RhoUtilities::CMS_RHO_RHOKT6PFJETS:
401 >    Rho = iPUEnergyDensity->At(0)->RhoKt6PFJets();
402 >    break;
403 >  default:
404 >    // use the old default
405 >    Rho = iPUEnergyDensity->At(0)->Rho();
406 >    break;
407 >  }
408 >    
409    const FourVectorM rawMom = iJet->RawMom();
410    iJetCorrector->setJetEta(rawMom.Eta());
411    iJetCorrector->setJetPt (rawMom.Pt());
412    iJetCorrector->setJetPhi(rawMom.Phi());
413    iJetCorrector->setJetE  (rawMom.E());
414 <  iJetCorrector->setRho   (iPUEnergyDensity->At(0)->RhoHighEta());
414 >  iJetCorrector->setRho   (Rho);
415    iJetCorrector->setJetA  (iJet->JetArea());
416    iJetCorrector->setJetEMF(-99.0);    
417    Double_t correction = iJetCorrector->getCorrection();

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines