diff --git a/DataFormats/HGCalReco/BuildFile.xml b/DataFormats/HGCalReco/BuildFile.xml
index ea59e3506573a..495243e3981c4 100644
--- a/DataFormats/HGCalReco/BuildFile.xml
+++ b/DataFormats/HGCalReco/BuildFile.xml
@@ -1,9 +1,12 @@
+
+
+
diff --git a/DataFormats/HGCalReco/interface/MtdHostCollection.h b/DataFormats/HGCalReco/interface/MtdHostCollection.h
new file mode 100644
index 0000000000000..49c137c53f498
--- /dev/null
+++ 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;
+using MtdHostCollectionView = PortableHostCollection::View;
+using MtdHostCollectionConstView = PortableHostCollection::ConstView;
+
+#endif
diff --git a/DataFormats/HGCalReco/interface/MtdSoA.h b/DataFormats/HGCalReco/interface/MtdSoA.h
new file mode 100644
index 0000000000000..9776224a11f8a
--- /dev/null
+++ 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
index 3d72d8ebf54c0..8dd610790ad41 100644
--- 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& 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 trackPtr, const edm::Ptr& 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 trackPtr() const { return trackPtr_; }
void setTrackPtr(const edm::Ptr& 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 idProbabilities_;
+ // vector of Ptr so Tracksters can come from different collections
+ // and there can be derived classes
std::vector > tracksters_;
edm::Ptr trackPtr_;
+ // Since it contains multiple tracksters, duplicate the probability interface
+ std::array 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
index 5591f5a0df38e..dd07d324e70e4 100644
--- 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
index d871bfb485a71..51d53d433e05a 100644
--- a/DataFormats/HGCalReco/src/classes.h
+++ b/DataFormats/HGCalReco/src/classes.h
@@ -1,5 +1,7 @@
#include
#include
+#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
index 7d4b3f16d2618..fe1e8bf89a68d 100644
--- a/DataFormats/HGCalReco/src/classes_def.xml
+++ b/DataFormats/HGCalReco/src/classes_def.xml
@@ -1,6 +1,7 @@
-
+
+
@@ -52,11 +53,30 @@
-
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/RecoHGCal/TICL/interface/TICLInterpretationAlgoBase.h b/RecoHGCal/TICL/interface/TICLInterpretationAlgoBase.h
index 784f1dcc4f751..8f016ef31ec8a 100644
--- a/RecoHGCal/TICL/interface/TICLInterpretationAlgoBase.h
+++ b/RecoHGCal/TICL/interface/TICLInterpretationAlgoBase.h
@@ -4,6 +4,7 @@
#include
#include
#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> tkTime_h;
const edm::Handle> tkTimeErr_h;
+ const edm::Handle> tkQuality_h;
const edm::Handle> tkBeta_h;
const edm::Handle> tkPath_h;
const edm::Handle> tkMtdPos_h;
TrackTimingInformation(const edm::Handle> tkT,
const edm::Handle> tkTE,
+ const edm::Handle> tkQ,
const edm::Handle> tkB,
const edm::Handle> tkP,
const edm::Handle> 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 inputTiming_h,
std::vector& resultTracksters,
std::vector& resultCandidate) = 0;
diff --git a/RecoHGCal/TICL/plugins/GeneralInterpretationAlgo.cc b/RecoHGCal/TICL/plugins/GeneralInterpretationAlgo.cc
index 7e166b8266b16..682bc6eb02dca 100644
--- 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("delta_tk_ts_layer1")),
del_tk_ts_int_(conf.getParameter("delta_tk_ts_interface")),
- del_ts_em_had_(conf.getParameter("delta_ts_em_had")),
- del_ts_had_had_(conf.getParameter("delta_ts_had_had")) {}
+ timing_quality_threshold_(conf.getParameter("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 inputTiming_h,
std::vector &resultTracksters,
std::vector &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> 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> 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 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("delta_tk_ts_layer1", 0.02);
desc.add("delta_tk_ts_interface", 0.03);
- desc.add("delta_ts_em_had", 0.03);
- desc.add("delta_ts_had_had", 0.03);
+ desc.add("timing_quality_threshold", 0.5);
TICLInterpretationAlgoBase::fillPSetDescription(desc);
}
diff --git a/RecoHGCal/TICL/plugins/GeneralInterpretationAlgo.h b/RecoHGCal/TICL/plugins/GeneralInterpretationAlgo.h
index 6ef3254ba46e9..f5346d4a35e79 100644
--- 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 inputTiming_h,
std::vector &resultTracksters,
std::vector &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_;
diff --git a/RecoHGCal/TICL/plugins/MTDSoAProducer.cc b/RecoHGCal/TICL/plugins/MTDSoAProducer.cc
new file mode 100644
index 0000000000000..f299069782617
--- /dev/null
+++ b/RecoHGCal/TICL/plugins/MTDSoAProducer.cc
@@ -0,0 +1,151 @@
+#include "FWCore/Framework/interface/Frameworkfwd.h"
+#include "FWCore/Framework/interface/stream/EDProducer.h"
+#include "FWCore/Framework/interface/Event.h"
+#include "FWCore/Framework/interface/EventSetup.h"
+#include "FWCore/Framework/interface/ESHandle.h"
+#include "FWCore/ParameterSet/interface/ParameterSet.h"
+#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
+#include "HeterogeneousCore/AlpakaInterface/interface/host.h"
+
+#include "DataFormats/Common/interface/ValueMap.h"
+#include "DataFormats/HGCalReco/interface/MtdHostCollection.h"
+#include "DataFormats/TrackReco/interface/TrackFwd.h"
+#include "DataFormats/TrackReco/interface/Track.h"
+//#include "DataFormats/TrackReco/interface/TrackBase.h"
+
+using namespace edm;
+
+class MTDSoAProducer : public edm::stream::EDProducer<> {
+public:
+ MTDSoAProducer(const ParameterSet& pset);
+
+ static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
+
+ void produce(edm::Event& ev, const edm::EventSetup& es) final;
+
+private:
+ edm::EDGetTokenT tracksToken_;
+ edm::EDGetTokenT> trackAssocToken_;
+ edm::EDGetTokenT> t0Token_;
+ edm::EDGetTokenT> sigmat0Token_;
+ edm::EDGetTokenT> tmtdToken_;
+ edm::EDGetTokenT> sigmatmtdToken_;
+ edm::EDGetTokenT> betaToken_;
+ edm::EDGetTokenT> pathToken_;
+ edm::EDGetTokenT> MVAQualityToken_;
+ edm::EDGetTokenT> posInMtdToken_;
+ edm::EDGetTokenT> momentumWithMTDToken_;
+ edm::EDGetTokenT> probPiToken_;
+ edm::EDGetTokenT> probKToken_;
+ edm::EDGetTokenT> probPToken_;
+};
+
+MTDSoAProducer::MTDSoAProducer(const ParameterSet& iConfig)
+ : tracksToken_(consumes(iConfig.getParameter("tracksSrc"))),
+ trackAssocToken_(consumes>(iConfig.getParameter("trackAssocSrc"))),
+ t0Token_(consumes>(iConfig.getParameter("t0Src"))),
+ sigmat0Token_(consumes>(iConfig.getParameter("sigmat0Src"))),
+ tmtdToken_(consumes>(iConfig.getParameter("tmtdSrc"))),
+ sigmatmtdToken_(consumes>(iConfig.getParameter("sigmatmtdSrc"))),
+ betaToken_(consumes>(iConfig.getParameter("betamtd"))),
+ pathToken_(consumes>(iConfig.getParameter("pathmtd"))),
+ MVAQualityToken_(consumes>(iConfig.getParameter("mvaquality"))),
+ posInMtdToken_(consumes>(iConfig.getParameter("posmtd"))),
+ momentumWithMTDToken_(consumes>(iConfig.getParameter("momentum"))),
+ probPiToken_(consumes>(iConfig.getParameter("probPi"))),
+ probKToken_(consumes>(iConfig.getParameter("probK"))),
+ probPToken_(consumes>(iConfig.getParameter("probP"))) {
+ produces();
+}
+
+// Configuration descriptions
+void MTDSoAProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
+ edm::ParameterSetDescription desc;
+ desc.add("tracksSrc", edm::InputTag("generalTracks"));
+ desc.add("trackAssocSrc", edm::InputTag("trackExtenderWithMTD:generalTrackassoc"));
+ desc.add("t0Src", edm::InputTag("tofPID:t0"));
+ desc.add("sigmat0Src", edm::InputTag("tofPID:sigmat0"));
+ desc.add("tmtdSrc", edm::InputTag("trackExtenderWithMTD:generalTracktmtd"));
+ desc.add("sigmatmtdSrc", edm::InputTag("trackExtenderWithMTD:generalTracksigmatmtd"));
+ desc.add("betamtd", edm::InputTag("trackExtenderWithMTD:generalTrackBeta"));
+ desc.add("pathmtd", edm::InputTag("trackExtenderWithMTD:generalTrackPathLength"));
+ desc.add("mvaquality", edm::InputTag("mtdTrackQualityMVA:mtdQualMVA"));
+ desc.add("posmtd", edm::InputTag("trackExtenderWithMTD:generalTrackmtdpos"));
+ desc.add("momentum", edm::InputTag("trackExtenderWithMTD:generalTrackp"));
+ desc.add("probPi", edm::InputTag("tofPID:probPi"));
+ desc.add("probK", edm::InputTag("tofPID:probK"));
+ desc.add("probP", edm::InputTag("tofPID:probP"));
+
+ descriptions.add("mtdSoAProducer", desc);
+}
+
+void MTDSoAProducer::produce(edm::Event& ev, const edm::EventSetup& es) {
+ edm::Handle tracksH;
+ ev.getByToken(tracksToken_, tracksH);
+ const auto& tracks = *tracksH;
+
+ const auto& trackAssoc = ev.get(trackAssocToken_);
+
+ const auto& t0 = ev.get(t0Token_);
+ const auto& sigmat0 = ev.get(sigmat0Token_);
+
+ const auto& tmtd = ev.get(tmtdToken_);
+ const auto& sigmatmtd = ev.get(sigmatmtdToken_);
+
+ const auto& beta = ev.get(betaToken_);
+ const auto& path = ev.get(pathToken_);
+ const auto& MVAquality = ev.get(MVAQualityToken_);
+ const auto& posInMTD = ev.get(posInMtdToken_);
+ const auto& momentum = ev.get(momentumWithMTDToken_);
+ const auto& probPi = ev.get(probPiToken_);
+ const auto& probK = ev.get(probKToken_);
+ const auto& probP = ev.get(probPToken_);
+
+ auto MtdInfo = std::make_unique(tracks.size(), cms::alpakatools::host());
+
+ auto& MtdInfoView = MtdInfo->view();
+ for (unsigned int iTrack = 0; iTrack < tracks.size(); ++iTrack) {
+ const reco::TrackRef trackref(tracksH, iTrack);
+
+ if (trackAssoc[trackref] == -1) {
+ MtdInfoView.trackAsocMTD()[iTrack] = -1;
+ MtdInfoView.time0()[iTrack] = 0.f;
+ MtdInfoView.time0Err()[iTrack] = -1.f;
+ MtdInfoView.time()[iTrack] = 0.f;
+ MtdInfoView.timeErr()[iTrack] = -1.f;
+ MtdInfoView.MVAquality()[iTrack] = 0.f;
+ MtdInfoView.pathLength()[iTrack] = 0.f;
+ MtdInfoView.beta()[iTrack] = 0.f;
+ MtdInfoView.posInMTD_x()[iTrack] = 0.f;
+ MtdInfoView.posInMTD_y()[iTrack] = 0.f;
+ MtdInfoView.posInMTD_z()[iTrack] = 0.f;
+ MtdInfoView.momentumWithMTD()[iTrack] = 0.f;
+ MtdInfoView.probPi()[iTrack] = 0.f;
+ MtdInfoView.probK()[iTrack] = 0.f;
+ MtdInfoView.probP()[iTrack] = 0.f;
+ continue;
+ }
+
+ MtdInfoView.trackAsocMTD()[iTrack] = trackAssoc[trackref];
+ MtdInfoView.time0()[iTrack] = t0[trackref];
+ MtdInfoView.time0Err()[iTrack] = sigmat0[trackref];
+ MtdInfoView.time()[iTrack] = tmtd[trackref];
+ MtdInfoView.timeErr()[iTrack] = sigmatmtd[trackref];
+ MtdInfoView.MVAquality()[iTrack] = MVAquality[trackref];
+ MtdInfoView.pathLength()[iTrack] = path[trackref];
+ MtdInfoView.beta()[iTrack] = beta[trackref];
+ MtdInfoView.posInMTD_x()[iTrack] = posInMTD[trackref].x();
+ MtdInfoView.posInMTD_y()[iTrack] = posInMTD[trackref].y();
+ MtdInfoView.posInMTD_z()[iTrack] = posInMTD[trackref].z();
+ MtdInfoView.momentumWithMTD()[iTrack] = momentum[trackref];
+ MtdInfoView.probPi()[iTrack] = probPi[trackref];
+ MtdInfoView.probK()[iTrack] = probK[trackref];
+ MtdInfoView.probP()[iTrack] = probP[trackref];
+ }
+
+ ev.put(std::move(MtdInfo));
+}
+
+//define this as a plug-in
+#include
+DEFINE_FWK_MODULE(MTDSoAProducer);
diff --git a/RecoHGCal/TICL/plugins/PFTICLProducerV5.cc b/RecoHGCal/TICL/plugins/PFTICLProducerV5.cc
new file mode 100644
index 0000000000000..fa70b50ef211c
--- /dev/null
+++ b/RecoHGCal/TICL/plugins/PFTICLProducerV5.cc
@@ -0,0 +1,136 @@
+// This producer converts a list of TICLCandidates to a list of PFCandidates.
+
+#include "FWCore/Framework/interface/Frameworkfwd.h"
+#include "FWCore/Framework/interface/stream/EDProducer.h"
+
+#include "FWCore/Framework/interface/Event.h"
+#include "FWCore/Framework/interface/MakerMacros.h"
+
+#include "DataFormats/Common/interface/View.h"
+
+#include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
+#include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
+
+#include "DataFormats/HGCalReco/interface/TICLCandidate.h"
+
+#include "RecoParticleFlow/PFProducer/interface/PFMuonAlgo.h"
+
+class PFTICLProducerV5 : public edm::stream::EDProducer<> {
+public:
+ PFTICLProducerV5(const edm::ParameterSet&);
+ ~PFTICLProducerV5() override {}
+
+ static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
+
+ void produce(edm::Event&, const edm::EventSetup&) override;
+
+private:
+ // parameters
+ const bool energy_from_regression_;
+ // inputs
+ const edm::EDGetTokenT> ticl_candidates_;
+ const edm::EDGetTokenT muons_;
+ // For PFMuonAlgo
+ std::unique_ptr pfmu_;
+};
+
+DEFINE_FWK_MODULE(PFTICLProducerV5);
+
+PFTICLProducerV5::PFTICLProducerV5(const edm::ParameterSet& conf)
+ : energy_from_regression_(conf.getParameter("energyFromRegression")),
+ ticl_candidates_(consumes>(conf.getParameter("ticlCandidateSrc"))),
+ muons_(consumes(conf.getParameter("muonSrc"))),
+ pfmu_(std::make_unique(conf.getParameterSet("pfMuonAlgoParameters"),
+ false)) { // postMuonCleaning = false
+ produces();
+}
+
+void PFTICLProducerV5::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
+ edm::ParameterSetDescription desc;
+ desc.add("ticlCandidateSrc", edm::InputTag("ticlCandidate"));
+ desc.add("energyFromRegression", true);
+ // For PFMuonAlgo
+ desc.add("muonSrc", edm::InputTag("muons1stStep"));
+ edm::ParameterSetDescription psd_PFMuonAlgo;
+ PFMuonAlgo::fillPSetDescription(psd_PFMuonAlgo);
+ desc.add("pfMuonAlgoParameters", psd_PFMuonAlgo);
+ //
+ descriptions.add("pfTICLProducerV5", desc);
+}
+
+void PFTICLProducerV5::produce(edm::Event& evt, const edm::EventSetup& es) {
+ //get TICLCandidates
+ edm::Handle> ticl_cand_h;
+ evt.getByToken(ticl_candidates_, ticl_cand_h);
+ const auto ticl_candidates = *ticl_cand_h;
+ const auto muonH = evt.getHandle(muons_);
+ const auto& muons = *muonH;
+
+ auto candidates = std::make_unique();
+
+ for (const auto& ticl_cand : ticl_candidates) {
+ const auto abs_pdg_id = std::abs(ticl_cand.pdgId());
+ const auto charge = ticl_cand.charge();
+ const auto& four_mom = ticl_cand.p4();
+ float total_raw_energy = 0.f;
+ float total_em_raw_energy = 0.f;
+ for (const auto& t : ticl_cand.tracksters()) {
+ total_raw_energy += t->raw_energy();
+ total_em_raw_energy += t->raw_em_energy();
+ }
+ float ecal_energy_fraction = total_em_raw_energy / total_raw_energy;
+ float ecal_energy = energy_from_regression_ ? ticl_cand.p4().energy() * ecal_energy_fraction
+ : ticl_cand.rawEnergy() * ecal_energy_fraction;
+ float hcal_energy =
+ energy_from_regression_ ? ticl_cand.p4().energy() - ecal_energy : ticl_cand.rawEnergy() - ecal_energy;
+ // fix for floating point rounding could go slightly below 0
+ hcal_energy = std::max(0.f, hcal_energy);
+ reco::PFCandidate::ParticleType part_type;
+ switch (abs_pdg_id) {
+ case 11:
+ part_type = reco::PFCandidate::e;
+ break;
+ case 13:
+ part_type = reco::PFCandidate::mu;
+ break;
+ case 22:
+ part_type = reco::PFCandidate::gamma;
+ break;
+ case 130:
+ part_type = reco::PFCandidate::h0;
+ break;
+ case 211:
+ part_type = reco::PFCandidate::h;
+ break;
+ // default also handles neutral pions (111) for the time being (not yet foreseen in PFCandidate)
+ default:
+ part_type = reco::PFCandidate::X;
+ }
+
+ candidates->emplace_back(charge, four_mom, part_type);
+
+ auto& candidate = candidates->back();
+ candidate.setEcalEnergy(ecal_energy, ecal_energy);
+ candidate.setHcalEnergy(hcal_energy, hcal_energy);
+ if (candidate.charge()) { // otherwise PFCandidate throws
+ // Construct edm::Ref from edm::Ptr. As of now, assumes type to be reco::Track. To be extended (either via
+ // dynamic type checking or configuration) if additional track types are needed.
+ reco::TrackRef trackref(ticl_cand.trackPtr().id(), int(ticl_cand.trackPtr().key()), &evt.productGetter());
+ candidate.setTrackRef(trackref);
+ // Utilize PFMuonAlgo
+ const int muId = PFMuonAlgo::muAssocToTrack(trackref, muons);
+ if (muId != -1) {
+ const reco::MuonRef muonref = reco::MuonRef(muonH, muId);
+ if ((PFMuonAlgo::isMuon(muonref) and not(*muonH)[muId].isTrackerMuon()) or
+ (ticl_cand.tracksters().empty() and muonref.isNonnull() and muonref->isGlobalMuon())) {
+ const bool allowLoose = (part_type == reco::PFCandidate::mu);
+ // Redefine pfmuon candidate kinematics and add muonref
+ pfmu_->reconstructMuon(candidate, muonref, allowLoose);
+ }
+ }
+ }
+ candidate.setTime(ticl_cand.time(), ticl_cand.timeError());
+ }
+
+ evt.put(std::move(candidates));
+}
diff --git a/RecoHGCal/TICL/plugins/SimTrackstersProducer.cc b/RecoHGCal/TICL/plugins/SimTrackstersProducer.cc
index ef9835e594631..4ff8fee7245ba 100644
--- a/RecoHGCal/TICL/plugins/SimTrackstersProducer.cc
+++ b/RecoHGCal/TICL/plugins/SimTrackstersProducer.cc
@@ -28,6 +28,8 @@
#include "SimDataFormats/CaloAnalysis/interface/CaloParticle.h"
#include "SimDataFormats/CaloAnalysis/interface/SimCluster.h"
+#include "SimDataFormats/CaloAnalysis/interface/MtdSimTrackster.h"
+#include "SimDataFormats/CaloAnalysis/interface/MtdSimTracksterFwd.h"
#include "RecoLocalCalo/HGCalRecAlgos/interface/RecHitTools.h"
#include "SimDataFormats/Vertex/interface/SimVertex.h"
@@ -85,6 +87,7 @@ class SimTrackstersProducer : public edm::stream::EDProducer<> {
const edm::EDGetTokenT> simclusters_token_;
const edm::EDGetTokenT> caloparticles_token_;
+ const edm::EDGetTokenT MTDSimTrackstersToken_;
const edm::EDGetTokenT associatorMapSimClusterToReco_token_;
const edm::EDGetTokenT associatorMapCaloParticleToReco_token_;
@@ -114,6 +117,7 @@ SimTrackstersProducer::SimTrackstersProducer(const edm::ParameterSet& ps)
filtered_layerclusters_mask_token_(consumes(ps.getParameter("filtered_mask"))),
simclusters_token_(consumes(ps.getParameter("simclusters"))),
caloparticles_token_(consumes(ps.getParameter("caloparticles"))),
+ MTDSimTrackstersToken_(consumes(ps.getParameter("MtdSimTracksters"))),
associatorMapSimClusterToReco_token_(
consumes(ps.getParameter("layerClusterSimClusterAssociator"))),
associatorMapCaloParticleToReco_token_(
@@ -146,6 +150,7 @@ void SimTrackstersProducer::fillDescriptions(edm::ConfigurationDescriptions& des
desc.add("filtered_mask", edm::InputTag("filteredLayerClustersSimTracksters", "ticlSimTracksters"));
desc.add("simclusters", edm::InputTag("mix", "MergedCaloTruth"));
desc.add("caloparticles", edm::InputTag("mix", "MergedCaloTruth"));
+ desc.add("MtdSimTracksters", edm::InputTag("mix", "MergedMtdTruthST"));
desc.add("layerClusterSimClusterAssociator",
edm::InputTag("layerClusterSimClusterAssociationProducer"));
desc.add("layerClusterCaloParticleAssociator",
@@ -249,6 +254,9 @@ void SimTrackstersProducer::produce(edm::Event& evt, const edm::EventSetup& es)
evt.getByToken(caloparticles_token_, caloParticles_h);
const auto& caloparticles = *caloParticles_h;
+ edm::Handle MTDSimTracksters_h;
+ evt.getByToken(MTDSimTrackstersToken_, MTDSimTracksters_h);
+
const auto& simClustersToRecoColl = evt.get(associatorMapSimClusterToReco_token_);
const auto& caloParticlesToRecoColl = evt.get(associatorMapCaloParticleToReco_token_);
const auto& simVertices = evt.get(simVerticesToken_);
@@ -422,6 +430,13 @@ void SimTrackstersProducer::produce(edm::Event& evt, const edm::EventSetup& es)
edm::OrphanHandle> simTracksters_h = evt.put(std::move(result));
+ // map between simTrack and Mtd SimTracksters to loop on them only one
+ std::unordered_map SimTrackToMtdST;
+ for (unsigned int i = 0; i < MTDSimTracksters_h->size(); ++i) {
+ const auto& simTrack = (*MTDSimTracksters_h)[i].g4Tracks()[0];
+ SimTrackToMtdST[simTrack.trackId()] = &((*MTDSimTracksters_h)[i]);
+ }
+
result_ticlCandidates->resize(result_fromCP->size());
std::vector toKeep;
for (size_t i = 0; i < simTracksters_h->size(); ++i) {
@@ -457,8 +472,15 @@ void SimTrackstersProducer::produce(edm::Event& evt, const edm::EventSetup& es)
float rawEnergy = 0.f;
float regressedEnergy = 0.f;
- cand.setTime(simVertices[cp.g4Tracks()[0].vertIndex()].position().t() * pow(10, 9));
- cand.setTimeError(0);
+ const auto& simTrack = cp.g4Tracks()[0];
+ auto pos = SimTrackToMtdST.find(simTrack.trackId());
+ if (pos != SimTrackToMtdST.end()) {
+ auto MTDst = pos->second;
+ // TODO: once the associators have been implemented check if the MTDst is associated with a reco before adding the MTD time
+ cand.setMTDTime(MTDst->time(), 0);
+ }
+
+ cand.setTime(simVertices[cp.g4Tracks()[0].vertIndex()].position().t() * pow(10, 9), 0);
for (const auto& trackster : cand.tracksters()) {
rawEnergy += trackster->raw_energy();
diff --git a/RecoHGCal/TICL/plugins/TICLCandidateProducer.cc b/RecoHGCal/TICL/plugins/TICLCandidateProducer.cc
index 442a8351fef63..80eb45caad9d1 100644
--- a/RecoHGCal/TICL/plugins/TICLCandidateProducer.cc
+++ b/RecoHGCal/TICL/plugins/TICLCandidateProducer.cc
@@ -18,6 +18,7 @@
#include "DataFormats/BeamSpot/interface/BeamSpot.h"
#include "DataFormats/CaloRecHit/interface/CaloCluster.h"
#include "DataFormats/HGCalReco/interface/Common.h"
+#include "DataFormats/HGCalReco/interface/MtdHostCollection.h"
#include "DataFormats/HGCalReco/interface/TICLLayerTile.h"
#include "DataFormats/HGCalReco/interface/Trackster.h"
#include "DataFormats/TrackReco/interface/Track.h"
@@ -81,7 +82,7 @@ class TICLCandidateProducer : public edm::stream::EDProducer<> {
template
void assignTimeToCandidates(std::vector &resultCandidates,
edm::Handle> track_h,
- TICLInterpretationAlgoBase::TrackTimingInformation inputTiming,
+ edm::Handle inputTiming_h,
TrajTrackAssociationCollection trjtrks,
F func) const;
@@ -98,15 +99,13 @@ class TICLCandidateProducer : public edm::stream::EDProducer<> {
std::vector>> original_masks_tokens_;
const edm::EDGetTokenT> tracks_token_;
- edm::EDGetTokenT> tracks_time_token_;
- edm::EDGetTokenT> tracks_time_err_token_;
- edm::EDGetTokenT> tracks_beta_token_;
- edm::EDGetTokenT> tracks_path_length_token_;
- edm::EDGetTokenT> tracks_glob_pos_token_;
const edm::EDGetTokenT trajTrackAssToken_;
+ edm::EDGetTokenT inputTimingToken_;
const edm::EDGetTokenT> muons_token_;
const bool useMTDTiming_;
+ const bool useTimingAverage_;
+ const float timingQualityThreshold_;
const edm::ESGetToken geometry_token_;
const edm::ESGetToken bfield_token_;
@@ -144,8 +143,11 @@ TICLCandidateProducer::TICLCandidateProducer(const edm::ParameterSet &ps)
consumes>>(ps.getParameter("layer_clustersTime"))),
tracks_token_(consumes>(ps.getParameter("tracks"))),
trajTrackAssToken_(consumes(ps.getParameter("trjtrkAss"))),
+ inputTimingToken_(consumes(ps.getParameter("timingSoA"))),
muons_token_(consumes>(ps.getParameter("muons"))),
useMTDTiming_(ps.getParameter("useMTDTiming")),
+ useTimingAverage_(ps.getParameter("useTimingAverage")),
+ timingQualityThreshold_(ps.getParameter("timingQualityThreshold")),
geometry_token_(esConsumes()),
bfield_token_(esConsumes()),
detector_(ps.getParameter("detector")),
@@ -195,16 +197,11 @@ TICLCandidateProducer::TICLCandidateProducer(const edm::ParameterSet &ps)
original_masks_tokens_.emplace_back(consumes>(tag));
}
+ std::string detectorName_ = (detector_ == "HFNose") ? "HGCalHFNoseSensitive" : "HGCalEESensitive";
+ hdc_token_ =
+ esConsumes(edm::ESInputTag("", detectorName_));
if (useMTDTiming_) {
- std::string detectorName_ = (detector_ == "HFNose") ? "HGCalHFNoseSensitive" : "HGCalEESensitive";
- hdc_token_ = esConsumes(
- edm::ESInputTag("", detectorName_));
- tracks_time_token_ = consumes>(ps.getParameter("tracksTime"));
- tracks_time_err_token_ = consumes>(ps.getParameter("tracksTimeErr"));
- tracks_beta_token_ = consumes>(ps.getParameter("tracksBeta"));
- tracks_path_length_token_ = consumes>(ps.getParameter("tracksPathLength"));
- tracks_glob_pos_token_ =
- consumes>(ps.getParameter("tracksGlobalPosition"));
+ inputTimingToken_ = consumes(ps.getParameter("timingSoA"));
}
produces>();
@@ -235,7 +232,7 @@ void TICLCandidateProducer::beginRun(edm::Run const &iEvent, edm::EventSetup con
};
void filterTracks(edm::Handle> tkH,
- const std::vector &muons,
+ const edm::Handle> &muons_h,
const StringCutObjectSelector cutTk_,
const float tkEnergyCut_,
std::vector &maskTracks) {
@@ -245,9 +242,10 @@ void filterTracks(edm::Handle> tkH,
reco::TrackRef trackref = reco::TrackRef(tkH, i);
// veto tracks associated to muons
- int muId = PFMuonAlgo::muAssocToTrack(trackref, muons);
+ int muId = PFMuonAlgo::muAssocToTrack(trackref, *muons_h);
+ const reco::MuonRef muonref = reco::MuonRef(muons_h, muId);
- if (!cutTk_((tk)) or muId != -1) {
+ if (!cutTk_((tk)) or (muId != -1 and PFMuonAlgo::isMuon(muonref) and not(*muons_h)[muId].isTrackerMuon())) {
maskTracks[i] = false;
continue;
}
@@ -271,24 +269,18 @@ void TICLCandidateProducer::produce(edm::Event &evt, const edm::EventSetup &es)
const auto &layerClusters = evt.get(clusters_token_);
const auto &layerClustersTimes = evt.get(clustersTime_token_);
- auto const &muons = evt.get(muons_token_);
+ edm::Handle muons_h;
+ evt.getByToken(muons_token_, muons_h);
edm::Handle> tracks_h;
- const auto &trjtrks = evt.get(trajTrackAssToken_);
-
- edm::Handle> trackTime_h;
- edm::Handle> trackTimeErr_h;
- edm::Handle> trackTimeBeta_h;
- edm::Handle> trackPathToMTD_h;
- edm::Handle> trackTimeGlobalPosition_h;
evt.getByToken(tracks_token_, tracks_h);
const auto &tracks = *tracks_h;
+
+ const auto &trjtrks = evt.get(trajTrackAssToken_);
+
+ edm::Handle inputTiming_h;
if (useMTDTiming_) {
- evt.getByToken(tracks_time_token_, trackTime_h);
- evt.getByToken(tracks_time_err_token_, trackTimeErr_h);
- evt.getByToken(tracks_beta_token_, trackTimeBeta_h);
- evt.getByToken(tracks_path_length_token_, trackPathToMTD_h);
- evt.getByToken(tracks_glob_pos_token_, trackTimeGlobalPosition_h);
+ evt.getByToken(inputTimingToken_, inputTiming_h);
}
const auto &bs = evt.get(bsToken_);
@@ -333,7 +325,7 @@ void TICLCandidateProducer::produce(edm::Event &evt, const edm::EventSetup &es)
std::vector maskTracks;
maskTracks.resize(tracks.size());
- filterTracks(tracks_h, muons, cutTk_, tkEnergyCut_, maskTracks);
+ filterTracks(tracks_h, muons_h, cutTk_, tkEnergyCut_, maskTracks);
const typename TICLInterpretationAlgoBase::Inputs input(evt,
es,
@@ -344,9 +336,6 @@ void TICLCandidateProducer::produce(edm::Event &evt, const edm::EventSetup &es)
tracks_h,
maskTracks);
- const typename TICLInterpretationAlgoBase::TrackTimingInformation inputTiming(
- trackTime_h, trackTimeErr_h, trackTimeBeta_h, trackPathToMTD_h, trackTimeGlobalPosition_h);
-
auto resultCandidates = std::make_unique>();
std::vector trackstersInTrackIndices(tracks.size(), -1);
@@ -354,7 +343,7 @@ void TICLCandidateProducer::produce(edm::Event &evt, const edm::EventSetup &es)
//egammaInterpretationAlg_->makecandidates(inputGSF, inputTiming, *resultTrackstersMerged, trackstersInGSFTrackIndices)
// mask generalTracks associated to GSFTrack linked in egammaInterpretationAlgo_
- generalInterpretationAlgo_->makeCandidates(input, inputTiming, *resultTracksters, trackstersInTrackIndices);
+ generalInterpretationAlgo_->makeCandidates(input, inputTiming_h, *resultTracksters, trackstersInTrackIndices);
assignPCAtoTracksters(
*resultTracksters, layerClusters, layerClustersTimes, rhtools_.getPositionLayer(rhtools_.lastLayerEE()).z(), true);
@@ -376,6 +365,13 @@ void TICLCandidateProducer::produce(edm::Event &evt, const edm::EventSetup &es)
//charged candidates track only
edm::Ptr tracksterPtr;
TICLCandidate chargedCandidate(trackPtr, tracksterPtr);
+ auto trackRef = edm::Ref(tracks_h, iTrack);
+ const int muId = PFMuonAlgo::muAssocToTrack(trackRef, *muons_h);
+ const reco::MuonRef muonRef = reco::MuonRef(muons_h, muId);
+ if (muonRef.isNonnull() and muonRef->isGlobalMuon()) {
+ // create muon candidate
+ chargedCandidate.setPdgId(13 * trackPtr.get()->charge());
+ }
resultCandidates->push_back(chargedCandidate);
}
}
@@ -442,7 +438,7 @@ void TICLCandidateProducer::produce(edm::Event &evt, const edm::EventSetup &es)
return 0.f;
};
- assignTimeToCandidates(*resultCandidates, tracks_h, inputTiming, trjtrks, getPathLength);
+ assignTimeToCandidates(*resultCandidates, tracks_h, inputTiming_h, trjtrks, getPathLength);
evt.put(std::move(resultCandidates));
}
@@ -574,19 +570,17 @@ void TICLCandidateProducer::energyRegressionAndID(const std::vector
-void TICLCandidateProducer::assignTimeToCandidates(
- std::vector &resultCandidates,
- edm::Handle> track_h,
- TICLInterpretationAlgoBase::TrackTimingInformation inputTiming,
- TrajTrackAssociationCollection trjtrks,
- F func) const {
+void TICLCandidateProducer::assignTimeToCandidates(std::vector &resultCandidates,
+ edm::Handle> track_h,
+ edm::Handle inputTiming_h,
+ TrajTrackAssociationCollection trjtrks,
+ F func) const {
for (auto &cand : resultCandidates) {
float beta = 1;
float time = 0.f;
float invTimeErr = 0.f;
+ float timeErr = -1.f;
- // if (not cand.tracksters().size())
- // continue;
for (const auto &tr : cand.tracksters()) {
if (tr->timeError() > 0) {
const auto invTimeESq = pow(tr->timeError(), -2);
@@ -596,16 +590,17 @@ void TICLCandidateProducer::assignTimeToCandidates(
auto path = std::sqrt(x * x + y * y + z * z);
if (cand.trackPtr().get() != nullptr) {
const auto &trackIndex = cand.trackPtr().get() - (edm::Ptr(track_h, 0)).get();
- const auto &trackRef = edm::Ref>(track_h, trackIndex);
- if (useMTDTiming_ and (*inputTiming.tkTimeErr_h)[trackRef] > 0) {
- const auto &trackMtdPos = (*inputTiming.tkMtdPos_h);
- const auto xMtd = trackMtdPos[trackRef].x();
- const auto yMtd = trackMtdPos[trackRef].y();
- const auto zMtd = trackMtdPos[trackRef].z();
-
- beta = (*inputTiming.tkBeta_h)[trackRef];
- path = std::sqrt((x - xMtd) * (x - xMtd) + (y - yMtd) * (y - yMtd) + (z - zMtd) * (z - zMtd)) +
- (*inputTiming.tkPath_h)[trackRef];
+ if (useMTDTiming_) {
+ auto const &inputTimingView = (*inputTiming_h).const_view();
+ if (inputTimingView.timeErr()[trackIndex] > 0) {
+ const auto xMtd = inputTimingView.posInMTD_x()[trackIndex];
+ const auto yMtd = inputTimingView.posInMTD_y()[trackIndex];
+ const auto zMtd = inputTimingView.posInMTD_z()[trackIndex];
+
+ beta = inputTimingView.beta()[trackIndex];
+ path = std::sqrt((x - xMtd) * (x - xMtd) + (y - yMtd) * (y - yMtd) + (z - zMtd) * (z - zMtd)) +
+ inputTimingView.pathLength()[trackIndex];
+ }
} else {
const auto &trackIndex = cand.trackPtr().get() - (edm::Ptr(track_h, 0)).get();
for (const auto &trj : trjtrks) {
@@ -627,10 +622,37 @@ void TICLCandidateProducer::assignTimeToCandidates(
}
}
if (invTimeErr > 0) {
- cand.setTime(time / invTimeErr);
+ time = time / invTimeErr;
// FIXME_ set a liminf of 0.02 ns on the ts error (based on residuals)
- auto timeErr = sqrt(1.f / invTimeErr) > 0.02 ? sqrt(1.f / invTimeErr) : 0.02;
- cand.setTimeError(timeErr);
+ timeErr = sqrt(1.f / invTimeErr) > 0.02 ? sqrt(1.f / invTimeErr) : 0.02;
+ cand.setTime(time, timeErr);
+ }
+
+ if (useMTDTiming_ and cand.charge()) {
+ // Check MTD timing availability
+ auto const &inputTimingView = (*inputTiming_h).const_view();
+ const auto &trackIndex = cand.trackPtr().get() - (edm::Ptr(track_h, 0)).get();
+ const bool assocQuality = inputTimingView.MVAquality()[trackIndex] > timingQualityThreshold_;
+ if (assocQuality) {
+ const auto timeHGC = cand.time();
+ const auto timeEHGC = cand.timeError();
+ const auto timeMTD = inputTimingView.time0()[trackIndex];
+ const auto timeEMTD = inputTimingView.time0Err()[trackIndex];
+
+ if (useTimingAverage_ && (timeEMTD > 0 && timeEHGC > 0)) {
+ // Compute weighted average between HGCAL and MTD timing
+ const auto invTimeESqHGC = pow(timeEHGC, -2);
+ const auto invTimeESqMTD = pow(timeEMTD, -2);
+ timeErr = 1.f / (invTimeESqHGC + invTimeESqMTD);
+ time = (timeHGC * invTimeESqHGC + timeMTD * invTimeESqMTD) * timeErr;
+ timeErr = sqrt(timeErr);
+ } else if (timeEMTD > 0) {
+ time = timeMTD;
+ timeErr = timeEMTD;
+ }
+ }
+ cand.setTime(time, timeErr);
+ cand.setMTDTime(inputTimingView.time()[trackIndex], inputTimingView.timeErr()[trackIndex]);
}
}
}
@@ -650,16 +672,14 @@ void TICLCandidateProducer::fillDescriptions(edm::ConfigurationDescriptions &des
desc.add("layer_clustersTime", edm::InputTag("hgcalMergeLayerClusters", "timeLayerCluster"));
desc.add("tracks", edm::InputTag("generalTracks"));
desc.add("trjtrkAss", edm::InputTag("generalTracks"));
- desc.add("tracksTime", edm::InputTag("trackExtenderWithMTD:generalTracktmtd"));
- desc.add("tracksTimeErr", edm::InputTag("trackExtenderWithMTD:generalTracksigmatmtd"));
- desc.add("tracksBeta", edm::InputTag("trackExtenderWithMTD:generalTrackBeta"));
- desc.add("tracksGlobalPosition", edm::InputTag("trackExtenderWithMTD:generalTrackmtdpos"));
- desc.add("tracksPathLength", edm::InputTag("trackExtenderWithMTD:generalTrackPathLength"));
+ desc.add("timingSoA", edm::InputTag("mtdSoA"));
desc.add("muons", edm::InputTag("muons1stStep"));
desc.add("detector", "HGCAL");
desc.add("propagator", "PropagatorWithMaterial");
desc.add("beamspot", edm::InputTag("offlineBeamSpot"));
desc.add("useMTDTiming", true);
+ desc.add("useTimingAverage", true);
+ desc.add("timingQualityThreshold", 0.5);
desc.add("tfDnnLabel", "tracksterSelectionTf");
desc.add("eid_input_name", "input");
desc.add("eid_output_name_energy", "output/regressed_energy");
diff --git a/RecoHGCal/TICL/plugins/TICLDumper.cc b/RecoHGCal/TICL/plugins/TICLDumper.cc
index fe60e71d2ea8d..8dc2cb914e98a 100644
--- a/RecoHGCal/TICL/plugins/TICLDumper.cc
+++ b/RecoHGCal/TICL/plugins/TICLDumper.cc
@@ -24,6 +24,7 @@
#include "DataFormats/CaloRecHit/interface/CaloCluster.h"
#include "DataFormats/HGCalReco/interface/Trackster.h"
#include "DataFormats/HGCalReco/interface/TICLCandidate.h"
+#include "DataFormats/MuonReco/interface/Muon.h"
#include "DataFormats/TrackReco/interface/Track.h"
#include "DataFormats/TrackReco/interface/TrackFwd.h"
#include "DataFormats/DetId/interface/DetId.h"
@@ -34,6 +35,7 @@
#include "SimDataFormats/CaloAnalysis/interface/SimCluster.h"
#include "RecoLocalCalo/HGCalRecAlgos/interface/RecHitTools.h"
+#include "RecoParticleFlow/PFProducer/interface/PFMuonAlgo.h"
#include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
#include "TrackingTools/GeomPropagators/interface/Propagator.h"
#include "TrackingTools/Records/interface/TrackingComponentsRecord.h"
@@ -80,6 +82,7 @@ class TICLDumper : public edm::one::EDAnalyzer> tracksters_token_;
+ const edm::EDGetTokenT> tracksters_in_candidate_token_;
const edm::EDGetTokenT> layer_clusters_token_;
const edm::EDGetTokenT> ticl_candidates_token_;
const edm::EDGetTokenT> tracks_token_;
@@ -100,6 +103,7 @@ class TICLDumper : public edm::one::EDAnalyzer> hgcaltracks_py_token_;
const edm::EDGetTokenT