Added GetFirstPointInDetector and GetLastPointInDetector methods and...

Added GetFirstPointInDetector and GetLastPointInDetector methods and first_points and last_points private variables to EDEPTrajectory

examples/EDEPReader.cpp

@@ -23,7 +23,9 @@ int main(int argc, char *argv[])
   // auto fMc = TFile::Open("/storage/gpfs_data/neutrino/SAND/GRAIN-PHYSICS-CASE/prod/numu_LBNF_GV2_prod/data/numu_LBNF_GV2_prod_1.genie.edep-sim.root");
   // auto fMc = TFile::Open("/storage/gpfs_data/neutrino/users/dcasazza/sand-reco/e_pos2_1000.edep.root");
   // auto fMc = TFile::Open("/storage/gpfs_data/neutrino/users/vpia/STTTrackReco/tools/tracklet/DRIFT3.root");
-  auto fMc = TFile::Open("/storage/gpfs_data/neutrino/SAND/PRODUCTIONS/PROD/SAND_opt3_DRIFT1_reverse_volSAND/SAND_opt3_DRIFT1_reverse_volSAND_0/sand-drift-spill-revere-volSAND.0.edep.root");
+  // auto fMc = TFile::Open("/storage/gpfs_data/neutrino/users/rdamico/lotsofnumu/production_14/events-in-rockBox_lv.4.edep.root");
+  auto fMc = TFile::Open("/storage/gpfs_data/neutrino/users/dcasazza/lot-of-nue/production_merge30cm/events-in-STTtracker.33.edep.root");
+  // auto fMc = TFile::Open("/storage/gpfs_data/neutrino/SAND/PRODUCTIONS/PROD/SAND_opt3_DRIFT1_reverse_volSAND/SAND_opt3_DRIFT1_reverse_volSAND_0/sand-drift-spill-revere-volSAND.0.edep.root");
   // TFile *f = TFile::Open("/storage/gpfs_data/neutrino/users/dcasazza/sand-reco/e_pos2_1000_branch4.digit.Clusters.root");
   // TFile *f = TFile::Open("/storage/gpfs_data/neutrino/users/dcasazza/sand-reco/mu_pos2_1000_branch4.digit.Clusters.root");
   
@@ -37,8 +39,8 @@ int main(int argc, char *argv[])
   auto nev = tMc.GetEntries();
   // // auto ev_min = std::atoi(argv[1]);
   // // auto ev_max = std::atoi(argv[2]);
-  auto ev_min = 2;
-  auto ev_max = 3;
+  auto ev_min = 0;
+  auto ev_max = 250;
 
   // loop over events
 
@@ -55,36 +57,45 @@ int main(int argc, char *argv[])
     std::string print_string = "";
     tree.Print(print_string, 1);
   
-    auto start_point = tree.begin()->GetHitMap().rbegin()->second.rbegin()->GetStart();
-    auto stop_point = tree.begin()->GetHitMap().rbegin()->second.rbegin()->GetStop();
-    auto middle_point = (start_point + stop_point) * 0.5;
-    double distance = 10;
-    std::cout << tree.begin()->GetId() << std::endl;
-    std::vector<EDEPTrajectory> filteredTrj;
-    tree.Filter(std::back_insert_iterator<std::vector<EDEPTrajectory>>(filteredTrj), 
-      [middle_point, distance](const EDEPTrajectory& trj) { return trj.HasHitNearPoint(middle_point.Vect(), distance);} );
-    std::cout << "Number of filtered traj: " << filteredTrj.size() << std::endl;
-    for (auto& trj:filteredTrj) {
-      std::cout << trj.GetId() << std::endl;
-    }
-    filteredTrj.clear();
-    tree.Filter(std::back_insert_iterator<std::vector<EDEPTrajectory>>(filteredTrj), 
-      [middle_point, distance](const EDEPTrajectory& trj) { return trj.HasHitNear4DPoint(middle_point, distance, middle_point.T());} );
-    std::cout << "Number of filtered traj: " << filteredTrj.size() << std::endl;
-    for (auto& trj:filteredTrj) {
-      std::cout << trj.GetHitNear4DPoint(middle_point, distance, middle_point.T())->GetId() << std::endl;
-    }
+    // auto start_point = tree.begin()->GetHitMap().rbegin()->second.rbegin()->GetStart();
+    // auto stop_point = tree.begin()->GetHitMap().rbegin()->second.rbegin()->GetStop();
+    // auto middle_point = (start_point + stop_point) * 0.5;
+    // double distance = 10;
+    // std::cout << tree.begin()->GetId() << std::endl;
+    // std::vector<EDEPTrajectory> filteredTrj;
+    // tree.Filter(std::back_insert_iterator<std::vector<EDEPTrajectory>>(filteredTrj), 
+    //   [middle_point, distance](const EDEPTrajectory& trj) { return trj.HasHitNearPoint(middle_point.Vect(), distance);} );
+    // std::cout << "Number of filtered traj: " << filteredTrj.size() << std::endl;
+    // for (auto& trj:filteredTrj) {
+    //   std::cout << trj.GetId() << std::endl;
+    // }
+    // filteredTrj.clear();
+    // tree.Filter(std::back_insert_iterator<std::vector<EDEPTrajectory>>(filteredTrj), 
+    //   [middle_point, distance](const EDEPTrajectory& trj) { return trj.HasHitNear4DPoint(middle_point, distance, middle_point.T());} );
+    // std::cout << "Number of filtered traj: " << filteredTrj.size() << std::endl;
+    // for (auto& trj:filteredTrj) {
+    //   std::cout << trj.GetHitNear4DPoint(middle_point, distance, middle_point.T())->GetId() << std::endl;
+    // }
 
     // std::cout << "########### Printing the tree" << std::endl;
     // tree.Print(print_string);
 
-    // // LArHit, Straw, EMCalSci
-    // component comp = component::ECAL;
+    // // LArHit, Straw, EMCalSci, DriftVolume
+    // component comp = component::DRIFT;
     // std::cout << "########### Filtering all the trajectories Entering in " << component_to_string.at(comp) << std::endl;
     // std::vector<EDEPTrajectory> filteredTrj;
     // tree.Filter(std::back_insert_iterator<std::vector<EDEPTrajectory>>(filteredTrj), 
-    //   [comp](const EDEPTrajectory& trj) { return trj.IsEntering(comp);} );
+    //   [comp](const EDEPTrajectory& trj) { return trj.IsExiting(comp);} );
     // std::cout << "Number of filtered traj: " << filteredTrj.size() << std::endl;
+
+    // for (auto trj:filteredTrj) {
+    //   for (auto p:trj.GetFirstPointsInDetector(comp)) {
+    //     p.GetPosition().Print();
+    //   }
+    //   for (auto p:trj.GetLastPointsInDetector(comp)) {
+    //     p.GetPosition().Print();
+    //   }
+    // }
     
     
     // filteredTrj.clear();