Skip to content
Snippets Groups Projects
Commit 1762d39d authored by Roger Naranjo's avatar Roger Naranjo Committed by John Kenneth Anders
Browse files

Adding IQE plots to Egamma ART tests

parent cf2eefae
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Reconstruction/egamma/egammaValidation/scripts/EgammaARTmonitoring_plotsMaker.py
+ 34
− 11
View file @ 1762d39d
......@@ -145,6 +145,20 @@ def make_comparison_plots(type, f_base, f_nightly, result_file):
make_ratio_plot(h_base, h_nightly, folder['title'], result_file)
def makeIQEPlots(inHist, name):
outHist = inHist.QuantilesX(0.75, "EResolution_IQE_mu")
outHist.GetXaxis().SetTitle("<#mu>")
outHist.GetYaxis().SetTitle("IQE")
outHist25 = inHist.QuantilesX(0.25, "EResolutio_IQE_mu_25")
outHist.Add(outHist25, -1)
outHist.Scale(1/1.349)
return outHist.Clone(inHist.GetName() + "_"+ name)
def make_profile_plots(f_base, f_nightly, result_file, particle_type):
cluster_list_to_loop = cluster_list
......@@ -155,17 +169,26 @@ def make_profile_plots(f_base, f_nightly, result_file, particle_type):
for histo in get_key_names(f_nightly, folder['name']):
if '2D' not in histo:
continue
h_base = f_base.Get(folder['name'] + '/' + histo)
h_base_profile = h_base.ProfileX(histo+"_ProfileB")
h_nightly = f_nightly.Get(folder['name'] + '/' + histo)
h_nightly_profile = h_nightly.ProfileX(histo+"_Profile")
h_base_profile.SetDirectory(0)
h_nightly_profile.SetDirectory(0)
y_axis_label = "Mean %s" % (h_base_profile.GetTitle() )
h_base_profile.SetTitle("")
make_ratio_plot(h_base_profile, h_nightly_profile, folder['title'], result_file, y_axis_label)
if 'mu' in histo:
h_base = f_base.Get(folder['name'] + '/' + histo)
h_nightly = f_nightly.Get(folder['name'] + '/' + histo)
h_base = makeIQEPlots(h_base,'IQE')
h_nightly = makeIQEPlots(h_nightly,'IQE')
make_ratio_plot(h_base, h_nightly, folder['title'], result_file, 'IQE')
else:
h_base = f_base.Get(folder['name'] + '/' + histo)
h_base_profile = h_base.ProfileX(histo+"_ProfileB")
h_nightly = f_nightly.Get(folder['name'] + '/' + histo)
h_nightly_profile = h_nightly.ProfileX(histo+"_Profile")
h_base_profile.SetDirectory(0)
h_nightly_profile.SetDirectory(0)
y_axis_label = "Mean %s" % (h_base_profile.GetTitle() )
h_base_profile.SetTitle("")
make_ratio_plot(h_base_profile, h_nightly_profile, folder['title'], result_file, y_axis_label)
......
Reconstruction/egamma/egammaValidation/src/ClusterHistograms.cxx
+ 10
− 7
View file @ 1762d39d
......@@ -18,7 +18,7 @@ StatusCode ClusterHistograms::initializePlots() {
histo2DMap["number_cells_vs_et_in_layer_3_2D"] = (new TH2D(Form("%s_%s",m_name.c_str(),"number_cells_vs_et_in_layer_3_2D"), "Number of cells;E_{T}", 600,0,300,50, 0,100.));
histo2DMap["number_cell_in_layer"] = (new TH2D(Form("%s_%s",m_name.c_str(),"number_cell_in_layer"), ";Number of cells;Layer",100,0,200, 4,0,4));
histo2DMap["mu_energy_resolution_2D"] = (new TH2D(Form("%s_%s",m_name.c_str(),"mu_energy_resolution_2D"), ";<#mu>; Energy Resolution", 5,0,80,20,-1,1));
ATH_CHECK(m_rootHistSvc->regHist(m_folder+"Eraw_Etruth_vs_Etruth_2D", histo2DMap["Eraw_Etruth_vs_Etruth_2D"]));
ATH_CHECK(m_rootHistSvc->regHist(m_folder+"Eraw_Etruth_vs_eta_2D", histo2DMap["Eraw_Etruth_vs_eta_2D"]));
......@@ -28,9 +28,9 @@ StatusCode ClusterHistograms::initializePlots() {
ATH_CHECK(m_rootHistSvc->regHist(m_folder+"number_cells_vs_et_in_layer_2_2D", histo2DMap["number_cells_vs_et_in_layer_2_2D"]));
ATH_CHECK(m_rootHistSvc->regHist(m_folder+"number_cells_vs_et_in_layer_3_2D", histo2DMap["number_cells_vs_et_in_layer_3_2D"]));
ATH_CHECK(m_rootHistSvc->regHist(m_folder+"number_cell_in_layer", histo2DMap["number_cell_in_layer"]));
ATH_CHECK(m_rootHistSvc->regHist(m_folder+"mu_energy_resolution_2D", histo2DMap["mu_energy_resolution_2D"]));
......@@ -44,17 +44,20 @@ void ClusterHistograms::fill(const xAOD::Egamma& egamma) {
void ClusterHistograms::fill(const xAOD::Egamma& egamma, float mu = 0) {
(void)mu;
const xAOD::CaloCluster *cluster = egamma.caloCluster();
const xAOD::TruthParticle *truth_egamma = xAOD::TruthHelpers::getTruthParticle(egamma);
if ( !truth_egamma ) return;
histo2DMap["Eraw_Etruth_vs_Etruth_2D"]->Fill(truth_egamma->e()/1000,cluster->rawE()/truth_egamma->e());
histo2DMap["Eraw_Etruth_vs_eta_2D"]->Fill(truth_egamma->eta(),cluster->rawE()/truth_egamma->e());
const auto Ereco = cluster->rawE();
const auto Etruth = truth_egamma->e();
const auto Eres = (Ereco - Etruth)/Etruth;
histo2DMap["Eraw_Etruth_vs_Etruth_2D"]->Fill(truth_egamma->e()/1000,Ereco/Etruth);
histo2DMap["Eraw_Etruth_vs_eta_2D"]->Fill(truth_egamma->eta(),Ereco/Etruth);
histo2DMap["mu_energy_resolution_2D"]->Fill(mu,Eres);
const CaloClusterCellLink* cellLinks = cluster->getCellLinks();
......
Reconstruction/egamma/egammaValidation/src/EgammaMonitoring.cxx
+ 18
− 15
View file @ 1762d39d
......@@ -276,6 +276,9 @@ StatusCode EgammaMonitoring::execute() {
// Retrieve things from the event store
const xAOD::EventInfo *eventInfo = nullptr;
ANA_CHECK(evtStore()->retrieve(eventInfo, "EventInfo"));
const float mu = eventInfo->averageInteractionsPerCrossing();
// Retrieve egamma truth particles
const xAOD::TruthParticleContainer *egTruthParticles = nullptr;
......@@ -325,9 +328,9 @@ StatusCode EgammaMonitoring::execute() {
if (!electron) continue;
clusterPromptAll->fill(*electron);
clusterPromptAll->fill(*electron,mu);
if (egtruth->pt() > 10*Gaudi::Units::GeV) {
clusterPrompt10GeV->fill(*electron);
clusterPrompt10GeV->fill(*electron,mu);
}
......@@ -502,9 +505,9 @@ StatusCode EgammaMonitoring::execute() {
for (auto elrec : *RecoEl) {
if (!elrec) continue;
clusterAll->fill(*elrec);
clusterAll->fill(*elrec,mu);
if (elrec->pt() > 10*Gaudi::Units::GeV) {
cluster10GeV->fill(*elrec);
cluster10GeV->fill(*elrec,mu);
}
recoElectronAll->fill(*elrec);
showerShapesAll->fill(*elrec);
......@@ -579,9 +582,9 @@ StatusCode EgammaMonitoring::execute() {
recoPhotonAll->fill(*phrec);
isolationAll->fill(*phrec);
showerShapesAll->fill(*phrec);
clusterAll->fill(*phrec);
clusterAll->fill(*phrec,mu);
if (phrec->pt() > 10*Gaudi::Units::GeV) {
cluster10GeV->fill(*phrec);
cluster10GeV->fill(*phrec,mu);
}
if (phrec->pt() > 10*Gaudi::Units::GeV){
showerShapes10GeV->fill(*phrec);
......@@ -600,9 +603,9 @@ StatusCode EgammaMonitoring::execute() {
if (!photon) continue;
truthPhotonRecoPhoton->fill(*egtruth);
clusterPromptAll->fill(*photon);
clusterPromptAll->fill(*photon,mu);
if (egtruth->pt() > 10*Gaudi::Units::GeV) {
clusterPrompt10GeV->fill(*photon);
clusterPrompt10GeV->fill(*photon,mu);
}
bool isRecoConv = xAOD::EgammaHelpers::isConvertedPhoton(photon);
......@@ -616,23 +619,23 @@ StatusCode EgammaMonitoring::execute() {
truthPhotonConvRecoConv->fill(*egtruth);
clusterConvPhoton->fill(*photon);
clusterConvPhoton->fill(*photon,mu);
if (convType == xAOD::EgammaParameters::singleSi) {
truthPhotonConvRecoConv1Si->fill(*egtruth);
clusterConvPhotonSi->fill(*photon);
clusterConvPhotonSi->fill(*photon,mu);
} else if (convType == xAOD::EgammaParameters::singleTRT) {
truthPhotonConvRecoConv1TRT->fill(*egtruth);
clusterConvPhotonTRT->fill(*photon);
clusterConvPhotonTRT->fill(*photon,mu);
} else if (convType == xAOD::EgammaParameters::doubleSi) {
truthPhotonConvRecoConv2Si->fill(*egtruth);
clusterConvPhotonSiSi->fill(*photon);
clusterConvPhotonSiSi->fill(*photon,mu);
} else if (convType == xAOD::EgammaParameters::doubleTRT) {
truthPhotonConvRecoConv2TRT->fill(*egtruth);
clusterConvPhotonTRTTRT->fill(*photon);
clusterConvPhotonTRTTRT->fill(*photon,mu);
} else if (convType == xAOD::EgammaParameters::doubleSiTRT) {
truthPhotonConvRecoConv2SiTRT->fill(*egtruth);
clusterConvPhotonSiTRT->fill(*photon);
clusterConvPhotonSiTRT->fill(*photon,mu);
}
if (m_IsoFixedCutTight->accept(*photon)) recoPhotonConvIsoFixedCutTight->fill(*egtruth);
......@@ -641,7 +644,7 @@ StatusCode EgammaMonitoring::execute() {
} // isRecoConv
else {
truthPhotonConvRecoUnconv->fill(*egtruth);
clusterUnconvPhoton->fill(*photon);
clusterUnconvPhoton->fill(*photon,mu);
}
} //isTrueConv
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment