CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_5/src/Validation/RecoVertex/src/BSvsPVHistogramMaker.cc

Go to the documentation of this file.
00001 #include "Validation/RecoVertex/interface/BSvsPVHistogramMaker.h"
00002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00003 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00004 #include "FWCore/Framework/interface/Event.h"
00005 #include "FWCore/ServiceRegistry/interface/Service.h"
00006 #include "CommonTools/UtilAlgos/interface/TFileService.h"
00007 #include "DataFormats/VertexReco/interface/Vertex.h"
00008 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
00009 #include "TH1F.h"
00010 #include "TH2F.h"
00011 #include "TProfile.h"
00012 
00013 
00014 BSvsPVHistogramMaker::BSvsPVHistogramMaker():
00015   _currdir(0), m_maxLS(100), useSlope_(true), _trueOnly(true), 
00016   _runHisto(true), _runHistoProfile(true), _runHistoBXProfile(true), _runHistoBX2D(false), _histoParameters() { }
00017 
00018 BSvsPVHistogramMaker::BSvsPVHistogramMaker(const edm::ParameterSet& iConfig):
00019   _currdir(0),
00020   m_maxLS(iConfig.getParameter<unsigned int>("maxLSBeforeRebin")),
00021   useSlope_(iConfig.getParameter<bool>("useSlope")),
00022   _trueOnly(iConfig.getUntrackedParameter<bool>("trueOnly",true)),
00023   _runHisto(iConfig.getUntrackedParameter<bool>("runHisto",true)),
00024   _runHistoProfile(iConfig.getUntrackedParameter<bool>("runHistoProfile",true)),
00025   _runHistoBXProfile(iConfig.getUntrackedParameter<bool>("runHistoBXProfile",true)),
00026   _runHistoBX2D(iConfig.getUntrackedParameter<bool>("runHistoBX2D",false)),
00027   _histoParameters(iConfig.getUntrackedParameter<edm::ParameterSet>("histoParameters",edm::ParameterSet())),
00028   _rhm()
00029 { }
00030 
00031 
00032 BSvsPVHistogramMaker::~BSvsPVHistogramMaker() {
00033 
00034   delete _currdir;
00035 
00036 }
00037 
00038 
00039 void BSvsPVHistogramMaker::book(const std::string dirname) {
00040 
00041   edm::Service<TFileService> tfserv;
00042   TFileDirectory* currdir = &(*tfserv);
00043 
00044   if(dirname!="") {
00045     currdir = new TFileDirectory(tfserv->mkdir(dirname));
00046     _currdir = currdir;
00047   }
00048 
00049   edm::LogInfo("HistogramBooking") << "Vertex histogram booking in directory " << dirname;
00050 
00051   _hdeltax = currdir->make<TH1F>("deltax","(PV-BS) X position",
00052                                _histoParameters.getUntrackedParameter<unsigned int>("nBinX",200),
00053                                _histoParameters.getUntrackedParameter<double>("xMin",-1.),
00054                                _histoParameters.getUntrackedParameter<double>("xMax",1.)
00055                                );
00056   _hdeltax->GetXaxis()->SetTitle("#Delta(X) [cm]");   _hdeltax->GetYaxis()->SetTitle("Vertices"); 
00057 
00058   _hdeltay = currdir->make<TH1F>("deltay","(PV-BS) Y position",
00059                                _histoParameters.getUntrackedParameter<unsigned int>("nBinY",200),
00060                                _histoParameters.getUntrackedParameter<double>("yMin",-1.),
00061                                _histoParameters.getUntrackedParameter<double>("yMax",1.)
00062                                );
00063   _hdeltay->GetXaxis()->SetTitle("#Delta(Y) [cm]");   _hdeltay->GetYaxis()->SetTitle("Vertices"); 
00064 
00065   _hdeltaz = currdir->make<TH1F>("deltaz","(PV-BS) Z position",
00066                                _histoParameters.getUntrackedParameter<unsigned int>("nBinZ",200),
00067                                _histoParameters.getUntrackedParameter<double>("zMin",-20.),
00068                                _histoParameters.getUntrackedParameter<double>("zMax",20.)
00069                                );
00070   _hdeltaz->GetXaxis()->SetTitle("#Delta(Z) [cm]");   _hdeltaz->GetYaxis()->SetTitle("Vertices"); 
00071 
00072   _hdeltaxvsz = currdir->make<TProfile>("deltaxvsz","(PV-BS) X position vs Z",
00073                                _histoParameters.getUntrackedParameter<unsigned int>("nBinZProfile",40),
00074                                _histoParameters.getUntrackedParameter<double>("zMinProfile",-20.),
00075                                _histoParameters.getUntrackedParameter<double>("zMaxProfile",20.)
00076                                );
00077   _hdeltaxvsz->GetXaxis()->SetTitle("Z [cm]");   _hdeltaxvsz->GetYaxis()->SetTitle("#Delta(X) [cm]"); 
00078 
00079   _hdeltayvsz = currdir->make<TProfile>("deltayvsz","(PV-BS) Y position vs Z",
00080                                _histoParameters.getUntrackedParameter<unsigned int>("nBinZProfile",40),
00081                                _histoParameters.getUntrackedParameter<double>("zMinProfile",-20.),
00082                                _histoParameters.getUntrackedParameter<double>("zMaxProfile",20.)
00083                                );
00084   _hdeltayvsz->GetXaxis()->SetTitle("Z [cm]");   _hdeltayvsz->GetYaxis()->SetTitle("#Delta(Y) [cm]"); 
00085 
00086 
00087 
00088 
00089   if(_runHisto) {
00090     _hdeltaxrun = _rhm.makeTH1F("deltaxrun","(PV-BS) X position",
00091                               _histoParameters.getUntrackedParameter<unsigned int>("nBinX",200),
00092                               _histoParameters.getUntrackedParameter<double>("xMin",-1.),
00093                               _histoParameters.getUntrackedParameter<double>("xMax",1.));
00094     
00095     _hdeltayrun = _rhm.makeTH1F("deltayrun","(PV-BS) Y position",
00096                               _histoParameters.getUntrackedParameter<unsigned int>("nBinY",200),
00097                               _histoParameters.getUntrackedParameter<double>("yMin",-1.),
00098                               _histoParameters.getUntrackedParameter<double>("yMax",1.));
00099     
00100     _hdeltazrun = _rhm.makeTH1F("deltazrun","(PV-BS) Z position",
00101                               _histoParameters.getUntrackedParameter<unsigned int>("nBinZ",200),
00102                               _histoParameters.getUntrackedParameter<double>("zMin",-20.),
00103                               _histoParameters.getUntrackedParameter<double>("zMax",20.));
00104     
00105     _hdeltaxvszrun = _rhm.makeTProfile("deltaxvszrun","(PV-BS) X position vs Z",
00106                                        _histoParameters.getUntrackedParameter<unsigned int>("nBinZProfile",40),
00107                                        _histoParameters.getUntrackedParameter<double>("zMinProfile",-20.),
00108                                        _histoParameters.getUntrackedParameter<double>("zMaxProfile",20.)
00109                                        );
00110 
00111     _hdeltayvszrun = _rhm.makeTProfile("deltayvszrun","(PV-BS) Y position vs Z",
00112                                        _histoParameters.getUntrackedParameter<unsigned int>("nBinZProfile",40),
00113                                        _histoParameters.getUntrackedParameter<double>("zMinProfile",-20.),
00114                                        _histoParameters.getUntrackedParameter<double>("zMaxProfile",20.)
00115                                        );
00116 
00117     if(_runHistoProfile) {
00118       _hdeltaxvsorbrun = _rhm.makeTProfile("deltaxvsorbrun","(PV-BS) X position vs orbit number",4*m_maxLS,0.5,m_maxLS*262144+0.5);
00119       _hdeltayvsorbrun = _rhm.makeTProfile("deltayvsorbrun","(PV-BS) Y position vs orbit number",4*m_maxLS,0.5,m_maxLS*262144+0.5);
00120       _hdeltazvsorbrun = _rhm.makeTProfile("deltazvsorbrun","(PV-BS) Z position vs orbit number",4*m_maxLS,0.5,m_maxLS*262144+0.5);
00121     }
00122     if(_runHistoBXProfile) {
00123       _hdeltaxvsbxrun = _rhm.makeTProfile("deltaxvsbxrun","(PV-BS) X position vs BX number",3564,-0.5,3563.5);
00124       _hdeltayvsbxrun = _rhm.makeTProfile("deltayvsbxrun","(PV-BS) Y position vs BX number",3564,-0.5,3563.5);
00125       _hdeltazvsbxrun = _rhm.makeTProfile("deltazvsbxrun","(PV-BS) Z position vs BX number",3564,-0.5,3563.5);
00126       if(_runHistoBX2D) {
00127         _hdeltaxvsbx2drun = _rhm.makeTH2F("deltaxvsbx2drun","(PV-BS) X position vs BX number",3564,-0.5,3563.5,
00128                                           _histoParameters.getUntrackedParameter<unsigned int>("nBinX",200),
00129                                           _histoParameters.getUntrackedParameter<double>("xMin",-1.),
00130                                           _histoParameters.getUntrackedParameter<double>("xMax",1.));
00131         _hdeltayvsbx2drun = _rhm.makeTH2F("deltayvsbx2drun","(PV-BS) Y position vs BX number",3564,-0.5,3563.5,
00132                                           _histoParameters.getUntrackedParameter<unsigned int>("nBinY",200),
00133                                           _histoParameters.getUntrackedParameter<double>("yMin",-1.),
00134                                           _histoParameters.getUntrackedParameter<double>("yMax",1.));
00135         _hdeltazvsbx2drun = _rhm.makeTH2F("deltazvsbx2drun","(PV-BS) Z position vs BX number",3564,-0.5,3563.5,
00136                                           _histoParameters.getUntrackedParameter<unsigned int>("nBinZ",200),
00137                                           _histoParameters.getUntrackedParameter<double>("zMin",-20.),
00138                                           _histoParameters.getUntrackedParameter<double>("zMax",20.));
00139       }
00140     }
00141     
00142   }
00143 }
00144 
00145 void BSvsPVHistogramMaker::beginRun(const unsigned int nrun) {
00146 
00147   char runname[100];
00148   sprintf(runname,"run_%d",nrun);
00149 
00150   TFileDirectory* currdir = _currdir;
00151   if(currdir==0) {
00152     edm::Service<TFileService> tfserv;
00153     currdir = &(*tfserv);
00154   }
00155 
00156   _rhm.beginRun(nrun,*currdir);
00157 
00158   if(_runHisto) {
00159     (*_hdeltaxrun)->GetXaxis()->SetTitle("#Delta(X) [cm]");   (*_hdeltaxrun)->GetYaxis()->SetTitle("Vertices"); 
00160     (*_hdeltayrun)->GetXaxis()->SetTitle("#Delta(Y) [cm]");   (*_hdeltayrun)->GetYaxis()->SetTitle("Vertices"); 
00161     (*_hdeltazrun)->GetXaxis()->SetTitle("#Delta(Z) [cm]");   (*_hdeltazrun)->GetYaxis()->SetTitle("Vertices"); 
00162     (*_hdeltaxvszrun)->GetXaxis()->SetTitle("Z [cm]");   (*_hdeltaxvszrun)->GetYaxis()->SetTitle("#Delta(X) [cm]"); 
00163     (*_hdeltayvszrun)->GetXaxis()->SetTitle("Z [cm]");   (*_hdeltayvszrun)->GetYaxis()->SetTitle("#Delta(Y) [cm]"); 
00164     
00165     if(_runHistoProfile) {
00166       (*_hdeltaxvsorbrun)->GetXaxis()->SetTitle("time [orbit#]");   (*_hdeltaxvsorbrun)->GetYaxis()->SetTitle("#Delta(X) [cm]"); 
00167       (*_hdeltaxvsorbrun)->SetBit(TH1::kCanRebin);
00168       (*_hdeltayvsorbrun)->GetXaxis()->SetTitle("time [orbit#]");   (*_hdeltayvsorbrun)->GetYaxis()->SetTitle("#Delta(Y) [cm]"); 
00169       (*_hdeltayvsorbrun)->SetBit(TH1::kCanRebin);
00170       (*_hdeltazvsorbrun)->GetXaxis()->SetTitle("time [orbit#]");   (*_hdeltazvsorbrun)->GetYaxis()->SetTitle("#Delta(Z) [cm]"); 
00171       (*_hdeltazvsorbrun)->SetBit(TH1::kCanRebin);
00172     }
00173     if(_runHistoBXProfile) {
00174       (*_hdeltaxvsbxrun)->GetXaxis()->SetTitle("BX");   (*_hdeltaxvsbxrun)->GetYaxis()->SetTitle("#Delta(X) [cm]"); 
00175       (*_hdeltayvsbxrun)->GetXaxis()->SetTitle("BX");   (*_hdeltayvsbxrun)->GetYaxis()->SetTitle("#Delta(Y) [cm]"); 
00176       (*_hdeltazvsbxrun)->GetXaxis()->SetTitle("BX");   (*_hdeltazvsbxrun)->GetYaxis()->SetTitle("#Delta(Z) [cm]"); 
00177       if(_runHistoBX2D) {
00178         (*_hdeltaxvsbx2drun)->GetXaxis()->SetTitle("BX");   (*_hdeltaxvsbx2drun)->GetYaxis()->SetTitle("#Delta(X) [cm]"); 
00179         (*_hdeltayvsbx2drun)->GetXaxis()->SetTitle("BX");   (*_hdeltayvsbx2drun)->GetYaxis()->SetTitle("#Delta(Y) [cm]"); 
00180         (*_hdeltazvsbx2drun)->GetXaxis()->SetTitle("BX");   (*_hdeltazvsbx2drun)->GetYaxis()->SetTitle("#Delta(Z) [cm]"); 
00181       }
00182     }
00183     
00184   }
00185 }
00186 
00187 void BSvsPVHistogramMaker::fill(const unsigned int orbit, const int bx, const reco::VertexCollection& vertices, const reco::BeamSpot& bs) {
00188   
00189   for(reco::VertexCollection::const_iterator vtx=vertices.begin();vtx!=vertices.end();++vtx) {
00190 
00191     if(!(_trueOnly && vtx->isFake())) {
00192       
00193       /*
00194       double deltax = vtx->x()-bs.x0();
00195       double deltay = vtx->y()-bs.y0();
00196       double deltaz = vtx->z()-bs.z0();
00197       */
00198       double deltax = vtx->x()-x(bs,vtx->z());
00199       double deltay = vtx->y()-y(bs,vtx->z());
00200       double deltaz = vtx->z()-bs.z0();
00201 
00202       _hdeltax->Fill(deltax);
00203       _hdeltay->Fill(deltay);
00204       _hdeltaz->Fill(deltaz);
00205       _hdeltaxvsz->Fill(vtx->z(),deltax);
00206       _hdeltayvsz->Fill(vtx->z(),deltay);
00207 
00208       if(_runHisto) {
00209         if(_hdeltaxrun && *_hdeltaxrun )  (*_hdeltaxrun)->Fill(deltax);
00210         if(_hdeltayrun && *_hdeltayrun )  (*_hdeltayrun)->Fill(deltay);
00211         if(_hdeltazrun && *_hdeltazrun )  (*_hdeltazrun)->Fill(deltaz);
00212         if(_hdeltaxvszrun && *_hdeltaxvszrun )  (*_hdeltaxvszrun)->Fill(vtx->z(),deltax);
00213         if(_hdeltayvszrun && *_hdeltayvszrun )  (*_hdeltayvszrun)->Fill(vtx->z(),deltay);
00214         if(_runHistoProfile) {
00215           if(_hdeltaxvsorbrun && *_hdeltaxvsorbrun )  (*_hdeltaxvsorbrun)->Fill(orbit,deltax);
00216           if(_hdeltayvsorbrun && *_hdeltayvsorbrun )  (*_hdeltayvsorbrun)->Fill(orbit,deltay);
00217           if(_hdeltazvsorbrun && *_hdeltazvsorbrun )  (*_hdeltazvsorbrun)->Fill(orbit,deltaz);
00218         }
00219         if(_runHistoBXProfile) {
00220           if(_hdeltaxvsbxrun && *_hdeltaxvsbxrun )  (*_hdeltaxvsbxrun)->Fill(bx,deltax);
00221           if(_hdeltayvsbxrun && *_hdeltayvsbxrun )  (*_hdeltayvsbxrun)->Fill(bx,deltay);
00222           if(_hdeltazvsbxrun && *_hdeltazvsbxrun )  (*_hdeltazvsbxrun)->Fill(bx,deltaz);
00223           if(_runHistoBX2D) {
00224             if(_hdeltaxvsbx2drun && *_hdeltaxvsbx2drun )  (*_hdeltaxvsbx2drun)->Fill(bx,deltax);
00225             if(_hdeltayvsbx2drun && *_hdeltayvsbx2drun )  (*_hdeltayvsbx2drun)->Fill(bx,deltay);
00226             if(_hdeltazvsbx2drun && *_hdeltazvsbx2drun )  (*_hdeltazvsbx2drun)->Fill(bx,deltaz);
00227           }
00228         }
00229       }
00230     }
00231   }
00232 }
00233 
00234 void BSvsPVHistogramMaker::fill(const edm::Event& iEvent, const reco::VertexCollection& vertices, const reco::BeamSpot& bs) {
00235 
00236   fill(iEvent.orbitNumber(),iEvent.bunchCrossing(),vertices,bs);
00237 
00238 }
00239 
00240 double BSvsPVHistogramMaker::x(const reco::BeamSpot& bs, const double z) const {
00241 
00242   double x = bs.x0();
00243   
00244   //  if(useSlope_) x += bs.dxdz()*z;
00245   if(useSlope_) x += bs.dxdz()*(z-bs.z0());
00246 
00247   //  if(useSlope_) x = bs.x(z);
00248  
00249   return x;
00250 
00251 }
00252 
00253 double BSvsPVHistogramMaker::y(const reco::BeamSpot& bs, const double z) const {
00254 
00255   double y = bs.y0();
00256   
00257   //  if(useSlope_) y += bs.dydz()*z;
00258   if(useSlope_) y += bs.dydz()*(z-bs.z0());
00259 
00260     //  if(useSlope_) y = bs.y(z);
00261  
00262   return y;
00263 
00264 }
00265