File indexing completed on 2022-12-04 23:58:00
0001
0002
0003 #include "FWCore/Framework/interface/Frameworkfwd.h"
0004 #include "FWCore/Framework/interface/stream/EDProducer.h"
0005
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "FWCore/Framework/interface/MakerMacros.h"
0008
0009 #include "DataFormats/Common/interface/View.h"
0010
0011 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
0012 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0013
0014 #include "DataFormats/HGCalReco/interface/TICLCandidate.h"
0015
0016 #include "RecoParticleFlow/PFProducer/interface/PFMuonAlgo.h"
0017
0018 class PFTICLProducer : public edm::stream::EDProducer<> {
0019 public:
0020 PFTICLProducer(const edm::ParameterSet&);
0021 ~PFTICLProducer() override {}
0022
0023 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0024
0025 void produce(edm::Event&, const edm::EventSetup&) override;
0026
0027 private:
0028
0029 const bool useMTDTiming_;
0030 const bool useTimingAverage_;
0031 const float timingQualityThreshold_;
0032 const bool energy_from_regression_;
0033
0034 const edm::EDGetTokenT<edm::View<TICLCandidate>> ticl_candidates_;
0035 const edm::EDGetTokenT<edm::ValueMap<float>> srcTrackTime_, srcTrackTimeError_, srcTrackTimeQuality_;
0036 const edm::EDGetTokenT<reco::MuonCollection> muons_;
0037
0038 std::unique_ptr<PFMuonAlgo> pfmu_;
0039 };
0040
0041 DEFINE_FWK_MODULE(PFTICLProducer);
0042
0043 PFTICLProducer::PFTICLProducer(const edm::ParameterSet& conf)
0044 : useMTDTiming_(conf.getParameter<bool>("useMTDTiming")),
0045 useTimingAverage_(conf.getParameter<bool>("useTimingAverage")),
0046 timingQualityThreshold_(conf.getParameter<double>("timingQualityThreshold")),
0047 energy_from_regression_(conf.getParameter<bool>("energyFromRegression")),
0048 ticl_candidates_(consumes<edm::View<TICLCandidate>>(conf.getParameter<edm::InputTag>("ticlCandidateSrc"))),
0049 srcTrackTime_(consumes<edm::ValueMap<float>>(conf.getParameter<edm::InputTag>("trackTimeValueMap"))),
0050 srcTrackTimeError_(consumes<edm::ValueMap<float>>(conf.getParameter<edm::InputTag>("trackTimeErrorMap"))),
0051 srcTrackTimeQuality_(consumes<edm::ValueMap<float>>(conf.getParameter<edm::InputTag>("trackTimeQualityMap"))),
0052 muons_(consumes<reco::MuonCollection>(conf.getParameter<edm::InputTag>("muonSrc"))),
0053 pfmu_(std::make_unique<PFMuonAlgo>(conf.getParameterSet("pfMuonAlgoParameters"),
0054 false)) {
0055 produces<reco::PFCandidateCollection>();
0056 }
0057
0058 void PFTICLProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0059 edm::ParameterSetDescription desc;
0060 desc.add<edm::InputTag>("ticlCandidateSrc", edm::InputTag("ticlTrackstersMerge"));
0061 desc.add<edm::InputTag>("trackTimeValueMap", edm::InputTag("tofPID:t0"));
0062 desc.add<edm::InputTag>("trackTimeErrorMap", edm::InputTag("tofPID:sigmat0"));
0063 desc.add<edm::InputTag>("trackTimeQualityMap", edm::InputTag("mtdTrackQualityMVA:mtdQualMVA"));
0064 desc.add<bool>("energyFromRegression", true);
0065 desc.add<double>("timingQualityThreshold", 0.5);
0066 desc.add<bool>("useMTDTiming", true);
0067 desc.add<bool>("useTimingAverage", false);
0068
0069 desc.add<edm::InputTag>("muonSrc", edm::InputTag("muons1stStep"));
0070 edm::ParameterSetDescription psd_PFMuonAlgo;
0071 PFMuonAlgo::fillPSetDescription(psd_PFMuonAlgo);
0072 desc.add<edm::ParameterSetDescription>("pfMuonAlgoParameters", psd_PFMuonAlgo);
0073
0074 descriptions.add("pfTICLProducer", desc);
0075 }
0076
0077 void PFTICLProducer::produce(edm::Event& evt, const edm::EventSetup& es) {
0078
0079 edm::Handle<edm::View<TICLCandidate>> ticl_cand_h;
0080 evt.getByToken(ticl_candidates_, ticl_cand_h);
0081 const auto ticl_candidates = *ticl_cand_h;
0082 edm::Handle<edm::ValueMap<float>> trackTimeH, trackTimeErrH, trackTimeQualH;
0083 evt.getByToken(srcTrackTime_, trackTimeH);
0084 evt.getByToken(srcTrackTimeError_, trackTimeErrH);
0085 evt.getByToken(srcTrackTimeQuality_, trackTimeQualH);
0086 const auto muonH = evt.getHandle(muons_);
0087 const auto& muons = *muonH;
0088
0089 auto candidates = std::make_unique<reco::PFCandidateCollection>();
0090
0091 for (const auto& ticl_cand : ticl_candidates) {
0092 const auto abs_pdg_id = std::abs(ticl_cand.pdgId());
0093 const auto charge = ticl_cand.charge();
0094 const auto& four_mom = ticl_cand.p4();
0095 float total_raw_energy = 0.f;
0096 float total_em_raw_energy = 0.f;
0097 for (const auto& t : ticl_cand.tracksters()) {
0098 total_raw_energy += t->raw_energy();
0099 total_em_raw_energy += t->raw_em_energy();
0100 }
0101 float ecal_energy_fraction = total_em_raw_energy / total_raw_energy;
0102 float ecal_energy = energy_from_regression_ ? ticl_cand.p4().energy() * ecal_energy_fraction
0103 : ticl_cand.rawEnergy() * ecal_energy_fraction;
0104 float hcal_energy =
0105 energy_from_regression_ ? ticl_cand.p4().energy() - ecal_energy : ticl_cand.rawEnergy() - ecal_energy;
0106
0107 hcal_energy = std::max(0.f, hcal_energy);
0108 reco::PFCandidate::ParticleType part_type;
0109 switch (abs_pdg_id) {
0110 case 11:
0111 part_type = reco::PFCandidate::e;
0112 break;
0113 case 13:
0114 part_type = reco::PFCandidate::mu;
0115 break;
0116 case 22:
0117 part_type = reco::PFCandidate::gamma;
0118 break;
0119 case 130:
0120 part_type = reco::PFCandidate::h0;
0121 break;
0122 case 211:
0123 part_type = reco::PFCandidate::h;
0124 break;
0125
0126 default:
0127 part_type = reco::PFCandidate::X;
0128 }
0129
0130 candidates->emplace_back(charge, four_mom, part_type);
0131
0132 auto& candidate = candidates->back();
0133 candidate.setEcalEnergy(ecal_energy, ecal_energy);
0134 candidate.setHcalEnergy(hcal_energy, hcal_energy);
0135 if (candidate.charge()) {
0136
0137
0138 reco::TrackRef trackref(ticl_cand.trackPtr().id(), int(ticl_cand.trackPtr().key()), &evt.productGetter());
0139 candidate.setTrackRef(trackref);
0140
0141 const int muId = PFMuonAlgo::muAssocToTrack(trackref, muons);
0142 if (muId != -1) {
0143 const reco::MuonRef muonref = reco::MuonRef(muonH, muId);
0144 const bool allowLoose = (part_type == reco::PFCandidate::mu);
0145
0146 pfmu_->reconstructMuon(candidate, muonref, allowLoose);
0147 }
0148 }
0149
0150
0151 auto time = ticl_cand.time();
0152 auto timeE = ticl_cand.timeError();
0153
0154 if (useMTDTiming_ and candidate.charge()) {
0155
0156 time = -99.;
0157 timeE = -1.;
0158
0159 const bool assocQuality = (*trackTimeQualH)[candidate.trackRef()] > timingQualityThreshold_;
0160 if (assocQuality) {
0161 const auto timeHGC = time;
0162 const auto timeEHGC = timeE;
0163 const auto timeMTD = (*trackTimeH)[candidate.trackRef()];
0164 const auto timeEMTD = (*trackTimeErrH)[candidate.trackRef()];
0165
0166 if (useTimingAverage_ && (timeEMTD > 0 && timeEHGC > 0)) {
0167
0168 const auto invTimeESqHGC = pow(timeEHGC, -2);
0169 const auto invTimeESqMTD = pow(timeEMTD, -2);
0170 timeE = (invTimeESqHGC * invTimeESqMTD) / (invTimeESqHGC + invTimeESqMTD);
0171 time = (timeHGC * invTimeESqHGC + timeMTD * invTimeESqMTD) * timeE;
0172 timeE = sqrt(timeE);
0173 } else if (timeEMTD > 0) {
0174 time = timeMTD;
0175 timeE = timeEMTD;
0176 }
0177 }
0178 }
0179 candidate.setTime(time, timeE);
0180 }
0181
0182 evt.put(std::move(candidates));
0183 }