Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-06 08:16:17

0001 #include "../CommonTools.h"
0002 #include <sPhenixStyle.C>
0003 
0004 
0005 void VertexingQA(std::string reffile, std::string newfile, std::string outfile)
0006 {
0007   SetsPhenixStyle();
0008   TFile *qa_file_new = TFile::Open(newfile.c_str());
0009   TFile *qa_file_ref = TFile::Open(reffile.c_str());
0010 
0011   TFile *outfilef = new TFile(outfile.c_str(), "recreate");
0012 
0013 
0014   {
0015     const char *hist_name_prefix = "QAG4SimulationTracking";
0016     TString prefix = TString("h_") + hist_name_prefix + TString("_");
0017 
0018     
0019   // obtain normalization
0020   double Nevent_new = 1;
0021   double Nevent_ref = 1;
0022     
0023   TH2 *h_new = (TH2 *) qa_file_new->GetObjectChecked(
0024       prefix + TString("DCArPhi_pT"), "TH2");
0025   assert(h_new);
0026 
0027   //  h_new->Rebin(1, 2);
0028   //h_new->Sumw2();
0029   //  h_new->Scale(1. / Nevent_new);
0030 
0031   TH2 *h_ref = NULL;
0032   if (qa_file_ref)
0033   {
0034     h_ref = (TH2 *) qa_file_ref->GetObjectChecked(
0035         prefix + TString("DCArPhi_pT"), "TH2");
0036     assert(h_ref);
0037 
0038     //    h_ref->Rebin(1, 2);
0039     //h_ref->Sumw2();
0040     h_ref->Scale(Nevent_new / Nevent_ref);
0041   }
0042 
0043   TCanvas *c1 = new TCanvas(TString("QA_Draw_Tracking_DCArPhi") + TString("_") + hist_name_prefix,
0044                             TString("QA_Draw_Tracking_DCArPhi") + TString("_") + hist_name_prefix,
0045                             1800, 1000);
0046   c1->Divide(4, 2);
0047   int idx = 1;
0048   TPad *p;
0049 
0050   vector<pair<double, double>> gpt_ranges{
0051       {0, 0.5},
0052       {0.5, 1},
0053       {1, 1.5},
0054       {1.5, 2},
0055       {2, 4},
0056       {4, 16},
0057       {16, 40}};
0058 
0059   TF1 *f1 = nullptr;
0060   TF1 *fit = nullptr;
0061   Double_t sigma = 0;
0062   Double_t sigma_unc = 0;
0063   char resstr[500];
0064   TLatex *res = nullptr;
0065   for (auto pt_range : gpt_ranges)
0066   {
0067     //cout << __PRETTY_FUNCTION__ << " process " << pt_range.first << " - " << pt_range.second << " GeV/c";
0068 
0069     p = (TPad *) c1->cd(idx++);
0070     c1->Update();
0071     p->SetLogy();
0072 
0073     const double epsilon = 1e-6;
0074     const int bin_start = h_new->GetXaxis()->FindBin(pt_range.first + epsilon);
0075     const int bin_end = h_new->GetXaxis()->FindBin(pt_range.second - epsilon);
0076 
0077     TH1 *h_proj_new = h_new->ProjectionY(
0078         TString::Format(
0079             "%s_New_ProjX_%d_%d",
0080             h_new->GetName(), bin_start, bin_end),
0081         bin_start, bin_end);
0082     if (pt_range.first < 2.0)
0083     {
0084       h_proj_new->GetXaxis()->SetRangeUser(-.05,.05);
0085       h_proj_new->Rebin(5);
0086     }
0087     else
0088     {
0089       h_proj_new->GetXaxis()->SetRangeUser(-.01,.01);
0090     }
0091     h_proj_new->SetTitle(TString(hist_name_prefix) + TString::Format(
0092                                                          ": %.1f - %.1f GeV/c", pt_range.first, pt_range.second));
0093     h_proj_new->GetXaxis()->SetTitle(TString::Format(
0094         "DCA (r #phi) [cm]"));
0095     h_proj_new->GetXaxis()->SetNdivisions(5,5);
0096     f1 = new TF1("f1","gaus",-.01,.01);
0097     h_proj_new->Fit(f1,"qm");
0098     sigma = f1->GetParameter(2);
0099     sigma_unc = f1->GetParError(2);
0100 
0101     TH1 *h_proj_ref = nullptr;
0102     if (h_ref)
0103     {
0104       h_proj_ref =
0105           h_ref->ProjectionY(
0106               TString::Format(
0107                   "%s_Ref_ProjX_%d_%d",
0108                   h_new->GetName(), bin_start, bin_end),
0109               bin_start, bin_end);
0110       if (pt_range.first < 2.0)
0111       {
0112     //h_proj_ref->GetXaxis()->SetRangeUser(-.05,.05);
0113     h_proj_ref->Rebin(5);
0114       }
0115     }
0116     
0117     DrawReference(h_proj_new, h_proj_ref);
0118 
0119     sprintf(resstr,"#sigma = %.5f #pm %.5f cm", sigma, sigma_unc);
0120     res = new TLatex(0.325,0.825,resstr);
0121     res->SetNDC();
0122     res->SetTextSize(0.05);
0123     res->SetTextAlign(13);
0124     res->Draw();
0125   }
0126   p = (TPad *) c1->cd(idx++);
0127   c1->Update();
0128   TPaveText *pt = new TPaveText(.05,.1,.95,.8);
0129   pt->AddText("No cuts");
0130   pt->Draw();
0131 
0132   //SaveCanvas(c1, TString(qa_file_name_new) + TString("_") + TString(c1->GetName()), true);
0133 
0134 
0135 
0136     c1->Draw();
0137     outfilef->cd();
0138     c1->Write();
0139 }
0140 {
0141     const char *hist_name_prefix = "QAG4SimulationTracking";
0142     TString prefix = TString("h_") + hist_name_prefix + TString("_");
0143 
0144     
0145   // obtain normalization
0146   double Nevent_new = 1;
0147   double Nevent_ref = 1;
0148     
0149   TH2 *h_new2 = (TH2 *) qa_file_new->GetObjectChecked(
0150       prefix + TString("DCArPhi_pT_cuts"), "TH2");
0151   assert(h_new2);
0152 
0153   //  h_new->Rebin(1, 2);
0154   //h_new2->Sumw2();
0155   //  h_new->Scale(1. / Nevent_new);
0156 
0157   TH2 *h_ref2 = NULL;
0158   if (qa_file_ref)
0159   {
0160     h_ref2 = (TH2 *) qa_file_ref->GetObjectChecked(
0161         prefix + TString("DCArPhi_pT_cuts"), "TH2");
0162     assert(h_ref2);
0163 
0164     //    h_ref->Rebin(1, 2);
0165     //h_ref2->Sumw2();
0166     h_ref2->Scale(Nevent_new / Nevent_ref);
0167   }
0168 
0169   TCanvas *c2 = new TCanvas(TString("QA_Draw_Tracking_DCArPhi2") + TString("_") + hist_name_prefix,
0170                             TString("QA_Draw_Tracking_DCArPhi2") + TString("_") + hist_name_prefix,
0171                             1800, 1000);
0172   c2->Divide(4, 2);
0173   int idx2 = 1;
0174   TPad *p2;
0175 
0176   vector<pair<double, double>> gpt_ranges2{
0177       {0, 0.5},
0178       {0.5, 1},
0179       {1, 1.5},
0180       {1.5, 2},
0181       {2, 4},
0182       {4, 16},
0183       {16, 40}};
0184   TF1 *f2 = nullptr;
0185   TF1 *fit2 = nullptr;
0186   Double_t sigma2 = 0;
0187   Double_t sigma_unc2 = 0;
0188   char resstr2[500];
0189   TLatex *res2 = nullptr;
0190   for (auto pt_range : gpt_ranges2)
0191   {
0192     //cout << __PRETTY_FUNCTION__ << " process " << pt_range.first << " - " << pt_range.second << " GeV/c";
0193 
0194     p2 = (TPad *) c2->cd(idx2++);
0195     c2->Update();
0196     p2->SetLogy();
0197 
0198     const double epsilon = 1e-6;
0199     const int bin_start = h_new2->GetXaxis()->FindBin(pt_range.first + epsilon);
0200     const int bin_end = h_new2->GetXaxis()->FindBin(pt_range.second - epsilon);
0201 
0202     TH1 *h_proj_new2 = h_new2->ProjectionY(
0203         TString::Format(
0204             "%s_New_ProjX_%d_%d",
0205             h_new2->GetName(), bin_start, bin_end),
0206         bin_start, bin_end);
0207     if (pt_range.first < 2.0)
0208     {
0209       h_proj_new2->GetXaxis()->SetRangeUser(-.05,.05);
0210       h_proj_new2->Rebin(5);
0211     }
0212     else
0213     {
0214       h_proj_new2->GetXaxis()->SetRangeUser(-.01,.01);
0215     }
0216     h_proj_new2->SetTitle(TString(hist_name_prefix) + TString::Format(
0217                                                          ": %.1f - %.1f GeV/c", pt_range.first, pt_range.second));
0218     h_proj_new2->GetXaxis()->SetTitle(TString::Format(
0219         "DCA (r #phi) [cm]"));
0220     h_proj_new2->GetXaxis()->SetNdivisions(5,5);
0221     f2 = new TF1("f2","gaus",-.01,.01);
0222     h_proj_new2->Fit(f2 , "mq");
0223     fit2 = h_proj_new2->GetFunction("f2");
0224     sigma2 = fit2->GetParameter(2);
0225     sigma_unc2 = fit2->GetParError(2);
0226 
0227     TH1 *h_proj_ref2 = nullptr;
0228     if (h_ref2)
0229     {
0230       h_proj_ref2 =
0231           h_ref2->ProjectionY(
0232               TString::Format(
0233                   "%s_Ref_ProjX_%d_%d",
0234                   h_new2->GetName(), bin_start, bin_end),
0235               bin_start, bin_end);
0236       if (pt_range.first < 2.0)
0237       {
0238     //h_proj_ref->GetXaxis()->SetRangeUser(-.05,.05);
0239     h_proj_ref2->Rebin(5);
0240       }
0241     }
0242     DrawReference(h_proj_new2, h_proj_ref2);
0243 
0244     sprintf(resstr2,"#sigma = %.5f #pm %.5f cm", sigma2, sigma_unc2);
0245     res2 = new TLatex(0.325,0.825,resstr2);
0246     res2->SetNDC();
0247     res2->SetTextSize(0.05);
0248     res2->SetTextAlign(13);
0249     res2->Draw();
0250   }
0251     
0252     
0253   p2 = (TPad *) c2->cd(idx2++);
0254   c2->Update();
0255   TPaveText *pt2 = new TPaveText(.05,.1,.95,.8);
0256   pt2->AddText("Cuts: MVTX hits>=2, INTT hits>=1,");
0257   pt2->AddText("TPC hits>=20");
0258   pt2->Draw();
0259 
0260   //SaveCanvas(c2, TString(qa_file_name_new) + TString("_") + TString(c2->GetName()), true);
0261 
0262     c2->Draw();
0263     outfilef->cd();
0264     c2->Write();
0265  }
0266 
0267 {
0268     const char *hist_name_prefix = "QAG4SimulationTracking";
0269       TString prefix = TString("h_") + hist_name_prefix + TString("_");
0270 
0271     
0272   // obtain normalization
0273   double Nevent_new = 1;
0274   double Nevent_ref = 1;
0275 
0276   if (qa_file_new)
0277   {
0278     TH1 *h_norm = (TH1 *) qa_file_new->GetObjectChecked(
0279         prefix + TString("Normalization"), "TH1");
0280     assert(h_norm);
0281 
0282     Nevent_new = h_norm->GetBinContent(h_norm->GetXaxis()->FindBin("Truth Track"));
0283   }
0284   if (qa_file_ref)
0285   {
0286     TH1 *h_norm = (TH1 *) qa_file_ref->GetObjectChecked(
0287         prefix + TString("Normalization"), "TH1");
0288     assert(h_norm);
0289 
0290     Nevent_ref = h_norm->GetBinContent(h_norm->GetXaxis()->FindBin("Truth Track"));
0291   }
0292 
0293   TH2 *h_new = (TH2 *) qa_file_new->GetObjectChecked(
0294       prefix + TString("DCAZ_pT"), "TH2");
0295   assert(h_new);
0296 
0297   //  h_new->Rebin(1, 2);
0298   //h_new->Sumw2();
0299   //  h_new->Scale(1. / Nevent_new);
0300 
0301   TH2 *h_ref = NULL;
0302   if (qa_file_ref)
0303   {
0304     h_ref = (TH2 *) qa_file_ref->GetObjectChecked(
0305         prefix + TString("DCAZ_pT"), "TH2");
0306     assert(h_ref);
0307 
0308     //    h_ref->Rebin(1, 2);
0309    // h_ref->Sumw2();
0310     h_ref->Scale(Nevent_new / Nevent_ref);
0311   }
0312 
0313   TCanvas *c1 = new TCanvas(TString("QA_Draw_Tracking_DCAZ") + TString("_") + hist_name_prefix,
0314                             TString("QA_Draw_Tracking_DCAZ") + TString("_") + hist_name_prefix,
0315                             1800, 1000);
0316   c1->Divide(4, 2);
0317   int idx = 1;
0318   TPad *p;
0319 
0320   vector<pair<double, double>> gpt_ranges{
0321       {0, 0.5},
0322       {0.5, 1},
0323       {1, 1.5},
0324       {1.5, 2},
0325       {2, 4},
0326       {4, 16},
0327       {16, 40}};
0328   TF1 *f1 = nullptr;
0329   TF1 *fit = nullptr;
0330   Double_t sigma = 0;
0331   Double_t sigma_unc = 0;
0332   char resstr[500];
0333   TLatex *res = nullptr;
0334   for (auto pt_range : gpt_ranges)
0335   {
0336     //cout << __PRETTY_FUNCTION__ << " process " << pt_range.first << " - " << pt_range.second << " GeV/c";
0337 
0338     p = (TPad *) c1->cd(idx++);
0339     c1->Update();
0340     p->SetLogy();
0341 
0342     const double epsilon = 1e-6;
0343     const int bin_start = h_new->GetXaxis()->FindBin(pt_range.first + epsilon);
0344     const int bin_end = h_new->GetXaxis()->FindBin(pt_range.second - epsilon);
0345 
0346     TH1 *h_proj_new = h_new->ProjectionY(
0347         TString::Format(
0348             "%s_New_ProjX_%d_%d",
0349             h_new->GetName(), bin_start, bin_end),
0350         bin_start, bin_end);
0351     if (pt_range.first < 2.0)
0352     {
0353       h_proj_new->GetXaxis()->SetRangeUser(-.05, .05);
0354       h_proj_new->Rebin(5);
0355     }
0356     else
0357     {
0358       h_proj_new->GetXaxis()->SetRangeUser(-.01, .01);
0359     }
0360     h_proj_new->SetTitle(TString(hist_name_prefix) + TString::Format(
0361                                                          ": %.1f - %.1f GeV/c", pt_range.first, pt_range.second));
0362     h_proj_new->GetXaxis()->SetTitle(TString::Format(
0363         "DCA (Z) [cm]"));
0364     h_proj_new->GetXaxis()->SetNdivisions(5, 5);
0365 
0366     f1 = new TF1("f1", "gaus", -.01, .01);
0367     h_proj_new->Fit(f1,"mq");
0368     sigma = f1->GetParameter(2);
0369     sigma_unc = f1->GetParError(2);
0370 
0371     TH1 *h_proj_ref = nullptr;
0372     if (h_ref)
0373     {
0374       h_proj_ref =
0375           h_ref->ProjectionY(
0376               TString::Format(
0377                   "%s_Ref_ProjX_%d_%d",
0378                   h_new->GetName(), bin_start, bin_end),
0379               bin_start, bin_end);
0380       if (pt_range.first < 2.0)
0381       {
0382         //h_proj_ref->GetXaxis()->SetRangeUser(-.05,.05);
0383         h_proj_ref->Rebin(5);
0384       }
0385     }
0386     DrawReference(h_proj_new, h_proj_ref);
0387     sprintf(resstr, "#sigma = %.5f #pm %.5f cm", sigma, sigma_unc);
0388     res = new TLatex(0.325, 0.825, resstr);
0389     res->SetNDC();
0390     res->SetTextSize(0.05);
0391     res->SetTextAlign(13);
0392     res->Draw();
0393   }
0394   p = (TPad *) c1->cd(idx++);
0395   c1->Update();
0396   TPaveText *pt = new TPaveText(.05, .1, .95, .8);
0397   pt->AddText("No cuts");
0398   pt->Draw();
0399 
0400 //  SaveCanvas(c1, TString(qa_file_name_new) + TString("_") + TString(c1->GetName()), true);
0401     
0402     c1->Draw();
0403     outfilef->cd();
0404     c1->Write();
0405  }
0406 
0407 {
0408     
0409     const char *hist_name_prefix = "QAG4SimulationTracking";
0410       TString prefix = TString("h_") + hist_name_prefix + TString("_");
0411 
0412     
0413   // obtain normalization
0414   double Nevent_new = 1;
0415   double Nevent_ref = 1;
0416     
0417   // cuts plots
0418   TH2 *h_new2 = (TH2 *) qa_file_new->GetObjectChecked(
0419       prefix + TString("DCAZ_pT_cuts"), "TH2");
0420   assert(h_new2);
0421 
0422   //  h_new->Rebin(1, 2);
0423   //h_new2->Sumw2();
0424   //  h_new->Scale(1. / Nevent_new);
0425 
0426   TH2 *h_ref2 = NULL;
0427   if (qa_file_ref)
0428   {
0429     h_ref2 = (TH2 *) qa_file_ref->GetObjectChecked(
0430         prefix + TString("DCAZ_pT_cuts"), "TH2");
0431     assert(h_ref2);
0432 
0433     //    h_ref->Rebin(1, 2);
0434     //h_ref2->Sumw2();
0435     h_ref2->Scale(Nevent_new / Nevent_ref);
0436   }
0437 
0438   TCanvas *c2 = new TCanvas(TString("QA_Draw_Tracking_DCAZ2") + TString("_") + hist_name_prefix,
0439                             TString("QA_Draw_Tracking_DCAZ2") + TString("_") + hist_name_prefix,
0440                             1800, 1000);
0441   c2->Divide(4, 2);
0442   int idx2 = 1;
0443   TPad *p2;
0444 
0445   vector<pair<double, double>> gpt_ranges2{
0446       {0, 0.5},
0447       {0.5, 1},
0448       {1, 1.5},
0449       {1.5, 2},
0450       {2, 4},
0451       {4, 16},
0452       {16, 40}};
0453   TF1 *f2 = nullptr;
0454   TF1 *fit2 = nullptr;
0455   Double_t sigma2 = 0;
0456   Double_t sigma_unc2 = 0;
0457   char resstr2[500];
0458   TLatex *res2 = nullptr;
0459   for (auto pt_range : gpt_ranges2)
0460   {
0461     //cout << __PRETTY_FUNCTION__ << " process " << pt_range.first << " - " << pt_range.second << " GeV/c";
0462 
0463     p2 = (TPad *) c2->cd(idx2++);
0464     c2->Update();
0465     p2->SetLogy();
0466 
0467     const double epsilon = 1e-6;
0468     const int bin_start = h_new2->GetXaxis()->FindBin(pt_range.first + epsilon);
0469     const int bin_end = h_new2->GetXaxis()->FindBin(pt_range.second - epsilon);
0470 
0471     TH1 *h_proj_new2 = h_new2->ProjectionY(
0472         TString::Format(
0473             "%s_New_ProjX_%d_%d",
0474             h_new2->GetName(), bin_start, bin_end),
0475         bin_start, bin_end);
0476     if (pt_range.first < 2.0)
0477     {
0478       h_proj_new2->GetXaxis()->SetRangeUser(-.05, .05);
0479       h_proj_new2->Rebin(5);
0480     }
0481     else
0482     {
0483       h_proj_new2->GetXaxis()->SetRangeUser(-.01, .01);
0484     }
0485     h_proj_new2->SetTitle(TString(hist_name_prefix) + TString::Format(
0486                                                           ": %.1f - %.1f GeV/c", pt_range.first, pt_range.second));
0487     h_proj_new2->GetXaxis()->SetTitle(TString::Format(
0488         "DCA (Z) [cm]"));
0489     h_proj_new2->GetXaxis()->SetNdivisions(5, 5);
0490 
0491     f2 = new TF1("f2", "gaus", -.01, .01);
0492     h_proj_new2->Fit(f2,"mq");
0493     fit2 = h_proj_new2->GetFunction("f2");
0494     sigma2 = fit2->GetParameter(2);
0495     sigma_unc2 = fit2->GetParError(2);
0496 
0497     TH1 *h_proj_ref2 = nullptr;
0498     if (h_ref2)
0499     {
0500       h_proj_ref2 =
0501           h_ref2->ProjectionY(
0502               TString::Format(
0503                   "%s_Ref_ProjX_%d_%d",
0504                   h_new2->GetName(), bin_start, bin_end),
0505               bin_start, bin_end);
0506       if (pt_range.first < 2.0)
0507       {
0508         //h_proj_ref->GetXaxis()->SetRangeUser(-.05,.05);
0509         h_proj_ref2->Rebin(5);
0510       }
0511     }
0512     DrawReference(h_proj_new2, h_proj_ref2);
0513     sprintf(resstr2, "#sigma = %.5f #pm %.5f cm", sigma2, sigma_unc2);
0514     res2 = new TLatex(0.325, 0.825, resstr2);
0515     res2->SetNDC();
0516     res2->SetTextSize(0.05);
0517     res2->SetTextAlign(13);
0518     res2->Draw();
0519   }
0520   p2 = (TPad *) c2->cd(idx2++);
0521   c2->Update();
0522   TPaveText *pt2 = new TPaveText(.05, .1, .95, .8);
0523   pt2->AddText("Cuts: MVTX hits>=2, INTT hits>=1,");
0524   pt2->AddText("TPC hits>=20");
0525   pt2->Draw();
0526 
0527 //  SaveCanvas(c2, TString(qa_file_name_new) + TString("_") + TString(c2->GetName()), true);
0528     c2->Draw();
0529     outfilef->cd();
0530     c2->Write();
0531  }
0532 
0533 {
0534     const char *hist_name_prefix = "QAG4SimulationTracking";
0535     TString prefix = TString("h_") + hist_name_prefix + TString("_");
0536     
0537     
0538   // obtain normalization
0539   double Nevent_new = 1;
0540   double Nevent_ref = 1;
0541 
0542   if (qa_file_new)
0543   {
0544 
0545     TH1 *h_norm = (TH1 *) qa_file_new->GetObjectChecked(
0546         prefix + TString("Normalization"), "TH1");
0547     assert(h_norm);
0548 
0549     Nevent_new = h_norm->GetBinContent(h_norm->GetXaxis()->FindBin("Event"));
0550   }
0551   if (qa_file_ref)
0552   {
0553     TH1 *h_norm = (TH1 *) qa_file_ref->GetObjectChecked(
0554         prefix + TString("Normalization"), "TH1");
0555     assert(h_norm);
0556 
0557     Nevent_ref = h_norm->GetBinContent(h_norm->GetXaxis()->FindBin("Event"));
0558   }
0559 
0560   TCanvas *c1 = new TCanvas(TString("QA_Draw_Tracking_DCA_Resolution") + TString("_") + hist_name_prefix,
0561                             TString("QA_Draw_Tracking_DCA_Resolution") + TString("_") + hist_name_prefix,
0562                             1800, 1000);
0563   c1->Divide(2, 1);
0564   int idx = 1;
0565   TPad *p;
0566 
0567   {
0568     p = (TPad *) c1->cd(idx++);
0569     c1->Update();
0570     p->SetLogx();
0571     TH1 *frame = p->DrawFrame(0.1, -0.01, 50, 0.01,
0572                               ";Truth p_{T} [GeV/c];<DCA (r #phi)> #pm #sigma(DCA (r #phi)) [cm]");
0573     gPad->SetLeftMargin(.2);
0574     frame->GetYaxis()->SetTitleOffset(2);
0575     TLine *l = new TLine(0.1, 0, 50, 0);
0576     l->SetLineColor(kGray);
0577     l->Draw();
0578 
0579     TH2 *h_QAG4SimulationTracking_DCArPhi = (TH2 *) qa_file_new->GetObjectChecked(
0580         prefix + "DCArPhi_pT_cuts", "TH2");
0581     assert(h_QAG4SimulationTracking_DCArPhi);
0582 
0583     h_QAG4SimulationTracking_DCArPhi->Rebin2D(20, 1);
0584 
0585     // h_QAG4SimulationTracking_DCArPhi->Draw("colz");
0586     TGraphErrors *ge_QAG4SimulationTracking_DCArPhi = FitProfile(h_QAG4SimulationTracking_DCArPhi);
0587     ge_QAG4SimulationTracking_DCArPhi->Draw("pe");
0588 
0589     TGraphErrors *h_ratio_ref = NULL;
0590     if (qa_file_ref)
0591     {
0592       TH2 *h_QAG4SimulationTracking_DCArPhi = (TH2 *) qa_file_ref->GetObjectChecked(
0593           prefix + "DCArPhi_pT_cuts", "TH2");
0594       assert(h_QAG4SimulationTracking_DCArPhi);
0595 
0596       h_QAG4SimulationTracking_DCArPhi->Rebin2D(20, 1);
0597 
0598       h_ratio_ref = FitProfile(h_QAG4SimulationTracking_DCArPhi);
0599       ge_QAG4SimulationTracking_DCArPhi->Draw("pe");
0600     }
0601 
0602     ge_QAG4SimulationTracking_DCArPhi->SetTitle("DCA (r #phi, #geq 2MVTX, #geq 1INTT, #geq 20TPC) [cm]");
0603     DrawReference(ge_QAG4SimulationTracking_DCArPhi, h_ratio_ref, true);
0604   }
0605 
0606   {
0607     p = (TPad *) c1->cd(idx++);
0608     c1->Update();
0609     p->SetLogx();
0610     TH1 *frame = p->DrawFrame(0.1, -0.01, 50, 0.01,
0611                               "DCA (Z) [cm];Truth p_{T} [GeV/c];<DCA (Z)> #pm #sigma(DCA (Z)) [cm]");
0612     // gPad->SetLeftMargin(.2);
0613     gPad->SetTopMargin(-1);
0614     frame->GetYaxis()->SetTitleOffset(1.7);
0615     //TLine *l = new TLine(0.1, 0, 50, 0);
0616     //l->SetLineColor(kGray);
0617     //l->Draw();
0618     HorizontalLine(gPad, 1)->Draw();
0619 
0620     TH2 *h_QAG4SimulationTracking_DCAZ = (TH2 *) qa_file_new->GetObjectChecked(
0621         prefix + "DCAZ_pT_cuts", "TH2");
0622     assert(h_QAG4SimulationTracking_DCAZ);
0623 
0624     h_QAG4SimulationTracking_DCAZ->Rebin2D(40, 1);
0625 
0626     TGraphErrors *ge_QAG4SimulationTracking_DCAZ = FitProfile(h_QAG4SimulationTracking_DCAZ);
0627     ge_QAG4SimulationTracking_DCAZ->Draw("pe");
0628     ge_QAG4SimulationTracking_DCAZ->SetTitle("DCA (Z) [cm]");
0629 
0630     TGraphErrors *h_ratio_ref = NULL;
0631     if (qa_file_ref)
0632     {
0633       TH2 *h_QAG4SimulationTracking_DCAZ = (TH2 *) qa_file_ref->GetObjectChecked(
0634           prefix + "DCAZ_pT_cuts", "TH2");
0635       assert(h_QAG4SimulationTracking_DCAZ);
0636 
0637       h_QAG4SimulationTracking_DCAZ->Rebin2D(40, 1);
0638 
0639       h_ratio_ref = FitProfile(h_QAG4SimulationTracking_DCAZ);
0640       ge_QAG4SimulationTracking_DCAZ->Draw("pe");
0641     }
0642 
0643     DrawReference(ge_QAG4SimulationTracking_DCAZ, h_ratio_ref, true);
0644   }
0645 
0646   //SaveCanvas(c1, TString(qa_file_name_new) + TString("_") + TString(c1->GetName()), true);
0647     c1->Draw();
0648     outfilef->cd();
0649     c1->Write();
0650  }
0651 
0652  {
0653     const char *hist_name_prefix = "QAG4SimulationTracking";
0654   TString prefix = TString("h_") + hist_name_prefix + TString("_");
0655 
0656   // obtain normalization
0657   double Nevent_new = 1;
0658   double Nevent_ref = 1;
0659 
0660   if (qa_file_new)
0661   {
0662     TH1 *h_norm = (TH1 *) qa_file_new->GetObjectChecked(
0663         prefix + TString("Normalization"), "TH1");
0664     assert(h_norm);
0665 
0666     Nevent_new = h_norm->GetBinContent(h_norm->GetXaxis()->FindBin("Truth Track"));
0667   }
0668   if (qa_file_ref)
0669   {
0670     TH1 *h_norm = (TH1 *) qa_file_ref->GetObjectChecked(
0671         prefix + TString("Normalization"), "TH1");
0672     assert(h_norm);
0673 
0674     Nevent_ref = h_norm->GetBinContent(h_norm->GetXaxis()->FindBin("Truth Track"));
0675   }
0676     
0677     
0678     
0679   TH2 *h_new = (TH2 *) qa_file_new->GetObjectChecked(
0680       prefix + TString("SigmalizedDCArPhi_pT"), "TH2");
0681   assert(h_new);
0682 
0683   //  h_new->Rebin(1, 2);
0684   //h_new->Sumw2();
0685   //  h_new->Scale(1. / Nevent_new);
0686 
0687   TH2 *h_ref = NULL;
0688   if (qa_file_ref)
0689   {
0690     h_ref = (TH2 *) qa_file_ref->GetObjectChecked(
0691         prefix + TString("SigmalizedDCArPhi_pT"), "TH2");
0692     assert(h_ref);
0693 
0694     //    h_ref->Rebin(1, 2);
0695     //h_ref->Sumw2();
0696     h_ref->Scale(Nevent_new / Nevent_ref);
0697   }
0698 
0699   TCanvas *c1 = new TCanvas(TString("QA_Draw_Tracking_SigmalizedDCArPhi") + TString("_") + hist_name_prefix,
0700                             TString("QA_Draw_Tracking_SigmalizedDCArPhi") + TString("_") + hist_name_prefix,
0701                             1800, 1000);
0702   c1->Divide(4, 2);
0703   int idx = 1;
0704   TPad *p;
0705 
0706   vector<pair<double, double>> gpt_ranges{
0707       {0, 0.5},
0708       {0.5, 1},
0709       {1, 1.5},
0710       {1.5, 2},
0711       {2, 4},
0712       {4, 16},
0713       {16, 40}};
0714   TF1 *f1 = nullptr;
0715   TF1 *fit = nullptr;
0716   Double_t sigma = 0;
0717   Double_t sigma_unc = 0;
0718   char resstr[500];
0719   TLatex *res = nullptr;
0720   for (auto pt_range : gpt_ranges)
0721   {
0722     // cout << __PRETTY_FUNCTION__ << " process " << pt_range.first << " - " << pt_range.second << " GeV/c";
0723 
0724     p = (TPad *) c1->cd(idx++);
0725     c1->Update();
0726     p->SetLogy();
0727 
0728     const double epsilon = 1e-6;
0729     const int bin_start = h_new->GetXaxis()->FindBin(pt_range.first + epsilon);
0730     const int bin_end = h_new->GetXaxis()->FindBin(pt_range.second - epsilon);
0731 
0732     TH1 *h_proj_new = h_new->ProjectionY(
0733         TString::Format(
0734             "%s_New_ProjX_%d_%d",
0735             h_new->GetName(), bin_start, bin_end),
0736         bin_start, bin_end);
0737     h_proj_new->GetXaxis()->SetRangeUser(-5.,5.);
0738     h_proj_new->Rebin(5);
0739     h_proj_new->SetTitle(TString(hist_name_prefix) + TString::Format(
0740                                                          ": %.1f - %.1f GeV/c", pt_range.first, pt_range.second));
0741     h_proj_new->GetXaxis()->SetTitle(TString::Format(
0742         "Sigmalized DCA (r #phi)"));
0743     h_proj_new->GetXaxis()->SetNdivisions(5,5);
0744 
0745     f1 = new TF1("f1","gaus",-4.,4.);
0746     h_proj_new->Fit(f1, "mq");
0747     sigma = f1->GetParameter(2);
0748     sigma_unc = f1->GetParError(2);
0749     
0750     TH1 *h_proj_ref = nullptr;
0751     if (h_ref)
0752     {
0753       h_proj_ref =
0754           h_ref->ProjectionY(
0755               TString::Format(
0756                   "%s_Ref_ProjX_%d_%d",
0757                   h_new->GetName(), bin_start, bin_end),
0758               bin_start, bin_end);
0759       //h_proj_ref->GetXaxis()->SetRangeUser(-.05,.05);
0760       h_proj_ref->Rebin(5);
0761     }
0762     
0763     DrawReference(h_proj_new, h_proj_ref);
0764     sprintf(resstr,"#sigma = %.5f #pm %.5f", sigma, sigma_unc);
0765     res = new TLatex(0.325,0.825,resstr);
0766     res->SetNDC();
0767     res->SetTextSize(0.05);
0768     res->SetTextAlign(13);
0769     res->Draw();
0770   }
0771   p = (TPad *) c1->cd(idx++);
0772   c1->Update();
0773   TPaveText *pt = new TPaveText(.05,.1,.95,.8);
0774   pt->AddText("Cuts: MVTX hits>=2, INTT hits>=1,");
0775   pt->AddText("TPC hits>=20");
0776   pt->Draw();
0777 
0778   // SaveCanvas(c1, TString(qa_file_name_new) + TString("_") + TString(c1->GetName()), true);
0779   c1->Draw();
0780   outfilef->cd();
0781   c1->Write();
0782 }
0783 
0784  {
0785     
0786     const char *hist_name_prefix = "QAG4SimulationTracking";
0787   TString prefix = TString("h_") + hist_name_prefix + TString("_");
0788 
0789   // obtain normalization
0790   double Nevent_new = 1;
0791   double Nevent_ref = 1;
0792 
0793   if (qa_file_new)
0794   {
0795     TH1 *h_norm = (TH1 *) qa_file_new->GetObjectChecked(
0796         prefix + TString("Normalization"), "TH1");
0797     assert(h_norm);
0798 
0799     Nevent_new = h_norm->GetBinContent(h_norm->GetXaxis()->FindBin("Truth Track"));
0800   }
0801   if (qa_file_ref)
0802   {
0803     TH1 *h_norm = (TH1 *) qa_file_ref->GetObjectChecked(
0804         prefix + TString("Normalization"), "TH1");
0805     assert(h_norm);
0806 
0807     Nevent_ref = h_norm->GetBinContent(h_norm->GetXaxis()->FindBin("Truth Track"));
0808   }
0809     
0810   TH2 *h_new2 = (TH2 *) qa_file_new->GetObjectChecked(
0811       prefix + TString("SigmalizedDCAZ_pT"), "TH2");
0812   assert(h_new2);
0813 
0814   //  h_new->Rebin(1, 2);
0815   //h_new2->Sumw2();
0816   //  h_new->Scale(1. / Nevent_new);
0817 
0818   TH2 *h_ref2 = NULL;
0819   if (qa_file_ref)
0820   {
0821     h_ref2 = (TH2 *) qa_file_ref->GetObjectChecked(
0822         prefix + TString("SigmalizedDCAZ_pT"), "TH2");
0823     assert(h_ref2);
0824 
0825     //    h_ref->Rebin(1, 2);
0826     //h_ref2->Sumw2();
0827     h_ref2->Scale(Nevent_new / Nevent_ref);
0828   }
0829 
0830   TCanvas *c2 = new TCanvas(TString("QA_Draw_Tracking_SigmalizedDCAZ") + TString("_") + hist_name_prefix,
0831                             TString("QA_Draw_Tracking_SigmalizedDCAZ") + TString("_") + hist_name_prefix,
0832                             1800, 1000);
0833   c2->Divide(4, 2);
0834   int idx2 = 1;
0835   TPad *p2;
0836 
0837   vector<pair<double, double>> gpt_ranges2{
0838       {0, 0.5},
0839       {0.5, 1},
0840       {1, 1.5},
0841       {1.5, 2},
0842       {2, 4},
0843       {4, 16},
0844       {16, 40}};
0845   TF1 *f2 = nullptr;
0846   TF1 *fit2 = nullptr;
0847   Double_t sigma2 = 0;
0848   Double_t sigma_unc2 = 0;
0849   char resstr2[500];
0850   TLatex *res2 = nullptr;
0851   for (auto pt_range : gpt_ranges2)
0852   {
0853    // cout << __PRETTY_FUNCTION__ << " process " << pt_range.first << " - " << pt_range.second << " GeV/c";
0854 
0855     p2 = (TPad *) c2->cd(idx2++);
0856     c2->Update();
0857     p2->SetLogy();
0858 
0859     const double epsilon = 1e-6;
0860     const int bin_start = h_new2->GetXaxis()->FindBin(pt_range.first + epsilon);
0861     const int bin_end = h_new2->GetXaxis()->FindBin(pt_range.second - epsilon);
0862 
0863     TH1 *h_proj_new2 = h_new2->ProjectionY(
0864         TString::Format(
0865             "%s_New_ProjX_%d_%d",
0866             h_new2->GetName(), bin_start, bin_end),
0867         bin_start, bin_end);
0868     h_proj_new2->GetXaxis()->SetRangeUser(-5.,5.);
0869     h_proj_new2->Rebin(5);
0870     h_proj_new2->SetTitle(TString(hist_name_prefix) + TString::Format(
0871                                                          ": %.1f - %.1f GeV/c", pt_range.first, pt_range.second));
0872     h_proj_new2->GetXaxis()->SetTitle(TString::Format(
0873         "Sigmalized DCA (Z)"));
0874     h_proj_new2->GetXaxis()->SetNdivisions(5,5);
0875     
0876     f2 = new TF1("f2","gaus",-4.,4.);
0877     h_proj_new2->Fit(f2, "mq");
0878     fit2 = f2; //h_proj_new2->GetFunction("f2");
0879     sigma2 = fit2->GetParameter(2);
0880     sigma_unc2 = fit2->GetParError(2);
0881 
0882     TH1 *h_proj_ref2 = nullptr;
0883     if (h_ref2)
0884     {
0885       h_proj_ref2 =
0886           h_ref2->ProjectionY(
0887               TString::Format(
0888                   "%s_Ref_ProjX_%d_%d",
0889                   h_new2->GetName(), bin_start, bin_end),
0890               bin_start, bin_end);
0891       //h_proj_ref->GetXaxis()->SetRangeUser(-.05,.05);
0892       h_proj_ref2->Rebin(5);
0893     }
0894     DrawReference(h_proj_new2, h_proj_ref2);
0895     sprintf(resstr2,"#sigma = %.5f #pm %.5f", sigma2, sigma_unc2);
0896     res2 = new TLatex(0.325,0.825,resstr2);
0897     res2->SetNDC();
0898     res2->SetTextSize(0.05);
0899     res2->SetTextAlign(13);
0900     res2->Draw();
0901   }
0902   p2 = (TPad *) c2->cd(idx2++);
0903   c2->Update();
0904   TPaveText *pt2 = new TPaveText(.05,.1,.95,.8);
0905   pt2->AddText("Cuts: MVTX hits>=2, INTT hits>=1,");
0906   pt2->AddText("TPC hits>=20");
0907   pt2->Draw();
0908 
0909   //SaveCanvas(c2, TString(qa_file_name_new) + TString("_") + TString(c2->GetName()), true);
0910     c2->Draw();
0911     outfilef->cd();
0912     c2->Write();
0913 }
0914 
0915 {
0916     const char *hist_name_prefix = "QAG4SimulationTracking";
0917     TString prefix = TString("h_") + hist_name_prefix + TString("_");
0918     
0919     
0920   // obtain normalization
0921   double Nevent_new = 1;
0922   double Nevent_ref = 1;
0923 
0924   if (qa_file_new)
0925   {
0926 
0927     TH1 *h_norm = (TH1 *) qa_file_new->GetObjectChecked(
0928         prefix + TString("Normalization"), "TH1");
0929     assert(h_norm);
0930 
0931     Nevent_new = h_norm->GetBinContent(h_norm->GetXaxis()->FindBin("Event"));
0932   }
0933   if (qa_file_ref)
0934   {
0935     TH1 *h_norm = (TH1 *) qa_file_ref->GetObjectChecked(
0936         prefix + TString("Normalization"), "TH1");
0937     assert(h_norm);
0938 
0939     Nevent_ref = h_norm->GetBinContent(h_norm->GetXaxis()->FindBin("Event"));
0940   }
0941 
0942   TCanvas *c1 = new TCanvas(TString("QA_Draw_Tracking_SigmalizedDCA_Resolution") + TString("_") + hist_name_prefix,
0943                             TString("QA_Draw_Tracking_SigmalizedDCA_Resolution") + TString("_") + hist_name_prefix,
0944                             1800, 1000);
0945   c1->Divide(2, 1);
0946   int idx = 1;
0947   TPad *p;
0948 
0949   {
0950     p = (TPad *) c1->cd(idx++);
0951     c1->Update();
0952     p->SetLogx();
0953     TH1 *frame = p->DrawFrame(0.1, -2, 50, 2,
0954                               ";Truth p_{T} [GeV/c];<Sigmalized DCA (r #phi)> #pm #sigma(Sigmalized DCA (r #phi))");
0955     gPad->SetLeftMargin(.2);
0956     frame->GetYaxis()->SetTitleOffset(2);
0957     TLine *l = new TLine(0.1, 0, 50, 0);
0958     l->SetLineColor(kGray);
0959     l->Draw();
0960     HorizontalLine(gPad, 1)->Draw();
0961     HorizontalLine(gPad, -1)->Draw();
0962 
0963     TH2 *h_QAG4SimulationTracking_DCArPhi = (TH2 *) qa_file_new->GetObjectChecked(
0964         prefix + "SigmalizedDCArPhi_pT", "TH2");
0965     assert(h_QAG4SimulationTracking_DCArPhi);
0966 
0967     h_QAG4SimulationTracking_DCArPhi->Rebin2D(20, 1);
0968 
0969     // h_QAG4SimulationTracking_DCArPhi->Draw("colz");
0970     TGraphErrors *ge_QAG4SimulationTracking_DCArPhi = FitProfile(h_QAG4SimulationTracking_DCArPhi);
0971     ge_QAG4SimulationTracking_DCArPhi->Draw("pe");
0972 
0973     TGraphErrors *h_ratio_ref = NULL;
0974     if (qa_file_ref)
0975     {
0976       TH2 *h_QAG4SimulationTracking_DCArPhi = (TH2 *) qa_file_ref->GetObjectChecked(
0977           prefix + "SigmalizedDCArPhi_pT", "TH2");
0978       assert(h_QAG4SimulationTracking_DCArPhi);
0979 
0980       h_QAG4SimulationTracking_DCArPhi->Rebin2D(20, 1);
0981 
0982       h_ratio_ref = FitProfile(h_QAG4SimulationTracking_DCArPhi);
0983       ge_QAG4SimulationTracking_DCArPhi->Draw("pe");
0984     }
0985 
0986     ge_QAG4SimulationTracking_DCArPhi->SetTitle("DCA_{r#phi}/#sigma[DCA_{r#phi}]");
0987     DrawReference(ge_QAG4SimulationTracking_DCArPhi, h_ratio_ref, true);
0988   }
0989 
0990   {
0991     p = (TPad *) c1->cd(idx++);
0992     c1->Update();
0993     p->SetLogx();
0994     TH1 *frame = p->DrawFrame(0.1, -2, 50, 2,
0995                               "DCA_z/#sigma[DCA_z];Truth p_{T} [GeV/c];<Sigmalized DCA (Z)> #pm #sigma(Sigmalized DCA (Z))");
0996     // gPad->SetLeftMargin(.2);
0997     gPad->SetTopMargin(-1);
0998     frame->GetYaxis()->SetTitleOffset(1.7);
0999     TLine *l = new TLine(0.1, 0, 50, 0);
1000     l->SetLineColor(kGray);
1001     l->Draw();
1002     HorizontalLine(gPad, 1)->Draw();
1003     HorizontalLine(gPad, -1)->Draw();
1004 
1005     TH2 *h_QAG4SimulationTracking_DCAZ = (TH2 *) qa_file_new->GetObjectChecked(
1006         prefix + "SigmalizedDCAZ_pT", "TH2");
1007     assert(h_QAG4SimulationTracking_DCAZ);
1008 
1009     h_QAG4SimulationTracking_DCAZ->Rebin2D(40, 1);
1010 
1011     TGraphErrors *ge_QAG4SimulationTracking_DCAZ = FitProfile(h_QAG4SimulationTracking_DCAZ);
1012     ge_QAG4SimulationTracking_DCAZ->Draw("pe");
1013     ge_QAG4SimulationTracking_DCAZ->SetTitle("DCA_z/#sigma[DCA_z]");
1014 
1015     TGraphErrors *h_ratio_ref = NULL;
1016     if (qa_file_ref)
1017     {
1018       TH2 *h_QAG4SimulationTracking_DCAZ = (TH2 *) qa_file_ref->GetObjectChecked(
1019           prefix + "SigmalizedDCAZ_pT", "TH2");
1020       assert(h_QAG4SimulationTracking_DCAZ);
1021 
1022       h_QAG4SimulationTracking_DCAZ->Rebin2D(40, 1);
1023 
1024       h_ratio_ref = FitProfile(h_QAG4SimulationTracking_DCAZ);
1025       ge_QAG4SimulationTracking_DCAZ->Draw("pe");
1026     }
1027 
1028     DrawReference(ge_QAG4SimulationTracking_DCAZ, h_ratio_ref, true);
1029   }
1030 
1031   //SaveCanvas(c1, TString(qa_file_name_new) + TString("_") + TString(c1->GetName()), true);
1032     c1->Draw();
1033     
1034     outfilef->cd();
1035     c1->Write();
1036 }
1037 
1038 
1039 {
1040     const char *hist_name_prefix = "QAG4SimulationVertex_SvtxVertexMap";
1041     TString prefix = TString("h_") + hist_name_prefix + TString("_");
1042     
1043     
1044   TCanvas *c1 = new TCanvas(TString("QA_Draw_Vertex_nVertex") +
1045                                 TString("_") + hist_name_prefix,
1046                             TString("QA_Draw_Vertex_nVertex") +
1047                                 TString("_") + hist_name_prefix,
1048                             1800, 1000);
1049   c1->Divide(2, 2);
1050   int idx = 1;
1051   TPad *p;
1052 
1053   {
1054     static const int nrebin = 1;
1055 
1056     p = (TPad *)c1->cd(idx++);
1057     c1->Update();
1058     // p->SetLogx();
1059     p->SetGridy();
1060 
1061     TH1 *h_pass =
1062         (TH1 *)qa_file_new->GetObjectChecked(prefix + "gntracks", "TH1");
1063     assert(h_pass);
1064 
1065     h_pass->Rebin(nrebin);
1066 
1067     //    h_ratio->GetXaxis()->SetRangeUser(min_Et, max_Et);
1068     h_pass->GetYaxis()->SetTitle("Counts");
1069     // h_pass->GetYaxis()->SetRangeUser(-0, 1.);
1070 
1071     TH1 *h_ref = NULL;
1072     if (qa_file_ref) {
1073       h_ref =
1074           (TH1 *)qa_file_ref->GetObjectChecked(prefix + "gntracks", "TH1");
1075       assert(h_ref);
1076 
1077       h_ref->Rebin(nrebin);
1078     }
1079 
1080     h_pass->SetTitle(TString(hist_name_prefix) + ": gntracks");
1081 
1082     DrawReference(h_pass, h_ref, false);
1083   }
1084 
1085   {
1086     static const int nrebin = 1;
1087     
1088     p = (TPad *)c1->cd(idx++);
1089     c1->Update();
1090     // p->SetLogx();
1091     p->SetGridy();
1092 
1093     TH1 *h_pass =
1094         (TH1 *)qa_file_new->GetObjectChecked(prefix + "gntracksmaps", "TH1");
1095     assert(h_pass);
1096 
1097     h_pass->Rebin(nrebin);
1098 
1099     //    h_ratio->GetXaxis()->SetRangeUser(min_Et, max_Et);
1100     h_pass->GetYaxis()->SetTitle("Counts");
1101     // h_pass->GetYaxis()->SetRangeUser(-0, 1.);
1102 
1103     TH1 *h_ref = NULL;
1104     if (qa_file_ref) {
1105       h_ref =
1106           (TH1 *)qa_file_ref->GetObjectChecked(prefix + "gntracksmaps", "TH1");
1107       assert(h_pass);
1108 
1109       h_ref->Rebin(nrebin);
1110     }
1111 
1112     h_pass->SetTitle(TString(hist_name_prefix) + ": gntracksmaps");
1113 
1114     DrawReference(h_pass, h_ref, false);
1115   }
1116 
1117   {
1118     static const int nrebin = 1;
1119     
1120     p = (TPad *)c1->cd(idx++);
1121     c1->Update();
1122     // p->SetLogx();
1123     p->SetGridy();
1124 
1125     TH1 *h_pass =
1126         (TH1 *)qa_file_new->GetObjectChecked(prefix + "ntracks", "TH1");
1127     assert(h_pass);
1128 
1129     h_pass->Rebin(nrebin);
1130 
1131     //    h_ratio->GetXaxis()->SetRangeUser(min_Et, max_Et);
1132     h_pass->GetYaxis()->SetTitle("Counts");
1133     // h_pass->GetYaxis()->SetRangeUser(-0, 1.);
1134 
1135     TH1 *h_ref = NULL;
1136     if (qa_file_ref) {
1137       h_ref =
1138           (TH1 *)qa_file_ref->GetObjectChecked(prefix + "ntracks", "TH1");
1139       assert(h_pass);
1140 
1141       h_ref->Rebin(nrebin);
1142     }
1143 
1144     h_pass->SetTitle(TString(hist_name_prefix) + ": ntracks");
1145 
1146     DrawReference(h_pass, h_ref, false);
1147   }
1148 
1149   {
1150     static const int nrebin = 1;
1151     
1152     p = (TPad *)c1->cd(idx++);
1153     c1->Update();
1154     // p->SetLogx();
1155     p->SetGridy();
1156 
1157     TH1 *h_pass =
1158         (TH1 *)qa_file_new->GetObjectChecked(prefix + "ntracks_cuts", "TH1");
1159     assert(h_pass);
1160 
1161     h_pass->Rebin(nrebin);
1162 
1163     //    h_ratio->GetXaxis()->SetRangeUser(min_Et, max_Et);
1164     h_pass->GetYaxis()->SetTitle("Counts");
1165     // h_pass->GetYaxis()->SetRangeUser(-0, 1.);
1166 
1167     TH1 *h_ref = NULL;
1168     if (qa_file_ref) {
1169       h_ref =
1170           (TH1 *)qa_file_ref->GetObjectChecked(prefix + "ntracks_cuts", "TH1");
1171       assert(h_pass);
1172 
1173       h_ref->Rebin(nrebin);
1174     }
1175 
1176     h_pass->SetTitle(TString(hist_name_prefix) + ": ntracks (#geq 2 MVTX)");
1177 
1178     DrawReference(h_pass, h_ref, false);
1179   }
1180 
1181  // SaveCanvas(c1,
1182  //            TString(qa_file_name_new) + TString("_") + TString(c1->GetName()),
1183  //            true);
1184     
1185     c1->Draw();
1186     outfilef->cd();
1187     c1->Write();
1188 }
1189 
1190 {
1191     const char *hist_name_prefix = "QAG4SimulationVertex_SvtxVertexMap";
1192     TString prefix = TString("h_") + hist_name_prefix + TString("_");
1193     
1194     
1195   // X-direction
1196 
1197   TH2 *h_new = (TH2 *) qa_file_new->GetObjectChecked(
1198       prefix + TString("vxRes_gvz"), "TH2");
1199   assert(h_new);
1200 
1201   // h_new->Rebin2D(1, 5);
1202   //h_new->Sumw2();
1203   // h_new->GetXaxis()->SetRangeUser(-15,15);
1204   //  h_new->Scale(1. / Nevent_new);
1205 
1206   TH2 *h_ref = NULL;
1207   if (qa_file_ref)
1208   {
1209     h_ref = (TH2 *) qa_file_ref->GetObjectChecked(
1210         prefix + TString("vxRes_gvz"), "TH2");
1211     assert(h_ref);
1212 
1213     // h_ref->Rebin2D(1, 5);
1214     //h_ref->Sumw2();
1215     // h_ref->Scale(Nevent_new / Nevent_ref);
1216   }
1217 
1218   TCanvas *c1 = new TCanvas(TString("QA_Draw_Vertex_Resolution_x") + TString("_") + hist_name_prefix,
1219                             TString("QA_Draw_Vertex_Resolution_x") + TString("_") + hist_name_prefix,
1220                             1800, 1000);
1221   c1->Divide(2,1);
1222   int idx = 1;
1223   TPad *p;
1224 
1225   vector<pair<double, double>> gvz_ranges{
1226       {-10.0, 10.0}};
1227   TF1 *f1 = nullptr;
1228   TF1 *fit = nullptr;
1229   Double_t sigma = 0;
1230   Double_t sigma_unc = 0;
1231   char resstr[500];
1232   TLatex *res = nullptr;
1233   for (auto gvz_range : gvz_ranges)
1234   {
1235    // cout << __PRETTY_FUNCTION__ << " process " << gvz_range.first << " - " << gvz_range.second << " cm";
1236 
1237     p = (TPad *) c1->cd(idx++);
1238     c1->Update();
1239     // p->SetLogy();
1240 
1241     const double epsilon = 1e-6;
1242     const int bin_start = h_new->GetXaxis()->FindBin(gvz_range.first + epsilon);
1243     const int bin_end = h_new->GetXaxis()->FindBin(gvz_range.second - epsilon);
1244 
1245     TH1 *h_proj_new = h_new->ProjectionY(
1246         TString::Format(
1247             "%s_New_ProjX_%d_%d",
1248             h_new->GetName(), bin_start, bin_end));
1249     // bin_start, bin_end);
1250 
1251     h_proj_new->SetTitle(TString(hist_name_prefix) + TString::Format(
1252                                                          ": %.1f - %.1f cm - gvz", gvz_range.first, gvz_range.second));
1253     h_proj_new->GetXaxis()->SetTitle(TString::Format(
1254         "Vertex Resolution (x) [cm]"));
1255     h_proj_new->GetXaxis()->SetNdivisions(5,5);
1256     h_proj_new->GetXaxis()->SetRangeUser(-0.002,0.002);
1257 
1258     f1 = new TF1("f1","gaus",-.002,.002);
1259     h_proj_new->Fit(f1, "qm");
1260     sigma = f1->GetParameter(2);
1261     sigma_unc = f1->GetParError(2);
1262 
1263 
1264     TH1 *h_proj_ref = nullptr;
1265     if (h_ref)
1266     {
1267       h_proj_ref =
1268           h_ref->ProjectionY(
1269               TString::Format(
1270                   "%s_Ref_ProjX_%d_%d",
1271                   h_new->GetName(), bin_start, bin_end));
1272               // bin_start, bin_end);
1273       h_proj_ref->GetXaxis()->SetRangeUser(-10,10);
1274     }
1275     
1276     DrawReference(h_proj_new, h_proj_ref);
1277     sprintf(resstr,"#sigma = %.5f #pm %.5f cm", sigma, sigma_unc);
1278     res = new TLatex(0.325,0.825,resstr);
1279     res->SetNDC();
1280     res->SetTextSize(0.05);
1281     res->SetTextAlign(13);
1282     res->Draw();
1283   }
1284   p = (TPad *) c1->cd(idx++);
1285   c1->Update();
1286   gPad->SetLeftMargin(.2);
1287   //h_new->GetYaxis()->SetTitleOffset(2);
1288   h_new->Draw("colz");
1289   
1290 
1291 //  SaveCanvas(c1, TString(qa_file_name_new) + TString("_") + TString(c1->GetName()), true);
1292     c1->Draw();
1293 
1294   // Y-direction
1295 
1296   TH2 *h_new2 = (TH2 *) qa_file_new->GetObjectChecked(
1297       prefix + TString("vyRes_gvz"), "TH2");
1298   assert(h_new2);
1299 
1300   //  h_new->Rebin(1, 2);
1301   //h_new2->Sumw2();
1302   //  h_new->Scale(1. / Nevent_new);
1303 
1304   TH2 *h_ref2 = NULL;
1305   if (qa_file_ref)
1306   {
1307     h_ref2 = (TH2 *) qa_file_ref->GetObjectChecked(
1308         prefix + TString("vyRes_gvz"), "TH2");
1309     assert(h_ref2);
1310 
1311     // h_ref->Rebin(1, 2);
1312     //h_ref2->Sumw2();
1313     // h_ref->Scale(Nevent_new / Nevent_ref);
1314   }
1315 
1316   TCanvas *c2 = new TCanvas(TString("QA_Draw_Vertex_Resolution_y") + TString("_") + hist_name_prefix,
1317                             TString("QA_Draw_Vertex_Resolution_y") + TString("_") + hist_name_prefix,
1318                             1800, 1000);
1319   c2->Divide(2,1);
1320   int idx2 = 1;
1321   TPad *p2;
1322 
1323   vector<pair<double, double>> gvz_ranges2{
1324       {-10.0, 10.0}};
1325   TF1 *f2 = nullptr;
1326   TF1 *fit2 = nullptr;
1327   Double_t sigma2 = 0;
1328   Double_t sigma_unc2 = 0;
1329   char resstr2[500];
1330   TLatex *res2 = nullptr;
1331   for (auto gvz_range : gvz_ranges2)
1332   {
1333     //cout << __PRETTY_FUNCTION__ << " process " << gvz_range.first << " - " << gvz_range.second << " cm";
1334 
1335     p2 = (TPad *) c2->cd(idx2++);
1336     c2->Update();
1337     // p->SetLogy();
1338 
1339     const double epsilon = 1e-6;
1340     const int bin_start = h_new2->GetXaxis()->FindBin(gvz_range.first + epsilon);
1341     const int bin_end = h_new2->GetXaxis()->FindBin(gvz_range.second - epsilon);
1342 
1343     TH1 *h_proj_new2 = h_new2->ProjectionY(
1344         TString::Format(
1345             "%s_New_ProjX_%d_%d",
1346             h_new2->GetName(), bin_start, bin_end),
1347         bin_start, bin_end);
1348 
1349     h_proj_new2->SetTitle(TString(hist_name_prefix) + TString::Format(
1350                                                          ": %.1f - %.1f cm - gvz", gvz_range.first, gvz_range.second));
1351     h_proj_new2->GetXaxis()->SetTitle(TString::Format(
1352         "Vertex Resolution (y) [cm]"));
1353     h_proj_new2->GetXaxis()->SetNdivisions(5,5);
1354     h_proj_new2->GetXaxis()->SetRangeUser(-0.002,0.002);
1355     
1356     f2 = new TF1("f2","gaus",-.002,.002);
1357     h_proj_new2->Fit(f2,  "qm");
1358     fit2 = h_proj_new2->GetFunction("f2");
1359     sigma2 = fit2->GetParameter(2);
1360     sigma_unc2 = fit2->GetParError(2);
1361 
1362     TH1 *h_proj_ref2 = nullptr;
1363     if (h_ref2)
1364     {
1365       h_proj_ref2 =
1366           h_ref2->ProjectionY(
1367               TString::Format(
1368                   "%s_Ref_ProjX_%d_%d",
1369                   h_new2->GetName(), bin_start, bin_end),
1370               bin_start, bin_end);
1371     }
1372     
1373     DrawReference(h_proj_new2, h_proj_ref2);
1374     sprintf(resstr2,"#sigma = %.5f #pm %.5f cm", sigma2, sigma_unc2);
1375     res2 = new TLatex(0.325,0.825,resstr2);
1376     res2->SetNDC();
1377     res2->SetTextSize(0.05);
1378     res2->SetTextAlign(13);
1379     res2->Draw();
1380   }
1381   p2 = (TPad *) c2->cd(idx2++);
1382   c2->Update();
1383   gPad->SetLeftMargin(.2);
1384   //h_new2->GetYaxis()->SetTitleOffset(2);
1385   h_new2->Draw("colz");
1386 
1387   //SaveCanvas(c2, TString(qa_file_name_new) + TString("_") + TString(c2->GetName()), true);
1388     c2->Draw();
1389     
1390   // Z-direction
1391 
1392   TH2 *h_new3 = (TH2 *) qa_file_new->GetObjectChecked(
1393       prefix + TString("vzRes_gvz"), "TH2");
1394   assert(h_new3);
1395 
1396   //  h_new->Rebin(1, 2);
1397   //h_new3->Sumw2();
1398   //  h_new->Scale(1. / Nevent_new);
1399 
1400   TH2 *h_ref3 = NULL;
1401   if (qa_file_ref)
1402   {
1403     h_ref3 = (TH2 *) qa_file_ref->GetObjectChecked(
1404         prefix + TString("vzRes_gvz"), "TH2");
1405     assert(h_ref3);
1406 
1407     // h_ref->Rebin(1, 2);
1408     //h_ref3->Sumw2();
1409     // h_ref->Scale(Nevent_new / Nevent_ref);
1410   }
1411 
1412   TCanvas *c3 = new TCanvas(TString("QA_Draw_Vertex_Resolution_z") + TString("_") + hist_name_prefix,
1413                             TString("QA_Draw_Vertex_Resolution_z") + TString("_") + hist_name_prefix,
1414                             1800, 1000);
1415   c3->Divide(2,1);
1416   int idx3 = 1;
1417   TPad *p3;
1418 
1419   vector<pair<double, double>> gvz_ranges3{
1420       {-10.0, 10.0}};
1421   TF1 *f3 = nullptr;
1422   TF1 *fit3 = nullptr;
1423   Double_t sigma3 = 0;
1424   Double_t sigma_unc3 = 0;
1425   char resstr3[500];
1426   TLatex *res3 = nullptr;
1427   for (auto gvz_range : gvz_ranges3)
1428   {
1429    // cout << __PRETTY_FUNCTION__ << " process " << gvz_range.first << " - " << gvz_range.second << " cm";
1430 
1431     p3 = (TPad *) c3->cd(idx3++);
1432     c3->Update();
1433     // p->SetLogy();
1434 
1435     const double epsilon = 1e-6;
1436     const int bin_start = h_new3->GetXaxis()->FindBin(gvz_range.first + epsilon);
1437     const int bin_end = h_new3->GetXaxis()->FindBin(gvz_range.second - epsilon);
1438 
1439     TH1 *h_proj_new3 = h_new3->ProjectionY(
1440         TString::Format(
1441             "%s_New_ProjX_%d_%d",
1442             h_new3->GetName(), bin_start, bin_end),
1443         bin_start, bin_end);
1444 
1445     h_proj_new3->SetTitle(TString(hist_name_prefix) + TString::Format(
1446                                                          ": %.1f - %.1f cm -gvz", gvz_range.first, gvz_range.second));
1447     h_proj_new3->GetXaxis()->SetTitle(TString::Format(
1448         "Vertex Resolution (z) [cm]"));
1449     h_proj_new3->GetXaxis()->SetNdivisions(5,5);
1450     h_proj_new3->GetXaxis()->SetRangeUser(-0.002,0.002);
1451     
1452     f3 = new TF1("f3","gaus",-.002,.002);
1453     h_proj_new3->Fit(f3,  "qm");
1454     fit3 = h_proj_new3->GetFunction("f3");
1455     sigma3 = fit3->GetParameter(2);
1456     sigma_unc3 = fit3->GetParError(2);
1457 
1458     TH1 *h_proj_ref3 = nullptr;
1459     if (h_ref3)
1460     {
1461       h_proj_ref3 =
1462           h_ref3->ProjectionY(
1463               TString::Format(
1464                   "%s_Ref_ProjX_%d_%d",
1465                   h_new3->GetName(), bin_start, bin_end),
1466               bin_start, bin_end);
1467     }
1468     
1469     DrawReference(h_proj_new3, h_proj_ref3);
1470     sprintf(resstr3,"#sigma = %.5f #pm %.5f cm", sigma3, sigma_unc3);
1471     res3 = new TLatex(0.325,0.825,resstr3);
1472     res3->SetNDC();
1473     res3->SetTextSize(0.05);
1474     res3->SetTextAlign(13);
1475     res3->Draw();
1476   }
1477   p3 = (TPad *) c3->cd(idx3++);
1478   c3->Update();
1479   gPad->SetLeftMargin(.2);
1480   //h_new3->GetYaxis()->SetTitleOffset(2);
1481   h_new3->Draw("colz");
1482 
1483 //  SaveCanvas(c3, TString(qa_file_name_new) + TString("_") + TString(c3->GetName()), true);
1484     c3->Draw();
1485     outfilef->cd();
1486     c3->Write();
1487 }
1488 
1489 }