Browse Source

plot styling

pull/2/head
Marius Pfeiffer 7 months ago
parent
commit
aeba3d41cd
  1. 65
      basic_analysis.h
  2. 28
      new_analysis_b02hphmmumu.cpp
  3. 57
      new_analysis_b02kppimmumu.cpp
  4. 24
      new_analysis_bu2hpmumu.cpp
  5. 11
      new_analysis_bu2kpmumu.cpp

65
basic_analysis.h

@ -151,6 +151,23 @@ void DrawInDefaultCanvas(RooPlot *rooPlot, const char *folder, double margin_lef
c->SaveAs(TString::Format("output_files/analysis/%s/%s.pdf", folder, name.Data()).Data());
}
void PrintStyles(RooPlot* plt) {
auto Xaxis = plt->GetXaxis();
auto Yaxis = plt->GetYaxis();
std::cout << "## Styles for plot " << plt->GetName() << std::endl;
std::cout << "## X axis\n"
<< "Label Size: " << Xaxis->GetLabelSize() << ",\n"
<< "Label Offset: " << Xaxis->GetLabelOffset() << ",\n"
<< "Title Size: " << Xaxis->GetTitleSize() << ",\n"
<< "Title Offset: " << Xaxis->GetTitleOffset() << "." << std::endl;
std::cout << "## Y axis\n"
<< "Label Size: " << Yaxis->GetLabelSize() << ",\n"
<< "Label Offset: " << Yaxis->GetLabelOffset() << ",\n"
<< "Title Size: " << Yaxis->GetTitleSize() << ",\n"
<< "Title Offset: " << Yaxis->GetTitleOffset() << "." << std::endl;
}
void DrawInDefaultCanvas(RooFitSummary fitSummary, const char *folder)
{
std::filesystem::create_directory(TString::Format("output_files/analysis/%s", folder).Data());
@ -161,29 +178,26 @@ void DrawInDefaultCanvas(RooFitSummary fitSummary, const char *folder)
c->GetPad(0)->GetPadPar(xlow, ylow, xup, yup);
c->Divide(1, 2);
Double_t upPad_ylow = ylow + 0.30 * (yup - ylow);
Double_t dwPad_yup = ylow + 0.30 * (yup - ylow);
Double_t upPad_ylow = ylow + 0.25 * (yup - ylow);
Double_t dwPad_yup = ylow + 0.25 * (yup - ylow);
TVirtualPad *upPad = c->GetPad(1);
upPad->SetPad(xlow, upPad_ylow, xup, yup);
TVirtualPad *dwPad = c->GetPad(2);
dwPad->SetPad(xlow, ylow, xup, dwPad_yup);
dwPad->SetTopMargin(0);
Double_t upPad_area = (yup - upPad_ylow) * (xup - xlow);
Double_t dwPad_area = (dwPad_yup - ylow) * (xup - xlow);
std::cout << "### UP AREA: " << upPad_area << " LOW AREA: " << dwPad_area << std::endl;
// std::cout << "### UP AREA: " << upPad_area << " LOW AREA: " << dwPad_area << std::endl;
const Double_t pull_title_size = 0.09;
const Double_t pull_label_size = 0.09;
const Double_t pull_title_offset_x = 3.0;
const Double_t pull_title_offset_y = 2.4;
Double_t fit_title_size = (pull_title_size * dwPad_area) / upPad_area;
Double_t fit_label_size = (pull_label_size * dwPad_area) / upPad_area;
Double_t fit_title_offset_x = (pull_title_offset_x * dwPad_area) / upPad_area;
Double_t fit_title_offset_y = (pull_title_offset_y * dwPad_area) / upPad_area;
// canvas->Update();
c->cd(1);
@ -201,13 +215,15 @@ void DrawInDefaultCanvas(RooFitSummary fitSummary, const char *folder)
auto fit_Xaxis = fitSummary.fit_histogram->GetXaxis();
auto fit_Yaxis = fitSummary.fit_histogram->GetYaxis();
fit_Xaxis->SetLabelOffset(0.02);
fit_Xaxis->SetLabelSize(fit_label_size);
fit_Xaxis->SetTitleOffset(fit_title_offset_x);
fit_Xaxis->SetTitleOffset(1.4);
fit_Xaxis->SetTitleSize(fit_title_size);
fit_Yaxis->SetLabelOffset(0.01);
fit_Yaxis->SetLabelSize(fit_label_size);
fit_Yaxis->SetTitleSize(fit_title_size);
fit_Yaxis->SetTitleOffset(fit_title_offset_y);
fit_Yaxis->SetTitleOffset(0);
fitSummary.fit_histogram->Draw();
@ -232,12 +248,10 @@ void DrawInDefaultCanvas(RooFitSummary fitSummary, const char *folder)
c->cd(2);
RooHist *pull_hist = (RooHist *)fitSummary.pull_histogram->getObject(0);
Double_t pull_min = fitSummary.pull_histogram->GetMinimum();
Double_t pull_max = fitSummary.pull_histogram->GetMaximum();
Double_t pull_min = pull_hist->GetMinimum();
Double_t pull_max = pull_hist->GetMaximum();
std::cout << "### (" << fitSummary.pull_histogram->GetName() << ") PULL MIN: " << pull_min << " PULL MAX: " << pull_max << std::endl;
// std::cout << "### (" << fitSummary.pull_histogram->GetName() << ") PULL MIN: " << pull_min << " PULL MAX: " << pull_max << std::endl;
double bound = 0;
if (TMath::Abs(pull_min) > TMath::Abs(pull_max))
@ -251,26 +265,35 @@ void DrawInDefaultCanvas(RooFitSummary fitSummary, const char *folder)
auto pull_Xaxis = fitSummary.pull_histogram->GetXaxis();
auto pull_Yaxis = fitSummary.pull_histogram->GetYaxis();
pull_Xaxis->SetLabelOffset(0.02);
pull_Xaxis->SetLabelSize(pull_label_size);
pull_Xaxis->SetTitleOffset(pull_title_offset_x);
pull_Xaxis->SetTitleSize(pull_title_size);
pull_Xaxis->SetTitle("");
pull_Yaxis->SetLabelOffset(0.01);
pull_Yaxis->SetLabelSize(pull_label_size);
pull_Yaxis->SetTitleSize(pull_title_size);
pull_Yaxis->SetTitleOffset(pull_title_offset_y);
pull_Yaxis->SetTitleOffset(0.45);
pull_Yaxis->SetTitle("Pull Distribution");
fitSummary.pull_histogram->SetTitle("");
fitSummary.pull_histogram->SetMaximum(bound);
fitSummary.pull_histogram->SetMinimum(-bound);
fitSummary.pull_histogram->Draw();
auto line_mid = new TLine(B_MASS_VAR_MIN, 0, B_MASS_VAR_MAX, 0);
line_mid->SetLineColor(kOrange+4);
line_mid->Draw();
auto line_up = new TLine(B_MASS_VAR_MIN, bound, B_MASS_VAR_MAX, bound);
auto line_up = new TLine(B_MASS_VAR_MIN, bound / 2, B_MASS_VAR_MAX, bound / 2);
line_up->SetLineStyle(kDashed);
line_up->SetLineColor(kOrange+4);
line_up->Draw();
auto line_low = new TLine(B_MASS_VAR_MIN, -bound, B_MASS_VAR_MAX, -bound);
auto line_low = new TLine(B_MASS_VAR_MIN, -bound / 2, B_MASS_VAR_MAX, -bound / 2);
line_low->SetLineStyle(kDashed);
line_low->SetLineColor(kOrange+4);
line_low->Draw();
c->Draw();
@ -400,10 +423,10 @@ RooFitSummary CreateRooFitHistogramAndFitCB(TH1D *hist, bool hasExpBkg, bool use
RooHist *pull_hist = roo_frame_mass->pullHist(hist_name, pull_compare_name);
roo_frame_pull->addPlotable(pull_hist, "P");
Double_t pull_min = pull_hist->GetMinimum();
Double_t pull_max = pull_hist->GetMaximum();
// Double_t pull_min = pull_hist->GetMinimum();
// Double_t pull_max = pull_hist->GetMaximum();
std::cout << "##### (" << pull_hist->GetName() << ") PULL MIN: " << pull_min << " PULL MAX: " << pull_max << std::endl;
// std::cout << "##### (" << pull_hist->GetName() << ") PULL MIN: " << pull_min << " PULL MAX: " << pull_max << std::endl;
return RooFitSummary{
roo_frame_mass,

28
new_analysis_b02hphmmumu.cpp

@ -47,8 +47,8 @@ int new_analysis_b02hphmmumu()
const char *analysis_name = "B0ToHpHmMuMu";
const char *data_tree_name = "B0ToHpHmMuMu";
const char *sim_tree_name = "B0ToHpHmMuMu_noPID";
const char *end_state_mass_literal = "m(#pi^{+}_{(#rightarrow K^{+})}#pi^{-}_{(#rightarrow K^{-})}#mu^{+}#mu^{-})";
const bool retrain_bdt = true;
const char *end_state_mass_literal = "m(#pi^{+}#pi^{-}_{(#rightarrow K^{-})}#mu^{+}#mu^{-} & #pi^{+}_{(#rightarrow K^{+})}#pi^{-}#mu^{+}#mu^{-})";
const bool retrain_bdt = false;
TChain *data_chain = new TChain(TString::Format("%s/DecayTree", data_tree_name));
// data_chain->Add("/auto/data/pfeiffer/inclusive_detached_dilepton/data_samples/spruce_magdown_2023_v0_tuple_90000000_v0r0p6288631.root");
@ -88,6 +88,7 @@ int new_analysis_b02hphmmumu()
h1_B_Mass_jpsi->GetXaxis()->SetTitle(end_state_mass_literal);
h1_B_Mass_psi2s->GetXaxis()->SetTitle(end_state_mass_literal);
h1_B_Mass_sim->GetXaxis()->SetTitle(end_state_mass_literal);
h2_Hlt1_flags_B_Mass->GetXaxis()->SetTitle(end_state_mass_literal);
h2_Hlt1_flags_excl_B_Mass->GetXaxis()->SetTitle(end_state_mass_literal);
@ -119,9 +120,10 @@ int new_analysis_b02hphmmumu()
ConnectVarsToData(vars, data_chain, sim_chain, sig_tree, bkg_tree);
unsigned int data_entries = data_chain->GetEntries();
unsigned int sim_entries = 100000; // sim_chain->GetEntries();
unsigned int sim_entries = sim_chain->GetEntries();
unsigned int bkg_events = 0;
unsigned int sig_events = 0;
if (retrain_bdt)
{
@ -159,10 +161,12 @@ int new_analysis_b02hphmmumu()
PrintProgress(TString::Format("%s BKG Collection", analysis_name), data_entries, 10000, i);
}
std::cout << "# Added " << bkg_events << " background events." << std::endl;
}
else
{
std::cout << "# Starting without BDT retrain." << std::endl;
bkg_events = data_entries;
}
for (unsigned int i = 0; i < sim_entries; i++)
@ -188,15 +192,17 @@ int new_analysis_b02hphmmumu()
{
Double_t reconstructed_B_Mass = (reconstructed_Kstar + l14v_sim->LorentzVector() + l24v_sim->LorentzVector()).M();
if (std::all_of(vars.begin(), vars.end(), [](TV *v)
{ return v->IsMCFinite(); }))
if (sig_events < bkg_events)
{
if (retrain_bdt)
if (retrain_bdt && std::all_of(vars.begin(), vars.end(), [](TV *v)
{ return v->IsMCFinite(); }))
{
sig_tree->Fill();
sig_events++;
}
h1_B_Mass_sim->Fill(reconstructed_B_Mass);
}
h1_B_Mass_sim->Fill(reconstructed_B_Mass);
}
}
@ -214,7 +220,7 @@ int new_analysis_b02hphmmumu()
Float_t *train_vars = new Float_t[vars.size()];
auto reader = SetupReader(vars, train_vars, analysis_name);
const double mva_cut_value = -0.12;
const double mva_cut_value = -0.0508;
for (unsigned int i = 0; i < data_entries; i++)
{
@ -265,8 +271,8 @@ int new_analysis_b02hphmmumu()
if (mva_response > mva_cut_value)
{
// if (TMath::Abs(reconstructed_Kstar.M() - KSTAR_MASS) < 100.)
// {
if (TMath::Abs(reconstructed_Kstar.M() - KSTAR_MASS) < 50.)
{
if (TMath::Abs(dimuon.M() - JPSI_MASS) < 100.)
{
h1_B_Mass_jpsi->Fill(reconstructed_B_Mass);
@ -275,7 +281,7 @@ int new_analysis_b02hphmmumu()
{
h1_B_Mass_psi2s->Fill(reconstructed_B_Mass);
}
// }
}
}
}

57
new_analysis_b02kppimmumu.cpp

@ -53,10 +53,14 @@ int new_analysis_b02kppimmumu()
TChain *data_chain = new TChain(TString::Format("%s/DecayTree", data_tree_name));
data_chain->Add("/auto/data/pfeiffer/inclusive_detached_dilepton/data_samples/Collision23_Beam6800GeV-VeloClosed-MagDown-Excl-UT_RealData_SprucingPass23r1_94000000_RD.root");
Int_t Hp_Q, Hm_Q;
FourVect *l14v_data = FourVect::Init(data_chain, "muplus");
FourVect *l24v_data = FourVect::Init(data_chain, "muminus");
FourVect *hp4v_data = FourVect::Init(data_chain, "Kplus");
FourVect *hm4v_data = FourVect::Init(data_chain, "piminus");
data_chain->SetBranchAddress("Kplus_Q", &Hp_Q);
data_chain->SetBranchAddress("piminus_Q", &Hm_Q);
TChain *sim_chain = new TChain(TString::Format("%s/DecayTree", sim_tree_name));
sim_chain->Add("/auto/data/pfeiffer/inclusive_detached_dilepton/MC/rd_btoxll_simulation_turbo_v0r0p6657752_B0ToKpPimMuMu_11144002_magdown.root");
@ -82,6 +86,7 @@ int new_analysis_b02kppimmumu()
h1_B_Mass_jpsi->GetXaxis()->SetTitle(end_state_mass_literal);
h1_B_Mass_psi2s->GetXaxis()->SetTitle(end_state_mass_literal);
h1_B_Mass_sim->GetXaxis()->SetTitle(end_state_mass_literal);
h2_Hlt1_flags_B_Mass->GetXaxis()->SetTitle(end_state_mass_literal);
h2_Hlt1_flags_excl_B_Mass->GetXaxis()->SetTitle(end_state_mass_literal);
@ -106,6 +111,11 @@ int new_analysis_b02kppimmumu()
std::cout << "# Starting evaluation of data." << std::endl;
int kplus = 0;
int kminus = 0;
int piplus = 0;
int piminus = 0;
for (unsigned int i = 0; i < data_entries; i++)
{
data_chain->GetEntry(i);
@ -120,15 +130,36 @@ int new_analysis_b02kppimmumu()
FillHlt1DecisionHistos(hlt1_decision_histos, reconstructed_B_Mass);
}
if (TMath::Abs(reconstructed_Kstar.M() - KSTAR_MASS) < 50.)
if (B_MASS_VAR_MIN <= reconstructed_B_Mass && reconstructed_B_Mass <= B_MASS_VAR_MAX)
{
if (TMath::Abs(dimuon.M() - JPSI_MASS) < 100.)
if (TMath::Abs(reconstructed_Kstar.M() - KSTAR_MASS) < 50.)
{
h1_B_Mass_jpsi->Fill(reconstructed_B_Mass);
}
else if (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.)
{
h1_B_Mass_psi2s->Fill(reconstructed_B_Mass);
if (TMath::Abs(dimuon.M() - JPSI_MASS) < 100.)
{
if (Hp_Q == 1)
{
kplus++;
}
else if (Hp_Q == -1)
{
kminus++;
}
if (Hm_Q == 1)
{
piplus++;
}
else if (Hm_Q == -1)
{
piminus++;
}
h1_B_Mass_jpsi->Fill(reconstructed_B_Mass);
}
else if (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.)
{
h1_B_Mass_psi2s->Fill(reconstructed_B_Mass);
}
}
}
@ -141,10 +172,10 @@ int new_analysis_b02kppimmumu()
std::cout << line << ": " << hits << std::endl;
}
DrawInDefaultCanvas(h2_Hlt1_flags_B_Mass, analysis_name, 0.16, "COLZ");
DrawInDefaultCanvas(h2_Hlt1_flags_excl_B_Mass, analysis_name, 0.16, "COLZ");
DrawInDefaultCanvas(h1_B_Mass_jpsi, analysis_name, 0.1);
DrawInDefaultCanvas(h1_B_Mass_psi2s, analysis_name, 0.1);
// DrawInDefaultCanvas(h2_Hlt1_flags_B_Mass, analysis_name, 0.16, "COLZ");
// DrawInDefaultCanvas(h2_Hlt1_flags_excl_B_Mass, analysis_name, 0.16, "COLZ");
// DrawInDefaultCanvas(h1_B_Mass_jpsi, analysis_name, 0.1);
// DrawInDefaultCanvas(h1_B_Mass_psi2s, analysis_name, 0.1);
auto roofit_hist_sim = CreateRooFitHistogramAndFitCB(h1_B_Mass_sim, false, false, ShapeParamters{});
auto roofit_hist_jpsi_fitsum = CreateRooFitHistogramAndFitCB(h1_B_Mass_jpsi, true, true, roofit_hist_sim.shape_parameters);
@ -154,7 +185,9 @@ int new_analysis_b02kppimmumu()
DrawInDefaultCanvas(roofit_hist_psi2s_fitsum, analysis_name);
DrawInDefaultCanvas(roofit_hist_sim, analysis_name);
DrawHlt1DecisionHistos(analysis_name, hlt1_decision_histos);
// DrawHlt1DecisionHistos(analysis_name, hlt1_decision_histos);
std::cout << "hist entries: " << h1_B_Mass_jpsi->GetEntries() << ", kplus: " << kplus << ", kminus: " << kminus << ", piplus: " << piplus << ", piminus: " << piminus << std::endl;
return 0;
}

24
new_analysis_bu2hpmumu.cpp

@ -48,7 +48,7 @@ int new_analysis_bu2hpmumu()
const char *data_tree_name = "BuToHpMuMu";
const char *sim_tree_name = "BuToHpMuMu_noPID";
const char *end_state_mass_literal = "m(#pi^{+}_{(#rightarrow K^{+})}#mu^{+}#mu^{-})";
const bool retrain_bdt = false;
const bool retrain_bdt = true;
TChain *data_chain = new TChain(TString::Format("%s/DecayTree", data_tree_name));
// data_chain->Add("/auto/data/pfeiffer/inclusive_detached_dilepton/data_samples/spruce_magdown_2023_v0_tuple_90000000_v0r0p6288631.root");
@ -81,6 +81,7 @@ int new_analysis_bu2hpmumu()
h1_B_Mass_jpsi->GetXaxis()->SetTitle(end_state_mass_literal);
h1_B_Mass_psi2s->GetXaxis()->SetTitle(end_state_mass_literal);
h1_B_Mass_sim->GetXaxis()->SetTitle(end_state_mass_literal);
h2_Hlt1_flags_B_Mass->GetXaxis()->SetTitle(end_state_mass_literal);
h2_Hlt1_flags_excl_B_Mass->GetXaxis()->SetTitle(end_state_mass_literal);
@ -108,8 +109,8 @@ int new_analysis_bu2hpmumu()
ConnectVarsToData(vars, data_chain, sim_chain, sig_tree, bkg_tree);
unsigned int data_entries = 100000;// data_chain->GetEntries();
unsigned int sim_entries = 100000; //sim_chain->GetEntries();
unsigned int data_entries = data_chain->GetEntries();
unsigned int sim_entries = sim_chain->GetEntries();
std::cout << "# Got " << data_entries << " data and " << sim_entries << " simulated events." << std::endl;
@ -150,22 +151,19 @@ int new_analysis_bu2hpmumu()
Double_t reconstructed_B_Mass = (hp4v_sim->LorentzVector(K_MASS) + l14v_sim->LorentzVector() + l24v_sim->LorentzVector()).M();
if (B_BKGCAT == 30 && TMath::Abs(L1_TRUEID) == PID_MUON && L2_TRUEID == -L1_TRUEID && TMath::Abs(Hp_TRUEID) == PID_KAON)
{
if (retrain_bdt && std::all_of(vars.begin(), vars.end(), [](TV *v)
{ return v->IsMCFinite(); }))
if (sig_events < bkg_events)
{
sig_tree->Fill();
if (retrain_bdt && std::all_of(vars.begin(), vars.end(), [](TV *v)
{ return v->IsMCFinite(); }))
{
sig_tree->Fill();
sig_events++;
}
}
sig_events++;
h1_B_Mass_sim->Fill(reconstructed_B_Mass);
}
if (sig_events >= bkg_events)
{
break;
}
PrintProgress(TString::Format("%s SIG Collection", analysis_name), sim_entries, 10000, i);
}

11
new_analysis_bu2kpmumu.cpp

@ -79,6 +79,7 @@ int new_analysis_bu2kpmumu()
h1_B_Mass_jpsi->GetXaxis()->SetTitle(end_state_mass_literal);
h1_B_Mass_psi2s->GetXaxis()->SetTitle(end_state_mass_literal);
h1_B_Mass_sim->GetXaxis()->SetTitle(end_state_mass_literal);
h2_Hlt1_flags_B_Mass->GetXaxis()->SetTitle(end_state_mass_literal);
h2_Hlt1_flags_excl_B_Mass->GetXaxis()->SetTitle(end_state_mass_literal);
@ -134,10 +135,10 @@ int new_analysis_bu2kpmumu()
std::cout << line << ": " << hits << std::endl;
}
DrawInDefaultCanvas(h2_Hlt1_flags_B_Mass, analysis_name, 0.16, "COLZ");
DrawInDefaultCanvas(h2_Hlt1_flags_excl_B_Mass, analysis_name, 0.16, "COLZ");
DrawInDefaultCanvas(h1_B_Mass_jpsi, analysis_name, 0.1);
DrawInDefaultCanvas(h1_B_Mass_psi2s, analysis_name, 0.1);
// DrawInDefaultCanvas(h2_Hlt1_flags_B_Mass, analysis_name, 0.16, "COLZ");
// DrawInDefaultCanvas(h2_Hlt1_flags_excl_B_Mass, analysis_name, 0.16, "COLZ");
// DrawInDefaultCanvas(h1_B_Mass_jpsi, analysis_name, 0.1);
// DrawInDefaultCanvas(h1_B_Mass_psi2s, analysis_name, 0.1);
auto roofit_hist_sim = CreateRooFitHistogramAndFitCB(h1_B_Mass_sim, false, false, ShapeParamters{});
auto roofit_hist_jpsi_fitsum = CreateRooFitHistogramAndFitCB(h1_B_Mass_jpsi, true, true, roofit_hist_sim.shape_parameters);
@ -147,7 +148,7 @@ int new_analysis_bu2kpmumu()
DrawInDefaultCanvas(roofit_hist_psi2s_fitsum, analysis_name);
DrawInDefaultCanvas(roofit_hist_sim, analysis_name);
DrawHlt1DecisionHistos(analysis_name, hlt1_decision_histos);
// DrawHlt1DecisionHistos(analysis_name, hlt1_decision_histos);
return 0;
}
Loading…
Cancel
Save