Skip to content
Snippets Groups Projects
Commit f5336293 authored by Maria Zurek's avatar Maria Zurek
Browse files

Add E/p plots in energy scan analysis

parent 29e01ae9
No related branches found
No related tags found
1 merge request!153Draft: Resolve "Add energy scan for Barrel Ecal"
...@@ -30,22 +30,22 @@ imaging_ecal_energy_scan_e: ...@@ -30,22 +30,22 @@ imaging_ecal_energy_scan_e:
extends: .rec_benchmark extends: .rec_benchmark
stage: run stage: run
timeout: 24 hours timeout: 24 hours
rules: # rules:
- if: '$RUN_EXTENDED_RECO_BENCHMARK == "true"' # - if: '$RUN_EXTENDED_RECO_BENCHMARK == "true"'
script: script:
- export E_file="sim_output/emcal_barrel_energy_scan_points_${PARTICLE}.txt" - export E_file="sim_output/emcal_barrel_energy_scan_points_${PARTICLE}.txt"
- bash benchmarks/imaging_ecal/run_emcal_barrel.sh -t emcal_barrel_${PARTICLE}_${ENERGY} -n 100 -p "${PARTICLE}" -e "${ENERGY}" && echo "${ENERGY}" >> "$E_file" - bash benchmarks/imaging_ecal/run_emcal_barrel.sh -t emcal_barrel_${PARTICLE}_${ENERGY} -n 100 -p "${PARTICLE}" -e "${ENERGY}" && echo "${ENERGY}" >> "$E_file"
parallel: parallel:
matrix: matrix:
- ENERGY: ["0.8", "1", "2", "5", "10"] - ENERGY: ["2", "5", "10"]
PARTICLE: [electron] PARTICLE: [electron, pion-]
imaging_ecal_energy_scan_ph: imaging_ecal_energy_scan_ph:
extends: .rec_benchmark extends: .rec_benchmark
stage: run stage: run
timeout: 24 hours timeout: 24 hours
# rules: rules:
# - if: '$RUN_EXTENDED_RECO_BENCHMARK == "true"' - if: '$RUN_EXTENDED_RECO_BENCHMARK == "true"'
script: script:
- export E_file="sim_output/emcal_barrel_energy_scan_points_${PARTICLE}.txt" - export E_file="sim_output/emcal_barrel_energy_scan_points_${PARTICLE}.txt"
- bash benchmarks/imaging_ecal/run_emcal_barrel.sh -t emcal_barrel_${PARTICLE}_${ENERGY} -n 100 -p "${PARTICLE}" -e "${ENERGY}" && echo "${ENERGY}" >> "$E_file" - bash benchmarks/imaging_ecal/run_emcal_barrel.sh -t emcal_barrel_${PARTICLE}_${ENERGY} -n 100 -p "${PARTICLE}" -e "${ENERGY}" && echo "${ENERGY}" >> "$E_file"
...@@ -58,8 +58,8 @@ imaging_ecal_energy_scan_analysis_electrons: ...@@ -58,8 +58,8 @@ imaging_ecal_energy_scan_analysis_electrons:
extends: .rec_benchmark extends: .rec_benchmark
stage: process stage: process
timeout: 24 hours timeout: 24 hours
rules: # rules:
- if: '$RUN_EXTENDED_RECO_BENCHMARK == "true"' # - if: '$RUN_EXTENDED_RECO_BENCHMARK == "true"'
needs: needs:
- ["imaging_ecal_energy_scan_e"] - ["imaging_ecal_energy_scan_e"]
script: script:
...@@ -70,14 +70,25 @@ imaging_ecal_energy_scan_analysis_photons: ...@@ -70,14 +70,25 @@ imaging_ecal_energy_scan_analysis_photons:
extends: .rec_benchmark extends: .rec_benchmark
stage: process stage: process
timeout: 24 hours timeout: 24 hours
# rules: rules:
# - if: '$RUN_EXTENDED_RECO_BENCHMARK == "true"' - if: '$RUN_EXTENDED_RECO_BENCHMARK == "true"'
needs: needs:
- ["imaging_ecal_energy_scan_ph"] - ["imaging_ecal_energy_scan_ph"]
script: script:
- ls -lhtR sim_output/ - ls -lhtR sim_output/
- root -b -q 'benchmarks/imaging_ecal/scripts/emcal_barrel_energy_scan_analysis.cxx+("photon")' - root -b -q 'benchmarks/imaging_ecal/scripts/emcal_barrel_energy_scan_analysis.cxx+("photon")'
imaging_ecal_energy_scan_analysis_electrons_pions:
extends: .rec_benchmark
stage: process
timeout: 24 hours
# rules:
# - if: '$RUN_EXTENDED_RECO_BENCHMARK == "true"'
needs:
- ["imaging_ecal_energy_scan_e"]
script:
- ls -lhtR sim_output/
- root -b -q 'benchmarks/imaging_ecal/scripts/emcal_barrel_energy_scan_analysis_pi_e.cxx+'
imaging_ecal:results: imaging_ecal:results:
extends: .rec_benchmark extends: .rec_benchmark
......
////////////////////////////////////////
// Read reconstruction ROOT output file
// Plot variables
////////////////////////////////////////
#include "ROOT/RDataFrame.hxx"
#include <iostream>
#include <fmt/core.h>
#include "dd4pod/CalorimeterHitCollection.h"
#include "dd4pod/Geant4ParticleCollection.h"
#include "eicd/CalorimeterHitCollection.h"
#include "eicd/CalorimeterHitData.h"
#include "eicd/ClusterCollection.h"
#include "eicd/ClusterData.h"
#include "eicd/RawCalorimeterHitCollection.h"
#include "eicd/RawCalorimeterHitData.h"
#include "TCanvas.h"
#include "TStyle.h"
#include "TMath.h"
#include "TH1.h"
#include "TF1.h"
#include "TH1D.h"
#include "TRandom.h"
#include "TGraphErrors.h"
#include "DD4hep/Detector.h"
#include "DDG4/Geant4Data.h"
#include "DDRec/CellIDPositionConverter.h"
using ROOT::RDataFrame;
using namespace ROOT::VecOps;
void save_canvas(TCanvas* c, std::string label)
{
c->SaveAs(fmt::format("results/{}.png",label).c_str());
}
void save_canvas(TCanvas* c, std::string label, double E)
{
std::string label_with_E = fmt::format("{}_{}", E, label);
save_canvas(c, label_with_E);
}
void save_canvas(TCanvas* c, std::string label, std::string E_label)
{
std::string label_with_E = fmt::format("{}_{}", E_label, label);
save_canvas(c, label_with_E);
}
void save_canvas(TCanvas* c, std::string var_label, std::string E_label, std::string particle_label)
{
std::string label_with_E = fmt::format("{}_emcal_barrel_{}_{}", E_label, particle_label, var_label);
save_canvas(c, label_with_E);
}
void set_histo_range(TH1D* h)
{
double up_fit = h->GetMean() + 5*h->GetStdDev();
double down_fit = h->GetMean() - 5*h->GetStdDev();
h->GetXaxis()->SetRangeUser(down_fit,up_fit);
}
std::tuple <double, double, double, double, double, double, double, double> extract_sampling_fraction_parameters(std::string E_label)
{
std::string input_fname_electron = fmt::format("sim_output/rec_emcal_barrel_{}_{}.root", electron, E_label);
ROOT::EnableImplicitMT();
ROOT::RDataFrame d0("events", input_fname_electron);
std::string input_fname_pionm = fmt::format("sim_output/rec_emcal_barrel_{}_{}.root", pion-, E_label);
ROOT::EnableImplicitMT();
ROOT::RDataFrame d2("events", input_fname_pionm)
// Thrown Energy [GeV]
auto Ethr = [](std::vector<dd4pod::Geant4ParticleData> const& input) {
auto p = input[2];
auto energy = TMath::Sqrt(p.ps.x * p.ps.x + p.ps.y * p.ps.y + p.ps.z * p.ps.z + p.mass * p.mass);
return energy;
};
// Thrown Momentum [GeV/c]
auto momentum = [](std::vector<dd4pod::Geant4ParticleData> const& input) {
auto p = input[2];
auto momentum = TMath::Sqrt(p.ps.x * p.ps.x + p.ps.y * p.ps.y + p.ps.z * p.ps.z);
return momentum;
}
// Smeared Momentum [GeV/c]
auto momentum_smeared = [](const double momentum, const double sigma_relative) {
TRandom smear();
auto sigma = smear.Gauss(0, sigma_relative*momentum);
return momentum + sigma;
}
// Reconstructed Energy [GeV]
auto Erec = [](const std::vector<eic::CalorimeterHitData>& evt) {
auto total_edep = 0.0;
for (const auto& i: evt)
total_edep += i.energy;
return total_edep;
};
// Number of Clusters
auto ncluster = [](const std::vector<eic::ClusterData>& evt) { return (int)evt.size(); };
// Cluster Energy [GeV]
auto Ecluster = [](const std::vector<eic::ClusterData>& evt) {
auto total_edep = 0.0;
for (const auto& i: evt)
total_edep += i.energy;
return total_edep;
};
// E/p
auto E_over_mom = [](const double E, const double momentum) {
return E / momentum;
};
// Define variables
auto delectron = d0.Define("Ethr", Ethr, {"mcparticles2"})
.Define("mom", momentum, {"mcparticles2"})
.Define("mom_smeared", momentum_smeared, {"mom", 0.005})
.Define("ErecImg", Erec, {"RecoEcalBarrelImagingHits"})
.Define("ErecScFi", Erec, {"EcalBarrelScFiHitsReco"})
.Define("NClusterImg", ncluster, {"EcalBarrelImagingClusters"})
.Define("NClusterScFi", ncluster, {"EcalBarrelScFiClusters"})
.Define("EClusterImg", Ecluster, {"EcalBarrelImagingClusters"})
.Define("EClusterScFi", Ecluster, {"EcalBarrelScFiClusters"})
.Define("E_over_p_ClusterScFi", E_over_p, {"EClusterScFi, mom_smeared"})
.Define("E_over_p_RecoScFi", E_over_p, {"ERecScFi, mom_smeared"});
auto dpion = d2.Define("Ethr", Ethr, {"mcparticles2"})
.Define("mom", momentum, {"mcparticles2"})
.Define("mom_smeared", momentum_smeared, {"mom", 0.005})
.Define("ErecImg", Erec, {"RecoEcalBarrelImagingHits"})
.Define("ErecScFi", Erec, {"EcalBarrelScFiHitsReco"})
.Define("NClusterImg", ncluster, {"EcalBarrelImagingClusters"})
.Define("NClusterScFi", ncluster, {"EcalBarrelScFiClusters"})
.Define("EClusterImg", Ecluster, {"EcalBarrelImagingClusters"})
.Define("EClusterScFi", Ecluster, {"EcalBarrelScFiClusters"})
.Define("E_over_p_ClusterScFi", E_over_p, {"EClusterScFi, mom_smeared"})
.Define("E_over_p_RecoScFi", E_over_p, {"ERecScFi, mom_smeared"});
// Define Histograms
auto hEthr_el = delectron.Histo1D({"hEthr_el", "Thrown Energy; Thrown Energy [GeV]; Events", 100, 0.0, 25.0},"Ethr");
auto hMom_el = delectron.Histo1D({"hMom_el", "Thrown Momentum; Thrown Momentum [GeV/c]; Events", 100, 0.0, 25.0},"mom");
auto hMomSmeared_el = delectron.Histo1D({"hMom_el", "Thrown Momentum Smeared; Thrown Momentum Smeared [GeV/c]; Events", 100, 0.0, 25.0},"mom_smeared");
// auto hEClusterImg_el = delectron.Histo1D({"hEClusterImg", "Cluster Energy; Cluster Energy [GeV]; Events", 200, 0.0, 25.0},"EClusterImg").GetPtr();
// auto hNClusterImg_el = delectron.Histo1D({"hNClusterImg", "Number of Clusters; # of Clusters; Events", 100, 0.0, 100.0}, "NClusterImg").GetPtr();
// auto hEoverPImg_el = delectron.Histo1D({"hfsamImg", "Cluster Energy/E true; Cluster Energy/E true; Events", 100, 0.0, 1.5},"fsamClusterImg").GetPtr();
// auto hEoverPRecImg_el = delectron.Histo1D({"hfsamRecImg", "Reco Hits Energy/E true; Reco Hits Energy/E true; Events", 50, 0.0, 0.1},"fsamRecImg").GetPtr();
auto hEClusterScFi_el = delectron.Histo1D({"hEClusterScFi_el", "Cluster Energy; Cluster Energy [GeV]; Events", 500, 0.0, 25.0},"EClusterScFi").GetPtr();
auto hNClusterScFi_el = delectron.Histo1D({"hNClusterScFi_el", "Number of Clusters; # of Clusters; Events", 100, 0.0, 100.0}, "NClusterScFi").GetPtr();
auto hEoverPScFi_el = delectron.Histo1D({"hfsamScFi_el", "Cluster Energy/E true; Cluster Energy/E true; Events", 100, 0.8, 1.2},"fsamClusterScFi").GetPtr();
auto hEoverPRecScFi_el = delectron.Histo1D({"hfsamRecScFi_el", "Reco Hits Energy/E true; Reco Hits Energy/E true; Events", 50, 0.0, 0.25},"fsamRecScFi").GetPtr();
auto hEthr_pion = dpion.Histo1D({"hEthr_pion", "Thrown Energy; Thrown Energy [GeV]; Events", 100, 0.0, 25.0},"Ethr");
auto hMom_pion = dpion.Histo1D({"hMom_pion", "Thrown Momentum; Thrown Momentum [GeV/c]; Events", 100, 0.0, 25.0},"mom");
auto hMomSmeared_pion = dpion.Histo1D({"hMomSmeared_pion", "Thrown Momentum Smeared; Thrown Momentum Smeared [GeV/c]; Events", 100, 0.0, 25.0},"mom_smeared");
// auto hEClusterImg_pion = dpion.Histo1D({"hEClusterImg", "Cluster Energy; Cluster Energy [GeV]; Events", 200, 0.0, 25.0},"EClusterImg").GetPtr();
// auto hNClusterImg_pion = dpion.Histo1D({"hNClusterImg", "Number of Clusters; # of Clusters; Events", 100, 0.0, 100.0}, "NClusterImg").GetPtr();
// auto hEoverPImg_pion = dpion.Histo1D({"hfsamImg", "Cluster Energy/E true; Cluster Energy/E true; Events", 100, 0.0, 1.5},"fsamClusterImg").GetPtr();
// auto hEoverPRecImg_pion = dpion.Histo1D({"hfsamRecImg", "Reco Hits Energy/E true; Reco Hits Energy/E true; Events", 50, 0.0, 0.1},"fsamRecImg").GetPtr();
auto hEClusterScFi_pion = dpion.Histo1D({"hEClusterScFi_pion", "Cluster Energy; Cluster Energy [GeV]; Events", 500, 0.0, 25.0},"EClusterScFi").GetPtr();
auto hNClusterScFi_pion = dpion.Histo1D({"hNClusterScFi_pion", "Number of Clusters; # of Clusters; Events", 100, 0.0, 100.0}, "NClusterScFi").GetPtr();
auto hEoverPScFi_pion = dpion.Histo1D({"hfsamScFi_pion", "Cluster Energy/E true; Cluster Energy/E true; Events", 100, 0.8, 1.2},"fsamClusterScFi").GetPtr();
auto hEoverPRecScFi_pion = dpion.Histo1D({"hfsamRecScFi_pion", "Reco Hits Energy/E true; Reco Hits Energy/E true; Events", 50, 0.0, 0.25},"fsamRecScFi").GetPtr();
// Event Counts
auto nevents_thrown = delectron.Count();
std::cout << "Number of Thrown Electron Events: " << (*nevents_thrown) << "\n";
// Draw Histograms
TCanvas* c1 = new TCanvas("c1", "c1", 700, 500);
c1->Divide(2,2);
c1->SetLogy(1);
c1->cd(1);
hEthr_el->GetYaxis()->SetTitleOffset(1.4);
hEthr_el->SetLineWidth(2);
hEthr_el->SetLineColor(kBlue);
hEthr_el->DrawClone();
hEthr_pion->SetLineColor(kRed);
hEthr_pion->DrawClone("same");
c1->cd(2);
hMom_el->GetYaxis()->SetTitleOffset(1.4);
hMom_el->SetLineWidth(2);
hMom_el->SetLineColor(kBlue);
hMom_el->DrawClone();
hMom_pion->SetLineColor(kRed);
hMom_pion->DrawClone("same");
c1->cd(3);
hMomSmeared_el->GetYaxis()->SetTitleOffset(1.4);
hMomSmeared_el->SetLineWidth(2);
hMomSmeared_el->SetLineColor(kBlue);
hMomSmeared_el->DrawClone();
hMomSmeared_pion->SetLineColor(kRed);
hMomSmeared_pion->DrawClone("same");
save_canvas(c1, "gen", E_label, "electron_pion");
TCanvas* c2 = new TCanvas("c2", "c2", 1400, 1000);
c2->cd();
hNClusterScFi_el->GetYaxis()->SetTitleOffset(1.4);
hNClusterScFi_el->SetLineWidth(2);
hNClusterScFi_el->SetLineColor(kBlue);
hNClusterScFi_el->DrawClone();
hNClusterScFi_pion->SetLineColor(kRed);
hNClusterScFi_pion->DrawClone("same");
save_canvas(c2, "NCluster", E_label, "electron_pion");
TCanvas* c3 = new TCanvas("c3", "c3", 1400, 1000);
c3->cd();
hEoverPScFi_el->GetYaxis()->SetTitleOffset(1.4);
hEoverPScFi_el->SetLineWidth(2);
hEoverPScFi_el->SetLineColor(kBlue);
hEoverPScFi_el->DrawClone();
hEoverPScFi_pion->SetLineColor(kRed);
hEoverPScFi_pion->DrawClone("same");
save_canvas(c3, "EoverPClusterScFi", E_label, "electron_pion");
TCanvas* c4 = new TCanvas("c4", "c4", 1400, 1000);
c4->cd();
hEoverPRecScFi_el->GetYaxis()->SetTitleOffset(1.4);
hEoverPRecScFi_el->SetLineWidth(2);
hEoverPRecScFi_el->SetLineColor(kBlue);
hEoverPRecScFii_el->DrawClone();
hEoverPRecScFi_pion->SetLineColor(kRed);
hEoverPRecScFi_pion->DrawClone("same");
save_canvas(c4, "EoverPRecScFi", E_label, "electron_pion");
return std::make_tuple(meanScFi, sigmaScFi, meanScFi_err, sigmaScFi_err, meanImg, sigmaImg, meanImg_err, sigmaImg_err);
}
std::vector<std::string> read_scanned_energies(std::string input_energies_fname)
{
std::vector<std::string> scanned_energies;
std::string E_label;
ifstream E_file (input_energies_fname);
if (E_file.is_open())
{
while (E_file >> E_label)
{
scanned_energies.push_back(E_label);
}
E_file.close();
return scanned_energies;
}
else
{
std::cout << "Unable to open file " << input_energies_fname << std::endl;
abort();
}
}
void fill_graph_error(TGraphErrors gr_fsam, TGraphErrors gr_fsam_res, double fsam, double fsam_res, double fsam_err, double fsam_res_err, double E) {
gr_fsam.SetPoint(gr_fsam.GetN(),E,100*fsam);
gr_fsam.SetPointError(gr_fsam.GetN()-1,0., 100*fsam_err);
gr_fsam_res.SetPoint(gr_fsam_res.GetN(),E,100.0*(fsam_res/fsam));
auto fsam_res_rel_err = 100.0*(sqrt(pow((fsam_res_err/fsam),2)+pow((fsam_err*fsam_res)/(fsam*fsam),2)));
gr_fsam_res.SetPointError(gr_fsam_res.GetN()-1,0.,fsam_res_rel_err);
}
void emcal_barrel_energy_scan_analysis_pi_e(std::string particle_label = "electron")
{
// Setting for graphs
gROOT->SetStyle("Plain");
gStyle->SetOptFit(1);
gStyle->SetLineWidth(2);
gStyle->SetPadTickX(1);
gStyle->SetPadTickY(1);
gStyle->SetPadGridX(1);
gStyle->SetPadGridY(1);
gStyle->SetPadLeftMargin(0.14);
gStyle->SetPadRightMargin(0.14);
auto scanned_energies = read_scanned_energies(fmt::format("sim_output/emcal_barrel_energy_scan_points_{}.txt", particle_label));
// TGraphErrors gr_fsamImg(scanned_energies.size()-1);
// TGraphErrors gr_fsamImg_res(scanned_energies.size()-1);
// TGraphErrors gr_fsamScFi(scanned_energies.size()-1);
// TGraphErrors gr_fsamScFi_res(scanned_energies.size()-1);
// for (const auto& E_label : scanned_energies) {
// auto [fsamScFi, fsamScFi_res, fsamScFi_err, fsamScFi_res_err, fsamImg, fsamImg_res, fsamImg_err, fsamImg_res_err] = extract_sampling_fraction_parameters(particle_label, E_label);
// auto E = std::stod(E_label);
// fill_graph_error(gr_fsamScFi, gr_fsamScFi_res, fsamScFi, fsamScFi_res, fsamScFi_err, fsamScFi_res_err, E);
// fill_graph_error(gr_fsamImg, gr_fsamImg_res, fsamImg, fsamImg_res, fsamImg_err, fsamImg_res_err, E);
// }
// TCanvas* c10 = new TCanvas("c10", "c10", 1400, 500);
// c10->cd(1);
// gr_fsamScFi.SetTitle("ScFi Cluster Energy/E true Scan;True Energy [GeV];Cluster Energy/E true [%]");
// gr_fsamScFi.SetMarkerStyle(20);
// gr_fsamScFi.Fit("pol0", "", "", 2., 20.);
// gr_fsamScFi.Draw("APE");
// c10->cd(2);
// TF1* func_res = new TF1("func_res", "[0]/sqrt(x) + [1]", 0.25, 20.);
// func_res->SetLineWidth(2);
// func_res->SetLineColor(kRed);
// gr_fsamScFi_res.SetTitle("ScFi Energy Resolution;True Energy [GeV];#Delta E/E [%]");
// gr_fsamScFi_res.SetMarkerStyle(20);
// gr_fsamScFi_res.Fit(func_res,"R");
// gr_fsamScFi_res.Draw("APE");
// save_canvas(c10, fmt::format("emcal_barrel_{}_fsamScFi_scan", particle_label));
// TCanvas* c11 = new TCanvas("c11", "c11", 1400, 500);
// c11->cd(1);
// gr_fsamImg.SetTitle("Img Cluster Energy/E true Scan;True Energy [GeV];Cluster Energy/E true [%]");
// gr_fsamImg.SetMarkerStyle(20);
// gr_fsamImg.Fit("pol0", "", "", 2., 20.);
// gr_fsamImg.Draw("APE");
// c11->cd(2);
// TF1* funcImg_res = new TF1("funcImg_res", "[0]/sqrt(x) + [1]", 0.25, 20.);
// funcImg_res->SetLineWidth(2);
// funcImg_res->SetLineColor(kRed);
// gr_fsamImg_res.SetTitle("Energy Resolution;True Energy [GeV];#Delta E/E [%]");
// gr_fsamImg_res.SetMarkerStyle(20);
// gr_fsamImg_res.Fit(funcImg_res,"R");
// gr_fsamImg_res.Draw("APE");
// save_canvas(c11, fmt::format("emcal_barrel_{}_fsamImg_scan", particle_label));
}
...@@ -155,14 +155,14 @@ std::tuple <double, double, double, double, double, double, double, double> extr ...@@ -155,14 +155,14 @@ std::tuple <double, double, double, double, double, double, double, double> extr
hEdigiImg->GetYaxis()->SetTitleOffset(1.4); hEdigiImg->GetYaxis()->SetTitleOffset(1.4);
hEdigiImg->SetLineWidth(2); hEdigiImg->SetLineWidth(2);
hEdigiImg->SetLineColor(kBlue); hEdigiImg->SetLineColor(kBlue);
set_histo_range(hEdigiImg);
hEdigiImg->DrawClone(); hEdigiImg->DrawClone();
//set_histo_range(hEdigiImg);
c2->cd(2); c2->cd(2);
hEdigiScFi->GetYaxis()->SetTitleOffset(1.4); hEdigiScFi->GetYaxis()->SetTitleOffset(1.4);
hEdigiScFi->SetLineWidth(2); hEdigiScFi->SetLineWidth(2);
hEdigiScFi->SetLineColor(kBlue); hEdigiScFi->SetLineColor(kBlue);
set_histo_range(hEdigiScFi); //set_histo_range(hEdigiScFi);
hEdigiScFi->DrawClone(); hEdigiScFi->DrawClone();
c2->cd(3); c2->cd(3);
...@@ -175,7 +175,7 @@ std::tuple <double, double, double, double, double, double, double, double> extr ...@@ -175,7 +175,7 @@ std::tuple <double, double, double, double, double, double, double, double> extr
hErecScFi->GetYaxis()->SetTitleOffset(1.4); hErecScFi->GetYaxis()->SetTitleOffset(1.4);
hErecScFi->SetLineWidth(2); hErecScFi->SetLineWidth(2);
hErecScFi->SetLineColor(kBlue); hErecScFi->SetLineColor(kBlue);
set_histo_range(hErecScFi); //set_histo_range(hErecScFi);
hErecScFi->DrawClone(); hErecScFi->DrawClone();
save_canvas(c2, "E_digi_rec", E_label, particle_label); save_canvas(c2, "E_digi_rec", E_label, particle_label);
{ {
...@@ -185,14 +185,20 @@ std::tuple <double, double, double, double, double, double, double, double> extr ...@@ -185,14 +185,20 @@ std::tuple <double, double, double, double, double, double, double, double> extr
hfsamRecImg->GetYaxis()->SetTitleOffset(1.4); hfsamRecImg->GetYaxis()->SetTitleOffset(1.4);
hfsamRecImg->SetLineWidth(2); hfsamRecImg->SetLineWidth(2);
hfsamRecImg->SetLineColor(kBlue); hfsamRecImg->SetLineColor(kBlue);
set_histo_range(hfsamRecImg); auto up_fit = hfsamRecImgi->GetMean() + 5*hfsamRecImg->GetStdDev();
auto down_fit = hfsamRecImg->GetMean() - 5*hfsamRecImg->GetStdDev();
hfsamRecImg->GetXaxis()->SetRangeUser(down_fit,up_fit);
//set_histo_range(hfsamRecImg);
hfsamRecImg->DrawClone(); hfsamRecImg->DrawClone();
c3->cd(2); c3->cd(2);
hfsamRecScFi->GetYaxis()->SetTitleOffset(1.4); hfsamRecScFi->GetYaxis()->SetTitleOffset(1.4);
hfsamRecScFi->SetLineWidth(2); hfsamRecScFi->SetLineWidth(2);
hfsamRecScFi->SetLineColor(kBlue); hfsamRecScFi->SetLineColor(kBlue);
set_histo_range(hfsamRecScFi); auto up_fit = hfsamRecScFi->GetMean() + 5*hfsamRecScFi->GetStdDev();
auto down_fit = hfsamRecScFi->GetMean() - 5*hfsamRecScFi->GetStdDev();
hfsamRecScFi->GetXaxis()->SetRangeUser(down_fit,up_fit);
//set_histo_range(hfsamRecScFi);
hfsamRecScFi->DrawClone(); hfsamRecScFi->DrawClone();
save_canvas(c3, "fsam_digi_rec", E_label, particle_label); save_canvas(c3, "fsam_digi_rec", E_label, particle_label);
...@@ -202,7 +208,10 @@ std::tuple <double, double, double, double, double, double, double, double> extr ...@@ -202,7 +208,10 @@ std::tuple <double, double, double, double, double, double, double, double> extr
hEClusterScFi->GetYaxis()->SetTitleOffset(1.4); hEClusterScFi->GetYaxis()->SetTitleOffset(1.4);
hEClusterScFi->SetLineWidth(2); hEClusterScFi->SetLineWidth(2);
hEClusterScFi->SetLineColor(kBlue); hEClusterScFi->SetLineColor(kBlue);
set_histo_range(hEClusterScFi); //set_histo_range(hEClusterScFi);
auto up_fit = hEClusterScFi->GetMean() + 5*hEClusterScFi->GetStdDev();
auto down_fit = hEClusterScFi->GetMean() - 5*hEClusterScFi->GetStdDev();
hEClusterScFi->GetXaxis()->SetRangeUser(down_fit,up_fit)
hEClusterScFi->DrawClone(); hEClusterScFi->DrawClone();
save_canvas(c5, "EClusterSCFi", E_label, particle_label); save_canvas(c5, "EClusterSCFi", E_label, particle_label);
...@@ -213,7 +222,10 @@ std::tuple <double, double, double, double, double, double, double, double> extr ...@@ -213,7 +222,10 @@ std::tuple <double, double, double, double, double, double, double, double> extr
hEClusterImg->GetYaxis()->SetTitleOffset(1.4); hEClusterImg->GetYaxis()->SetTitleOffset(1.4);
hEClusterImg->SetLineWidth(2); hEClusterImg->SetLineWidth(2);
hEClusterImg->SetLineColor(kBlue); hEClusterImg->SetLineColor(kBlue);
set_histo_range(hEClusterImg); auto up_fit = hEClusterImg->GetMean() + 5*hEClusterImg->GetStdDev();
auto down_fit = hEClusterImg->GetMean() - 5*hEClusterImg->GetStdDev();
hEClusterImg->GetXaxis()->SetRangeUser(down_fit,up_fit);
//set_histo_range(hEClusterImg);
hEClusterImg->DrawClone(); hEClusterImg->DrawClone();
save_canvas(c6, "EClusterImg", E_label, particle_label); save_canvas(c6, "EClusterImg", E_label, particle_label);
...@@ -223,12 +235,13 @@ std::tuple <double, double, double, double, double, double, double, double> extr ...@@ -223,12 +235,13 @@ std::tuple <double, double, double, double, double, double, double, double> extr
hfsamImg->GetYaxis()->SetTitleOffset(1.4); hfsamImg->GetYaxis()->SetTitleOffset(1.4);
hfsamImg->SetLineWidth(2); hfsamImg->SetLineWidth(2);
hfsamImg->SetLineColor(kBlue); hfsamImg->SetLineColor(kBlue);
set_histo_range(hfsamImg); //set_histo_range(hfsamImg);
hfsamImg->DrawClone(); hfsamImg->DrawClone();
auto up_fit = hfsamImg->GetMean() + 5*hfsamImg->GetStdDev(); auto up_fit = hfsamImg->GetMean() + 5*hfsamImg->GetStdDev();
auto down_fit = hfsamImg->GetMean() - 5*hfsamImg->GetStdDev(); auto down_fit = hfsamImg->GetMean() - 5*hfsamImg->GetStdDev();
if(down_fit <=0 ) down_fit = hfsamImg->GetXaxis()->GetBinUpEdge(1); if(down_fit <=0 ) down_fit = hfsamImg->GetXaxis()->GetBinUpEdge(1);
hfsamImg->Fit("gaus", "", "", down_fit, up_fit); hfsamImg->Fit("gaus", "", "", down_fit, up_fit);
hfsamImg->GetXaxis()->SetRangeUser(down_fit,up_fit);
TF1 *gausImg = hfsamImg->GetFunction("gaus"); TF1 *gausImg = hfsamImg->GetFunction("gaus");
gausImg->SetLineWidth(2); gausImg->SetLineWidth(2);
gausImg->SetLineColor(kRed); gausImg->SetLineColor(kRed);
...@@ -243,13 +256,14 @@ std::tuple <double, double, double, double, double, double, double, double> extr ...@@ -243,13 +256,14 @@ std::tuple <double, double, double, double, double, double, double, double> extr
hfsamScFi->GetYaxis()->SetTitleOffset(1.4); hfsamScFi->GetYaxis()->SetTitleOffset(1.4);
hfsamScFi->SetLineWidth(2); hfsamScFi->SetLineWidth(2);
hfsamScFi->SetLineColor(kBlue); hfsamScFi->SetLineColor(kBlue);
set_histo_range(hfsamScFi); //set_histo_range(hfsamScFi);
hfsamScFi->DrawClone(); hfsamScFi->DrawClone();
hfsamScFi->Fit("gaus", "", "", 0.01, 0.1); hfsamScFi->Fit("gaus", "", "", 0.01, 0.1);
up_fit = hfsamScFi->GetMean() + 5*hfsamScFi->GetStdDev(); up_fit = hfsamScFi->GetMean() + 5*hfsamScFi->GetStdDev();
down_fit = hfsamScFi->GetMean() - 5*hfsamScFi->GetStdDev(); down_fit = hfsamScFi->GetMean() - 5*hfsamScFi->GetStdDev();
if(down_fit <=0 ) down_fit = hfsamScFi->GetXaxis()->GetBinUpEdge(1); if(down_fit <=0 ) down_fit = hfsamScFi->GetXaxis()->GetBinUpEdge(1);
hfsamScFi->Fit("gaus", "", "", down_fit, up_fit); hfsamScFi->Fit("gaus", "", "", down_fit, up_fit);
hfsamScFi->GetXaxis()->SetRangeUser(down_fit,up_fit);
TF1 *gaus = hfsamScFi->GetFunction("gaus"); TF1 *gaus = hfsamScFi->GetFunction("gaus");
gaus->SetLineWidth(2); gaus->SetLineWidth(2);
gaus->SetLineColor(kRed); gaus->SetLineColor(kRed);
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment