Reference Manual

Go to the documentation of this file.
00001 #include "DQMOffline/RecoB/interface/BTagDifferentialPlot.h"
00002 #include "DQMOffline/RecoB/interface/EffPurFromHistos.h"
00003 #include "DQMOffline/RecoB/interface/Tools.h"
00004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00005 #include "TF1.h"
00006 
00007 #include "DQMServices/Core/interface/DQMStore.h"
00008 #include "DQMServices/Core/interface/MonitorElement.h"
00009 #include "FWCore/ServiceRegistry/interface/Service.h" 
00010 
00011 #include "DQMOffline/RecoB/interface/HistoProviderDQM.h"
00012 
00013 #include <algorithm>
00014 #include <iostream>
00015 #include <sstream>
00016 using namespace std ;
00017 
00018 
00019 
00020 BTagDifferentialPlot::BTagDifferentialPlot (const double& bEff, const ConstVarType& constVariable,
00021         const std::string & tagName) :
00022         fixedBEfficiency     ( bEff )  ,
00023         noProcessing         ( false ) , processed(false), constVar(constVariable),
00024         constVariableName    ( "" )    , diffVariableName     ( "" )    ,
00025         constVariableValue   ( 999.9 , 999.9 ) , commonName( "MisidForBEff_" + tagName+"_") ,
00026         theDifferentialHistoB_d    ( 0 ) ,
00027         theDifferentialHistoB_u    ( 0 ) ,
00028         theDifferentialHistoB_s    ( 0 ) ,
00029         theDifferentialHistoB_c    ( 0 ) ,
00030         theDifferentialHistoB_b    ( 0 ) ,
00031         theDifferentialHistoB_g    ( 0 ) ,
00032         theDifferentialHistoB_ni   ( 0 ) ,
00033         theDifferentialHistoB_dus  ( 0 ) ,
00034         theDifferentialHistoB_dusg ( 0 )  {}
00035 
00036 
00037 BTagDifferentialPlot::~BTagDifferentialPlot () {
00038 }
00039 
00040 
00041 
00042 
00043 
00044 
00045 void BTagDifferentialPlot::plot (TCanvas & thePlotCanvas ) {
00046 
00047 //   thePlotCanvas = new TCanvas(  commonName ,
00048 //                              commonName ,
00049 //                              btppXCanvas , btppYCanvas ) ;
00050 //
00051 //   if ( !btppTitle ) gStyle->SetOptTitle ( 0 ) ;
00052 
00053   if (!processed) return;
00054 //fixme:
00055   bool btppNI = false;
00056   bool btppColour = true;
00057 
00058   thePlotCanvas.SetFillColor ( 0 ) ;
00059   thePlotCanvas.cd ( 1 ) ;
00060   gPad->SetLogy  ( 1 ) ;
00061   gPad->SetGridx ( 1 ) ;
00062   gPad->SetGridy ( 1 ) ;
00063 
00064   int col_b   ;
00065   int col_c   ;
00066   int col_g   ;
00067   int col_dus ;
00068   int col_ni  ;
00069 
00070   int mStyle_b   ;
00071   int mStyle_c   ;
00072   int mStyle_g   ;
00073   int mStyle_dus ;
00074   int mStyle_ni  ;
00075 
00076   // marker size (same for all)
00077   float mSize = 1.5 ;
00078 
00079   if ( btppColour ) {
00080     col_b    = 2 ;
00081     col_c    = 6 ;
00082     col_g    = 3 ;
00083     col_dus  = 4 ;
00084     col_ni   = 5 ;
00085     mStyle_b   = 20 ;
00086     mStyle_c   = 20 ;
00087     mStyle_g   = 20 ;
00088     mStyle_dus = 20 ;
00089     mStyle_ni  = 20 ;
00090   }
00091   else {
00092     col_b    = 1 ;
00093     col_c    = 1 ;
00094     col_g    = 1 ;
00095     col_dus  = 1 ;
00096     col_ni   = 1 ;
00097     mStyle_b   = 12 ;
00098     mStyle_c   = 22 ;
00099     mStyle_g   = 29 ;
00100     mStyle_dus = 20 ;
00101     mStyle_ni  = 27 ;
00102   }
00103 
00104   // for the moment: plot b (to see what the constant b-efficiency is), c, g, uds
00105   // b in red
00106   // No, do not plot b (because only visible for the soft leptons)
00107   // theDifferentialHistoB_b   -> GetXaxis()->SetTitle ( diffVariableName ) ;
00108   // theDifferentialHistoB_b   -> GetYaxis()->SetTitle ( "non b-jet efficiency" ) ;
00109   // theDifferentialHistoB_b   -> GetYaxis()->SetTitleOffset ( 1.25 ) ;
00110   // theDifferentialHistoB_b   -> SetMaximum ( 0.4 )  ;
00111   // theDifferentialHistoB_b   -> SetMinimum ( 1.e-4 )  ;
00112   // theDifferentialHistoB_b   -> SetMarkerColor ( col_b ) ;
00113   // theDifferentialHistoB_b   -> SetLineColor   ( col_b ) ;
00114   // theDifferentialHistoB_b   -> SetMarkerSize  ( mSize ) ;
00115   // theDifferentialHistoB_b   -> SetMarkerStyle ( mStyle_b ) ;
00116   // theDifferentialHistoB_b   -> SetStats ( false ) ;
00117   // theDifferentialHistoB_b   -> Draw ( "pe" ) ;
00118   // c in magenta
00119   theDifferentialHistoB_c ->getTH1F()  -> SetMaximum ( 0.4 )  ;
00120   theDifferentialHistoB_c ->getTH1F()  -> SetMinimum ( 1.e-4 )  ;
00121   theDifferentialHistoB_c ->getTH1F()  -> SetMarkerColor ( col_c ) ;
00122   theDifferentialHistoB_c ->getTH1F()  -> SetLineColor   ( col_c ) ;
00123   theDifferentialHistoB_c ->getTH1F()  -> SetMarkerSize  ( mSize ) ;
00124   theDifferentialHistoB_c ->getTH1F()  -> SetMarkerStyle ( mStyle_c ) ;
00125   theDifferentialHistoB_c ->getTH1F()  -> SetStats     ( false ) ;
00126   //  theDifferentialHistoB_c   -> Draw("peSame") ;
00127   theDifferentialHistoB_c   ->getTH1F()-> Draw("pe") ;
00128   // uds in blue
00129   theDifferentialHistoB_dus ->getTH1F()-> SetMarkerColor ( col_dus ) ;
00130   theDifferentialHistoB_dus ->getTH1F()-> SetLineColor   ( col_dus ) ;
00131   theDifferentialHistoB_dus ->getTH1F()-> SetMarkerSize  ( mSize ) ;
00132   theDifferentialHistoB_dus ->getTH1F()-> SetMarkerStyle ( mStyle_dus ) ;
00133   theDifferentialHistoB_dus ->getTH1F()-> SetStats     ( false ) ;
00134   theDifferentialHistoB_dus ->getTH1F()-> Draw("peSame") ;
00135   // g in green
00136   // only uds not to confuse
00137   theDifferentialHistoB_g   ->getTH1F()-> SetMarkerColor ( col_g ) ;
00138   theDifferentialHistoB_g   ->getTH1F()-> SetLineColor   ( col_g ) ;
00139   theDifferentialHistoB_g   ->getTH1F()-> SetMarkerSize  ( mSize ) ;
00140   theDifferentialHistoB_g   ->getTH1F()-> SetMarkerStyle ( mStyle_g ) ;
00141   theDifferentialHistoB_g   ->getTH1F()-> SetStats     ( false ) ;
00142   theDifferentialHistoB_g   ->getTH1F()-> Draw("peSame") ;
00143 
00144   // NI if wanted
00145   if ( btppNI ) {
00146     theDifferentialHistoB_ni ->getTH1F()-> SetMarkerColor ( col_ni ) ;
00147     theDifferentialHistoB_ni ->getTH1F()-> SetLineColor   ( col_ni ) ;
00148     theDifferentialHistoB_ni ->getTH1F()-> SetMarkerSize  ( mSize ) ;
00149     theDifferentialHistoB_ni ->getTH1F()-> SetMarkerStyle ( mStyle_ni ) ;
00150     theDifferentialHistoB_ni ->getTH1F()-> SetStats     ( false ) ;
00151     theDifferentialHistoB_ni ->getTH1F()-> Draw("peSame") ;
00152   }
00153 }
00154 
00155 void BTagDifferentialPlot::epsPlot(const std::string & name)
00156 {
00157   plot(name, ".eps");
00158 }
00159 
00160 void BTagDifferentialPlot::psPlot(const std::string & name)
00161 {
00162   plot(name, ".ps");
00163 }
00164 
00165 void BTagDifferentialPlot::plot(const std::string & name, const std::string & ext)
00166 {
00167   if (!processed) return;
00168    TCanvas tc(commonName.c_str(), commonName.c_str());
00169    plot(tc);
00170    tc.Print((name + commonName + ext).c_str());
00171 }
00172 
00173 
00174 void BTagDifferentialPlot::process () {
00175   setVariableName () ; // also sets noProcessing if not OK
00176   if ( noProcessing ) return ;
00177   bookHisto () ;
00178   fillHisto () ;
00179   processed = true;
00180 }
00181 
00182 
00183 void BTagDifferentialPlot::setVariableName ()
00184 {
00185   if ( constVar==constETA ) {
00186     constVariableName  = "eta" ;
00187     diffVariableName   = "pt"  ;
00188     constVariableValue = make_pair ( theBinPlotters[0]->etaPtBin().getEtaMin() , theBinPlotters[0]->etaPtBin().getEtaMax() ) ;
00189   }
00190   if ( constVar==constPT  ) {
00191     constVariableName = "pt"  ;
00192     diffVariableName  = "eta" ;
00193     constVariableValue = make_pair ( theBinPlotters[0]->etaPtBin().getPtMin() , theBinPlotters[0]->etaPtBin().getPtMax() ) ;
00194   }
00195 
00196   /*  std::cout
00197      << "====>>>> BTagDifferentialPlot::setVariableName() : set const/diffVariableName to : "
00198      << constVariableName << " / " << diffVariableName << endl
00199      << "====>>>>                                            constant value interval : "
00200      << constVariableValue.first  << " - " << constVariableValue.second << endl ;
00201   */
00202 }
00203 
00204 
00205 
00206 void BTagDifferentialPlot::bookHisto () {
00207 
00208   // vector with ranges
00209   vector<float> variableRanges ;
00210 
00211   for ( vector<JetTagPlotter *>::const_iterator iP = theBinPlotters.begin() ; 
00212         iP != theBinPlotters.end() ; ++iP ) {
00213     const EtaPtBin & currentBin = (*iP)->etaPtBin()  ;
00214     if ( diffVariableName == "eta" ) {
00215       // only active bins in the variable on x-axis
00216       if ( currentBin.getEtaActive() ) {
00217         variableRanges.push_back ( currentBin.getEtaMin() ) ;
00218         // also max if last one
00219         if ( iP == --theBinPlotters.end() ) variableRanges.push_back ( currentBin.getEtaMax() ) ;
00220       }
00221     }
00222     if ( diffVariableName == "pt" ) {
00223       // only active bins in the variable on x-axis
00224       if ( currentBin.getPtActive() ) {
00225         variableRanges.push_back ( currentBin.getPtMin() ) ;
00226         // also max if last one
00227         if ( iP == --theBinPlotters.end() ) variableRanges.push_back ( currentBin.getPtMax() ) ;
00228       }
00229     }
00230   }
00231 
00232   // to book histo with variable binning -> put into array
00233   int      nBins    = variableRanges.size() - 1 ;
00234   float * binArray = &variableRanges[0];
00235   //float * binArray = new float [nBins+1] ;
00236 
00237   //for ( int i = 0 ; i < nBins + 1 ; i++ ) {
00238   //  binArray[i] = variableRanges[i] ;
00239   //}
00240 
00241 
00242   // part of the name common to all flavours
00243   std::stringstream stream("");
00244   stream << fixedBEfficiency << "_Const_" << constVariableName << "_" << constVariableValue.first << "-" ;
00245   stream << constVariableValue.second << "_" << "_Vs_" << diffVariableName ;
00246   commonName += stream.str();
00247   std::remove(commonName.begin(), commonName.end(), ' ') ;
00248   std::replace(commonName.begin(), commonName.end(), '.' , 'v' ) ;
00249 
00250   std::string label(commonName);
00251   HistoProviderDQM prov ("Btag",label);
00252 
00253   theDifferentialHistoB_d    = (prov.book1D ( "D_"    + commonName , "D_"    + commonName , nBins , binArray )) ;
00254   theDifferentialHistoB_u    = (prov.book1D ( "U_"    + commonName , "U_"    + commonName , nBins , binArray )) ;
00255   theDifferentialHistoB_s    = (prov.book1D ( "S_"    + commonName , "S_"    + commonName , nBins , binArray )) ;
00256   theDifferentialHistoB_c    = (prov.book1D ( "C_"    + commonName , "C_"    + commonName , nBins , binArray )) ;
00257   theDifferentialHistoB_b    = (prov.book1D ( "B_"    + commonName , "B_"    + commonName , nBins , binArray )) ;
00258   theDifferentialHistoB_g    = (prov.book1D ( "G_"    + commonName , "G_"    + commonName , nBins , binArray )) ;
00259   theDifferentialHistoB_ni   = (prov.book1D ( "NI_"   + commonName , "NI_"   + commonName , nBins , binArray )) ;
00260   theDifferentialHistoB_dus  = (prov.book1D ( "DUS_"  + commonName , "DUS_"  + commonName , nBins , binArray )) ;
00261   theDifferentialHistoB_dusg = (prov.book1D ( "DUSG_" + commonName , "DUSG_" + commonName , nBins , binArray )) ;
00262 }
00263 
00264 
00265 void BTagDifferentialPlot::fillHisto () {
00266   // loop over bins and find corresponding misid. in the MisIdVs..... histo
00267   for ( vector<JetTagPlotter *>::const_iterator iP = theBinPlotters.begin() ;
00268         iP != theBinPlotters.end() ; ++iP ) {
00269     const EtaPtBin   & currentBin              = (*iP)->etaPtBin() ;
00270     EffPurFromHistos * currentEffPurFromHistos = (*iP)->getEffPurFromHistos() ;
00271     //
00272     bool   isActive   = true ;
00273     double valueXAxis = -999.99 ;
00274     // find right bin based on middle of the interval
00275     if ( diffVariableName == "eta" ) {
00276       isActive = currentBin.getEtaActive() ;
00277       valueXAxis = 0.5 * ( currentBin.getEtaMin() + currentBin.getEtaMax() ) ;
00278     } else if ( diffVariableName == "pt"  ) {
00279       isActive = currentBin.getPtActive() ;
00280       valueXAxis = 0.5 * ( currentBin.getPtMin() + currentBin.getPtMax() ) ;
00281     } else {
00282       throw cms::Exception("Configuration")
00283         << "====>>>> BTagDifferentialPlot::fillHisto() : illegal diffVariableName = " << diffVariableName << endl;
00284     }
00285 
00286     // for the moment: ignore inactive bins
00287     // (maybe later: if a Bin is inactive -> set value to fill well below left edge of histogram to have it in the underflow)
00288 
00289     if ( !isActive ) continue ;
00290 
00291     // to have less lines of code ....
00292     vector< pair<TH1F*,TH1F*> > effPurDifferentialPairs ;
00293 
00294     // all flavours (b is a good cross check! must be constant and = fixed b-efficiency)
00295     // get histo; find the bin of the fixed b-efficiency in the histo and get misid; fill
00296 
00297 
00298     effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_d()    , theDifferentialHistoB_d ->getTH1F()   ) ) ;
00299     effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_u()    , theDifferentialHistoB_u ->getTH1F()   ) ) ;
00300     effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_s()    , theDifferentialHistoB_s ->getTH1F()   ) ) ;
00301     effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_c()    , theDifferentialHistoB_c  ->getTH1F()  ) ) ;
00302     effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_b()    , theDifferentialHistoB_b  ->getTH1F()  ) ) ;
00303     effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_g()    , theDifferentialHistoB_g  ->getTH1F()  ) ) ;
00304     effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_ni()   , theDifferentialHistoB_ni ->getTH1F()  ) ) ;
00305     effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_dus()  , theDifferentialHistoB_dus->getTH1F()  ) ) ;
00306     effPurDifferentialPairs.push_back ( make_pair ( currentEffPurFromHistos->getEffFlavVsBEff_dusg() , theDifferentialHistoB_dusg->getTH1F() ) ) ;
00307 
00308     for ( vector< pair<TH1F*,TH1F*> >::const_iterator itP  = effPurDifferentialPairs.begin() ;
00309                                                       itP != effPurDifferentialPairs.end()   ; ++itP ) {
00310       TH1F * effPurHist = itP->first  ;
00311       TH1F * diffHist   = itP->second ;
00312       pair<double, double> mistag = getMistag(fixedBEfficiency, effPurHist);
00313       int iBinSet = diffHist->FindBin(valueXAxis) ;
00314       diffHist->SetBinContent(iBinSet, mistag.first);
00315       diffHist->SetBinError(iBinSet, mistag.second);
00316     }
00317   }
00318 
00319 }
00320 
00321 pair<double, double>
00322 BTagDifferentialPlot::getMistag(const double& fixedBEfficiency, TH1F * effPurHist)
00323 {
00324   int iBinGet = effPurHist->FindBin ( fixedBEfficiency ) ;
00325   double effForBEff    = effPurHist->GetBinContent ( iBinGet ) ;
00326   double effForBEffErr = effPurHist->GetBinError   ( iBinGet ) ;
00327 
00328   if (effForBEff==0. && effForBEffErr==0.) {
00329     // The bin was empty. Could be that it was not filled, as the scan-plot
00330     //  did not have an entry at the requested value, or that the curve
00331     // is above or below.
00332     // Fit a plynomial, and evaluate the mistag at the requested value.
00333     int fitStatus;
00334     try {
00335     fitStatus = effPurHist->Fit("pol4", "q");
00336     }catch (...){
00337   return pair<double, double>(effForBEff, effForBEffErr);
00338 }
00339     if (fitStatus != 0) {
00340       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.";
00341       //    } else {
00342       //      edm::LogInfo("BTagDifferentialPlot")<<"Fit OK to hisogram " << effPurHist->GetTitle() << " entries = " << effPurHist->GetEntries();
00343       return pair<double, double>(effForBEff, effForBEffErr);
00344     }
00345     TF1 *myfunc = effPurHist->GetFunction("pol4");
00346     effForBEff = myfunc->Eval(fixedBEfficiency);
00347     effPurHist->RecursiveRemove(myfunc);
00348     //Error: first non-empty bin on the right and take its error
00349     for (int i = iBinGet+1; i< effPurHist->GetNbinsX(); ++i) {
00350       if (effPurHist->GetBinContent(i)!=0) {
00351         effForBEffErr = effPurHist->GetBinError(i);
00352         break;
00353       }
00354     }
00355   }
00356 
00357   return pair<double, double>(effForBEff, effForBEffErr);
00358 }
00359
/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/DQMOffline/RecoB/src/BTagDifferentialPlot.cc
