CRT Clustering

sbncode/CAFMaker/CAFMakerParams.h

@@ -252,6 +252,18 @@ namespace caf
       "pandoraTrackCRTTrack"
     };
 
+    Atom<string> CRTSpacePointMatchLabel {
+      Name("CRTSpacePointMatchLabel"),
+      Comment("Base label of track to CRT spacepoint matching producer."),
+      "crtspacepointmatching"
+    };
+
+    Atom<string> SBNDCRTTrackMatchLabel {
+      Name("SBNDCRTTrackMatchLabel"),
+      Comment("Base label of track to SBND CRT track matching producer."),
+      "crttrackmatching"
+    };
+
     Atom<string> TrackMCSLabel {
       Name("TrackMCSLabel"),
       Comment("Base label of track MCS momentum calculation producer."),
@@ -267,15 +279,27 @@ namespace caf
     Atom<string> CRTHitLabel {
       Name("CRTHitLabel"),
       Comment("Label of sbn CRT hits."),
-      "crthit" // same for icarus and sbnd
+      "crthit" // icarus
     };
 
     Atom<string> CRTTrackLabel {
       Name("CRTTrackLabel"),
       Comment("Label of sbn CRT tracks."),
-      "crttrack" // same for icarus and sbnd
+      "crttrack" // icarus
     };
-
+
+    Atom<string> CRTSpacePointLabel {
+      Name("CRTSpacePointLabel"),
+      Comment("Label of sbnd CRT spacepoints."),
+      "crtspacepoints" // sbnd
+    };
+
+    Atom<string> SBNDCRTTrackLabel {
+      Name("SBNDCRTTrackLabel"),
+      Comment("Label of sbnd CRT tracks."),
+      "crttracks" // sbnd
+    };
+
     Atom<string> CRTPMTLabel {
       Name("CRTPMTLabel"),
       Comment("Label for the CRTPMT Matched variables from the crtpmt data product"),

sbncode/CAFMaker/CAFMaker_module.cc

@@ -276,6 +276,10 @@ class CAFMaker : public art::EDProducer {
   art::FindOneP<T> FindOnePStrict(const U& from, const art::Event& evt,
 				  const art::InputTag& label) const;
 
+  template <class T, class D, class U>
+  art::FindOneP<T, D> FindOnePDStrict(const U& from,
+                                      const art::Event& evt,
+                                      const art::InputTag& tag) const;
 
   /// \brief Retrieve an object from an association, with error handling
   ///
@@ -1005,6 +1009,25 @@ art::FindOneP<T> CAFMaker::FindOnePStrict(const U& from,
   return ret;
 }
 
+//......................................................................
+template <class T, class D, class U>
+art::FindOneP<T, D> CAFMaker::FindOnePDStrict(const U& from,
+                                              const art::Event& evt,
+                                              const art::InputTag& tag) const {
+  art::FindOneP<T, D> ret(from, evt, tag);
+
+  if (!tag.label().empty() && !ret.isValid() && fParams.StrictMode()) {
+    std::cout << "CAFMaker: No Assn from '"
+              << cet::demangle_symbol(typeid(from).name()) << "' to '"
+              << cet::demangle_symbol(typeid(T).name())
+              << "' found under label '" << tag << "'. "
+              << "Set 'StrictMode: false' to continue anyway." << std::endl;
+    abort();
+  }
+
+  return ret;
+}
+
 //......................................................................
 template <class T>
 bool CAFMaker::GetAssociatedProduct(const art::FindManyP<T>& fm, int idx,
@@ -1311,40 +1334,64 @@ void CAFMaker::produce(art::Event& evt) noexcept {
     pass_flash_trig = *flashtrig_handle;
   }
 
-  // Fill various detector information associated with the event
-  //
-  // Get all of the CRT hits
-  std::vector<caf::SRCRTHit> srcrthits;
-
-  art::Handle<std::vector<sbn::crt::CRTHit>> crthits_handle;
-  GetByLabelStrict(evt, fParams.CRTHitLabel(), crthits_handle);
-  // fill into event
-  //int64_t CRT_T0_reference_time = fParams.ReferenceCRTT0ToBeam() ? -srtrigger.beam_gate_time_abs : 0; // ns, signed
-  //double CRT_T1_reference_time = fParams.ReferenceCRTT1FromTriggerToBeam() ? srtrigger.trigger_within_gate : 0.;
   int64_t CRT_T0_reference_time = isRealData ?  -srtrigger.beam_gate_time_abs : -fParams.CRTSimT0Offset();
   double CRT_T1_reference_time = isRealData ? srtrigger.trigger_within_gate : -fParams.CRTSimT0Offset();
 
-  if (crthits_handle.isValid()) {
-    const std::vector<sbn::crt::CRTHit> &crthits = *crthits_handle;
-    for (unsigned i = 0; i < crthits.size(); i++) {
-      srcrthits.emplace_back();
-      FillCRTHit(crthits[i], fParams.CRTUseTS0(), CRT_T0_reference_time, CRT_T1_reference_time, srcrthits.back());
-    }
-  }
+  // Fill various detector information associated with the event
 
-  // Get all of the CRT Tracks
+  std::vector<caf::SRCRTHit> srcrthits;
   std::vector<caf::SRCRTTrack> srcrttracks;
+  std::vector<caf::SRCRTSpacePoint> srcrtspacepoints;
+  std::vector<caf::SRSBNDCRTTrack> srsbndcrttracks;
+
+  if(fDet == kICARUS)
+    {
+      art::Handle<std::vector<sbn::crt::CRTHit>> crthits_handle;
+      GetByLabelStrict(evt, fParams.CRTHitLabel(), crthits_handle);
+      // fill into event
+      if (crthits_handle.isValid()) {
+        const std::vector<sbn::crt::CRTHit> &crthits = *crthits_handle;
+        for (unsigned i = 0; i < crthits.size(); i++) {
+          srcrthits.emplace_back();
+          FillCRTHit(crthits[i], fParams.CRTUseTS0(), CRT_T0_reference_time, CRT_T1_reference_time, srcrthits.back());
+        }
+      }
 
-  art::Handle<std::vector<sbn::crt::CRTTrack>> crttracks_handle;
-  GetByLabelStrict(evt, fParams.CRTTrackLabel(), crttracks_handle);
-  // fill into event
-  if (crttracks_handle.isValid()) {
-    const std::vector<sbn::crt::CRTTrack> &crttracks = *crttracks_handle;
-    for (unsigned i = 0; i < crttracks.size(); i++) {
-      srcrttracks.emplace_back();
-      FillCRTTrack(crttracks[i], fParams.CRTUseTS0(), srcrttracks.back());
+      art::Handle<std::vector<sbn::crt::CRTTrack>> crttracks_handle;
+      GetByLabelStrict(evt, fParams.CRTTrackLabel(), crttracks_handle);
+      // fill into event
+      if (crttracks_handle.isValid()) {
+        const std::vector<sbn::crt::CRTTrack> &crttracks = *crttracks_handle;
+        for (unsigned i = 0; i < crttracks.size(); i++) {
+          srcrttracks.emplace_back();
+          FillCRTTrack(crttracks[i], fParams.CRTUseTS0(), srcrttracks.back());
+        }
+      }
+    }
+  else if(fDet == kSBND)
+    {
+      art::Handle<std::vector<sbnd::crt::CRTSpacePoint>> crtspacepoints_handle;
+      GetByLabelStrict(evt, fParams.CRTSpacePointLabel(), crtspacepoints_handle);
+
+      if (crtspacepoints_handle.isValid()) {
+        const std::vector<sbnd::crt::CRTSpacePoint> &crtspacepoints = *crtspacepoints_handle;
+        for (unsigned i = 0; i < crtspacepoints.size(); i++) {
+          srcrtspacepoints.emplace_back();
+          FillCRTSpacePoint(crtspacepoints[i], srcrtspacepoints.back());
+        }
+      }
+
+      art::Handle<std::vector<sbnd::crt::CRTTrack>> sbndcrttracks_handle;
+      GetByLabelStrict(evt, fParams.SBNDCRTTrackLabel(), sbndcrttracks_handle);
+      // fill into event
+      if (sbndcrttracks_handle.isValid()) {
+        const std::vector<sbnd::crt::CRTTrack> &sbndcrttracks = *sbndcrttracks_handle;
+        for (unsigned i = 0; i < sbndcrttracks.size(); i++) {
+          srsbndcrttracks.emplace_back();
+          FillSBNDCRTTrack(sbndcrttracks[i], srsbndcrttracks.back());
+        }
+      }
     }
-  }
 
   // Get all of the CRTPMT Matches .. 
   std::vector<caf::SRCRTPMTMatch> srcrtpmtmatches;
@@ -1622,6 +1669,14 @@ void CAFMaker::produce(art::Event& evt) noexcept {
       FindManyPStrict<anab::T0>(slcTracks, evt,
                fParams.CRTTrackMatchLabel() + slice_tag_suff);
 
+    art::FindOneP<sbnd::crt::CRTSpacePoint, anab::T0> foCRTSpacePointMatch =
+      FindOnePDStrict<sbnd::crt::CRTSpacePoint, anab::T0>(slcTracks, evt,
+               fParams.CRTSpacePointMatchLabel() + slice_tag_suff);
+
+    art::FindOneP<sbnd::crt::CRTTrack, anab::T0> foSBNDCRTTrackMatch =
+      FindOnePDStrict<sbnd::crt::CRTTrack, anab::T0>(slcTracks, evt,
+               fParams.SBNDCRTTrackMatchLabel() + slice_tag_suff);
+
     std::vector<art::FindManyP<recob::MCSFitResult>> fmMCSs;
     static const std::vector<std::string> PIDnames {"muon", "pion", "kaon", "proton"};
     for (std::string pid: PIDnames) {
@@ -1836,7 +1891,7 @@ void CAFMaker::produce(art::Event& evt) noexcept {
               fParams.TrackHitFillRRStartCut(), fParams.TrackHitFillRREndCut(),
               lar::providerFrom<geo::Geometry>(), dprop, trk);
         }
-        if (fmCRTHitMatch.isValid()) {
+        if (fmCRTHitMatch.isValid() && fDet == kICARUS) {
           art::FindManyP<sbn::crt::CRTHit> CRTT02Hit = FindManyPStrict<sbn::crt::CRTHit>
               (fmCRTHitMatch.at(iPart), evt, fParams.CRTHitMatchLabel() + slice_tag_suff);
 
@@ -1846,9 +1901,25 @@ void CAFMaker::produce(art::Event& evt) noexcept {
           FillTrackCRTHit(fmCRTHitMatch.at(iPart), crthitmatch, fParams.CRTUseTS0(), CRT_T0_reference_time, CRT_T1_reference_time, trk);
         }
         // NOTE: SEE TODO AT fmCRTTrackMatch
-        if (fmCRTTrackMatch.isValid()) {
+        if (fmCRTTrackMatch.isValid() && fDet == kICARUS) {
           FillTrackCRTTrack(fmCRTTrackMatch.at(iPart), trk);
         }
+
+        if(foCRTSpacePointMatch.isValid() && fDet == kSBND)
+          {
+            const art::Ptr<sbnd::crt::CRTSpacePoint> crtspacepoint = foCRTSpacePointMatch.at(iPart);
+
+            if(crtspacepoint.isNonnull())
+              FillTrackCRTSpacePoint(foCRTSpacePointMatch.data(iPart).ref(), crtspacepoint, trk);
+          }
+        if(foSBNDCRTTrackMatch.isValid() && fDet == kSBND)
+          {
+            const art::Ptr<sbnd::crt::CRTTrack> sbndcrttrack = foSBNDCRTTrackMatch.at(iPart);
+
+            if(sbndcrttrack.isNonnull())
+              FillTrackSBNDCRTTrack(foSBNDCRTTrackMatch.data(iPart).ref(), sbndcrttrack, trk);
+          }
+
         // Truth matching
         if (fmTrackHit.isValid()) {
           if ( !isRealData ) {
@@ -1912,17 +1983,21 @@ void CAFMaker::produce(art::Event& evt) noexcept {
   //#######################################################
   //  Fill rec Tree
   //#######################################################
-  rec.nslc            = rec.slc.size();
-  rec.mc              = srtruthbranch;
-  rec.fake_reco       = srfakereco;
-  rec.nfake_reco      = srfakereco.size();
-  rec.pass_flashtrig  = pass_flash_trig;  // trigger result
-  rec.crt_hits        = srcrthits;
-  rec.ncrt_hits       = srcrthits.size();
-  rec.crt_tracks        = srcrttracks;
-  rec.ncrt_tracks       = srcrttracks.size();
-  rec.opflashes       = srflashes;
-  rec.nopflashes      = srflashes.size();
+  rec.nslc             = rec.slc.size();
+  rec.mc               = srtruthbranch;
+  rec.fake_reco        = srfakereco;
+  rec.nfake_reco       = srfakereco.size();
+  rec.pass_flashtrig   = pass_flash_trig;  // trigger result
+  rec.crt_hits         = srcrthits;
+  rec.ncrt_hits        = srcrthits.size();
+  rec.crt_tracks       = srcrttracks;
+  rec.ncrt_tracks      = srcrttracks.size();
+  rec.crt_spacepoints  = srcrtspacepoints;
+  rec.ncrt_spacepoints = srcrtspacepoints.size();
+  rec.sbnd_crt_tracks  = srsbndcrttracks;
+  rec.nsbnd_crt_tracks = srsbndcrttracks.size();
+  rec.opflashes        = srflashes;
+  rec.nopflashes       = srflashes.size();
   if (fParams.FillTrueParticles()) {
     rec.true_particles  = true_particles;
   }

sbncode/CAFMaker/CMakeLists.txt

@@ -34,6 +34,7 @@ art_make_library( LIBRARY_NAME sbncode_CAFMaker
                   sbnobj::Common_CRT
                   sbnobj::Common_Reco
                   sbnobj::Common_Analysis
+                  sbnobj::SBND_CRT
                   lardataalg::DetectorInfo
                   art::Framework_Services_System_TriggerNamesService_service
                   sbncode_Metadata_MetadataSBN_service

sbncode/CAFMaker/FillReco.cxx

@@ -112,6 +112,31 @@ namespace caf
     srtrack.hitb.plane = track.plane2;
   }
 
+  void FillCRTSpacePoint(const sbnd::crt::CRTSpacePoint &spacepoint,
+                         caf::SRCRTSpacePoint &srspacepoint,
+                         bool allowEmpty)
+  {
+    srspacepoint.position     = SRVector3D(spacepoint.X(), spacepoint.Y(), spacepoint.Z());
+    srspacepoint.position_err = SRVector3D(spacepoint.XErr(), spacepoint.YErr(), spacepoint.ZErr());
+    srspacepoint.pe           = spacepoint.PE();
+    srspacepoint.time         = spacepoint.Time();
+    srspacepoint.time_err     = spacepoint.TimeErr();
+    srspacepoint.complete     = spacepoint.Complete();
+  }
+
+  void FillSBNDCRTTrack(const sbnd::crt::CRTTrack &track,
+                        caf::SRSBNDCRTTrack &srsbndcrttrack,
+                        bool allowEmpty)
+  {
+    for(auto const& point : track.Points())
+      srsbndcrttrack.points.emplace_back(point.X(), point.Y(), point.Z());
+
+    srsbndcrttrack.time    = track.Time();
+    srsbndcrttrack.time_err = track.TimeErr();
+    srsbndcrttrack.pe       = track.PE();
+    srsbndcrttrack.tof      = track.ToF();
+  }
+
   void FillCRTPMTMatch(const sbn::crt::CRTPMTMatching &match,
 		       caf::SRCRTPMTMatch &srmatch,
 		       bool allowEmpty){
@@ -540,6 +565,24 @@ namespace caf
     }
   }
 
+  void FillTrackCRTSpacePoint(const anab::T0 &t0match,
+                              const art::Ptr<sbnd::crt::CRTSpacePoint> &spacepointmatch,
+                              caf::SRTrack &srtrack,
+                              bool allowEmpty)
+  {
+    srtrack.crtspacepoint.score = t0match.fTriggerConfidence;
+    FillCRTSpacePoint(*spacepointmatch, srtrack.crtspacepoint.spacepoint);
+  }
+
+  void FillTrackSBNDCRTTrack(const anab::T0 &t0match,
+                             const art::Ptr<sbnd::crt::CRTTrack> &trackmatch,
+                             caf::SRTrack &srtrack,
+                             bool allowEmpty)
+  {
+    srtrack.crtsbndtrack.score = t0match.fTriggerConfidence;
+    FillSBNDCRTTrack(*trackmatch, srtrack.crtsbndtrack.track);
+  }
+
   void FillTrackMCS(const recob::Track& track,
                     const std::array<std::vector<art::Ptr<recob::MCSFitResult>>, 4> &mcs_results,
                     caf::SRTrack& srtrack,

sbncode/CAFMaker/FillReco.h

@@ -37,6 +37,8 @@
 #include "sbnobj/Common/Reco/TPCPMTBarycenterMatch.h"
 #include "sbnobj/Common/CRT/CRTHit.hh"
 #include "sbnobj/Common/CRT/CRTTrack.hh"
+#include "sbnobj/SBND/CRT/CRTSpacePoint.hh"
+#include "sbnobj/SBND/CRT/CRTTrack.hh"
 #include "sbnobj/Common/CRT/CRTPMTMatching.hh"
 #include "nusimdata/SimulationBase/MCParticle.h"
 #include "nusimdata/SimulationBase/MCTruth.h"
@@ -135,6 +137,16 @@ namespace caf
                        caf::SRTrack &srtrack,
                        bool allowEmpty = false);
 
+  void FillTrackCRTSpacePoint(const anab::T0 &t0match,
+                              const art::Ptr<sbnd::crt::CRTSpacePoint> &spacepointmatch,
+                              caf::SRTrack &srtrack,
+                              bool allowEmpty = false);
+
+  void FillTrackSBNDCRTTrack(const anab::T0 &t0match,
+                             const art::Ptr<sbnd::crt::CRTTrack> &trackmatch,
+                             caf::SRTrack &srtrack,
+                             bool allowEmpty = false);
+
   void FillTrackMCS(const recob::Track& track,
                     const std::array<std::vector<art::Ptr<recob::MCSFitResult>>, 4> &mcs_results,
                     caf::SRTrack& srtrack,
@@ -193,6 +205,14 @@ namespace caf
                   caf::SRCRTTrack &srtrack,
                   bool allowEmpty = false);
 
+  void FillCRTSpacePoint(const sbnd::crt::CRTSpacePoint &spacepoint,
+                         caf::SRCRTSpacePoint &srspacepoint,
+                         bool allowEmpty = false);
+
+  void FillSBNDCRTTrack(const sbnd::crt::CRTTrack &track,
+                        caf::SRSBNDCRTTrack &srsbndcrttrack,
+                        bool allowEmpty = false);
+
   void FillOpFlash(const recob::OpFlash &flash,
                   std::vector<recob::OpHit const*> const& hits,
                   int cryo,

sbncode/CosmicID/CMakeLists.txt

@@ -4,6 +4,7 @@ cet_build_plugin( CRUMBS art::module
           fhiclcpp::fhiclcpp
           art::Persistency_Provenance canvas::canvas
           sbnobj::Common_Reco
+          sbnobj::SBND_CRT
           sbncode_GeoWrappers
           sbncode_LArRecoProducer
           lardataobj::RecoBase_AssnsDicts_dict