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

Squashed commit of the following:

commit cbee19db
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Mon Aug 16 15:23:55 2021 -0400

    This is really not the place to do this test. Should be done with final output file after reco.

commit 582f72b6
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Mon Aug 16 14:52:58 2021 -0400

    Something got messed up in the rebase. Fixed and now testing original functionality.

commit 510cb6cb
Merge: 0d80f473 548e97b5
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Mon Aug 16 14:40:25 2021 -0400

    Merge branch '50-b0tracker-acceptance-tests' of https://eicweb.phy.anl.gov/EIC/benchmarks/detector_benchmarks into 50-b0tracker-acceptance-tests
    Still trying to add MC vs. accepted information - rebase issue.

commit 0d80f473
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Mon Aug 16 14:37:54 2021 -0400

    Trying to grab MC data so I can make some very simple acceptance plots.

commit 548e97b5
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 11:42:19 2021 -0400

    Add crossing angle to gen_particle script.

commit d60f1fa5
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 11:18:21 2021 -0400

    Things seem to work now. Need more events to test operation.

commit 21bab2a4
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 11:09:30 2021 -0400

    Another issue.

commit cc86cb9a
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 11:03:21 2021 -0400

    Fixed typo(s).

commit e24c2728
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 10:49:10 2021 -0400

    Trying to add occupancy plots for the B0.

commit 726efc4f
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 09:11:30 2021 -0400

    Fixed typo - still getting used to syntax.

commit a2e0d295
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 09:06:41 2021 -0400

    Working on updating particle gun and hits reader for B0.

commit 0bcbf163
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 11:42:19 2021 -0400

    Add crossing angle to gen_particle script.

commit 2b677af9
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 11:18:21 2021 -0400

    Things seem to work now. Need more events to test operation.

commit cab2b05c
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 11:09:30 2021 -0400

    Another issue.

commit 21eeb20d
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 11:03:21 2021 -0400

    Fixed typo(s).

commit d3a75da2
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 10:49:10 2021 -0400

    Trying to add occupancy plots for the B0.

commit 8984a51f
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 09:11:30 2021 -0400

    Fixed typo - still getting used to syntax.

commit cf469f0c
Author: Alexander Jentsch <ajentsch@bnl.gov>
Date:   Wed Aug 11 09:06:41 2021 -0400

    Working on updating particle gun and hits reader for B0.
parent 53e627d2
No related branches found
No related tags found
1 merge request!77Resolve "B0Tracker acceptance tests"
...@@ -7,7 +7,7 @@ if [[ ! -n "${JUGGLER_DETECTOR}" ]] ; then ...@@ -7,7 +7,7 @@ if [[ ! -n "${JUGGLER_DETECTOR}" ]] ; then
fi fi
if [[ ! -n "${JUGGLER_N_EVENTS}" ]] ; then if [[ ! -n "${JUGGLER_N_EVENTS}" ]] ; then
export JUGGLER_N_EVENTS=100 export JUGGLER_N_EVENTS=1000
fi fi
export FILE_NAME_TAG="forward_protons" export FILE_NAME_TAG="forward_protons"
......
...@@ -36,25 +36,67 @@ void b0_tracker_hits(const char* fname = "./sim_output/sim_forward_protons.root" ...@@ -36,25 +36,67 @@ void b0_tracker_hits(const char* fname = "./sim_output/sim_forward_protons.root"
ROOT::RDataFrame d0(*t); ROOT::RDataFrame d0(*t);
auto hits_eta = [&](const std::vector<dd4pod::TrackerHitData>& hits) { auto hits_theta = [&](const std::vector<dd4pod::TrackerHitData>& hits) {
std::vector<double> result; std::vector<double> result;
for (const auto& h : hits) { for (const auto& h : hits) {
ROOT::Math::XYZVector vec(h.position.x,h.position.y,h.position.z); ROOT::Math::XYZVector vec(h.position.x,h.position.y,h.position.z);
result.push_back(vec.eta()); result.push_back(1000*vec.theta());
std::cout << vec.eta() << "\n"; std::cout << 1000*vec.theta() << "\n";
} }
return result; return result;
}; };
auto d1 = d0.Define("hits_eta", hits_eta, {"B0TrackerHits"}); auto local_position = [&](const std::vector<dd4pod::TrackerHitData>& hits) {
std::vector<std::array<double, 2>> result;
for (const auto& h : hits) {
auto pos0 = (h.position);
result.push_back({pos0.x , pos0.y});
}
return result;
};
auto x_pos = [&](const std::vector<std::array<double, 2>>& xypos) {
std::vector<double> result;
for (const auto& h : xypos) {
result.push_back(h.at(0));
}
return result;
};
auto h1 = d1.Histo1D({"h1", "hits_eta", 100, 0,20}, "hits_eta"); auto y_pos = [&](const std::vector<std::array<double, 2>>& xypos) {
std::vector<double> result;
for (const auto& h : xypos) {
result.push_back(h.at(1));
}
return result;
};
auto d1 = d0.Define("nhits", hits_theta, {"B0TrackerHits"})
.Define("xy_hit_pos", local_position, {"B0TrackerHits"})
.Define("x_pos", x_pos, {"xy_hit_pos"})
.Define("y_pos", y_pos, {"xy_hit_pos"});
auto h_local_pos = d1.Histo2D({"h_local_pos", ";x [mm]; y [mm] ", 100, -100.0, -200.0, 100, -100.0, 100.0}, "x_pos", "y_pos");
auto d2 = d0.Define("hits_theta", hits_theta, {"B0TrackerHits"});
auto h1 = d2.Histo1D({"h1", "hits_theta", 100, 0,20}, "hits_theta");
TCanvas* c = new TCanvas(); TCanvas* c = new TCanvas();
h1->DrawCopy(); h1->DrawCopy();
c->SaveAs("results/b0_tracker_hits_eta.png"); c->SaveAs("results/b0_tracker_hits_theta.png");
c->SaveAs("results/b0_tracker_hits_eta.pdf"); c->SaveAs("results/b0_tracker_hits_theta.pdf");
h_local_pos->DrawCopy("colz");
c->SaveAs("results/b0_tracker_hits_occupancy_disk_1.png");
c->SaveAs("results/b0_tracker_hits_occupancy_disk_1.pdf");
auto n1 = h1->GetMean(); auto n1 = h1->GetMean();
std::cout << "Pseudorapidity of hits: " << n1 << std::endl; std::cout << "Polar angle of hits: " << n1 << std::endl;
//if (n1 < 5) { //if (n1 < 5) {
// std::quick_exit(1); // std::quick_exit(1);
......
...@@ -20,12 +20,19 @@ using namespace HepMC3; ...@@ -20,12 +20,19 @@ using namespace HepMC3;
/** Generate electrons in the central region. /** Generate electrons in the central region.
* This is for testing detectors in the "barrel" region. * This is for testing detectors in the "barrel" region.
*/ */
void gen_forward_protons(int n_events = 100, void gen_forward_protons(int n_events = 1000,
const char* out_fname = "forward_protons.hepmc") const char* out_fname = "forward_protons.hepmc")
{ {
double cos_theta_min = std::cos(0.5*(M_PI/180.0));
double crossingAngle = -0.025; //radiansc
// generate protons in B0 acceptance - roughly 5 - 20 mrad
double cos_theta_min = std::cos(0.02);
double cos_theta_max = std::cos(0.0*(M_PI/180.0)); double cos_theta_max = std::cos(0.0*(M_PI/180.0));
double partEnergyMin = 270.0; // xL 0.98
double partEnergyMax = 275.0; // top beam energy
const double M_p = common_bench::particleMap.at(2212).mass; const double M_p = common_bench::particleMap.at(2212).mass;
WriterAscii hepmc_output(out_fname); WriterAscii hepmc_output(out_fname);
...@@ -42,12 +49,12 @@ void gen_forward_protons(int n_events = 100, ...@@ -42,12 +49,12 @@ void gen_forward_protons(int n_events = 100,
// pdgid 111 - pi0 // pdgid 111 - pi0
// pdgid 2212 - proton // pdgid 2212 - proton
GenParticlePtr p1 = GenParticlePtr p1 =
std::make_shared<GenParticle>(FourVector(0.0, 0.0, 10.0, 10.0), 2212, 4); std::make_shared<GenParticle>(FourVector(0.0, 0.0, partEnergyMax, partEnergyMax), 2212, 4);
GenParticlePtr p2 = std::make_shared<GenParticle>( GenParticlePtr p2 = std::make_shared<GenParticle>(
FourVector(0.0, 0.0, 0.0, M_p), 2212, 4); FourVector(0.0, 0.0, 0.0, M_p), 2212, 4);
// Define momentum // Define momentum
Double_t p = r1->Uniform(200.0, 275.0); Double_t p = r1->Uniform(partEnergyMin, partEnergyMax);
Double_t phi = r1->Uniform(0.0, 2.0 * M_PI); Double_t phi = r1->Uniform(0.0, 2.0 * M_PI);
Double_t costh = r1->Uniform(cos_theta_min, cos_theta_max); Double_t costh = r1->Uniform(cos_theta_min, cos_theta_max);
Double_t th = std::acos(costh); Double_t th = std::acos(costh);
...@@ -55,20 +62,12 @@ void gen_forward_protons(int n_events = 100, ...@@ -55,20 +62,12 @@ void gen_forward_protons(int n_events = 100,
Double_t py = p * std::sin(phi) * std::sin(th); Double_t py = p * std::sin(phi) * std::sin(th);
Double_t pz = p * std::cos(th); Double_t pz = p * std::cos(th);
ROOT::Math::XYZVector p0 = {px,py,pz}; ROOT::Math::XYZVector p0 = {px,py,pz};
//ROOT::Math::Rotation3D r = (-0.025); //ROOT::Math::Rotation3D r = (-0.025);
ROOT::Math::RotationY r(-0.025); ROOT::Math::RotationY r(-0.025);
auto p_rot = r*p0; auto p_rot = r*p0;
// Generates random vectors, uniformly distributed over the surface of a
// sphere of given radius, in this case momentum.
// r1->Sphere(px, py, pz, p);
//std::cout << std::sqrt(px*px + py*py + pz*pz) - p << " is zero? \n";
// 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>( GenParticlePtr p3 = std::make_shared<GenParticle>(
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment