CMS 3D CMS Logo

List of all members | Public Member Functions | Protected Attributes | Private Member Functions | Private Attributes
PFTauElecRejectionBenchmark Class Reference

#include <PFTauElecRejectionBenchmark.h>

Public Member Functions

 PFTauElecRejectionBenchmark ()
 
void process (edm::Handle< edm::HepMCProduct > mcevt, edm::Handle< reco::PFTauCollection > pfTaus, edm::Handle< reco::PFTauDiscriminator > pfTauIsoDiscr, edm::Handle< reco::PFTauDiscriminator > pfTauElecDiscr)
 
void setup (std::string Filename, std::string benchmarkLabel, double maxDeltaR, double minRecoPt, double maxRecoAbsEta, double minMCPt, double maxMCAbsEta, std::string sGenMatchObjectLabel, bool applyEcalCrackCut, DQMStore *db_store)
 
void write ()
 
virtual ~PFTauElecRejectionBenchmark ()
 

Protected Attributes

DQMStoredb_
 

Private Member Functions

bool isInEcalCrack (double eta) const
 

Private Attributes

std::vector< TLorentzVector > _GenObjects
 
bool applyEcalCrackCut_
 
std::string benchmarkLabel_
 
TFile * file_
 
TH1F * hElecMVA
 
TH1F * hElecPreID
 
TH1F * hEmfrac
 
TH1F * hEmfrac_barrel
 
TH1F * hEmfrac_endcap
 
TH1F * hEmfrac_preid0
 
TH1F * hEmfrac_preid1
 
TH2F * hEmfracvsEoP
 
TH2F * hEmfracvsEoP_preid0
 
TH2F * hEmfracvsEoP_preid1
 
TH1F * hEoverP
 
TH1F * hEoverP_barrel
 
TH1F * hEoverP_endcap
 
TH1F * hEoverP_preid0
 
TH1F * hEoverP_preid1
 
TH2F * hHoPvsEoP
 
TH2F * hHoPvsEoP_preid0
 
TH2F * hHoPvsEoP_preid1
 
TH1F * hHoverP
 
TH1F * hHoverP_barrel
 
TH1F * hHoverP_endcap
 
TH1F * hHoverP_preid0
 
TH1F * hHoverP_preid1
 
TH1F * hleadGsfTk_eta
 
TH1F * hleadGsfTk_phi
 
TH1F * hleadGsfTk_pt
 
TH1F * hleadTk_eta
 
TH1F * hleadTk_phi
 
TH1F * hleadTk_pt
 
TH1F * hpfcand_deltaEta
 
TH1F * hpfcand_deltaEta_weightE
 
TH1F * hpfcand_deltaPhiOverQ
 
TH1F * hpfcand_deltaPhiOverQ_weightE
 
TH1F * hTauElecDiscriminant
 
double maxDeltaR_
 
double maxMCAbsEta_
 
double maxRecoAbsEta_
 
double minMCPt_
 
double minRecoPt_
 
std::string outputFile_
 
std::string sGenMatchObjectLabel_
 

Detailed Description

Definition at line 34 of file PFTauElecRejectionBenchmark.h.

Constructor & Destructor Documentation

PFTauElecRejectionBenchmark::PFTauElecRejectionBenchmark ( )

Definition at line 24 of file PFTauElecRejectionBenchmark.cc.

PFTauElecRejectionBenchmark::~PFTauElecRejectionBenchmark ( )
virtual

Definition at line 26 of file PFTauElecRejectionBenchmark.cc.

References file_.

26  {
27  if(file_) file_->Close();
28 }

Member Function Documentation

bool PFTauElecRejectionBenchmark::isInEcalCrack ( double  eta) const
private

Definition at line 362 of file PFTauElecRejectionBenchmark.cc.

Referenced by process().

362  {
363  return (eta < 0.018 ||
364  (eta>0.423 && eta<0.461) ||
365  (eta>0.770 && eta<0.806) ||
366  (eta>1.127 && eta<1.163) ||
367  (eta>1.460 && eta<1.558));
368 }
void PFTauElecRejectionBenchmark::process ( edm::Handle< edm::HepMCProduct mcevt,
edm::Handle< reco::PFTauCollection pfTaus,
edm::Handle< reco::PFTauDiscriminator pfTauIsoDiscr,
edm::Handle< reco::PFTauDiscriminator pfTauElecDiscr 
)

Definition at line 217 of file PFTauElecRejectionBenchmark.cc.

References _GenObjects, funct::abs(), applyEcalCrackCut_, spr::deltaEta, allConversions_cfi::DeltaPhi, hiPixelPairStep_cff::deltaPhi, reco::PFCandidate::ecalEnergy(), reco::tau::disc::Eta(), edm::HepMCProduct::GetEvent(), hElecMVA, hEmfrac, hEmfrac_barrel, hEmfrac_endcap, hEmfrac_preid0, hEmfrac_preid1, hEmfracvsEoP, hEmfracvsEoP_preid0, hEmfracvsEoP_preid1, hEoverP, hEoverP_barrel, hEoverP_endcap, hEoverP_preid0, hEoverP_preid1, hHoPvsEoP, hHoPvsEoP_preid0, hHoPvsEoP_preid1, hHoverP, hHoverP_barrel, hHoverP_endcap, hHoverP_preid0, hHoverP_preid1, hleadTk_eta, hleadTk_phi, hleadTk_pt, hpfcand_deltaEta, hpfcand_deltaEta_weightE, hpfcand_deltaPhiOverQ, hpfcand_deltaPhiOverQ_weightE, hTauElecDiscriminant, mps_fire::i, createfilelist::int, isInEcalCrack(), edm::Ptr< T >::isNonnull(), edm::Ref< C, T, F >::isNonnull(), maxDeltaR_, maxMCAbsEta_, maxRecoAbsEta_, minMCPt_, minRecoPt_, AlCaHLTBitMon_ParallelJobs::p, reco::LeafCandidate::pdgId(), reco::PFCandidate::positionAtECALEntrance(), reco::LeafCandidate::px(), reco::LeafCandidate::py(), reco::LeafCandidate::pz(), sGenMatchObjectLabel_, metsig::tau, reco::PFCandidate::trackRef(), and z.

219  {
220 
221 
222  // Find Gen Objects to be matched with
223  HepMC::GenEvent * generated_event = new HepMC::GenEvent(*(mcevt->GetEvent()));
224  _GenObjects.clear();
225 
226  TLorentzVector taunet;
227  HepMC::GenEvent::particle_iterator p;
228  for (p = generated_event->particles_begin(); p != generated_event->particles_end(); p++) {
229  if(std::abs((*p)->pdg_id()) == 15&&(*p)->status()==2) {
230  bool lept_decay = false;
231  TLorentzVector tau((*p)->momentum().px(),(*p)->momentum().py(),(*p)->momentum().pz(),(*p)->momentum().e());
232  HepMC::GenVertex::particle_iterator z = (*p)->end_vertex()->particles_begin(HepMC::descendants);
233  for(; z != (*p)->end_vertex()->particles_end(HepMC::descendants); z++) {
234  if(std::abs((*z)->pdg_id()) == 11 || std::abs((*z)->pdg_id()) == 13) lept_decay=true;
235  if(std::abs((*z)->pdg_id()) == 16)
236  taunet.SetPxPyPzE((*z)->momentum().px(),(*z)->momentum().py(),(*z)->momentum().pz(),(*z)->momentum().e());
237 
238  }
239  if(lept_decay==false) {
240  TLorentzVector jetMom=tau-taunet;
241  if (sGenMatchObjectLabel_=="tau") _GenObjects.push_back(jetMom);
242  }
243  } else if(std::abs((*p)->pdg_id()) == 11&&(*p)->status()==1) {
244  TLorentzVector elec((*p)->momentum().px(),(*p)->momentum().py(),(*p)->momentum().pz(),(*p)->momentum().e());
245  if (sGenMatchObjectLabel_=="e") _GenObjects.push_back(elec);
246  }
247  }
248 
250 
251 
252  // Loop over all PFTaus
253  math::XYZPointF myleadTkEcalPos;
254  for (PFTauCollection::size_type iPFTau=0;iPFTau<pfTaus->size();iPFTau++) {
255  PFTauRef thePFTau(pfTaus,iPFTau);
256  if ((*pfTauIsoDiscr)[thePFTau] == 1) {
257  if ((*thePFTau).et() > minRecoPt_ && std::abs((*thePFTau).eta()) < maxRecoAbsEta_) {
258 
259  // Check if track goes to Ecal crack
260  reco::TrackRef myleadTk;
261  const reco::PFCandidatePtr& pflch = thePFTau->leadPFChargedHadrCand();
262  if(pflch.isNonnull()){
263  myleadTk=pflch->trackRef();
264  myleadTkEcalPos = pflch->positionAtECALEntrance();
265 
266  if(myleadTk.isNonnull()){
267  if (applyEcalCrackCut_ && isInEcalCrack(std::abs((double)myleadTkEcalPos.eta()))) {
268  continue; // do nothing
269  } else {
270 
271  // Match with gen object
272  for (unsigned int i = 0; i<_GenObjects.size();i++) {
273  if (_GenObjects[i].Et() >= minMCPt_ && std::abs(_GenObjects[i].Eta()) < maxMCAbsEta_ ) {
274  TLorentzVector pftau((*thePFTau).px(),(*thePFTau).py(),(*thePFTau).pz(),(*thePFTau).energy());
275  double GenDeltaR = pftau.DeltaR(_GenObjects[i]);
276  if (GenDeltaR<maxDeltaR_) {
277 
278  hleadTk_pt->Fill((float)myleadTk->pt());
279  hleadTk_eta->Fill((float)myleadTk->eta());
280  hleadTk_phi->Fill((float)myleadTk->phi());
281 
282  hEoverP->Fill((*thePFTau).ecalStripSumEOverPLead());
283  hHoverP->Fill((*thePFTau).hcal3x3OverPLead());
284  hEmfrac->Fill((*thePFTau).emFraction());
285 
286  if (std::abs(myleadTk->eta())<1.5) {
287  hEoverP_barrel->Fill((*thePFTau).ecalStripSumEOverPLead());
288  hHoverP_barrel->Fill((*thePFTau).hcal3x3OverPLead());
289  hEmfrac_barrel->Fill((*thePFTau).emFraction());
290  } else if (std::abs(myleadTk->eta())>1.5 && std::abs(myleadTk->eta())<2.5) {
291  hEoverP_endcap->Fill((*thePFTau).ecalStripSumEOverPLead());
292  hHoverP_endcap->Fill((*thePFTau).hcal3x3OverPLead());
293  hEmfrac_endcap->Fill((*thePFTau).emFraction());
294  }
295 
296  // if -999 fill in -1 bin!
297  if ((*thePFTau).electronPreIDOutput()<-1)
298  hElecMVA->Fill(-1);
299  else
300  hElecMVA->Fill((*thePFTau).electronPreIDOutput());
301 
302  hTauElecDiscriminant->Fill((*pfTauElecDiscr)[thePFTau]);
303 
304  hHoPvsEoP->Fill((*thePFTau).ecalStripSumEOverPLead(),(*thePFTau).hcal3x3OverPLead());
305  hEmfracvsEoP->Fill((*thePFTau).emFraction(),(*thePFTau).hcal3x3OverPLead());
306 
307  if ((*thePFTau).electronPreIDDecision()==1) {
308  hEoverP_preid1->Fill((*thePFTau).ecalStripSumEOverPLead());
309  hHoverP_preid1->Fill((*thePFTau).hcal3x3OverPLead());
310  hEmfrac_preid1->Fill((*thePFTau).emFraction());
311  hHoPvsEoP_preid1->Fill((*thePFTau).ecalStripSumEOverPLead(),(*thePFTau).hcal3x3OverPLead());
312  hEmfracvsEoP_preid1->Fill((*thePFTau).emFraction(),(*thePFTau).hcal3x3OverPLead());
313  } else {
314  hEoverP_preid0->Fill((*thePFTau).ecalStripSumEOverPLead());
315  hHoverP_preid0->Fill((*thePFTau).hcal3x3OverPLead());
316  hEmfrac_preid0->Fill((*thePFTau).emFraction());
317  hHoPvsEoP_preid0->Fill((*thePFTau).ecalStripSumEOverPLead(),(*thePFTau).hcal3x3OverPLead());
318  hEmfracvsEoP_preid0->Fill((*thePFTau).emFraction(),(*thePFTau).hcal3x3OverPLead());
319  }
320 
321 
322  }
323  }
324 
325  // Loop over all PFCands for cluster plots
326  std::vector<CandidatePtr> myPFCands=(*thePFTau).pfTauTagInfoRef()->PFCands();
327  for(int i=0;i<(int)myPFCands.size();i++){
328  const reco::PFCandidate* pfCand = dynamic_cast<const reco::PFCandidate*>(myPFCands[i].get());
329  if (pfCand == nullptr)
330  continue;
331 
332  math::XYZPointF candPos;
333  if (std::abs(pfCand->pdgId()) == 211 || std::abs(pfCand->pdgId()) == 11) // if charged hadron or electron
334  candPos = pfCand->positionAtECALEntrance();
335  else
336  candPos = math::XYZPointF(pfCand->px(),pfCand->py(),pfCand->pz());
337 
338  //double deltaR = ROOT::Math::VectorUtil::DeltaR(myleadTkEcalPos,candPos);
339  double deltaPhi = ROOT::Math::VectorUtil::DeltaPhi(myleadTkEcalPos,candPos);
340  double deltaEta = std::abs(myleadTkEcalPos.eta()-myPFCands[i]->eta());
341  double deltaPhiOverQ = deltaPhi/(double)myleadTk->charge();
342 
343  hpfcand_deltaEta->Fill(deltaEta);
344  hpfcand_deltaEta_weightE->Fill(deltaEta*pfCand->ecalEnergy());
345  hpfcand_deltaPhiOverQ->Fill(deltaPhiOverQ);
346  hpfcand_deltaPhiOverQ_weightE->Fill(deltaPhiOverQ*pfCand->ecalEnergy());
347 
348  }
349 
350  }
351 
352  }
353  }
354  }
355 
356  }
357  }
358  }
359 }
double ecalEnergy() const
return corrected Ecal energy
Definition: PFCandidate.h:222
int pdgId() const final
PDG identifier.
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:251
double px() const final
x coordinate of momentum vector
ROOT::Math::PositionVector3D< ROOT::Math::Cartesian3D< float > > XYZPointF
point in space with cartesian internal representation
Definition: Point3D.h:10
uint16_t size_type
const math::XYZPointF & positionAtECALEntrance() const
Definition: PFCandidate.h:368
static const double deltaEta
Definition: CaloConstants.h:8
reco::TrackRef trackRef() const
Definition: PFCandidate.cc:442
std::vector< TLorentzVector > _GenObjects
double pz() const final
z coordinate of momentum vector
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool isNonnull() const
Checks for non-null.
Definition: Ptr.h:168
const HepMC::GenEvent * GetEvent() const
Definition: HepMCProduct.h:38
double py() const final
y coordinate of momentum vector
Particle reconstructed by the particle flow algorithm.
Definition: PFCandidate.h:40
void PFTauElecRejectionBenchmark::setup ( std::string  Filename,
std::string  benchmarkLabel,
double  maxDeltaR,
double  minRecoPt,
double  maxRecoAbsEta,
double  minMCPt,
double  maxMCAbsEta,
std::string  sGenMatchObjectLabel,
bool  applyEcalCrackCut,
DQMStore db_store 
)

