00001 #ifndef Validation_RecoMuon_Histograms_H
00002 #define Validation_RecoMuon_Histograms_H
00003
00012 #include "TH1F.h"
00013 #include "TH2F.h"
00014 #include "TString.h"
00015 #include "TFile.h"
00016
00017 #include "DQMServices/Core/interface/DQMStore.h"
00018 #include "DQMServices/Core/interface/MonitorElement.h"
00019 #include "FWCore/ServiceRegistry/interface/Service.h"
00020
00021 #include "DataFormats/GeometryVector/interface/Pi.h"
00022 #include <iostream>
00023 #include <vector>
00024 #include <math.h>
00025
00026 class HTrackVariables{
00027 public:
00028
00029 HTrackVariables(std::string dirName_, std::string name,std::string whereIs =""):theName(name.c_str()),where(whereIs.c_str()){
00030 dbe_ = edm::Service<DQMStore>().operator->();
00031 dbe_->cd();
00032 std::string dirName=dirName_;
00033
00034
00035
00036 dbe_->setCurrentFolder(dirName.c_str());
00037
00038 hEta = dbe_->book1D(theName+"_Eta_"+where,"#eta at "+where,120,-3.,3.);
00039 hPhi = dbe_->book1D(theName+"_Phi_"+where,"#phi at "+where,100,-Geom::pi(),Geom::pi());
00040 hP = dbe_->book1D(theName+"_P_"+where,"p_{t} at "+where,1000,0,2000);
00041 hPt = dbe_->book1D(theName+"_Pt_"+where,"p_{t} at "+where,1000,0,2000);
00042 hCharge = dbe_->book1D(theName+"_charge_"+where,"Charge at "+where,4,-2.,2.);
00043
00044 hEtaVsGen = dbe_->book1D(theName+"_EtaVsGen_"+where,"#eta at "+where,120,-3.,3.);
00045 hPhiVsGen = dbe_->book1D(theName+"_PhiVsGen_"+where,"#phi at "+where,100,-Geom::pi(),Geom::pi());
00046 hPtVsGen = dbe_->book1D(theName+"_PtVsGen_"+where,"p_{t} at "+where,1000,0,2000);
00047
00048 hDeltaR = dbe_->book1D(theName+"_DeltaR_"+where,"Delta R w.r.t. sim track for "+where,1000,0,20);
00049
00050 theEntries = 0;
00051 }
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064 ~HTrackVariables(){}
00065
00066 MonitorElement *eta() {return hEta;}
00067 MonitorElement *phi() {return hPhi;}
00068 MonitorElement *p() {return hP;}
00069 MonitorElement *pt() {return hPt;}
00070 MonitorElement *charge() {return hCharge;}
00071 int entries() {return theEntries;}
00072
00073 void Fill(double p, double pt, double eta, double phi, double charge){
00074 hEta->Fill(eta);
00075 hPhi->Fill(phi);
00076 hP->Fill(p);
00077 hPt->Fill(pt);
00078 hCharge->Fill(charge);
00079 ++theEntries;
00080 }
00081
00082 void Fill(double pt, double eta, double phi){
00083 hEtaVsGen->Fill(eta);
00084 hPhiVsGen->Fill(phi);
00085 hPtVsGen->Fill(pt);
00086 }
00087
00088 void FillDeltaR(double deltaR){
00089 hDeltaR->Fill(deltaR);
00090 }
00091
00092 double computeEfficiency(HTrackVariables *sim){
00093
00094 efficiencyHistos.push_back(computeEfficiency(hEtaVsGen,sim->eta()));
00095 efficiencyHistos.push_back(computeEfficiency(hPhiVsGen,sim->phi()));
00096
00097 efficiencyHistos.push_back(computeEfficiency(hPtVsGen,sim->pt()));
00098
00099
00100 double efficiency = 100*entries()/sim->entries();
00101 return efficiency;
00102 }
00103
00104 MonitorElement* computeEfficiency(MonitorElement *reco, MonitorElement *sim){
00105
00106 TH1F* hReco = reco->getTH1F();
00107 TH1F* hSim = sim->getTH1F();
00108
00109 std::string name = hReco->GetName();
00110 std::string title = hReco->GetTitle();
00111
00112 MonitorElement * me = dbe_->book1D(
00113 "Eff_"+name,
00114 "Efficiecny as a function of "+title,
00115 hSim->GetNbinsX(),
00116 hSim->GetXaxis()->GetXmin(),
00117 hSim->GetXaxis()->GetXmax()
00118 );
00119
00120 me->getTH1F()->Divide(hReco,hSim,1.,1.,"b");
00121
00122
00123 int nBinsEta = me->getTH1F()->GetNbinsX();
00124 for(int bin = 1; bin <= nBinsEta; bin++) {
00125 float nSimHit = hSim->GetBinContent(bin);
00126 float eff = me->getTH1F()->GetBinContent(bin);
00127 float error = 0;
00128 if(nSimHit != 0 && eff <= 1) {
00129 error = sqrt(eff*(1-eff)/nSimHit);
00130 }
00131 me->getTH1F()->SetBinError(bin, error);
00132 }
00133
00134 return me;
00135 }
00136
00137
00138 private:
00139 DQMStore* dbe_;
00140
00141 std::string theName;
00142 std::string where;
00143
00144 int theEntries;
00145
00146 MonitorElement *hEta;
00147 MonitorElement *hPhi;
00148 MonitorElement *hP;
00149 MonitorElement *hPt;
00150 MonitorElement *hCharge;
00151
00152 MonitorElement *hEtaVsGen;
00153 MonitorElement *hPhiVsGen;
00154 MonitorElement *hPtVsGen;
00155
00156 MonitorElement* hDeltaR;
00157
00158 std::vector<MonitorElement*> efficiencyHistos;
00159
00160
00161 };
00162
00163
00164 class HResolution {
00165 public:
00166
00167 HResolution(std::string dirName_,std::string name,std::string whereIs):theName(name.c_str()),where(whereIs.c_str()){
00168
00169 dbe_ = edm::Service<DQMStore>().operator->();
00170 dbe_->cd();
00171 std::string dirName=dirName_;
00172
00173
00174
00175 dbe_->setCurrentFolder(dirName.c_str());
00176
00177 double eta = 15.; int neta = 800;
00178 double phi = 12.; int nphi = 400;
00179 double pt = 60.; int npt = 2000;
00180
00181 hEta = dbe_->book1D(theName+"_Eta_"+where,"#eta "+theName,neta,-eta,eta);
00182 hPhi = dbe_->book1D(theName+"_Phi_"+where,"#phi "+theName,nphi,-phi,phi);
00183
00184 hP = dbe_->book1D(theName+"_P_"+where,"P "+theName,400,-4,4);
00185 hPt = dbe_->book1D(theName+"_Pt_"+where,"P_{t} "+theName,npt,-pt,pt);
00186
00187 hCharge = dbe_->book1D(theName+"_charge_"+where,"Charge "+theName,4,-2.,2.);
00188
00189
00190 h2Eta = dbe_->book2D(theName+"_Eta_vs_Eta"+where,"#eta "+theName+" as a function of #eta",200,-2.5,2.5,neta,-eta,eta);
00191 h2Phi = dbe_->book2D(theName+"_Phi_vs_Phi"+where,"#phi "+theName+" as a function of #phi",100,-Geom::pi(),Geom::pi(),nphi,-phi,phi);
00192
00193 h2P = dbe_->book2D(theName+"_P_vs_P"+where,"P "+theName+" as a function of P",1000,0,2000,400,-4,4);
00194 h2Pt = dbe_->book2D(theName+"_Pt_vs_Pt"+where,"P_{t} "+theName+" as a function of P_{t}",1000,0,2000,npt,-pt,pt);
00195
00196 h2PtVsEta = dbe_->book2D(theName+"_Pt_vs_Eta"+where,"P_{t} "+theName+" as a function of #eta",200,-2.5,2.5,npt,-pt,pt);
00197 h2PtVsPhi = dbe_->book2D(theName+"_Pt_vs_Phi"+where,"P_{t} "+theName+" as a function of #phi",100,-Geom::pi(),Geom::pi(),npt,-pt,pt);
00198
00199 h2EtaVsPt = dbe_->book2D(theName+"_Eta_vs_Pt"+where,"#eta "+theName+" as a function of P_{t}",1000,0,2000,neta,-eta,eta);
00200 h2EtaVsPhi = dbe_->book2D(theName+"_Eta_vs_Phi"+where,"#eta "+theName+" as a function of #phi",100,-Geom::pi(),Geom::pi(),neta,-eta,eta);
00201
00202 h2PhiVsPt = dbe_->book2D(theName+"_Phi_vs_Pt"+where,"#phi "+theName+" as a function of P_{t}",1000,0,2000,nphi,-phi,phi);
00203 h2PhiVsEta = dbe_->book2D(theName+"_Phi_vs_Eta"+where,"#phi "+theName+" as a function of #eta",200,-2.5,2.5,nphi,-phi,phi);
00204 }
00205
00206 HResolution(std::string name, TFile* file):theName(name.c_str()){
00207
00208 }
00209
00210 ~HResolution(){}
00211
00212
00213 void Fill(double p, double pt, double eta, double phi,
00214 double rp, double rpt,double reta, double rphi, double rcharge){
00215
00216 Fill(rp, rpt, reta, rphi, rcharge);
00217
00218
00219 h2Eta->Fill(eta,reta);
00220 h2Phi->Fill(phi,rphi);
00221
00222 h2P->Fill(p,rp);
00223 h2Pt->Fill(pt,rpt);
00224
00225 h2PtVsEta->Fill(eta,rpt);
00226 h2PtVsPhi->Fill(phi,rpt);
00227
00228 h2EtaVsPt ->Fill(pt,reta);
00229 h2EtaVsPhi->Fill(phi,reta);
00230
00231 h2PhiVsPt ->Fill(pt,rphi);
00232 h2PhiVsEta->Fill(eta,rphi);
00233 }
00234
00235
00236 void Fill(double p, double pt, double eta, double phi,
00237 double rp, double rpt){
00238
00239 hP->Fill(rp);
00240 hPt->Fill(rpt);
00241
00242 h2P->Fill(p,rp);
00243
00244
00245
00246 h2Pt->Fill(pt,rpt);
00247 h2PtVsEta->Fill(eta,rpt);
00248 h2PtVsPhi->Fill(phi,rpt);
00249 }
00250
00251 void Fill(double rp, double rpt,
00252 double reta, double rphi, double rcharge){
00253
00254 hEta->Fill(reta);
00255 hPhi->Fill(rphi);
00256
00257 hP->Fill(rp);
00258 hPt->Fill(rpt);
00259
00260 hCharge->Fill(rcharge);
00261 }
00262
00263
00264
00265 private:
00266 DQMStore* dbe_;
00267
00268 std::string theName;
00269 std::string where;
00270
00271 MonitorElement *hEta;
00272 MonitorElement *hPhi;
00273
00274 MonitorElement *hP;
00275 MonitorElement *hPt;
00276
00277 MonitorElement *hCharge;
00278
00279 MonitorElement *h2Eta;
00280 MonitorElement *h2Phi;
00281
00282 MonitorElement *h2P;
00283 MonitorElement *h2Pt;
00284
00285 MonitorElement *h2PtVsEta;
00286 MonitorElement *h2PtVsPhi;
00287
00288 MonitorElement *h2EtaVsPt;
00289 MonitorElement *h2EtaVsPhi;
00290
00291 MonitorElement *h2PhiVsPt;
00292 MonitorElement *h2PhiVsEta;
00293
00294 };
00295
00296
00297 class HResolution1DRecHit{
00298 public:
00299 HResolution1DRecHit(std::string name):theName(name.c_str()){
00300
00301
00302 hResX = dbe_->book1D (theName+"_X_Res", "X residual", 5000, -0.5,0.5);
00303 hResY = dbe_->book1D (theName+"_Y_Res", "Y residual", 5000, -1.,1.);
00304 hResZ = dbe_->book1D (theName+"_Z_Res", "Z residual", 5000, -1.5,1.5);
00305
00306 hResXVsEta = dbe_->book2D(theName+"_XResVsEta", "X residual vs eta",
00307 200, -2.5,2.5, 5000, -1.5,1.5);
00308 hResYVsEta = dbe_->book2D(theName+"_YResVsEta", "Y residual vs eta",
00309 200, -2.5,2.5, 5000, -1.,1.);
00310 hResZVsEta = dbe_->book2D(theName+"_ZResVsEta", "Z residual vs eta",
00311 200, -2.5,2.5, 5000, -1.5,1.5);
00312
00313 hXResVsPhi = dbe_->book2D(theName+"_XResVsPhi", "X residual vs phi",
00314 100,-Geom::pi(),Geom::pi(), 5000, -0.5,0.5);
00315 hYResVsPhi = dbe_->book2D(theName+"_YResVsPhi", "Y residual vs phi",
00316 100,-Geom::pi(),Geom::pi(), 5000, -1.,1.);
00317 hZResVsPhi = dbe_->book2D(theName+"_ZResVsPhi", "Z residual vs phi",
00318 100,-Geom::pi(),Geom::pi(), 5000, -1.5,1.5);
00319
00320 hXResVsPos = dbe_->book2D(theName+"_XResVsPos", "X residual vs position",
00321 10000, -750,750, 5000, -0.5,0.5);
00322 hYResVsPos = dbe_->book2D(theName+"_YResVsPos", "Y residual vs position",
00323 10000, -740,740, 5000, -1.,1.);
00324 hZResVsPos = dbe_->book2D(theName+"_ZResVsPos", "Z residual vs position",
00325 10000, -1100,1100, 5000, -1.5,1.5);
00326
00327 hXPull = dbe_->book1D (theName+"_XPull", "X pull", 600, -2,2);
00328 hYPull = dbe_->book1D (theName+"_YPull", "Y pull", 600, -2,2);
00329 hZPull = dbe_->book1D (theName+"_ZPull", "Z pull", 600, -2,2);
00330
00331 hXPullVsPos = dbe_->book2D (theName+"_XPullVsPos", "X pull vs position",10000, -750,750, 600, -2,2);
00332 hYPullVsPos = dbe_->book2D (theName+"_YPullVsPos", "Y pull vs position",10000, -740,740, 600, -2,2);
00333 hZPullVsPos = dbe_->book2D (theName+"_ZPullVsPos", "Z pull vs position",10000, -1100,1100, 600, -2,2);
00334 }
00335
00336 HResolution1DRecHit(TString name_, TFile* file){}
00337
00338 ~HResolution1DRecHit(){}
00339
00340 void Fill(double x, double y, double z,
00341 double dx, double dy, double dz,
00342 double errx, double erry, double errz,
00343 double eta, double phi) {
00344
00345 double rx = dx/x, ry = dy/y, rz = dz/z;
00346
00347 hResX->Fill(rx);
00348 hResY->Fill(ry);
00349 hResZ->Fill(rz);
00350
00351 hResXVsEta->Fill(eta,rx);
00352 hResYVsEta->Fill(eta,ry);
00353 hResZVsEta->Fill(eta,rz);
00354
00355 hXResVsPhi->Fill(phi,rx);
00356 hYResVsPhi->Fill(phi,ry);
00357 hZResVsPhi->Fill(phi,rz);
00358
00359 hXResVsPos->Fill(x,rx);
00360 hYResVsPos->Fill(y,ry);
00361 hZResVsPos->Fill(z,rz);
00362
00363 if(errx < 1e-6)
00364 std::cout << "NO proper error set for X: "<<errx<<std::endl;
00365 else{
00366 hXPull->Fill(dx/errx);
00367 hXPullVsPos->Fill(x,dx/errx);
00368 }
00369
00370 if(erry < 1e-6)
00371 std::cout << "NO proper error set for Y: "<<erry<<std::endl;
00372 else{
00373 hYPull->Fill(dy/erry);
00374 hYPullVsPos->Fill(y,dy/erry);
00375 }
00376 if(errz < 1e-6)
00377 std::cout << "NO proper error set for Z: "<<errz<<std::endl;
00378 else{
00379 hZPull->Fill(dz/errz);
00380 hZPullVsPos->Fill(z,dz/errz);
00381 }
00382 }
00383
00384
00385
00386 public:
00387 std::string theName;
00388
00389 MonitorElement *hResX;
00390 MonitorElement *hResY;
00391 MonitorElement *hResZ;
00392
00393 MonitorElement *hResXVsEta;
00394 MonitorElement *hResYVsEta;
00395 MonitorElement *hResZVsEta;
00396
00397 MonitorElement *hXResVsPhi;
00398 MonitorElement *hYResVsPhi;
00399 MonitorElement *hZResVsPhi;
00400
00401 MonitorElement *hXResVsPos;
00402 MonitorElement *hYResVsPos;
00403 MonitorElement *hZResVsPos;
00404
00405 MonitorElement *hXPull;
00406 MonitorElement *hYPull;
00407 MonitorElement *hZPull;
00408
00409 MonitorElement *hXPullVsPos;
00410 MonitorElement *hYPullVsPos;
00411 MonitorElement *hZPullVsPos;
00412
00413 private:
00414 DQMStore* dbe_;
00415 };
00416 #endif
00417