[ticlv5] additional changes concerning timing and linking (cms-sw#24)

* ticlv5: use quality in time compatibility and move combination of HGCAL and MTD time in the TICLCandidateProducer * ticlv5: fix comments * ticlv5: consider trackerMuon for linking - new PFTICLProducer for timing and muons * ticlv5: add SoA for MTD info (NOTE: use PortableHostCollectionReadRules when available --> pre3) * ticlv5: use the MTD SoA in the TICLCandidate producer * ticlv5: add ETL time to candidates * code checks and code format * ticlv5: fix when MTD timing is not used * ticlv5: fix validation plots
waredjeb · Apr 19, 2024 · b7339f7 · b7339f7
1 parent 14ad18b
commit b7339f7
Show file tree

Hide file tree

Showing 20 changed files with 569 additions and 112 deletions.
diff --git a/DataFormats/HGCalReco/BuildFile.xml b/DataFormats/HGCalReco/BuildFile.xml
@@ -1,9 +1,12 @@
 <use name="DataFormats/Math"/>
 <use name="DataFormats/Provenance"/>
+<use name="DataFormats/Portable"/>
+<use name="DataFormats/SoATemplate"/>
 <use name="DataFormats/Candidate"/>
 <use name="DataFormats/Common"/>
 <use name="DataFormats/TrackReco"/>
 <use name="DataFormats/GeometryVector"/>
+<use name="HeterogeneousCore/AlpakaInterface"/>
 <use name="eigen"/>
 <export>
   <lib name="1"/>

diff --git a/DataFormats/HGCalReco/interface/MtdHostCollection.h b/DataFormats/HGCalReco/interface/MtdHostCollection.h
@@ -0,0 +1,12 @@
+#ifndef DataFormats_HGCalReco_MtdHostCollection_h
+#define DataFormats_HGCalReco_MtdHostCollection_h
+
+#include "DataFormats/Portable/interface/PortableHostCollection.h"
+#include "DataFormats/HGCalReco/interface/MtdSoA.h"
+
+// MtdSoA in host memory
+using MtdHostCollection = PortableHostCollection<MtdSoA>;
+using MtdHostCollectionView = PortableHostCollection<MtdSoA>::View;
+using MtdHostCollectionConstView = PortableHostCollection<MtdSoA>::ConstView;
+
+#endif
diff --git a/DataFormats/HGCalReco/interface/MtdSoA.h b/DataFormats/HGCalReco/interface/MtdSoA.h
@@ -0,0 +1,27 @@
+#ifndef DataFormats_HGCalReco_MtdSoA_h
+#define DataFormats_HGCalReco_MtdSoA_h
+
+#include "DataFormats/SoATemplate/interface/SoALayout.h"
+
+GENERATE_SOA_LAYOUT(MtdSoALayout,
+                    SOA_COLUMN(int32_t, trackAsocMTD),
+                    SOA_COLUMN(float, time0),
+                    SOA_COLUMN(float, time0Err),
+                    SOA_COLUMN(float, time),
+                    SOA_COLUMN(float, timeErr),
+                    SOA_COLUMN(float, MVAquality),
+                    SOA_COLUMN(float, pathLength),
+                    SOA_COLUMN(float, beta),
+                    SOA_COLUMN(float, posInMTD_x),
+                    SOA_COLUMN(float, posInMTD_y),
+                    SOA_COLUMN(float, posInMTD_z),
+                    SOA_COLUMN(float, momentumWithMTD),
+                    SOA_COLUMN(float, probPi),
+                    SOA_COLUMN(float, probK),
+                    SOA_COLUMN(float, probP))
+
+using MtdSoA = MtdSoALayout<>;
+using MtdSoAView = MtdSoA::View;
+using MtdSoAConstView = MtdSoA::ConstView;
+
+#endif
diff --git a/DataFormats/HGCalReco/interface/TICLCandidate.h b/DataFormats/HGCalReco/interface/TICLCandidate.h
@@ -24,10 +24,12 @@ class TICLCandidate : public reco::LeafCandidate {
 
   TICLCandidate(const edm::Ptr<ticl::Trackster>& trackster)
       : LeafCandidate(),
-        idProbabilities_{},
         tracksters_({trackster}),
+        idProbabilities_{},
         time_(trackster->time()),
         timeError_(trackster->timeError()),
+        MTDtime_{0.f},
+        MTDtimeError_{-1.f},
         rawEnergy_(0.f) {}
 
   TICLCandidate(const edm::Ptr<reco::Track> trackPtr, const edm::Ptr<ticl::Trackster>& tracksterPtr)
@@ -36,7 +38,7 @@ class TICLCandidate : public reco::LeafCandidate {
       edm::LogError("TICLCandidate") << "At least one between track and trackster must be valid\n";
 
     if (tracksterPtr.isNonnull()) {
-      tracksters_.push_back(std::move(tracksterPtr));
+      tracksters_.push_back(tracksterPtr);
       auto const& trackster = tracksters_[0].get();
       idProbabilities_ = trackster->id_probabilities();
       if (trackPtr_.isNonnull()) {
@@ -79,8 +81,18 @@ class TICLCandidate : public reco::LeafCandidate {
   inline float time() const { return time_; }
   inline float timeError() const { return timeError_; }
 
-  void setTime(float time) { time_ = time; };
-  void setTimeError(float timeError) { timeError_ = timeError; }
+  void setTime(float time, float timeError) {
+    time_ = time;
+    timeError_ = timeError;
+  };
+
+  inline float MTDtime() const { return MTDtime_; }
+  inline float MTDtimeError() const { return MTDtimeError_; }
+
+  void setMTDTime(float time, float timeError) {
+    MTDtime_ = time;
+    MTDtimeError_ = timeError;
+  };
 
   inline const edm::Ptr<reco::Track> trackPtr() const { return trackPtr_; }
   void setTrackPtr(const edm::Ptr<reco::Track>& trackPtr) { trackPtr_ = trackPtr; }
@@ -114,17 +126,17 @@ class TICLCandidate : public reco::LeafCandidate {
   inline void setIdProbability(ParticleType type, float value) { idProbabilities_[int(type)] = value; }
 
 private:
-  std::array<float, 8> idProbabilities_;
+  // vector of Ptr so Tracksters can come from different collections
+  // and there can be derived classes
   std::vector<edm::Ptr<ticl::Trackster> > tracksters_;
   edm::Ptr<reco::Track> trackPtr_;
+  // Since it contains multiple tracksters, duplicate the probability interface
+  std::array<float, 8> idProbabilities_;
 
   float time_;
   float timeError_;
+  float MTDtime_;
+  float MTDtimeError_;
   float rawEnergy_;
-
-  // vector of Ptr so Tracksters can come from different collections
-  // and there can be derived classes
-
-  // Since it contains multiple tracksters, duplicate the probability interface
 };
 #endif
diff --git a/DataFormats/HGCalReco/interface/Trackster.h b/DataFormats/HGCalReco/interface/Trackster.h
@@ -40,8 +40,8 @@ namespace ticl {
         : barycenter_({0.f, 0.f, 0.f}),
           regressed_energy_(0.f),
           raw_energy_(0.f),
-          time_(0.f),
           boundTime_(0.f),
+          time_(0.f),
           timeError_(-1.f),
           id_probabilities_{},
           raw_pt_(0.f),
@@ -176,8 +176,8 @@ namespace ticl {
     float regressed_energy_;
     float raw_energy_;
     // -99, -1 if not available. ns units otherwise
-    float time_;
     float boundTime_;
+    float time_;
     float timeError_;
 
     // trackster ID probabilities

diff --git a/DataFormats/HGCalReco/src/classes.h b/DataFormats/HGCalReco/src/classes.h
@@ -1,5 +1,7 @@
 #include <vector>
 #include <array>
+#include "DataFormats/HGCalReco/interface/MtdSoA.h"
+#include "DataFormats/HGCalReco/interface/MtdHostCollection.h"
 #include "DataFormats/HGCalReco/interface/Trackster.h"
 #include "DataFormats/HGCalReco/interface/TICLLayerTile.h"
 #include "DataFormats/HGCalReco/interface/TICLSeedingRegion.h"

diff --git a/DataFormats/HGCalReco/src/classes_def.xml b/DataFormats/HGCalReco/src/classes_def.xml
@@ -1,6 +1,7 @@
 
 <lcgdict>
-  <class name="ticl::Trackster" ClassVersion="11">
+  <class name="ticl::Trackster" ClassVersion="12">
+    <version ClassVersion="12" checksum="2111726404"/>
     <version ClassVersion="11" checksum="2995118904"/>
     <version ClassVersion="10" checksum="556627704"/>
     <version ClassVersion="9" checksum="1001808235"/>
@@ -52,11 +53,30 @@
   <class name="edm::Wrapper<TICLSeedingRegion>" />
   <class name="edm::Wrapper<std::vector<TICLSeedingRegion> >" />
 
-  <class name="TICLCandidate" ClassVersion="4">
+  <class name="TICLCandidate" ClassVersion="5">
+   <version ClassVersion="5" checksum="2260471800"/>
    <version ClassVersion="4" checksum="3155265136"/>
    <version ClassVersion="3" checksum="450468662"/>
   </class>
   <class name="std::vector<TICLCandidate>" />
   <class name="edm::Wrapper<TICLCandidate>" />
   <class name="edm::Wrapper<std::vector<TICLCandidate> >" />
+
+  <class name="MtdSoA"/>
+  <class name="MtdSoA::View"/>
+
+  <class name="MtdHostCollection"/>
+  <read
+    sourceClass="MtdHostCollection"
+    targetClass="MtdHostCollection"
+    version="[1-]"
+    source="MtdSoA layout_;"
+    target="buffer_,layout_,view_"
+    embed="false">
+  <![CDATA[
+    MtdHostCollection::ROOTReadStreamer(newObj, onfile.layout_);
+  ]]>
+  </read>
+  <class name="edm::Wrapper<MtdHostCollection>" splitLevel="0"/> 
+
 </lcgdict>
diff --git a/RecoHGCal/TICL/interface/TICLInterpretationAlgoBase.h b/RecoHGCal/TICL/interface/TICLInterpretationAlgoBase.h
@@ -4,6 +4,7 @@
 #include <memory>
 #include <vector>
 #include "DataFormats/CaloRecHit/interface/CaloCluster.h"
+#include "DataFormats/HGCalReco/interface/MtdHostCollection.h"
 #include "DataFormats/HGCalReco/interface/Trackster.h"
 #include "DataFormats/HGCalReco/interface/TICLCandidate.h"
 #include "DataFormats/HGCalReco/interface/TICLLayerTile.h"
@@ -65,20 +66,22 @@ namespace ticl {
     struct TrackTimingInformation {
       const edm::Handle<edm::ValueMap<float>> tkTime_h;
       const edm::Handle<edm::ValueMap<float>> tkTimeErr_h;
+      const edm::Handle<edm::ValueMap<float>> tkQuality_h;
       const edm::Handle<edm::ValueMap<float>> tkBeta_h;
       const edm::Handle<edm::ValueMap<float>> tkPath_h;
       const edm::Handle<edm::ValueMap<GlobalPoint>> tkMtdPos_h;
 
       TrackTimingInformation(const edm::Handle<edm::ValueMap<float>> tkT,
                              const edm::Handle<edm::ValueMap<float>> tkTE,
+                             const edm::Handle<edm::ValueMap<float>> tkQ,
                              const edm::Handle<edm::ValueMap<float>> tkB,
                              const edm::Handle<edm::ValueMap<float>> tkP,
                              const edm::Handle<edm::ValueMap<GlobalPoint>> mtdPos)
-          : tkTime_h(tkT), tkTimeErr_h(tkTE), tkBeta_h(tkB), tkPath_h(tkP), tkMtdPos_h(mtdPos) {}
+          : tkTime_h(tkT), tkTimeErr_h(tkTE), tkQuality_h(tkQ), tkBeta_h(tkB), tkPath_h(tkP), tkMtdPos_h(mtdPos) {}
     };
 
     virtual void makeCandidates(const Inputs& input,
-                                const TrackTimingInformation& inputTiming,
+                                edm::Handle<MtdHostCollection> inputTiming_h,
                                 std::vector<Trackster>& resultTracksters,
                                 std::vector<int>& resultCandidate) = 0;
 

diff --git a/RecoHGCal/TICL/plugins/GeneralInterpretationAlgo.cc b/RecoHGCal/TICL/plugins/GeneralInterpretationAlgo.cc
@@ -13,8 +13,7 @@ GeneralInterpretationAlgo::GeneralInterpretationAlgo(const edm::ParameterSet &co
     : TICLInterpretationAlgoBase(conf, cc),
       del_tk_ts_layer1_(conf.getParameter<double>("delta_tk_ts_layer1")),
       del_tk_ts_int_(conf.getParameter<double>("delta_tk_ts_interface")),
-      del_ts_em_had_(conf.getParameter<double>("delta_ts_em_had")),
-      del_ts_had_had_(conf.getParameter<double>("delta_ts_had_had")) {}
+      timing_quality_threshold_(conf.getParameter<double>("timing_quality_threshold")) {}
 
 void GeneralInterpretationAlgo::initialize(const HGCalDDDConstants *hgcons,
                                            const hgcal::RecHitTools rhtools,
@@ -154,6 +153,7 @@ bool GeneralInterpretationAlgo::timeAndEnergyCompatible(float &total_raw_energy,
                                                         const Trackster &trackster,
                                                         const float &tkT,
                                                         const float &tkTErr,
+                                                        const float &tkQual,
                                                         const float &tkBeta,
                                                         const GlobalPoint &tkMtdPos,
                                                         bool useMTDTiming) {
@@ -171,7 +171,7 @@ bool GeneralInterpretationAlgo::timeAndEnergyCompatible(float &total_raw_energy,
   float tsTErr = trackster.timeError();
 
   bool timeCompatible = false;
-  if (tsT == -99. or tkTErr == -1)
+  if (tsT == -99. or tkTErr == -1 or tkQual < timing_quality_threshold_)
     timeCompatible = true;
   else {
     const auto &barycenter = trackster.barycenter();
@@ -203,10 +203,10 @@ bool GeneralInterpretationAlgo::timeAndEnergyCompatible(float &total_raw_energy,
 }
 
 void GeneralInterpretationAlgo::makeCandidates(const Inputs &input,
-                                               const TrackTimingInformation &inputTiming,
+                                               edm::Handle<MtdHostCollection> inputTiming_h,
                                                std::vector<Trackster> &resultTracksters,
                                                std::vector<int> &resultCandidate) {
-  bool useMTDTiming = inputTiming.tkTime_h.isValid();
+  bool useMTDTiming = inputTiming_h.isValid();
   std::cout << "GeneralInterpretationAlgo " << std::endl;
   const auto tkH = input.tracksHandle;
   const auto maskTracks = input.maskedTracks;
@@ -295,11 +295,12 @@ void GeneralInterpretationAlgo::makeCandidates(const Inputs &input,
 
   }  // TS
 
+  // step 1: tracks -> all tracksters, at firstLayerEE
   std::vector<std::vector<unsigned>> tsNearTk(tracks.size());
   findTrackstersInWindow(
       tracksters, trackPColl, tracksterPropTiles, tsAllProp, del_tk_ts_layer1_, tracksters.size(), tsNearTk);
-  // step 4: tracks -> all tracksters, at lastLayerEE
 
+  // step 2: tracks -> all tracksters, at lastLayerEE
   std::vector<std::vector<unsigned>> tsNearTkAtInt(tracks.size());
   findTrackstersInWindow(
       tracksters, tkPropIntColl, tsPropIntTiles, tsAllPropInt, del_tk_ts_int_, tracksters.size(), tsNearTkAtInt);
@@ -317,16 +318,19 @@ void GeneralInterpretationAlgo::makeCandidates(const Inputs &input,
     std::vector<unsigned int> chargedCandidate;
     float total_raw_energy = 0.f;
 
-    auto tkRef = reco::TrackRef(tkH, i);
     float track_time = 0.f;
     float track_timeErr = 0.f;
+    float track_quality = 0.f;
     float track_beta = 0.f;
     GlobalPoint track_MtdPos{0.f, 0.f, 0.f};
     if (useMTDTiming) {
-      track_time = (*inputTiming.tkTime_h)[tkRef];
-      track_timeErr = (*inputTiming.tkTimeErr_h)[tkRef];
-      track_beta = (*inputTiming.tkBeta_h)[tkRef];
-      track_MtdPos = (*inputTiming.tkMtdPos_h)[tkRef];
+      auto const &inputTimingView = (*inputTiming_h).const_view();
+      track_time = inputTimingView.time()[i];
+      track_timeErr = inputTimingView.timeErr()[i];
+      track_quality = inputTimingView.MVAquality()[i];
+      track_beta = inputTimingView.beta()[i];
+      track_MtdPos = {
+          inputTimingView.posInMTD_x()[i], inputTimingView.posInMTD_y()[i], inputTimingView.posInMTD_z()[i]};
     }
 
     for (auto const tsIdx : tsNearTk[i]) {
@@ -335,6 +339,7 @@ void GeneralInterpretationAlgo::makeCandidates(const Inputs &input,
                                                         tracksters[tsIdx],
                                                         track_time,
                                                         track_timeErr,
+                                                        track_quality,
                                                         track_beta,
                                                         track_MtdPos,
                                                         useMTDTiming)) {
@@ -348,6 +353,7 @@ void GeneralInterpretationAlgo::makeCandidates(const Inputs &input,
                                                         tracksters[tsIdx],
                                                         track_time,
                                                         track_timeErr,
+                                                        track_quality,
                                                         track_beta,
                                                         track_MtdPos,
                                                         useMTDTiming)) {
@@ -383,7 +389,6 @@ void GeneralInterpretationAlgo::fillPSetDescription(edm::ParameterSetDescription
                         "hitPattern().numberOfLostHits(\"MISSING_OUTER_HITS\") < 5");
   desc.add<double>("delta_tk_ts_layer1", 0.02);
   desc.add<double>("delta_tk_ts_interface", 0.03);
-  desc.add<double>("delta_ts_em_had", 0.03);
-  desc.add<double>("delta_ts_had_had", 0.03);
+  desc.add<double>("timing_quality_threshold", 0.5);
   TICLInterpretationAlgoBase::fillPSetDescription(desc);
 }
diff --git a/RecoHGCal/TICL/plugins/GeneralInterpretationAlgo.h b/RecoHGCal/TICL/plugins/GeneralInterpretationAlgo.h
@@ -18,7 +18,7 @@ namespace ticl {
     ~GeneralInterpretationAlgo() override;
 
     void makeCandidates(const Inputs &input,
-                        const TrackTimingInformation &inputTiming,
+                        edm::Handle<MtdHostCollection> inputTiming_h,
                         std::vector<Trackster> &resultTracksters,
                         std::vector<int> &resultCandidate) override;
 
@@ -51,6 +51,7 @@ namespace ticl {
                                  const Trackster &trackster,
                                  const float &tkTime,
                                  const float &tkTimeErr,
+                                 const float &tkQual,
                                  const float &tkBeta,
                                  const GlobalPoint &tkMtdPos,
                                  bool useMTDTiming);
@@ -59,8 +60,7 @@ namespace ticl {
     const float maxDeltaT_ = 3.0f;
     const float del_tk_ts_layer1_;
     const float del_tk_ts_int_;
-    const float del_ts_em_had_;
-    const float del_ts_had_had_;
+    const float timing_quality_threshold_;
 
     const HGCalDDDConstants *hgcons_;