From 2ce5edf86ae3cd93f973d293de11db7cc9fa5932 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 17:00:26 -0400 Subject: [PATCH 01/46] added muon benchmarks for zdc and insert --- .gitlab-ci.yml | 6 +- Snakefile | 2 + benchmarks/insert_muon/Snakefile | 57 ++++++++ benchmarks/insert_muon/Snakefile~ | 58 ++++++++ .../insert_muon/analysis/gen_particles.cxx | 127 ++++++++++++++++++ benchmarks/insert_muon/analysis/muon_plots.py | 0 benchmarks/insert_muon/config.yml | 34 +++++ benchmarks/insert_muon/config.yml~ | 34 +++++ benchmarks/zdc_muon/Snakefile | 57 ++++++++ benchmarks/zdc_muon/Snakefile~ | 57 ++++++++ .../zdc_muon/analysis/gen_particles.cxx | 127 ++++++++++++++++++ benchmarks/zdc_muon/analysis/muon_plots.py | 0 benchmarks/zdc_muon/config.yml | 34 +++++ benchmarks/zdc_muon/config.yml~ | 34 +++++ 14 files changed, 626 insertions(+), 1 deletion(-) create mode 100644 benchmarks/insert_muon/Snakefile create mode 100644 benchmarks/insert_muon/Snakefile~ create mode 100644 benchmarks/insert_muon/analysis/gen_particles.cxx create mode 100644 benchmarks/insert_muon/analysis/muon_plots.py create mode 100644 benchmarks/insert_muon/config.yml create mode 100644 benchmarks/insert_muon/config.yml~ create mode 100644 benchmarks/zdc_muon/Snakefile create mode 100644 benchmarks/zdc_muon/Snakefile~ create mode 100644 benchmarks/zdc_muon/analysis/gen_particles.cxx create mode 100644 benchmarks/zdc_muon/analysis/muon_plots.py create mode 100644 benchmarks/zdc_muon/config.yml create mode 100644 benchmarks/zdc_muon/config.yml~ diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index 4e9d681c..1248cd6d 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -44,6 +44,7 @@ stages: - simulate - calibrate - benchmarks + - analyze - collect - deploy - status-report @@ -145,7 +146,8 @@ include: - local: 'benchmarks/pid/config.yml' - local: 'benchmarks/timing/config.yml' - local: 'benchmarks/b0_tracker/config.yml' - + - local: 'benchmarks/insert_muon/config.yml' + - local: 'benchmarks/zdc_muon/config.yml' deploy_results: allow_failure: true stage: deploy @@ -161,6 +163,8 @@ deploy_results: - "collect_results:tracking_performance_campaigns" - "collect_results:zdc" - "collect_results:zdc_lyso" + - "collect_results:insert_muon" + - "collect_results:zdc_muon" script: - echo "deploy results!" - find results -print | sort | tee summary.txt diff --git a/Snakefile b/Snakefile index d153297c..1069d32b 100644 --- a/Snakefile +++ b/Snakefile @@ -7,6 +7,8 @@ include: "benchmarks/ecal_gaps/Snakefile" include: "benchmarks/material_scan/Snakefile" include: "benchmarks/tracking_performances/Snakefile" include: "benchmarks/zdc_lyso/Snakefile" +include: "benchmarks/insert_muon/Snakefile" +include: "benchmarks/zdc_muon/Snakefile" use_s3 = config["remote_provider"].lower() == "s3" use_xrootd = config["remote_provider"].lower() == "xrootd" diff --git a/benchmarks/insert_muon/Snakefile b/benchmarks/insert_muon/Snakefile new file mode 100644 index 00000000..b5b989d7 --- /dev/null +++ b/benchmarks/insert_muon/Snakefile @@ -0,0 +1,57 @@ +rule insert_muon_generate: + input: + script="benchmarks/insert_muon/analysis/gen_particles.cxx", + params: + NEVENTS_GEN=1000, + th_max=5.7, + th_min=2.0 + output: + GEN_FILE="sim_output/insert_muon/mu-_{P}GeV.hepmc" + shell: + """ +root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "mu-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' +""" + +rule insert_muon_simulate: + input: + GEN_FILE="sim_output/insert_muon/mu-_{P}GeV.hepmc" + params: + PHYSICS_LIST="FTFP_BERT" + output: + SIM_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_SIM=1000 +# Running simulation +npsim \ + --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ + --numberOfEvents $NEVENTS_SIM \ + --physicsList {params.PHYSICS_LIST} \ + --inputFiles {input.GEN_FILE} \ + --outputFile {output.SIM_FILE} +""" + +rule insert_muon_recon: + input: + SIM_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_sim_muon_{P}GeV.edm4hep.root" + output: + REC_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_REC=1000 +eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters -Pjana:nevents=$NEVENTS_REC +""" + +rule insert_muon_analysis: + input: + expand("sim_output/insert_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root", + P=[20, 30, 40, 50, 60, 70, 80], + DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), + script="benchmarks/insert_muon/analysis/muon_plots.py", + output: + results_dir=directory("results/{DETECTOR_CONFIG}/insert_muon"), + shell: + """ +mkdir -p {output.results_dir} +python {input.script} {output.results_dir} +""" diff --git a/benchmarks/insert_muon/Snakefile~ b/benchmarks/insert_muon/Snakefile~ new file mode 100644 index 00000000..3f24d738 --- /dev/null +++ b/benchmarks/insert_muon/Snakefile~ @@ -0,0 +1,58 @@ +rule neutron_generate: + input: + script="benchmarks/neutron/analysis/gen_particles.cxx", + params: + NEVENTS_GEN=1000, + th_max=5.7, + th_min=2.0 + output: + GEN_FILE="sim_output/neutron/neutron_{P}GeV.hepmc" + shell: + """ +mkdir -p sim_output/neutron +root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "neutron", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' +""" + +rule neutron_simulate: + input: + GEN_FILE="sim_output/neutron/neutron_{P}GeV.hepmc" + params: + PHYSICS_LIST="FTFP_BERT" + output: + SIM_FILE="sim_output/neutron/{DETECTOR_CONFIG}_sim_neutron_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_SIM=1000 +# Running simulation +npsim \ + --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ + --numberOfEvents $NEVENTS_SIM \ + --physicsList {params.PHYSICS_LIST} \ + --inputFiles {input.GEN_FILE} \ + --outputFile {output.SIM_FILE} +""" + +rule neutron_recon: + input: + SIM_FILE="sim_output/neutron/{DETECTOR_CONFIG}_sim_neutron_{P}GeV.edm4hep.root" + output: + REC_FILE="sim_output/neutron/{DETECTOR_CONFIG}_rec_neutron_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_REC=1000 +eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters -Pjana:nevents=$NEVENTS_REC +""" + +rule neutron_analysis: + input: + expand("sim_output/neutron/{DETECTOR_CONFIG}_rec_neutron_{P}GeV.edm4hep.root", + P=[20, 30, 40, 50, 60, 70, 80], + DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), + script="benchmarks/neutron/analysis/neutron_plots.py", + output: + results_dir=directory("results/{DETECTOR_CONFIG}/neutron"), + shell: + """ +mkdir -p {output.results_dir} +python {input.script} {output.results_dir} +""" diff --git a/benchmarks/insert_muon/analysis/gen_particles.cxx b/benchmarks/insert_muon/analysis/gen_particles.cxx new file mode 100644 index 00000000..0877a852 --- /dev/null +++ b/benchmarks/insert_muon/analysis/gen_particles.cxx @@ -0,0 +1,127 @@ +#include "HepMC3/GenEvent.h" +#include "HepMC3/ReaderAscii.h" +#include "HepMC3/WriterAscii.h" +#include "HepMC3/Print.h" + +#include "TRandom3.h" +#include "TVector3.h" + +#include +#include +#include +#include +#include +#include +#include + +using namespace HepMC3; + +// Generate single electron as input to the Insert simulation. +// -- +// We generate events with a constant polar angle with respect to +// the proton direction and then rotate so that the events are given +// in normal lab coordinate system +// -- +void gen_particles( + int n_events = 1000, + const char* out_fname = "gen_particles.hepmc", + TString particle_name = "e-", + double th_min = 3., // Minimum polar angle, in degrees + double th_max = 3., // Maximum polar angle, in degrees + double phi_min = 0., // Minimum azimuthal angle, in degrees + double phi_max = 360., // Maximum azimuthal angle, in degrees + double p = 10., // Momentum in GeV/c + int dist = 0, //Momentum distribution: 0=fixed, 1=uniform, 2=Gaussian + int useCrossingAngle=1 // 0= no rotation, 1 = -25 mrad + ) +{ + WriterAscii hepmc_output(out_fname); + int events_parsed = 0; + GenEvent evt(Units::GEV, Units::MM); + + // Random number generator + TRandom3 *r1 = new TRandom3(0); //Use time as random seed + + // Getting generated particle information + TDatabasePDG *pdg = new TDatabasePDG(); + TParticlePDG *particle = pdg->GetParticle(particle_name); + const double mass = particle->Mass(); + const int pdgID = particle->PdgCode(); + + for (events_parsed = 0; events_parsed < n_events; events_parsed++) { + + //Set the event number + evt.set_event_number(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(FourVector(0.0, 0.0, 10.0, 10.0), 11, 4); + GenParticlePtr p2 = std::make_shared( + FourVector(0.0, 0.0, 0.0, 0.938), 2212, 4); + + // Define momentum with respect to proton direction + double phi = r1->Uniform(phi_min*TMath::DegToRad(),phi_max*TMath::DegToRad()); + double th = r1->Uniform(th_min*TMath::DegToRad(),th_max*TMath::DegToRad()); + + //Total momentum distribution + double pevent = -1; + if(dist==0){ //fixed + pevent = p; + } + else if(dist==1){ //Uniform: +-50% variation + pevent = p*(1. + r1->Uniform(-0.5,0.5) ); + } + else if(dist==2){ //Gaussian: Sigma = 0.1*mean + while(pevent<0) //Avoid negative values + pevent = r1->Gaus(p,0.1*p); + } + + double px = pevent * std::cos(phi) * std::sin(th); + double py = pevent * std::sin(phi) * std::sin(th); + double pz = pevent * std::cos(th); + TVector3 pvec(px,py,pz); + + //Rotate to lab coordinate system + double cross_angle = -25./1000.*useCrossingAngle; //in Rad + TVector3 pbeam_dir(sin(cross_angle),0,cos(cross_angle)); //proton beam direction + pvec.RotateY(-pbeam_dir.Theta()); // Theta is returned positive, beam in negative X + // type 1 is final state + // pdgid 11 - electron 0.510 MeV/c^2 + GenParticlePtr p3 = std::make_shared( + FourVector( + pvec.X(), pvec.Y(), pvec.Z(), + sqrt(pevent*pevent + (mass * mass))), + pdgID, 1); + + //If wanted, set non-zero vertex + double vx = 0.; + double vy = 0.; + double vz = 0.; + double vt = 0.; + + GenVertexPtr v1 = std::make_shared(); + evt.shift_position_by(FourVector(vx, vy, vz, vt)); + 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 % 100 == 0) { + std::cout << "Event: " << events_parsed << std::endl; + } + evt.clear(); + } + hepmc_output.close(); + std::cout << "Events parsed and written: " << events_parsed << std::endl; +} diff --git a/benchmarks/insert_muon/analysis/muon_plots.py b/benchmarks/insert_muon/analysis/muon_plots.py new file mode 100644 index 00000000..e69de29b diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml new file mode 100644 index 00000000..15de8aa7 --- /dev/null +++ b/benchmarks/insert_muon/config.yml @@ -0,0 +1,34 @@ +simulate:insert_muon: + stage: simulate + extends: .phy_benchmark + needs: ["common:detector"] + parallel: + matrix: + - P: 20 + - P: 30 + - P: 40 + - P: 50 + - P: 60 + - P: 70 + - P: 80 + timeout: 1 hours + script: + - snakemake --cores 1 sim_output/insert_muon/epic_craterlake_rec_mu-_${P}GeV.edm4hep.root + retry: + max: 2 + when: + - runner_system_failure + +analyze:insert_muon: + stage: analyze + extends: .phy_benchmark + needs: ["simulate:insert_muon"] + script: + - mkdir -p results/epic_craterlake + - python benchmarks/insert_muon/analysis/muon_plots.py results/epic_craterlake/insert_muon + +collect_results:insert_muon: + stage: collect_results + needs: ["analyze:insert_muon"] + script: + - ls -al diff --git a/benchmarks/insert_muon/config.yml~ b/benchmarks/insert_muon/config.yml~ new file mode 100644 index 00000000..1510d771 --- /dev/null +++ b/benchmarks/insert_muon/config.yml~ @@ -0,0 +1,34 @@ +neutron:simulate: + stage: simulate + extends: .phy_benchmark + needs: ["common:detector"] + parallel: + matrix: + - P: 20 + - P: 30 + - P: 40 + - P: 50 + - P: 60 + - P: 70 + - P: 80 + timeout: 1 hours + script: + - snakemake --cores 1 sim_output/neutron/epic_craterlake_rec_neutron_${P}GeV.edm4hep.root + retry: + max: 2 + when: + - runner_system_failure + +neutron:analyze: + stage: analyze + extends: .phy_benchmark + needs: ["neutron:simulate"] + script: + - mkdir -p results/epic_craterlake + - python benchmarks/neutron/analysis/neutron_plots.py results/epic_craterlake/neutron + +neutron:results: + stage: collect_results + needs: ["neutron:analyze"] + script: + - ls -al diff --git a/benchmarks/zdc_muon/Snakefile b/benchmarks/zdc_muon/Snakefile new file mode 100644 index 00000000..8547a132 --- /dev/null +++ b/benchmarks/zdc_muon/Snakefile @@ -0,0 +1,57 @@ +rule zdc_muon_generate: + input: + script="benchmarks/zdc_muon/analysis/gen_particles.cxx", + params: + NEVENTS_GEN=1000, + th_max=5.7, + th_min=2.0 + output: + GEN_FILE="sim_output/zdc_muon/mu-_{P}GeV.hepmc" + shell: + """ +root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "mu-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' +""" + +rule zdc_muon_simulate: + input: + GEN_FILE="sim_output/zdc_muon/mu-_{P}GeV.hepmc" + params: + PHYSICS_LIST="FTFP_BERT" + output: + SIM_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_SIM=1000 +# Running simulation +npsim \ + --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ + --numberOfEvents $NEVENTS_SIM \ + --physicsList {params.PHYSICS_LIST} \ + --inputFiles {input.GEN_FILE} \ + --outputFile {output.SIM_FILE} +""" + +rule zdc_muon_recon: + input: + SIM_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_sim_muon_{P}GeV.edm4hep.root" + output: + REC_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_REC=1000 +eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalFarForwardZDCRecHits,HcalFarForwardZDCClusters,HcalFarForwardZDCSubcellHits, -Pjana:nevents=$NEVENTS_REC +""" + +rule zdc_muon_analysis: + input: + expand("sim_output/zdc_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root", + P=[20, 30, 40, 50, 60, 70, 80], + DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), + script="benchmarks/zdc_muon/analysis/muon_plots.py", + output: + results_dir=directory("results/{DETECTOR_CONFIG}/zdc_muon"), + shell: + """ +mkdir -p {output.results_dir} +python {input.script} {output.results_dir} +""" diff --git a/benchmarks/zdc_muon/Snakefile~ b/benchmarks/zdc_muon/Snakefile~ new file mode 100644 index 00000000..ab859b44 --- /dev/null +++ b/benchmarks/zdc_muon/Snakefile~ @@ -0,0 +1,57 @@ +rule insert_muon_generate: + input: + script="benchmarks/insert_muon/analysis/gen_particles.cxx", + params: + NEVENTS_GEN=1000, + th_max=5.7, + th_min=2.0 + output: + GEN_FILE="sim_output/insert_muon/mu-_{P}GeV.hepmc" + shell: + """ +root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "mu-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' +""" + +rule insert_muon_simulate: + input: + GEN_FILE="sim_output/insert_muon/mu-_{P}GeV.hepmc" + params: + PHYSICS_LIST="FTFP_BERT" + output: + SIM_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_SIM=1000 +# Running simulation +npsim \ + --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ + --numberOfEvents $NEVENTS_SIM \ + --physicsList {params.PHYSICS_LIST} \ + --inputFiles {input.GEN_FILE} \ + --outputFile {output.SIM_FILE} +""" + +rule insert_muon_recon: + input: + SIM_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_sim_muon_{P}GeV.edm4hep.root" + output: + REC_FILE="sim_output/inset_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_REC=1000 +eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters -Pjana:nevents=$NEVENTS_REC +""" + +rule insert_muon_analysis: + input: + expand("sim_output/insert_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root", + P=[20, 30, 40, 50, 60, 70, 80], + DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), + script="benchmarks/insert_muon/analysis/muon_plots.py", + output: + results_dir=directory("results/{DETECTOR_CONFIG}/insert_muon"), + shell: + """ +mkdir -p {output.results_dir} +python {input.script} {output.results_dir} +""" diff --git a/benchmarks/zdc_muon/analysis/gen_particles.cxx b/benchmarks/zdc_muon/analysis/gen_particles.cxx new file mode 100644 index 00000000..0877a852 --- /dev/null +++ b/benchmarks/zdc_muon/analysis/gen_particles.cxx @@ -0,0 +1,127 @@ +#include "HepMC3/GenEvent.h" +#include "HepMC3/ReaderAscii.h" +#include "HepMC3/WriterAscii.h" +#include "HepMC3/Print.h" + +#include "TRandom3.h" +#include "TVector3.h" + +#include +#include +#include +#include +#include +#include +#include + +using namespace HepMC3; + +// Generate single electron as input to the Insert simulation. +// -- +// We generate events with a constant polar angle with respect to +// the proton direction and then rotate so that the events are given +// in normal lab coordinate system +// -- +void gen_particles( + int n_events = 1000, + const char* out_fname = "gen_particles.hepmc", + TString particle_name = "e-", + double th_min = 3., // Minimum polar angle, in degrees + double th_max = 3., // Maximum polar angle, in degrees + double phi_min = 0., // Minimum azimuthal angle, in degrees + double phi_max = 360., // Maximum azimuthal angle, in degrees + double p = 10., // Momentum in GeV/c + int dist = 0, //Momentum distribution: 0=fixed, 1=uniform, 2=Gaussian + int useCrossingAngle=1 // 0= no rotation, 1 = -25 mrad + ) +{ + WriterAscii hepmc_output(out_fname); + int events_parsed = 0; + GenEvent evt(Units::GEV, Units::MM); + + // Random number generator + TRandom3 *r1 = new TRandom3(0); //Use time as random seed + + // Getting generated particle information + TDatabasePDG *pdg = new TDatabasePDG(); + TParticlePDG *particle = pdg->GetParticle(particle_name); + const double mass = particle->Mass(); + const int pdgID = particle->PdgCode(); + + for (events_parsed = 0; events_parsed < n_events; events_parsed++) { + + //Set the event number + evt.set_event_number(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(FourVector(0.0, 0.0, 10.0, 10.0), 11, 4); + GenParticlePtr p2 = std::make_shared( + FourVector(0.0, 0.0, 0.0, 0.938), 2212, 4); + + // Define momentum with respect to proton direction + double phi = r1->Uniform(phi_min*TMath::DegToRad(),phi_max*TMath::DegToRad()); + double th = r1->Uniform(th_min*TMath::DegToRad(),th_max*TMath::DegToRad()); + + //Total momentum distribution + double pevent = -1; + if(dist==0){ //fixed + pevent = p; + } + else if(dist==1){ //Uniform: +-50% variation + pevent = p*(1. + r1->Uniform(-0.5,0.5) ); + } + else if(dist==2){ //Gaussian: Sigma = 0.1*mean + while(pevent<0) //Avoid negative values + pevent = r1->Gaus(p,0.1*p); + } + + double px = pevent * std::cos(phi) * std::sin(th); + double py = pevent * std::sin(phi) * std::sin(th); + double pz = pevent * std::cos(th); + TVector3 pvec(px,py,pz); + + //Rotate to lab coordinate system + double cross_angle = -25./1000.*useCrossingAngle; //in Rad + TVector3 pbeam_dir(sin(cross_angle),0,cos(cross_angle)); //proton beam direction + pvec.RotateY(-pbeam_dir.Theta()); // Theta is returned positive, beam in negative X + // type 1 is final state + // pdgid 11 - electron 0.510 MeV/c^2 + GenParticlePtr p3 = std::make_shared( + FourVector( + pvec.X(), pvec.Y(), pvec.Z(), + sqrt(pevent*pevent + (mass * mass))), + pdgID, 1); + + //If wanted, set non-zero vertex + double vx = 0.; + double vy = 0.; + double vz = 0.; + double vt = 0.; + + GenVertexPtr v1 = std::make_shared(); + evt.shift_position_by(FourVector(vx, vy, vz, vt)); + 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 % 100 == 0) { + std::cout << "Event: " << events_parsed << std::endl; + } + evt.clear(); + } + hepmc_output.close(); + std::cout << "Events parsed and written: " << events_parsed << std::endl; +} diff --git a/benchmarks/zdc_muon/analysis/muon_plots.py b/benchmarks/zdc_muon/analysis/muon_plots.py new file mode 100644 index 00000000..e69de29b diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml new file mode 100644 index 00000000..f26cea79 --- /dev/null +++ b/benchmarks/zdc_muon/config.yml @@ -0,0 +1,34 @@ +simulate:zdc_muon: + stage: simulate + extends: .phy_benchmark + needs: ["common:detector"] + parallel: + matrix: + - P: 20 + - P: 30 + - P: 40 + - P: 50 + - P: 60 + - P: 70 + - P: 80 + timeout: 1 hours + script: + - snakemake --cores 1 sim_output/zdc_muon/epic_craterlake_rec_mu-_${P}GeV.edm4hep.root + retry: + max: 2 + when: + - runner_system_failure + +analyze:zdc_muon: + stage: analyze + extends: .phy_benchmark + needs: ["simulate:zdc_muon"] + script: + - mkdir -p results/epic_craterlake + - python benchmarks/zdc_muon/analysis/muon_plots.py results/epic_craterlake/zdc_muon + +collect_results:zdc_muon: + stage: collect_results + needs: ["analyze:zdc_muon"] + script: + - ls -al diff --git a/benchmarks/zdc_muon/config.yml~ b/benchmarks/zdc_muon/config.yml~ new file mode 100644 index 00000000..15de8aa7 --- /dev/null +++ b/benchmarks/zdc_muon/config.yml~ @@ -0,0 +1,34 @@ +simulate:insert_muon: + stage: simulate + extends: .phy_benchmark + needs: ["common:detector"] + parallel: + matrix: + - P: 20 + - P: 30 + - P: 40 + - P: 50 + - P: 60 + - P: 70 + - P: 80 + timeout: 1 hours + script: + - snakemake --cores 1 sim_output/insert_muon/epic_craterlake_rec_mu-_${P}GeV.edm4hep.root + retry: + max: 2 + when: + - runner_system_failure + +analyze:insert_muon: + stage: analyze + extends: .phy_benchmark + needs: ["simulate:insert_muon"] + script: + - mkdir -p results/epic_craterlake + - python benchmarks/insert_muon/analysis/muon_plots.py results/epic_craterlake/insert_muon + +collect_results:insert_muon: + stage: collect_results + needs: ["analyze:insert_muon"] + script: + - ls -al From 2f4ed37d1bffd8ba7eeadbf7e64a7fca042a0154 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 17:00:41 -0400 Subject: [PATCH 02/46] removed ~ files --- benchmarks/insert_muon/Snakefile~ | 58 ------------------------------ benchmarks/insert_muon/config.yml~ | 34 ------------------ benchmarks/zdc_muon/Snakefile~ | 57 ----------------------------- benchmarks/zdc_muon/config.yml~ | 34 ------------------ 4 files changed, 183 deletions(-) delete mode 100644 benchmarks/insert_muon/Snakefile~ delete mode 100644 benchmarks/insert_muon/config.yml~ delete mode 100644 benchmarks/zdc_muon/Snakefile~ delete mode 100644 benchmarks/zdc_muon/config.yml~ diff --git a/benchmarks/insert_muon/Snakefile~ b/benchmarks/insert_muon/Snakefile~ deleted file mode 100644 index 3f24d738..00000000 --- a/benchmarks/insert_muon/Snakefile~ +++ /dev/null @@ -1,58 +0,0 @@ -rule neutron_generate: - input: - script="benchmarks/neutron/analysis/gen_particles.cxx", - params: - NEVENTS_GEN=1000, - th_max=5.7, - th_min=2.0 - output: - GEN_FILE="sim_output/neutron/neutron_{P}GeV.hepmc" - shell: - """ -mkdir -p sim_output/neutron -root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "neutron", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' -""" - -rule neutron_simulate: - input: - GEN_FILE="sim_output/neutron/neutron_{P}GeV.hepmc" - params: - PHYSICS_LIST="FTFP_BERT" - output: - SIM_FILE="sim_output/neutron/{DETECTOR_CONFIG}_sim_neutron_{P}GeV.edm4hep.root" - shell: - """ -NEVENTS_SIM=1000 -# Running simulation -npsim \ - --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ - --numberOfEvents $NEVENTS_SIM \ - --physicsList {params.PHYSICS_LIST} \ - --inputFiles {input.GEN_FILE} \ - --outputFile {output.SIM_FILE} -""" - -rule neutron_recon: - input: - SIM_FILE="sim_output/neutron/{DETECTOR_CONFIG}_sim_neutron_{P}GeV.edm4hep.root" - output: - REC_FILE="sim_output/neutron/{DETECTOR_CONFIG}_rec_neutron_{P}GeV.edm4hep.root" - shell: - """ -NEVENTS_REC=1000 -eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters -Pjana:nevents=$NEVENTS_REC -""" - -rule neutron_analysis: - input: - expand("sim_output/neutron/{DETECTOR_CONFIG}_rec_neutron_{P}GeV.edm4hep.root", - P=[20, 30, 40, 50, 60, 70, 80], - DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), - script="benchmarks/neutron/analysis/neutron_plots.py", - output: - results_dir=directory("results/{DETECTOR_CONFIG}/neutron"), - shell: - """ -mkdir -p {output.results_dir} -python {input.script} {output.results_dir} -""" diff --git a/benchmarks/insert_muon/config.yml~ b/benchmarks/insert_muon/config.yml~ deleted file mode 100644 index 1510d771..00000000 --- a/benchmarks/insert_muon/config.yml~ +++ /dev/null @@ -1,34 +0,0 @@ -neutron:simulate: - stage: simulate - extends: .phy_benchmark - needs: ["common:detector"] - parallel: - matrix: - - P: 20 - - P: 30 - - P: 40 - - P: 50 - - P: 60 - - P: 70 - - P: 80 - timeout: 1 hours - script: - - snakemake --cores 1 sim_output/neutron/epic_craterlake_rec_neutron_${P}GeV.edm4hep.root - retry: - max: 2 - when: - - runner_system_failure - -neutron:analyze: - stage: analyze - extends: .phy_benchmark - needs: ["neutron:simulate"] - script: - - mkdir -p results/epic_craterlake - - python benchmarks/neutron/analysis/neutron_plots.py results/epic_craterlake/neutron - -neutron:results: - stage: collect_results - needs: ["neutron:analyze"] - script: - - ls -al diff --git a/benchmarks/zdc_muon/Snakefile~ b/benchmarks/zdc_muon/Snakefile~ deleted file mode 100644 index ab859b44..00000000 --- a/benchmarks/zdc_muon/Snakefile~ +++ /dev/null @@ -1,57 +0,0 @@ -rule insert_muon_generate: - input: - script="benchmarks/insert_muon/analysis/gen_particles.cxx", - params: - NEVENTS_GEN=1000, - th_max=5.7, - th_min=2.0 - output: - GEN_FILE="sim_output/insert_muon/mu-_{P}GeV.hepmc" - shell: - """ -root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "mu-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' -""" - -rule insert_muon_simulate: - input: - GEN_FILE="sim_output/insert_muon/mu-_{P}GeV.hepmc" - params: - PHYSICS_LIST="FTFP_BERT" - output: - SIM_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root" - shell: - """ -NEVENTS_SIM=1000 -# Running simulation -npsim \ - --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ - --numberOfEvents $NEVENTS_SIM \ - --physicsList {params.PHYSICS_LIST} \ - --inputFiles {input.GEN_FILE} \ - --outputFile {output.SIM_FILE} -""" - -rule insert_muon_recon: - input: - SIM_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_sim_muon_{P}GeV.edm4hep.root" - output: - REC_FILE="sim_output/inset_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root" - shell: - """ -NEVENTS_REC=1000 -eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters -Pjana:nevents=$NEVENTS_REC -""" - -rule insert_muon_analysis: - input: - expand("sim_output/insert_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root", - P=[20, 30, 40, 50, 60, 70, 80], - DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), - script="benchmarks/insert_muon/analysis/muon_plots.py", - output: - results_dir=directory("results/{DETECTOR_CONFIG}/insert_muon"), - shell: - """ -mkdir -p {output.results_dir} -python {input.script} {output.results_dir} -""" diff --git a/benchmarks/zdc_muon/config.yml~ b/benchmarks/zdc_muon/config.yml~ deleted file mode 100644 index 15de8aa7..00000000 --- a/benchmarks/zdc_muon/config.yml~ +++ /dev/null @@ -1,34 +0,0 @@ -simulate:insert_muon: - stage: simulate - extends: .phy_benchmark - needs: ["common:detector"] - parallel: - matrix: - - P: 20 - - P: 30 - - P: 40 - - P: 50 - - P: 60 - - P: 70 - - P: 80 - timeout: 1 hours - script: - - snakemake --cores 1 sim_output/insert_muon/epic_craterlake_rec_mu-_${P}GeV.edm4hep.root - retry: - max: 2 - when: - - runner_system_failure - -analyze:insert_muon: - stage: analyze - extends: .phy_benchmark - needs: ["simulate:insert_muon"] - script: - - mkdir -p results/epic_craterlake - - python benchmarks/insert_muon/analysis/muon_plots.py results/epic_craterlake/insert_muon - -collect_results:insert_muon: - stage: collect_results - needs: ["analyze:insert_muon"] - script: - - ls -al From 7f6014d2754ddf02cfbfd160c0d48c96f2e4dfc3 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 17:08:05 -0400 Subject: [PATCH 03/46] fixed naming issue --- benchmarks/insert_muon/Snakefile | 6 +++--- benchmarks/zdc_muon/Snakefile | 6 +++--- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/benchmarks/insert_muon/Snakefile b/benchmarks/insert_muon/Snakefile index b5b989d7..cee43378 100644 --- a/benchmarks/insert_muon/Snakefile +++ b/benchmarks/insert_muon/Snakefile @@ -33,9 +33,9 @@ npsim \ rule insert_muon_recon: input: - SIM_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_sim_muon_{P}GeV.edm4hep.root" + SIM_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root" output: - REC_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root" + REC_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root" shell: """ NEVENTS_REC=1000 @@ -44,7 +44,7 @@ eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_fil rule insert_muon_analysis: input: - expand("sim_output/insert_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root", + expand("sim_output/insert_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root", P=[20, 30, 40, 50, 60, 70, 80], DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), script="benchmarks/insert_muon/analysis/muon_plots.py", diff --git a/benchmarks/zdc_muon/Snakefile b/benchmarks/zdc_muon/Snakefile index 8547a132..73fc4c3a 100644 --- a/benchmarks/zdc_muon/Snakefile +++ b/benchmarks/zdc_muon/Snakefile @@ -33,9 +33,9 @@ npsim \ rule zdc_muon_recon: input: - SIM_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_sim_muon_{P}GeV.edm4hep.root" + SIM_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root" output: - REC_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root" + REC_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root" shell: """ NEVENTS_REC=1000 @@ -44,7 +44,7 @@ eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_fil rule zdc_muon_analysis: input: - expand("sim_output/zdc_muon/{DETECTOR_CONFIG}_rec_muon_{P}GeV.edm4hep.root", + expand("sim_output/zdc_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root", P=[20, 30, 40, 50, 60, 70, 80], DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), script="benchmarks/zdc_muon/analysis/muon_plots.py", From 671460e1faf24cdffcd4450bb0197ed95699adef Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 19:10:17 -0400 Subject: [PATCH 04/46] Update benchmarks/zdc_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/zdc_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml index f26cea79..1f95ffe2 100644 --- a/benchmarks/zdc_muon/config.yml +++ b/benchmarks/zdc_muon/config.yml @@ -1,6 +1,6 @@ simulate:zdc_muon: stage: simulate - extends: .phy_benchmark + extends: .det_benchmark needs: ["common:detector"] parallel: matrix: From b737e106617b70aadb38f538abbe7a7b76190137 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 19:10:24 -0400 Subject: [PATCH 05/46] Update benchmarks/insert_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/insert_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index 15de8aa7..ab9c62f5 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -21,7 +21,7 @@ simulate:insert_muon: analyze:insert_muon: stage: analyze - extends: .phy_benchmark + extends: .det_benchmark needs: ["simulate:insert_muon"] script: - mkdir -p results/epic_craterlake From 623bea0b991efd32c1bda177d8554276747b10fb Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 19:10:31 -0400 Subject: [PATCH 06/46] Update benchmarks/insert_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/insert_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index ab9c62f5..146910d0 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -1,6 +1,6 @@ simulate:insert_muon: stage: simulate - extends: .phy_benchmark + extends: .det_benchmark needs: ["common:detector"] parallel: matrix: From fb5dbb63fe9b96b8aa488bd13634f658b7266929 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 19:10:40 -0400 Subject: [PATCH 07/46] Update benchmarks/zdc_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/zdc_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml index 1f95ffe2..496bb97a 100644 --- a/benchmarks/zdc_muon/config.yml +++ b/benchmarks/zdc_muon/config.yml @@ -21,7 +21,7 @@ simulate:zdc_muon: analyze:zdc_muon: stage: analyze - extends: .phy_benchmark + extends: .det_benchmark needs: ["simulate:zdc_muon"] script: - mkdir -p results/epic_craterlake From 14828bd3e994aff8610e919ac550408e887263ca Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 19:10:51 -0400 Subject: [PATCH 08/46] Update benchmarks/zdc_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/zdc_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml index 496bb97a..9621a10f 100644 --- a/benchmarks/zdc_muon/config.yml +++ b/benchmarks/zdc_muon/config.yml @@ -22,7 +22,7 @@ simulate:zdc_muon: analyze:zdc_muon: stage: analyze extends: .det_benchmark - needs: ["simulate:zdc_muon"] + needs: ["sim:zdc_muon"] script: - mkdir -p results/epic_craterlake - python benchmarks/zdc_muon/analysis/muon_plots.py results/epic_craterlake/zdc_muon From 9e5cdc225b222061878fa98ee90bab543fc75aa8 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 19:11:02 -0400 Subject: [PATCH 09/46] Update benchmarks/zdc_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/zdc_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml index 9621a10f..6021cdc1 100644 --- a/benchmarks/zdc_muon/config.yml +++ b/benchmarks/zdc_muon/config.yml @@ -1,4 +1,4 @@ -simulate:zdc_muon: +sim:zdc_muon: stage: simulate extends: .det_benchmark needs: ["common:detector"] From 1e7efe5e1280b88cc29f655ddeb4fe7f6950a803 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 19:11:08 -0400 Subject: [PATCH 10/46] Update benchmarks/insert_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/insert_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index 146910d0..83c86ef9 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -22,7 +22,7 @@ simulate:insert_muon: analyze:insert_muon: stage: analyze extends: .det_benchmark - needs: ["simulate:insert_muon"] + needs: ["sim:insert_muon"] script: - mkdir -p results/epic_craterlake - python benchmarks/insert_muon/analysis/muon_plots.py results/epic_craterlake/insert_muon From 8e76e1d261f4657c777234fd51f5c760701e6d83 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 19:11:16 -0400 Subject: [PATCH 11/46] Update benchmarks/insert_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/insert_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index 83c86ef9..99046b64 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -1,4 +1,4 @@ -simulate:insert_muon: +sim:insert_muon: stage: simulate extends: .det_benchmark needs: ["common:detector"] From db0fe2ac44c12fa6ca1655cdc4d53d9d6a762ce4 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 23:09:43 -0400 Subject: [PATCH 12/46] Update benchmarks/zdc_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/zdc_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml index 6021cdc1..d5bbd226 100644 --- a/benchmarks/zdc_muon/config.yml +++ b/benchmarks/zdc_muon/config.yml @@ -20,7 +20,7 @@ sim:zdc_muon: - runner_system_failure analyze:zdc_muon: - stage: analyze + stage: benchmarks extends: .det_benchmark needs: ["sim:zdc_muon"] script: From 262ae4aff351a78a0977ec9e16d3087a5532c1d6 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 23:09:55 -0400 Subject: [PATCH 13/46] Update benchmarks/zdc_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/zdc_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml index d5bbd226..3cb43519 100644 --- a/benchmarks/zdc_muon/config.yml +++ b/benchmarks/zdc_muon/config.yml @@ -28,7 +28,7 @@ analyze:zdc_muon: - python benchmarks/zdc_muon/analysis/muon_plots.py results/epic_craterlake/zdc_muon collect_results:zdc_muon: - stage: collect_results + stage: collect needs: ["analyze:zdc_muon"] script: - ls -al From 4d7f4058e0ee4bbd7711a61cc991310a6877468e Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 23:10:02 -0400 Subject: [PATCH 14/46] Update benchmarks/insert_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/insert_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index 99046b64..9f8a8c46 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -28,7 +28,7 @@ analyze:insert_muon: - python benchmarks/insert_muon/analysis/muon_plots.py results/epic_craterlake/insert_muon collect_results:insert_muon: - stage: collect_results + stage: collect needs: ["analyze:insert_muon"] script: - ls -al From 85918a787298411f79c147d6d6a55fc4371680ae Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 23:10:10 -0400 Subject: [PATCH 15/46] Update benchmarks/insert_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/insert_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index 9f8a8c46..a4ff06ff 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -20,7 +20,7 @@ sim:insert_muon: - runner_system_failure analyze:insert_muon: - stage: analyze + stage: benchmarks extends: .det_benchmark needs: ["sim:insert_muon"] script: From 893b45666f2370bb770307645b7b2b9fa79389f0 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Thu, 19 Sep 2024 23:10:17 -0400 Subject: [PATCH 16/46] Update .gitlab-ci.yml Co-authored-by: Dmitry Kalinkin --- .gitlab-ci.yml | 1 - 1 file changed, 1 deletion(-) diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index 95179a3b..82489844 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -44,7 +44,6 @@ stages: - simulate - calibrate - benchmarks - - analyze - collect - deploy - status-report From 7c4129d8f025d28ead908482316071db1e0eacd4 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 08:57:06 -0400 Subject: [PATCH 17/46] Update benchmarks/insert_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/insert_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index a4ff06ff..c26fe62f 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -29,6 +29,6 @@ analyze:insert_muon: collect_results:insert_muon: stage: collect - needs: ["analyze:insert_muon"] + needs: ["bench:insert_muon"] script: - ls -al From 2ef75907fb3de64b4c285b8e8c940bc77b073c70 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 08:57:13 -0400 Subject: [PATCH 18/46] Update benchmarks/zdc_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/zdc_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml index 3cb43519..e56282d4 100644 --- a/benchmarks/zdc_muon/config.yml +++ b/benchmarks/zdc_muon/config.yml @@ -29,6 +29,6 @@ analyze:zdc_muon: collect_results:zdc_muon: stage: collect - needs: ["analyze:zdc_muon"] + needs: ["bench:zdc_muon"] script: - ls -al From 8045593f65e613dc50ca88b333aa5753858a63b1 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 08:57:20 -0400 Subject: [PATCH 19/46] Update benchmarks/zdc_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/zdc_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml index e56282d4..c9e391c6 100644 --- a/benchmarks/zdc_muon/config.yml +++ b/benchmarks/zdc_muon/config.yml @@ -19,7 +19,7 @@ sim:zdc_muon: when: - runner_system_failure -analyze:zdc_muon: +bench:zdc_muon: stage: benchmarks extends: .det_benchmark needs: ["sim:zdc_muon"] From fd3c784ee153424f81f7569b5d1f4f99b2f61690 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 08:57:33 -0400 Subject: [PATCH 20/46] Update benchmarks/insert_muon/config.yml Co-authored-by: Dmitry Kalinkin --- benchmarks/insert_muon/config.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index c26fe62f..945a0e4a 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -19,7 +19,7 @@ sim:insert_muon: when: - runner_system_failure -analyze:insert_muon: +bench:insert_muon: stage: benchmarks extends: .det_benchmark needs: ["sim:insert_muon"] From bd612d992ade6fa311a6addc479321a8a361c1e7 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 12:08:00 -0400 Subject: [PATCH 21/46] Update muon_plots.py --- benchmarks/insert_muon/analysis/muon_plots.py | 71 +++++++++++++++++++ 1 file changed, 71 insertions(+) diff --git a/benchmarks/insert_muon/analysis/muon_plots.py b/benchmarks/insert_muon/analysis/muon_plots.py index e69de29b..8fad7628 100644 --- a/benchmarks/insert_muon/analysis/muon_plots.py +++ b/benchmarks/insert_muon/analysis/muon_plots.py @@ -0,0 +1,71 @@ +import numpy as np, pandas as pd, matplotlib.pyplot as plt, matplotlib as mpl, awkward as ak, sys +import mplhep as hep +hep.style.use("CMS") + +plt.rcParams['figure.facecolor']='white' +plt.rcParams['savefig.facecolor']='white' +plt.rcParams['savefig.bbox']='tight' + +plt.rcParams["figure.figsize"] = (7, 7) + +config=sys.argv[1].split("/")[1] #results/{config}/{benchmark_name} +outdir=sys.argv[1]+"/" +try: + import os + os.mkdir(outdir[:-1]) +except: + pass + +def Landau(x, normalization,location,stdev): + #print(type(x)) + u=(x-location)*3.591/stdev/2.355 + renormalization = 1.64872*normalization + return renormalization * np.exp(-u/2 - np.exp(-u)/2) + +import uproot as ur +arrays_sim={} +momenta=20, 30, 40, 50, 60,70, 80 +for p in momenta: + filename=f'sim_output/insert_muon/{config}_sim_mu-_{p}GeV.edm4hep.root' + print("opening file", filename) + events = ur.open(filename+':events') + arrays_sim[p] = events.arrays() + import gc + gc.collect() + print("read", filename) + +for p in 50,: + E=arrays_sim[p]["HcalEndcapPInsertHits.energy"] + y, x,_=plt.hist(1e3*ak.flatten(E),bins=100, range=(0, 1.2), histtype='step') + bc=(x[1:]+x[:-1])/2 + from scipy.optimize import curve_fit + slc=abs(bc-.48)<.15 + fnc=Landau + p0=[100, .5, .05] + #print(list(y), list(x)) + coeff, var_matrix = curve_fit(fnc, list(bc[slc]), list(y[slc]), p0=p0, + sigma=list(np.sqrt(y[slc])+(y[slc]==0))) + print(coeff) + xx=np.linspace(0,.7, 100) + plt.plot(xx, fnc(xx,*coeff), label=f'Landau fit:\nMIP={coeff[1]*1e3:.0f}$\\pm${1e3*np.sqrt(var_matrix[1][1]):.0f} keV') + plt.xlabel("hit energy [MeV]") + plt.ylabel("hits") + plt.title(f"$E_{{\\mu^-}}=${p} GeV") + plt.legend(fontsize=20) + plt.savefig(outdir+"/MIP.pdf") + +for p in 50,: + array=arrays_sim[p] + bins=30 + selection=np.sum(array["HcalEndcapPInsertHits.energy"],axis=-1)>0 + h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['eta_truth']), bins=bins) + h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['eta_truth']), bins=bins) + + h = h1 / h2 + pc=plt.pcolor(xedges, yedges, h.T) + plt.xlabel("$\\phi^*$ [rad]") + plt.ylabel("$\\eta^*$ [rad]") + cb = plt.colorbar(pc) + cb.set_label("acceptance") + plt.title(f"$E_{{\\mu^-}}=${p} GeV") + plt.savefig(outdir+"/acceptance.pdf") From 882ad876439f14689ff86448b0fc79e3405e334a Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 13:32:41 -0400 Subject: [PATCH 22/46] use only 50 GeV muons (no need for momentum dependence) --- benchmarks/insert_muon/Snakefile | 12 ++++++------ benchmarks/insert_muon/analysis/muon_plots.py | 2 +- benchmarks/insert_muon/config.yml | 6 ------ benchmarks/zdc_muon/Snakefile | 6 +++--- benchmarks/zdc_muon/config.yml | 10 ++-------- 5 files changed, 12 insertions(+), 24 deletions(-) diff --git a/benchmarks/insert_muon/Snakefile b/benchmarks/insert_muon/Snakefile index cee43378..66b7850f 100644 --- a/benchmarks/insert_muon/Snakefile +++ b/benchmarks/insert_muon/Snakefile @@ -2,9 +2,9 @@ rule insert_muon_generate: input: script="benchmarks/insert_muon/analysis/gen_particles.cxx", params: - NEVENTS_GEN=1000, - th_max=5.7, - th_min=2.0 + NEVENTS_GEN=5000, + th_max=7.0, + th_min=1.7 output: GEN_FILE="sim_output/insert_muon/mu-_{P}GeV.hepmc" shell: @@ -21,7 +21,7 @@ rule insert_muon_simulate: SIM_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root" shell: """ -NEVENTS_SIM=1000 +NEVENTS_SIM=5000 # Running simulation npsim \ --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ @@ -38,14 +38,14 @@ rule insert_muon_recon: REC_FILE="sim_output/insert_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root" shell: """ -NEVENTS_REC=1000 +NEVENTS_REC=5000 eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters -Pjana:nevents=$NEVENTS_REC """ rule insert_muon_analysis: input: expand("sim_output/insert_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root", - P=[20, 30, 40, 50, 60, 70, 80], + P=[50], DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), script="benchmarks/insert_muon/analysis/muon_plots.py", output: diff --git a/benchmarks/insert_muon/analysis/muon_plots.py b/benchmarks/insert_muon/analysis/muon_plots.py index 8fad7628..c63b56a9 100644 --- a/benchmarks/insert_muon/analysis/muon_plots.py +++ b/benchmarks/insert_muon/analysis/muon_plots.py @@ -24,7 +24,7 @@ def Landau(x, normalization,location,stdev): import uproot as ur arrays_sim={} -momenta=20, 30, 40, 50, 60,70, 80 +momenta=50, for p in momenta: filename=f'sim_output/insert_muon/{config}_sim_mu-_{p}GeV.edm4hep.root' print("opening file", filename) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index 945a0e4a..e2edfdbf 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -4,13 +4,7 @@ sim:insert_muon: needs: ["common:detector"] parallel: matrix: - - P: 20 - - P: 30 - - P: 40 - P: 50 - - P: 60 - - P: 70 - - P: 80 timeout: 1 hours script: - snakemake --cores 1 sim_output/insert_muon/epic_craterlake_rec_mu-_${P}GeV.edm4hep.root diff --git a/benchmarks/zdc_muon/Snakefile b/benchmarks/zdc_muon/Snakefile index 73fc4c3a..e38677c0 100644 --- a/benchmarks/zdc_muon/Snakefile +++ b/benchmarks/zdc_muon/Snakefile @@ -3,8 +3,8 @@ rule zdc_muon_generate: script="benchmarks/zdc_muon/analysis/gen_particles.cxx", params: NEVENTS_GEN=1000, - th_max=5.7, - th_min=2.0 + th_max=0.23, + th_min=0 output: GEN_FILE="sim_output/zdc_muon/mu-_{P}GeV.hepmc" shell: @@ -45,7 +45,7 @@ eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_fil rule zdc_muon_analysis: input: expand("sim_output/zdc_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root", - P=[20, 30, 40, 50, 60, 70, 80], + P=[50], DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), script="benchmarks/zdc_muon/analysis/muon_plots.py", output: diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml index c9e391c6..d79338dc 100644 --- a/benchmarks/zdc_muon/config.yml +++ b/benchmarks/zdc_muon/config.yml @@ -4,16 +4,10 @@ sim:zdc_muon: needs: ["common:detector"] parallel: matrix: - - P: 20 - - P: 30 - - P: 40 - P: 50 - - P: 60 - - P: 70 - - P: 80 timeout: 1 hours script: - - snakemake --cores 1 sim_output/zdc_muon/epic_craterlake_rec_mu-_${P}GeV.edm4hep.root + - snakemake --cores 1 sim_output/zdc_muon/epic_zdc_sipm_on_tile_only_rec_mu-_${P}GeV.edm4hep.root retry: max: 2 when: @@ -25,7 +19,7 @@ bench:zdc_muon: needs: ["sim:zdc_muon"] script: - mkdir -p results/epic_craterlake - - python benchmarks/zdc_muon/analysis/muon_plots.py results/epic_craterlake/zdc_muon + - python benchmarks/zdc_muon/analysis/muon_plots.py results/epic_zdc_sipm_on_tile_only/zdc_muon collect_results:zdc_muon: stage: collect From f79b3237f1e500093c376001c0d6b5b209c193a3 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 14:25:35 -0400 Subject: [PATCH 23/46] updated insert muon plots --- benchmarks/insert_muon/analysis/muon_plots.py | 27 ++++++++++++++----- 1 file changed, 21 insertions(+), 6 deletions(-) diff --git a/benchmarks/insert_muon/analysis/muon_plots.py b/benchmarks/insert_muon/analysis/muon_plots.py index c63b56a9..e1096136 100644 --- a/benchmarks/insert_muon/analysis/muon_plots.py +++ b/benchmarks/insert_muon/analysis/muon_plots.py @@ -34,6 +34,20 @@ def Landau(x, normalization,location,stdev): gc.collect() print("read", filename) +for array in arrays_sim.values(): + tilt=-0.025 + px=array['MCParticles.momentum.x'][:,2] + py=array['MCParticles.momentum.y'][:,2] + pz=array['MCParticles.momentum.z'][:,2] + p=np.sqrt(px**2+py**2+pz**2) + + pxp=px*np.cos(tilt)-pz*np.sin(tilt) + pyp=py + pzp=pz*np.cos(tilt)+px*np.sin(tilt) + + array['eta_truth']=1/2*np.log((p+pzp)/(p-pzp)) + array['phi_truth']=np.arctan2(pyp,pxp) + for p in 50,: E=arrays_sim[p]["HcalEndcapPInsertHits.energy"] y, x,_=plt.hist(1e3*ak.flatten(E),bins=100, range=(0, 1.2), histtype='step') @@ -47,6 +61,7 @@ def Landau(x, normalization,location,stdev): sigma=list(np.sqrt(y[slc])+(y[slc]==0))) print(coeff) xx=np.linspace(0,.7, 100) + MIP=coeff[1]/1000 plt.plot(xx, fnc(xx,*coeff), label=f'Landau fit:\nMIP={coeff[1]*1e3:.0f}$\\pm${1e3*np.sqrt(var_matrix[1][1]):.0f} keV') plt.xlabel("hit energy [MeV]") plt.ylabel("hits") @@ -54,15 +69,15 @@ def Landau(x, normalization,location,stdev): plt.legend(fontsize=20) plt.savefig(outdir+"/MIP.pdf") -for p in 50,: + plt.figure(figsize=(10,7)) array=arrays_sim[p] - bins=30 - selection=np.sum(array["HcalEndcapPInsertHits.energy"],axis=-1)>0 - h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['eta_truth']), bins=bins) - h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['eta_truth']), bins=bins) + bins=30; r=((-np.pi, np.pi),(2.8, 4.2)) + selection=np.sum(array["HcalEndcapPInsertHits.energy"]>0.5*MIP,axis=-1)>0 + h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['eta_truth']), bins=bins, range=r) + h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['eta_truth']), bins=bins, range=r) h = h1 / h2 - pc=plt.pcolor(xedges, yedges, h.T) + pc=plt.pcolor(xedges, yedges, h.T,linewidth=0) plt.xlabel("$\\phi^*$ [rad]") plt.ylabel("$\\eta^*$ [rad]") cb = plt.colorbar(pc) From 656d99a39f43b8b9477aae0dc0afd1d7b289a683 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 14:31:00 -0400 Subject: [PATCH 24/46] new muon plots for the zdc muon benchmark --- benchmarks/zdc_muon/analysis/muon_plots.py | 86 ++++++++++++++++++++++ 1 file changed, 86 insertions(+) diff --git a/benchmarks/zdc_muon/analysis/muon_plots.py b/benchmarks/zdc_muon/analysis/muon_plots.py index e69de29b..24a340fb 100644 --- a/benchmarks/zdc_muon/analysis/muon_plots.py +++ b/benchmarks/zdc_muon/analysis/muon_plots.py @@ -0,0 +1,86 @@ +import numpy as np, pandas as pd, matplotlib.pyplot as plt, matplotlib as mpl, awkward as ak, sys +import mplhep as hep +hep.style.use("CMS") + +plt.rcParams['figure.facecolor']='white' +plt.rcParams['savefig.facecolor']='white' +plt.rcParams['savefig.bbox']='tight' + +plt.rcParams["figure.figsize"] = (7, 7) + +config=sys.argv[1].split("/")[1] #results/{config}/{benchmark_name} +outdir=sys.argv[1]+"/" +try: + import os + os.mkdir(outdir[:-1]) +except: + pass + +def Landau(x, normalization,location,stdev): + #print(type(x)) + u=(x-location)*3.591/stdev/2.355 + renormalization = 1.64872*normalization + return renormalization * np.exp(-u/2 - np.exp(-u)/2) + +import uproot as ur +arrays_sim={} +momenta=50, +for p in momenta: + filename=f'sim_output/zdc_muon/{config}_sim_mu-_{p}GeV.edm4hep.root' + print("opening file", filename) + events = ur.open(filename+':events') + arrays_sim[p] = events.arrays() + import gc + gc.collect() + print("read", filename) + +for array in arrays_sim.values(): + tilt=-0.025 + px=array['MCParticles.momentum.x'][:,2] + py=array['MCParticles.momentum.y'][:,2] + pz=array['MCParticles.momentum.z'][:,2] + p=np.sqrt(px**2+py**2+pz**2) + + pxp=px*np.cos(tilt)-pz*np.sin(tilt) + pyp=py + pzp=pz*np.cos(tilt)+px*np.sin(tilt) + + array['eta_truth']=1/2*np.log((p+pzp)/(p-pzp)) + array['phi_truth']=np.arctan2(pyp,pxp) + +for p in 50,: + E=arrays_sim[p]["HcalFarForwardZDCHits.energy"] + y, x,_=plt.hist(1e3*ak.flatten(E),bins=100, range=(0, 1.2), histtype='step') + bc=(x[1:]+x[:-1])/2 + from scipy.optimize import curve_fit + slc=abs(bc-.48)<.15 + fnc=Landau + p0=[100, .5, .05] + #print(list(y), list(x)) + coeff, var_matrix = curve_fit(fnc, list(bc[slc]), list(y[slc]), p0=p0, + sigma=list(np.sqrt(y[slc])+(y[slc]==0))) + print(coeff) + xx=np.linspace(0,.7, 100) + MIP=coeff[1]/1000 + plt.plot(xx, fnc(xx,*coeff), label=f'Landau fit:\nMIP={coeff[1]*1e3:.0f}$\\pm${1e3*np.sqrt(var_matrix[1][1]):.0f} keV') + plt.xlabel("hit energy [MeV]") + plt.ylabel("hits") + plt.title(f"$E_{{\\mu^-}}=${p} GeV") + plt.legend(fontsize=20) + plt.savefig(outdir+"/MIP.pdf") + + plt.figure(figsize=(10,7)) + array=arrays_sim[p] + bins=30; r=((-np.pi, np.pi),(2.8, 4.2)) + selection=np.sum(array["HcalFarForwardZDCHits.energy"]>0.5*MIP,axis=-1)>0 + h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['eta_truth']), bins=bins, range=r) + h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['eta_truth']), bins=bins, range=r) + + h = h1 / h2 + pc=plt.pcolor(xedges, yedges, h.T,linewidth=0) + plt.xlabel("$\\phi^*$ [rad]") + plt.ylabel("$\\eta^*$ [rad]") + cb = plt.colorbar(pc) + cb.set_label("acceptance") + plt.title(f"$E_{{\\mu^-}}=${p} GeV") + plt.savefig(outdir+"/acceptance.pdf") From 4df9c2a34ee223ff42e3da3b04321ca12ef02c08 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 16:02:06 -0400 Subject: [PATCH 25/46] more edits --- benchmarks/insert_muon/Snakefile | 2 +- benchmarks/insert_muon/analysis/muon_plots.py | 2 +- benchmarks/zdc_muon/Snakefile | 10 +++++----- benchmarks/zdc_muon/analysis/muon_plots.py | 10 +++++----- 4 files changed, 12 insertions(+), 12 deletions(-) diff --git a/benchmarks/insert_muon/Snakefile b/benchmarks/insert_muon/Snakefile index 66b7850f..b36026b0 100644 --- a/benchmarks/insert_muon/Snakefile +++ b/benchmarks/insert_muon/Snakefile @@ -44,7 +44,7 @@ eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_fil rule insert_muon_analysis: input: - expand("sim_output/insert_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root", + expand("sim_output/insert_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root", P=[50], DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), script="benchmarks/insert_muon/analysis/muon_plots.py", diff --git a/benchmarks/insert_muon/analysis/muon_plots.py b/benchmarks/insert_muon/analysis/muon_plots.py index e1096136..4c45aac3 100644 --- a/benchmarks/insert_muon/analysis/muon_plots.py +++ b/benchmarks/insert_muon/analysis/muon_plots.py @@ -79,7 +79,7 @@ def Landau(x, normalization,location,stdev): h = h1 / h2 pc=plt.pcolor(xedges, yedges, h.T,linewidth=0) plt.xlabel("$\\phi^*$ [rad]") - plt.ylabel("$\\eta^*$ [rad]") + plt.ylabel("$\\eta^*$") cb = plt.colorbar(pc) cb.set_label("acceptance") plt.title(f"$E_{{\\mu^-}}=${p} GeV") diff --git a/benchmarks/zdc_muon/Snakefile b/benchmarks/zdc_muon/Snakefile index e38677c0..b115df42 100644 --- a/benchmarks/zdc_muon/Snakefile +++ b/benchmarks/zdc_muon/Snakefile @@ -2,8 +2,8 @@ rule zdc_muon_generate: input: script="benchmarks/zdc_muon/analysis/gen_particles.cxx", params: - NEVENTS_GEN=1000, - th_max=0.23, + NEVENTS_GEN=5000, + th_max=0.6, th_min=0 output: GEN_FILE="sim_output/zdc_muon/mu-_{P}GeV.hepmc" @@ -21,7 +21,7 @@ rule zdc_muon_simulate: SIM_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root" shell: """ -NEVENTS_SIM=1000 +NEVENTS_SIM=5000 # Running simulation npsim \ --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ @@ -38,13 +38,13 @@ rule zdc_muon_recon: REC_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root" shell: """ -NEVENTS_REC=1000 +NEVENTS_REC=5000 eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalFarForwardZDCRecHits,HcalFarForwardZDCClusters,HcalFarForwardZDCSubcellHits, -Pjana:nevents=$NEVENTS_REC """ rule zdc_muon_analysis: input: - expand("sim_output/zdc_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root", + expand("sim_output/zdc_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root", P=[50], DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), script="benchmarks/zdc_muon/analysis/muon_plots.py", diff --git a/benchmarks/zdc_muon/analysis/muon_plots.py b/benchmarks/zdc_muon/analysis/muon_plots.py index 24a340fb..ced74d00 100644 --- a/benchmarks/zdc_muon/analysis/muon_plots.py +++ b/benchmarks/zdc_muon/analysis/muon_plots.py @@ -45,7 +45,7 @@ def Landau(x, normalization,location,stdev): pyp=py pzp=pz*np.cos(tilt)+px*np.sin(tilt) - array['eta_truth']=1/2*np.log((p+pzp)/(p-pzp)) + array['theta_truth']=np.arctan(np.hypot(pxp,pyp)/pzp) array['phi_truth']=np.arctan2(pyp,pxp) for p in 50,: @@ -71,15 +71,15 @@ def Landau(x, normalization,location,stdev): plt.figure(figsize=(10,7)) array=arrays_sim[p] - bins=30; r=((-np.pi, np.pi),(2.8, 4.2)) + bins=30; r=((-np.pi, np.pi),(0, 10)) selection=np.sum(array["HcalFarForwardZDCHits.energy"]>0.5*MIP,axis=-1)>0 - h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['eta_truth']), bins=bins, range=r) - h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['eta_truth']), bins=bins, range=r) + h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['theta_truth']*1000), bins=bins, range=r) + h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['theta_truth']*1000), bins=bins, range=r) h = h1 / h2 pc=plt.pcolor(xedges, yedges, h.T,linewidth=0) plt.xlabel("$\\phi^*$ [rad]") - plt.ylabel("$\\eta^*$ [rad]") + plt.ylabel("$\\theta^*$ [mrad]") cb = plt.colorbar(pc) cb.set_label("acceptance") plt.title(f"$E_{{\\mu^-}}=${p} GeV") From 9f94b244c9384dda542419e07913ebac42abc2b4 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 16:07:47 -0400 Subject: [PATCH 26/46] don't run the rec in muon benchmarks --- benchmarks/insert_muon/config.yml | 2 +- benchmarks/zdc_muon/config.yml | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index e2edfdbf..3894e352 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -7,7 +7,7 @@ sim:insert_muon: - P: 50 timeout: 1 hours script: - - snakemake --cores 1 sim_output/insert_muon/epic_craterlake_rec_mu-_${P}GeV.edm4hep.root + - snakemake --cores 1 sim_output/insert_muon/epic_craterlake_sim_mu-_${P}GeV.edm4hep.root retry: max: 2 when: diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml index d79338dc..8e639b8c 100644 --- a/benchmarks/zdc_muon/config.yml +++ b/benchmarks/zdc_muon/config.yml @@ -7,7 +7,7 @@ sim:zdc_muon: - P: 50 timeout: 1 hours script: - - snakemake --cores 1 sim_output/zdc_muon/epic_zdc_sipm_on_tile_only_rec_mu-_${P}GeV.edm4hep.root + - snakemake --cores 1 sim_output/zdc_muon/epic_zdc_sipm_on_tile_only_sim_mu-_${P}GeV.edm4hep.root retry: max: 2 when: From c4f6055784b2c7fc28b6eaed7467d7b8d8e31bf4 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 20 Sep 2024 16:11:40 -0400 Subject: [PATCH 27/46] removed zdc muon benchmarks (may be added as a separate pull request) --- .gitlab-ci.yml | 2 - Snakefile | 2 +- benchmarks/zdc_muon/Snakefile | 57 -------- .../zdc_muon/analysis/gen_particles.cxx | 127 ------------------ benchmarks/zdc_muon/analysis/muon_plots.py | 86 ------------ benchmarks/zdc_muon/config.yml | 28 ---- 6 files changed, 1 insertion(+), 301 deletions(-) delete mode 100644 benchmarks/zdc_muon/Snakefile delete mode 100644 benchmarks/zdc_muon/analysis/gen_particles.cxx delete mode 100644 benchmarks/zdc_muon/analysis/muon_plots.py delete mode 100644 benchmarks/zdc_muon/config.yml diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index 813ba8b0..6c616e7a 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -146,7 +146,6 @@ include: - local: 'benchmarks/timing/config.yml' - local: 'benchmarks/b0_tracker/config.yml' - local: 'benchmarks/insert_muon/config.yml' - - local: 'benchmarks/zdc_muon/config.yml' deploy_results: allow_failure: true stage: deploy @@ -164,7 +163,6 @@ deploy_results: - "collect_results:zdc" - "collect_results:zdc_lyso" - "collect_results:insert_muon" - - "collect_results:zdc_muon" script: - echo "deploy results!" - find results -print | sort | tee summary.txt diff --git a/Snakefile b/Snakefile index b1d5e272..9a25a3f1 100644 --- a/Snakefile +++ b/Snakefile @@ -9,7 +9,7 @@ include: "benchmarks/tracking_performances/Snakefile" include: "benchmarks/tracking_performances_dis/Snakefile" include: "benchmarks/zdc_lyso/Snakefile" include: "benchmarks/insert_muon/Snakefile" -include: "benchmarks/zdc_muon/Snakefile" + use_s3 = config["remote_provider"].lower() == "s3" use_xrootd = config["remote_provider"].lower() == "xrootd" diff --git a/benchmarks/zdc_muon/Snakefile b/benchmarks/zdc_muon/Snakefile deleted file mode 100644 index b115df42..00000000 --- a/benchmarks/zdc_muon/Snakefile +++ /dev/null @@ -1,57 +0,0 @@ -rule zdc_muon_generate: - input: - script="benchmarks/zdc_muon/analysis/gen_particles.cxx", - params: - NEVENTS_GEN=5000, - th_max=0.6, - th_min=0 - output: - GEN_FILE="sim_output/zdc_muon/mu-_{P}GeV.hepmc" - shell: - """ -root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "mu-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' -""" - -rule zdc_muon_simulate: - input: - GEN_FILE="sim_output/zdc_muon/mu-_{P}GeV.hepmc" - params: - PHYSICS_LIST="FTFP_BERT" - output: - SIM_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root" - shell: - """ -NEVENTS_SIM=5000 -# Running simulation -npsim \ - --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ - --numberOfEvents $NEVENTS_SIM \ - --physicsList {params.PHYSICS_LIST} \ - --inputFiles {input.GEN_FILE} \ - --outputFile {output.SIM_FILE} -""" - -rule zdc_muon_recon: - input: - SIM_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root" - output: - REC_FILE="sim_output/zdc_muon/{DETECTOR_CONFIG}_rec_mu-_{P}GeV.edm4hep.root" - shell: - """ -NEVENTS_REC=5000 -eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalFarForwardZDCRecHits,HcalFarForwardZDCClusters,HcalFarForwardZDCSubcellHits, -Pjana:nevents=$NEVENTS_REC -""" - -rule zdc_muon_analysis: - input: - expand("sim_output/zdc_muon/{DETECTOR_CONFIG}_sim_mu-_{P}GeV.edm4hep.root", - P=[50], - DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), - script="benchmarks/zdc_muon/analysis/muon_plots.py", - output: - results_dir=directory("results/{DETECTOR_CONFIG}/zdc_muon"), - shell: - """ -mkdir -p {output.results_dir} -python {input.script} {output.results_dir} -""" diff --git a/benchmarks/zdc_muon/analysis/gen_particles.cxx b/benchmarks/zdc_muon/analysis/gen_particles.cxx deleted file mode 100644 index 0877a852..00000000 --- a/benchmarks/zdc_muon/analysis/gen_particles.cxx +++ /dev/null @@ -1,127 +0,0 @@ -#include "HepMC3/GenEvent.h" -#include "HepMC3/ReaderAscii.h" -#include "HepMC3/WriterAscii.h" -#include "HepMC3/Print.h" - -#include "TRandom3.h" -#include "TVector3.h" - -#include -#include -#include -#include -#include -#include -#include - -using namespace HepMC3; - -// Generate single electron as input to the Insert simulation. -// -- -// We generate events with a constant polar angle with respect to -// the proton direction and then rotate so that the events are given -// in normal lab coordinate system -// -- -void gen_particles( - int n_events = 1000, - const char* out_fname = "gen_particles.hepmc", - TString particle_name = "e-", - double th_min = 3., // Minimum polar angle, in degrees - double th_max = 3., // Maximum polar angle, in degrees - double phi_min = 0., // Minimum azimuthal angle, in degrees - double phi_max = 360., // Maximum azimuthal angle, in degrees - double p = 10., // Momentum in GeV/c - int dist = 0, //Momentum distribution: 0=fixed, 1=uniform, 2=Gaussian - int useCrossingAngle=1 // 0= no rotation, 1 = -25 mrad - ) -{ - WriterAscii hepmc_output(out_fname); - int events_parsed = 0; - GenEvent evt(Units::GEV, Units::MM); - - // Random number generator - TRandom3 *r1 = new TRandom3(0); //Use time as random seed - - // Getting generated particle information - TDatabasePDG *pdg = new TDatabasePDG(); - TParticlePDG *particle = pdg->GetParticle(particle_name); - const double mass = particle->Mass(); - const int pdgID = particle->PdgCode(); - - for (events_parsed = 0; events_parsed < n_events; events_parsed++) { - - //Set the event number - evt.set_event_number(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(FourVector(0.0, 0.0, 10.0, 10.0), 11, 4); - GenParticlePtr p2 = std::make_shared( - FourVector(0.0, 0.0, 0.0, 0.938), 2212, 4); - - // Define momentum with respect to proton direction - double phi = r1->Uniform(phi_min*TMath::DegToRad(),phi_max*TMath::DegToRad()); - double th = r1->Uniform(th_min*TMath::DegToRad(),th_max*TMath::DegToRad()); - - //Total momentum distribution - double pevent = -1; - if(dist==0){ //fixed - pevent = p; - } - else if(dist==1){ //Uniform: +-50% variation - pevent = p*(1. + r1->Uniform(-0.5,0.5) ); - } - else if(dist==2){ //Gaussian: Sigma = 0.1*mean - while(pevent<0) //Avoid negative values - pevent = r1->Gaus(p,0.1*p); - } - - double px = pevent * std::cos(phi) * std::sin(th); - double py = pevent * std::sin(phi) * std::sin(th); - double pz = pevent * std::cos(th); - TVector3 pvec(px,py,pz); - - //Rotate to lab coordinate system - double cross_angle = -25./1000.*useCrossingAngle; //in Rad - TVector3 pbeam_dir(sin(cross_angle),0,cos(cross_angle)); //proton beam direction - pvec.RotateY(-pbeam_dir.Theta()); // Theta is returned positive, beam in negative X - // type 1 is final state - // pdgid 11 - electron 0.510 MeV/c^2 - GenParticlePtr p3 = std::make_shared( - FourVector( - pvec.X(), pvec.Y(), pvec.Z(), - sqrt(pevent*pevent + (mass * mass))), - pdgID, 1); - - //If wanted, set non-zero vertex - double vx = 0.; - double vy = 0.; - double vz = 0.; - double vt = 0.; - - GenVertexPtr v1 = std::make_shared(); - evt.shift_position_by(FourVector(vx, vy, vz, vt)); - 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 % 100 == 0) { - std::cout << "Event: " << events_parsed << std::endl; - } - evt.clear(); - } - hepmc_output.close(); - std::cout << "Events parsed and written: " << events_parsed << std::endl; -} diff --git a/benchmarks/zdc_muon/analysis/muon_plots.py b/benchmarks/zdc_muon/analysis/muon_plots.py deleted file mode 100644 index ced74d00..00000000 --- a/benchmarks/zdc_muon/analysis/muon_plots.py +++ /dev/null @@ -1,86 +0,0 @@ -import numpy as np, pandas as pd, matplotlib.pyplot as plt, matplotlib as mpl, awkward as ak, sys -import mplhep as hep -hep.style.use("CMS") - -plt.rcParams['figure.facecolor']='white' -plt.rcParams['savefig.facecolor']='white' -plt.rcParams['savefig.bbox']='tight' - -plt.rcParams["figure.figsize"] = (7, 7) - -config=sys.argv[1].split("/")[1] #results/{config}/{benchmark_name} -outdir=sys.argv[1]+"/" -try: - import os - os.mkdir(outdir[:-1]) -except: - pass - -def Landau(x, normalization,location,stdev): - #print(type(x)) - u=(x-location)*3.591/stdev/2.355 - renormalization = 1.64872*normalization - return renormalization * np.exp(-u/2 - np.exp(-u)/2) - -import uproot as ur -arrays_sim={} -momenta=50, -for p in momenta: - filename=f'sim_output/zdc_muon/{config}_sim_mu-_{p}GeV.edm4hep.root' - print("opening file", filename) - events = ur.open(filename+':events') - arrays_sim[p] = events.arrays() - import gc - gc.collect() - print("read", filename) - -for array in arrays_sim.values(): - tilt=-0.025 - px=array['MCParticles.momentum.x'][:,2] - py=array['MCParticles.momentum.y'][:,2] - pz=array['MCParticles.momentum.z'][:,2] - p=np.sqrt(px**2+py**2+pz**2) - - pxp=px*np.cos(tilt)-pz*np.sin(tilt) - pyp=py - pzp=pz*np.cos(tilt)+px*np.sin(tilt) - - array['theta_truth']=np.arctan(np.hypot(pxp,pyp)/pzp) - array['phi_truth']=np.arctan2(pyp,pxp) - -for p in 50,: - E=arrays_sim[p]["HcalFarForwardZDCHits.energy"] - y, x,_=plt.hist(1e3*ak.flatten(E),bins=100, range=(0, 1.2), histtype='step') - bc=(x[1:]+x[:-1])/2 - from scipy.optimize import curve_fit - slc=abs(bc-.48)<.15 - fnc=Landau - p0=[100, .5, .05] - #print(list(y), list(x)) - coeff, var_matrix = curve_fit(fnc, list(bc[slc]), list(y[slc]), p0=p0, - sigma=list(np.sqrt(y[slc])+(y[slc]==0))) - print(coeff) - xx=np.linspace(0,.7, 100) - MIP=coeff[1]/1000 - plt.plot(xx, fnc(xx,*coeff), label=f'Landau fit:\nMIP={coeff[1]*1e3:.0f}$\\pm${1e3*np.sqrt(var_matrix[1][1]):.0f} keV') - plt.xlabel("hit energy [MeV]") - plt.ylabel("hits") - plt.title(f"$E_{{\\mu^-}}=${p} GeV") - plt.legend(fontsize=20) - plt.savefig(outdir+"/MIP.pdf") - - plt.figure(figsize=(10,7)) - array=arrays_sim[p] - bins=30; r=((-np.pi, np.pi),(0, 10)) - selection=np.sum(array["HcalFarForwardZDCHits.energy"]>0.5*MIP,axis=-1)>0 - h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['theta_truth']*1000), bins=bins, range=r) - h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['theta_truth']*1000), bins=bins, range=r) - - h = h1 / h2 - pc=plt.pcolor(xedges, yedges, h.T,linewidth=0) - plt.xlabel("$\\phi^*$ [rad]") - plt.ylabel("$\\theta^*$ [mrad]") - cb = plt.colorbar(pc) - cb.set_label("acceptance") - plt.title(f"$E_{{\\mu^-}}=${p} GeV") - plt.savefig(outdir+"/acceptance.pdf") diff --git a/benchmarks/zdc_muon/config.yml b/benchmarks/zdc_muon/config.yml deleted file mode 100644 index 8e639b8c..00000000 --- a/benchmarks/zdc_muon/config.yml +++ /dev/null @@ -1,28 +0,0 @@ -sim:zdc_muon: - stage: simulate - extends: .det_benchmark - needs: ["common:detector"] - parallel: - matrix: - - P: 50 - timeout: 1 hours - script: - - snakemake --cores 1 sim_output/zdc_muon/epic_zdc_sipm_on_tile_only_sim_mu-_${P}GeV.edm4hep.root - retry: - max: 2 - when: - - runner_system_failure - -bench:zdc_muon: - stage: benchmarks - extends: .det_benchmark - needs: ["sim:zdc_muon"] - script: - - mkdir -p results/epic_craterlake - - python benchmarks/zdc_muon/analysis/muon_plots.py results/epic_zdc_sipm_on_tile_only/zdc_muon - -collect_results:zdc_muon: - stage: collect - needs: ["bench:zdc_muon"] - script: - - ls -al From d22c8b84a6b5fa4b60107b8a1c2817d3ca0f548d Mon Sep 17 00:00:00 2001 From: Dmitry Kalinkin Date: Mon, 23 Sep 2024 02:34:28 -0400 Subject: [PATCH 28/46] Update config.yml --- benchmarks/insert_muon/config.yml | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index 3894e352..c8bb3b85 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -18,11 +18,13 @@ bench:insert_muon: extends: .det_benchmark needs: ["sim:insert_muon"] script: - - mkdir -p results/epic_craterlake - - python benchmarks/insert_muon/analysis/muon_plots.py results/epic_craterlake/insert_muon + - snakemake --cores 1 results/epic_zdc_sipm_on_tile_only/insert_muon collect_results:insert_muon: stage: collect needs: ["bench:insert_muon"] script: - ls -al + - mv results{,_save}/ # move results directory out of the way to preserve it + - snakemake --cores 1 --delete-all-output results/epic_zdc_sipm_on_tile_only/insert_muon + - mv results{_save,}/ From e43eb9466eccd23b6a84b160546345fe88410272 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Mon, 23 Sep 2024 12:51:57 -0400 Subject: [PATCH 29/46] Update config.yml use craterlake configuration, not zdc_sipm_on_tile_only --- benchmarks/insert_muon/config.yml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/benchmarks/insert_muon/config.yml b/benchmarks/insert_muon/config.yml index c8bb3b85..bf69f555 100644 --- a/benchmarks/insert_muon/config.yml +++ b/benchmarks/insert_muon/config.yml @@ -18,7 +18,7 @@ bench:insert_muon: extends: .det_benchmark needs: ["sim:insert_muon"] script: - - snakemake --cores 1 results/epic_zdc_sipm_on_tile_only/insert_muon + - snakemake --cores 1 results/epic_craterlake/insert_muon collect_results:insert_muon: stage: collect @@ -26,5 +26,5 @@ collect_results:insert_muon: script: - ls -al - mv results{,_save}/ # move results directory out of the way to preserve it - - snakemake --cores 1 --delete-all-output results/epic_zdc_sipm_on_tile_only/insert_muon + - snakemake --cores 1 --delete-all-output results/epic_craterlake/insert_muon - mv results{_save,}/ From f08596ff51fa811164a3de7ad7604b2432e7608f Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Wed, 25 Sep 2024 09:17:54 -0400 Subject: [PATCH 30/46] first commit for branch --- .gitlab-ci.yml | 2 + Snakefile | 3 +- benchmarks/insert_tau/Snakefile | 57 ++++++++ benchmarks/insert_tau/Snakefile~ | 57 ++++++++ .../insert_tau/analysis/gen_particles.cxx | 127 ++++++++++++++++++ benchmarks/insert_tau/analysis/muon_plots.py | 86 ++++++++++++ benchmarks/insert_tau/analysis/tau_plots.py | 0 benchmarks/insert_tau/config.yml | 30 +++++ 8 files changed, 360 insertions(+), 2 deletions(-) create mode 100644 benchmarks/insert_tau/Snakefile create mode 100644 benchmarks/insert_tau/Snakefile~ create mode 100644 benchmarks/insert_tau/analysis/gen_particles.cxx create mode 100644 benchmarks/insert_tau/analysis/muon_plots.py create mode 100644 benchmarks/insert_tau/analysis/tau_plots.py create mode 100644 benchmarks/insert_tau/config.yml diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index fb5b57e3..e59b9ed1 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -148,6 +148,7 @@ include: - local: 'benchmarks/timing/config.yml' - local: 'benchmarks/b0_tracker/config.yml' - local: 'benchmarks/insert_muon/config.yml' + - local: 'benchmarks/insert_tau/config.yml' - local: 'benchmarks/zdc_sigma/config.yml' - local: 'benchmarks/zdc_lambda/config.yml' - local: 'benchmarks/insert_neutron/config.yml' @@ -172,6 +173,7 @@ deploy_results: - "collect_results:zdc" - "collect_results:zdc_lyso" - "collect_results:insert_muon" + - "collect_results:insert_tau" - "collect_results:zdc_photon" - "collect_results:zdc_pi0" script: diff --git a/Snakefile b/Snakefile index bce37f8d..8d1de7a6 100644 --- a/Snakefile +++ b/Snakefile @@ -10,12 +10,11 @@ include: "benchmarks/tracking_performances_dis/Snakefile" include: "benchmarks/zdc_lyso/Snakefile" include: "benchmarks/insert_muon/Snakefile" include: "benchmarks/zdc_lambda/Snakefile" -include: "benchmarks/zdc_lyso/Snakefile" include: "benchmarks/zdc_photon/Snakefile" include: "benchmarks/zdc_pi0/Snakefile" include: "benchmarks/zdc_sigma/Snakefile" include: "benchmarks/insert_neutron/Snakefile" - +include: "benchmarks/insert_tau/Snakefile" use_s3 = config["remote_provider"].lower() == "s3" use_xrootd = config["remote_provider"].lower() == "xrootd" diff --git a/benchmarks/insert_tau/Snakefile b/benchmarks/insert_tau/Snakefile new file mode 100644 index 00000000..e0f9a77c --- /dev/null +++ b/benchmarks/insert_tau/Snakefile @@ -0,0 +1,57 @@ +rule insert_tau_generate: + input: + script="benchmarks/insert_tau/analysis/gen_particles.cxx", + params: + NEVENTS_GEN=1000, + th_max=7.0, + th_min=1.7 + output: + GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" + shell: + """ +root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "tau-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' +""" + +rule insert_tau_simulate: + input: + GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" + params: + PHYSICS_LIST="FTFP_BERT" + output: + SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_SIM=1000 +# Running simulation +npsim \ + --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ + --numberOfEvents $NEVENTS_SIM \ + --physicsList {params.PHYSICS_LIST} \ + --inputFiles {input.GEN_FILE} \ + --outputFile {output.SIM_FILE} +""" + +rule insert_tau_recon: + input: + SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV.edm4hep.root" + output: + REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_REC=1000 +eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters -Pjana:nevents=$NEVENTS_REC +""" + +rule insert_tau_analysis: + input: + expand("sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep.root", + P=[20,30, 40,50,60, 80, 100], + DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), + script="benchmarks/insert_tau/analysis/tau_plots.py", + output: + results_dir=directory("results/{DETECTOR_CONFIG}/insert_tau"), + shell: + """ +mkdir -p {output.results_dir} +python {input.script} {output.results_dir} +""" diff --git a/benchmarks/insert_tau/Snakefile~ b/benchmarks/insert_tau/Snakefile~ new file mode 100644 index 00000000..e7b0cb34 --- /dev/null +++ b/benchmarks/insert_tau/Snakefile~ @@ -0,0 +1,57 @@ +rule insert_tau_generate: + input: + script="benchmarks/insert_tau/analysis/gen_particles.cxx", + params: + NEVENTS_GEN=1000, + th_max=7.0, + th_min=1.7 + output: + GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" + shell: + """ +root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "tau-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' +""" + +rule insert_tau_simulate: + input: + GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" + params: + PHYSICS_LIST="FTFP_BERT" + output: + SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_SIM=1000 +# Running simulation +npsim \ + --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ + --numberOfEvents $NEVENTS_SIM \ + --physicsList {params.PHYSICS_LIST} \ + --inputFiles {input.GEN_FILE} \ + --outputFile {output.SIM_FILE} +""" + +rule insert_tau_recon: + input: + SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV.edm4hep.root" + output: + REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep.root" + shell: + """ +NEVENTS_REC=1000 +eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters -Pjana:nevents=$NEVENTS_REC +""" + +rule insert_tau_analysis: + input: + expand("sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep.root", + P=[20,50, 100, 200], + DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), + script="benchmarks/insert_tau/analysis/tau_plots.py", + output: + results_dir=directory("results/{DETECTOR_CONFIG}/insert_tau"), + shell: + """ +mkdir -p {output.results_dir} +python {input.script} {output.results_dir} +""" diff --git a/benchmarks/insert_tau/analysis/gen_particles.cxx b/benchmarks/insert_tau/analysis/gen_particles.cxx new file mode 100644 index 00000000..0877a852 --- /dev/null +++ b/benchmarks/insert_tau/analysis/gen_particles.cxx @@ -0,0 +1,127 @@ +#include "HepMC3/GenEvent.h" +#include "HepMC3/ReaderAscii.h" +#include "HepMC3/WriterAscii.h" +#include "HepMC3/Print.h" + +#include "TRandom3.h" +#include "TVector3.h" + +#include +#include +#include +#include +#include +#include +#include + +using namespace HepMC3; + +// Generate single electron as input to the Insert simulation. +// -- +// We generate events with a constant polar angle with respect to +// the proton direction and then rotate so that the events are given +// in normal lab coordinate system +// -- +void gen_particles( + int n_events = 1000, + const char* out_fname = "gen_particles.hepmc", + TString particle_name = "e-", + double th_min = 3., // Minimum polar angle, in degrees + double th_max = 3., // Maximum polar angle, in degrees + double phi_min = 0., // Minimum azimuthal angle, in degrees + double phi_max = 360., // Maximum azimuthal angle, in degrees + double p = 10., // Momentum in GeV/c + int dist = 0, //Momentum distribution: 0=fixed, 1=uniform, 2=Gaussian + int useCrossingAngle=1 // 0= no rotation, 1 = -25 mrad + ) +{ + WriterAscii hepmc_output(out_fname); + int events_parsed = 0; + GenEvent evt(Units::GEV, Units::MM); + + // Random number generator + TRandom3 *r1 = new TRandom3(0); //Use time as random seed + + // Getting generated particle information + TDatabasePDG *pdg = new TDatabasePDG(); + TParticlePDG *particle = pdg->GetParticle(particle_name); + const double mass = particle->Mass(); + const int pdgID = particle->PdgCode(); + + for (events_parsed = 0; events_parsed < n_events; events_parsed++) { + + //Set the event number + evt.set_event_number(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(FourVector(0.0, 0.0, 10.0, 10.0), 11, 4); + GenParticlePtr p2 = std::make_shared( + FourVector(0.0, 0.0, 0.0, 0.938), 2212, 4); + + // Define momentum with respect to proton direction + double phi = r1->Uniform(phi_min*TMath::DegToRad(),phi_max*TMath::DegToRad()); + double th = r1->Uniform(th_min*TMath::DegToRad(),th_max*TMath::DegToRad()); + + //Total momentum distribution + double pevent = -1; + if(dist==0){ //fixed + pevent = p; + } + else if(dist==1){ //Uniform: +-50% variation + pevent = p*(1. + r1->Uniform(-0.5,0.5) ); + } + else if(dist==2){ //Gaussian: Sigma = 0.1*mean + while(pevent<0) //Avoid negative values + pevent = r1->Gaus(p,0.1*p); + } + + double px = pevent * std::cos(phi) * std::sin(th); + double py = pevent * std::sin(phi) * std::sin(th); + double pz = pevent * std::cos(th); + TVector3 pvec(px,py,pz); + + //Rotate to lab coordinate system + double cross_angle = -25./1000.*useCrossingAngle; //in Rad + TVector3 pbeam_dir(sin(cross_angle),0,cos(cross_angle)); //proton beam direction + pvec.RotateY(-pbeam_dir.Theta()); // Theta is returned positive, beam in negative X + // type 1 is final state + // pdgid 11 - electron 0.510 MeV/c^2 + GenParticlePtr p3 = std::make_shared( + FourVector( + pvec.X(), pvec.Y(), pvec.Z(), + sqrt(pevent*pevent + (mass * mass))), + pdgID, 1); + + //If wanted, set non-zero vertex + double vx = 0.; + double vy = 0.; + double vz = 0.; + double vt = 0.; + + GenVertexPtr v1 = std::make_shared(); + evt.shift_position_by(FourVector(vx, vy, vz, vt)); + 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 % 100 == 0) { + std::cout << "Event: " << events_parsed << std::endl; + } + evt.clear(); + } + hepmc_output.close(); + std::cout << "Events parsed and written: " << events_parsed << std::endl; +} diff --git a/benchmarks/insert_tau/analysis/muon_plots.py b/benchmarks/insert_tau/analysis/muon_plots.py new file mode 100644 index 00000000..4c45aac3 --- /dev/null +++ b/benchmarks/insert_tau/analysis/muon_plots.py @@ -0,0 +1,86 @@ +import numpy as np, pandas as pd, matplotlib.pyplot as plt, matplotlib as mpl, awkward as ak, sys +import mplhep as hep +hep.style.use("CMS") + +plt.rcParams['figure.facecolor']='white' +plt.rcParams['savefig.facecolor']='white' +plt.rcParams['savefig.bbox']='tight' + +plt.rcParams["figure.figsize"] = (7, 7) + +config=sys.argv[1].split("/")[1] #results/{config}/{benchmark_name} +outdir=sys.argv[1]+"/" +try: + import os + os.mkdir(outdir[:-1]) +except: + pass + +def Landau(x, normalization,location,stdev): + #print(type(x)) + u=(x-location)*3.591/stdev/2.355 + renormalization = 1.64872*normalization + return renormalization * np.exp(-u/2 - np.exp(-u)/2) + +import uproot as ur +arrays_sim={} +momenta=50, +for p in momenta: + filename=f'sim_output/insert_muon/{config}_sim_mu-_{p}GeV.edm4hep.root' + print("opening file", filename) + events = ur.open(filename+':events') + arrays_sim[p] = events.arrays() + import gc + gc.collect() + print("read", filename) + +for array in arrays_sim.values(): + tilt=-0.025 + px=array['MCParticles.momentum.x'][:,2] + py=array['MCParticles.momentum.y'][:,2] + pz=array['MCParticles.momentum.z'][:,2] + p=np.sqrt(px**2+py**2+pz**2) + + pxp=px*np.cos(tilt)-pz*np.sin(tilt) + pyp=py + pzp=pz*np.cos(tilt)+px*np.sin(tilt) + + array['eta_truth']=1/2*np.log((p+pzp)/(p-pzp)) + array['phi_truth']=np.arctan2(pyp,pxp) + +for p in 50,: + E=arrays_sim[p]["HcalEndcapPInsertHits.energy"] + y, x,_=plt.hist(1e3*ak.flatten(E),bins=100, range=(0, 1.2), histtype='step') + bc=(x[1:]+x[:-1])/2 + from scipy.optimize import curve_fit + slc=abs(bc-.48)<.15 + fnc=Landau + p0=[100, .5, .05] + #print(list(y), list(x)) + coeff, var_matrix = curve_fit(fnc, list(bc[slc]), list(y[slc]), p0=p0, + sigma=list(np.sqrt(y[slc])+(y[slc]==0))) + print(coeff) + xx=np.linspace(0,.7, 100) + MIP=coeff[1]/1000 + plt.plot(xx, fnc(xx,*coeff), label=f'Landau fit:\nMIP={coeff[1]*1e3:.0f}$\\pm${1e3*np.sqrt(var_matrix[1][1]):.0f} keV') + plt.xlabel("hit energy [MeV]") + plt.ylabel("hits") + plt.title(f"$E_{{\\mu^-}}=${p} GeV") + plt.legend(fontsize=20) + plt.savefig(outdir+"/MIP.pdf") + + plt.figure(figsize=(10,7)) + array=arrays_sim[p] + bins=30; r=((-np.pi, np.pi),(2.8, 4.2)) + selection=np.sum(array["HcalEndcapPInsertHits.energy"]>0.5*MIP,axis=-1)>0 + h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['eta_truth']), bins=bins, range=r) + h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['eta_truth']), bins=bins, range=r) + + h = h1 / h2 + pc=plt.pcolor(xedges, yedges, h.T,linewidth=0) + plt.xlabel("$\\phi^*$ [rad]") + plt.ylabel("$\\eta^*$") + cb = plt.colorbar(pc) + cb.set_label("acceptance") + plt.title(f"$E_{{\\mu^-}}=${p} GeV") + plt.savefig(outdir+"/acceptance.pdf") diff --git a/benchmarks/insert_tau/analysis/tau_plots.py b/benchmarks/insert_tau/analysis/tau_plots.py new file mode 100644 index 00000000..e69de29b diff --git a/benchmarks/insert_tau/config.yml b/benchmarks/insert_tau/config.yml new file mode 100644 index 00000000..bf69f555 --- /dev/null +++ b/benchmarks/insert_tau/config.yml @@ -0,0 +1,30 @@ +sim:insert_muon: + stage: simulate + extends: .det_benchmark + needs: ["common:detector"] + parallel: + matrix: + - P: 50 + timeout: 1 hours + script: + - snakemake --cores 1 sim_output/insert_muon/epic_craterlake_sim_mu-_${P}GeV.edm4hep.root + retry: + max: 2 + when: + - runner_system_failure + +bench:insert_muon: + stage: benchmarks + extends: .det_benchmark + needs: ["sim:insert_muon"] + script: + - snakemake --cores 1 results/epic_craterlake/insert_muon + +collect_results:insert_muon: + stage: collect + needs: ["bench:insert_muon"] + script: + - ls -al + - mv results{,_save}/ # move results directory out of the way to preserve it + - snakemake --cores 1 --delete-all-output results/epic_craterlake/insert_muon + - mv results{_save,}/ From 693b08f06b6f7a5b8c505d8250d63e9fe5b0a94b Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Wed, 23 Oct 2024 12:50:48 -0400 Subject: [PATCH 31/46] changed which banks are included --- benchmarks/insert_tau/Snakefile | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/insert_tau/Snakefile b/benchmarks/insert_tau/Snakefile index e0f9a77c..59dfe56b 100644 --- a/benchmarks/insert_tau/Snakefile +++ b/benchmarks/insert_tau/Snakefile @@ -39,7 +39,7 @@ rule insert_tau_recon: shell: """ NEVENTS_REC=1000 -eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters -Pjana:nevents=$NEVENTS_REC +eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters,EcalEndcapPClusters,LFHCALClusters -Pjana:nevents=$NEVENTS_REC """ rule insert_tau_analysis: From 956ac1776bde8e09eb0f378939c7d892f58c4616 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Wed, 23 Oct 2024 13:14:17 -0400 Subject: [PATCH 32/46] updated tau snakefile and config.yml to do multiple processes simultaneously --- benchmarks/insert_tau/Snakefile | 14 +++++++++----- benchmarks/insert_tau/config.yml | 22 ++++++++++++++-------- 2 files changed, 23 insertions(+), 13 deletions(-) diff --git a/benchmarks/insert_tau/Snakefile b/benchmarks/insert_tau/Snakefile index 59dfe56b..d8aef969 100644 --- a/benchmarks/insert_tau/Snakefile +++ b/benchmarks/insert_tau/Snakefile @@ -15,10 +15,11 @@ root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "tau-", rule insert_tau_simulate: input: GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" + warmup="warmup/{DETECTOR_CONFIG}.edm4hep.root" params: PHYSICS_LIST="FTFP_BERT" output: - SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV.edm4hep.root" + SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV_{INDEX}.edm4hep.root" shell: """ NEVENTS_SIM=1000 @@ -26,6 +27,7 @@ NEVENTS_SIM=1000 npsim \ --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ --numberOfEvents $NEVENTS_SIM \ + --skipNEvents $(( $NEVENTS_SIM * {wildcards.INDEX} )) \ --physicsList {params.PHYSICS_LIST} \ --inputFiles {input.GEN_FILE} \ --outputFile {output.SIM_FILE} @@ -33,9 +35,9 @@ npsim \ rule insert_tau_recon: input: - SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV.edm4hep.root" + SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV_{INDEX}.edm4hep.root" output: - REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep.root" + REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV_{INDEX}.edm4hep.root" shell: """ NEVENTS_REC=1000 @@ -44,9 +46,11 @@ eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_fil rule insert_tau_analysis: input: - expand("sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep.root", + expand("sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep_{INDEX}.root", P=[20,30, 40,50,60, 80, 100], - DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), + DETECTOR_CONFIG=["{DETECTOR_CONFIG}"], + INDEX=range(5), + ), script="benchmarks/insert_tau/analysis/tau_plots.py", output: results_dir=directory("results/{DETECTOR_CONFIG}/insert_tau"), diff --git a/benchmarks/insert_tau/config.yml b/benchmarks/insert_tau/config.yml index bf69f555..14ffa34a 100644 --- a/benchmarks/insert_tau/config.yml +++ b/benchmarks/insert_tau/config.yml @@ -1,30 +1,36 @@ -sim:insert_muon: +sim:insert_tau: stage: simulate extends: .det_benchmark needs: ["common:detector"] parallel: matrix: + - P: 20 + - P: 30 + - P: 40 - P: 50 + - P: 60 + - P: 80 + - P: 100 timeout: 1 hours script: - - snakemake --cores 1 sim_output/insert_muon/epic_craterlake_sim_mu-_${P}GeV.edm4hep.root + - snakemake $SNAKEMAKE_FLAGS --cores 5 sim_output/insert_tau/epic_craterlake_sim_tau-_${P}GeV_{0,1,2,3,4}.edm4hep.root retry: max: 2 when: - runner_system_failure -bench:insert_muon: +bench:insert_tau: stage: benchmarks extends: .det_benchmark - needs: ["sim:insert_muon"] + needs: ["sim:insert_tau"] script: - - snakemake --cores 1 results/epic_craterlake/insert_muon + - snakemake $SNAKEMAKE_FLAGS --cores 1 results/epic_craterlake/insert_tau -collect_results:insert_muon: +collect_results:insert_tau: stage: collect - needs: ["bench:insert_muon"] + needs: ["bench:insert_tau"] script: - ls -al - mv results{,_save}/ # move results directory out of the way to preserve it - - snakemake --cores 1 --delete-all-output results/epic_craterlake/insert_muon + - snakemake $SNAKEMAKE_FLAGS --cores 1 --delete-all-output results/epic_craterlake/insert_tau - mv results{_save,}/ From 04bf9d2a99331d07a6b38d860c3bdd46050237a2 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Wed, 23 Oct 2024 13:36:18 -0400 Subject: [PATCH 33/46] updated snakefiles --- Snakefile | 2 -- benchmarks/insert_tau/Snakefile | 2 +- 2 files changed, 1 insertion(+), 3 deletions(-) diff --git a/Snakefile b/Snakefile index 9d7a6127..9c125483 100644 --- a/Snakefile +++ b/Snakefile @@ -8,8 +8,6 @@ include: "benchmarks/material_scan/Snakefile" include: "benchmarks/tracking_performances/Snakefile" include: "benchmarks/tracking_performances_dis/Snakefile" include: "benchmarks/lfhcal/Snakefile" -include: "benchmarks/insert_muon/Snakefile" -include: "benchmarks/zdc_lambda/Snakefile" include: "benchmarks/zdc_lyso/Snakefile" include: "benchmarks/insert_muon/Snakefile" include: "benchmarks/zdc_lambda/Snakefile" diff --git a/benchmarks/insert_tau/Snakefile b/benchmarks/insert_tau/Snakefile index d8aef969..2e6b120e 100644 --- a/benchmarks/insert_tau/Snakefile +++ b/benchmarks/insert_tau/Snakefile @@ -6,7 +6,7 @@ rule insert_tau_generate: th_max=7.0, th_min=1.7 output: - GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" + GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc", shell: """ root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "tau-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' From 9bf8ce26794f91c4a57b65dc0176221ac3af56ec Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Wed, 23 Oct 2024 17:13:26 -0400 Subject: [PATCH 34/46] Delete benchmarks/insert_tau/Snakefile~ --- benchmarks/insert_tau/Snakefile~ | 57 -------------------------------- 1 file changed, 57 deletions(-) delete mode 100644 benchmarks/insert_tau/Snakefile~ diff --git a/benchmarks/insert_tau/Snakefile~ b/benchmarks/insert_tau/Snakefile~ deleted file mode 100644 index e7b0cb34..00000000 --- a/benchmarks/insert_tau/Snakefile~ +++ /dev/null @@ -1,57 +0,0 @@ -rule insert_tau_generate: - input: - script="benchmarks/insert_tau/analysis/gen_particles.cxx", - params: - NEVENTS_GEN=1000, - th_max=7.0, - th_min=1.7 - output: - GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" - shell: - """ -root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "tau-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' -""" - -rule insert_tau_simulate: - input: - GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" - params: - PHYSICS_LIST="FTFP_BERT" - output: - SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV.edm4hep.root" - shell: - """ -NEVENTS_SIM=1000 -# Running simulation -npsim \ - --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ - --numberOfEvents $NEVENTS_SIM \ - --physicsList {params.PHYSICS_LIST} \ - --inputFiles {input.GEN_FILE} \ - --outputFile {output.SIM_FILE} -""" - -rule insert_tau_recon: - input: - SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV.edm4hep.root" - output: - REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep.root" - shell: - """ -NEVENTS_REC=1000 -eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters -Pjana:nevents=$NEVENTS_REC -""" - -rule insert_tau_analysis: - input: - expand("sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep.root", - P=[20,50, 100, 200], - DETECTOR_CONFIG=["{DETECTOR_CONFIG}"]), - script="benchmarks/insert_tau/analysis/tau_plots.py", - output: - results_dir=directory("results/{DETECTOR_CONFIG}/insert_tau"), - shell: - """ -mkdir -p {output.results_dir} -python {input.script} {output.results_dir} -""" From 1589f4ab8e7f1d3bbf49e6310a72f18519a08cc8 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Wed, 23 Oct 2024 17:15:24 -0400 Subject: [PATCH 35/46] revert muon_plots.py --- benchmarks/insert_muon/analysis/muon_plots.py | 166 +++++++++--------- 1 file changed, 83 insertions(+), 83 deletions(-) diff --git a/benchmarks/insert_muon/analysis/muon_plots.py b/benchmarks/insert_muon/analysis/muon_plots.py index 117ad6f5..ac2aa100 100644 --- a/benchmarks/insert_muon/analysis/muon_plots.py +++ b/benchmarks/insert_muon/analysis/muon_plots.py @@ -1,83 +1,83 @@ -import numpy as np, pandas as pd, matplotlib.pyplot as plt, matplotlib as mpl, awkward as ak, sys -import mplhep as hep -hep.style.use("CMS") - -plt.rcParams['figure.facecolor']='white' -plt.rcParams['savefig.facecolor']='white' -plt.rcParams['savefig.bbox']='tight' - -plt.rcParams["figure.figsize"] = (7, 7) - -config=sys.argv[1].split("/")[1] #results/{config}/{benchmark_name} -outdir=sys.argv[1]+"/" -try: - import os - os.mkdir(outdir[:-1]) -except: - pass - -def Landau(x, normalization,location,stdev): - #print(type(x)) - u=(x-location)*3.591/stdev/2.355 - renormalization = 1.64872*normalization - return renormalization * np.exp(-u/2 - np.exp(-u)/2) - -import uproot as ur -arrays_sim={} -momenta=50, -for p in momenta: - arrays_sim[p] = ur.concatenate({ - f'sim_output/insert_muon/{config}_sim_mu-_{p}GeV_{index}.edm4hep.root': 'events' - for index in range(5) - }) - -for array in arrays_sim.values(): - tilt=-0.025 - px=array['MCParticles.momentum.x'][:,2] - py=array['MCParticles.momentum.y'][:,2] - pz=array['MCParticles.momentum.z'][:,2] - p=np.sqrt(px**2+py**2+pz**2) - - pxp=px*np.cos(tilt)-pz*np.sin(tilt) - pyp=py - pzp=pz*np.cos(tilt)+px*np.sin(tilt) - - array['eta_truth']=1/2*np.log((p+pzp)/(p-pzp)) - array['phi_truth']=np.arctan2(pyp,pxp) - -for p in 50,: - E=arrays_sim[p]["HcalEndcapPInsertHits.energy"] - y, x,_=plt.hist(1e3*ak.flatten(E),bins=100, range=(0, 1.2), histtype='step') - bc=(x[1:]+x[:-1])/2 - from scipy.optimize import curve_fit - slc=abs(bc-.48)<.15 - fnc=Landau - p0=[100, .5, .05] - #print(list(y), list(x)) - coeff, var_matrix = curve_fit(fnc, list(bc[slc]), list(y[slc]), p0=p0, - sigma=list(np.sqrt(y[slc])+(y[slc]==0)), maxfev=10000) - print(coeff) - xx=np.linspace(0,.7, 100) - MIP=coeff[1]/1000 - plt.plot(xx, fnc(xx,*coeff), label=f'Landau fit:\nMIP={coeff[1]*1e3:.0f}$\\pm${1e3*np.sqrt(var_matrix[1][1]):.0f} keV') - plt.xlabel("hit energy [MeV]") - plt.ylabel("hits") - plt.title(f"$E_{{\\mu^-}}=${p} GeV") - plt.legend(fontsize=20) - plt.savefig(outdir+"/MIP.pdf") - - plt.figure(figsize=(10,7)) - array=arrays_sim[p] - bins=30; r=((-np.pi, np.pi),(2.8, 4.2)) - selection=np.sum(array["HcalEndcapPInsertHits.energy"]>0.5*MIP,axis=-1)>0 - h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['eta_truth']), bins=bins, range=r) - h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['eta_truth']), bins=bins, range=r) - - h = h1 / h2 - pc=plt.pcolor(xedges, yedges, h.T,linewidth=0) - plt.xlabel("$\\phi^*$ [rad]") - plt.ylabel("$\\eta^*$") - cb = plt.colorbar(pc) - cb.set_label("acceptance") - plt.title(f"$E_{{\\mu^-}}=${p} GeV") - plt.savefig(outdir+"/acceptance.pdf") +import numpy as np, pandas as pd, matplotlib.pyplot as plt, matplotlib as mpl, awkward as ak, sys +import mplhep as hep +hep.style.use("CMS") + +plt.rcParams['figure.facecolor']='white' +plt.rcParams['savefig.facecolor']='white' +plt.rcParams['savefig.bbox']='tight' + +plt.rcParams["figure.figsize"] = (7, 7) + +config=sys.argv[1].split("/")[1] #results/{config}/{benchmark_name} +outdir=sys.argv[1]+"/" +try: + import os + os.mkdir(outdir[:-1]) +except: + pass + +def Landau(x, normalization,location,stdev): + #print(type(x)) + u=(x-location)*3.591/stdev/2.355 + renormalization = 1.64872*normalization + return renormalization * np.exp(-u/2 - np.exp(-u)/2) + +import uproot as ur +arrays_sim={} +momenta=50, +for p in momenta: + arrays_sim[p] = ur.concatenate({ + f'sim_output/insert_muon/{config}_sim_mu-_{p}GeV_{index}.edm4hep.root': 'events' + for index in range(5) + }) + +for array in arrays_sim.values(): + tilt=-0.025 + px=array['MCParticles.momentum.x'][:,2] + py=array['MCParticles.momentum.y'][:,2] + pz=array['MCParticles.momentum.z'][:,2] + p=np.sqrt(px**2+py**2+pz**2) + + pxp=px*np.cos(tilt)-pz*np.sin(tilt) + pyp=py + pzp=pz*np.cos(tilt)+px*np.sin(tilt) + + array['eta_truth']=1/2*np.log((p+pzp)/(p-pzp)) + array['phi_truth']=np.arctan2(pyp,pxp) + +for p in 50,: + E=arrays_sim[p]["HcalEndcapPInsertHits.energy"] + y, x,_=plt.hist(1e3*ak.flatten(E),bins=100, range=(0, 1.2), histtype='step') + bc=(x[1:]+x[:-1])/2 + from scipy.optimize import curve_fit + slc=abs(bc-.48)<.15 + fnc=Landau + p0=[100, .5, .05] + #print(list(y), list(x)) + coeff, var_matrix = curve_fit(fnc, list(bc[slc]), list(y[slc]), p0=p0, + sigma=list(np.sqrt(y[slc])+(y[slc]==0)), maxfev=10000) + print(coeff) + xx=np.linspace(0,.7, 100) + MIP=coeff[1]/1000 + plt.plot(xx, fnc(xx,*coeff), label=f'Landau fit:\nMIP={coeff[1]*1e3:.0f}$\\pm${1e3*np.sqrt(var_matrix[1][1]):.0f} keV') + plt.xlabel("hit energy [MeV]") + plt.ylabel("hits") + plt.title(f"$E_{{\\mu^-}}=${p} GeV") + plt.legend(fontsize=20) + plt.savefig(outdir+"/MIP.pdf") + + plt.figure(figsize=(10,7)) + array=arrays_sim[p] + bins=30; r=((-np.pi, np.pi),(2.8, 4.2)) + selection=np.sum(array["HcalEndcapPInsertHits.energy"]>0.5*MIP,axis=-1)>0 + h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['eta_truth']), bins=bins, range=r) + h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['eta_truth']), bins=bins, range=r) + + h = h1 / h2 + pc=plt.pcolor(xedges, yedges, h.T,linewidth=0) + plt.xlabel("$\\phi^*$ [rad]") + plt.ylabel("$\\eta^*$") + cb = plt.colorbar(pc) + cb.set_label("acceptance") + plt.title(f"$E_{{\\mu^-}}=${p} GeV") + plt.savefig(outdir+"/acceptance.pdf") From 996e71e462555c300ce84103e1c9bf4ce923067e Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Wed, 23 Oct 2024 17:25:45 -0400 Subject: [PATCH 36/46] Delete benchmarks/insert_tau/analysis/muon_plots.py --- benchmarks/insert_tau/analysis/muon_plots.py | 86 -------------------- 1 file changed, 86 deletions(-) delete mode 100644 benchmarks/insert_tau/analysis/muon_plots.py diff --git a/benchmarks/insert_tau/analysis/muon_plots.py b/benchmarks/insert_tau/analysis/muon_plots.py deleted file mode 100644 index 4c45aac3..00000000 --- a/benchmarks/insert_tau/analysis/muon_plots.py +++ /dev/null @@ -1,86 +0,0 @@ -import numpy as np, pandas as pd, matplotlib.pyplot as plt, matplotlib as mpl, awkward as ak, sys -import mplhep as hep -hep.style.use("CMS") - -plt.rcParams['figure.facecolor']='white' -plt.rcParams['savefig.facecolor']='white' -plt.rcParams['savefig.bbox']='tight' - -plt.rcParams["figure.figsize"] = (7, 7) - -config=sys.argv[1].split("/")[1] #results/{config}/{benchmark_name} -outdir=sys.argv[1]+"/" -try: - import os - os.mkdir(outdir[:-1]) -except: - pass - -def Landau(x, normalization,location,stdev): - #print(type(x)) - u=(x-location)*3.591/stdev/2.355 - renormalization = 1.64872*normalization - return renormalization * np.exp(-u/2 - np.exp(-u)/2) - -import uproot as ur -arrays_sim={} -momenta=50, -for p in momenta: - filename=f'sim_output/insert_muon/{config}_sim_mu-_{p}GeV.edm4hep.root' - print("opening file", filename) - events = ur.open(filename+':events') - arrays_sim[p] = events.arrays() - import gc - gc.collect() - print("read", filename) - -for array in arrays_sim.values(): - tilt=-0.025 - px=array['MCParticles.momentum.x'][:,2] - py=array['MCParticles.momentum.y'][:,2] - pz=array['MCParticles.momentum.z'][:,2] - p=np.sqrt(px**2+py**2+pz**2) - - pxp=px*np.cos(tilt)-pz*np.sin(tilt) - pyp=py - pzp=pz*np.cos(tilt)+px*np.sin(tilt) - - array['eta_truth']=1/2*np.log((p+pzp)/(p-pzp)) - array['phi_truth']=np.arctan2(pyp,pxp) - -for p in 50,: - E=arrays_sim[p]["HcalEndcapPInsertHits.energy"] - y, x,_=plt.hist(1e3*ak.flatten(E),bins=100, range=(0, 1.2), histtype='step') - bc=(x[1:]+x[:-1])/2 - from scipy.optimize import curve_fit - slc=abs(bc-.48)<.15 - fnc=Landau - p0=[100, .5, .05] - #print(list(y), list(x)) - coeff, var_matrix = curve_fit(fnc, list(bc[slc]), list(y[slc]), p0=p0, - sigma=list(np.sqrt(y[slc])+(y[slc]==0))) - print(coeff) - xx=np.linspace(0,.7, 100) - MIP=coeff[1]/1000 - plt.plot(xx, fnc(xx,*coeff), label=f'Landau fit:\nMIP={coeff[1]*1e3:.0f}$\\pm${1e3*np.sqrt(var_matrix[1][1]):.0f} keV') - plt.xlabel("hit energy [MeV]") - plt.ylabel("hits") - plt.title(f"$E_{{\\mu^-}}=${p} GeV") - plt.legend(fontsize=20) - plt.savefig(outdir+"/MIP.pdf") - - plt.figure(figsize=(10,7)) - array=arrays_sim[p] - bins=30; r=((-np.pi, np.pi),(2.8, 4.2)) - selection=np.sum(array["HcalEndcapPInsertHits.energy"]>0.5*MIP,axis=-1)>0 - h1, xedges, yedges = np.histogram2d(list(array[selection]['phi_truth']),list(array[selection]['eta_truth']), bins=bins, range=r) - h2, xedges, yedges = np.histogram2d(list(array['phi_truth']),list(array['eta_truth']), bins=bins, range=r) - - h = h1 / h2 - pc=plt.pcolor(xedges, yedges, h.T,linewidth=0) - plt.xlabel("$\\phi^*$ [rad]") - plt.ylabel("$\\eta^*$") - cb = plt.colorbar(pc) - cb.set_label("acceptance") - plt.title(f"$E_{{\\mu^-}}=${p} GeV") - plt.savefig(outdir+"/acceptance.pdf") From cef20c9cf1e63c6faa853a6e0e0ee8f50162fc34 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 25 Oct 2024 10:34:45 -0400 Subject: [PATCH 37/46] Update Snakefile incorporated changes similar to those from pull request #88 --- benchmarks/insert_tau/Snakefile | 79 ++++++++++++++++++--------------- 1 file changed, 43 insertions(+), 36 deletions(-) diff --git a/benchmarks/insert_tau/Snakefile b/benchmarks/insert_tau/Snakefile index 2e6b120e..37cdef34 100644 --- a/benchmarks/insert_tau/Snakefile +++ b/benchmarks/insert_tau/Snakefile @@ -1,33 +1,39 @@ +def get_n_events(wildcards): + energy = float(wildcards.P) + n_events = 1000 + n_events = int(n_events // ((energy / 20) ** 0.5)) + return n_events + rule insert_tau_generate: - input: - script="benchmarks/insert_tau/analysis/gen_particles.cxx", - params: - NEVENTS_GEN=1000, - th_max=7.0, - th_min=1.7 - output: - GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc", - shell: - """ -root -l -b -q '{input.script}({params.NEVENTS_GEN},"{output.GEN_FILE}", "tau-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' + input: + script="benchmarks/insert_tau/analysis/gen_particles.cxx", + params: + N_EVENTS=get_n_events, + th_max=7.0, + th_min=1.7 + output: + GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc", + shell: + """ +root -l -b -q '{input.script}({params.N_EVENTS},"{output.GEN_FILE}", "tau-", {params.th_min}, {params.th_max}, 0., 360., {wildcards.P})' """ rule insert_tau_simulate: - input: - GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" - warmup="warmup/{DETECTOR_CONFIG}.edm4hep.root" - params: - PHYSICS_LIST="FTFP_BERT" - output: - SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV_{INDEX}.edm4hep.root" - shell: - """ -NEVENTS_SIM=1000 + input: + GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" + warmup="warmup/{DETECTOR_CONFIG}.edm4hep.root" + params: + N_EVENTS=get_n_events, + PHYSICS_LIST="FTFP_BERT" + output: + SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV_{INDEX}.edm4hep.root" + shell: + """ # Running simulation npsim \ --compactFile $DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml \ - --numberOfEvents $NEVENTS_SIM \ - --skipNEvents $(( $NEVENTS_SIM * {wildcards.INDEX} )) \ + --numberOfEvents {params.N_EVENTS} \ + --skipNEvents $(( {params.N_EVENTS} * {wildcards.INDEX} )) \ --physicsList {params.PHYSICS_LIST} \ --inputFiles {input.GEN_FILE} \ --outputFile {output.SIM_FILE} @@ -36,26 +42,27 @@ npsim \ rule insert_tau_recon: input: SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV_{INDEX}.edm4hep.root" + params: + N_EVENTS=get_n_events, output: REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV_{INDEX}.edm4hep.root" shell: """ -NEVENTS_REC=1000 -eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters,EcalEndcapPClusters,LFHCALClusters -Pjana:nevents=$NEVENTS_REC +eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters,EcalEndcapPClusters,LFHCALClusters -Pjana:nevents={params.N_EVENTS} """ rule insert_tau_analysis: - input: - expand("sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep_{INDEX}.root", - P=[20,30, 40,50,60, 80, 100], - DETECTOR_CONFIG=["{DETECTOR_CONFIG}"], - INDEX=range(5), - ), - script="benchmarks/insert_tau/analysis/tau_plots.py", - output: - results_dir=directory("results/{DETECTOR_CONFIG}/insert_tau"), - shell: - """ + input: + expand("sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep_{INDEX}.root", + P=[20,30, 40,50,60, 80, 100], + DETECTOR_CONFIG=["{DETECTOR_CONFIG}"], + INDEX=range(5), + ), + script="benchmarks/insert_tau/analysis/tau_plots.py", + output: + results_dir=directory("results/{DETECTOR_CONFIG}/insert_tau"), + shell: + """ mkdir -p {output.results_dir} python {input.script} {output.results_dir} """ From 5ec80d07aca71bd381e41ce8279738e547591f38 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 25 Oct 2024 10:42:40 -0400 Subject: [PATCH 38/46] Update tau_plots.py created tau_plots.py --- benchmarks/insert_tau/analysis/tau_plots.py | 155 ++++++++++++++++++++ 1 file changed, 155 insertions(+) diff --git a/benchmarks/insert_tau/analysis/tau_plots.py b/benchmarks/insert_tau/analysis/tau_plots.py index e69de29b..e2653833 100644 --- a/benchmarks/insert_tau/analysis/tau_plots.py +++ b/benchmarks/insert_tau/analysis/tau_plots.py @@ -0,0 +1,155 @@ +import numpy as np, pandas as pd, matplotlib.pyplot as plt, matplotlib as mpl, awkward as ak, sys +import mplhep as hep +hep.style.use("CMS") + +plt.rcParams['figure.facecolor']='white' +plt.rcParams['savefig.facecolor']='white' +plt.rcParams['savefig.bbox']='tight' + +plt.rcParams["figure.figsize"] = (7, 7) + +config=sys.argv[1].split("/")[1] #results/{config}/{benchmark_name} +outdir=sys.argv[1]+"/" +try: + import os + os.mkdir(outdir[:-1]) +except: + pass + +import uproot as ur +arrays_sim={p:ur.open(f'sim_output/insert_tau/{config}_rec_tau-_{p}GeV.edm4eic.root:events').arrays() for p in (20, 30, 40, 50, 60,80,100)} + +for a in arrays_sim.values(): + #recon + Etot=0 + px=0 + py=0 + pz=0 + for det in "HcalEndcapPInsert", "EcalEndcapPInsert", "EcalEndcapP", "LFHCAL": + + E=a[f'{det}Clusters.energy'] + + + #if det=="EcalEndcapPInsert": + # E=E/1.08 + E=E*(-0*E+1.2) + + + #uncorr=(e/w+h) + #s=-0.0064*uncorr+1.80 + #r=uncorr*s #reconstructed energy with correction + + x=a[f'{det}Clusters.position.x'] + y=a[f'{det}Clusters.position.y'] + z=a[f'{det}Clusters.position.z'] + r=np.sqrt(x**2+y**2+z**2) + Etot=Etot+np.sum(E, axis=-1) + px=px+np.sum(E*x/r,axis=-1) + py=py+np.sum(E*y/r,axis=-1) + pz=pz+np.sum(E*z/r,axis=-1) + + a['jet_p_recon']=np.sqrt(px**2+py**2+pz**2) + a['jet_E_recon']=Etot + + a['jet_theta_recon']=np.arctan2(np.hypot(px*np.cos(-.025)-pz*np.sin(-.025),py), + pz*np.cos(-.025)+px*np.sin(-.025)) + + #truth + m=a['MCParticles.mass'][::,2:] + px=a['MCParticles.momentum.x'][::,2:] + py=a['MCParticles.momentum.y'][::,2:] + pz=a['MCParticles.momentum.z'][::,2:] + E=np.sqrt(m**2+px**2+py**2+pz**2) + status=a['MCParticles.simulatorStatus'][::,2:] + PDG=a['MCParticles.PDG'][::,2:] + + #find the hadronic part: initial-state tau - all leptons + selection=1*(PDG==15)-1*(np.abs(PDG)==16) + is_hadronic=1*(np.sum((PDG==-14)+(PDG==-12), axis=-1)==0) + + px_hfs, py_hfs, pz_hfs= np.sum(px*selection,axis=-1)*is_hadronic,np.sum(py*selection,axis=-1)*is_hadronic, np.sum(pz*selection,axis=-1)*is_hadronic + + a['hfs_p_truth']=np.sqrt(px_hfs**2+py_hfs**2+pz_hfs**2) + a['hfs_E_truth']=np.sum(E*selection,axis=-1)*is_hadronic + + + a['hfs_theta_truth']=np.arctan2(np.hypot(px_hfs*np.cos(-.025)-pz_hfs*np.sin(-.025),py_hfs), + pz_hfs*np.cos(-.025)+px_hfs*np.sin(-.025)) + a['hfs_eta_truth']=-np.log(np.tan(a['hfs_theta_truth']/2)) + a['n_mu']=np.sum(np.abs(PDG)==13, axis=-1) + a['n_e']=np.sum(np.abs(PDG)==13, axis=-1) + a['hfs_mass_truth']=np.sqrt(a['hfs_E_truth']**2-a['hfs_p_truth']**2) + +for a in arrays_sim.values(): + selection=(a['hfs_eta_truth']>3.1) & (a['hfs_eta_truth']<3.8)\ + &(a['n_mu']==0)&(a['n_e']==0)&(a['hfs_mass_truth']>.140)&(a['jet_E_recon']>0) + +Etruth=[] +Erecon=[] + +theta_truth=[] +theta_recon=[] + +eta_max=3.7 +eta_min=3.3 +for a in arrays_sim.values(): + selection=(a['hfs_eta_truth']>eta_min) & (a['hfs_eta_truth'].140)&(a['jet_E_recon']>1) + theta_truth=np.concatenate((theta_truth,a['hfs_theta_truth'][selection])) + theta_recon=np.concatenate((theta_recon,a['jet_theta_recon'][selection])) + Etruth=np.concatenate((Etruth,a['hfs_E_truth'][selection])) + Erecon=np.concatenate((Erecon,a['jet_E_recon'][selection])) + +plt.figure() +plt.scatter(theta_truth, theta_recon, 1) +plt.xlabel("$\\theta^{hfs}_{\\rm truth}$ [rad]") +plt.ylabel("$\\theta^{hfs}_{\\rm rec}$ [rad]") +plt.title(f"$E_{{\\tau}}$=20-100 GeV, ${eta_min}<\\eta_{{hfs}}<{eta_max}$") +plt.plot([0.04,0.1], [0.04, 0.1], color='tab:orange') +plt.ylim(0, 0.15) +plt.savefig(outdir+"/theta_scatter.pdf") + +plt.figure() +plt.scatter(Etruth, Erecon, 1) +plt.xlabel("$E^{hfs}_{\\rm truth}$ [GeV]") +plt.ylabel("$E^{hfs}_{\\rm rec}$ [GeV]") +plt.title(f"$E_{{\\tau}}$=20-100 GeV, ${eta_min}<\\eta_{{hfs}}<{eta_max}$") +plt.plot((0,100), (0, 100), color='tab:orange') +plt.savefig(outdir+"/energy_scatter.pdf") + +def gauss(x, A,mu, sigma): + return A * np.exp(-(x-mu)**2/(2*sigma**2)) +from scipy.optimize import curve_fit + +res=[] +dres=[] +emid=[] +ew=[] +partitions=(20,30, 40, 60,80,100) +for emin, emax in zip(partitions[:-1], partitions[1:]): + + y,x = np.histogram((theta_recon-theta_truth)[(Etruth>emin)&(Etruth Date: Fri, 25 Oct 2024 10:46:39 -0400 Subject: [PATCH 39/46] Update tau_plots.py read from the correct file locations --- benchmarks/insert_tau/analysis/tau_plots.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/benchmarks/insert_tau/analysis/tau_plots.py b/benchmarks/insert_tau/analysis/tau_plots.py index e2653833..a944a2b5 100644 --- a/benchmarks/insert_tau/analysis/tau_plots.py +++ b/benchmarks/insert_tau/analysis/tau_plots.py @@ -17,7 +17,14 @@ pass import uproot as ur -arrays_sim={p:ur.open(f'sim_output/insert_tau/{config}_rec_tau-_{p}GeV.edm4eic.root:events').arrays() for p in (20, 30, 40, 50, 60,80,100)} +arrays_sim={} +momenta=20, 30, 40, 50, 60,80,100 +for p in momenta: + arrays_sim[p] = ur.concatenate({ + f'sim_output/insert_tau/{config}_rec_tau-_{p}GeV_{index}.edm4eic.root': 'events' + for index in range(5) + }) + for a in arrays_sim.values(): #recon From 0ce0843d21a63128ec7e5705905eded47343854e Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Fri, 25 Oct 2024 10:50:58 -0400 Subject: [PATCH 40/46] Update tau_plots.py --- benchmarks/insert_tau/analysis/tau_plots.py | 10 +--------- 1 file changed, 1 insertion(+), 9 deletions(-) diff --git a/benchmarks/insert_tau/analysis/tau_plots.py b/benchmarks/insert_tau/analysis/tau_plots.py index a944a2b5..f9f5c862 100644 --- a/benchmarks/insert_tau/analysis/tau_plots.py +++ b/benchmarks/insert_tau/analysis/tau_plots.py @@ -36,15 +36,7 @@ E=a[f'{det}Clusters.energy'] - - #if det=="EcalEndcapPInsert": - # E=E/1.08 - E=E*(-0*E+1.2) - - - #uncorr=(e/w+h) - #s=-0.0064*uncorr+1.80 - #r=uncorr*s #reconstructed energy with correction + #todo apply corrections depending on whether this is an electromagnetic or hadronic shower. x=a[f'{det}Clusters.position.x'] y=a[f'{det}Clusters.position.y'] From e6c0006abbb5444bf9b8e91ba0dd6313e8b6a415 Mon Sep 17 00:00:00 2001 From: Dmitry Kalinkin Date: Sun, 27 Oct 2024 20:56:53 -0400 Subject: [PATCH 41/46] config.yml: remove reference to `common:detector` Not a thing since https://github.com/eic/detector_benchmarks/pull/69 --- benchmarks/insert_tau/config.yml | 1 - 1 file changed, 1 deletion(-) diff --git a/benchmarks/insert_tau/config.yml b/benchmarks/insert_tau/config.yml index 14ffa34a..719c0705 100644 --- a/benchmarks/insert_tau/config.yml +++ b/benchmarks/insert_tau/config.yml @@ -1,7 +1,6 @@ sim:insert_tau: stage: simulate extends: .det_benchmark - needs: ["common:detector"] parallel: matrix: - P: 20 From 0b024338474d47905d6a6c3720a80922eba9bb52 Mon Sep 17 00:00:00 2001 From: Dmitry Kalinkin Date: Mon, 28 Oct 2024 03:00:46 -0400 Subject: [PATCH 42/46] Snakefile: trying to fix what breaks snakemake's parser --- benchmarks/insert_tau/Snakefile | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/benchmarks/insert_tau/Snakefile b/benchmarks/insert_tau/Snakefile index 37cdef34..b200c0ad 100644 --- a/benchmarks/insert_tau/Snakefile +++ b/benchmarks/insert_tau/Snakefile @@ -10,7 +10,7 @@ rule insert_tau_generate: params: N_EVENTS=get_n_events, th_max=7.0, - th_min=1.7 + th_min=1.7, output: GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc", shell: @@ -24,7 +24,7 @@ rule insert_tau_simulate: warmup="warmup/{DETECTOR_CONFIG}.edm4hep.root" params: N_EVENTS=get_n_events, - PHYSICS_LIST="FTFP_BERT" + PHYSICS_LIST="FTFP_BERT", output: SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV_{INDEX}.edm4hep.root" shell: @@ -40,21 +40,21 @@ npsim \ """ rule insert_tau_recon: - input: - SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV_{INDEX}.edm4hep.root" + input: + SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV_{INDEX}.edm4hep.root" params: N_EVENTS=get_n_events, - output: - REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV_{INDEX}.edm4hep.root" - shell: - """ + output: + REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV_{INDEX}.edm4hep.root" + shell: + """ eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters,EcalEndcapPClusters,LFHCALClusters -Pjana:nevents={params.N_EVENTS} """ rule insert_tau_analysis: input: expand("sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep_{INDEX}.root", - P=[20,30, 40,50,60, 80, 100], + P=[20, 30, 40, 50, 60, 80, 100], DETECTOR_CONFIG=["{DETECTOR_CONFIG}"], INDEX=range(5), ), From 5a9065c64900f048dceded259860b192cde335bc Mon Sep 17 00:00:00 2001 From: Dmitry Kalinkin Date: Mon, 28 Oct 2024 03:02:11 -0400 Subject: [PATCH 43/46] tau_plots.py: maxfev=10000 --- benchmarks/insert_tau/analysis/tau_plots.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/benchmarks/insert_tau/analysis/tau_plots.py b/benchmarks/insert_tau/analysis/tau_plots.py index f9f5c862..77c2dd32 100644 --- a/benchmarks/insert_tau/analysis/tau_plots.py +++ b/benchmarks/insert_tau/analysis/tau_plots.py @@ -135,7 +135,7 @@ def gauss(x, A,mu, sigma): fnc=gauss sigma=np.sqrt(y[slc])+(y[slc]==0) - coeff, var_matrix = curve_fit(fnc, list(bc[slc]), list(y[slc]), p0=p0,sigma=list(sigma)) + coeff, var_matrix = curve_fit(fnc, list(bc[slc]), list(y[slc]), p0=p0, sigma=list(sigma), maxfev=10000) res.append(abs(coeff[2])) dres.append(np.sqrt(var_matrix[2][2])) emid.append((emin+emax)/2) @@ -147,7 +147,7 @@ def gauss(x, A,mu, sigma): fnc=lambda E,B:B/E p0=[1,] -coeff, var_matrix = curve_fit(fnc, emid, res, p0=p0,sigma=list(dres)) +coeff, var_matrix = curve_fit(fnc, emid, res, p0=p0, sigma=list(dres), maxfev=10000) xx=np.linspace(10, 100, 100) plt.plot(xx, 1000*fnc(xx, *coeff), label=f"fit: ${coeff[0]:.2f}/E$ mrad") plt.legend() From 07b0b1e0f8125614865ce460b67ee9fb57b71542 Mon Sep 17 00:00:00 2001 From: Dmitry Kalinkin Date: Mon, 28 Oct 2024 03:14:24 -0400 Subject: [PATCH 44/46] Snakefile: add more missing commas --- benchmarks/insert_tau/Snakefile | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/benchmarks/insert_tau/Snakefile b/benchmarks/insert_tau/Snakefile index b200c0ad..4e81284f 100644 --- a/benchmarks/insert_tau/Snakefile +++ b/benchmarks/insert_tau/Snakefile @@ -20,8 +20,8 @@ root -l -b -q '{input.script}({params.N_EVENTS},"{output.GEN_FILE}", "tau-", {pa rule insert_tau_simulate: input: - GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc" - warmup="warmup/{DETECTOR_CONFIG}.edm4hep.root" + GEN_FILE="sim_output/insert_tau/tau-_{P}GeV.hepmc", + warmup="warmup/{DETECTOR_CONFIG}.edm4hep.root", params: N_EVENTS=get_n_events, PHYSICS_LIST="FTFP_BERT", @@ -41,11 +41,11 @@ npsim \ rule insert_tau_recon: input: - SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV_{INDEX}.edm4hep.root" + SIM_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_sim_tau-_{P}GeV_{INDEX}.edm4hep.root", params: N_EVENTS=get_n_events, output: - REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV_{INDEX}.edm4hep.root" + REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV_{INDEX}.edm4hep.root", shell: """ eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters,EcalEndcapPClusters,LFHCALClusters -Pjana:nevents={params.N_EVENTS} From 7ec42945088881fad03e3ffa3fbbc2cfef75d97f Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Mon, 4 Nov 2024 10:43:07 -0500 Subject: [PATCH 45/46] Update Snakefile fixed a filename bug --- benchmarks/insert_tau/Snakefile | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/benchmarks/insert_tau/Snakefile b/benchmarks/insert_tau/Snakefile index 4e81284f..57e4390e 100644 --- a/benchmarks/insert_tau/Snakefile +++ b/benchmarks/insert_tau/Snakefile @@ -45,7 +45,7 @@ rule insert_tau_recon: params: N_EVENTS=get_n_events, output: - REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV_{INDEX}.edm4hep.root", + REC_FILE="sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV_{INDEX}.edm4eic.root", shell: """ eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_files=$DETECTOR_PATH/{wildcards.DETECTOR_CONFIG}.xml -Ppodio:output_collections=MCParticles,HcalEndcapPInsertRecHits,HcalEndcapPInsertClusters,HcalEndcapPInsertSubcellHits,EcalEndcapPInsertRecHits,EcalEndcapPInsertClusters,EcalEndcapPClusters,LFHCALClusters -Pjana:nevents={params.N_EVENTS} @@ -53,7 +53,7 @@ eicrecon {input.SIM_FILE} -Ppodio:output_file={output.REC_FILE} -Pdd4hep:xml_fil rule insert_tau_analysis: input: - expand("sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV.edm4hep_{INDEX}.root", + expand("sim_output/insert_tau/{DETECTOR_CONFIG}_rec_tau-_{P}GeV_{INDEX}.edm4eic.root", P=[20, 30, 40, 50, 60, 80, 100], DETECTOR_CONFIG=["{DETECTOR_CONFIG}"], INDEX=range(5), From 22a69de3af826494a18765019da3323e73969951 Mon Sep 17 00:00:00 2001 From: Sebouh Paul Date: Mon, 4 Nov 2024 16:17:50 -0500 Subject: [PATCH 46/46] Update config.yml extends: .det_benchmark --- benchmarks/insert_tau/config.yml | 1 + 1 file changed, 1 insertion(+) diff --git a/benchmarks/insert_tau/config.yml b/benchmarks/insert_tau/config.yml index 719c0705..1c43d0e2 100644 --- a/benchmarks/insert_tau/config.yml +++ b/benchmarks/insert_tau/config.yml @@ -26,6 +26,7 @@ bench:insert_tau: - snakemake $SNAKEMAKE_FLAGS --cores 1 results/epic_craterlake/insert_tau collect_results:insert_tau: + extends: .det_benchmark stage: collect needs: ["bench:insert_tau"] script: