CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Types | Public Member Functions | Private Member Functions | Private Attributes
BTagDifferentialPlot Class Reference

#include <BTagDifferentialPlot.h>

Public Types

enum  ConstVarType { constPT, constETA }
 

Public Member Functions

void addBinPlotter (JetTagPlotter *aPlotter)
 
 BTagDifferentialPlot (const double &bEff, const ConstVarType &constVariable, const std::string &tagName)
 
void epsPlot (const std::string &name)
 
TH1F * getDifferentialHistoB_b ()
 
TH1F * getDifferentialHistoB_c ()
 
TH1F * getDifferentialHistoB_d ()
 
TH1F * getDifferentialHistoB_dus ()
 
TH1F * getDifferentialHistoB_dusg ()
 
TH1F * getDifferentialHistoB_g ()
 
TH1F * getDifferentialHistoB_ni ()
 
TH1F * getDifferentialHistoB_s ()
 
TH1F * getDifferentialHistoB_u ()
 
void plot (TCanvas &theCanvas)
 
void plot (const std::string &name, const std::string &ext)
 
void process ()
 
void psPlot (const std::string &name)
 
 ~BTagDifferentialPlot ()
 

Private Member Functions

void bookHisto ()
 
void fillHisto ()
 
std::pair< double, double > getMistag (const double &fixedBEfficiency, TH1F *effPurHist)
 
void setVariableName ()
 

Private Attributes

std::string commonName
 
ConstVarType constVar
 
std::string constVariableName
 
std::pair< double, double > constVariableValue
 
std::string diffVariableName
 
double fixedBEfficiency
 
bool noProcessing
 
bool processed
 
std::vector< JetTagPlotter * > theBinPlotters
 
MonitorElementtheDifferentialHistoB_b
 
MonitorElementtheDifferentialHistoB_c
 
MonitorElementtheDifferentialHistoB_d
 
MonitorElementtheDifferentialHistoB_dus
 
MonitorElementtheDifferentialHistoB_dusg
 
MonitorElementtheDifferentialHistoB_g
 
MonitorElementtheDifferentialHistoB_ni
 
MonitorElementtheDifferentialHistoB_s
 
MonitorElementtheDifferentialHistoB_u
 

Detailed Description

Definition at line 15 of file BTagDifferentialPlot.h.

Member Enumeration Documentation

Enumerator
constPT 
constETA 

Definition at line 19 of file BTagDifferentialPlot.h.

Constructor & Destructor Documentation

BTagDifferentialPlot::BTagDifferentialPlot ( const double &  bEff,
const ConstVarType constVariable,
const std::string &  tagName 
)

Definition at line 20 of file BTagDifferentialPlot.cc.

21  :
22  fixedBEfficiency ( bEff ) ,
23  noProcessing ( false ) , processed(false), constVar(constVariable),
24  constVariableName ( "" ) , diffVariableName ( "" ) ,
25  constVariableValue ( 999.9 , 999.9 ) , commonName( "MisidForBEff_" + tagName+"_") ,
MonitorElement * theDifferentialHistoB_b
MonitorElement * theDifferentialHistoB_g
MonitorElement * theDifferentialHistoB_s
MonitorElement * theDifferentialHistoB_dusg
MonitorElement * theDifferentialHistoB_dus
std::pair< double, double > constVariableValue
MonitorElement * theDifferentialHistoB_u
MonitorElement * theDifferentialHistoB_d
MonitorElement * theDifferentialHistoB_c
MonitorElement * theDifferentialHistoB_ni
BTagDifferentialPlot::~BTagDifferentialPlot ( )

Definition at line 37 of file BTagDifferentialPlot.cc.

37  {
38 }

Member Function Documentation

void BTagDifferentialPlot::addBinPlotter ( JetTagPlotter aPlotter)
inline

Definition at line 26 of file BTagDifferentialPlot.h.

References theBinPlotters.

26 { theBinPlotters.push_back ( aPlotter ) ; }
std::vector< JetTagPlotter * > theBinPlotters
void BTagDifferentialPlot::bookHisto ( )
private

Definition at line 206 of file BTagDifferentialPlot.cc.

References HistoProviderDQM::book1D(), commonName, constVariableName, constVariableValue, diffVariableName, fixedBEfficiency, EtaPtBin::getEtaActive(), EtaPtBin::getEtaMax(), EtaPtBin::getEtaMin(), EtaPtBin::getPtActive(), EtaPtBin::getPtMax(), EtaPtBin::getPtMin(), diffTwoXMLs::label, python.multivaluedict::remove(), python.rootplot.root2matplotlib::replace(), theBinPlotters, theDifferentialHistoB_b, theDifferentialHistoB_c, theDifferentialHistoB_d, theDifferentialHistoB_dus, theDifferentialHistoB_dusg, theDifferentialHistoB_g, theDifferentialHistoB_ni, theDifferentialHistoB_s, and theDifferentialHistoB_u.

Referenced by process().

206  {
207 
208  // vector with ranges
209  vector<float> variableRanges ;
210 
211  for ( vector<JetTagPlotter *>::const_iterator iP = theBinPlotters.begin() ;
212  iP != theBinPlotters.end() ; ++iP ) {
213  const EtaPtBin & currentBin = (*iP)->etaPtBin() ;
214  if ( diffVariableName == "eta" ) {
215  // only active bins in the variable on x-axis
216  if ( currentBin.getEtaActive() ) {
217  variableRanges.push_back ( currentBin.getEtaMin() ) ;
218  // also max if last one
219  if ( iP == --theBinPlotters.end() ) variableRanges.push_back ( currentBin.getEtaMax() ) ;
220  }
221  }
222  if ( diffVariableName == "pt" ) {
223  // only active bins in the variable on x-axis
224  if ( currentBin.getPtActive() ) {
225  variableRanges.push_back ( currentBin.getPtMin() ) ;
226  // also max if last one
227  if ( iP == --theBinPlotters.end() ) variableRanges.push_back ( currentBin.getPtMax() ) ;
228  }
229  }
230  }
231 
232  // to book histo with variable binning -> put into array
233  int nBins = variableRanges.size() - 1 ;
234  float * binArray = &variableRanges[0];
235  //float * binArray = new float [nBins+1] ;
236 
237  //for ( int i = 0 ; i < nBins + 1 ; i++ ) {
238  // binArray[i] = variableRanges[i] ;
239  //}
240 
241 
242  // part of the name common to all flavours
243  std::stringstream stream("");
244  stream << fixedBEfficiency << "_Const_" << constVariableName << "_" << constVariableValue.first << "-" ;
245  stream << constVariableValue.second << "_" << "_Vs_" << diffVariableName ;
246  commonName += stream.str();
247  std::remove(commonName.begin(), commonName.end(), ' ') ;
248  std::replace(commonName.begin(), commonName.end(), '.' , 'v' ) ;
249 
250  std::string label(commonName);
251  HistoProviderDQM prov ("Btag",label);
252 
253  theDifferentialHistoB_d = (prov.book1D ( "D_" + commonName , "D_" + commonName , nBins , binArray )) ;
254  theDifferentialHistoB_u = (prov.book1D ( "U_" + commonName , "U_" + commonName , nBins , binArray )) ;
255  theDifferentialHistoB_s = (prov.book1D ( "S_" + commonName , "S_" + commonName , nBins , binArray )) ;
256  theDifferentialHistoB_c = (prov.book1D ( "C_" + commonName , "C_" + commonName , nBins , binArray )) ;
257  theDifferentialHistoB_b = (prov.book1D ( "B_" + commonName , "B_" + commonName , nBins , binArray )) ;
258  theDifferentialHistoB_g = (prov.book1D ( "G_" + commonName , "G_" + commonName , nBins , binArray )) ;
259  theDifferentialHistoB_ni = (prov.book1D ( "NI_" + commonName , "NI_" + commonName , nBins , binArray )) ;
260  theDifferentialHistoB_dus = (prov.book1D ( "DUS_" + commonName , "DUS_" + commonName , nBins , binArray )) ;
261  theDifferentialHistoB_dusg = (prov.book1D ( "DUSG_" + commonName , "DUSG_" + commonName , nBins , binArray )) ;
262 }
MonitorElement * theDifferentialHistoB_b
MonitorElement * theDifferentialHistoB_g
MonitorElement * theDifferentialHistoB_s
MonitorElement * theDifferentialHistoB_dusg
double getEtaMax() const
Definition: EtaPtBin.h:37
double getEtaMin() const
Definition: EtaPtBin.h:36
double getPtMin() const
Definition: EtaPtBin.h:40
bool getEtaActive() const
Get rapidity/pt ranges and check whether rapidity/pt cuts are active.
Definition: EtaPtBin.h:35
MonitorElement * theDifferentialHistoB_dus
bool getPtActive() const
Definition: EtaPtBin.h:39
std::pair< double, double > constVariableValue
double getPtMax() const
Definition: EtaPtBin.h:41
std::vector< JetTagPlotter * > theBinPlotters
MonitorElement * theDifferentialHistoB_u
MonitorElement * theDifferentialHistoB_d
MonitorElement * theDifferentialHistoB_c
MonitorElement * theDifferentialHistoB_ni
void BTagDifferentialPlot::epsPlot ( const std::string &  name)

Definition at line 155 of file BTagDifferentialPlot.cc.

References plot().

156 {
157  plot(name, ".eps");
158 }
void plot(TCanvas &theCanvas)
void BTagDifferentialPlot::fillHisto ( )
private

Definition at line 265 of file BTagDifferentialPlot.cc.

References diffVariableName, edm::hlt::Exception, fixedBEfficiency, EffPurFromHistos::getEffFlavVsBEff_b(), EffPurFromHistos::getEffFlavVsBEff_c(), EffPurFromHistos::getEffFlavVsBEff_d(), EffPurFromHistos::getEffFlavVsBEff_dus(), EffPurFromHistos::getEffFlavVsBEff_dusg(), EffPurFromHistos::getEffFlavVsBEff_g(), EffPurFromHistos::getEffFlavVsBEff_ni(), EffPurFromHistos::getEffFlavVsBEff_s(), EffPurFromHistos::getEffFlavVsBEff_u(), EtaPtBin::getEtaActive(), EtaPtBin::getEtaMax(), EtaPtBin::getEtaMin(), getMistag(), EtaPtBin::getPtActive(), EtaPtBin::getPtMax(), EtaPtBin::getPtMin(), MonitorElement::getTH1F(), theBinPlotters, theDifferentialHistoB_b, theDifferentialHistoB_c, theDifferentialHistoB_d, theDifferentialHistoB_dus, theDifferentialHistoB_dusg, theDifferentialHistoB_g, theDifferentialHistoB_ni, theDifferentialHistoB_s, theDifferentialHistoB_u, and funct::true.

Referenced by process().

265  {
266  // loop over bins and find corresponding misid. in the MisIdVs..... histo
267  for ( vector<JetTagPlotter *>::const_iterator iP = theBinPlotters.begin() ;
268  iP != theBinPlotters.end() ; ++iP ) {
269  const EtaPtBin & currentBin = (*iP)->etaPtBin() ;
270  EffPurFromHistos * currentEffPurFromHistos = (*iP)->getEffPurFromHistos() ;
271  //
272  bool isActive = true ;
273  double valueXAxis = -999.99 ;
274  // find right bin based on middle of the interval
275  if ( diffVariableName == "eta" ) {
276  isActive = currentBin.getEtaActive() ;
277  valueXAxis = 0.5 * ( currentBin.getEtaMin() + currentBin.getEtaMax() ) ;
278  } else if ( diffVariableName == "pt" ) {
279  isActive = currentBin.getPtActive() ;
280  valueXAxis = 0.5 * ( currentBin.getPtMin() + currentBin.getPtMax() ) ;
281  } else {
282  throw cms::Exception("Configuration")
283  << "====>>>> BTagDifferentialPlot::fillHisto() : illegal diffVariableName = " << diffVariableName << endl;
284  }
285 
286  // for the moment: ignore inactive bins
287  // (maybe later: if a Bin is inactive -> set value to fill well below left edge of histogram to have it in the underflow)
288 
289  if ( !isActive ) continue ;
290 
291  // to have less lines of code ....
292  vector< pair<TH1F*,TH1F*> > effPurDifferentialPairs ;
293 
294  // all flavours (b is a good cross check! must be constant and = fixed b-efficiency)
295  // get histo; find the bin of the fixed b-efficiency in the histo and get misid; fill
296 
297 
298  effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_d() , theDifferentialHistoB_d ->getTH1F() ) ) ;
299  effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_u() , theDifferentialHistoB_u ->getTH1F() ) ) ;
300  effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_s() , theDifferentialHistoB_s ->getTH1F() ) ) ;
301  effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_c() , theDifferentialHistoB_c ->getTH1F() ) ) ;
302  effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_b() , theDifferentialHistoB_b ->getTH1F() ) ) ;
303  effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_g() , theDifferentialHistoB_g ->getTH1F() ) ) ;
304  effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_ni() , theDifferentialHistoB_ni ->getTH1F() ) ) ;
305  effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_dus() , theDifferentialHistoB_dus->getTH1F() ) ) ;
306  effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_dusg() , theDifferentialHistoB_dusg->getTH1F() ) ) ;
307 
308  for ( vector< pair<TH1F*,TH1F*> >::const_iterator itP = effPurDifferentialPairs.begin() ;
309  itP != effPurDifferentialPairs.end() ; ++itP ) {
310  TH1F * effPurHist = itP->first ;
311  TH1F * diffHist = itP->second ;
312  pair<double, double> mistag = getMistag(fixedBEfficiency, effPurHist);
313  int iBinSet = diffHist->FindBin(valueXAxis) ;
314  diffHist->SetBinContent(iBinSet, mistag.first);
315  diffHist->SetBinError(iBinSet, mistag.second);
316  }
317  }
318 
319 }
MonitorElement * theDifferentialHistoB_b
MonitorElement * theDifferentialHistoB_g
TH1F * getEffFlavVsBEff_g()
MonitorElement * theDifferentialHistoB_s
TH1F * getEffFlavVsBEff_dusg()
TH1F * getEffFlavVsBEff_c()
MonitorElement * theDifferentialHistoB_dusg
double getEtaMax() const
Definition: EtaPtBin.h:37
double getEtaMin() const
Definition: EtaPtBin.h:36
double getPtMin() const
Definition: EtaPtBin.h:40
bool getEtaActive() const
Get rapidity/pt ranges and check whether rapidity/pt cuts are active.
Definition: EtaPtBin.h:35
MonitorElement * theDifferentialHistoB_dus
bool getPtActive() const
Definition: EtaPtBin.h:39
TH1F * getEffFlavVsBEff_u()
std::pair< double, double > getMistag(const double &fixedBEfficiency, TH1F *effPurHist)
double getPtMax() const
Definition: EtaPtBin.h:41
std::vector< JetTagPlotter * > theBinPlotters
TH1F * getEffFlavVsBEff_dus()
TH1F * getEffFlavVsBEff_d()
MonitorElement * theDifferentialHistoB_u
TH1F * getTH1F(void) const
TH1F * getEffFlavVsBEff_ni()
MonitorElement * theDifferentialHistoB_d
TH1F * getEffFlavVsBEff_b()
TH1F * getEffFlavVsBEff_s()
MonitorElement * theDifferentialHistoB_c
MonitorElement * theDifferentialHistoB_ni
TH1F* BTagDifferentialPlot::getDifferentialHistoB_b ( )
inline

Definition at line 47 of file BTagDifferentialPlot.h.

References MonitorElement::getTH1F(), and theDifferentialHistoB_b.

47 { return theDifferentialHistoB_b ->getTH1F() ; }
MonitorElement * theDifferentialHistoB_b
TH1F * getTH1F(void) const
TH1F* BTagDifferentialPlot::getDifferentialHistoB_c ( )
inline

Definition at line 46 of file BTagDifferentialPlot.h.

References MonitorElement::getTH1F(), and theDifferentialHistoB_c.

46 { return theDifferentialHistoB_c ->getTH1F() ; }
TH1F * getTH1F(void) const
MonitorElement * theDifferentialHistoB_c
TH1F* BTagDifferentialPlot::getDifferentialHistoB_d ( )
inline

Definition at line 43 of file BTagDifferentialPlot.h.

References MonitorElement::getTH1F(), and theDifferentialHistoB_d.

43 { return theDifferentialHistoB_d ->getTH1F() ; }
TH1F * getTH1F(void) const
MonitorElement * theDifferentialHistoB_d
TH1F* BTagDifferentialPlot::getDifferentialHistoB_dus ( )
inline

Definition at line 50 of file BTagDifferentialPlot.h.

References MonitorElement::getTH1F(), and theDifferentialHistoB_dus.

50 { return theDifferentialHistoB_dus->getTH1F() ; }
MonitorElement * theDifferentialHistoB_dus
TH1F * getTH1F(void) const
TH1F* BTagDifferentialPlot::getDifferentialHistoB_dusg ( )
inline

Definition at line 51 of file BTagDifferentialPlot.h.

References MonitorElement::getTH1F(), and theDifferentialHistoB_dusg.

51 { return theDifferentialHistoB_dusg->getTH1F() ; }
MonitorElement * theDifferentialHistoB_dusg
TH1F * getTH1F(void) const
TH1F* BTagDifferentialPlot::getDifferentialHistoB_g ( )
inline

Definition at line 48 of file BTagDifferentialPlot.h.

References MonitorElement::getTH1F(), and theDifferentialHistoB_g.

48 { return theDifferentialHistoB_g ->getTH1F() ; }
MonitorElement * theDifferentialHistoB_g
TH1F * getTH1F(void) const
TH1F* BTagDifferentialPlot::getDifferentialHistoB_ni ( )
inline

Definition at line 49 of file BTagDifferentialPlot.h.

References MonitorElement::getTH1F(), and theDifferentialHistoB_ni.

49 { return theDifferentialHistoB_ni->getTH1F() ; }
TH1F * getTH1F(void) const
MonitorElement * theDifferentialHistoB_ni
TH1F* BTagDifferentialPlot::getDifferentialHistoB_s ( )
inline

Definition at line 45 of file BTagDifferentialPlot.h.

References MonitorElement::getTH1F(), and theDifferentialHistoB_s.

45 { return theDifferentialHistoB_s ->getTH1F() ; }
MonitorElement * theDifferentialHistoB_s
TH1F * getTH1F(void) const
TH1F* BTagDifferentialPlot::getDifferentialHistoB_u ( )
inline

Definition at line 44 of file BTagDifferentialPlot.h.

References MonitorElement::getTH1F(), and theDifferentialHistoB_u.

44 { return theDifferentialHistoB_u ->getTH1F() ; }
MonitorElement * theDifferentialHistoB_u
TH1F * getTH1F(void) const
pair< double, double > BTagDifferentialPlot::getMistag ( const double &  fixedBEfficiency,
TH1F *  effPurHist 
)
private

Definition at line 322 of file BTagDifferentialPlot.cc.

References i.

Referenced by fillHisto().

323 {
324  int iBinGet = effPurHist->FindBin ( fixedBEfficiency ) ;
325  double effForBEff = effPurHist->GetBinContent ( iBinGet ) ;
326  double effForBEffErr = effPurHist->GetBinError ( iBinGet ) ;
327 
328  if (effForBEff==0. && effForBEffErr==0.) {
329  // The bin was empty. Could be that it was not filled, as the scan-plot
330  // did not have an entry at the requested value, or that the curve
331  // is above or below.
332  // Fit a plynomial, and evaluate the mistag at the requested value.
333  int fitStatus;
334  try {
335  fitStatus = effPurHist->Fit("pol4", "q");
336  }catch (...){
337  return pair<double, double>(effForBEff, effForBEffErr);
338 }
339  if (fitStatus != 0) {
340  edm::LogWarning("BTagDifferentialPlot")<<"Fit failed to hisogram " << effPurHist->GetTitle() << " , perhaps because too few entries = " << effPurHist->GetEntries() <<". This bin will be missing in plots at fixed b efficiency.";
341  // } else {
342  // edm::LogInfo("BTagDifferentialPlot")<<"Fit OK to hisogram " << effPurHist->GetTitle() << " entries = " << effPurHist->GetEntries();
343  return pair<double, double>(effForBEff, effForBEffErr);
344  }
345  TF1 *myfunc = effPurHist->GetFunction("pol4");
346  effForBEff = myfunc->Eval(fixedBEfficiency);
347  effPurHist->RecursiveRemove(myfunc);
348  //Error: first non-empty bin on the right and take its error
349  for (int i = iBinGet+1; i< effPurHist->GetNbinsX(); ++i) {
350  if (effPurHist->GetBinContent(i)!=0) {
351  effForBEffErr = effPurHist->GetBinError(i);
352  break;
353  }
354  }
355  }
356 
357  return pair<double, double>(effForBEff, effForBEffErr);
358 }
int i
Definition: DBlmapReader.cc:9
void BTagDifferentialPlot::plot ( TCanvas &  theCanvas)

Definition at line 45 of file BTagDifferentialPlot.cc.

References alignCSCRings::e, MonitorElement::getTH1F(), processed, theDifferentialHistoB_c, theDifferentialHistoB_dus, theDifferentialHistoB_g, and theDifferentialHistoB_ni.

Referenced by epsPlot(), plot(), and psPlot().

45  {
46 
47 // thePlotCanvas = new TCanvas( commonName ,
48 // commonName ,
49 // btppXCanvas , btppYCanvas ) ;
50 //
51 // if ( !btppTitle ) gStyle->SetOptTitle ( 0 ) ;
52 
53  if (!processed) return;
54 //fixme:
55  bool btppNI = false;
56  bool btppColour = true;
57 
58  thePlotCanvas.SetFillColor ( 0 ) ;
59  thePlotCanvas.cd ( 1 ) ;
60  gPad->SetLogy ( 1 ) ;
61  gPad->SetGridx ( 1 ) ;
62  gPad->SetGridy ( 1 ) ;
63 
64 // int col_b ;
65  int col_c ;
66  int col_g ;
67  int col_dus ;
68  int col_ni ;
69 
70 // int mStyle_b ;
71  int mStyle_c ;
72  int mStyle_g ;
73  int mStyle_dus ;
74  int mStyle_ni ;
75 
76  // marker size (same for all)
77  float mSize = 1.5 ;
78 
79  if ( btppColour ) {
80 // col_b = 2 ;
81  col_c = 6 ;
82  col_g = 3 ;
83  col_dus = 4 ;
84  col_ni = 5 ;
85 // mStyle_b = 20 ;
86  mStyle_c = 20 ;
87  mStyle_g = 20 ;
88  mStyle_dus = 20 ;
89  mStyle_ni = 20 ;
90  }
91  else {
92 // col_b = 1 ;
93  col_c = 1 ;
94  col_g = 1 ;
95  col_dus = 1 ;
96  col_ni = 1 ;
97 // mStyle_b = 12 ;
98  mStyle_c = 22 ;
99  mStyle_g = 29 ;
100  mStyle_dus = 20 ;
101  mStyle_ni = 27 ;
102  }
103 
104  // for the moment: plot b (to see what the constant b-efficiency is), c, g, uds
105  // b in red
106  // No, do not plot b (because only visible for the soft leptons)
107  // theDifferentialHistoB_b -> GetXaxis()->SetTitle ( diffVariableName ) ;
108  // theDifferentialHistoB_b -> GetYaxis()->SetTitle ( "non b-jet efficiency" ) ;
109  // theDifferentialHistoB_b -> GetYaxis()->SetTitleOffset ( 1.25 ) ;
110  // theDifferentialHistoB_b -> SetMaximum ( 0.4 ) ;
111  // theDifferentialHistoB_b -> SetMinimum ( 1.e-4 ) ;
112  // theDifferentialHistoB_b -> SetMarkerColor ( col_b ) ;
113  // theDifferentialHistoB_b -> SetLineColor ( col_b ) ;
114  // theDifferentialHistoB_b -> SetMarkerSize ( mSize ) ;
115  // theDifferentialHistoB_b -> SetMarkerStyle ( mStyle_b ) ;
116  // theDifferentialHistoB_b -> SetStats ( false ) ;
117  // theDifferentialHistoB_b -> Draw ( "pe" ) ;
118  // c in magenta
119  theDifferentialHistoB_c ->getTH1F() -> SetMaximum ( 0.4 ) ;
120  theDifferentialHistoB_c ->getTH1F() -> SetMinimum ( 1.e-4 ) ;
121  theDifferentialHistoB_c ->getTH1F() -> SetMarkerColor ( col_c ) ;
122  theDifferentialHistoB_c ->getTH1F() -> SetLineColor ( col_c ) ;
123  theDifferentialHistoB_c ->getTH1F() -> SetMarkerSize ( mSize ) ;
124  theDifferentialHistoB_c ->getTH1F() -> SetMarkerStyle ( mStyle_c ) ;
125  theDifferentialHistoB_c ->getTH1F() -> SetStats ( false ) ;
126  // theDifferentialHistoB_c -> Draw("peSame") ;
127  theDifferentialHistoB_c ->getTH1F()-> Draw("pe") ;
128  // uds in blue
129  theDifferentialHistoB_dus ->getTH1F()-> SetMarkerColor ( col_dus ) ;
130  theDifferentialHistoB_dus ->getTH1F()-> SetLineColor ( col_dus ) ;
131  theDifferentialHistoB_dus ->getTH1F()-> SetMarkerSize ( mSize ) ;
132  theDifferentialHistoB_dus ->getTH1F()-> SetMarkerStyle ( mStyle_dus ) ;
133  theDifferentialHistoB_dus ->getTH1F()-> SetStats ( false ) ;
134  theDifferentialHistoB_dus ->getTH1F()-> Draw("peSame") ;
135  // g in green
136  // only uds not to confuse
137  theDifferentialHistoB_g ->getTH1F()-> SetMarkerColor ( col_g ) ;
138  theDifferentialHistoB_g ->getTH1F()-> SetLineColor ( col_g ) ;
139  theDifferentialHistoB_g ->getTH1F()-> SetMarkerSize ( mSize ) ;
140  theDifferentialHistoB_g ->getTH1F()-> SetMarkerStyle ( mStyle_g ) ;
141  theDifferentialHistoB_g ->getTH1F()-> SetStats ( false ) ;
142  theDifferentialHistoB_g ->getTH1F()-> Draw("peSame") ;
143 
144  // NI if wanted
145  if ( btppNI ) {
146  theDifferentialHistoB_ni ->getTH1F()-> SetMarkerColor ( col_ni ) ;
147  theDifferentialHistoB_ni ->getTH1F()-> SetLineColor ( col_ni ) ;
148  theDifferentialHistoB_ni ->getTH1F()-> SetMarkerSize ( mSize ) ;
149  theDifferentialHistoB_ni ->getTH1F()-> SetMarkerStyle ( mStyle_ni ) ;
150  theDifferentialHistoB_ni ->getTH1F()-> SetStats ( false ) ;
151  theDifferentialHistoB_ni ->getTH1F()-> Draw("peSame") ;
152  }
153 }
MonitorElement * theDifferentialHistoB_g
MonitorElement * theDifferentialHistoB_dus
TH1F * getTH1F(void) const
MonitorElement * theDifferentialHistoB_c
MonitorElement * theDifferentialHistoB_ni
void BTagDifferentialPlot::plot ( const std::string &  name,
const std::string &  ext 
)

Definition at line 165 of file BTagDifferentialPlot.cc.

References commonName, plot(), and processed.

166 {
167  if (!processed) return;
168  TCanvas tc(commonName.c_str(), commonName.c_str());
169  plot(tc);
170  tc.Print((name + commonName + ext).c_str());
171 }
void plot(TCanvas &theCanvas)
void BTagDifferentialPlot::process ( )
void BTagDifferentialPlot::psPlot ( const std::string &  name)

Definition at line 160 of file BTagDifferentialPlot.cc.

References plot().

161 {
162  plot(name, ".ps");
163 }
void plot(TCanvas &theCanvas)
void BTagDifferentialPlot::setVariableName ( )
private

Definition at line 183 of file BTagDifferentialPlot.cc.

References constETA, constPT, constVar, constVariableName, constVariableValue, diffVariableName, and theBinPlotters.

Referenced by process().

184 {
185  if ( constVar==constETA ) {
186  constVariableName = "eta" ;
187  diffVariableName = "pt" ;
188  constVariableValue = make_pair ( theBinPlotters[0]->etaPtBin().getEtaMin() , theBinPlotters[0]->etaPtBin().getEtaMax() ) ;
189  }
190  if ( constVar==constPT ) {
191  constVariableName = "pt" ;
192  diffVariableName = "eta" ;
193  constVariableValue = make_pair ( theBinPlotters[0]->etaPtBin().getPtMin() , theBinPlotters[0]->etaPtBin().getPtMax() ) ;
194  }
195 
196  /* std::cout
197  << "====>>>> BTagDifferentialPlot::setVariableName() : set const/diffVariableName to : "
198  << constVariableName << " / " << diffVariableName << endl
199  << "====>>>> constant value interval : "
200  << constVariableValue.first << " - " << constVariableValue.second << endl ;
201  */
202 }
std::pair< double, double > constVariableValue
std::vector< JetTagPlotter * > theBinPlotters

Member Data Documentation

std::string BTagDifferentialPlot::commonName
private

Definition at line 85 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), and plot().

ConstVarType BTagDifferentialPlot::constVar
private

Definition at line 75 of file BTagDifferentialPlot.h.

Referenced by setVariableName().

std::string BTagDifferentialPlot::constVariableName
private

Definition at line 77 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), and setVariableName().

std::pair<double,double> BTagDifferentialPlot::constVariableValue
private

Definition at line 82 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), and setVariableName().

std::string BTagDifferentialPlot::diffVariableName
private

Definition at line 79 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), fillHisto(), and setVariableName().

double BTagDifferentialPlot::fixedBEfficiency
private

Definition at line 69 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), and fillHisto().

bool BTagDifferentialPlot::noProcessing
private

Definition at line 72 of file BTagDifferentialPlot.h.

Referenced by process().

bool BTagDifferentialPlot::processed
private

Definition at line 73 of file BTagDifferentialPlot.h.

Referenced by plot(), and process().

std::vector<JetTagPlotter *> BTagDifferentialPlot::theBinPlotters
private

Definition at line 89 of file BTagDifferentialPlot.h.

Referenced by addBinPlotter(), bookHisto(), fillHisto(), and setVariableName().

MonitorElement* BTagDifferentialPlot::theDifferentialHistoB_b
private

Definition at line 96 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), fillHisto(), and getDifferentialHistoB_b().

MonitorElement* BTagDifferentialPlot::theDifferentialHistoB_c
private

Definition at line 95 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), fillHisto(), getDifferentialHistoB_c(), and plot().

MonitorElement* BTagDifferentialPlot::theDifferentialHistoB_d
private

Definition at line 92 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), fillHisto(), and getDifferentialHistoB_d().

MonitorElement* BTagDifferentialPlot::theDifferentialHistoB_dus
private

Definition at line 99 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), fillHisto(), getDifferentialHistoB_dus(), and plot().

MonitorElement* BTagDifferentialPlot::theDifferentialHistoB_dusg
private

Definition at line 100 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), fillHisto(), and getDifferentialHistoB_dusg().

MonitorElement* BTagDifferentialPlot::theDifferentialHistoB_g
private

Definition at line 97 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), fillHisto(), getDifferentialHistoB_g(), and plot().

MonitorElement* BTagDifferentialPlot::theDifferentialHistoB_ni
private

Definition at line 98 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), fillHisto(), getDifferentialHistoB_ni(), and plot().

MonitorElement* BTagDifferentialPlot::theDifferentialHistoB_s
private

Definition at line 94 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), fillHisto(), and getDifferentialHistoB_s().

MonitorElement* BTagDifferentialPlot::theDifferentialHistoB_u
private

Definition at line 93 of file BTagDifferentialPlot.h.

Referenced by bookHisto(), fillHisto(), and getDifferentialHistoB_u().