Classes | Public Member Functions | Private Member Functions | Private Attributes | List of all members
trkf::TrackAnaCT Class Reference
Inheritance diagram for trkf::TrackAnaCT:
art::EDAnalyzer art::detail::Analyzer art::detail::LegacyModule art::Observer art::ModuleBase

Classes

struct  MCHists
 
struct  RecoHists
 

Public Member Functions

 TrackAnaCT (fhicl::ParameterSet const &pset)
 
virtual ~TrackAnaCT ()
 
void analyze (const art::Event &evt)
 
void anaStitch (detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp, const art::Event &evt)
 
void endJob ()
 
- Public Member Functions inherited from art::EDAnalyzer
 EDAnalyzer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 EDAnalyzer (Table< Config > const &config)
 
std::string workerType () const
 
- Public Member Functions inherited from art::detail::Analyzer
virtual ~Analyzer () noexcept
 
 Analyzer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 Analyzer (Table< Config > const &config)
 
void doBeginJob (SharedResources const &resources)
 
void doEndJob ()
 
void doRespondToOpenInputFile (FileBlock const &fb)
 
void doRespondToCloseInputFile (FileBlock const &fb)
 
void doRespondToOpenOutputFiles (FileBlock const &fb)
 
void doRespondToCloseOutputFiles (FileBlock const &fb)
 
bool doBeginRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doEndRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doBeginSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEndSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEvent (EventPrincipal &ep, ModuleContext const &mc, std::atomic< std::size_t > &counts_run, std::atomic< std::size_t > &counts_passed, std::atomic< std::size_t > &counts_failed)
 
- Public Member Functions inherited from art::Observer
 ~Observer () noexcept
 
 Observer (Observer const &)=delete
 
 Observer (Observer &&)=delete
 
Observer & operator= (Observer const &)=delete
 
Observer & operator= (Observer &&)=delete
 
void registerProducts (ProductDescriptions &, ModuleDescription const &)
 
void fillDescriptions (ModuleDescription const &)
 
fhicl::ParameterSetID selectorConfig () const
 
- Public Member Functions inherited from art::ModuleBase
virtual ~ModuleBase () noexcept
 
 ModuleBase ()
 
ModuleDescription const & moduleDescription () const
 
void setModuleDescription (ModuleDescription const &)
 
std::array< std::vector< ProductInfo >, NumBranchTypes > const & getConsumables () const
 
void sortConsumables (std::string const &current_process_name)
 
template<typename T , BranchType BT>
ViewToken< T > consumesView (InputTag const &tag)
 
template<typename T , BranchType BT>
ViewToken< T > mayConsumeView (InputTag const &tag)
 

Private Member Functions

template<typename T >
std::vector< size_t > fsort_indexes (const std::vector< T > &v)
 

Private Attributes

std::string fTrackModuleLabel
 
std::string fSpacepointModuleLabel
 
std::string fStitchModuleLabel
 
std::string fTrkSpptAssocModuleLabel
 
std::string fHitSpptAssocModuleLabel
 
int fDump
 
double fMinMCKE
 
double fMinMCLen
 
double fMatchColinearity
 
double fMatchDisp
 
double fWMatchDisp
 
bool fIgnoreSign
 
bool fStitchedAnalysis
 
std::map< int, MCHists > fMCHistMap
 
std::map< int, RecoHists > fRecoHistMap
 
int fNumEvent
 

Additional Inherited Members

- Public Types inherited from art::EDAnalyzer
using WorkerType = WorkerT< EDAnalyzer >
 
using ModuleType = EDAnalyzer
 
- Protected Member Functions inherited from art::Observer
std::string const & processName () const
 
bool wantAllEvents () const noexcept
 
bool wantEvent (ScheduleID id, Event const &e) const
 
Handle< TriggerResults > getTriggerResults (Event const &e) const
 
 Observer (fhicl::ParameterSet const &config)
 
 Observer (std::vector< std::string > const &select_paths, std::vector< std::string > const &reject_paths, fhicl::ParameterSet const &config)
 
- Protected Member Functions inherited from art::ModuleBase
ConsumesCollector & consumesCollector ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > consumes (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > consumesView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void consumesMany ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > mayConsume (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > mayConsumeView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void mayConsumeMany ()
 

Detailed Description

Definition at line 247 of file TrackAnaCT_module.cc.

Constructor & Destructor Documentation

trkf::TrackAnaCT::TrackAnaCT ( fhicl::ParameterSet const &  pset)
explicit
Todo:
Move this module to DUNE code and remove it from larreco

Definition at line 691 of file TrackAnaCT_module.cc.

697  : EDAnalyzer(pset)
698  , fTrackModuleLabel(pset.get<std::string>("TrackModuleLabel"))
699  , fSpacepointModuleLabel(pset.get<std::string>("SpacepointModuleLabel"))
700  , fStitchModuleLabel(pset.get<std::string>("StitchModuleLabel"))
701  , fTrkSpptAssocModuleLabel(pset.get<std::string>("TrkSpptAssocModuleLabel"))
702  , fHitSpptAssocModuleLabel(pset.get<std::string>("HitSpptAssocModuleLabel"))
703  , fDump(pset.get<int>("Dump"))
704  , fMinMCKE(pset.get<double>("MinMCKE"))
705  , fMinMCLen(pset.get<double>("MinMCLen"))
706  , fMatchColinearity(pset.get<double>("MatchColinearity"))
707  , fMatchDisp(pset.get<double>("MatchDisp"))
708  , fWMatchDisp(pset.get<double>("WMatchDisp"))
709  , fIgnoreSign(pset.get<bool>("IgnoreSign"))
710  , fStitchedAnalysis(pset.get<bool>("StitchedAnalysis",false))
711  , fNumEvent(0)
712  {
713 
714  ///\todo Move this module to DUNE code and remove it from larreco
715 
716  art::ServiceHandle<geo::Geometry> geom;
717  if(geom->DetectorName().find("dune") == std::string::npos)
718  throw cet::exception("TrackAnaCT") << "TrackAnaCT should only be used with DUNE "
719  << "geometries, the name for this detector, "
720  << geom->DetectorName() << ", does not contain "
721  << "dune, bail.";
722 
723  // Report.
724 
725  mf::LogInfo("TrackAnaCT")
726  << "TrackAnaCT configured with the following parameters:\n"
727  << " TrackModuleLabel = " << fTrackModuleLabel << "\n"
728  << " StitchModuleLabel = " << fStitchModuleLabel << "\n"
729  << " TrkSpptAssocModuleLabel = " << fTrkSpptAssocModuleLabel << "\n"
730  << " HitSpptAssocModuleLabel = " << fHitSpptAssocModuleLabel << "\n"
731  << " Dump = " << fDump << "\n"
732  << " MinMCKE = " << fMinMCKE << "\n"
733  << " MinMCLen = " << fMinMCLen;
734  }
string
std::string string
Definition: nybbler.cc:12
mf::LogInfo
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
Definition: MessageLogger.h:209
art::EDAnalyzer::EDAnalyzer
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.h:25
trkf::TrackAnaCT::fStitchModuleLabel
std::string fStitchModuleLabel
Definition: TrackAnaCT_module.cc:400
geo::GeometryCore::DetectorName
std::string DetectorName() const
Returns a string with the name of the detector, as configured.
Definition: GeometryCore.h:1761
trkf::TrackAnaCT::fHitSpptAssocModuleLabel
std::string fHitSpptAssocModuleLabel
Definition: TrackAnaCT_module.cc:402
trkf::TrackAnaCT::fTrackModuleLabel
std::string fTrackModuleLabel
Definition: TrackAnaCT_module.cc:398
trkf::TrackAnaCT::fTrkSpptAssocModuleLabel
std::string fTrkSpptAssocModuleLabel
Definition: TrackAnaCT_module.cc:401
trkf::TrackAnaCT::fSpacepointModuleLabel
std::string fSpacepointModuleLabel
Definition: TrackAnaCT_module.cc:399
fhicl::exception
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
trkf::TrackAnaCT::~TrackAnaCT ( )
virtual

Definition at line 736 of file TrackAnaCT_module.cc.

740  {}

Member Function Documentation

void trkf::TrackAnaCT::analyze ( const art::Event &  evt)

figuring out which TPC

Definition at line 742 of file TrackAnaCT_module.cc.

748  {
749  art::ServiceHandle<geo::Geometry> geom;
750 
751 
752  ++fNumEvent;
753 
754  // Optional dump stream.
755 
756  std::unique_ptr<mf::LogInfo> pdump;
757  if(fDump > 0) {
758  --fDump;
759  pdump = std::unique_ptr<mf::LogInfo>(new mf::LogInfo("TrackAnaCT"));
760  }
761 
762  // Make sure histograms are booked.
763 
764  bool mc = !evt.isRealData();
765 
766  // Get mc particles.
767 
768  auto const clockData = art::ServiceHandle<detinfo::DetectorClocksService const>()->DataFor(evt);
769  auto const detProp = art::ServiceHandle<detinfo::DetectorPropertiesService const>()->DataFor(evt, clockData);
770 
771  std::vector<const simb::MCParticle*> plist2;
772  if(mc) {
773 
774 // art::ServiceHandle<cheat::BackTrackerService> bt_serv;
775  art::ServiceHandle<cheat::ParticleInventoryService> pi_serv;
776  sim::ParticleList const& plist = pi_serv->ParticleList();
777  plist2.reserve(plist.size());
778 
779 
780  if(pdump) {
781  *pdump << "MC Particles\n"
782  << " Id pdg x y z dx dy dz p\n"
783  << "-------------------------------------------------------------------------------------------\n";
784  }
785 
786  // Loop over mc particles, and fill histograms that depend only
787  // on mc particles. Also, fill a secondary list of mc particles
788  // that pass various selection criteria.
789 
790  for(sim::ParticleList::const_iterator ipart = plist.begin();
791  ipart != plist.end(); ++ipart) {
792  const simb::MCParticle* part = (*ipart).second;
793  if (!part)
794  throw cet::exception("SeedAna") << "no particle!\n";
795  int pdg = part->PdgCode();
796  if(fIgnoreSign)
797  pdg = std::abs(pdg);
798 
799  // Ignore everything except stable charged nonshowering particles.
800 
801  int apdg = std::abs(pdg);
802  if(apdg == 13 || // Muon
803  apdg == 211 || // Charged pion
804  apdg == 321 || // Charged kaon
805  apdg == 2212) { // (Anti)proton
806 
807  // Apply minimum energy cut.
808 
809  if(part->E() >= 0.001*part->Mass() + fMinMCKE) {
810 
811  // Calculate the x offset due to nonzero mc particle time.
812 
813  double mctime = part->T(); // nsec
814  double mcdx = mctime * 1.e-3 * detProp.DriftVelocity(); // cm
815 
816  // Calculate the length of this mc particle inside the fiducial volume.
817 
818  TVector3 mcstart;
819  TVector3 mcend;
820  TVector3 mcstartmom;
821  TVector3 mcendmom;
822  double plen = length(detProp, *part, mcdx, mcstart, mcend, mcstartmom, mcendmom);
823 
824  // Apply minimum fiducial length cut. Always reject particles that have
825  // zero length in the tpc regardless of the configured cut.
826 
827  if(plen > 0. && plen > fMinMCLen) {
828 
829  // This is a good mc particle (capable of making a track).
830 
831  plist2.push_back(part);
832 
833  // Dump MC particle information here.
834 
835  if(pdump) {
836  double pstart = mcstartmom.Mag();
837  double pend = mcendmom.Mag();
838  *pdump << "\nOffset"
839  << std::setw(3) << part->TrackId()
840  << std::setw(6) << part->PdgCode()
841  << " "
842  << std::fixed << std::setprecision(2)
843  << std::setw(10) << mcdx
844  << "\nStart "
845  << std::setw(3) << part->TrackId()
846  << std::setw(6) << part->PdgCode()
847  << " "
848  << std::fixed << std::setprecision(2)
849  << std::setw(10) << mcstart[0]
850  << std::setw(10) << mcstart[1]
851  << std::setw(10) << mcstart[2];
852  if(pstart > 0.) {
853  *pdump << " "
854  << std::fixed << std::setprecision(3)
855  << std::setw(10) << mcstartmom[0]/pstart
856  << std::setw(10) << mcstartmom[1]/pstart
857  << std::setw(10) << mcstartmom[2]/pstart;
858  }
859  else
860  *pdump << std::setw(32) << " ";
861  *pdump << std::setw(12) << std::fixed << std::setprecision(3) << pstart;
862  *pdump << "\nEnd "
863  << std::setw(3) << part->TrackId()
864  << std::setw(6) << part->PdgCode()
865  << " "
866  << std::fixed << std::setprecision(2)
867  << std::setw(10) << mcend[0]
868  << std::setw(10) << mcend[1]
869  << std::setw(10) << mcend[2];
870  if(pend > 0.01) {
871  *pdump << " "
872  << std::fixed << std::setprecision(3)
873  << std::setw(10) << mcendmom[0]/pend
874  << std::setw(10) << mcendmom[1]/pend
875  << std::setw(10) << mcendmom[2]/pend;
876  }
877  else
878  *pdump << std::setw(32)<< " ";
879  *pdump << std::setw(12) << std::fixed << std::setprecision(3) << pend << "\n";
880  }
881 
882  // Fill histograms.
883 
884  if(fMCHistMap.count(pdg) == 0) {
885  std::ostringstream ostr;
886  ostr << "MC" << (fIgnoreSign ? "All" : (pdg > 0 ? "Pos" : "Neg")) << std::abs(pdg);
887  fMCHistMap[pdg] = MCHists(ostr.str());
888  }
889  const MCHists& mchists = fMCHistMap[pdg];
890 
891  double mctheta_xz = std::atan2(mcstartmom.X(), mcstartmom.Z());
892  double mctheta_yz = std::atan2(mcstartmom.Y(), mcstartmom.Z());
893 
894  mchists.fHmcstartx->Fill(mcstart.X());
895  mchists.fHmcstarty->Fill(mcstart.Y());
896  mchists.fHmcstartz->Fill(mcstart.Z());
897  mchists.fHmcendx->Fill(mcend.X());
898  mchists.fHmcendy->Fill(mcend.Y());
899  mchists.fHmcendz->Fill(mcend.Z());
900  mchists.fHmctheta->Fill(mcstartmom.Theta());
901  mchists.fHmcphi->Fill(mcstartmom.Phi());
902  mchists.fHmctheta_xz->Fill(mctheta_xz);
903  mchists.fHmctheta_yz->Fill(mctheta_yz);
904  mchists.fHmcmom->Fill(mcstartmom.Mag());
905  mchists.fHmclen->Fill(plen);
906  }
907  }
908  }
909  }
910  } //mc
911 
912 
913  // Get tracks and spacepoints and hits
914 
915  auto trackh = evt.getHandle< std::vector<recob::Track> >(fTrackModuleLabel);
916  auto trackvh = evt.getHandle< std::vector< art::PtrVector < recob::Track > > >(fStitchModuleLabel);
917 
918  // This new top part of TrackAnaCT between two long lines of ************s
919  // is particular to analyzing Stitched Tracks.
920  // *******************************************************************//
921 
922  if (trackvh.isValid() && fStitchedAnalysis)
923  {
924  mf::LogDebug("TrackAnaCT")
925  << "TrackAnaCT read " << trackvh->size()
926  << " vectors of Stitched PtrVectorsof tracks.";
927  anaStitch(clockData, detProp, evt);
928  }
929 
930  if(trackh.isValid()) {
931 
932  if(pdump) {
933  *pdump << "\nReconstructed Tracks\n"
934  << " Id MCid x y z dx dy dz p\n"
935  << "-------------------------------------------------------------------------------------------\n";
936  }
937 
938  // Loop over tracks.
939 
940 
941  int ntrack = trackh->size();
942  for(int i = 0; i < ntrack; ++i) {
943  art::Ptr<recob::Track> ptrack(trackh, i);
944  const recob::Track& track = *ptrack;
945  art::FindManyP<recob::Hit> fh(trackh, evt, fTrkSpptAssocModuleLabel);
946 
947  ////
948  /// figuring out which TPC
949  ///
950  ///
951  //
952  // auto pcoll { ptrack };
953  //art::FindManyP<recob::SpacePoint> fs( pcoll, evt, fTrkSpptAssocModuleLabel);
954  // auto sppt = fs.at(0);//.at(is);
955  // art::FindManyP<recob::Hit> fh( sppt, evt, fHitSpptAssocModuleLabel);
956  auto hit = fh.at(0).at(0);
957  geo::WireID tmpWireid=hit->WireID();
958  int hit_tpc=tmpWireid.TPC;
959  ///
960  ///
961  //
962  //
963  //
964  //
965  //
966  //
967 
968 
969 
970 
971  /*
972  std::vector<art::Ptr<recob::Hit> > hitlist;
973  auto hitListHandle = evt.getHandle< std::vector<recob::Hit> >(fHitsModuleLabel);
974  if (hitListHandle) art::fill_ptr_vector(hitlist, hitListHandle);
975  */
976 
977  // Calculate the x offset due to nonzero reconstructed time.
978 
979  //double recotime = track.Time() * sampling_rate(clockData); // nsec
980  double recotime = 0.;
981  double trackdx = recotime * 1.e-3 * detProp.DriftVelocity(); // cm
982 
983  // Fill histograms involving reco tracks only.
984 
985  int ntraj = track.NumberTrajectoryPoints();
986  if(ntraj > 0) {
987  TVector3 pos = track.Vertex<TVector3>();
988  TVector3 dir = track.VertexDirection<TVector3>();
989  TVector3 end = track.End<TVector3>();
990  pos[0] += trackdx;
991  end[0] += trackdx;
992 
993  double dpos = bdist(pos,hit_tpc);
994  double dend = bdist(end,hit_tpc);
995  double tlen = length(track);
996  double theta_xz = std::atan2(dir.X(), dir.Z());
997  double theta_yz = std::atan2(dir.Y(), dir.Z());
998 
999  if(fRecoHistMap.count(0) == 0)
1000  fRecoHistMap[0] = RecoHists("Reco");
1001  const RecoHists& rhists = fRecoHistMap[0];
1002 
1003  rhists.fHstartx->Fill(pos.X());
1004  rhists.fHstarty->Fill(pos.Y());
1005  rhists.fHstartz->Fill(pos.Z());
1006  rhists.fHstartd->Fill(dpos);
1007  rhists.fHendx->Fill(end.X());
1008  rhists.fHendy->Fill(end.Y());
1009  rhists.fHendz->Fill(end.Z());
1010  rhists.fHendd->Fill(dend);
1011  rhists.fHtheta->Fill(dir.Theta());
1012  rhists.fHphi->Fill(dir.Phi());
1013  rhists.fHtheta_xz->Fill(theta_xz);
1014  rhists.fHtheta_yz->Fill(theta_yz);
1015 
1016  double mom = 0.;
1017  if(track.NumberTrajectoryPoints() > 0)
1018  mom = track.VertexMomentum();
1019  rhists.fHmom->Fill(mom);
1020  rhists.fHlen->Fill(tlen);
1021 
1022  // Id of matching mc particle.
1023 
1024  int mcid = -1;
1025 
1026  // Loop over direction.
1027 
1028  for(int swap=0; swap<2; ++swap) {
1029 
1030  // Analyze reversed tracks only if start momentum = end momentum.
1031 
1032  if(swap != 0 && track.NumberTrajectoryPoints() > 0 &&
1033  std::abs(track.VertexMomentum() - track.EndMomentum()) > 1.e-3)
1034  continue;
1035 
1036  // Calculate the global-to-local rotation matrix.
1037 
1038  int start_point = (swap == 0 ? 0 : ntraj-1);
1039  TMatrixD rot = track.GlobalToLocalRotationAtPoint<TMatrixD>(start_point);
1040 
1041  // Update track data for reversed case.
1042 
1043  if(swap != 0) {
1044  rot(1, 0) = -rot(1, 0);
1045  rot(2, 0) = -rot(2, 0);
1046  rot(1, 1) = -rot(1, 1);
1047  rot(2, 1) = -rot(2, 1);
1048  rot(1, 2) = -rot(1, 2);
1049  rot(2, 2) = -rot(2, 2);
1050 
1051  pos = track.End<TVector3>();
1052  dir = -track.EndDirection<TVector3>();
1053  end = track.Vertex<TVector3>();
1054  pos[0] += trackdx;
1055  end[0] += trackdx;
1056 
1057  dpos = bdist(pos,hit_tpc);
1058  dend = bdist(end,hit_tpc);
1059  theta_xz = std::atan2(dir.X(), dir.Z());
1060  theta_yz = std::atan2(dir.Y(), dir.Z());
1061 
1062  if(track.NumberTrajectoryPoints() > 0)
1063  mom = track.EndMomentum();
1064  }
1065 
1066  // Get covariance matrix.
1067 
1068  // const TMatrixD& cov = (swap == 0 ? track.VertexCovariance() : track.EndCovariance());
1069 
1070  // Loop over track-like mc particles.
1071 
1072  for(auto ipart = plist2.begin(); ipart != plist2.end(); ++ipart) {
1073  const simb::MCParticle* part = *ipart;
1074  if (!part)
1075  throw cet::exception("SeedAna") << "no particle! [II]\n";
1076  int pdg = part->PdgCode();
1077  if(fIgnoreSign) pdg = std::abs(pdg);
1078  auto iMCHistMap = fMCHistMap.find(pdg);
1079  if (iMCHistMap == fMCHistMap.end())
1080  throw cet::exception("SeedAna") << "no particle with ID=" << pdg << "\n";
1081  const MCHists& mchists = iMCHistMap->second;
1082 
1083  // Calculate the x offset due to nonzero mc particle time.
1084 
1085  double mctime = part->T(); // nsec
1086  double mcdx = mctime * 1.e-3 * detProp.DriftVelocity(); // cm
1087 
1088  // Calculate the points where this mc particle enters and leaves the
1089  // fiducial volume, and the length in the fiducial volume.
1090 
1091  TVector3 mcstart;
1092  TVector3 mcend;
1093  TVector3 mcstartmom;
1094  TVector3 mcendmom;
1095  double plen = length(detProp, *part, mcdx, mcstart, mcend, mcstartmom, mcendmom);
1096 
1097  // Get the displacement of this mc particle in the global coordinate system.
1098 
1099  TVector3 mcpos = mcstart - pos;
1100 
1101  // Rotate the momentum and position to the
1102  // track-local coordinate system.
1103 
1104  TVector3 mcmoml = rot * mcstartmom;
1105  TVector3 mcposl = rot * mcpos;
1106 
1107  double colinearity = mcmoml.Z() / mcmoml.Mag();
1108 
1109  double u = mcposl.X();
1110  double v = mcposl.Y();
1111  double w = mcposl.Z();
1112 
1113  double pu = mcmoml.X();
1114  double pv = mcmoml.Y();
1115  double pw = mcmoml.Z();
1116 
1117  double dudw = pu / pw;
1118  double dvdw = pv / pw;
1119 
1120  double u0 = u - w * dudw;
1121  double v0 = v - w * dvdw;
1122  double uv0 = std::sqrt(u0*u0 + v0*v0);
1123 
1124  mchists.fHduvcosth->Fill(colinearity, uv0);
1125  if(std::abs(uv0) < fMatchDisp) {
1126 
1127  // Fill slope matching histograms.
1128 
1129  mchists.fHmcdudw->Fill(dudw);
1130  mchists.fHmcdvdw->Fill(dvdw);
1131  // mchists.fHdudwpull->Fill(dudw / std::sqrt(cov(2,2)));
1132  // mchists.fHdvdwpull->Fill(dvdw / std::sqrt(cov(3,3)));
1133  }
1134  mchists.fHcosth->Fill(colinearity);
1135  if(colinearity > fMatchColinearity) {
1136 
1137  // Fill displacement matching histograms.
1138 
1139  mchists.fHmcu->Fill(u0);
1140  mchists.fHmcv->Fill(v0);
1141  mchists.fHmcw->Fill(w);
1142  // mchists.fHupull->Fill(u0 / std::sqrt(cov(0,0)));
1143  // mchists.fHvpull->Fill(v0 / std::sqrt(cov(1,1)));
1144 
1145  if(std::abs(uv0) < fMatchDisp) {
1146 
1147  // Fill matching histograms.
1148 
1149  double mctheta_xz = std::atan2(mcstartmom.X(), mcstartmom.Z());
1150  double mctheta_yz = std::atan2(mcstartmom.Y(), mcstartmom.Z());
1151 
1152  mchists.fHstartdx->Fill(pos.X() - mcstart.X());
1153  mchists.fHstartdy->Fill(pos.Y() - mcstart.Y());
1154  mchists.fHstartdz->Fill(pos.Z() - mcstart.Z());
1155  mchists.fHenddx->Fill(end.X() - mcend.X());
1156  mchists.fHenddy->Fill(end.Y() - mcend.Y());
1157  mchists.fHenddz->Fill(end.Z() - mcend.Z());
1158  mchists.fHlvsl->Fill(plen, tlen);
1159  mchists.fHdl->Fill(tlen - plen);
1160  mchists.fHpvsp->Fill(mcstartmom.Mag(), mom);
1161  double dp = mom - mcstartmom.Mag();
1162  mchists.fHdp->Fill(dp);
1163  // mchists.fHppull->Fill(dp / std::sqrt(cov(4,4)));
1164  if(std::abs(dpos) >= 5. && std::abs(dend) >= 5.) {
1165  mchists.fHpvspc->Fill(mcstartmom.Mag(), mom);
1166  mchists.fHdpc->Fill(dp);
1167  // mchists.fHppullc->Fill(dp / std::sqrt(cov(4,4)));
1168  }
1169 
1170  // Count this track as well-reconstructed if it is matched to an
1171  // mc particle (which it is if get here), and if in addition the
1172  // starting position (w) matches and the reconstructed track length
1173  // is more than 0.5 of the mc particle trajectory length.
1174 
1175  bool good = std::abs(w) <= fWMatchDisp &&
1176  tlen > 0.5 * plen;
1177  if(good) {
1178  mcid = part->TrackId();
1179  mchists.fHgstartx->Fill(mcstart.X());
1180  mchists.fHgstarty->Fill(mcstart.Y());
1181  mchists.fHgstartz->Fill(mcstart.Z());
1182  mchists.fHgendx->Fill(mcend.X());
1183  mchists.fHgendy->Fill(mcend.Y());
1184  mchists.fHgendz->Fill(mcend.Z());
1185  mchists.fHgtheta->Fill(mcstartmom.Theta());
1186  mchists.fHgphi->Fill(mcstartmom.Phi());
1187  mchists.fHgtheta_xz->Fill(mctheta_xz);
1188  mchists.fHgtheta_yz->Fill(mctheta_yz);
1189  mchists.fHgmom->Fill(mcstartmom.Mag());
1190  mchists.fHglen->Fill(plen);
1191  }
1192  }
1193  }
1194  }
1195  }
1196 
1197  // Dump track information here.
1198 
1199  if(pdump) {
1200  TVector3 pos = track.Vertex<TVector3>();
1201  TVector3 dir = track.VertexDirection<TVector3>();
1202  TVector3 end = track.End<TVector3>();
1203  pos[0] += trackdx;
1204  end[0] += trackdx;
1205  TVector3 enddir = track.EndDirection<TVector3>();
1206  double pstart = track.VertexMomentum();
1207  double pend = track.EndMomentum();
1208  *pdump << "\nOffset"
1209  << std::setw(3) << track.ID()
1210  << std::setw(6) << mcid
1211  << " "
1212  << std::fixed << std::setprecision(2)
1213  << std::setw(10) << trackdx
1214  << "\nStart "
1215  << std::setw(3) << track.ID()
1216  << std::setw(6) << mcid
1217  << " "
1218  << std::fixed << std::setprecision(2)
1219  << std::setw(10) << pos[0]
1220  << std::setw(10) << pos[1]
1221  << std::setw(10) << pos[2];
1222  if(pstart > 0.) {
1223  *pdump << " "
1224  << std::fixed << std::setprecision(3)
1225  << std::setw(10) << dir[0]
1226  << std::setw(10) << dir[1]
1227  << std::setw(10) << dir[2];
1228  }
1229  else
1230  *pdump << std::setw(32) << " ";
1231  *pdump << std::setw(12) << std::fixed << std::setprecision(3) << pstart;
1232  *pdump << "\nEnd "
1233  << std::setw(3) << track.ID()
1234  << std::setw(6) << mcid
1235  << " "
1236  << std::fixed << std::setprecision(2)
1237  << std::setw(10) << end[0]
1238  << std::setw(10) << end[1]
1239  << std::setw(10) << end[2];
1240  if(pend > 0.01) {
1241  *pdump << " "
1242  << std::fixed << std::setprecision(3)
1243  << std::setw(10) << enddir[0]
1244  << std::setw(10) << enddir[1]
1245  << std::setw(10) << enddir[2];
1246  }
1247  else
1248  *pdump << std::setw(32)<< " ";
1249  *pdump << std::setw(12) << std::fixed << std::setprecision(3) << pend << "\n";
1250  }
1251  }
1252  }
1253  } // i
1254 
1255  }
simb::MCParticle::E
double E(const int i=0) const
Definition: MCParticle.h:233
ValidateOpDetReco.end
end
while True: pbar.update(maxval-len(onlies[E][S])) #print iS, "/", len(onlies[E][S]) found = False for...
Definition: ValidateOpDetReco.py:548
recob::Track::EndMomentum
double EndMomentum() const
Definition: Track.h:144
simb::MCParticle::PdgCode
int PdgCode() const
Definition: MCParticle.h:212
recob::Track::VertexMomentum
double VertexMomentum() const
Definition: Track.h:142
trkf::TrackAnaCT::anaStitch
void anaStitch(detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp, const art::Event &evt)
Definition: TrackAnaCT_module.cc:1257
art::DataViewImpl::getHandle
Handle< PROD > getHandle(SelectorBase const &) const
Definition: DataViewImpl.h:382
mf::LogInfo
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
Definition: MessageLogger.h:209
simb::MCParticle::Mass
double Mass() const
Definition: MCParticle.h:239
recob::Track::NumberTrajectoryPoints
size_t NumberTrajectoryPoints() const
Various functions related to the presence and the number of (valid) points.
Definition: Track.h:102
recob::Track::VertexDirection
Vector_t VertexDirection() const
Definition: Track.h:132
const_iterator
intermediate_table::const_iterator const_iterator
Definition: intermediate_table.cc:28
dir
string dir
Definition: INukeNucleonCorr.cxx:61
setprecision
Q_EXPORT QTSManip setprecision(int p)
Definition: qtextstream.h:343
simb::MCParticle::TrackId
int TrackId() const
Definition: MCParticle.h:210
recob::Track::GlobalToLocalRotationAtPoint
Rotation_t GlobalToLocalRotationAtPoint(size_t p) const
Definition: Track.h:190
trkf::TrackAnaCT::fStitchModuleLabel
std::string fStitchModuleLabel
Definition: TrackAnaCT_module.cc:400
art::DataViewImpl::isRealData
bool isRealData() const
Definition: DataViewImpl.cc:113
abs
T abs(T value)
Definition: sparse_vector_test.cc:30
e
const double e
Definition: gUpMuFluxGen.cxx:165
art::swap
void swap(Handle< T > &a, Handle< T > &b)
recob::Track::Vertex
Point_t const & Vertex() const
Definition: Track.h:124
simb::MCParticle::T
double T(const int i=0) const
Definition: MCParticle.h:224
trkf::TrackAnaCT::fTrackModuleLabel
std::string fTrackModuleLabel
Definition: TrackAnaCT_module.cc:398
hit
Detector simulation of raw signals on wires.
Definition: NumberOfHitsFilter_module.cc:25
setw
Q_EXPORT QTSManip setw(int w)
Definition: qtextstream.h:331
cheat::ParticleInventoryService::ParticleList
const sim::ParticleList & ParticleList() const
Definition: ParticleInventoryService_service.cc:114
recob::Track::ID
int ID() const
Definition: Track.h:198
mf::LogDebug
MaybeLogger_< ELseverityLevel::ELsev_success, false > LogDebug
Definition: MessageLogger.h:208
recob::Track::EndDirection
Vector_t EndDirection() const
Definition: Track.h:133
trkf::TrackAnaCT::fTrkSpptAssocModuleLabel
std::string fTrkSpptAssocModuleLabel
Definition: TrackAnaCT_module.cc:401
trkf::TrackAnaCT::fMCHistMap
std::map< int, MCHists > fMCHistMap
Definition: TrackAnaCT_module.cc:415
run_cnn_1class.pdg
pdg
Definition: run_cnn_1class.py:54
recob::Track::End
Point_t const & End() const
Definition: Track.h:125
geo::TPCID::TPC
TPCID_t TPC
Index of the TPC within its cryostat.
Definition: geo_types.h:406
trkf::TrackAnaCT::fRecoHistMap
std::map< int, RecoHists > fRecoHistMap
Definition: TrackAnaCT_module.cc:416
recob::Track
Track from a non-cascading particle.A recob::Track consists of a recob::TrackTrajectory, plus additional members relevant for a "fitted" track:
Definition: Track.h:49
fhicl::exception
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
void trkf::TrackAnaCT::anaStitch ( detinfo::DetectorClocksData const &  clockData,
detinfo::DetectorPropertiesData const &  detProp,
const art::Event &  evt 
)

art::FindManyP<recob::Hit> fh(sppth, evt, fHitSpptAssocModuleLabel);

Definition at line 1257 of file TrackAnaCT_module.cc.

1260  {
1261 
1262  art::ServiceHandle<cheat::BackTrackerService> bt_serv;
1263  art::ServiceHandle<cheat::ParticleInventoryService> pi_serv;
1264  art::ServiceHandle<geo::Geometry> geom;
1265 
1266  std::map<int, std::map<int, art::PtrVector<recob::Hit>> > hitmap; // trkID, otrk, hitvec
1267  std::map<int, int > KEmap; // length traveled in det [cm]?, trkID want to sort by KE
1268  bool mc = !evt.isRealData();
1269 
1270  auto trackh = evt.getHandle< std::vector<recob::Track> >(fTrackModuleLabel);
1271  auto sppth = evt.getHandle< std::vector< recob::SpacePoint> >(fSpacepointModuleLabel);
1272  auto trackvh = evt.getHandle< std::vector< art::PtrVector < recob::Track > > >(fStitchModuleLabel);
1273 
1274  int ntv(trackvh->size());
1275 
1276  std::vector < art::PtrVector<recob::Track> >::const_iterator cti = trackvh->begin();
1277  /// art::FindManyP<recob::Hit> fh(sppth, evt, fHitSpptAssocModuleLabel);
1278 
1279  if(trackh.isValid()) {
1280  art::FindManyP<recob::SpacePoint> fswhole(trackh, evt, fTrkSpptAssocModuleLabel);
1281  int nsppts_assnwhole = fswhole.size();
1282  std::cout << "TrackAnaCT: Number of clumps of Spacepoints from Assn for all Tracks: " << nsppts_assnwhole << std::endl;
1283  }
1284 
1285  if(fRecoHistMap.count(0) == 0)
1286  {
1287  fRecoHistMap[0] = RecoHists("Reco");
1288  std::cout << "\n" << "\t\t TrkAna: Fresh fRecoHistMap[0] ******* \n" << std::endl;
1289  }
1290  const RecoHists& rhistsStitched = fRecoHistMap[0];
1291 
1292  std::vector < std::vector <unsigned int> > NtrkIdsAll;
1293  std::vector < double > ntvsorted;
1294  hitmap.clear();
1295  KEmap.clear();
1296 
1297 
1298  for (int o = 0; o < ntv; ++o) // o for outer
1299  {
1300 
1301  const art::PtrVector<recob::Track> pvtrack(*(cti++));
1302  auto it = pvtrack.begin();
1303  int ntrack = pvtrack.size();
1304  // if (ntrack>1) std::cout << "\t\t TrkAna: New Stitched Track ******* " << std::endl;
1305  std::vector< std::vector <unsigned int> > NtrkId_Hit; // hit IDs in inner tracks
1306  std::vector<unsigned int> vecMode;
1307 
1308  for(int i = 0; i < ntrack; ++i) {
1309 
1310  const art::Ptr<recob::Track> ptrack(*(it++));
1311  // const recob::Track& track = *ptrack;
1312  auto pcoll={ ptrack };
1313  art::FindManyP<recob::SpacePoint> fs( pcoll, evt, fTrkSpptAssocModuleLabel);
1314  // From gdb> ptype fs, the vector of Ptr<SpacePoint>s it appears is grabbed after fs.at(0)
1315  bool assns(true);
1316  try {
1317  // Get Spacepoints from this Track, get Hits from those Spacepoints.
1318  // int nsppts = ptrack->NumberTrajectoryPoints();
1319 
1320  int nsppts_assn = fs.at(0).size();
1321  // if (ntrack>1) std::cout << "\t\tTrackAnaCT: Number of Spacepoints from Track.NumTrajPts(): " << nsppts << std::endl;
1322  // if (ntrack>1) std::cout << "\t\tTrackAnaCT: Number of Spacepoints from Assns for this Track: " << nsppts_assn << std::endl;
1323  //assert (nsppts_assn == nsppts);
1324  auto sppt = fs.at(0);//.at(is);
1325  art::FindManyP<recob::Hit> fh( sppt, evt, fHitSpptAssocModuleLabel);
1326  // Importantly, loop on all sppts, though they don't all contribute to the track.
1327  // As opposed to looping on the trajectory pts, which is a lower number.
1328  // Also, important, in job in whch this runs I set TrackKal3DSPS parameter MaxPass=1,
1329  // cuz I don't want merely the sparse set of sppts as follows from the uncontained
1330  // MS-measurement in 2nd pass.
1331  std::vector <unsigned int> vecNtrkIds;
1332  for(int is = 0; is < nsppts_assn; ++is) {
1333  int nhits = fh.at(is).size(); // should be 2 or 3: number of planes.
1334  for(int ih = 0; ih < nhits; ++ih) {
1335  auto hit = fh.at(is).at(ih); // Our vector is after the .at(is) this time.
1336  if (hit->SignalType()!=geo::kCollection) continue;
1337  rhistsStitched.fHHitChg->Fill(hit->Integral());
1338  rhistsStitched.fHHitWidth->Fill(hit->RMS() * 2.);
1339  if (mc)
1340  {
1341  std::vector<sim::TrackIDE> tids = bt_serv->HitToTrackIDEs(clockData, hit);
1342  // more here.
1343  // Loop over track ids.
1344  bool justOne(true); // Only take first trk that contributed to this hit
1345  // std::cout << "\t\t TrkAna: TrkId tids.size() ******* " << tids.size() <<std::endl;
1346  for(std::vector<sim::TrackIDE>::const_iterator itid = tids.begin();itid != tids.end(); ++itid) {
1347  int trackID = std::abs(itid->trackID);
1348  hitmap[trackID][o].push_back(hit);
1349 
1350  if (justOne) { vecNtrkIds.push_back(trackID); justOne=false; }
1351  // Add hit to PtrVector corresponding to this track id.
1352  rhistsStitched.fHHitTrkId->Fill(trackID);
1353  const simb::MCParticle* part = pi_serv->TrackIdToParticle_P(trackID);
1354  rhistsStitched.fHHitPdg->Fill(part->PdgCode());
1355  // This really needs to be indexed as KE deposited in volTPC, not just KE. EC, 24-July-2014.
1356 
1357  TVector3 mcstart;
1358  TVector3 mcend;
1359  TVector3 mcstartmom;
1360  TVector3 mcendmom;
1361  double mctime = part->T(); // nsec
1362  double mcdx = mctime * 1.e-3 * detProp.DriftVelocity(); // cm
1363 
1364  double plen = length(detProp, *part, mcdx, mcstart, mcend, mcstartmom, mcendmom);
1365 
1366  KEmap[(int)(1e6*plen)] = trackID; // multiple assignment but always the same, so fine.
1367  // std::cout << "\t\t TrkAna: TrkId trackID, KE [MeV] ******* " << trackID << ", " << (int)(1e3*(part->E()-part->Mass())) <<std::endl;
1368  }
1369 
1370  } // mc
1371  } // hits
1372 
1373  } // spacepoints
1374 
1375  if (mc)
1376  {
1377  NtrkId_Hit.push_back(vecNtrkIds);
1378  // Find the trkID mode for this i^th track
1379  unsigned int ii(1);
1380  int max(-12), n(1), ind(0);
1381  std::sort(vecNtrkIds.begin(),vecNtrkIds.end());
1382  std::vector<unsigned int> strkIds(vecNtrkIds);
1383  while ( ii < vecNtrkIds.size() )
1384  {
1385  if (strkIds.at(ii) != strkIds.at(ii-1))
1386  {
1387  n=1;
1388  }
1389  else
1390  {
1391  n++;
1392  }
1393  if (n>max) { max = n; ind = ii;}
1394  ii++;
1395  }
1396  // std::cout << "\t\t TrkAna: TrkId ind for this track is ******* " << ind <<std::endl;
1397  unsigned int mode(sim::NoParticleId);
1398  if (strkIds.begin()!=strkIds.end())
1399  mode = strkIds.at(ind);
1400  vecMode.push_back(mode);
1401 
1402  // if (ntrack>1) std::cout << "\t\t TrkAna: TrkId mode for this component track is ******* " << mode <<std::endl;
1403  if (strkIds.size()!=0)
1404  rhistsStitched.fModeFrac->Fill((double)max/(double)strkIds.size());
1405  else
1406  rhistsStitched.fModeFrac->Fill(-1.0);
1407  } // mc
1408 
1409  } // end try
1410  catch (cet::exception& x) {
1411  assns = false;
1412  }
1413  if (!assns) throw cet::exception("TrackAnaCT") << "Bad Associations. \n";
1414 
1415  } // i
1416 
1417  if (mc)
1418  {
1419  // one vector per o trk, for all modes of stitched i trks
1420  NtrkIdsAll.push_back(vecMode);
1421 
1422  std::unique(NtrkIdsAll.back().begin(),NtrkIdsAll.back().end());
1423  double sum(0.0);
1424  for (auto const val : NtrkIdsAll.back())
1425  {
1426  // rhistsStitched.fNTrkIdTrks3->Fill(o,val%100,hitmap[val][o].size());
1427  sum += hitmap[val][o].size();
1428  }
1429  ntvsorted.push_back(sum);
1430 
1431  }
1432 
1433  //
1434  } // o
1435 
1436  int vtmp(0);
1437  // get KEmap indices by most energetic first, least last.
1438  for (auto it = KEmap.rbegin(); it!=KEmap.rend(); ++it)
1439  {
1440  // int tval = it->second; // grab trkIDs in order, since they're sorted by KE
1441  // int ke = it->first; // grab trkIDs in order, since they're sorted by KE
1442  // const simb::MCParticle* part = pi_serv->TrackIDToParticle(tval);
1443 
1444  // std::cout << "TrackAnaCTStitch: KEmap cntr vtmp, Length ke, tval, pdg : " << vtmp << ", " << ke <<", " << tval <<", " << part->PdgCode() << ", " << std::endl;
1445 
1446  vtmp++;
1447  }
1448 
1449  // get o trk indices by their hits. Most populous track first, least last.
1450  for (auto const o : fsort_indexes(ntvsorted))
1451  {
1452  int v(0);
1453  // get KEmap indices by longest trajectory first, least last.
1454  for (auto it = KEmap.rbegin(); it!=KEmap.rend(); ++it)
1455  {
1456  int val = it->second; // grab trkIDs in order, since they're sorted by KE
1457  // const simb::MCParticle* part = pi_serv->TrackIDToParticle(val);
1458  // std::cout << "TrackAnaCTStitch: trk o, KEmap cntr v, KE val, pdg hitmap[val][o].size(): " << o <<", " << v << ", " << val <<", " << part->PdgCode() << ", " << hitmap[val][o].size() << std::endl;
1459  rhistsStitched.fNTrkIdTrks3->Fill(o,v,hitmap[val][o].size());
1460  v++;
1461  }
1462 
1463  }
1464 
1465  // In how many o tracks did each trkId appear? Histo it. Would like it to be precisely 1.
1466  // Histo it vs. particle KE.
1467  flattener flat(NtrkIdsAll);
1468  std::vector <unsigned int> &v = flat;
1469  //auto const it ( std::unique(v.begin(),v.end()) ); // never use this it, perhaps.
1470  for (auto const val : v)
1471  {
1472  if (val != (unsigned int)sim::NoParticleId)
1473  {
1474  const simb::MCParticle* part = pi_serv->TrackIdToParticle_P( val );
1475  double T(part->E() - 0.001*part->Mass());
1476  rhistsStitched.fNTrkIdTrks->Fill(std::count(v.begin(),v.end(),val));
1477  rhistsStitched.fNTrkIdTrks2->Fill(std::count(v.begin(),v.end(),val),T);
1478  }
1479  else
1480  {
1481  rhistsStitched.fNTrkIdTrks2->Fill(-1.0,0.0);
1482  }
1483  }
1484 
1485  }
simb::MCParticle::E
double E(const int i=0) const
Definition: MCParticle.h:233
simb::MCParticle::PdgCode
int PdgCode() const
Definition: MCParticle.h:212
cheat::BackTrackerService::HitToTrackIDEs
std::vector< sim::TrackIDE > HitToTrackIDEs(detinfo::DetectorClocksData const &clockData, recob::Hit const &hit) const
Definition: BackTrackerService_service.cc:173
art::DataViewImpl::getHandle
Handle< PROD > getHandle(SelectorBase const &) const
Definition: DataViewImpl.h:382
cheat::ParticleInventoryService::TrackIdToParticle_P
const simb::MCParticle * TrackIdToParticle_P(int id) const
Definition: ParticleInventoryService_service.cc:128
simb::MCParticle::Mass
double Mass() const
Definition: MCParticle.h:239
genie::units::fs
static constexpr double fs
Definition: Units.h:100
const_iterator
intermediate_table::const_iterator const_iterator
Definition: intermediate_table.cc:28
util::size
decltype(auto) constexpr size(T &&obj)
ADL-aware version of std::size.
Definition: StdUtils.h:92
trkf::TrackAnaCT::fsort_indexes
std::vector< size_t > fsort_indexes(const std::vector< T > &v)
Definition: TrackAnaCT_module.cc:1520
trkf::TrackAnaCT::fStitchModuleLabel
std::string fStitchModuleLabel
Definition: TrackAnaCT_module.cc:400
art::DataViewImpl::isRealData
bool isRealData() const
Definition: DataViewImpl.cc:113
abs
T abs(T value)
Definition: sparse_vector_test.cc:30
sim::NoParticleId
static const int NoParticleId
Definition: sim.h:28
cet::n
std::void_t< T > n
Definition: metaprogramming.h:28
simb::MCParticle::T
double T(const int i=0) const
Definition: MCParticle.h:224
max
static int max(int a, int b)
Definition: fortranscanner.cpp:60366
trkf::TrackAnaCT::fHitSpptAssocModuleLabel
std::string fHitSpptAssocModuleLabel
Definition: TrackAnaCT_module.cc:402
trkf::TrackAnaCT::fTrackModuleLabel
std::string fTrackModuleLabel
Definition: TrackAnaCT_module.cc:398
hit
Detector simulation of raw signals on wires.
Definition: NumberOfHitsFilter_module.cc:25
keras_to_tensorflow.int
int
Definition: keras_to_tensorflow.py:69
trkf::TrackAnaCT::fTrkSpptAssocModuleLabel
std::string fTrkSpptAssocModuleLabel
Definition: TrackAnaCT_module.cc:401
trkf::TrackAnaCT::fSpacepointModuleLabel
std::string fSpacepointModuleLabel
Definition: TrackAnaCT_module.cc:399
reco_momentum_tuples.count
int count
Definition: reco_momentum_tuples.py:202
train.x
list x
Definition: train.py:276
pandora_metadata.ih
ih
Definition: pandora_metadata.py:81
make_pion_angle.flat
def flat(h_in)
Definition: make_pion_angle.py:38
trkf::TrackAnaCT::fRecoHistMap
std::map< int, RecoHists > fRecoHistMap
Definition: TrackAnaCT_module.cc:416
fhicl::exception
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
endl
QTextStream & endl(QTextStream &s)
Definition: qtextstream.cpp:2030
geo::kCollection
Signal from collection planes.
Definition: geo_types.h:146
void trkf::TrackAnaCT::endJob ( )
virtual

Reimplemented from art::EDAnalyzer.

Definition at line 1487 of file TrackAnaCT_module.cc.

1491  {
1492  // Print summary.
1493 
1494  mf::LogInfo("TrackAnaCT")
1495  << "TrackAnaCT statistics:\n"
1496  << " Number of events = " << fNumEvent;
1497 
1498  // Fill efficiency histograms.
1499 
1500  for(std::map<int, MCHists>::const_iterator i = fMCHistMap.begin();
1501  i != fMCHistMap.end(); ++i) {
1502  const MCHists& mchists = i->second;
1503  effcalc(mchists.fHgstartx, mchists.fHmcstartx, mchists.fHestartx);
1504  effcalc(mchists.fHgstarty, mchists.fHmcstarty, mchists.fHestarty);
1505  effcalc(mchists.fHgstartz, mchists.fHmcstartz, mchists.fHestartz);
1506  effcalc(mchists.fHgendx, mchists.fHmcendx, mchists.fHeendx);
1507  effcalc(mchists.fHgendy, mchists.fHmcendy, mchists.fHeendy);
1508  effcalc(mchists.fHgendz, mchists.fHmcendz, mchists.fHeendz);
1509  effcalc(mchists.fHgtheta, mchists.fHmctheta, mchists.fHetheta);
1510  effcalc(mchists.fHgphi, mchists.fHmcphi, mchists.fHephi);
1511  effcalc(mchists.fHgtheta_xz, mchists.fHmctheta_xz, mchists.fHetheta_xz);
1512  effcalc(mchists.fHgtheta_yz, mchists.fHmctheta_yz, mchists.fHetheta_yz);
1513  effcalc(mchists.fHgmom, mchists.fHmcmom, mchists.fHemom);
1514  effcalc(mchists.fHglen, mchists.fHmclen, mchists.fHelen);
1515  }
1516  }
mf::LogInfo
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
Definition: MessageLogger.h:209
const_iterator
intermediate_table::const_iterator const_iterator
Definition: intermediate_table.cc:28
trkf::TrackAnaCT::fMCHistMap
std::map< int, MCHists > fMCHistMap
Definition: TrackAnaCT_module.cc:415
template<typename T >
std::vector< size_t > trkf::TrackAnaCT::fsort_indexes ( const std::vector< T > &  v)
private

Definition at line 1520 of file TrackAnaCT_module.cc.

1520  {
1521  // initialize original index locations
1522  std::vector<size_t> idx(v.size());
1523  for (size_t i = 0; i != idx.size(); ++i) idx[i] = i;
1524  // sort indexes based on comparing values in v
1525  std::sort(idx.begin(), idx.end(),
1526  [&v](size_t i1, size_t i2) {return v[i1] > v[i2];}); // Want most occupied trks first. was <, EC, 21-July.
1527  return idx;
1528  }

Member Data Documentation

int trkf::TrackAnaCT::fDump
private

Definition at line 404 of file TrackAnaCT_module.cc.

std::string trkf::TrackAnaCT::fHitSpptAssocModuleLabel
private

Definition at line 402 of file TrackAnaCT_module.cc.

bool trkf::TrackAnaCT::fIgnoreSign
private

Definition at line 410 of file TrackAnaCT_module.cc.

double trkf::TrackAnaCT::fMatchColinearity
private

Definition at line 407 of file TrackAnaCT_module.cc.

double trkf::TrackAnaCT::fMatchDisp
private

Definition at line 408 of file TrackAnaCT_module.cc.

std::map<int, MCHists> trkf::TrackAnaCT::fMCHistMap
private

Definition at line 415 of file TrackAnaCT_module.cc.

double trkf::TrackAnaCT::fMinMCKE
private

Definition at line 405 of file TrackAnaCT_module.cc.

double trkf::TrackAnaCT::fMinMCLen
private

Definition at line 406 of file TrackAnaCT_module.cc.

int trkf::TrackAnaCT::fNumEvent
private

Definition at line 420 of file TrackAnaCT_module.cc.

std::map<int, RecoHists> trkf::TrackAnaCT::fRecoHistMap
private

Definition at line 416 of file TrackAnaCT_module.cc.

std::string trkf::TrackAnaCT::fSpacepointModuleLabel
private

Definition at line 399 of file TrackAnaCT_module.cc.

bool trkf::TrackAnaCT::fStitchedAnalysis
private

Definition at line 411 of file TrackAnaCT_module.cc.

std::string trkf::TrackAnaCT::fStitchModuleLabel
private

Definition at line 400 of file TrackAnaCT_module.cc.

std::string trkf::TrackAnaCT::fTrackModuleLabel
private

Definition at line 398 of file TrackAnaCT_module.cc.

std::string trkf::TrackAnaCT::fTrkSpptAssocModuleLabel
private

Definition at line 401 of file TrackAnaCT_module.cc.

double trkf::TrackAnaCT::fWMatchDisp
private

Definition at line 409 of file TrackAnaCT_module.cc.

The documentation for this class was generated from the following file:
Generated by   doxygen 1.8.11