Definition at line 48 of file PFTauElecRejectionBenchmark.cc.

References applyEcalCrackCut_, benchmarkLabel_, BOOK1D, BOOK2D, gather_cfg::cout, db_, file_, electrons_cff::maxDeltaR, maxDeltaR_, pfTauBenchmarkElecRejection_cfi::maxMCAbsEta, maxMCAbsEta_, pfTauBenchmarkElecRejection_cfi::maxRecoAbsEta, maxRecoAbsEta_, pfTauBenchmarkElecRejection_cfi::minMCPt, minMCPt_, pfTauBenchmarkElecRejection_cfi::minRecoPt, minRecoPt_, outputFile_, callgraph::path, SETAXES, and sGenMatchObjectLabel_.

58  {
60  benchmarkLabel_ = benchmarkLabel;
61  outputFile_=Filename;
64  minMCPt_ = minMCPt;
66  sGenMatchObjectLabel_ = sGenMatchObjectLabel;
67  applyEcalCrackCut_= applyEcalCrackCut;
68 
69  file_ = nullptr;
70  db_ = db_store;
71 
72  // print parameters
73  cout<< "PFTauElecRejectionBenchmark Setup parameters =============================================="<<endl;
74  cout << "Filename to write histograms " << Filename<<endl;
75  cout << "Benchmark label name " << benchmarkLabel_<<endl;
76  cout << "maxDeltaRMax " << maxDeltaR << endl;
77  cout << "minRecoPt " << minRecoPt_ << endl;
78  cout << "maxRecoAbsEta " << maxRecoAbsEta_ << endl;
79  cout << "minMCPt " << minMCPt_ << endl;
80  cout << "maxMCAbsEta " << maxMCAbsEta_ << endl;
81  cout << "applyEcalCrackCut " << applyEcalCrackCut_ << endl;
82  cout << "sGenMatchObjectLabel " << sGenMatchObjectLabel_ << endl;
83 
84  // Book histogram
85 
86  // Establish DQM Store
87  string path = "PFTask/Benchmarks/"+ benchmarkLabel_ + "/";
88  path += "Gen";
89  if (db_) {
90  db_->setCurrentFolder(path);
91  }
92  else {
93  file_ = new TFile(outputFile_.c_str(), "recreate");
94  cout << "Info: DQM is not available to provide data storage service. Using TFile to save histograms. "<<endl;
95  }
96 
97  // E/p
98  BOOK1D(EoverP,"E/p",100, 0., 4.);
99  SETAXES(EoverP, "E/p", "Entries");
100 
101  BOOK1D(EoverP_barrel,"E/p barrel",100, 0., 4.);
102  SETAXES(EoverP_barrel, "E/p barrel", "Entries");
103 
104  BOOK1D(EoverP_endcap,"E/p endcap",100, 0., 4.);
105  SETAXES(EoverP_endcap, "E/p endcap", "Entries");
106 
107  BOOK1D(EoverP_preid0,"E/p (preid=0)",100, 0., 4.);
108  SETAXES(EoverP_preid0, "E/p", "Entries");
109 
110  BOOK1D(EoverP_preid1,"E/p (preid=1)",100, 0., 4.);
111  SETAXES(EoverP_preid1, "E/p", "Entries");
112 
113  // H/p
114  BOOK1D(HoverP,"H/p",100, 0., 2.);
115  SETAXES(HoverP, "H/p", "Entries");
116 
117  BOOK1D(HoverP_barrel,"H/p barrel",100, 0., 2.);
118  SETAXES(HoverP_barrel, "H/p barrel", "Entries");
119 
120  BOOK1D(HoverP_endcap,"H/p endcap",100, 0., 2.);
121  SETAXES(HoverP_endcap, "H/p endcap", "Entries");
122 
123  BOOK1D(HoverP_preid0,"H/p (preid=0)",100, 0., 2.);
124  SETAXES(HoverP_preid0, "H/p", "Entries");
125 
126  BOOK1D(HoverP_preid1,"H/p (preid=1)",100, 0., 2.);
127  SETAXES(HoverP_preid1, "H/p", "Entries");
128 
129  // emfrac
130  BOOK1D(Emfrac,"EM fraction",100, 0., 1.01);
131  SETAXES(Emfrac, "em fraction", "Entries");
132 
133  BOOK1D(Emfrac_barrel,"EM fraction barrel",100, 0., 1.01);
134  SETAXES(Emfrac_barrel, "em fraction barrel", "Entries");
135 
136  BOOK1D(Emfrac_endcap,"EM fraction endcap",100, 0., 1.01);
137  SETAXES(Emfrac_endcap, "em fraction endcap", "Entries");
138 
139  BOOK1D(Emfrac_preid0,"EM fraction (preid=0)",100, 0., 1.01);
140  SETAXES(Emfrac_preid0, "em fraction", "Entries");
141 
142  BOOK1D(Emfrac_preid1,"EM fraction (preid=1)",100, 0., 1.01);
143  SETAXES(Emfrac_preid1, "em fraction", "Entries");
144 
145  // PreID
146  BOOK1D(ElecPreID,"PFElectron PreID decision",6, 0., 1.01);
147  SETAXES(ElecPreID, "PFElectron PreID decision", "Entries");
148 
149  // MVA
150  BOOK1D(ElecMVA,"PFElectron MVA",100, -1.01, 1.01);
151  SETAXES(ElecMVA, "PFElectron MVA", "Entries");
152 
153  // Discriminant
154  BOOK1D(TauElecDiscriminant,"PFTau-Electron Discriminant",6, 0., 1.01);
155  SETAXES(TauElecDiscriminant, "PFTau-Electron Discriminant", "Entries");
156 
157 
158  // PFCand clusters
159  BOOK1D(pfcand_deltaEta,"PFCand cluster dEta",100, 0., 0.8);
160  SETAXES(pfcand_deltaEta, "PFCand cluster #Delta(#eta)", "Entries");
161 
162  BOOK1D(pfcand_deltaEta_weightE,"PFCand cluster dEta, energy weighted",100, 0., 0.8);
163  SETAXES(pfcand_deltaEta_weightE, "PFCand cluster #Delta(#eta)", "Entries");
164 
165  BOOK1D(pfcand_deltaPhiOverQ,"PFCand cluster dPhi/q",100, -0.8, 0.8);
166  SETAXES(pfcand_deltaPhiOverQ, "PFCand cluster #Delta(#phi)/q", "Entries");
167 
168  BOOK1D(pfcand_deltaPhiOverQ_weightE,"PFCand cluster dEta/q, energy weighted",100, -0.8, 0.8);
169  SETAXES(pfcand_deltaPhiOverQ_weightE, "PFCand cluster #Delta(#phi)/q", "Entries");
170 
171 
172  // Leading KF track
173  BOOK1D(leadTk_pt,"leading KF track pt",100, 0., 80.);
174  SETAXES(leadTk_pt, "leading KF track p_{T} (GeV)", "Entries");
175 
176  BOOK1D(leadTk_eta,"leading KF track eta",100, -4., 4.);
177  SETAXES(leadTk_eta, "leading KF track #eta", "Entries");
178 
179  BOOK1D(leadTk_phi,"leading KF track phi",100, -3.2, 3.2);
180  SETAXES(leadTk_phi, "leading KF track #phi", "Entries");
181 
182  // Leading Gsf track
183  BOOK1D(leadGsfTk_pt,"leading Gsf track pt",100, 0., 80.);
184  SETAXES(leadGsfTk_pt, "leading Gsf track p_{T} (GeV)", "Entries");
185 
186  BOOK1D(leadGsfTk_eta,"leading Gsf track eta",100, -4., 4.);
187  SETAXES(leadGsfTk_eta, "leading Gsf track #eta", "Entries");
188 
189  BOOK1D(leadGsfTk_phi,"leading Gsf track phi",100, -3.2, 3.2);
190  SETAXES(leadGsfTk_phi, "leading Gsf track #phi", "Entries");
191 
192 
193  // H/p vs E/p
194  BOOK2D(HoPvsEoP,"H/p vs. E/p",100, 0., 2., 100, 0., 2.);
195  SETAXES(HoPvsEoP, "E/p", "H/p");
196 
197  BOOK2D(HoPvsEoP_preid0,"H/p vs. E/p (preid=0)",100, 0., 2., 100, 0., 2.);
198  SETAXES(HoPvsEoP_preid0, "E/p", "H/p");
199 
200  BOOK2D(HoPvsEoP_preid1,"H/p vs. E/p (preid=0)",100, 0., 2., 100, 0., 2.);
201  SETAXES(HoPvsEoP_preid1, "E/p", "H/p");
202 
203  // em fraction vs E/p
204  BOOK2D(EmfracvsEoP,"emfrac vs. E/p",100, 0., 2., 100, 0., 1.01);
205  SETAXES(EmfracvsEoP, "E/p", "em fraction");
206 
207  BOOK2D(EmfracvsEoP_preid0,"emfrac vs. E/p (preid=0)",100, 0., 2., 100, 0., 1.01);
208  SETAXES(EmfracvsEoP_preid0, "E/p", "em fraction");
209 
210  BOOK2D(EmfracvsEoP_preid1,"emfrac vs. E/p (preid=0)",100, 0., 2., 100, 0., 1.01);
211  SETAXES(EmfracvsEoP_preid1, "E/p", "em fraction");
212 
213 
214 }
#define BOOK1D(name, title, nbinsx, lowx, highx)
#define BOOK2D(name, title, nbinsx, lowx, highx, nbinsy, lowy, highy)
#define SETAXES(name, xtitle, ytitle)
void PFTauElecRejectionBenchmark::write ( )

Definition at line 30 of file PFTauElecRejectionBenchmark.cc.

References gather_cfg::cout, db_, file_, and outputFile_.

30  {
31  // Store the DAQ Histograms
32  if (!outputFile_.empty()) {
33  if (db_)
34  db_->save(outputFile_);
35  // use bare Root if no DQM (FWLite applications)
36  else if (file_) {
37  file_->Write(outputFile_.c_str());
38  cout << "Benchmark output written to file " << outputFile_.c_str() << endl;
39  file_->Close();
40  }
41  }
42  else
43  cout << "No output file specified ("<<outputFile_<<"). Results will not be saved!" << endl;
44 
45 }

Member Data Documentation

std::vector<TLorentzVector> PFTauElecRejectionBenchmark::_GenObjects
private

Definition at line 120 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

bool PFTauElecRejectionBenchmark::applyEcalCrackCut_
private

Definition at line 70 of file PFTauElecRejectionBenchmark.h.

Referenced by process(), and setup().

std::string PFTauElecRejectionBenchmark::benchmarkLabel_
private

Definition at line 63 of file PFTauElecRejectionBenchmark.h.

Referenced by setup().

DQMStore* PFTauElecRejectionBenchmark::db_
protected

Definition at line 124 of file PFTauElecRejectionBenchmark.h.

Referenced by setup(), and write().

TFile* PFTauElecRejectionBenchmark::file_
private

Definition at line 61 of file PFTauElecRejectionBenchmark.h.

Referenced by setup(), write(), and ~PFTauElecRejectionBenchmark().

TH1F* PFTauElecRejectionBenchmark::hElecMVA
private

Definition at line 94 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hElecPreID
private

Definition at line 93 of file PFTauElecRejectionBenchmark.h.

TH1F* PFTauElecRejectionBenchmark::hEmfrac
private

Definition at line 75 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hEmfrac_barrel
private

Definition at line 79 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hEmfrac_endcap
private

Definition at line 83 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hEmfrac_preid0
private

Definition at line 87 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hEmfrac_preid1
private

Definition at line 91 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH2F* PFTauElecRejectionBenchmark::hEmfracvsEoP
private

Definition at line 101 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH2F* PFTauElecRejectionBenchmark::hEmfracvsEoP_preid0
private

Definition at line 102 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH2F* PFTauElecRejectionBenchmark::hEmfracvsEoP_preid1
private

Definition at line 103 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hEoverP
private

Definition at line 73 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hEoverP_barrel
private

Definition at line 77 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hEoverP_endcap
private

Definition at line 81 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hEoverP_preid0
private

Definition at line 85 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hEoverP_preid1
private

Definition at line 89 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH2F* PFTauElecRejectionBenchmark::hHoPvsEoP
private

Definition at line 97 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH2F* PFTauElecRejectionBenchmark::hHoPvsEoP_preid0
private

Definition at line 98 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH2F* PFTauElecRejectionBenchmark::hHoPvsEoP_preid1
private

Definition at line 99 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hHoverP
private

Definition at line 74 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hHoverP_barrel
private

Definition at line 78 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hHoverP_endcap
private

Definition at line 82 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hHoverP_preid0
private

Definition at line 86 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hHoverP_preid1
private

Definition at line 90 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hleadGsfTk_eta
private

Definition at line 116 of file PFTauElecRejectionBenchmark.h.

TH1F* PFTauElecRejectionBenchmark::hleadGsfTk_phi
private

Definition at line 117 of file PFTauElecRejectionBenchmark.h.

TH1F* PFTauElecRejectionBenchmark::hleadGsfTk_pt
private

Definition at line 115 of file PFTauElecRejectionBenchmark.h.

TH1F* PFTauElecRejectionBenchmark::hleadTk_eta
private

Definition at line 111 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hleadTk_phi
private

Definition at line 112 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hleadTk_pt
private

Definition at line 110 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hpfcand_deltaEta
private

Definition at line 105 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hpfcand_deltaEta_weightE
private

Definition at line 106 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hpfcand_deltaPhiOverQ
private

Definition at line 107 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hpfcand_deltaPhiOverQ_weightE
private

Definition at line 108 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

TH1F* PFTauElecRejectionBenchmark::hTauElecDiscriminant
private

Definition at line 95 of file PFTauElecRejectionBenchmark.h.

Referenced by process().

double PFTauElecRejectionBenchmark::maxDeltaR_
private

Definition at line 64 of file PFTauElecRejectionBenchmark.h.

Referenced by process(), and setup().

double PFTauElecRejectionBenchmark::maxMCAbsEta_
private

Definition at line 66 of file PFTauElecRejectionBenchmark.h.

Referenced by process(), and setup().

double PFTauElecRejectionBenchmark::maxRecoAbsEta_
private

Definition at line 68 of file PFTauElecRejectionBenchmark.h.

Referenced by process(), and setup().

double PFTauElecRejectionBenchmark::minMCPt_
private

Definition at line 65 of file PFTauElecRejectionBenchmark.h.

Referenced by process(), and setup().

double PFTauElecRejectionBenchmark::minRecoPt_
private

Definition at line 67 of file PFTauElecRejectionBenchmark.h.

Referenced by process(), and setup().

std::string PFTauElecRejectionBenchmark::outputFile_
private

Definition at line 62 of file PFTauElecRejectionBenchmark.h.

Referenced by setup(), and write().

std::string PFTauElecRejectionBenchmark::sGenMatchObjectLabel_
private

Definition at line 69 of file PFTauElecRejectionBenchmark.h.

Referenced by process(), and setup().