Skip to content

Commit

Permalink
Refactor validation code (#431)
Browse files Browse the repository at this point in the history 
* Refactor utilities

* Add copyright. Add tuple header. Fix path.

* Add header guards

* Move compare.py

* Factor out plotting part of validation

* Fix plotting factorisation

* Factor out MakePlots into a header file.

* Fix casting

* Break line

* Rename headers

* Fix header guards

* Improve logging

* Improve legends in validation plots

* Fix JE validation plots and factor out code

* Remove link to utils_plot.h

* Improve y-range calculation and plotting, add protections

* Improve AliPhysics task macros

* Add documentation comments in AliPhysics macros
  • Loading branch information
@vkucera
vkucera authored Nov 22, 2023
1 parent 390f6f1 commit a01425c
Show file tree
Hide file tree
Showing 13 changed files with 493 additions and 523 deletions.
184 changes: 15 additions & 169 deletions codeHF/Compare.C
View file
Original file line number Diff line number Diff line change
@@ -1,42 +1,21 @@
// Comparison of AliPhysics and O2 histograms

#include "utils_plot.h"
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.

// vectors of histogram specifications
using VecSpecHis = std::vector<std::tuple<TString, TString, TString, int, bool, bool, TString>>;
// Comparison of AliPhysics and O2 histograms

// Add histogram specification in the vector.
void AddHistogram(VecSpecHis& vec, TString label, TString nameAli, TString nameO2, int rebin, bool logH, bool logR, TString proj = "x")
{
vec.push_back(std::make_tuple(label, nameAli, nameO2, rebin, logH, logR, proj));
}
#include "../exec/utilitiesValidation.h"

Int_t Compare(TString fileO2 = "AnalysisResults_O2.root", TString fileAli = "AnalysisResults_ALI.root", TString options = "", bool doRatio = false)
Int_t Compare(TString pathFileO2 = "AnalysisResults_O2.root", TString pathFileAli = "AnalysisResults_ALI.root", TString options = "", bool doRatio = false)
{
gStyle->SetOptStat(0);
gStyle->SetPalette(0);
gStyle->SetCanvasColor(0);
gStyle->SetFrameFillColor(0);

TFile* fO2 = new TFile(fileO2.Data());
if (fO2->IsZombie()) {
printf("Failed to open file %s\n", fileO2.Data());
return 1;
}
TFile* fAli = new TFile(fileAli.Data());
if (fAli->IsZombie()) {
printf("Failed to open file %s\n", fileAli.Data());
return 1;
}

TString pathListAli = "HFVertices/clistHFVertices";
TList* lAli = nullptr;
fAli->GetObject(pathListAli.Data(), lAli);
if (!lAli) {
printf("Failed to load list %s from %s\n", pathListAli.Data(), fileAli.Data());
return 1;
}

TString labelParticle = "";

// Histogram specification: axis label, AliPhysics name, O2Physics path/name, rebin, log scale histogram, log scale ratio, projection axis
Expand Down Expand Up @@ -251,7 +230,7 @@ Int_t Compare(TString fileO2 = "AnalysisResults_O2.root", TString fileAli = "Ana
AddHistogram(vecHisJetSubstructureMC, "#it{n}_{SD,gen}", "fHistJetNsd_Part", "jet-substructure-hf-mcp/h_jet_nsd", 1, 0, 0);

// vector of specifications of vectors: name, VecSpecHis, pads X, pads Y
std::vector<std::tuple<TString, VecSpecHis, int, int>> vecSpecVecSpec;
VecSpecVecSpec vecSpecVecSpec;

// Add vector specifications in the vector.
if (options.Contains(" events "))
Expand Down Expand Up @@ -295,138 +274,5 @@ Int_t Compare(TString fileO2 = "AnalysisResults_O2.root", TString fileAli = "Ana
if (options.Contains(" jets-substructure-mc "))
vecSpecVecSpec.push_back(std::make_tuple("jets-substructure-mc", vecHisJetSubstructureMC, 5, 3));

// Histogram plot vertical margins
Float_t marginHigh = 0.05;
Float_t marginLow = 0.05;
bool logScaleH = false;
// Ratio plot vertical margins
Float_t marginRHigh = 0.05;
Float_t marginRLow = 0.05;
bool logScaleR = false;
Float_t yMin, yMax;
Int_t nAli, nO2, rebin;

TH1F* hAli = nullptr;
TH1D* hO2 = nullptr;
TH1F* hRatio = nullptr;
TString labelAxis = "";
TString nameHisAli = "";
TString nameHisO2 = "";
TString projAx = "";
TCanvas* canHis = nullptr;
TCanvas* canRat = nullptr;

// loop over lists
for (const auto& specVecSpec : vecSpecVecSpec) {
auto nameSpec = std::get<0>(specVecSpec); // list name
auto vecSpec = std::get<1>(specVecSpec); // list of histogram specs.
int nPadsX = std::get<2>(specVecSpec); // number of horizontal pads
int nPadsY = std::get<3>(specVecSpec); // number of vertical pads
Printf("\nProcessing histogram list: %s (%d)", nameSpec.Data(), (int)vecSpec.size());
if (nPadsX * nPadsY < vecSpec.size()) {
Printf("Not enough pads (%d)", nPadsX * nPadsY);
return 1;
}

canHis = new TCanvas(Form("canHis_%s", nameSpec.Data()), Form("Histos_%s", nameSpec.Data()), 3000, 1600);
SetCanvas(canHis, nPadsX, nPadsY);
if (doRatio) {
canRat = new TCanvas(Form("canRat_%s", nameSpec.Data()), Form("Ratios_%s", nameSpec.Data()), 3000, 1600);
SetCanvas(canRat, nPadsX, nPadsY);
}

// loop over histograms
for (int index = 0; index < vecSpec.size(); index++) {
auto spec = vecSpec[index];
labelAxis = std::get<0>(spec);
nameHisAli = std::get<1>(spec);
nameHisO2 = std::get<2>(spec);
rebin = std::get<3>(spec);
logScaleH = std::get<4>(spec);
logScaleR = std::get<5>(spec);
projAx = std::get<6>(spec);

// Get AliPhysics histogram.
hAli = (TH1F*)lAli->FindObject(nameHisAli.Data());
if (!hAli) {
printf("Failed to load %s from %s\n", nameHisAli.Data(), fileAli.Data());
return 1;
}

// Get O2 histogram.
auto oO2 = fO2->Get(nameHisO2.Data());
if (!oO2) {
printf("Failed to load %s from %s\n", nameHisO2.Data(), fileO2.Data());
return 1;
}

if (oO2->InheritsFrom("TH3")) {
if (projAx == "x") {
hO2 = ((TH3D*)oO2)->ProjectionX();
} else if (projAx == "y") {
hO2 = ((TH3D*)oO2)->ProjectionY();
}
} else if (oO2->InheritsFrom("TH2")) {
if (projAx == "x") {
hO2 = ((TH2D*)oO2)->ProjectionX();
} else if (projAx == "y") {
hO2 = ((TH2D*)oO2)->ProjectionY();
}
} else {
hO2 = (TH1D*)oO2;
}

Printf("%d (%s, %s): bins: %d, %d, ranges: %g-%g, %g-%g",
index, nameHisAli.Data(), nameHisO2.Data(),
hAli->GetNbinsX(), hO2->GetNbinsX(),
hAli->GetXaxis()->GetBinLowEdge(1), hAli->GetXaxis()->GetBinUpEdge(hAli->GetNbinsX()),
hO2->GetXaxis()->GetBinLowEdge(1), hO2->GetXaxis()->GetBinUpEdge(hO2->GetNbinsX()));

nAli = hAli->GetEntries();
nO2 = hO2->GetEntries();

// Histograms
auto padH = canHis->cd(index + 1);
hAli->Rebin(rebin);
hO2->Rebin(rebin);
hAli->SetLineColor(1);
hAli->SetLineWidth(2);
hO2->SetLineColor(2);
hO2->SetLineWidth(1);
hAli->SetTitle(Form("Entries: Ali: %d, O^{2}: %d;%s;Entries", nAli, nO2, labelAxis.Data()));
hAli->GetYaxis()->SetMaxDigits(3);
yMin = TMath::Min(hO2->GetMinimum(0), hAli->GetMinimum(0));
yMax = TMath::Max(hO2->GetMaximum(), hAli->GetMaximum());
SetHistogram(hAli, yMin, yMax, marginLow, marginHigh, logScaleH);
SetPad(padH, logScaleH);
hAli->Draw();
hO2->Draw("same");
TLegend* legend = new TLegend(0.8, 0.72, 1., 0.92);
legend->AddEntry(hAli, "Ali", "L");
legend->AddEntry(hO2, "O^{2}", "L");
legend->Draw();

// Ratio
if (doRatio) {
auto padR = canRat->cd(index + 1);
hRatio = (TH1F*)hO2->Clone(Form("hRatio%d", index));
hRatio->Divide(hAli);
hRatio->SetTitle(Form("Entries ratio: %g;%s;O^{2}/Ali", (double)nO2 / (double)nAli, labelAxis.Data()));
yMin = hRatio->GetMinimum(0);
yMax = hRatio->GetMaximum();
SetHistogram(hRatio, yMin, yMax, marginRLow, marginRHigh, logScaleR);
SetPad(padR, logScaleR);
hRatio->Draw();
}
}
canHis->SaveAs(Form("comparison_histos_%s.pdf", nameSpec.Data()));
canHis->SaveAs(Form("comparison_histos_%s.png", nameSpec.Data()));
if (doRatio) {
canRat->SaveAs(Form("comparison_ratios_%s.pdf", nameSpec.Data()));
canRat->SaveAs(Form("comparison_ratios_%s.png", nameSpec.Data()));
}
delete canHis;
delete canRat;
}
return 0;
return MakePlots(vecSpecVecSpec, pathFileO2, pathFileAli, pathListAli, doRatio);
}
21 changes: 16 additions & 5 deletions codeHF/PlotEfficiency.C
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,17 @@
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.

// Plotting of reconstruction efficiency

#include "utils_plot.h"
#include "../exec/utilitiesPlot.h"

Int_t PlotEfficiency(TString pathFile = "AnalysisResults.root", TString particles = "d0")
{
Expand Down Expand Up @@ -88,12 +99,12 @@ Int_t PlotEfficiency(TString pathFile = "AnalysisResults.root", TString particle
bool okPrompt = true;
TH1F* hPtRecPrompt = (TH1F*)file->Get(nameHistRec.Data());
if (!hPtRecPrompt) {
Printf("Warning: Failed to load %s from %s", nameHistRec.Data(), pathFile.Data());
Warning("PlotEfficiency", "Failed to load %s from %s", nameHistRec.Data(), pathFile.Data());
okPrompt = false;
}
TH1F* hPtGenPrompt = (TH1F*)file->Get(nameHistgen.Data());
if (!hPtGenPrompt) {
Printf("Warning: Failed to load %s from %s", nameHistgen.Data(), pathFile.Data());
Warning("PlotEfficiency", "Failed to load %s from %s", nameHistgen.Data(), pathFile.Data());
okPrompt = false;
}

Expand All @@ -109,12 +120,12 @@ Int_t PlotEfficiency(TString pathFile = "AnalysisResults.root", TString particle
bool okNonPrompt = true;
TH1F* hPtRecNonPrompt = (TH1F*)file->Get(nameHistRec.Data());
if (!hPtRecNonPrompt) {
Printf("Warning: Failed to load %s from %s", nameHistRec.Data(), pathFile.Data());
Warning("PlotEfficiency", "Failed to load %s from %s", nameHistRec.Data(), pathFile.Data());
okNonPrompt = false;
}
TH1F* hPtGenNonPrompt = (TH1F*)file->Get(nameHistgen.Data());
if (!hPtGenNonPrompt) {
Printf("Warning: Failed to load %s from %s", nameHistgen.Data(), pathFile.Data());
Warning("PlotEfficiency", "Failed to load %s from %s", nameHistgen.Data(), pathFile.Data());
okNonPrompt = false;
}

Expand Down
13 changes: 12 additions & 1 deletion codeHF/PlotEfficiencyRecoStep.C
View file
Original file line number Diff line number Diff line change
@@ -1,3 +1,14 @@
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.

// Computation and plotting of reconstruction efficiency step-by-step
// Four steps defined: RecoHFFlag, RecoTopol, RecoCand, RecoPID
// RecoHFFlag: candidates properly flagged (e.g. in HFD0CandidateSelector --> hfflag() is D0ToPiK)
Expand All @@ -10,7 +21,7 @@
// .L PlotEfficiencyRecoStep.C
// PlotEfficiencyRecoStep("InputName.root","particlename",true);

#include "utils_plot.h"
#include "../exec/utilitiesPlot.h"

void SetProperAxisRange(TH1F** histo, int NIteration, float marginHigh, float marginLow, bool logScaleH);

Expand Down
48 changes: 15 additions & 33 deletions codeHF/RunHFTaskLocal.C
View file
Original file line number Diff line number Diff line change
@@ -1,4 +1,17 @@
TChain* CreateLocalChain(const char* txtfile);
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.

// Macro to run the HF AliPhysics task that produces validation histograms

#include "../exec/utilitiesAli.h"

Long64_t RunHFTaskLocal(TString txtfile = "./list_ali.txt",
TString jsonfilename = "dpl-config_std.json",
Expand Down Expand Up @@ -26,7 +39,7 @@ Long64_t RunHFTaskLocal(TString txtfile = "./list_ali.txt",

TChain* chainESD = CreateLocalChain(txtfile.Data());
if (!chainESD) {
Error("CreateLocalChain", "Failed to create chain from file %s", txtfile.Data());
Fatal("CreateLocalChain", "Failed to create chain from file %s", txtfile.Data());
return -1;
}

Expand Down Expand Up @@ -60,34 +73,3 @@ Long64_t RunHFTaskLocal(TString txtfile = "./list_ali.txt",
mgr->PrintStatus();
return mgr->StartAnalysis("local", chainESD);
};

TChain* CreateLocalChain(const char* txtfile)
{
// Open the file
ifstream in;
in.open(txtfile);
Int_t count = 0;
// Read the input list of files and add them to the chain
TString line;
TChain* chain = new TChain("esdTree");
while (in.good()) {
in >> line;
if (line.IsNull() || line.BeginsWith("#"))
continue;
TString esdFile(line);
TFile* file = TFile::Open(esdFile);
if (file && !file->IsZombie()) {
chain->Add(esdFile);
file->Close();
} else {
Error("CreateLocalChain", "Skipping un-openable file: %s", esdFile.Data());
}
}
in.close();
if (!chain->GetListOfFiles()->GetEntries()) {
Error("CreateLocalChain", "No file from %s could be opened", txtfile);
delete chain;
return nullptr;
}
return chain;
}
45 changes: 0 additions & 45 deletions codeHF/utils_plot.h
View file

This file was deleted.

Loading

0 comments on commit a01425c

Please sign in to comment.