Skip to content
Snippets Groups Projects
Commit 3cde54b7 authored by Jihee Kim's avatar Jihee Kim Committed by Whitney Armstrong
Browse files

Included pi0 script in the benchmarks CI

Moved datasets scripts over benchmarks ecal/scripts directory
parent 27388635
Branches
No related tags found
1 merge request!4pi0 data set
...@@ -12,3 +12,15 @@ ecal_1_emcal_electrons: ...@@ -12,3 +12,15 @@ ecal_1_emcal_electrons:
script: script:
- bash ecal/emcal_electrons.sh - bash ecal/emcal_electrons.sh
ecal_1_emcal_pi0s:
image: eicweb.phy.anl.gov:4567/eic/juggler/juggler:latest
tags:
- silicon
timeout: 12 hours 30 minutes
artifacts:
expire_in: 20 weeks
paths:
- results/
stage: run
script:
- bash ecal/emcal_pi0s.sh
#!/bin/bash
# Based on emcal_electrons.sh script
# Detector
if [[ ! -n "${JUGGLER_DETECTOR}" ]] ; then
export JUGGLER_DETECTOR="topside"
fi
# Number of events
if [[ ! -n "${JUGGLER_N_EVENTS}" ]] ; then
export JUGGLER_N_EVENTS=100
fi
# File names
export JUGGLER_FILE_NAME_TAG="emcal_uniform_pi0s"
export JUGGLER_GEN_FILE="${JUGGLER_FILE_NAME_TAG}.hepmc"
export JUGGLER_SIM_FILE="sim_${JUGGLER_FILE_NAME_TAG}.root"
export JUGGLER_REC_FILE="rec_${JUGGLER_FILE_NAME_TAG}.root"
echo "JUGGLER_N_EVENTS = ${JUGGLER_N_EVENTS}"
echo "JUGGLER_DETECTOR = ${JUGGLER_DETECTOR}"
# Datasets
#git clone https://eicweb.phy.anl.gov/EIC/datasets.git datasets
root -b -q "ecal/scripts/emcal_pi0.cxx(${JUGGLER_N_EVENTS}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
root -b -q "ecal/scripts/emcal_pi0_reader.cxx(\"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
# Detector TOPSiDE
git clone https://eicweb.phy.anl.gov/EIC/detectors/${JUGGLER_DETECTOR}.git
mkdir ${JUGGLER_DETECTOR}/build
pushd ${JUGGLER_DETECTOR}/build
cmake ../. -DCMAKE_INSTALL_PREFIX=/usr/local && make -j30 install
popd
echo "CHECK POINT FOR GEANT4 SIMULATION"
pwd
pushd ${JUGGLER_DETECTOR}
pwd
ls -l
echo "JUGGLER_N_EVENTS = ${JUGGLER_N_EVENTS}"
echo "JUGGLER_DETECTOR = ${JUGGLER_DETECTOR}"
echo "JUGGLER_FILE_NAME_TAG = ${JUGGLER_FILE_NAME_TAG}"
echo "JUGGLER_SIM_FILE = ${JUGGLER_SIM_FILE}"
# run geant4 simulations
npsim --runType batch \
-v WARNING \
--numberOfEvents ${JUGGLER_N_EVENTS} \
--compactFile ${JUGGLER_DETECTOR}.xml \
--inputFiles ../${JUGGLER_FILE_NAME_TAG}.hepmc \
--outputFile ${JUGGLER_SIM_FILE}
# Need to figure out how to pass file name to juggler from the commandline
xenv -x /usr/local/Juggler.xenv gaudirun.py ../ecal/options/example_crystal.py
ls -l
popd
ls -l
pwd
mkdir -p results
root -b -q "ecal/scripts/makeplot_pi0.C(\"${JUGGLER_DETECTOR}/${JUGGLER_REC_FILE}\")"
mkdir -p sim_output
cp "${JUGGLER_DETECTOR}/${JUGGLER_REC_FILE}" sim_output/.
cp "${JUGGLER_DETECTOR}/${JUGGLER_SIM_FILE}" sim_output/.
//////////////////////////////////////////////////////////////
// Crystal EMCAL detector
// Single Pion dataset
// J.KIM 10/18/2020
//////////////////////////////////////////////////////////////
#include "HepMC3/GenEvent.h"
#include "HepMC3/ReaderAscii.h"
#include "HepMC3/WriterAscii.h"
#include "HepMC3/Print.h"
#include <iostream>
#include<random>
#include<cmath>
#include <math.h>
#include <TMath.h>
using namespace HepMC3;
void emcal_pi0(int n_events = 1e6, const char* out_fname = "./data/emcal_pi0_0GeVto30GeV_100kEvt.hepmc")
{
WriterAscii hepmc_output(out_fname);
int events_parsed = 0;
GenEvent evt(Units::GEV, Units::MM);
// Random number generator
TRandom *r1 = new TRandom();
for (events_parsed = 0; events_parsed < n_events; events_parsed++) {
// FourVector(px,py,pz,e,pdgid,status)
// type 4 is beam
// pdgid 11 - electron
// pdgid 111 - pi0
// pdgid 2212 - proton
GenParticlePtr p1 =
std::make_shared<GenParticle>(FourVector(0.0, 0.0, 10.0, 10.0), 11, 4);
GenParticlePtr p2 = std::make_shared<GenParticle>(
FourVector(0.0, 0.0, 0.0, 0.938), 2212, 4);
// Define momentum
Double_t p = r1->Uniform(0.0,30.0);
Double_t px;
Double_t py;
Double_t pz;
// Generates random vectors, uniformly distributed over the surface of a
// sphere of given radius, in this case momentum.
r1->Sphere(px, py, pz, p);
// type 1 is final state
// pdgid 111 - pi0 135 MeV/c^2
GenParticlePtr p3 = std::make_shared<GenParticle>(
FourVector(
px, py, pz,
sqrt(p*p + (0.134976*0.134976))),
111, 1);
GenVertexPtr v1 = std::make_shared<GenVertex>();
v1->add_particle_in(p1);
v1->add_particle_in(p2);
v1->add_particle_out(p3);
evt.add_vertex(v1);
if (events_parsed == 0) {
std::cout << "First event: " << std::endl;
Print::listing(evt);
}
hepmc_output.write_event(evt);
if (events_parsed % 10000 == 0) {
std::cout << "Event: " << events_parsed << std::endl;
}
evt.clear();
}
hepmc_output.close();
std::cout << "Events parsed and written: " << events_parsed << std::endl;
}
//////////////////////////
// Crystal EMCAL detector
// Singe Pion dataset
// J.KIM 10/18/2020
//////////////////////////
#include "HepMC3/GenEvent.h"
#include "HepMC3/ReaderAscii.h"
#include "HepMC3/WriterAscii.h"
#include "HepMC3/Print.h"
#include "TH1F.h"
#include <iostream>
#include "TStyle.h"
using namespace HepMC3;
void emcal_pi0_reader(const char* in_fname = "./data/emcal_electron_0GeVto30GeV_100kEvt.hepmc"){
// 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);
ReaderAscii hepmc_input(in_fname);
int events_parsed = 0;
GenEvent evt(Units::GEV, Units::MM);
// Histograms
TH1F* h_pi0_energy = new TH1F("h_pi0_energy", "pi0 energy;E [GeV];Events", 100,-0.5,30.5);
TH1F* h_pi0_eta = new TH1F("h_pi0_eta", "pi0 #eta;#eta;Events", 100,-10.0,10.0);
TH1F* h_pi0_theta = new TH1F("h_pi0_theta", "pi0 #theta;#theta [degree];Events", 100,-0.5,180.5);
TH1F* h_pi0_phi = new TH1F("h_pi0_phi", "pi0 #phi;#phi [degree];Events", 100,-180.5,180.5);
TH2F* h_pi0_pzpt = new TH2F("h_pi0_pzpt", "pi0 pt vs pz;pt [GeV];pz [GeV]", 100,-0.5,30.5,100,-30.5,30.5);
TH2F* h_pi0_pxpy = new TH2F("h_pi0_pxpy", "pi0 px vs py;px [GeV];py [GeV]", 100,-30.5,30.5,100,-30.5,30.5);
TH3F* h_pi0_p = new TH3F("h_pi0_p", "pi0 p;px [GeV];py [GeV];pz [GeV]", 100,-30.5,30.5,100,-30.5,30.5,100,-30.5,30.5);
while(!hepmc_input.failed()) {
// Read event from input file
hepmc_input.read_event(evt);
// If reading failed - exit loop
if( hepmc_input.failed() ) break;
for(const auto& v : evt.vertices() ) {
for(const auto& p : v->particles_out() ) {
if(p->pid() == 111) {
h_pi0_energy->Fill(p->momentum().e());
h_pi0_eta->Fill(p->momentum().eta());
h_pi0_theta->Fill(p->momentum().theta()*TMath::RadToDeg());
h_pi0_phi->Fill(p->momentum().phi()*TMath::RadToDeg());
h_pi0_pzpt->Fill(TMath::Sqrt(p->momentum().px()*p->momentum().px() + p->momentum().py()*p->momentum().py()),p->momentum().pz());
h_pi0_pxpy->Fill(p->momentum().px(),p->momentum().py());
h_pi0_p->Fill(p->momentum().px(),p->momentum().py(),p->momentum().pz());
}
}
}
evt.clear();
events_parsed++;
}
std::cout << "Events parsed and written: " << events_parsed << std::endl;
TCanvas* c = new TCanvas("c","c", 500, 500);
h_pi0_energy->GetYaxis()->SetTitleOffset(1.8);
h_pi0_energy->SetLineWidth(2);
h_pi0_energy->SetLineColor(kBlue);
h_pi0_energy->DrawClone();
c->SaveAs("results/input_pi0_energy_0GeVto30GeV.png");
c->SaveAs("results/input_pi0_energy_0GeVto30GeV.pdf");
TCanvas* c1 = new TCanvas("c1","c1", 500, 500);
h_pi0_eta->GetYaxis()->SetTitleOffset(1.9);
h_pi0_eta->SetLineWidth(2);
h_pi0_eta->SetLineColor(kBlue);
h_pi0_eta->DrawClone();
c1->SaveAs("results/input_pi0_eta_0GeVto30GeV.png");
c1->SaveAs("results/input_pi0_eta_0GeVto30GeV.pdf");
TCanvas* c2 = new TCanvas("c2","c2", 500, 500);
h_pi0_theta->GetYaxis()->SetTitleOffset(1.8);
h_pi0_theta->SetLineWidth(2);
h_pi0_theta->SetLineColor(kBlue);
h_pi0_theta->DrawClone();
c2->SaveAs("results/input_pi0_theta_0GeVto30GeV.png");
c2->SaveAs("results/input_pi0_theta_0GeVto30GeV.pdf");
TCanvas* c3 = new TCanvas("c3","c3", 500, 500);
h_pi0_phi->GetYaxis()->SetTitleOffset(1.8);
h_pi0_phi->SetLineWidth(2);
h_pi0_phi->GetYaxis()->SetRangeUser(0.0,h_pi0_phi->GetMaximum()+100.0);
h_pi0_phi->SetLineColor(kBlue);
h_pi0_phi->DrawClone();
c3->SaveAs("results/input_pi0_phi_0GeVto30GeV.png");
c3->SaveAs("results/input_pi0_phi_0GeVto30GeV.pdf");
TCanvas* c4 = new TCanvas("c4","c4", 500, 500);
h_pi0_pzpt->GetYaxis()->SetTitleOffset(1.4);
h_pi0_pzpt->SetLineWidth(2);
h_pi0_pzpt->SetLineColor(kBlue);
h_pi0_pzpt->DrawClone("COLZ");
c4->SaveAs("results/input_pi0_pzpt_0GeVto30GeV.png");
c4->SaveAs("results/input_pi0_pzpt_0GeVto30GeV.pdf");
TCanvas* c5 = new TCanvas("c5","c5", 500, 500);
h_pi0_pxpy->GetYaxis()->SetTitleOffset(1.4);
h_pi0_pxpy->SetLineWidth(2);
h_pi0_pxpy->SetLineColor(kBlue);
h_pi0_pxpy->DrawClone("COLZ");
c5->SaveAs("results/input_pi0_pxpy_0GeVto30GeV.png");
c5->SaveAs("results/input_pi0_pxpy_0GeVto30GeV.pdf");
TCanvas* c6 = new TCanvas("c6","c6", 500, 500);
h_pi0_p->GetYaxis()->SetTitleOffset(1.8);
h_pi0_p->GetXaxis()->SetTitleOffset(1.6);
h_pi0_p->GetZaxis()->SetTitleOffset(1.6);
h_pi0_p->SetLineWidth(2);
h_pi0_p->SetLineColor(kBlue);
h_pi0_p->DrawClone();
c6->SaveAs("results/input_pi0_p_0GeVto30GeV.png");
c6->SaveAs("results/input_pi0_p_0GeVto30GeV.pdf");
}
////////////////////////////////////////
// Read reconstruction ROOT output file
// Plot variables
////////////////////////////////////////
int makeplot_pi0(const char* input_fname = "../sim_output/sim_emcal_pi0s_output.root")
{
// Setting figures
gROOT->SetStyle("Plain");
gStyle->SetLineWidth(3);
gStyle->SetOptStat("nem");
gStyle->SetPadTickX(1);
gStyle->SetPadTickY(1);
gStyle->SetPadGridX(1);
gStyle->SetPadGridY(1);
gStyle->SetPadLeftMargin(0.14);
gStyle->SetPadRightMargin(0.14);
// Input ROOT file
TFile *f = new TFile(input_fname,"READ");
TTree *t = (TTree *)f->Get("events");
// Set Branch status and addressed
t->SetMakeClass(1);
t->SetBranchStatus("*", 0);
Int_t EcalClusters_;
t->SetBranchStatus("EcalClusters", 1);
t->SetBranchAddress("EcalClusters", &EcalClusters_);
Int_t RecoEcalHits_;
t->SetBranchStatus("RecoEcalHits", 1);
t->SetBranchAddress("RecoEcalHits", &RecoEcalHits_);
const Int_t kMaxEcalClusters = 100000;
Double_t cluster_x_pos[kMaxEcalClusters];
Double_t cluster_y_pos[kMaxEcalClusters];
Double_t cluster_z_pos[kMaxEcalClusters];
Float_t cluster_energy[kMaxEcalClusters];
t->SetBranchStatus("EcalClusters.position.x",1);
t->SetBranchStatus("EcalClusters.position.y",1);
t->SetBranchStatus("EcalClusters.position.z",1);
t->SetBranchStatus("EcalClusters.energy",1);
t->SetBranchAddress("EcalClusters.position.x",cluster_x_pos);
t->SetBranchAddress("EcalClusters.position.y",cluster_y_pos);
t->SetBranchAddress("EcalClusters.position.z",cluster_z_pos);
t->SetBranchAddress("EcalClusters.energy",cluster_energy);
const Int_t kMaxRecoEcalHits = 100000;
Double_t rec_x_pos[kMaxRecoEcalHits];
Double_t rec_y_pos[kMaxRecoEcalHits];
Double_t rec_energy[kMaxRecoEcalHits];
t->SetBranchStatus("RecoEcalHits.position.x",1);
t->SetBranchStatus("RecoEcalHits.position.y",1);
t->SetBranchStatus("RecoEcalHits.energy",1);
t->SetBranchAddress("RecoEcalHits.position.x",rec_x_pos);
t->SetBranchAddress("RecoEcalHits.position.y",rec_y_pos);
t->SetBranchAddress("RecoEcalHits.energy",rec_energy);
// Setting for Canvas
TCanvas *c1 = new TCanvas("c1", "c1", 500, 500);
TCanvas *c2 = new TCanvas("c2", "c2", 500, 500);
TCanvas *c3 = new TCanvas("c3", "c3", 500, 500);
TCanvas *c4 = new TCanvas("c4", "c4", 500, 500);
TCanvas *c5 = new TCanvas("c5", "c5", 500, 500);
TCanvas *c6 = new TCanvas("c6", "c6", 500, 500);
TCanvas *c7 = new TCanvas("c7", "c7", 500, 500);
TCanvas *c8 = new TCanvas("c8", "c8", 500, 500);
TCanvas *c9 = new TCanvas("c9", "c9", 500, 500);
TCanvas *c10 = new TCanvas("c10","c10", 500, 500);
TCanvas *c11 = new TCanvas("c11","c11", 500, 500);
TCanvas *c12 = new TCanvas("c12","c12", 500, 500);
TCanvas *c13 = new TCanvas("c13","c13", 500, 500);
// Declare histograms
TH1D *h1 = new TH1D("h1", "Scattering Angle(#theta)", 100,130.0,180.0);
TH1D *h2 = new TH1D("h2", "Pseudo-rapidity(#eta)", 100,-5.0,0.0);
TH2D *h3 = new TH2D("h3", "Cluster E vs Pseudo-rapidity", 100,-0.5,30.5,100,-5.0,0.0);
TH1D *h4 = new TH1D("h4", "Reconstructed energy per event", 100,-0.5,30.5);
TH1D *h5 = new TH1D("h5", "Number of Clusters per event", 5,-0.5,4.5);
TH1D *h6 = new TH1D("h6", "Scattering Angle(#theta) with CUT", 100,130.0,180.0);
TH1D *h7 = new TH1D("h7", "Pseudo-rapidity(#eta) with CUT", 100,-5.0,0.0);
TH2D *h8 = new TH2D("h8", "Cluster Hit Position", 62,-62.0,62.0,62,-62.0,62.0);
TH2D *h9 = new TH2D("h9", "All Hit Position", 62,-62.0,62.0,62,-62.0,62.0);
TH1D *h10 = new TH1D("h10","Invariant mass", 60,0.0,300.0);
TH1D *h11 = new TH1D("h11","E1", 100,-0.5,30.5);
TH1D *h12 = new TH1D("h12","E2", 100,-0.5,30.5);
TH1D *h13 = new TH1D("h13","angle", 100,0.0,90.0);
// Declare ellipse for boundary of crystal calorimeter
TEllipse *ell1 = new TEllipse(0.0, 0.0, 60.0, 60.0);
ell1->SetFillStyle(4000);
TEllipse *ell2 = new TEllipse(0.0, 0.0, 12.0, 12.0);
ell2->SetFillStyle(4000);
// Total number of entries
Int_t nentries = t->GetEntries();
// Variables are used in calculation
Double_t r; // Radius [cm]
Double_t phi; // Azimuth [degree]
Double_t theta; // Inclination [degree]
Double_t eta; // Pseudo-rapidity [unitless]
Float_t cluster_e; // Cluster energy [GeV]
Float_t total_cluster_e; // Add up clusters per event [GeV]
Double_t dot_product_pos_clusters; // dot product of positions of two photons
Double_t mag_pos2_cluster_1; // squared magnitude of position
Double_t mag_pos2_cluster_2; // squared magnitude of position
Double_t cosine_clusters; // cos(theta_photons)
Double_t theta_photons; // angle between two photons
Double_t invariant_mass; // M^2 = 2 * p_1 * p_2 * (1 - cos(theta_photons))
// Loop over event by event
for (int ievent = 0; ievent < nentries; ievent++)
{
t->GetEntry(ievent);
Int_t ncluster = EcalClusters_;
Int_t nreconhits = RecoEcalHits_;
total_cluster_e = 0.0;
h5->Fill(ncluster, 1.0);
// Loop over cluster by cluster
for (int icluster=0; icluster < ncluster; icluster++)
{
r = TMath::Sqrt((cluster_x_pos[icluster]*cluster_x_pos[icluster]) +
(cluster_y_pos[icluster]*cluster_y_pos[icluster]) +
(cluster_z_pos[icluster]*cluster_z_pos[icluster]));
phi = TMath::ATan(cluster_y_pos[icluster]/cluster_x_pos[icluster]) * TMath::RadToDeg();
theta = TMath::ACos(cluster_z_pos[icluster] / r) * TMath::RadToDeg();
eta = -1.0 * TMath::Log(TMath::Tan((theta*TMath::DegToRad())/2.0));
cluster_e = cluster_energy[icluster] / 1.e+3;
total_cluster_e += cluster_e;
}
// Select events with two clusters
// To calculate invariant mass
// M^2 = 2p1p2(1-cos(theta))
// p1 = E1
// p2 = E2
// theta: angle between two photons
if(ncluster == 2)
{
h1->Fill(theta, 1.0);
h2->Fill(eta, 1.0);
h3->Fill(cluster_e, eta, 1.0);
h4->Fill(total_cluster_e, 1.0);
h8->Fill(cluster_x_pos[0],cluster_y_pos[0], 1.0);
if(total_cluster_e > 0.5)
{
h6->Fill(theta, 1.0);
h7->Fill(eta, 1.0);
}
// Calculate invariant mass
dot_product_pos_clusters = cluster_x_pos[0]*cluster_x_pos[1] + cluster_y_pos[0]*cluster_y_pos[1] + cluster_z_pos[0]*cluster_z_pos[1];
mag_pos2_cluster_1 = (cluster_x_pos[0]*cluster_x_pos[0]) + (cluster_y_pos[0]*cluster_y_pos[0]) + (cluster_z_pos[0]*cluster_z_pos[0]);
mag_pos2_cluster_2 = (cluster_x_pos[1]*cluster_x_pos[1]) + (cluster_y_pos[1]*cluster_y_pos[1]) + (cluster_z_pos[1]*cluster_z_pos[1]);
cosine_clusters = (dot_product_pos_clusters/TMath::Sqrt(mag_pos2_cluster_1*mag_pos2_cluster_2));
theta_photons = TMath::ACos(cosine_clusters)*TMath::RadToDeg();
invariant_mass = TMath::Sqrt(2.0*cluster_energy[0]*cluster_energy[1]*(1.0 - cosine_clusters));
// Fill histograms
h10->Fill(invariant_mass, 1.0);
h11->Fill(cluster_energy[0], 1.0);
h12->Fill(cluster_energy[1], 1.0);
h13->Fill(theta_photons, 1.0);
}
// Loop over hit by hit
for(int ireconhit=0; ireconhit < nreconhits; ireconhit++)
h9->Fill(rec_x_pos[ireconhit],rec_y_pos[ireconhit], 1.0);
}
// Drawing and Saving figures
c1->cd();
h1->SetLineColor(kBlue);
h1->SetLineWidth(2);
h1->GetYaxis()->SetRangeUser(0.0,h1->GetMaximum()+10.0);
h1->GetXaxis()->SetTitle("#theta [degree]");
h1->GetYaxis()->SetTitle("events");
h1->GetYaxis()->SetTitleOffset(1.4);
h1->DrawClone();
c1->SaveAs("results/pi0_theta_hist_0GeVto30GeV.png");
c1->SaveAs("results/pi0_theta_hist_0GeVto30GeV.pdf");
c2->cd();
h2->SetLineColor(kBlue);
h2->SetLineWidth(2);
h2->GetYaxis()->SetRangeUser(0.0,h2->GetMaximum()+10.0);
h2->GetXaxis()->SetTitle("#eta");
h2->GetYaxis()->SetTitle("events");
h2->GetYaxis()->SetTitleOffset(1.4);
h2->DrawClone();
c2->SaveAs("results/pi0_eta_hist_0GeVto30GeV.png");
c2->SaveAs("results/pi0_eta_hist_0GeVto30GeV.pdf");
c3->cd();
h3->GetXaxis()->SetTitle("Cluster energy [GeV]");
h3->GetYaxis()->SetTitle("#eta");
h3->GetYaxis()->SetTitleOffset(1.4);
h3->DrawClone("COLZ");
c3->SaveAs("results/pi0_E_vs_eta_hist_0GeVto30GeV.png");
c3->SaveAs("results/pi0_E_vs_eta_hist_0GeVto30GeV.pdf");
c4->cd();
c4->SetLogy(1);
h4->SetLineColor(kBlue);
h4->SetLineWidth(2);
h4->GetXaxis()->SetTitle("reconstructed energy [GeV]");
h4->GetYaxis()->SetTitle("events");
h4->GetYaxis()->SetTitleOffset(1.4);
h4->DrawClone();
c4->SaveAs("results/pi0_Erec_hist_0GeVto30GeV.png");
c4->SaveAs("results/pi0_Erec_hist_0GeVto30GeV.pdf");
c5->cd();
c5->SetLogy(1);
h5->SetLineColor(kBlue);
h5->SetLineWidth(2);
h5->GetXaxis()->SetTitle("Number of Clusters");
h5->GetYaxis()->SetTitle("events");
h5->GetYaxis()->SetTitleOffset(1.4);
h5->DrawClone();
c5->SaveAs("results/pi0_ncluster_hist_0GeVto30GeV.png");
c5->SaveAs("results/pi0_ncluster_hist_0GeVto30GeV.pdf");
c6->cd();
h6->SetLineColor(kBlue);
h6->SetLineWidth(2);
h6->GetYaxis()->SetRangeUser(0.0,h1->GetMaximum()+10.0);
h6->GetXaxis()->SetTitle("#theta [degree]");
h6->GetYaxis()->SetTitle("events");
h6->GetYaxis()->SetTitleOffset(1.4);
h6->DrawClone();
c6->SaveAs("results/pi0_theta_hist_CUT_0GeVto30GeV.png");
c6->SaveAs("results/pi0_theta_hist_CUT_0GeVto30GeV.pdf");
c7->cd();
h7->SetLineColor(kBlue);
h7->SetLineWidth(2);
h7->GetYaxis()->SetRangeUser(0.0,h2->GetMaximum()+10.0);
h7->GetXaxis()->SetTitle("#eta");
h7->GetYaxis()->SetTitle("events");
h7->GetYaxis()->SetTitleOffset(1.4);
h7->DrawClone();
c7->SaveAs("results/pi0_eta_hist_CUT_0GeVto30GeV.png");
c7->SaveAs("results/pi0_eta_hist_CUT_0GeVto30GeV.pdf");
c8->cd();
h8->GetXaxis()->SetTitle("Hit position X [cm]");
h8->GetYaxis()->SetTitle("Hit position Y [cm]");
h8->GetYaxis()->SetTitleOffset(1.4);
h8->DrawClone("COLZ");
ell1->Draw("same");
ell2->Draw("same");
c8->SaveAs("results/pi0_hit_pos_cluster_0GeVto30GeV.png");
c8->SaveAs("results/pi0_hit_pos_cluster_0GeVto30GeV.pdf");
c9->cd();
h9->GetXaxis()->SetTitle("Hit position X [cm]");
h9->GetYaxis()->SetTitle("Hit position Y [cm]");
h9->GetYaxis()->SetTitleOffset(1.4);
h9->DrawClone("COLZ");
ell1->Draw("same");
ell2->Draw("same");
c9->SaveAs("results/pi0_hit_pos_all_0GeVto30GeV.png");
c9->SaveAs("results/pi0_hit_pos_all_0GeVto30GeV.pdf");
c10->cd();
h10->SetLineColor(kBlue);
h10->SetLineWidth(2);
h10->GetXaxis()->SetTitle("Invariant mass [MeV]");
h10->GetYaxis()->SetTitle("events");
h10->GetYaxis()->SetTitleOffset(1.4);
h10->DrawClone();
c10->SaveAs("results/pi0_invariant_mass_hist_0GeVto30GeV.png");
c10->SaveAs("results/pi0_invariant_mass_hist_0GeVto30GeV.pdf");
c11->cd();
h11->SetLineColor(kBlue);
h11->SetLineWidth(2);
h11->GetXaxis()->SetTitle("Cluster energy 1 [MeV]");
h11->GetYaxis()->SetTitle("events");
h11->GetYaxis()->SetTitleOffset(1.4);
h11->DrawClone();
c11->SaveAs("results/pi0_E1_hist_0GeVto30GeV.png");
c11->SaveAs("results/pi0_E1_hist_0GeVto30GeV.pdf");
c12->cd();
h12->SetLineColor(kBlue);
h12->SetLineWidth(2);
h12->GetXaxis()->SetTitle("Cluster energy 2 [MeV]");
h12->GetYaxis()->SetTitle("events");
h12->GetYaxis()->SetTitleOffset(1.4);
h12->DrawClone();
c12->SaveAs("results/pi0_E2_hist_0GeVto30GeV.png");
c12->SaveAs("results/pi0_E2_hist_0GeVto30GeV.pdf");
c13->cd();
h13->SetLineColor(kBlue);
h13->SetLineWidth(2);
h13->GetXaxis()->SetTitle("Angle between two photons [degree]");
h13->GetYaxis()->SetTitle("events");
h13->GetYaxis()->SetTitleOffset(1.4);
h13->DrawClone();
c13->SaveAs("results/pi0_angle_2photons_0GeVto30GeV.png");
c13->SaveAs("results/pi0_angle_2photons_0GeVto30GeV.pdf");
return 0;
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment