diff --git a/data/emcal_electron_0GeVto30GeV_1000kEvt.hepmc b/data/emcal_electron_0GeVto30GeV_1000kEvt.hepmc new file mode 100644 index 0000000000000000000000000000000000000000..96d6c5786cace9a68bbf4b8efdb4362797ea2333 Binary files /dev/null and b/data/emcal_electron_0GeVto30GeV_1000kEvt.hepmc differ diff --git a/emcal_electrons.cxx b/emcal_electrons.cxx index ca45a209de64a25df72c5b0c628d3e2acb936090..80b21dddbb4f5c7ac2acd816b8814326948ffe46 100644 --- a/emcal_electrons.cxx +++ b/emcal_electrons.cxx @@ -19,7 +19,7 @@ using namespace HepMC3; -void emcal_electrons(int n_events = 1e6, const char* out_fname = "./data/emcal_electron_0GeVto30GeV_100kEvt.hepmc") +void emcal_electrons(int n_events = 1e6, const char* out_fname = "./data/emcal_electron_0GeVto30GeV_1000kEvt.hepmc") { WriterAscii hepmc_output(out_fname); int events_parsed = 0; @@ -28,6 +28,12 @@ void emcal_electrons(int n_events = 1e6, const char* out_fname = "./data/emcal_e // Random number generator TRandom *r1 = new TRandom(); + // Constraining the solid angle, but larger than that subtended by the detector + // https://indico.bnl.gov/event/7449/contributions/35966/attachments/27177/41430/EIC-DWG-Calo-03192020.pdf + // See a figure on slide 26 + double cos_theta_min = std::cos(M_PI * (120.0 / 180.0)); + double cos_theta_max = std::cos(M_PI); + for (events_parsed = 0; events_parsed < n_events; events_parsed++) { // FourVector(px,py,pz,e,pdgid,status) // type 4 is beam @@ -40,13 +46,18 @@ void emcal_electrons(int n_events = 1e6, const char* out_fname = "./data/emcal_e FourVector(0.0, 0.0, 0.0, 0.938), 2212, 4); // Define momentum - Double_t p = r1->Uniform(0.0, 30.0); - Double_t px; - Double_t py; - Double_t pz; + // Momentum starting at 0 GeV/c to 30 GeV/c + Double_t p = r1->Uniform(0.0, 30.0); + Double_t phi = r1->Uniform(0.0, 2.0 * M_PI); + Double_t costheta = r1->Uniform(cos_theta_min, cos_theta_max); + Double_t theta = std::acos(costheta); + Double_t px = p * std::cos(phi) * std::sin(theta); + Double_t py = p * std::sin(phi) * std::sin(theta); + Double_t pz = p * std::cos(theta); + // Generates random vectors, uniformly distributed over the surface of a // sphere of given radius, in this case momentum. - r1->Sphere(px, py, pz, p); + //r1->Sphere(px, py, pz, p); // type 1 is final state // pdgid 11 - electron 0.510 MeV/c^2 diff --git a/emcal_electrons_reader.cxx b/emcal_electrons_reader.cxx index 5537e3b8783be44e464f714157bb1621a5835c73..928f055e609a0892cd57395f4e9dddf495e6b940 100644 --- a/emcal_electrons_reader.cxx +++ b/emcal_electrons_reader.cxx @@ -25,9 +25,9 @@ void emcal_electrons_reader(){ gStyle->SetPadGridX(1); gStyle->SetPadGridY(1); gStyle->SetPadLeftMargin(0.14); - gStyle->SetPadRightMargin(0.14); + gStyle->SetPadRightMargin(0.17); - ReaderAscii hepmc_input("./data/emcal_electron_0GeVto30GeV_100kEvt.hepmc"); + ReaderAscii hepmc_input("./data/emcal_electron_0GeVto30GeV_1000kEvt.hepmc"); int events_parsed = 0; GenEvent evt(Units::GEV, Units::MM); @@ -91,7 +91,7 @@ void emcal_electrons_reader(){ TCanvas* c3 = new TCanvas("c3","c3", 500, 500); h_electrons_phi->GetYaxis()->SetTitleOffset(1.8); h_electrons_phi->SetLineWidth(2); - h_electrons_phi->GetYaxis()->SetRangeUser(0.0,h_electrons_phi->GetMaximum()+100.0); + h_electrons_phi->GetYaxis()->SetRangeUser(0.0,h_electrons_phi->GetMaximum()+300.0); h_electrons_phi->SetLineColor(kBlue); h_electrons_phi->DrawClone(); c3->SaveAs("results/emcal_electrons_phi_0GeVto30GeV.png");