CMSSW Reference Manual

00001 // -*- C++ -*-
00002 //
00003 // Package:    ImpactParameterCalibration
00004 // Class:      ImpactParameterCalibration
00005 // 
00013 //
00014 // Original Author:  Jeremy Andrea/Andrea Rizzi
00015 //         Created:  Mon Aug  6 16:10:38 CEST 2007
00016 // $Id: ImpactParameterCalibration.cc,v 1.9 2008/03/06 08:18:49 tboccali Exp $
00017 //
00018 //
00019 // system include files
00020 
00021 #include "FWCore/Framework/interface/Frameworkfwd.h"
00022 #include "FWCore/Framework/interface/EDAnalyzer.h"
00023 #include "FWCore/Framework/interface/Event.h"
00024 #include "FWCore/Framework/interface/MakerMacros.h"
00025 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00026 #include "FWCore/ParameterSet/interface/InputTag.h"
00027 #include "FWCore/Framework/interface/ESHandle.h"
00028 #include "FWCore/Framework/interface/EventSetup.h"
00029 #include "Utilities/General/interface/FileInPath.h"
00030 
00031 #include "FWCore/Framework/interface/IOVSyncValue.h"
00032 #include "FWCore/ServiceRegistry/interface/Service.h"
00033 #include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
00034 
00035 
00036 #include "DataFormats/BTauReco/interface/TrackIPTagInfo.h"
00037 #include "CondFormats/BTauObjects/interface/TrackProbabilityCalibration.h"
00038 #include "CondFormats/BTauObjects/interface/CalibratedHistogram.h"
00039 
00040 #include "CondFormats/DataRecord/interface/BTagTrackProbability2DRcd.h"
00041 #include "CondFormats/DataRecord/interface/BTagTrackProbability3DRcd.h"
00042 
00043 #include <TClass.h>
00044 
00045 
00046 #include "RVersion.h"
00047 #if ROOT_VERSION_CODE >= ROOT_VERSION(5,15,0)
00048 #include "TBufferFile.h"
00049 typedef TBufferFile MyTBuffer;
00050 #else
00051 #include "TBuffer.h"
00052 typedef TBuffer MyTBuffer;
00053 #endif
00054 
00055 
00056 #include <TBufferXML.h>
00057 #include <iostream>
00058 #include <fstream>
00059 #include <sstream>
00060 #include <string>
00061 #include <vector>
00062 #include <memory>
00063 
00064 #include "TrackClassMatch.h"
00065 
00066 
00067 using namespace reco;
00068 using namespace std;
00069 //
00070 // class decleration
00071 //
00072 
00073 
00074 
00075 class ImpactParameterCalibration : public edm::EDAnalyzer {
00076    public:
00077       explicit ImpactParameterCalibration(const edm::ParameterSet&);
00078       ~ImpactParameterCalibration();
00079 
00080 
00081    private:
00082       virtual void beginJob(const edm::EventSetup&) ;
00083       virtual void analyze(const edm::Event&, const edm::EventSetup&);
00084       virtual void endJob() ;
00085       edm::ParameterSet config;
00086 
00087    TrackProbabilityCalibration * fromXml(edm::FileInPath xmlCalibration);
00088 
00089    static TrackProbabilityCategoryData createCategory(double  pmin,double  pmax,
00090                  double  etamin,  double  etamax,
00091                   int  nhitmin,  int  nhitmax,
00092                  int  npixelhitsmin, int  npixelhitsmax,
00093                   double cmin, double cmax, int withFirst)
00094   {
00095   TrackProbabilityCategoryData c;
00096   c.pMin=pmin;
00097   c.pMax=pmax;
00098   c.etaMin=etamin;
00099   c.etaMax=etamax;
00100   c.nHitsMin=nhitmin;
00101   c.nHitsMax=nhitmax;
00102   c.nPixelHitsMin=npixelhitsmin;
00103   c.nPixelHitsMax=npixelhitsmax;
00104   c.chiMin=cmin;
00105   c.chiMax=cmax;
00106   c.withFirstPixel=withFirst;
00107   return c;
00108  }
00109 
00110    TrackProbabilityCalibration * m_calibration[2];
00111    edm::InputTag m_iptaginfo;
00112    edm::InputTag m_pv;
00113   unsigned int minLoop, maxLoop;
00114 
00115 };
00116 
00117 ImpactParameterCalibration::ImpactParameterCalibration(const edm::ParameterSet& iConfig):config(iConfig)
00118 {
00119   m_iptaginfo = iConfig.getParameter<edm::InputTag>("tagInfoSrc");
00120   m_pv = iConfig.getParameter<edm::InputTag>("primaryVertexSrc");
00121   bool createOnlyOne = iConfig.getUntrackedParameter<bool>("createOnlyOneCalibration", false);
00122   minLoop=0;
00123   maxLoop=1;
00124   if (createOnlyOne == true){
00125     int whichCalib = iConfig.getUntrackedParameter<int>("dimension", 2);
00126     if (whichCalib==2){
00127       std::cout <<" Writing only 2D calibrations"<<std::endl;
00128       minLoop=1;
00129       maxLoop=1;
00130     }else if (whichCalib==3){
00131       std::cout <<" Writing only 3D calibrations"<<std::endl;
00132       minLoop=0;
00133       maxLoop=0;
00134     }else {
00135       std::cout <<" Dimension not found: "<<whichCalib<<"; it must be either 2 or 3"<<std::endl;
00136     }
00137   }
00138 
00139 }
00140 
00141 
00142 ImpactParameterCalibration::~ImpactParameterCalibration()
00143 {
00144 }
00145 
00146 
00147 void
00148 ImpactParameterCalibration::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00149 {
00150   using namespace edm;
00151   using namespace reco;
00152 
00153   Handle<TrackIPTagInfoCollection> ipHandle;
00154   iEvent.getByLabel(m_iptaginfo, ipHandle);
00155   const TrackIPTagInfoCollection & ip = *(ipHandle.product());
00156 
00157 //  cout << "Found " << ip.size() << " TagInfo" << endl;
00158 
00159   Handle<reco::VertexCollection> primaryVertex;
00160   iEvent.getByLabel(m_pv,primaryVertex);
00161 
00162   vector<TrackProbabilityCalibration::Entry>::iterator found;
00163   vector<TrackProbabilityCalibration::Entry>::iterator it_begin;
00164   vector<TrackProbabilityCalibration::Entry>::iterator it_end;
00165 
00166       
00167    TrackIPTagInfoCollection::const_iterator it = ip.begin();
00168    for(; it != ip.end(); it++)
00169      {
00170       TrackRefVector selTracks=it->selectedTracks();
00171 //      if(it->primaryVertex().isNull()) continue;
00172       if(primaryVertex.product()->size() == 0) 
00173        {
00174        std::cout << "No PV in the event!!" << std::endl;
00175          continue;
00176         }
00177        const Vertex & pv = *(primaryVertex.product()->begin());
00178            
00179       for(unsigned int i=minLoop; i <= maxLoop;i++)
00180       { 
00181         it_begin=m_calibration[i]->data.begin();
00182         it_end=m_calibration[i]->data.end();
00183   
00184       for(unsigned int j=0;j<ip.size(); j++)
00185         {
00186           double ipsig;
00187           if (i==0) ipsig  = it->impactParameterData()[j].ip3d.significance();
00188           else  ipsig  = it->impactParameterData()[j].ip2d.significance();
00189           TrackClassMatch::Input input(*selTracks[j],*it->jet(),pv);
00190           if(ipsig < 0) 
00191            {
00192             found = std::find_if(it_begin,it_end,bind1st(TrackClassMatch(),input));
00193 //            std::cout << ip[j].significance() << std::endl; 
00194             if(found!=it_end) 
00195               found->histogram.fill(-ipsig);
00196             else
00197               std::cout << "No category for this track!!" << std::endl;
00198            }
00199          }
00200       } 
00201      }  
00202       
00203       
00204   
00205 }
00206 
00207 
00208 
00209 
00210 
00211 
00212 
00213 
00214 
00215 
00216 // ------------ method called once each job just before starting event loop  ------------
00217 void 
00218 ImpactParameterCalibration::beginJob(const edm::EventSetup & iSetup)
00219 {
00220   using namespace edm;
00221   m_calibration[0] =   new TrackProbabilityCalibration();
00222   m_calibration[1] =   new TrackProbabilityCalibration();
00223 
00224   CalibratedHistogram hist(config.getParameter<int>("nBins"),0,config.getParameter<double>("maxSignificance"));
00225 
00226   std::string categories = config.getParameter<std::string>("inputCategories");
00227 
00228   if(categories == "HardCoded")
00229   {
00230    vector<TrackProbabilityCategoryData> v;
00231   v.push_back(createCategory(0,5000,0,2.4,8,50,1,1,0,5,0));
00232   v.push_back(createCategory(0,5000,0,2.4,8,50,2,5,2.5,5,0));
00233   v.push_back(createCategory(0,8,0,0.8,8,50,3,5,0,2.5,0));
00234   v.push_back(createCategory(0,8,0.8,1.6,8,50,3,5,0,2.5,0));
00235   v.push_back(createCategory(0,8,1.6,2.4,8,50,3,5,0,2.5,0));
00236   v.push_back(createCategory(0,8,0,2.4,8,50,2,2,0,2.5,0));
00237   v.push_back(createCategory(8,5000,0,0.8,8,50,3,5,0,2.5,0));
00238   v.push_back(createCategory(8,5000,0.8,1.6,8,50,3,5,0,2.5,0));
00239   v.push_back(createCategory(8,5000,1.6,2.4,8,50,3,5,0,2.5,0));
00240   v.push_back(createCategory(8,5000,0,2.4,8,50,2,2,0,2.5,0));
00241   for(unsigned int i=minLoop;i <=maxLoop ;i++)
00242    for(unsigned int j=0;j<v.size() ; j++)
00243     {
00244      TrackProbabilityCalibration::Entry e;
00245      e.category=v[j];
00246      e.histogram=hist;
00247      m_calibration[i]->data.push_back(e);
00248     }
00249 
00250   }
00251   if(categories == "RootXML")
00252   {
00253     bool resetHistogram = config.getParameter<bool>("resetHistograms");
00254     const TrackProbabilityCalibration * ca[2];
00255     ca[0]  = fromXml(config.getParameter<edm::FileInPath>("calibFile3d"));
00256     ca[1]  = fromXml(config.getParameter<edm::FileInPath>("calibFile2d"));
00257   
00258     for(unsigned int i=minLoop;i <=maxLoop ;i++)
00259      for(unsigned int j=0;j<ca[i]->data.size() ; j++)
00260      {
00261       TrackProbabilityCalibration::Entry e;
00262       e.category=ca[i]->data[j].category;
00263 
00264       if(resetHistogram)
00265        e.histogram=hist;
00266       else
00267        e.histogram=ca[i]->data[j].histogram;
00268 
00269       m_calibration[i]->data.push_back(e);
00270      }
00271 
00272    delete ca[0];
00273    delete ca[1];
00274 
00275    }
00276   if(categories == "EventSetup")
00277    {
00278     bool resetHistogram = config.getParameter<bool>("resetHistograms");
00279     ESHandle<TrackProbabilityCalibration> calib2DHandle;
00280     iSetup.get<BTagTrackProbability2DRcd>().get(calib2DHandle);
00281     ESHandle<TrackProbabilityCalibration> calib3DHandle;
00282     iSetup.get<BTagTrackProbability3DRcd>().get(calib3DHandle);
00283     const TrackProbabilityCalibration * ca[2];
00284     ca[0]  = calib3DHandle.product();
00285     ca[1]  = calib2DHandle.product();
00286     for(unsigned int i=minLoop;i <=maxLoop ;i++)
00287     for(unsigned int j=0;j<ca[i]->data.size() ; j++)
00288     {
00289      TrackProbabilityCalibration::Entry e;
00290      e.category=ca[i]->data[j].category;
00291 
00292      if(resetHistogram)
00293       e.histogram=hist;
00294      else 
00295       e.histogram=ca[i]->data[j].histogram;
00296 
00297      m_calibration[i]->data.push_back(e);
00298     }
00299 
00300   }
00301 
00302 
00303 
00304 /*  edm::FileInPath f2d(m_xmlfilename2D);
00305     edm::FileInPath f3d(m_xmlfilename3D);
00306     calibrationNew   =  new AlgorithmCalibration<TrackClassFilterCategory,CalibratedHistogramXML>((f3d.fullPath()).c_str());
00307     calibration2dNew =  new AlgorithmCalibration<TrackClassFilterCategory,CalibratedHistogramXML>((f2d.fullPath()).c_str());
00308     vector<float> * bins =0;
00309     if(m_resetData)
00310       {
00311         if(m_newBinning)  bins = new  vector<float>(CalibratedHistogram::constantBinning(m_nBin,0,m_range));
00312         vector<pair<TrackClassFilterCategory, CalibratedHistogramXML> > data = calibrationNew->categoriesWithData();
00313         vector<pair<TrackClassFilterCategory, CalibratedHistogramXML> > data2d = calibration2dNew->categoriesWithData();
00314         std::cout <<  data.size() <<  std::endl;
00315         for(unsigned int i = 0 ; i < data.size();i++)
00316           {
00317             data[i].second.reset();
00318             if(bins)  data[i].second.setUpperLimits(*bins);
00319           }
00320         for(unsigned int i = 0 ; i < data2d.size();i++)
00321           {
00322             data2d[i].second.reset();
00323             if(bins)  data2d[i].second.setUpperLimits(*bins);
00324           }
00325         
00326       }
00327     if(bins) delete bins;
00328     
00329 */
00330 
00331 
00332 }
00333 
00334 TrackProbabilityCalibration * ImpactParameterCalibration::fromXml(edm::FileInPath xmlCalibration)   
00335 {
00336      std::ifstream xmlFile(xmlCalibration.fullPath().c_str());
00337         if (!xmlFile.good())
00338                 throw cms::Exception("BTauFakeMVAJetTagConditions")
00339                         << "File \"" << xmlCalibration.fullPath()
00340                         << "\" could not be opened for reading."
00341                         << std::endl;
00342         std::ostringstream ss;
00343         ss << xmlFile.rdbuf();
00344         xmlFile.close();
00345         TClass *classType = 0;
00346         void *ptr = TBufferXML(TBuffer::kRead).ConvertFromXMLAny(
00347                                 ss.str().c_str(), &classType, kTRUE, kFALSE);
00348         if (!ptr)
00349                 throw cms::Exception("ImpactParameterCalibration")
00350                         << "Unknown error parsing XML serialization"
00351                         << std::endl;
00352 
00353         if (std::strcmp(classType->GetName(),
00354                 "TrackProbabilityCalibration")) {
00355                 classType->Destructor(ptr);
00356                 throw cms::Exception("ImpactParameterCalibration")
00357                         << "Serialized object has wrong C++ type."
00358                         << std::endl;
00359         }
00360 
00361         return static_cast<TrackProbabilityCalibration*>(ptr);
00362 }
00363 
00364 
00365 
00366 
00367 // ------------ method called once each job just after ending the event loop  ------------
00368 void 
00369 ImpactParameterCalibration::endJob() {
00370 
00371   if(config.getParameter<bool>("writeToDB"))
00372   {
00373     edm::Service<cond::service::PoolDBOutputService> mydbservice;
00374     if( !mydbservice.isAvailable() ) return;
00375     if(minLoop == 0 )    mydbservice->createNewIOV<TrackProbabilityCalibration>(m_calibration[0], mydbservice->beginOfTime(), mydbservice->endOfTime(),"BTagTrackProbability3DRcd");
00376     if(maxLoop == 1)   mydbservice->createNewIOV<TrackProbabilityCalibration>(m_calibration[1],  mydbservice->beginOfTime(), mydbservice->endOfTime(),"BTagTrackProbability2DRcd");
00377   } 
00378     
00379 
00380   if(config.getParameter<bool>("writeToRootXML"))
00381   {
00382     if(maxLoop == 1 ){
00383       std::ofstream of2("2d.xml");
00384       TBufferXML b2(TBuffer::kWrite);
00385       of2 << b2.ConvertToXML(const_cast<void*>(static_cast<const void*>(m_calibration[1])),
00386                              TClass::GetClass("TrackProbabilityCalibration"),
00387                              kTRUE, kFALSE);
00388       of2.close();
00389     }
00390     if(minLoop == 0 ){
00391       std::ofstream of3("3d.xml");
00392       TBufferXML b3(TBuffer::kWrite);
00393       of3 << b3.ConvertToXML(const_cast<void*>(static_cast<const void*>(m_calibration[0])),
00394                              TClass::GetClass("TrackProbabilityCalibration"),
00395                              kTRUE, kFALSE);
00396       of3.close();
00397     }
00398   }
00399 
00400  
00401   if(config.getParameter<bool>("writeToBinary"))
00402   {
00403     if(maxLoop == 1 ){
00404       
00405       std::ofstream ofile("2d.dat");
00406       MyTBuffer buffer(TBuffer::kWrite);
00407       buffer.StreamObject(const_cast<void*>(static_cast<const void*>(m_calibration[1])),
00408                           TClass::GetClass("TrackProbabilityCalibration"));
00409       Int_t size = buffer.Length();
00410       ofile.write(buffer.Buffer(),size);
00411       ofile.close();
00412     }
00413     if(minLoop == 0 ){
00414       std::ofstream ofile3("3d.dat");
00415       MyTBuffer buffer3(TBuffer::kWrite);
00416       buffer3.StreamObject(const_cast<void*>(static_cast<const void*>(m_calibration[0])),
00417                            TClass::GetClass("TrackProbabilityCalibration"));
00418       Int_t size3 = buffer3.Length();
00419       ofile3.write(buffer3.Buffer(),size3);
00420       ofile3.close();
00421     }
00422   }
00423   
00424   
00425 
00426 
00427 
00428 }
00429 
00430 //define this as a plug-in
00431 DEFINE_FWK_MODULE(ImpactParameterCalibration);