-
Eoin Clerkin authored
Decision to not use doxygen for licence headers. Removes doxygen formatting and file tag.
Eoin Clerkin authoredDecision to not use doxygen for licence headers. Removes doxygen formatting and file tag.
Draw_L1_histo.C 14.43 KiB
/* Copyright (C) 2010-2013 Frankfurt Institute for Advanced Studies, Goethe-Universität Frankfurt, Frankfurt
SPDX-License-Identifier: GPL-3.0-only
Authors: Igor Kulakov [committer] */
#include <unistd.h> // for dir navigation
const int textFont = 22; // TNewRoman
const bool divide = 0; // = 0 - each histo in separate file\screen. = 1 - all in one file\screen
void FitHistoGaus(TH1* hist, float* sigma = 0)
{
if (hist && (hist->GetEntries() != 0)) {
TF1* fit = new TF1("fit", "gaus");
fit->SetLineColor(2);
fit->SetLineWidth(3);
hist->Fit("fit", "", "", hist->GetXaxis()->GetXmin(), hist->GetXaxis()->GetXmax());
if (sigma) *sigma = fit->GetParameter(2);
}
else {
std::cout << " E: Read hists error! " << std::endl;
}
}
void FitHistoLine(TH1* hist)
{
if (hist && (hist->GetEntries() != 0)) {
TF1* fit = new TF1("fit", "[0]");
fit->SetLineColor(2);
fit->SetLineWidth(3);
hist->Fit("fit", "", "", hist->GetXaxis()->GetXmin(), hist->GetXaxis()->GetXmax());
}
else {
std::cout << " E: Read hists error! " << std::endl;
}
}
void makeUpHisto(TH1* hist, TString title)
{
if (hist && (hist->GetEntries() != 0)) {
hist->GetXaxis()->SetLabelFont(textFont);
hist->GetXaxis()->SetTitleFont(textFont);
hist->GetYaxis()->SetLabelFont(textFont);
hist->GetYaxis()->SetTitleFont(textFont);
hist->GetXaxis()->SetTitle(title);
hist->GetXaxis()->SetTitleOffset(1);
hist->GetYaxis()->SetTitle("Entries");
hist->GetYaxis()->SetTitleOffset(1.05); // good then entries per bit <= 9999
}
else {
std::cout << " E: Read hists error! " << std::endl;
}
}
void makeUpHisto(TH2* hist, TString title)
{
if (hist && (hist->GetEntries() != 0)) {
// hist->GetXaxis()->SetLabelFont(textFont);
// hist->GetXaxis()->SetTitleFont(textFont);
// hist->GetYaxis()->SetLabelFont(textFont);
// hist->GetYaxis()->SetTitleFont(textFont);
// hist->GetXaxis()->SetTitle(title);
// hist->GetXaxis()->SetTitleOffset(1);
// hist->GetYaxis()->SetTitle("Entries");
// hist->GetYaxis()->SetTitleOffset(1.05); // good then entries per bit <= 9999
hist->SetOption("CONTZ");
}
else {
std::cout << " E: Read hists error! " << std::endl;
}
}
void makeUpProfile(TProfile* prof, TString titleX, TString titleY)
{
if (prof && (prof->GetEntries() != 0)) {
prof->SetMarkerColor(kRed);
prof->SetLineColor(kRed);
prof->GetXaxis()->SetTitle(titleX);
prof->GetXaxis()->SetTitleOffset(1);
prof->GetYaxis()->SetTitle(titleY);
prof->GetYaxis()->SetTitleOffset(0.85);
prof->GetXaxis()->SetLabelFont(textFont);
prof->GetXaxis()->SetTitleFont(textFont);
prof->GetYaxis()->SetLabelFont(textFont);
prof->GetYaxis()->SetTitleFont(textFont);
}
else {
std::cout << " E: Read profile error! " << std::endl;
}
}
void Draw_L1_histo()
{
// ============================ Set Styles ============================
TStyle* histoStyle = new TStyle("histoStyle", "Plain Style(no colors/fill areas)");
// TStyle *histoStyle = gStyle;
histoStyle->SetTextFont(textFont);
histoStyle->SetPadColor(0);
histoStyle->SetCanvasColor(0);
histoStyle->SetTitleColor(0);
histoStyle->SetStatColor(0);
histoStyle->SetOptTitle(0); // without main up title
histoStyle->SetOptStat(1000111010);
// The parameter mode can be = ksiourmen (default = 000001111)
// k = 1; kurtosis printed
// k = 2; kurtosis and kurtosis error printed
// s = 1; skewness printed
// s = 2; skewness and skewness error printed
// i = 1; integral of bins printed
// o = 1; number of overflows printed
// u = 1; number of underflows printed
// r = 1; rms printed
// r = 2; rms and rms error printed
// m = 1; mean value printed
// m = 2; mean and mean error values printed
// e = 1; number of entries printed
// n = 1; name of histogram is printed
histoStyle->SetOptFit(10001);
// The parameter mode can be = pcev (default = 0111)
// p = 1; print Probability
// c = 1; print Chisquare/Number of degress of freedom
// e = 1; print errors (if e=1, v must be 1)
// v = 1; print name/values of parameters
histoStyle->SetStatW(0.175);
histoStyle->SetStatH(0.02);
histoStyle->SetStatX(0.95);
histoStyle->SetStatY(0.97);
histoStyle->SetStatFontSize(0.05);
histoStyle->SetStatFont(textFont);
histoStyle->SetPadBorderMode(0);
histoStyle->SetFrameBorderMode(0);
histoStyle->SetPalette(1); // for 2D
TStyle* profStyle = new TStyle("profStyle", "Plain Style(no colors/fill areas)");
// TStyle *profStyle = gStyle;
profStyle->SetTextFont(textFont);
profStyle->SetPadColor(0);
profStyle->SetCanvasColor(0);
profStyle->SetTitleColor(0);
profStyle->SetStatColor(0);
profStyle->SetOptTitle(0); // without main up title
profStyle->SetOptStat(0);
profStyle->SetOptFit(0);
profStyle->SetGridStyle(2);
profStyle->SetPadBorderMode(0);
profStyle->SetFrameBorderMode(0);
// ============================ Open file ============================
const TString fileName = "L1_histo.root";
TFile* fileIn = new TFile(fileName.Data(), "read");
FILE *ress, *pulls;
ress = fopen("resolution.txt", "w");
pulls = fopen("pull.txt", "w");
// ------------ make histos -----------
// ------------ tracks -----------
TString dirFitName = "/L1/Fit";
dirFitName = fileName + ":" + dirFitName;
TDirectory* dirFit = fileIn->Get(dirFitName);
const int nHisto = 40;
struct THistoData {
char* name; // one, which used in input file
char* title; // for X-axis
};
const THistoData histoData[nHisto] = {
{"fst_x", "Residual (x^{reco}-x^{mc}) [#mum]"} {"fst_y", "Residual (y^{reco}-y^{mc}) [#mum]"} {
"fst_tx", "Residual (t_{x}^{reco}-t_{x}^{mc}) [10^{-3}]"} {
"fst_ty", "Residual (t_{y}^{reco}-t_{y}^{mc}) [10^{-3}]"} {"fst_P", "Resolution (p^{reco}-p^{mc})/p^{mc}"} {
"fst_px", "Pull x"} {"fst_py", "Pull y"} {"fst_ptx", "Pull t_{x}"} {"fst_pty", "Pull t_{y}"} {"fst_pQP",
"Pull q/p"}
{"lst_pQP", "Pull q/p"} {"lst_ptx", "Pull t_{x}"} {"lst_pty", "Pull t_{y}"} {"lst_px", "Pull x"} {
"lst_py", "Pull y"} {"lst_P", "Resolution (p^{reco}-p^{mc})/p^{mc}"} {
"lst_tx", "Residual (t_{x}^{reco}-t_{x}^{mc}) [mrad]"} {"lst_ty", "Residual (t_{y}^{reco}-t_{y}^{mc}) [mrad]"} {
"lst_x", "Residual (x^{reco}-x^{mc}) [#mum]"} {"lst_y", "Residual (y^{reco}-y^{mc}) [#mum]"} {
"svrt_pQP", "Pull q/p"} {"svrt_ptx", "Pull t_{x}"} {"svrt_pty", "Pull t_{y}"} {"svrt_px", "Pull x"} {
"svrt_py", "Pull y"} {"svrt_P", "Resolution (p^{reco}-p^{mc})/p^{mc}"} {
"svrt_tx", "Residual (t_{x}^{reco}-t_{x}^{mc}) [mrad]"} {"svrt_ty", "Residual (t_{y}^{reco}-t_{y}^{mc}) [mrad]"} {
"svrt_x", "Residual (x^{reco}-x^{mc}) [#mum]"} {"svrt_y", "Residual (y^{reco}-y^{mc}) [#mum]"} {
"pvrt_pQP", "Pull q/p"} {"pvrt_ptx", "Pull t_{x}"} {"pvrt_pty", "Pull t_{y}"} {"pvrt_px", "Pull x"} {
"pvrt_py", "Pull y"} {"pvrt_P", "Resolution (p^{reco}-p^{mc})/p^{mc}"} {
"pvrt_tx", "Residual (t_{x}^{reco}-t_{x}^{mc}) [mrad]"} {"pvrt_ty", "Residual (t_{y}^{reco}-t_{y}^{mc}) [mrad]"} {
"pvrt_x", "Residual (x^{reco}-x^{mc}) [#mum]"} {"pvrt_y", "Residual (y^{reco}-y^{mc}) [#mum]"}};
TH1F* histo[nHisto];
float charat[nHisto][2]; // 0 - sigma, 1 - RMS
for (int i = 0; i < nHisto; i++) {
histo[i] = (TH1F*) dirFit->Get(histoData[i].name);
makeUpHisto(histo[i], histoData[i].title);
FitHistoGaus(histo[i], &(charat[i][0]));
charat[i][1] = histo[i]->GetRMS();
}
cout.setf(std::ios::scientific, std::ios::floatfield);
cout.precision(2);
cout << endl
<< "Name "
<< " \t"
<< " Sigma "
<< "\t"
<< " RMS " << endl;
for (int i = 0; i < nHisto; i++) {
cout << histoData[i].name << " \t" << charat[i][0] << "\t" << charat[i][1] << endl;
}
// ------------ hits, strips -----------
TString dirInputName = "/L1/Input";
dirInputName = fileName + ":" + dirInputName;
TDirectory* dirInput = fileIn->Get(dirInputName);
const int nHistoInput = 4;
const THistoData histoDataInput[nHisto] = {{"Px", "Pull x"} {"Py", "Pull y"} {
"x", "Residual (x^{hit}-x^{mc}) [#mum]"} {"y", "Residual (y^{hit}-y^{mc}) [#mum]"}};
TH1F* histoInput[nHistoInput];
for (int i = 0; i < nHistoInput; i++) {
histoInput[i] = (TH1F*) dirInput->Get(histoDataInput[i].name);
makeUpHisto(histoInput[i], histoDataInput[i].title);
FitHistoGaus(histoInput[i]);
}
// ------------ other histo -----------
TString dirOtherName = "/L1";
dirOtherName = fileName + ":" + dirOtherName;
TDirectory* dirOther = fileIn->Get(dirOtherName);
const int nHistoOther = 15;
const THistoData histoDataOther[nHisto] = {
{"h_acc_MCmom", "Momentum of tracks"} {"h_reco_MCmom", "Momentum of tracks"} {
"h_unreco_MCmom", "Momentum of tracks"} {"h_ghost_Rmom", "Momentum of tracks"} {
"h_acc_prim_MCmom", "Momentum of tracks"} {"h_reco_prim_MCmom", "Momentum of tracks"} {
"h_unreco_prim_MCmom", "Momentum of tracks"} {"h_ghost_Rmom", "Momentum of tracks"} {
"h_acc_sec_MCmom", "Momentum of tracks"} {"h_reco_sec_MCmom", "Momentum of tracks"} {
"h_unreco_sec_MCmom", "Momentum of tracks"} {"h_ghost_Rmom", "Momentum of tracks"}
// reg,acc,reco,ghost distribution always have to be first
{"h_reco_prob", "Prob of reco track"} {"h_ghost_prob", "Prob of ghost track"} {"h_rest_prob",
"Prob of reco rest track"}};
TH1F* histoOther[nHistoOther];
for (int i = 0; i < 12; i++) {
if (i % 4 != 2) histoOther[i] = (TH1F*) dirOther->Get(histoDataOther[i].name);
else {
histoOther[i] = new TH1F(histoDataOther[i].name, "Not reconstructed tracks", 100, 0.0, 5.0);
for (j = 0; j < histoOther[i]->GetNbinsX(); j++) {
histoOther[i]->SetBinContent(j, histoOther[i - 2]->GetBinContent(j) - histoOther[i - 1]->GetBinContent(j));
}
}
makeUpHisto(histoOther[i], histoDataOther[i].title); }
for (int i = 12; i < nHistoOther; i++) {
histoOther[i] = (TH1F*) dirOther->Get(histoDataOther[i].name);
makeUpHisto(histoOther[i], histoDataOther[i].title);
FitHistoLine(histoOther[i]);
}
// ------------ other 2D histo -----------
TString dirOther2DName = "/L1";
dirOther2DName = fileName + ":" + dirOther2DName;
TDirectory* dirOther2D = fileIn->Get(dirOther2DName);
const int nHistoOther2D = 1;
const THistoData histoDataOther2D[nHisto] = {
{"h2_reco_nhits_vs_mom_prim", "NHits and momentum for primary"}
// {"h2_reco_prob_ndf", "Prob and NDF" }
// {"h2_reco_prob_mom", "Prob and momentum" }
// {"h2_reco_prob_teta", "Prob and teta" }
// {"h2_teta_prob", "Prob" }
// {"h2_teta_probI", "ProbI" }
// {"h2_teta_probO", "ProbO" }
};
TH2F* histoOther2D[nHistoOther2D];
for (int i = 0; i < nHistoOther2D; i++) {
histoOther2D[i] = (TH2F*) dirOther2D->Get(histoDataOther2D[i].name);
makeUpHisto(histoOther2D[i], histoDataOther2D[i].title);
}
// ------------ make profiles -----------
TString dirName = "/L1";
dirName = fileName + ":" + dirName;
TDirectory* dir = fileIn->Get(dirName);
const int nProf = 4;
struct TProfData {
char* name; // one, which used in input file
char* titleX; // for X-axis
char* titleY; // for Y-axis
};
const TProfData profData[nProf] = {
{"p_eff_all_vs_mom", "Momentum [GeV/c]", "Efficiency [%]"} {"p_eff_all_vs_nhits", "NMCHits", "Efficiency [%]"} {
"p_eff_prim_vs_nhits", "NMCHits", "Efficiency [%]"} {"p_eff_sec_vs_nhits", "NMCHits", "Efficiency [%]"}};
TProfile* profile[nProf];
for (int i = 0; i < nProf; i++) {
profile[i] = (TProfile*) dir->Get(profData[i].name);
makeUpProfile(profile[i], profData[i].titleX, profData[i].titleY);
}
// ============================ Draw ============================
// TCanvas *c1;
// c1 = new TCanvas("c1","c1",0,0,1200,800);
// c1 -> Divide(nHisto);
//
// c1->UseCurrentStyle();
//
// for (int i = 0; i < nHisto; i++){
// c1->cd(i+1);
// histo[i]->Draw();
// }
// c1->SaveAs("histo.png");
system("mkdir L1_histo -p");
chdir("L1_histo");
histoStyle->cd();
TCanvas *c2, *c3;
c2 = new TCanvas("c2", "c2", 0, 0, 900, 600);
c3 = new TCanvas("c3", "c3", 0, 0, 900, 600);
c3->Divide(5, 2);
int iPart = 1;
for (int i = 0; i < 10; i++) {
c3->cd(iPart++);
histo[i]->Draw();
}
c3->SaveAs("fst.pdf");
c2->cd();
for (int i = 0; i < nHisto; i++) {
histo[i]->Draw();
TString name = TString(histoData[i].name) + ".pdf";
c2->SaveAs(name);
}
for (int i = 0; i < nHistoInput; i++) {
histoInput[i]->Draw();
TString name = TString(histoDataInput[i].name) + ".pdf";
c2->SaveAs(name);
}
// draw reg,acc,reco,ghost distribution
histoStyle->SetOptStat(000000000);
for (int j = 0, iH = 0; j < 3; j++) {
TLegend* leg1 = new TLegend(0.2 + 0.3, .62, 0.55 + 0.3, .75);
leg1->SetFillColor(0);
leg1->SetLineWidth(0);
leg1->SetBorderSize(0);
leg1->SetTextFont(22);
for (int i = 0; i < 4; i++, iH++) {
switch (i) {
case 0: histoOther[iH]->SetLineColor(kBlue); break;
case 1: histoOther[iH]->SetLineColor(kRed); break;
case 2:
histoOther[iH]->SetLineColor(kBlue);
histoOther[iH]->SetLineStyle(2);
break;
case 3: histoOther[iH]->SetLineColor(kBlack); break;
}
if (i == 0) histoOther[iH]->Draw();
else
histoOther[iH]->Draw("same");
leg1->AddEntry(histoOther[iH], histoOther[iH]->GetTitle());
}
leg1->Draw();
switch (j) {
case 0: c2->SaveAs("MultiplicityMom.pdf"); break;
case 1: c2->SaveAs("MultiplicityMomPrim.pdf"); break;
case 2: c2->SaveAs("MultiplicityMomSec.pdf"); break;
}
}
histoStyle->SetOptStat(1000111010);
for (int i = 12; i < nHistoOther; i++) {
histoOther[i]->Draw();
TString name = TString(histoDataOther[i].name) + ".pdf";
c2->SaveAs(name);
}
histoStyle->SetOptStat(1000000010);
for (int i = 0; i < nHistoOther2D; i++) { histoOther2D[i]->Draw();
TString name = TString(histoDataOther2D[i].name) + ".pdf";
c2->SaveAs(name);
}
histoStyle->SetOptStat(1000111010);
profStyle->cd();
c2->cd();
for (int i = 0; i < nProf; i++) {
c2->SetGridx();
c2->SetGridy();
profile[i]->Draw();
TString name = TString(profData[i].name) + ".pdf";
c2->SaveAs(name);
}
chdir("..");
}