CMS 3D CMS Logo
AntiElectronIDMVA6.cc
Go to the documentation of this file.
1 #include "RecoTauTag/RecoTau/interface/AntiElectronIDMVA6.h"
2 
3 #include "FWCore/Utilities/interface/Exception.h"
4 
5 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
6 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
7 #include "DataFormats/TrackReco/interface/Track.h"
8 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
9 #include "DataFormats/MuonReco/interface/Muon.h"
10 
11 #include "CondFormats/DataRecord/interface/GBRWrapperRcd.h"
12 #include "FWCore/Framework/interface/ESHandle.h"
13 
14 #include "MagneticField/Engine/interface/MagneticField.h"
15 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
16 
17 #include <TMath.h>
18 #include <TFile.h>
19 #include <array>
20 
21 AntiElectronIDMVA6::AntiElectronIDMVA6(const edm::ParameterSet& cfg)
22  : isInitialized_(false),
23  mva_NoEleMatch_woGwoGSF_BL_(nullptr),
24  mva_NoEleMatch_wGwoGSF_BL_(nullptr),
25  mva_woGwGSF_BL_(nullptr),
26  mva_wGwGSF_BL_(nullptr),
27  mva_NoEleMatch_woGwoGSF_EC_(nullptr),
28  mva_NoEleMatch_wGwoGSF_EC_(nullptr),
29  mva_woGwGSF_EC_(nullptr),
30  mva_wGwGSF_EC_(nullptr)
31 {
32  loadMVAfromDB_ = cfg.exists("loadMVAfromDB") ? cfg.getParameter<bool>("loadMVAfromDB"): false;
33  if ( !loadMVAfromDB_ ) {
34  if(cfg.exists("inputFileName")){
35  inputFileName_ = cfg.getParameter<edm::FileInPath>("inputFileName");
36  }else throw cms::Exception("MVA input not defined") << "Requested to load tau MVA input from ROOT file but no file provided in cfg file";
37 
38  }
39 
40  mvaName_NoEleMatch_woGwoGSF_BL_ = cfg.getParameter<std::string>("mvaName_NoEleMatch_woGwoGSF_BL");
41  mvaName_NoEleMatch_wGwoGSF_BL_ = cfg.getParameter<std::string>("mvaName_NoEleMatch_wGwoGSF_BL");
42  mvaName_woGwGSF_BL_ = cfg.getParameter<std::string>("mvaName_woGwGSF_BL");
43  mvaName_wGwGSF_BL_ = cfg.getParameter<std::string>("mvaName_wGwGSF_BL");
44  mvaName_NoEleMatch_woGwoGSF_EC_ = cfg.getParameter<std::string>("mvaName_NoEleMatch_woGwoGSF_EC");
45  mvaName_NoEleMatch_wGwoGSF_EC_ = cfg.getParameter<std::string>("mvaName_NoEleMatch_wGwoGSF_EC");
46  mvaName_woGwGSF_EC_ = cfg.getParameter<std::string>("mvaName_woGwGSF_EC");
47  mvaName_wGwGSF_EC_ = cfg.getParameter<std::string>("mvaName_wGwGSF_EC");
48 
49  Var_NoEleMatch_woGwoGSF_Barrel_ = new Float_t[10];
50  Var_NoEleMatch_wGwoGSF_Barrel_ = new Float_t[18];
51  Var_woGwGSF_Barrel_ = new Float_t[24];
52  Var_wGwGSF_Barrel_ = new Float_t[32];
53  Var_NoEleMatch_woGwoGSF_Endcap_ = new Float_t[9];
54  Var_NoEleMatch_wGwoGSF_Endcap_ = new Float_t[17];
55  Var_woGwGSF_Endcap_ = new Float_t[23];
56  Var_wGwGSF_Endcap_ = new Float_t[31];
57 
58  bField_ = 0;
59  verbosity_ = 0;
60 }
61 
62 AntiElectronIDMVA6::~AntiElectronIDMVA6()
63 {
64  delete [] Var_NoEleMatch_woGwoGSF_Barrel_;
65  delete [] Var_NoEleMatch_wGwoGSF_Barrel_;
66  delete [] Var_woGwGSF_Barrel_;
67  delete [] Var_wGwGSF_Barrel_;
68  delete [] Var_NoEleMatch_woGwoGSF_Endcap_;
69  delete [] Var_NoEleMatch_wGwoGSF_Endcap_;
70  delete [] Var_woGwGSF_Endcap_;
71  delete [] Var_wGwGSF_Endcap_;
72 
73  if ( !loadMVAfromDB_ ){
74  delete mva_NoEleMatch_woGwoGSF_BL_;
75  delete mva_NoEleMatch_wGwoGSF_BL_;
76  delete mva_woGwGSF_BL_;
77  delete mva_wGwGSF_BL_;
78  delete mva_NoEleMatch_woGwoGSF_EC_;
79  delete mva_NoEleMatch_wGwoGSF_EC_;
80  delete mva_woGwGSF_EC_;
81  delete mva_wGwGSF_EC_;
82  }
83 
84  for ( std::vector<TFile*>::iterator it = inputFilesToDelete_.begin();
85  it != inputFilesToDelete_.end(); ++it ) {
86  delete (*it);
87  }
88 }
89 
90 namespace
91 {
92  const GBRForest* loadMVAfromFile(TFile* inputFile, const std::string& mvaName)
93  {
94  const GBRForest* mva = (GBRForest*)inputFile->Get(mvaName.data());
95  if ( !mva )
96  throw cms::Exception("PFRecoTauDiscriminationAgainstElectronMVA6::loadMVA")
97  << " Failed to load MVA = " << mvaName.data() << " from file " << " !!\n";
98 
99  return mva;
100  }
101 
102  const GBRForest* loadMVAfromDB(const edm::EventSetup& es, const std::string& mvaName)
103  {
104  edm::ESHandle<GBRForest> mva;
105  es.get<GBRWrapperRcd>().get(mvaName, mva);
106  return mva.product();
107  }
108 }
109 
110 void AntiElectronIDMVA6::beginEvent(const edm::Event& evt, const edm::EventSetup& es)
111 {
112  if ( !isInitialized_ ) {
113  if ( loadMVAfromDB_ ) {
114  mva_NoEleMatch_woGwoGSF_BL_ = loadMVAfromDB(es, mvaName_NoEleMatch_woGwoGSF_BL_);
115  mva_NoEleMatch_wGwoGSF_BL_ = loadMVAfromDB(es, mvaName_NoEleMatch_wGwoGSF_BL_);
116  mva_woGwGSF_BL_ = loadMVAfromDB(es, mvaName_woGwGSF_BL_);
117  mva_wGwGSF_BL_ = loadMVAfromDB(es, mvaName_wGwGSF_BL_);
118  mva_NoEleMatch_woGwoGSF_EC_ = loadMVAfromDB(es, mvaName_NoEleMatch_woGwoGSF_EC_);
119  mva_NoEleMatch_wGwoGSF_EC_ = loadMVAfromDB(es, mvaName_NoEleMatch_wGwoGSF_EC_);
120  mva_woGwGSF_EC_ = loadMVAfromDB(es, mvaName_woGwGSF_EC_);
121  mva_wGwGSF_EC_ = loadMVAfromDB(es, mvaName_wGwGSF_EC_);
122  } else {
123  if ( inputFileName_.location() == edm::FileInPath::Unknown ) throw cms::Exception("PFRecoTauDiscriminationAgainstElectronMVA6::loadMVA")
124  << " Failed to find File = " << inputFileName_ << " !!\n";
125  TFile* inputFile = new TFile(inputFileName_.fullPath().data());
126 
127  mva_NoEleMatch_woGwoGSF_BL_ = loadMVAfromFile(inputFile, mvaName_NoEleMatch_woGwoGSF_BL_);
128  mva_NoEleMatch_wGwoGSF_BL_ = loadMVAfromFile(inputFile, mvaName_NoEleMatch_wGwoGSF_BL_);
129  mva_woGwGSF_BL_ = loadMVAfromFile(inputFile, mvaName_woGwGSF_BL_);
130  mva_wGwGSF_BL_ = loadMVAfromFile(inputFile, mvaName_wGwGSF_BL_);
131  mva_NoEleMatch_woGwoGSF_EC_ = loadMVAfromFile(inputFile, mvaName_NoEleMatch_woGwoGSF_EC_);
132  mva_NoEleMatch_wGwoGSF_EC_ = loadMVAfromFile(inputFile, mvaName_NoEleMatch_wGwoGSF_EC_);
133  mva_woGwGSF_EC_ = loadMVAfromFile(inputFile, mvaName_woGwGSF_EC_);
134  mva_wGwGSF_EC_ = loadMVAfromFile(inputFile, mvaName_wGwGSF_EC_);
135  inputFilesToDelete_.push_back(inputFile);
136  }
137  isInitialized_ = true;
138  }
139 
140  edm::ESHandle<MagneticField> pSetup;
141  es.get<IdealMagneticFieldRecord>().get(pSetup);
142  bField_ = pSetup->inTesla(GlobalPoint(0,0,0)).z();
143 }
144 
145 double AntiElectronIDMVA6::MVAValue(Float_t TauPt,
146  Float_t TauEtaAtEcalEntrance,
147  Float_t TauPhi,
148  Float_t TauLeadChargedPFCandPt,
149  Float_t TauLeadChargedPFCandEtaAtEcalEntrance,
150  Float_t TauEmFraction,
151  Float_t TauLeadPFChargedHadrHoP,
152  Float_t TauLeadPFChargedHadrEoP,
153  Float_t TauVisMassIn,
154  Float_t TaudCrackEta,
155  Float_t TaudCrackPhi,
156  Float_t TauHasGsf,
157  Int_t TauSignalPFGammaCandsIn,
158  Int_t TauSignalPFGammaCandsOut,
159  const std::vector<Float_t>& GammasdEtaInSigCone,
160  const std::vector<Float_t>& GammasdPhiInSigCone,
161  const std::vector<Float_t>& GammasPtInSigCone,
162  const std::vector<Float_t>& GammasdEtaOutSigCone,
163  const std::vector<Float_t>& GammasdPhiOutSigCone,
164  const std::vector<Float_t>& GammasPtOutSigCone,
165  Float_t ElecEta,
166  Float_t ElecPhi,
167  Float_t ElecEtotOverPin,
168  Float_t ElecChi2NormGSF,
169  Float_t ElecChi2NormKF,
170  Float_t ElecGSFNumHits,
171  Float_t ElecKFNumHits,
172  Float_t ElecGSFTrackResol,
173  Float_t ElecGSFTracklnPt,
174  Float_t ElecPin,
175  Float_t ElecPout,
176  Float_t ElecEecal,
177  Float_t ElecDeltaEta,
178  Float_t ElecDeltaPhi,
179  Float_t ElecMvaInSigmaEtaEta,
180  Float_t ElecMvaInHadEnergy,
181  Float_t ElecMvaInDeltaEta)
182 {
183  double sumPt = 0.;
184  double dEta2 = 0.;
185  double dPhi2 = 0.;
186  double sumPt2 = 0.;
187  for ( unsigned int i = 0 ; i < GammasPtInSigCone.size() ; ++i ) {
188  double pt_i = GammasPtInSigCone[i];
189  double phi_i = GammasdPhiInSigCone[i];
190  if ( GammasdPhiInSigCone[i] > M_PI ) phi_i = GammasdPhiInSigCone[i] - 2*M_PI;
191  else if ( GammasdPhiInSigCone[i] < -M_PI ) phi_i = GammasdPhiInSigCone[i] + 2*M_PI;
192  double eta_i = GammasdEtaInSigCone[i];
193  sumPt += pt_i;
194  sumPt2 += (pt_i*pt_i);
195  dEta2 += (pt_i*eta_i*eta_i);
196  dPhi2 += (pt_i*phi_i*phi_i);
197  }
198  Float_t TauGammaEnFracIn = -99.;
199  if ( TauPt > 0. ) {
200  TauGammaEnFracIn = sumPt/TauPt;
201  }
202  if ( sumPt > 0. ) {
203  dEta2 /= sumPt;
204  dPhi2 /= sumPt;
205  }
206  Float_t TauGammaEtaMomIn = std::sqrt(dEta2)*std::sqrt(TauGammaEnFracIn)*TauPt;
207  Float_t TauGammaPhiMomIn = std::sqrt(dPhi2)*std::sqrt(TauGammaEnFracIn)*TauPt;
208 
209  sumPt = 0.;
210  dEta2 = 0.;
211  dPhi2 = 0.;
212  sumPt2 = 0.;
213  for ( unsigned int i = 0 ; i < GammasPtOutSigCone.size() ; ++i ) {
214  double pt_i = GammasPtOutSigCone[i];
215  double phi_i = GammasdPhiOutSigCone[i];
216  if ( GammasdPhiOutSigCone[i] > M_PI ) phi_i = GammasdPhiOutSigCone[i] - 2*M_PI;
217  else if ( GammasdPhiOutSigCone[i] < -M_PI ) phi_i = GammasdPhiOutSigCone[i] + 2*M_PI;
218  double eta_i = GammasdEtaOutSigCone[i];
219  sumPt += pt_i;
220  sumPt2 += (pt_i*pt_i);
221  dEta2 += (pt_i*eta_i*eta_i);
222  dPhi2 += (pt_i*phi_i*phi_i);
223  }
224  Float_t TauGammaEnFracOut = sumPt/TauPt;
225  if ( sumPt > 0. ) {
226  dEta2 /= sumPt;
227  dPhi2 /= sumPt;
228  }
229  Float_t TauGammaEtaMomOut = std::sqrt(dEta2)*std::sqrt(TauGammaEnFracOut)*TauPt;
230  Float_t TauGammaPhiMomOut = std::sqrt(dPhi2)*std::sqrt(TauGammaEnFracOut)*TauPt;
231 
232  return MVAValue(TauPt,
233  TauEtaAtEcalEntrance,
234  TauPhi,
235  TauLeadChargedPFCandPt,
236  TauLeadChargedPFCandEtaAtEcalEntrance,
237  TauEmFraction,
238  TauLeadPFChargedHadrHoP,
239  TauLeadPFChargedHadrEoP,
240  TauVisMassIn,
241  TaudCrackEta,
242  TaudCrackPhi,
243  TauHasGsf,
244  TauSignalPFGammaCandsIn,
245  TauSignalPFGammaCandsOut,
246  TauGammaEtaMomIn,
247  TauGammaEtaMomOut,
248  TauGammaPhiMomIn,
249  TauGammaPhiMomOut,
250  TauGammaEnFracIn,
251  TauGammaEnFracOut,
252  ElecEta,
253  ElecPhi,
254  ElecEtotOverPin,
255  ElecChi2NormGSF,
256  ElecChi2NormKF,
257  ElecGSFNumHits,
258  ElecKFNumHits,
259  ElecGSFTrackResol,
260  ElecGSFTracklnPt,
261  ElecPin,
262  ElecPout,
263  ElecEecal,
264  ElecDeltaEta,
265  ElecDeltaPhi,
266  ElecMvaInSigmaEtaEta,
267  ElecMvaInHadEnergy,
268  ElecMvaInDeltaEta);
269 }
270 
271 double AntiElectronIDMVA6::MVAValue(Float_t TauPt,
272  Float_t TauEtaAtEcalEntrance,
273  Float_t TauPhi,
274  Float_t TauLeadChargedPFCandPt,
275  Float_t TauLeadChargedPFCandEtaAtEcalEntrance,
276  Float_t TauEmFraction,
277  Float_t TauLeadPFChargedHadrHoP,
278  Float_t TauLeadPFChargedHadrEoP,
279  Float_t TauVisMassIn,
280  Float_t TaudCrackEta,
281  Float_t TaudCrackPhi,
282  Float_t TauHasGsf,
283  Int_t TauSignalPFGammaCandsIn,
284  Int_t TauSignalPFGammaCandsOut,
285  Float_t TauGammaEtaMomIn,
286  Float_t TauGammaEtaMomOut,
287  Float_t TauGammaPhiMomIn,
288  Float_t TauGammaPhiMomOut,
289  Float_t TauGammaEnFracIn,
290  Float_t TauGammaEnFracOut,
291  Float_t ElecEta,
292  Float_t ElecPhi,
293  Float_t ElecEtotOverPin,
294  Float_t ElecChi2NormGSF,
295  Float_t ElecChi2NormKF,
296  Float_t ElecGSFNumHits,
297  Float_t ElecKFNumHits,
298  Float_t ElecGSFTrackResol,
299  Float_t ElecGSFTracklnPt,
300  Float_t ElecPin,
301  Float_t ElecPout,
302  Float_t ElecEecal,
303  Float_t ElecDeltaEta,
304  Float_t ElecDeltaPhi,
305  Float_t ElecMvaInSigmaEtaEta,
306  Float_t ElecMvaInHadEnergy,
307  Float_t ElecMvaInDeltaEta)
308 {
309 
310  if ( !isInitialized_ ) {
311  throw cms::Exception("ClassNotInitialized")
312  << " AntiElectronMVA not properly initialized !!\n";
313  }
314 
315  double mvaValue = -99.;
316 
317  const float ECALBarrelEndcapEtaBorder = 1.479;
318  float ElecDeltaPinPoutOverPin = (ElecPin > 0.0) ? (std::abs(ElecPin - ElecPout)/ElecPin) : 1.0;
319  float ElecEecalOverPout = (ElecPout > 0.0) ? (ElecEecal/ElecPout) : 20.0;
320  float ElecNumHitsDiffOverSum = ((ElecGSFNumHits + ElecKFNumHits) > 0.0) ?
321  ((ElecGSFNumHits - ElecKFNumHits)/(ElecGSFNumHits + ElecKFNumHits)) : 1.0;
322 
323  if ( deltaR(TauEtaAtEcalEntrance, TauPhi, ElecEta, ElecPhi) > 0.3 && TauSignalPFGammaCandsIn == 0 && TauHasGsf < 0.5) {
324  if ( std::abs(TauEtaAtEcalEntrance) < ECALBarrelEndcapEtaBorder ){
325  Var_NoEleMatch_woGwoGSF_Barrel_[0] = TauEtaAtEcalEntrance;
326  Var_NoEleMatch_woGwoGSF_Barrel_[1] = TauLeadChargedPFCandEtaAtEcalEntrance;
327  Var_NoEleMatch_woGwoGSF_Barrel_[2] = std::min(float(2.), TauLeadChargedPFCandPt/std::max(float(1.), TauPt));
328  Var_NoEleMatch_woGwoGSF_Barrel_[3] = std::log(std::max(float(1.), TauPt));
329  Var_NoEleMatch_woGwoGSF_Barrel_[4] = TauEmFraction;
330  Var_NoEleMatch_woGwoGSF_Barrel_[5] = TauLeadPFChargedHadrHoP;
331  Var_NoEleMatch_woGwoGSF_Barrel_[6] = TauLeadPFChargedHadrEoP;
332  Var_NoEleMatch_woGwoGSF_Barrel_[7] = TauVisMassIn;
333  Var_NoEleMatch_woGwoGSF_Barrel_[8] = TaudCrackEta;
334  Var_NoEleMatch_woGwoGSF_Barrel_[9] = TaudCrackPhi;
335  mvaValue = mva_NoEleMatch_woGwoGSF_BL_->GetClassifier(Var_NoEleMatch_woGwoGSF_Barrel_);
336  } else {
337  Var_NoEleMatch_woGwoGSF_Endcap_[0] = TauEtaAtEcalEntrance;
338  Var_NoEleMatch_woGwoGSF_Endcap_[1] = TauLeadChargedPFCandEtaAtEcalEntrance;
339  Var_NoEleMatch_woGwoGSF_Endcap_[2] = std::min(float(2.), TauLeadChargedPFCandPt/std::max(float(1.), TauPt));
340  Var_NoEleMatch_woGwoGSF_Endcap_[3] = std::log(std::max(float(1.), TauPt));
341  Var_NoEleMatch_woGwoGSF_Endcap_[4] = TauEmFraction;
342  Var_NoEleMatch_woGwoGSF_Endcap_[5] = TauLeadPFChargedHadrHoP;
343  Var_NoEleMatch_woGwoGSF_Endcap_[6] = TauLeadPFChargedHadrEoP;
344  Var_NoEleMatch_woGwoGSF_Endcap_[7] = TauVisMassIn;
345  Var_NoEleMatch_woGwoGSF_Endcap_[8] = TaudCrackEta;
346  mvaValue = mva_NoEleMatch_woGwoGSF_EC_->GetClassifier(Var_NoEleMatch_woGwoGSF_Endcap_);
347  }
348  }
349  else if ( deltaR(TauEtaAtEcalEntrance, TauPhi, ElecEta, ElecPhi) > 0.3 && TauSignalPFGammaCandsIn > 0 && TauHasGsf < 0.5 ) {
350  if ( std::abs(TauEtaAtEcalEntrance) < ECALBarrelEndcapEtaBorder ){
351  Var_NoEleMatch_wGwoGSF_Barrel_[0] = TauEtaAtEcalEntrance;
352  Var_NoEleMatch_wGwoGSF_Barrel_[1] = TauLeadChargedPFCandEtaAtEcalEntrance;
353  Var_NoEleMatch_wGwoGSF_Barrel_[2] = std::min(float(2.), TauLeadChargedPFCandPt/std::max(float(1.), TauPt));
354  Var_NoEleMatch_wGwoGSF_Barrel_[3] = std::log(std::max(float(1.), TauPt));
355  Var_NoEleMatch_wGwoGSF_Barrel_[4] = TauEmFraction;
356  Var_NoEleMatch_wGwoGSF_Barrel_[5] = TauSignalPFGammaCandsIn;
357  Var_NoEleMatch_wGwoGSF_Barrel_[6] = TauSignalPFGammaCandsOut;
358  Var_NoEleMatch_wGwoGSF_Barrel_[7] = TauLeadPFChargedHadrHoP;
359  Var_NoEleMatch_wGwoGSF_Barrel_[8] = TauLeadPFChargedHadrEoP;
360  Var_NoEleMatch_wGwoGSF_Barrel_[9] = TauVisMassIn;
361  Var_NoEleMatch_wGwoGSF_Barrel_[10] = TauGammaEtaMomIn;
362  Var_NoEleMatch_wGwoGSF_Barrel_[11] = TauGammaEtaMomOut;
363  Var_NoEleMatch_wGwoGSF_Barrel_[12] = TauGammaPhiMomIn;
364  Var_NoEleMatch_wGwoGSF_Barrel_[13] = TauGammaPhiMomOut;
365  Var_NoEleMatch_wGwoGSF_Barrel_[14] = TauGammaEnFracIn;
366  Var_NoEleMatch_wGwoGSF_Barrel_[15] = TauGammaEnFracOut;
367  Var_NoEleMatch_wGwoGSF_Barrel_[16] = TaudCrackEta;
368  Var_NoEleMatch_wGwoGSF_Barrel_[17] = TaudCrackPhi;
369  mvaValue = mva_NoEleMatch_wGwoGSF_BL_->GetClassifier(Var_NoEleMatch_wGwoGSF_Barrel_);
370  } else {
371  Var_NoEleMatch_wGwoGSF_Endcap_[0] = TauEtaAtEcalEntrance;
372  Var_NoEleMatch_wGwoGSF_Endcap_[1] = TauLeadChargedPFCandEtaAtEcalEntrance;
373  Var_NoEleMatch_wGwoGSF_Endcap_[2] = std::min(float(2.), TauLeadChargedPFCandPt/std::max(float(1.), TauPt));
374  Var_NoEleMatch_wGwoGSF_Endcap_[3] = std::log(std::max(float(1.), TauPt));
375  Var_NoEleMatch_wGwoGSF_Endcap_[4] = TauEmFraction;
376  Var_NoEleMatch_wGwoGSF_Endcap_[5] = TauSignalPFGammaCandsIn;
377  Var_NoEleMatch_wGwoGSF_Endcap_[6] = TauSignalPFGammaCandsOut;
378  Var_NoEleMatch_wGwoGSF_Endcap_[7] = TauLeadPFChargedHadrHoP;
379  Var_NoEleMatch_wGwoGSF_Endcap_[8] = TauLeadPFChargedHadrEoP;
380  Var_NoEleMatch_wGwoGSF_Endcap_[9] = TauVisMassIn;
381  Var_NoEleMatch_wGwoGSF_Endcap_[10] = TauGammaEtaMomIn;
382  Var_NoEleMatch_wGwoGSF_Endcap_[11] = TauGammaEtaMomOut;
383  Var_NoEleMatch_wGwoGSF_Endcap_[12] = TauGammaPhiMomIn;
384  Var_NoEleMatch_wGwoGSF_Endcap_[13] = TauGammaPhiMomOut;
385  Var_NoEleMatch_wGwoGSF_Endcap_[14] = TauGammaEnFracIn;
386  Var_NoEleMatch_wGwoGSF_Endcap_[15] = TauGammaEnFracOut;
387  Var_NoEleMatch_wGwoGSF_Endcap_[16] = TaudCrackEta;
388  mvaValue = mva_NoEleMatch_wGwoGSF_EC_->GetClassifier(Var_NoEleMatch_wGwoGSF_Endcap_);
389  }
390  }
391  else if ( TauSignalPFGammaCandsIn == 0 && TauHasGsf > 0.5 ) {
392  if ( std::abs(TauEtaAtEcalEntrance) < ECALBarrelEndcapEtaBorder ) {
393  Var_woGwGSF_Barrel_[0] = std::max(float(-0.1), ElecEtotOverPin);
394  Var_woGwGSF_Barrel_[1] = std::log(ElecChi2NormGSF);
395  Var_woGwGSF_Barrel_[2] = ElecGSFNumHits;
396  Var_woGwGSF_Barrel_[3] = std::log(ElecGSFTrackResol);
397  Var_woGwGSF_Barrel_[4] = ElecGSFTracklnPt;
398  Var_woGwGSF_Barrel_[5] = ElecNumHitsDiffOverSum;
399  Var_woGwGSF_Barrel_[6] = std::log(ElecChi2NormKF);
400  Var_woGwGSF_Barrel_[7] = std::min(ElecDeltaPinPoutOverPin, float(1.));
401  Var_woGwGSF_Barrel_[8] = std::min(ElecEecalOverPout, float(20.));
402  Var_woGwGSF_Barrel_[9] = ElecDeltaEta;
403  Var_woGwGSF_Barrel_[10] = ElecDeltaPhi;
404  Var_woGwGSF_Barrel_[11] = std::min(ElecMvaInSigmaEtaEta, float(0.01));
405  Var_woGwGSF_Barrel_[12] = std::min(ElecMvaInHadEnergy, float(20.));
406  Var_woGwGSF_Barrel_[13] = std::min(ElecMvaInDeltaEta, float(0.1));
407  Var_woGwGSF_Barrel_[14] = TauEtaAtEcalEntrance;
408  Var_woGwGSF_Barrel_[15] = TauLeadChargedPFCandEtaAtEcalEntrance;
409  Var_woGwGSF_Barrel_[16] = std::min(float(2.), TauLeadChargedPFCandPt/std::max(float(1.), TauPt));
410  Var_woGwGSF_Barrel_[17] = std::log(std::max(float(1.), TauPt));
411  Var_woGwGSF_Barrel_[18] = TauEmFraction;
412  Var_woGwGSF_Barrel_[19] = TauLeadPFChargedHadrHoP;
413  Var_woGwGSF_Barrel_[20] = TauLeadPFChargedHadrEoP;
414  Var_woGwGSF_Barrel_[21] = TauVisMassIn;
415  Var_woGwGSF_Barrel_[22] = TaudCrackEta;
416  Var_woGwGSF_Barrel_[23] = TaudCrackPhi;
417  mvaValue = mva_woGwGSF_BL_->GetClassifier(Var_woGwGSF_Barrel_);
418  } else {
419  Var_woGwGSF_Endcap_[0] = std::max(float(-0.1), ElecEtotOverPin);
420  Var_woGwGSF_Endcap_[1] = std::log(ElecChi2NormGSF);
421  Var_woGwGSF_Endcap_[2] = ElecGSFNumHits;
422  Var_woGwGSF_Endcap_[3] = std::log(ElecGSFTrackResol);
423  Var_woGwGSF_Endcap_[4] = ElecGSFTracklnPt;
424  Var_woGwGSF_Endcap_[5] = ElecNumHitsDiffOverSum;
425  Var_woGwGSF_Endcap_[6] = std::log(ElecChi2NormKF);
426  Var_woGwGSF_Endcap_[7] = std::min(ElecDeltaPinPoutOverPin, float(1.));
427  Var_woGwGSF_Endcap_[8] = std::min(ElecEecalOverPout, float(20.));
428  Var_woGwGSF_Endcap_[9] = ElecDeltaEta;
429  Var_woGwGSF_Endcap_[10] = ElecDeltaPhi;
430  Var_woGwGSF_Endcap_[11] = std::min(ElecMvaInSigmaEtaEta, float(0.01));
431  Var_woGwGSF_Endcap_[12] = std::min(ElecMvaInHadEnergy, float(20.));
432  Var_woGwGSF_Endcap_[13] = std::min(ElecMvaInDeltaEta, float(0.1));
433  Var_woGwGSF_Endcap_[14] = TauEtaAtEcalEntrance;
434  Var_woGwGSF_Endcap_[15] = TauLeadChargedPFCandEtaAtEcalEntrance;
435  Var_woGwGSF_Endcap_[16] = std::min(float(2.), TauLeadChargedPFCandPt/std::max(float(1.), TauPt));
436  Var_woGwGSF_Endcap_[17] = std::log(std::max(float(1.), TauPt));
437  Var_woGwGSF_Endcap_[18] = TauEmFraction;
438  Var_woGwGSF_Endcap_[19] = TauLeadPFChargedHadrHoP;
439  Var_woGwGSF_Endcap_[20] = TauLeadPFChargedHadrEoP;
440  Var_woGwGSF_Endcap_[21] = TauVisMassIn;
441  Var_woGwGSF_Endcap_[22] = TaudCrackEta;
442  mvaValue = mva_woGwGSF_EC_->GetClassifier(Var_woGwGSF_Endcap_);
443  }
444  }
445  else if ( TauSignalPFGammaCandsIn > 0 && TauHasGsf > 0.5 ) {
446  if ( std::abs(TauEtaAtEcalEntrance) < ECALBarrelEndcapEtaBorder ) {
447  Var_wGwGSF_Barrel_[0] = std::max(float(-0.1), ElecEtotOverPin);
448  Var_wGwGSF_Barrel_[1] = std::log(ElecChi2NormGSF);
449  Var_wGwGSF_Barrel_[2] = ElecGSFNumHits;
450  Var_wGwGSF_Barrel_[3] = std::log(ElecGSFTrackResol);
451  Var_wGwGSF_Barrel_[4] = ElecGSFTracklnPt;
452  Var_wGwGSF_Barrel_[5] = ElecNumHitsDiffOverSum;
453  Var_wGwGSF_Barrel_[6] = std::log(ElecChi2NormKF);
454  Var_wGwGSF_Barrel_[7] = std::min(ElecDeltaPinPoutOverPin, float(1.));
455  Var_wGwGSF_Barrel_[8] = std::min(ElecEecalOverPout, float(20.));
456  Var_wGwGSF_Barrel_[9] = ElecDeltaEta;
457  Var_wGwGSF_Barrel_[10] = ElecDeltaPhi;
458  Var_wGwGSF_Barrel_[11] = std::min(ElecMvaInSigmaEtaEta, float(0.01));
459  Var_wGwGSF_Barrel_[12] = std::min(ElecMvaInHadEnergy, float(20.));
460  Var_wGwGSF_Barrel_[13] = std::min(ElecMvaInDeltaEta, float(0.1));
461  Var_wGwGSF_Barrel_[14] = TauEtaAtEcalEntrance;
462  Var_wGwGSF_Barrel_[15] = TauLeadChargedPFCandEtaAtEcalEntrance;
463  Var_wGwGSF_Barrel_[16] = std::min(float(2.), TauLeadChargedPFCandPt/std::max(float(1.), TauPt));
464  Var_wGwGSF_Barrel_[17] = std::log(std::max(float(1.), TauPt));
465  Var_wGwGSF_Barrel_[18] = TauEmFraction;
466  Var_wGwGSF_Barrel_[19] = TauSignalPFGammaCandsIn;
467  Var_wGwGSF_Barrel_[20] = TauSignalPFGammaCandsOut;
468  Var_wGwGSF_Barrel_[21] = TauLeadPFChargedHadrHoP;
469  Var_wGwGSF_Barrel_[22] = TauLeadPFChargedHadrEoP;
470  Var_wGwGSF_Barrel_[23] = TauVisMassIn;
471  Var_wGwGSF_Barrel_[24] = TauGammaEtaMomIn;
472  Var_wGwGSF_Barrel_[25] = TauGammaEtaMomOut;
473  Var_wGwGSF_Barrel_[26] = TauGammaPhiMomIn;
474  Var_wGwGSF_Barrel_[27] = TauGammaPhiMomOut;
475  Var_wGwGSF_Barrel_[28] = TauGammaEnFracIn;
476  Var_wGwGSF_Barrel_[29] = TauGammaEnFracOut;
477  Var_wGwGSF_Barrel_[30] = TaudCrackEta;
478  Var_wGwGSF_Barrel_[31] = TaudCrackPhi;
479  mvaValue = mva_wGwGSF_BL_->GetClassifier(Var_wGwGSF_Barrel_);
480  } else {
481  Var_wGwGSF_Endcap_[0] = std::max(float(-0.1), ElecEtotOverPin);
482  Var_wGwGSF_Endcap_[1] = std::log(ElecChi2NormGSF);
483  Var_wGwGSF_Endcap_[2] = ElecGSFNumHits;
484  Var_wGwGSF_Endcap_[3] = std::log(ElecGSFTrackResol);
485  Var_wGwGSF_Endcap_[4] = ElecGSFTracklnPt;
486  Var_wGwGSF_Endcap_[5] = ElecNumHitsDiffOverSum;
487  Var_wGwGSF_Endcap_[6] = std::log(ElecChi2NormKF);
488  Var_wGwGSF_Endcap_[7] = std::min(ElecDeltaPinPoutOverPin, float(1.));
489  Var_wGwGSF_Endcap_[8] = std::min(ElecEecalOverPout, float(20.));
490  Var_wGwGSF_Endcap_[9] = ElecDeltaEta;
491  Var_wGwGSF_Endcap_[10] = ElecDeltaPhi;
492  Var_wGwGSF_Endcap_[11] = std::min(ElecMvaInSigmaEtaEta, float(0.01));
493  Var_wGwGSF_Endcap_[12] = std::min(ElecMvaInHadEnergy, float(20.));
494  Var_wGwGSF_Endcap_[13] = std::min(ElecMvaInDeltaEta, float(0.1));
495  Var_wGwGSF_Endcap_[14] = TauEtaAtEcalEntrance;
496  Var_wGwGSF_Endcap_[15] = TauLeadChargedPFCandEtaAtEcalEntrance;
497  Var_wGwGSF_Endcap_[16] = std::min(float(2.), TauLeadChargedPFCandPt/std::max(float(1.), TauPt));
498  Var_wGwGSF_Endcap_[17] = std::log(std::max(float(1.), TauPt));
499  Var_wGwGSF_Endcap_[18] = TauEmFraction;
500  Var_wGwGSF_Endcap_[19] = TauSignalPFGammaCandsIn;
501  Var_wGwGSF_Endcap_[20] = TauSignalPFGammaCandsOut;
502  Var_wGwGSF_Endcap_[21] = TauLeadPFChargedHadrHoP;
503  Var_wGwGSF_Endcap_[22] = TauLeadPFChargedHadrEoP;
504  Var_wGwGSF_Endcap_[23] = TauVisMassIn;
505  Var_wGwGSF_Endcap_[24] = TauGammaEtaMomIn;
506  Var_wGwGSF_Endcap_[25] = TauGammaEtaMomOut;
507  Var_wGwGSF_Endcap_[26] = TauGammaPhiMomIn;
508  Var_wGwGSF_Endcap_[27] = TauGammaPhiMomOut;
509  Var_wGwGSF_Endcap_[28] = TauGammaEnFracIn;
510  Var_wGwGSF_Endcap_[29] = TauGammaEnFracOut;
511  Var_wGwGSF_Endcap_[30] = TaudCrackEta;
512  mvaValue = mva_wGwGSF_EC_->GetClassifier(Var_wGwGSF_Endcap_);
513  }
514  }
515  return mvaValue;
516 }
517 
518 double AntiElectronIDMVA6::MVAValue(const reco::PFTau& thePFTau,
519  const reco::GsfElectron& theGsfEle, bool usePhiAtEcalEntranceExtrapolation)
520 
521 {
522  // === tau variables ===
523  float TauEtaAtEcalEntrance = -99.;
524  float sumEtaTimesEnergy = 0.;
525  float sumEnergy = 0.;
526  const std::vector<reco::PFCandidatePtr>& signalPFCands = thePFTau.signalPFCands();
527  for ( std::vector<reco::PFCandidatePtr>::const_iterator pfCandidate = signalPFCands.begin();
528  pfCandidate != signalPFCands.end(); ++pfCandidate ) {
529  sumEtaTimesEnergy += (*pfCandidate)->positionAtECALEntrance().eta()*(*pfCandidate)->energy();
530  sumEnergy += (*pfCandidate)->energy();
531  }
532  if ( sumEnergy > 0. ) {
533  TauEtaAtEcalEntrance = sumEtaTimesEnergy/sumEnergy;
534  }
535 
536  float TauLeadChargedPFCandEtaAtEcalEntrance = -99.;
537  float TauLeadChargedPFCandPt = -99.;
538  for ( std::vector<reco::PFCandidatePtr>::const_iterator pfCandidate = signalPFCands.begin();
539  pfCandidate != signalPFCands.end(); ++pfCandidate ) {
540  const reco::Track* track = nullptr;
541  if ( (*pfCandidate)->trackRef().isNonnull() ) track = (*pfCandidate)->trackRef().get();
542  else if ( (*pfCandidate)->muonRef().isNonnull() && (*pfCandidate)->muonRef()->innerTrack().isNonnull() ) track = (*pfCandidate)->muonRef()->innerTrack().get();
543  else if ( (*pfCandidate)->muonRef().isNonnull() && (*pfCandidate)->muonRef()->globalTrack().isNonnull() ) track = (*pfCandidate)->muonRef()->globalTrack().get();
544  else if ( (*pfCandidate)->muonRef().isNonnull() && (*pfCandidate)->muonRef()->outerTrack().isNonnull() ) track = (*pfCandidate)->muonRef()->outerTrack().get();
545  else if ( (*pfCandidate)->gsfTrackRef().isNonnull() ) track = (*pfCandidate)->gsfTrackRef().get();
546  if ( track ) {
547  if ( track->pt() > TauLeadChargedPFCandPt ) {
548  TauLeadChargedPFCandEtaAtEcalEntrance = (*pfCandidate)->positionAtECALEntrance().eta();
549  TauLeadChargedPFCandPt = track->pt();
550  }
551  }
552  }
553 
554  Float_t TauPt = thePFTau.pt();
555  Float_t TauEmFraction = std::max(thePFTau.emFraction(), (Float_t)0.);
556  Float_t TauLeadPFChargedHadrHoP = 0.;
557  Float_t TauLeadPFChargedHadrEoP = 0.;
558  if ( thePFTau.leadPFChargedHadrCand()->p() > 0. ) {
559  TauLeadPFChargedHadrHoP = thePFTau.leadPFChargedHadrCand()->hcalEnergy()/thePFTau.leadPFChargedHadrCand()->p();
560  TauLeadPFChargedHadrEoP = thePFTau.leadPFChargedHadrCand()->ecalEnergy()/thePFTau.leadPFChargedHadrCand()->p();
561  }
562 
563  std::vector<Float_t> GammasdEtaInSigCone;
564  std::vector<Float_t> GammasdPhiInSigCone;
565  std::vector<Float_t> GammasPtInSigCone;
566  std::vector<Float_t> GammasdEtaOutSigCone;
567  std::vector<Float_t> GammasdPhiOutSigCone;
568  std::vector<Float_t> GammasPtOutSigCone;
569  reco::Candidate::LorentzVector pfGammaSum(0,0,0,0);
570  reco::Candidate::LorentzVector pfChargedSum(0,0,0,0);
571 
572  for ( unsigned i = 0 ; i < thePFTau.signalPFGammaCands().size(); ++i ) {
573  reco::PFCandidatePtr gamma = thePFTau.signalPFGammaCands().at(i);
574  float dR = deltaR(gamma->p4(), thePFTau.leadPFChargedHadrCand()->p4());
575  float signalrad = std::max(0.05, std::min(0.10, 3.0/std::max(1.0, thePFTau.pt())));
576 
577  // pfGammas inside the tau signal cone
578  if (dR < signalrad) {
579  if ( thePFTau.leadPFChargedHadrCand().isNonnull() ) {
580  GammasdEtaInSigCone.push_back(gamma->eta() - thePFTau.leadPFChargedHadrCand()->eta());
581  GammasdPhiInSigCone.push_back(gamma->phi() - thePFTau.leadPFChargedHadrCand()->phi());
582  }
583  else {
584  GammasdEtaInSigCone.push_back(gamma->eta() - thePFTau.eta());
585  GammasdPhiInSigCone.push_back(gamma->phi() - thePFTau.phi());
586  }
587  GammasPtInSigCone.push_back(gamma->pt());
588  pfGammaSum += gamma->p4();
589  }
590  // pfGammas outside the tau signal cone
591  else {
592  if ( thePFTau.leadPFChargedHadrCand().isNonnull() ) {
593  GammasdEtaOutSigCone.push_back(gamma->eta() - thePFTau.leadPFChargedHadrCand()->eta());
594  GammasdPhiOutSigCone.push_back(gamma->phi() - thePFTau.leadPFChargedHadrCand()->phi());
595  }
596  else {
597  GammasdEtaOutSigCone.push_back(gamma->eta() - thePFTau.eta());
598  GammasdPhiOutSigCone.push_back(gamma->phi() - thePFTau.phi());
599  }
600  GammasPtOutSigCone.push_back(gamma->pt());
601  }
602  }
603 
604  for ( unsigned i = 0 ; i < thePFTau.signalPFChargedHadrCands().size(); ++i ) {
605  reco::PFCandidatePtr charged = thePFTau.signalPFChargedHadrCands().at(i);
606  float dR = deltaR(charged->p4(), thePFTau.leadPFChargedHadrCand()->p4());
607  float signalrad = std::max(0.05, std::min(0.10, 3.0/std::max(1.0, thePFTau.pt())));
608 
609  // charged particles inside the tau signal cone
610  if (dR < signalrad) {
611  pfChargedSum += charged->p4();
612  }
613  }
614 
615  Int_t TauSignalPFGammaCandsIn = GammasPtInSigCone.size();
616  Int_t TauSignalPFGammaCandsOut = GammasPtOutSigCone.size();
617  Float_t TauVisMassIn = (pfGammaSum + pfChargedSum).mass();
618 
619  Float_t TauPhi = thePFTau.phi();
620  float sumPhiTimesEnergy = 0.;
621  float sumEnergyPhi = 0.;
622  if ( !usePhiAtEcalEntranceExtrapolation ){
623  for (auto const& pfc : signalPFCands){
624  sumPhiTimesEnergy += pfc->positionAtECALEntrance().phi()*pfc->energy();
625  sumEnergyPhi += pfc->energy();
626  }
627  }
628  else{
629  TauPhi= -99.;
630  for (unsigned int o = 0; o < signalPFCands.size(); ++o ) {
631  reco::Candidate const* signalCand = signalPFCands[o].get();
632  float phi = thePFTau.phi();
633  math::XYZPoint aPos;
634  if ( atECalEntrance(signalCand, aPos) ) phi = aPos.Phi();
635  sumPhiTimesEnergy += phi*signalCand->energy();
636  sumEnergy += signalCand->energy();
637  }
638  }
639  if ( sumEnergyPhi > 0. ) {
640  TauPhi = sumPhiTimesEnergy/sumEnergyPhi;
641  }
642  Float_t TaudCrackPhi = dCrackPhi(TauPhi, TauEtaAtEcalEntrance);
643  Float_t TaudCrackEta = dCrackEta(TauEtaAtEcalEntrance);
644  Float_t TauHasGsf = thePFTau.leadPFChargedHadrCand()->gsfTrackRef().isNonnull();
645 
646 
647  // === electron variables ===
648  Float_t ElecEta = theGsfEle.eta();
649  Float_t ElecPhi = theGsfEle.phi();
650 
651  //Variables related to the electron Cluster
652  Float_t ElecEe = 0.;
653  Float_t ElecEgamma = 0.;
654  reco::SuperClusterRef pfSuperCluster = theGsfEle.superCluster();
655  if ( pfSuperCluster.isNonnull() && pfSuperCluster.isAvailable() ) {
656  for ( reco::CaloCluster_iterator pfCluster = pfSuperCluster->clustersBegin();
657  pfCluster != pfSuperCluster->clustersEnd(); ++pfCluster ) {
658  double pfClusterEn = (*pfCluster)->energy();
659  if ( pfCluster == pfSuperCluster->clustersBegin() ) ElecEe += pfClusterEn;
660  else ElecEgamma += pfClusterEn;
661  }
662  }
663 
664  Float_t ElecPin = std::sqrt(theGsfEle.trackMomentumAtVtx().Mag2());
665  Float_t ElecPout = std::sqrt(theGsfEle.trackMomentumOut().Mag2());
666  Float_t ElecEtotOverPin = (ElecPin > 0.0) ? ((ElecEe + ElecEgamma)/ElecPin) : -0.1;
667  Float_t ElecEecal = theGsfEle.ecalEnergy();
668  Float_t ElecDeltaEta = theGsfEle.deltaEtaSeedClusterTrackAtCalo();
669  Float_t ElecDeltaPhi = theGsfEle.deltaPhiSeedClusterTrackAtCalo();
670  Float_t ElecMvaInSigmaEtaEta = (theGsfEle).mvaInput().sigmaEtaEta;
671  Float_t ElecMvaInHadEnergy = (theGsfEle).mvaInput().hadEnergy;
672  Float_t ElecMvaInDeltaEta = (theGsfEle).mvaInput().deltaEta;
673 
674  //Variables related to the GsfTrack
675  Float_t ElecChi2NormGSF = -99.;
676  Float_t ElecGSFNumHits = -99.;
677  Float_t ElecGSFTrackResol = -99.;
678  Float_t ElecGSFTracklnPt = -99.;
679  if ( theGsfEle.gsfTrack().isNonnull() ) {
680  ElecChi2NormGSF = (theGsfEle).gsfTrack()->normalizedChi2();
681  ElecGSFNumHits = (theGsfEle).gsfTrack()->numberOfValidHits();
682  if ( theGsfEle.gsfTrack()->pt() > 0. ) {
683  ElecGSFTrackResol = theGsfEle.gsfTrack()->ptError()/theGsfEle.gsfTrack()->pt();
684  ElecGSFTracklnPt = log(theGsfEle.gsfTrack()->pt())*M_LN10;
685  }
686  }
687 
688  //Variables related to the CtfTrack
689  Float_t ElecChi2NormKF = -99.;
690  Float_t ElecKFNumHits = -99.;
691  if ( theGsfEle.closestCtfTrackRef().isNonnull() ) {
692  ElecChi2NormKF = (theGsfEle).closestCtfTrackRef()->normalizedChi2();
693  ElecKFNumHits = (theGsfEle).closestCtfTrackRef()->numberOfValidHits();
694  }
695 
696  return MVAValue(TauPt,
697  TauEtaAtEcalEntrance,
698  TauPhi,
699  TauLeadChargedPFCandPt,
700  TauLeadChargedPFCandEtaAtEcalEntrance,
701  TauEmFraction,
702  TauLeadPFChargedHadrHoP,
703  TauLeadPFChargedHadrEoP,
704  TauVisMassIn,
705  TaudCrackEta,
706  TaudCrackPhi,
707  TauHasGsf,
708  TauSignalPFGammaCandsIn,
709  TauSignalPFGammaCandsOut,
710  GammasdEtaInSigCone,
711  GammasdPhiInSigCone,
712  GammasPtInSigCone,
713  GammasdEtaOutSigCone,
714  GammasdPhiOutSigCone,
715  GammasPtOutSigCone,
716  ElecEta,
717  ElecPhi,
718  ElecEtotOverPin,
719  ElecChi2NormGSF,
720  ElecChi2NormKF,
721  ElecGSFNumHits,
722  ElecKFNumHits,
723  ElecGSFTrackResol,
724  ElecGSFTracklnPt,
725  ElecPin,
726  ElecPout,
727  ElecEecal,
728  ElecDeltaEta,
729  ElecDeltaPhi,
730  ElecMvaInSigmaEtaEta,
731  ElecMvaInHadEnergy,
732  ElecMvaInDeltaEta);
733 }
734 
735 double AntiElectronIDMVA6::MVAValue(const reco::PFTau& thePFTau, bool usePhiAtEcalEntranceExtrapolation)
736 {
737  // === tau variables ===
738  float TauEtaAtEcalEntrance = -99.;
739  float sumEtaTimesEnergy = 0.;
740  float sumEnergy = 0.;
741  const std::vector<reco::PFCandidatePtr>& signalPFCands = thePFTau.signalPFCands();
742  for ( std::vector<reco::PFCandidatePtr>::const_iterator pfCandidate = signalPFCands.begin();
743  pfCandidate != signalPFCands.end(); ++pfCandidate ) {
744  sumEtaTimesEnergy += (*pfCandidate)->positionAtECALEntrance().eta()*(*pfCandidate)->energy();
745  sumEnergy += (*pfCandidate)->energy();
746  }
747  if ( sumEnergy > 0. ) {
748  TauEtaAtEcalEntrance = sumEtaTimesEnergy/sumEnergy;
749  }
750 
751  float TauLeadChargedPFCandEtaAtEcalEntrance = -99.;
752  float TauLeadChargedPFCandPt = -99.;
753  for ( std::vector<reco::PFCandidatePtr>::const_iterator pfCandidate = signalPFCands.begin();
754  pfCandidate != signalPFCands.end(); ++pfCandidate ) {
755  const reco::Track* track = nullptr;
756  if ( (*pfCandidate)->trackRef().isNonnull() ) track = (*pfCandidate)->trackRef().get();
757  else if ( (*pfCandidate)->muonRef().isNonnull() && (*pfCandidate)->muonRef()->innerTrack().isNonnull() ) track = (*pfCandidate)->muonRef()->innerTrack().get();
758  else if ( (*pfCandidate)->muonRef().isNonnull() && (*pfCandidate)->muonRef()->globalTrack().isNonnull() ) track = (*pfCandidate)->muonRef()->globalTrack().get();
759  else if ( (*pfCandidate)->muonRef().isNonnull() && (*pfCandidate)->muonRef()->outerTrack().isNonnull() ) track = (*pfCandidate)->muonRef()->outerTrack().get();
760  else if ( (*pfCandidate)->gsfTrackRef().isNonnull() ) track = (*pfCandidate)->gsfTrackRef().get();
761  if ( track ) {
762  if ( track->pt() > TauLeadChargedPFCandPt ) {
763  TauLeadChargedPFCandEtaAtEcalEntrance = (*pfCandidate)->positionAtECALEntrance().eta();
764  TauLeadChargedPFCandPt = track->pt();
765  }
766  }
767  }
768 
769  Float_t TauPt = thePFTau.pt();
770  Float_t TauEmFraction = std::max(thePFTau.emFraction(), (Float_t)0.);
771  Float_t TauLeadPFChargedHadrHoP = 0.;
772  Float_t TauLeadPFChargedHadrEoP = 0.;
773  if ( thePFTau.leadPFChargedHadrCand()->p() > 0. ) {
774  TauLeadPFChargedHadrHoP = thePFTau.leadPFChargedHadrCand()->hcalEnergy()/thePFTau.leadPFChargedHadrCand()->p();
775  TauLeadPFChargedHadrEoP = thePFTau.leadPFChargedHadrCand()->ecalEnergy()/thePFTau.leadPFChargedHadrCand()->p();
776  }
777 
778  std::vector<Float_t> GammasdEtaInSigCone;
779  std::vector<Float_t> GammasdPhiInSigCone;
780  std::vector<Float_t> GammasPtInSigCone;
781  std::vector<Float_t> GammasdEtaOutSigCone;
782  std::vector<Float_t> GammasdPhiOutSigCone;
783  std::vector<Float_t> GammasPtOutSigCone;
784  reco::Candidate::LorentzVector pfGammaSum(0,0,0,0);
785  reco::Candidate::LorentzVector pfChargedSum(0,0,0,0);
786 
787  for ( unsigned i = 0 ; i < thePFTau.signalPFGammaCands().size(); ++i ) {
788  reco::PFCandidatePtr gamma = thePFTau.signalPFGammaCands().at(i);
789  float dR = deltaR(gamma->p4(), thePFTau.leadPFChargedHadrCand()->p4());
790  float signalrad = std::max(0.05, std::min(0.10, 3.0/std::max(1.0, thePFTau.pt())));
791 
792  // pfGammas inside the tau signal cone
793  if (dR < signalrad) {
794  if ( thePFTau.leadPFChargedHadrCand().isNonnull() ) {
795  GammasdEtaInSigCone.push_back(gamma->eta() - thePFTau.leadPFChargedHadrCand()->eta());
796  GammasdPhiInSigCone.push_back(gamma->phi() - thePFTau.leadPFChargedHadrCand()->phi());
797  }
798  else {
799  GammasdEtaInSigCone.push_back(gamma->eta() - thePFTau.eta());
800  GammasdPhiInSigCone.push_back(gamma->phi() - thePFTau.phi());
801  }
802  GammasPtInSigCone.push_back(gamma->pt());
803  pfGammaSum += gamma->p4();
804  }
805  // pfGammas outside the tau signal cone
806  else {
807  if ( thePFTau.leadPFChargedHadrCand().isNonnull() ) {
808  GammasdEtaOutSigCone.push_back(gamma->eta() - thePFTau.leadPFChargedHadrCand()->eta());
809  GammasdPhiOutSigCone.push_back(gamma->phi() - thePFTau.leadPFChargedHadrCand()->phi());
810  }
811  else {
812  GammasdEtaOutSigCone.push_back(gamma->eta() - thePFTau.eta());
813  GammasdPhiOutSigCone.push_back(gamma->phi() - thePFTau.phi());
814  }
815  GammasPtOutSigCone.push_back(gamma->pt());
816  }
817  }
818 
819  for ( unsigned i = 0 ; i < thePFTau.signalPFChargedHadrCands().size(); ++i ) {
820  reco::PFCandidatePtr charged = thePFTau.signalPFChargedHadrCands().at(i);
821  float dR = deltaR(charged->p4(), thePFTau.leadPFChargedHadrCand()->p4());
822  float signalrad = std::max(0.05, std::min(0.10, 3.0/std::max(1.0, thePFTau.pt())));
823 
824  // charged particles inside the tau signal cone
825  if (dR < signalrad) {
826  pfChargedSum += charged->p4();
827  }
828  }
829 
830  Int_t TauSignalPFGammaCandsIn = GammasPtInSigCone.size();
831  Int_t TauSignalPFGammaCandsOut = GammasPtOutSigCone.size();
832  Float_t TauVisMassIn = (pfGammaSum + pfChargedSum).mass();
833 
834  Float_t TauPhi = thePFTau.phi();
835  float sumPhiTimesEnergy = 0.;
836  float sumEnergyPhi = 0.;
837  if ( !usePhiAtEcalEntranceExtrapolation ){
838  for ( std::vector<reco::PFCandidatePtr>::const_iterator pfCandidate = signalPFCands.begin();
839  pfCandidate != signalPFCands.end(); ++pfCandidate ) {
840  sumPhiTimesEnergy += (*pfCandidate)->positionAtECALEntrance().phi()*(*pfCandidate)->energy();
841  sumEnergyPhi += (*pfCandidate)->energy();
842  }
843  }
844  else{
845  TauPhi= -99.;
846  for (unsigned int o = 0; o < signalPFCands.size(); o++ ) {
847  reco::Candidate const* signalCand = signalPFCands[o].get();
848  float phi = thePFTau.phi();
849  math::XYZPoint aPos;
850  if ( atECalEntrance(signalCand, aPos) == true ) phi = aPos.Phi();
851  sumPhiTimesEnergy += phi*signalCand->energy();
852  sumEnergy += signalCand->energy();
853  }
854  }
855  if ( sumEnergyPhi > 0. ) {
856  TauPhi = sumPhiTimesEnergy/sumEnergyPhi;
857  }
858  Float_t TaudCrackPhi = dCrackPhi(TauPhi, TauEtaAtEcalEntrance);
859  Float_t TaudCrackEta = dCrackEta(TauEtaAtEcalEntrance);
860  Float_t TauHasGsf = thePFTau.leadPFChargedHadrCand()->gsfTrackRef().isNonnull();
861 
862 
863  // === electron variables ===
864  Float_t dummyElecEta = 9.9;
865 
866  return MVAValue(TauPt,
867  TauEtaAtEcalEntrance,
868  TauPhi,
869  TauLeadChargedPFCandPt,
870  TauLeadChargedPFCandEtaAtEcalEntrance,
871  TauEmFraction,
872  TauLeadPFChargedHadrHoP,
873  TauLeadPFChargedHadrEoP,
874  TauVisMassIn,
875  TaudCrackEta,
876  TaudCrackPhi,
877  TauHasGsf,
878  TauSignalPFGammaCandsIn,
879  TauSignalPFGammaCandsOut,
880  GammasdEtaInSigCone,
881  GammasdPhiInSigCone,
882  GammasPtInSigCone,
883  GammasdEtaOutSigCone,
884  GammasdPhiOutSigCone,
885  GammasPtOutSigCone,
886  dummyElecEta,
887  0.,
888  0.,
889  0.,
890  0.,
891  0.,
892  0.,
893  0.,
894  0.,
895  0.,
896  0.,
897  0.,
898  0.,
899  0.,
900  0.,
901  0.,
902  0.);
903 }
904 
905 double AntiElectronIDMVA6::MVAValue(const pat::Tau& theTau, const pat::Electron& theEle, bool usePhiAtEcalEntranceExtrapolation)
906 {
907  // === tau variables ===
908  float TauEtaAtEcalEntrance = theTau.etaAtEcalEntrance();
909 
910  float TauLeadChargedPFCandEtaAtEcalEntrance = theTau.etaAtEcalEntranceLeadChargedCand();
911  float TauLeadChargedPFCandPt = theTau.ptLeadChargedCand();
912 
913  Float_t TauPt = theTau.pt();
914  //Float_t TauEmFraction = std::max(theTau.ecalEnergy()/(theTau.ecalEnergy()+theTau.hcalEnergy()), (Float_t)0.);
915  Float_t TauEmFraction = std::max(theTau.emFraction_MVA(), (Float_t)0.);
916  Float_t TauLeadPFChargedHadrHoP = 0.;
917  Float_t TauLeadPFChargedHadrEoP = 0.;
918  if ( theTau.leadChargedHadrCand()->p() > 0. ) {
919  TauLeadPFChargedHadrHoP = theTau.hcalEnergyLeadChargedHadrCand()/theTau.leadChargedHadrCand()->p();
920  TauLeadPFChargedHadrEoP = theTau.ecalEnergyLeadChargedHadrCand()/theTau.leadChargedHadrCand()->p();
921  }
922 
923  std::vector<Float_t> GammasdEtaInSigCone;
924  std::vector<Float_t> GammasdPhiInSigCone;
925  std::vector<Float_t> GammasPtInSigCone;
926  std::vector<Float_t> GammasdEtaOutSigCone;
927  std::vector<Float_t> GammasdPhiOutSigCone;
928  std::vector<Float_t> GammasPtOutSigCone;
929  reco::Candidate::LorentzVector pfGammaSum(0,0,0,0);
930  reco::Candidate::LorentzVector pfChargedSum(0,0,0,0);
931 
932  const reco::CandidatePtrVector signalGammaCands = theTau.signalGammaCands();
933  for ( reco::CandidatePtrVector::const_iterator gamma = signalGammaCands.begin(); gamma != signalGammaCands.end(); ++gamma ){
934  float dR = deltaR((*gamma)->p4(), theTau.leadChargedHadrCand()->p4());
935  float signalrad = std::max(0.05, std::min(0.10, 3.0/std::max(1.0, theTau.pt())));
936 
937  // pfGammas inside the tau signal cone
938  if (dR < signalrad) {
939  if ( theTau.leadChargedHadrCand().isNonnull() ) {
940  GammasdEtaInSigCone.push_back((*gamma)->eta() - theTau.leadChargedHadrCand()->eta());
941  GammasdPhiInSigCone.push_back((*gamma)->phi() - theTau.leadChargedHadrCand()->phi());
942  //A.-C. please check whether this change is safe against future trainings
943  //GammasdPhiInSigCone.push_back(deltaPhi((*gamma)->phi(), theTau.leadChargedHadrCand()->phi()));
944  }
945  else {
946  GammasdEtaInSigCone.push_back((*gamma)->eta() - theTau.eta());
947  GammasdPhiInSigCone.push_back((*gamma)->phi() - theTau.phi());
948  //A.-C. please check whether this change is safe against future trainings
949  //GammasdPhiInSigCone.push_back(deltaPhi((*gamma)->phi(), theTau.phi()));
950  }
951  GammasPtInSigCone.push_back((*gamma)->pt());
952  pfGammaSum += (*gamma)->p4();
953  }
954  // pfGammas outside the tau signal cone
955  else {
956  if ( theTau.leadChargedHadrCand().isNonnull() ) {
957  GammasdEtaOutSigCone.push_back((*gamma)->eta() - theTau.leadChargedHadrCand()->eta());
958  GammasdPhiOutSigCone.push_back((*gamma)->phi() - theTau.leadChargedHadrCand()->phi());
959  //A.-C. please check whether this change is safe against future trainings
960  //GammasdPhiOutSigCone.push_back(deltaPhi((*gamma)->phi(), theTau.leadChargedHadrCand()->phi()));
961  }
962  else {
963  GammasdEtaOutSigCone.push_back((*gamma)->eta() - theTau.eta());
964  GammasdPhiOutSigCone.push_back((*gamma)->phi() - theTau.phi());
965  //A.-C. please chaekc whether this change is safe against future trainings
966  //GammasdPhiOutSigCone.push_back(deltaPhi((*gamma)->phi(), theTau.phi()));
967  }
968  GammasPtOutSigCone.push_back((*gamma)->pt());
969  }
970  }
971 
972  const reco::CandidatePtrVector signalChargedCands = theTau.signalChargedHadrCands();
973  for ( reco::CandidatePtrVector::const_iterator charged = signalChargedCands.begin(); charged != signalChargedCands.end(); ++charged ){
974  float dR = deltaR((*charged)->p4(), theTau.leadChargedHadrCand()->p4());
975  float signalrad = std::max(0.05, std::min(0.10, 3.0/std::max(1.0, theTau.pt())));
976 
977  // charged particles inside the tau signal cone
978  if (dR < signalrad) {
979  pfChargedSum += (*charged)->p4();
980  }
981  }
982 
983  Int_t TauSignalPFGammaCandsIn = GammasPtInSigCone.size();
984  Int_t TauSignalPFGammaCandsOut = GammasPtOutSigCone.size();
985  Float_t TauVisMassIn = (pfGammaSum + pfChargedSum).mass();
986  Float_t TauPhi = -99.;
987  if ( usePhiAtEcalEntranceExtrapolation ) {
988  float sumPhiTimesEnergy = 0.;
989  float sumEnergy = 0.;
990  const reco::CandidatePtrVector signalCands = theTau.signalCands();
991  for (unsigned int o = 0; o < signalCands.size(); o++ ) {
992  reco::Candidate const* signalCand = signalCands[o].get();
993  float phi = theTau.phi();
994  math::XYZPoint aPos;
995  if ( atECalEntrance(signalCand, aPos) == true ) phi = aPos.Phi();
996  sumPhiTimesEnergy += phi*signalCand->energy();
997  sumEnergy += signalCand->energy();
998  }
999  if ( sumEnergy > 0. ) {
1000  TauPhi = sumPhiTimesEnergy/sumEnergy;
1001  }
1002  }
1003  else {
1004  TauPhi = theTau.phiAtEcalEntrance();
1005  }
1006 
1007  Float_t TaudCrackPhi = dCrackPhi(TauPhi, TauEtaAtEcalEntrance);
1008  Float_t TaudCrackEta = dCrackEta(TauEtaAtEcalEntrance);
1009 
1010  Float_t TauHasGsf = 0;
1011  pat::PackedCandidate const* packedLeadTauCand = dynamic_cast<pat::PackedCandidate const*>(theTau.leadChargedHadrCand().get());
1012  if( abs(packedLeadTauCand->pdgId()) == 11 ) TauHasGsf = 1;
1013 
1014  // === electron variables ===
1015  Float_t ElecEta = theEle.eta();
1016  Float_t ElecPhi = theEle.phi();
1017 
1018  //Variables related to the electron Cluster
1019  Float_t ElecEe = 0.;
1020  Float_t ElecEgamma = 0.;
1021  reco::SuperClusterRef pfSuperCluster = theEle.superCluster();
1022  if ( pfSuperCluster.isNonnull() && pfSuperCluster.isAvailable() ) {
1023  for ( reco::CaloCluster_iterator pfCluster = pfSuperCluster->clustersBegin(); pfCluster != pfSuperCluster->clustersEnd(); ++pfCluster ) {
1024  double pfClusterEn = (*pfCluster)->energy();
1025  if ( pfCluster == pfSuperCluster->clustersBegin() ) ElecEe += pfClusterEn;
1026  else ElecEgamma += pfClusterEn;
1027  }
1028  }
1029 
1030  Float_t ElecPin = std::sqrt(theEle.trackMomentumAtVtx().Mag2());
1031  Float_t ElecPout = std::sqrt(theEle.trackMomentumOut().Mag2());
1032  Float_t ElecEtotOverPin = (ElecPin > 0.0) ? ((ElecEe + ElecEgamma)/ElecPin) : -0.1;
1033  Float_t ElecEecal = theEle.ecalEnergy();
1034  Float_t ElecDeltaEta = theEle.deltaEtaSeedClusterTrackAtCalo();
1035  Float_t ElecDeltaPhi = theEle.deltaPhiSeedClusterTrackAtCalo();
1036  Float_t ElecMvaInSigmaEtaEta = (theEle).mvaInput().sigmaEtaEta;
1037  Float_t ElecMvaInHadEnergy = (theEle).mvaInput().hadEnergy;
1038  Float_t ElecMvaInDeltaEta = (theEle).mvaInput().deltaEta;
1039 
1040  //Variables related to the GsfTrack
1041  Float_t ElecChi2NormGSF = -99.;
1042  Float_t ElecGSFNumHits = -99.;
1043  Float_t ElecGSFTrackResol = -99.;
1044  Float_t ElecGSFTracklnPt = -99.;
1045  if ( theEle.gsfTrack().isNonnull() ) {
1046  ElecChi2NormGSF = (theEle).gsfTrack()->normalizedChi2();
1047  ElecGSFNumHits = (theEle).gsfTrack()->numberOfValidHits();
1048  if ( theEle.gsfTrack()->pt() > 0. ) {
1049  ElecGSFTrackResol = theEle.gsfTrack()->ptError()/theEle.gsfTrack()->pt();
1050  ElecGSFTracklnPt = log(theEle.gsfTrack()->pt())*M_LN10;
1051  }
1052  }
1053 
1054  //Variables related to the CtfTrack
1055  Float_t ElecChi2NormKF = -99.;
1056  Float_t ElecKFNumHits = -99.;
1057  if ( theEle.closestCtfTrackRef().isNonnull() ) {
1058  ElecChi2NormKF = (theEle).closestCtfTrackRef()->normalizedChi2();
1059  ElecKFNumHits = (theEle).closestCtfTrackRef()->numberOfValidHits();
1060  }
1061 
1062  return MVAValue(TauPt,
1063  TauEtaAtEcalEntrance,
1064  TauPhi,
1065  TauLeadChargedPFCandPt,
1066  TauLeadChargedPFCandEtaAtEcalEntrance,
1067  TauEmFraction,
1068  TauLeadPFChargedHadrHoP,
1069  TauLeadPFChargedHadrEoP,
1070  TauVisMassIn,
1071  TaudCrackEta,
1072  TaudCrackPhi,
1073  TauHasGsf,
1074  TauSignalPFGammaCandsIn,
1075  TauSignalPFGammaCandsOut,
1076  GammasdEtaInSigCone,
1077  GammasdPhiInSigCone,
1078  GammasPtInSigCone,
1079  GammasdEtaOutSigCone,
1080  GammasdPhiOutSigCone,
1081  GammasPtOutSigCone,
1082  ElecEta,
1083  ElecPhi,
1084  ElecEtotOverPin,
1085  ElecChi2NormGSF,
1086  ElecChi2NormKF,
1087  ElecGSFNumHits,
1088  ElecKFNumHits,
1089  ElecGSFTrackResol,
1090  ElecGSFTracklnPt,
1091  ElecPin,
1092  ElecPout,
1093  ElecEecal,
1094  ElecDeltaEta,
1095  ElecDeltaPhi,
1096  ElecMvaInSigmaEtaEta,
1097  ElecMvaInHadEnergy,
1098  ElecMvaInDeltaEta);
1099 }
1100 
1101 double AntiElectronIDMVA6::MVAValue(const pat::Tau& theTau, bool usePhiAtEcalEntranceExtrapolation)
1102 {
1103  // === tau variables ===
1104  float TauEtaAtEcalEntrance = theTau.etaAtEcalEntrance();
1105 
1106  float TauLeadChargedPFCandEtaAtEcalEntrance = theTau.etaAtEcalEntranceLeadChargedCand();
1107  float TauLeadChargedPFCandPt = theTau.ptLeadChargedCand();
1108 
1109  Float_t TauPt = theTau.pt();
1110  //Float_t TauEmFraction = std::max(theTau.ecalEnergy()/(theTau.ecalEnergy()+theTau.hcalEnergy()), (Float_t)0.);
1111  Float_t TauEmFraction = std::max(theTau.emFraction_MVA(), (Float_t)0.);
1112  Float_t TauLeadPFChargedHadrHoP = 0.;
1113  Float_t TauLeadPFChargedHadrEoP = 0.;
1114  if ( theTau.leadChargedHadrCand()->p() > 0. ) {
1115  TauLeadPFChargedHadrHoP = theTau.hcalEnergyLeadChargedHadrCand()/theTau.leadChargedHadrCand()->p();
1116  TauLeadPFChargedHadrEoP = theTau.ecalEnergyLeadChargedHadrCand()/theTau.leadChargedHadrCand()->p();
1117  }
1118 
1119  std::vector<Float_t> GammasdEtaInSigCone;
1120  std::vector<Float_t> GammasdPhiInSigCone;
1121  std::vector<Float_t> GammasPtInSigCone;
1122  std::vector<Float_t> GammasdEtaOutSigCone;
1123  std::vector<Float_t> GammasdPhiOutSigCone;
1124  std::vector<Float_t> GammasPtOutSigCone;
1125  reco::Candidate::LorentzVector pfGammaSum(0,0,0,0);
1126  reco::Candidate::LorentzVector pfChargedSum(0,0,0,0);
1127 
1128  const reco::CandidatePtrVector signalGammaCands = theTau.signalGammaCands();
1129  for ( reco::CandidatePtrVector::const_iterator gamma = signalGammaCands.begin(); gamma != signalGammaCands.end(); ++gamma ) {
1130  float dR = deltaR((*gamma)->p4(), theTau.leadChargedHadrCand()->p4());
1131  float signalrad = std::max(0.05, std::min(0.10, 3.0/std::max(1.0, theTau.pt())));
1132 
1133  // pfGammas inside the tau signal cone
1134  if (dR < signalrad) {
1135  if ( theTau.leadChargedHadrCand().isNonnull() ) {
1136  GammasdEtaInSigCone.push_back((*gamma)->eta() - theTau.leadChargedHadrCand()->eta());
1137  GammasdPhiInSigCone.push_back((*gamma)->phi() - theTau.leadChargedHadrCand()->phi());
1138  }
1139  else {
1140  GammasdEtaInSigCone.push_back((*gamma)->eta() - theTau.eta());
1141  GammasdPhiInSigCone.push_back((*gamma)->phi() - theTau.phi());
1142  }
1143  GammasPtInSigCone.push_back((*gamma)->pt());
1144  pfGammaSum += (*gamma)->p4();
1145  }
1146  // pfGammas outside the tau signal cone
1147  else {
1148  if ( theTau.leadChargedHadrCand().isNonnull() ) {
1149  GammasdEtaOutSigCone.push_back((*gamma)->eta() - theTau.leadChargedHadrCand()->eta());
1150  GammasdPhiOutSigCone.push_back((*gamma)->phi() - theTau.leadChargedHadrCand()->phi());
1151  }
1152  else {
1153  GammasdEtaOutSigCone.push_back((*gamma)->eta() - theTau.eta());
1154  GammasdPhiOutSigCone.push_back((*gamma)->phi() - theTau.phi());
1155  }
1156  GammasPtOutSigCone.push_back((*gamma)->pt());
1157  }
1158  }
1159 
1160  const reco::CandidatePtrVector signalChargedCands = theTau.signalChargedHadrCands();
1161  for ( reco::CandidatePtrVector::const_iterator charged = signalChargedCands.begin(); charged != signalChargedCands.end(); ++charged ) {
1162  float dR = deltaR((*charged)->p4(), theTau.leadChargedHadrCand()->p4());
1163  float signalrad = std::max(0.05, std::min(0.10, 3.0/std::max(1.0, theTau.pt())));
1164 
1165  // charged particles inside the tau signal cone
1166  if (dR < signalrad) {
1167  pfChargedSum += (*charged)->p4();
1168  }
1169  }
1170 
1171  Int_t TauSignalPFGammaCandsIn = GammasPtInSigCone.size();
1172  Int_t TauSignalPFGammaCandsOut = GammasPtOutSigCone.size();
1173  Float_t TauVisMassIn = (pfGammaSum + pfChargedSum).mass();
1174  Float_t TauPhi = -99.;
1175  if ( usePhiAtEcalEntranceExtrapolation ) {
1176  float sumPhiTimesEnergy = 0.;
1177  float sumEnergy = 0.;
1178  const reco::CandidatePtrVector signalCands = theTau.signalCands();
1179  for (unsigned int o = 0; o < signalCands.size(); o++ ) {
1180  reco::Candidate const* signalCand = signalCands[o].get();
1181  float phi = theTau.phi();
1182  math::XYZPoint aPos;
1183  if ( atECalEntrance(signalCand, aPos) == true ) phi = aPos.Phi();
1184  sumPhiTimesEnergy += phi*signalCand->energy();
1185  sumEnergy += signalCand->energy();
1186  }
1187  if ( sumEnergy > 0. ) {
1188  TauPhi = sumPhiTimesEnergy/sumEnergy;
1189  }
1190  }
1191  else {
1192  TauPhi = theTau.phiAtEcalEntrance();
1193  }
1194 
1195  Float_t TaudCrackPhi = dCrackPhi(TauPhi, TauEtaAtEcalEntrance);
1196  Float_t TaudCrackEta = dCrackEta(TauEtaAtEcalEntrance);
1197 
1198  Float_t TauHasGsf = 0;
1199  pat::PackedCandidate const* packedLeadTauCand = dynamic_cast<pat::PackedCandidate const*>(theTau.leadChargedHadrCand().get());
1200  //const reco::Track & pseudoTrack = packedLeadTauCand->pseudoTrack();
1201  if( abs(packedLeadTauCand->pdgId()) == 11 ) TauHasGsf = 1;
1202 
1203  // === electron variables ===
1204  Float_t dummyElecEta = 9.9;
1205 
1206  return MVAValue(TauPt,
1207  TauEtaAtEcalEntrance,
1208  TauPhi,
1209  TauLeadChargedPFCandPt,
1210  TauLeadChargedPFCandEtaAtEcalEntrance,
1211  TauEmFraction,
1212  TauLeadPFChargedHadrHoP,
1213  TauLeadPFChargedHadrEoP,
1214  TauVisMassIn,
1215  TaudCrackEta,
1216  TaudCrackPhi,
1217  TauHasGsf,
1218  TauSignalPFGammaCandsIn,
1219  TauSignalPFGammaCandsOut,
1220  GammasdEtaInSigCone,
1221  GammasdPhiInSigCone,
1222  GammasPtInSigCone,
1223  GammasdEtaOutSigCone,
1224  GammasdPhiOutSigCone,
1225  GammasPtOutSigCone,
1226  dummyElecEta,
1227  0.,
1228  0.,
1229  0.,
1230  0.,
1231  0.,
1232  0.,
1233  0.,
1234  0.,
1235  0.,
1236  0.,
1237  0.,
1238  0.,
1239  0.,
1240  0.,
1241  0.,
1242  0.);
1243 }
1244 
1245 double AntiElectronIDMVA6::minimum(double a, double b)
1246 {
1247  if ( std::abs(b) < std::abs(a) ) return b;
1248  else return a;
1249 }
1250 
1251 namespace {
1252 
1253  // IN: define locations of the 18 phi-cracks
1254  std::array<double,18> fill_cPhi() {
1255  constexpr double pi = M_PI; // 3.14159265358979323846;
1256  std::array<double,18> cPhi;
1257  // IN: define locations of the 18 phi-cracks
1258  cPhi[0] = 2.97025;
1259  for ( unsigned iCrack = 1; iCrack <= 17; ++iCrack )
1260  cPhi[iCrack] = cPhi[0] - 2.*iCrack*pi/18;
1261  return cPhi;
1262  }
1263 
1264  const std::array<double,18> cPhi = fill_cPhi();
1265 
1266 }
1267 
1268 double AntiElectronIDMVA6::dCrackPhi(double phi, double eta)
1269 {
1270 //--- compute the (unsigned) distance to the closest phi-crack in the ECAL barrel
1271 
1272  constexpr double pi = M_PI; // 3.14159265358979323846;
1273 
1274  // IN: shift of this location if eta < 0
1275  constexpr double delta_cPhi = 0.00638;
1276 
1277  double retVal = 99.;
1278 
1279  if ( eta >= -1.47464 && eta <= 1.47464 ) {
1280 
1281  // the location is shifted
1282  if ( eta < 0. ) phi += delta_cPhi;
1283 
1284  // CV: need to bring-back phi into interval [-pi,+pi]
1285  if ( phi > pi ) phi -= 2.*pi;
1286  if ( phi < -pi ) phi += 2.*pi;
1287 
1288  if ( phi >= -pi && phi <= pi ) {
1289 
1290  // the problem of the extrema:
1291  if ( phi < cPhi[17] || phi >= cPhi[0] ) {
1292  if ( phi < 0. ) phi += 2.*pi;
1293  retVal = minimum(phi - cPhi[0], phi - cPhi[17] - 2.*pi);
1294  } else {
1295  // between these extrema...
1296  bool OK = false;
1297  unsigned iCrack = 16;
1298  while( !OK ) {
1299  if ( phi < cPhi[iCrack] ) {
1300  retVal = minimum(phi - cPhi[iCrack + 1], phi - cPhi[iCrack]);
1301  OK = true;
1302  } else {
1303  iCrack -= 1;
1304  }
1305  }
1306  }
1307  } else {
1308  retVal = 0.; // IN: if there is a problem, we assume that we are in a crack
1309  }
1310  } else {
1311  return -99.;
1312  }
1313 
1314  return std::abs(retVal);
1315 }
1316 
1317 double AntiElectronIDMVA6::dCrackEta(double eta)
1318 {
1319 //--- compute the (unsigned) distance to the closest eta-crack in the ECAL barrel
1320 
1321  // IN: define locations of the eta-cracks
1322  double cracks[5] = { 0., 4.44747e-01, 7.92824e-01, 1.14090e+00, 1.47464e+00 };
1323 
1324  double retVal = 99.;
1325 
1326  for ( int iCrack = 0; iCrack < 5 ; ++iCrack ) {
1327  double d = minimum(eta - cracks[iCrack], eta + cracks[iCrack]);
1328  if ( std::abs(d) < std::abs(retVal) ) {
1329  retVal = d;
1330  }
1331  }
1332 
1333  return std::abs(retVal);
1334 }
1335 
1336 bool AntiElectronIDMVA6::atECalEntrance(const reco::Candidate* part, math::XYZPoint &pos)
1337 {
1338  bool result = false;
1339  BaseParticlePropagator theParticle =
1340  BaseParticlePropagator(RawParticle(math::XYZTLorentzVector(part->px(),
1341  part->py(),
1342  part->pz(),
1343  part->energy()),
1344  math::XYZTLorentzVector(part->vertex().x(),
1345  part->vertex().y(),
1346  part->vertex().z(),
1347  0.)),
1348  0.,0.,bField_);
1349  theParticle.setCharge(part->charge());
1350  theParticle.propagateToEcalEntrance(false);
1351  if(theParticle.getSuccess()!=0){
1352  pos = math::XYZPoint(theParticle.vertex());
1353  result = true;
1354  }
1355  else {
1356  result = false;
1357  }
1358  return result;
1359 }
edm::Ref::isAvailable
bool isAvailable() const
Definition: Ref.h:577
RawParticle::setCharge
void setCharge(float q)
set the MEASURED charge
Definition: RawParticle.cc:139
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
reco::GsfElectron::gsfTrack
GsfTrackRef gsfTrack() const override
reference to a GsfTrack
Definition: GsfElectron.h:185
AntiElectronIDMVA6::minimum
double minimum(double a, double b)
Definition: AntiElectronIDMVA6.cc:1245
PATTauDiscriminationAgainstElectronMVA6_cfi.usePhiAtEcalEntranceExtrapolation
usePhiAtEcalEntranceExtrapolation
Definition: PATTauDiscriminationAgainstElectronMVA6_cfi.py:37
AntiElectronIDMVA6::mvaName_woGwGSF_BL_
std::string mvaName_woGwGSF_BL_
Definition: AntiElectronIDMVA6.h:150
reco::Candidate::pz
virtual double pz() const =0
z coordinate of momentum vector
AntiElectronIDMVA6::mvaName_NoEleMatch_woGwoGSF_BL_
std::string mvaName_NoEleMatch_woGwoGSF_BL_
Definition: AntiElectronIDMVA6.h:148
reco::GsfElectron::closestCtfTrackRef
virtual TrackRef closestCtfTrackRef() const
Definition: GsfElectron.h:201
edm::Ref::isNonnull
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:253
mps_fire.i
i
Definition: mps_fire.py:269
edm::PtrVectorBase::size
size_type size() const
Size of the RefVector.
Definition: PtrVectorBase.h:74
AntiElectronIDMVA6::atECalEntrance
bool atECalEntrance(const reco::Candidate *part, math::XYZPoint &pos)
Definition: AntiElectronIDMVA6.cc:1336
reco::LeafCandidate::eta
double eta() const final
momentum pseudorapidity
Definition: LeafCandidate.h:137
reco::PFTau::leadPFChargedHadrCand
const PFCandidatePtr & leadPFChargedHadrCand() const
Definition: PFTau.cc:67
mps_fire.result
result
Definition: mps_fire.py:226
pat::PackedCandidate::pdgId
int pdgId() const override
PDG identifier.
Definition: PackedCandidate.h:565
AntiElectronIDMVA6::mvaName_NoEleMatch_woGwoGSF_EC_
std::string mvaName_NoEleMatch_woGwoGSF_EC_
Definition: AntiElectronIDMVA6.h:152
pat::Tau::hcalEnergyLeadChargedHadrCand
float hcalEnergyLeadChargedHadrCand() const
return hcal energy from LeadChargedHadrCand
Definition: Tau.h:345
AlCaHLTBitMon_QueryRunRegistry.string
string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:255
edm::Ptr::get
T const * get() const
Returns C++ pointer to the item.
Definition: Ptr.h:159
AntiElectronIDMVA6::inputFileName_
edm::FileInPath inputFileName_
Definition: AntiElectronIDMVA6.h:146
AntiElectronIDMVA6::AntiElectronIDMVA6
AntiElectronIDMVA6(const edm::ParameterSet &)
Definition: AntiElectronIDMVA6.cc:21
AntiElectronIDMVA6::mva_wGwGSF_EC_
const GBRForest * mva_wGwGSF_EC_
Definition: AntiElectronIDMVA6.h:173
reco::GsfElectron::trackMomentumAtVtx
math::XYZVectorF trackMomentumAtVtx() const
Definition: GsfElectron.h:291
reco::PFTau::signalPFGammaCands
const std::vector< reco::PFCandidatePtr > & signalPFGammaCands() const
Gamma candidates in signal region.
Definition: PFTau.cc:84
GlobalPoint
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
edm::ParameterSet::exists
bool exists(std::string const &parameterName) const
checks if a parameter exists
Definition: ParameterSet.cc:760
AntiElectronIDMVA6::Var_NoEleMatch_wGwoGSF_Barrel_
Float_t * Var_NoEleMatch_wGwoGSF_Barrel_
Definition: AntiElectronIDMVA6.h:158
BaseParticlePropagator::getSuccess
int getSuccess() const
Has propagation been performed and was barrel or endcap reached ?
Definition: BaseParticlePropagator.h:297
AntiElectronIDMVA6::dCrackPhi
double dCrackPhi(double phi, double eta)
Definition: AntiElectronIDMVA6.cc:1268
reco::LeafCandidate::pt
double pt() const final
transverse momentum
Definition: LeafCandidate.h:131
reco::PFTau::emFraction
float emFraction() const
Definition: PFTau.cc:227
AntiElectronIDMVA6::mva_NoEleMatch_woGwoGSF_EC_
const GBRForest * mva_NoEleMatch_woGwoGSF_EC_
Definition: AntiElectronIDMVA6.h:170
pat::Tau::signalCands
reco::CandidatePtrVector signalCands() const
AntiElectronIDMVA6::mvaName_NoEleMatch_wGwoGSF_EC_
std::string mvaName_NoEleMatch_wGwoGSF_EC_
Definition: AntiElectronIDMVA6.h:153
nullptr
#define nullptr
taus_updatedMVAIds_cff.loadMVAfromDB
loadMVAfromDB
Definition: taus_updatedMVAIds_cff.py:15
AntiElectronIDMVA6::beginEvent
void beginEvent(const edm::Event &, const edm::EventSetup &)
Definition: AntiElectronIDMVA6.cc:110
AntiElectronIDMVA6::mva_NoEleMatch_wGwoGSF_BL_
const GBRForest * mva_NoEleMatch_wGwoGSF_BL_
Definition: AntiElectronIDMVA6.h:167
reco::PFTau::signalPFCands
const std::vector< reco::PFCandidatePtr > & signalPFCands() const
PFCandidates in signal region.
Definition: PFTau.cc:78
pat::Tau::signalGammaCands
reco::CandidatePtrVector signalGammaCands() const
constexpr
#define constexpr
AntiElectronIDMVA6::Var_NoEleMatch_woGwoGSF_Endcap_
Float_t * Var_NoEleMatch_woGwoGSF_Endcap_
Definition: AntiElectronIDMVA6.h:161
pi
const Double_t pi
Definition: trackSplitPlot.h:36
edm::PtrVector::begin
const_iterator begin() const
Definition: PtrVector.h:129
math::XYZTLorentzVector
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
AntiElectronIDMVA6::inputFilesToDelete_
std::vector< TFile * > inputFilesToDelete_
Definition: AntiElectronIDMVA6.h:175
reco::Candidate::energy
virtual double energy() const =0
energy
reco::Candidate::py
virtual double py() const =0
y coordinate of momentum vector
AntiElectronIDMVA6::mva_woGwGSF_EC_
const GBRForest * mva_woGwGSF_EC_
Definition: AntiElectronIDMVA6.h:172
CustomPhysics_cfi.gamma
gamma
Definition: CustomPhysics_cfi.py:17
reco::TrackBase::eta
double eta() const
pseudorapidity of momentum vector
Definition: TrackBase.h:651
reco::GsfElectron::deltaPhiSeedClusterTrackAtCalo
float deltaPhiSeedClusterTrackAtCalo() const
Definition: GsfElectron.h:253
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:18
reco::GsfElectron::trackMomentumOut
math::XYZVectorF trackMomentumOut() const
Definition: GsfElectron.h:293
reco::TrackBase::pt
double pt() const
track transverse momentum
Definition: TrackBase.h:621
DetachedQuadStep_cff.mva
mva
Definition: DetachedQuadStep_cff.py:209
PV3DBase::z
T z() const
Definition: PV3DBase.h:64
pat::Tau::ecalEnergyLeadChargedHadrCand
float ecalEnergyLeadChargedHadrCand() const
Definition: Tau.h:343
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
edm::PtrVector::end
const_iterator end() const
Definition: PtrVector.h:134
min
T min(T a, T b)
Definition: MathUtil.h:58
pat::Tau::phiAtEcalEntrance
float phiAtEcalEntrance() const
return phiAtEcalEntrance
Definition: Tau.h:347
pat::Tau::leadChargedHadrCand
const reco::CandidatePtr leadChargedHadrCand() const
edm::FileInPath::location
LocationCode location() const
Where was the file found?
Definition: FileInPath.cc:191
analyzePatCleaning_cfg.inputFile
inputFile
Definition: analyzePatCleaning_cfg.py:6
taus_updatedMVAIds_cff.mvaName
mvaName
Definition: taus_updatedMVAIds_cff.py:16
pat::Tau
Analysis-level tau class.
Definition: Tau.h:55
AntiElectronIDMVA6::mvaName_wGwGSF_EC_
std::string mvaName_wGwGSF_EC_
Definition: AntiElectronIDMVA6.h:155
MagneticField::inTesla
virtual GlobalVector inTesla(const GlobalPoint &gp) const =0
Field value ad specified global point, in Tesla.
BaseParticlePropagator::propagateToEcalEntrance
bool propagateToEcalEntrance(bool first=true)
Definition: BaseParticlePropagator.cc:487
edm::Ptr::isNonnull
bool isNonnull() const
Checks for non-null.
Definition: Ptr.h:168
M_PI
#define M_PI
Definition: BXVectorInputProducer.cc:51
deltaR
double deltaR(double eta1, double eta2, double phi1, double phi2)
Definition: TreeUtility.cc:17
pat::Tau::ptLeadChargedCand
float ptLeadChargedCand() const
return pt from LeadChargedCand
Definition: Tau.h:353
PFRecoTauDiscriminationAgainstElectronDeadECAL_cfi.dR
dR
Definition: PFRecoTauDiscriminationAgainstElectronDeadECAL_cfi.py:19
RawParticle::vertex
const XYZTLorentzVector & vertex() const
the vertex fourvector
Definition: RawParticle.h:285
reco::OK
std::pair< int, edm::FunctionWithDict > OK
Definition: findMethod.cc:136
AntiElectronIDMVA6::MVAValue
double MVAValue(Float_t TauPt, Float_t TauEtaAtEcalEntrance, Float_t TauPhi, Float_t TauLeadChargedPFCandPt, Float_t TauLeadChargedPFCandEtaAtEcalEntrance, Float_t TauEmFraction, Float_t TauLeadPFChargedHadrHoP, Float_t TauLeadPFChargedHadrEoP, Float_t TauVisMassIn, Float_t TaudCrackEta, Float_t TaudCrackPhi, Float_t TauHasGsf, Int_t TauSignalPFGammaCandsIn, Int_t TauSignalPFGammaCandsOut, const std::vector< Float_t > &GammasdEtaInSigCone, const std::vector< Float_t > &GammasdPhiInSigCone, const std::vector< Float_t > &GammasPtInSigCone, const std::vector< Float_t > &GammasdEtaOutSigCone, const std::vector< Float_t > &GammasdPhiOutSigCone, const std::vector< Float_t > &GammasPtOutSigCone, Float_t ElecEta, Float_t ElecPhi, Float_t ElecEtotOverPin, Float_t ElecChi2NormGSF, Float_t ElecChi2NormKF, Float_t ElecGSFNumHits, Float_t ElecKFNumHits, Float_t ElecGSFTrackResol, Float_t ElecGSFTracklnPt, Float_t ElecPin, Float_t ElecPout, Float_t ElecEecal, Float_t ElecDeltaEta, Float_t ElecDeltaPhi, Float_t ElecMvaInSigmaEtaEta, Float_t ElecMvaInHadEnergy, Float_t ElecMvaInDeltaEta)
Definition: AntiElectronIDMVA6.cc:145
HiIsolationCommonParameters_cff.track
track
Definition: HiIsolationCommonParameters_cff.py:9
math::XYZPoint
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
part
part
Definition: HCALResponse.h:20
AntiElectronIDMVA6::Var_NoEleMatch_wGwoGSF_Endcap_
Float_t * Var_NoEleMatch_wGwoGSF_Endcap_
Definition: AntiElectronIDMVA6.h:162
b
double b
Definition: hdecay.h:120
pat::Electron
Analysis-level electron class.
Definition: Electron.h:52
pat::Electron::closestCtfTrackRef
reco::TrackRef closestCtfTrackRef() const override
override the reco::GsfElectron::closestCtfTrackRef method, to access the internal storage of the trac...
AntiElectronIDMVA6::mva_wGwGSF_BL_
const GBRForest * mva_wGwGSF_BL_
Definition: AntiElectronIDMVA6.h:169
AntiElectronIDMVA6::mva_woGwGSF_BL_
const GBRForest * mva_woGwGSF_BL_
Definition: AntiElectronIDMVA6.h:168
reco::GsfElectron::ecalEnergy
float ecalEnergy() const
Definition: GsfElectron.h:848
reco::Candidate::LorentzVector
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Candidate.h:37
pat::Tau::emFraction_MVA
float emFraction_MVA() const
return emFraction_MVA
Definition: Tau.h:355
reco::Candidate::charge
virtual int charge() const =0
electric charge
reco::GsfElectron::deltaEtaSeedClusterTrackAtCalo
float deltaEtaSeedClusterTrackAtCalo() const
Definition: GsfElectron.h:250
AntiElectronIDMVA6::mva_NoEleMatch_wGwoGSF_EC_
const GBRForest * mva_NoEleMatch_wGwoGSF_EC_
Definition: AntiElectronIDMVA6.h:171
a
double a
Definition: hdecay.h:121
edm::EventSetup::get
T get() const
Definition: EventSetup.h:63
reco::GsfElectron::superCluster
SuperClusterRef superCluster() const override
reference to a SuperCluster
Definition: GsfElectron.h:184
reco::Candidate::get
T get() const
get a component
Definition: Candidate.h:217
AntiElectronIDMVA6::mvaName_wGwGSF_BL_
std::string mvaName_wGwGSF_BL_
Definition: AntiElectronIDMVA6.h:151
AntiElectronIDMVA6::mvaName_NoEleMatch_wGwoGSF_BL_
std::string mvaName_NoEleMatch_wGwoGSF_BL_
Definition: AntiElectronIDMVA6.h:149
edm::FileInPath::fullPath
std::string fullPath() const
Definition: FileInPath.cc:197
reco::Candidate::vertex
virtual const Point & vertex() const =0
vertex position
pat::Electron::gsfTrack
reco::GsfTrackRef gsfTrack() const override
override the reco::GsfElectron::gsfTrack method, to access the internal storage of the supercluster ...
pat::Tau::signalChargedHadrCands
reco::CandidatePtrVector signalChargedHadrCands() const
reco::Candidate::px
virtual double px() const =0
x coordinate of momentum vector
pat::Tau::etaAtEcalEntrance
float etaAtEcalEntrance() const
return etaAtEcalEntrance
Definition: Tau.h:349
pat::Tau::etaAtEcalEntranceLeadChargedCand
float etaAtEcalEntranceLeadChargedCand() const
return etaAtEcalEntrance from LeadChargedCand
Definition: Tau.h:351
AntiElectronIDMVA6::dCrackEta
double dCrackEta(double eta)
Definition: AntiElectronIDMVA6.cc:1317
AntiElectronIDMVA6::mva_NoEleMatch_woGwoGSF_BL_
const GBRForest * mva_NoEleMatch_woGwoGSF_BL_
Definition: AntiElectronIDMVA6.h:166
GBRForest::GetClassifier
double GetClassifier(const float *vector) const
Definition: GBRForest.h:43
reco::LeafCandidate::phi
double phi() const final
momentum azimuthal angle
Definition: LeafCandidate.h:133
edm::ESHandle::product
T const * product() const
Definition: ESHandle.h:86
AntiElectronIDMVA6::mvaName_woGwGSF_EC_
std::string mvaName_woGwGSF_EC_
Definition: AntiElectronIDMVA6.h:154
pat::Electron::superCluster
reco::SuperClusterRef superCluster() const override
override the reco::GsfElectron::superCluster method, to access the internal storage of the superclust...
TtFullHadEvtBuilder_cfi.sumPt
sumPt
Definition: TtFullHadEvtBuilder_cfi.py:38
reco::PFTau::signalPFChargedHadrCands
const std::vector< reco::PFCandidatePtr > & signalPFChargedHadrCands() const
Charged hadrons in signal region.
Definition: PFTau.cc:80
AntiElectronIDMVA6::Var_NoEleMatch_woGwoGSF_Barrel_
Float_t * Var_NoEleMatch_woGwoGSF_Barrel_
Definition: AntiElectronIDMVA6.h:157
Generated for CMSSW Reference Manual by doxygen 1.8.11