Skip to content
Snippets Groups Projects
Commit 5b3fecbb authored by Whitney Armstrong's avatar Whitney Armstrong
Browse files

Updated script. emcal_electron

- Script can take # of events and file name as arguments
- `root -b -q "emcal_electron.cxx(100,\"file_name.root\")" `
parent 7bcdcb38
No related branches found
No related tags found
1 merge request!14Updated script. emcal_electron
...@@ -19,52 +19,61 @@ ...@@ -19,52 +19,61 @@
using namespace HepMC3; using namespace HepMC3;
void emcal_electrons(){ void emcal_electrons(int n_events = 1e6, const char* out_fname = "./data/emcal_electron_0GeVto30GeV_100kEvt.hepmc")
WriterAscii hepmc_output("./data/emcal_electron_0GeVto30GeV_100kEvt.hepmc"); {
int events_parsed = 0; WriterAscii hepmc_output(out_fname);
GenEvent evt(Units::GEV, Units::MM); int events_parsed = 0;
GenEvent evt(Units::GEV, Units::MM);
// Random number generator // Random number generator
TRandom *r1 = new TRandom(); TRandom *r1 = new TRandom();
for (events_parsed = 0; events_parsed < 100000; events_parsed++) { for (events_parsed = 0; events_parsed < n_events; events_parsed++) {
// FourVector(px,py,pz,e,pdgid,status) // FourVector(px,py,pz,e,pdgid,status)
// type 4 is beam // type 4 is beam
// pdgid 11 - electron // pdgid 11 - electron
// pdgid 2212 - proton // pdgid 111 - pi0
GenParticlePtr p1 = std::make_shared<GenParticle>(FourVector(0.0, 0.0, 10.0, 10.0), 11, 4); // pdgid 2212 - proton
GenParticlePtr p2 = std::make_shared<GenParticle>(FourVector(0.0, 0.0, 0.0, 0.938), 2212, 4); GenParticlePtr p1 =
std::make_shared<GenParticle>(FourVector(0.0, 0.0, 10.0, 10.0), 11, 4);
GenParticlePtr p2 = std::make_shared<GenParticle>(
FourVector(0.0, 0.0, 0.0, 0.938), 2212, 4);
// Define momentum // Define momentum
Double_t p = 0.0 + events_parsed * 3.e-4; Double_t p = 0.0 + events_parsed * 3.e-4;
Double_t px; Double_t px;
Double_t py; Double_t py;
Double_t pz; Double_t pz;
// Generates random vectors, uniformly distributed over the surface of a sphere of given radius, in this case momentum. // Generates random vectors, uniformly distributed over the surface of a
r1->Sphere(px,py,pz,p); // sphere of given radius, in this case momentum.
r1->Sphere(px, py, pz, p);
// type 1 is final state // type 1 is final state
// pdgid 11 - electron 0.510 MeV/c^2 // pdgid 11 - electron 0.510 MeV/c^2
GenParticlePtr p3 = std::make_shared<GenParticle>(FourVector(px,py,pz,sqrt((px*px)+(py*py)+(pz*pz)+(0.000511*0.000511))), 11, 1); GenParticlePtr p3 = std::make_shared<GenParticle>(
FourVector(
px, py, pz,
sqrt((px * px) + (py * py) + (pz * pz) + (0.000511 * 0.000511))),
11, 1);
GenVertexPtr v1 = std::make_shared<GenVertex>(); GenVertexPtr v1 = std::make_shared<GenVertex>();
v1->add_particle_in(p1); v1->add_particle_in(p1);
v1->add_particle_in(p2); v1->add_particle_in(p2);
v1->add_particle_out(p3); v1->add_particle_out(p3);
evt.add_vertex(v1); evt.add_vertex(v1);
if (events_parsed == 0) { if (events_parsed == 0) {
std::cout << "First event: " << std::endl; std::cout << "First event: " << std::endl;
Print::listing(evt); Print::listing(evt);
} }
hepmc_output.write_event(evt); hepmc_output.write_event(evt);
if (events_parsed % 10 == 0) { if (events_parsed % 10 == 0) {
std::cout << "Event: " << events_parsed << std::endl; std::cout << "Event: " << events_parsed << std::endl;
} }
evt.clear(); evt.clear();
} }
hepmc_output.close(); hepmc_output.close();
std::cout << "Events parsed and written: " << events_parsed << std::endl; std::cout << "Events parsed and written: " << events_parsed << std::endl;
} }
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment