File indexing completed on 2021-06-16 03:20:14
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "DataFormats/EgammaReco/interface/ElectronSeed.h"
0010 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0011 #include "DataFormats/MuonReco/interface/Muon.h"
0012 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0013 #include "DataFormats/ParticleFlowReco/interface/PFRecTrack.h"
0014 #include "DataFormats/TrackReco/interface/Track.h"
0015 #include "DataFormats/VertexReco/interface/Vertex.h"
0016 #include "FWCore/Framework/interface/ESHandle.h"
0017 #include "FWCore/Framework/interface/Event.h"
0018 #include "FWCore/Framework/interface/stream/EDProducer.h"
0019 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0020 #include "MagneticField/Engine/interface/MagneticField.h"
0021 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0022 #include "RecoParticleFlow/PFTracking/interface/PFTrackTransformer.h"
0023 #include "TrackingTools/IPTools/interface/IPTools.h"
0024 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0025 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0026 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0027
0028 class PFTrackProducer : public edm::stream::EDProducer<> {
0029 public:
0030
0031 explicit PFTrackProducer(const edm::ParameterSet&);
0032
0033 private:
0034 void beginRun(const edm::Run&, const edm::EventSetup&) override;
0035 void endRun(const edm::Run&, const edm::EventSetup&) override;
0036
0037
0038 void produce(edm::Event&, const edm::EventSetup&) override;
0039
0040 const edm::ESGetToken<TransientTrackBuilder, TransientTrackRecord> transientTrackToken_;
0041 const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magneticFieldToken_;
0042
0043
0044 std::unique_ptr<PFTrackTransformer> pfTransformer_;
0045 std::vector<edm::EDGetTokenT<reco::TrackCollection>> tracksContainers_;
0046 std::vector<edm::EDGetTokenT<std::vector<Trajectory>>> trajContainers_;
0047 edm::EDGetTokenT<reco::GsfTrackCollection> gsfTrackLabel_;
0048 edm::EDGetTokenT<reco::MuonCollection> muonColl_;
0049 edm::EDGetTokenT<reco::VertexCollection> vtx_h;
0050
0051 bool useQuality_;
0052 reco::TrackBase::TrackQuality trackQuality_;
0053 bool trajinev_;
0054 bool gsfinev_;
0055 };
0056
0057 #include "FWCore/Framework/interface/MakerMacros.h"
0058 DEFINE_FWK_MODULE(PFTrackProducer);
0059
0060 using namespace std;
0061 using namespace edm;
0062 using namespace reco;
0063 PFTrackProducer::PFTrackProducer(const ParameterSet& iConfig)
0064 : transientTrackToken_(esConsumes(edm::ESInputTag("", "TransientTrackBuilder"))),
0065 magneticFieldToken_(esConsumes<edm::Transition::BeginRun>()),
0066 pfTransformer_() {
0067 produces<reco::PFRecTrackCollection>();
0068
0069 std::vector<InputTag> tags = iConfig.getParameter<vector<InputTag>>("TkColList");
0070 trajinev_ = iConfig.getParameter<bool>("TrajInEvents");
0071 tracksContainers_.reserve(tags.size());
0072 if (trajinev_) {
0073 trajContainers_.reserve(tags.size());
0074 }
0075 for (auto const& tag : tags) {
0076 tracksContainers_.push_back(consumes<reco::TrackCollection>(tag));
0077 if (trajinev_) {
0078 trajContainers_.push_back(consumes<std::vector<Trajectory>>(tag));
0079 }
0080 }
0081
0082 useQuality_ = iConfig.getParameter<bool>("UseQuality");
0083
0084 gsfinev_ = iConfig.getParameter<bool>("GsfTracksInEvents");
0085 if (gsfinev_) {
0086 gsfTrackLabel_ = consumes<reco::GsfTrackCollection>(iConfig.getParameter<InputTag>("GsfTrackModuleLabel"));
0087 }
0088
0089 trackQuality_ = reco::TrackBase::qualityByName(iConfig.getParameter<std::string>("TrackQuality"));
0090
0091 muonColl_ = consumes<reco::MuonCollection>(iConfig.getParameter<InputTag>("MuColl"));
0092
0093 vtx_h = consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("PrimaryVertexLabel"));
0094 }
0095
0096 void PFTrackProducer::produce(Event& iEvent, const EventSetup& iSetup) {
0097
0098 auto PfTrColl = std::make_unique<reco::PFRecTrackCollection>();
0099
0100
0101 Handle<GsfTrackCollection> gsftrackcoll;
0102 bool foundgsf = false;
0103 if (gsfinev_) {
0104 foundgsf = iEvent.getByToken(gsfTrackLabel_, gsftrackcoll);
0105 }
0106
0107
0108 Handle<reco::VertexCollection> vertex;
0109 iEvent.getByToken(vtx_h, vertex);
0110 reco::Vertex dummy;
0111 const reco::Vertex* pv = &dummy;
0112 if (vertex.isValid()) {
0113 pv = &*vertex->begin();
0114 } else {
0115 reco::Vertex::Error e;
0116 e(0, 0) = 0.0015 * 0.0015;
0117 e(1, 1) = 0.0015 * 0.0015;
0118 e(2, 2) = 15. * 15.;
0119 reco::Vertex::Point p(0, 0, 0);
0120 dummy = reco::Vertex(p, e, 0, 0, 0);
0121 }
0122
0123
0124 TransientTrackBuilder const& thebuilder = iSetup.getData(transientTrackToken_);
0125
0126
0127 Handle<reco::MuonCollection> recMuons;
0128 iEvent.getByToken(muonColl_, recMuons);
0129
0130
0131 const vector<Trajectory> dummyTj(0);
0132
0133 for (unsigned int istr = 0; istr < tracksContainers_.size(); istr++) {
0134
0135 Handle<reco::TrackCollection> tkRefCollection;
0136 iEvent.getByToken(tracksContainers_[istr], tkRefCollection);
0137 reco::TrackCollection Tk = *(tkRefCollection.product());
0138
0139
0140 const vector<Trajectory>* Tj = &dummyTj;
0141 if (trajinev_) {
0142
0143 Handle<vector<Trajectory>> tjCollection;
0144 iEvent.getByToken(trajContainers_[istr], tjCollection);
0145
0146 Tj = tjCollection.product();
0147 }
0148
0149 for (unsigned int i = 0; i < Tk.size(); i++) {
0150 reco::TrackRef trackRef(tkRefCollection, i);
0151
0152 if (useQuality_ && (!(Tk[i].quality(trackQuality_)))) {
0153 bool isMuCandidate = false;
0154
0155
0156
0157 if (recMuons.isValid()) {
0158 for (unsigned j = 0; j < recMuons->size(); j++) {
0159 reco::MuonRef muonref(recMuons, j);
0160 if (muonref->track().isNonnull())
0161 if (muonref->track() == trackRef && muonref->isGlobalMuon()) {
0162 isMuCandidate = true;
0163
0164 break;
0165 }
0166 }
0167 }
0168 if (!isMuCandidate) {
0169 continue;
0170 }
0171 }
0172
0173
0174 bool preId = false;
0175 if (foundgsf) {
0176
0177 for (auto const& gsfTrack : *gsftrackcoll) {
0178 if (gsfTrack.seedRef().isNull())
0179 continue;
0180 auto const& seed = *(gsfTrack.extra()->seedRef());
0181 auto const& ElSeed = dynamic_cast<ElectronSeed const&>(seed);
0182 if (ElSeed.ctfTrack().isNonnull()) {
0183 if (ElSeed.ctfTrack() == trackRef) {
0184 preId = true;
0185 break;
0186 }
0187 }
0188 }
0189 }
0190 if (preId) {
0191
0192 reco::PFRecTrack pftrack(trackRef->charge(), reco::PFRecTrack::KF_ELCAND, i, trackRef);
0193
0194 bool valid = false;
0195 if (trajinev_) {
0196 valid = pfTransformer_->addPoints(pftrack, *trackRef, (*Tj)[i]);
0197 } else {
0198 Trajectory FakeTraj;
0199 valid = pfTransformer_->addPoints(pftrack, *trackRef, FakeTraj);
0200 }
0201 if (valid) {
0202
0203 double stip = -999;
0204 const reco::PFTrajectoryPoint& atECAL = pftrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALEntrance);
0205 if (atECAL.isValid())
0206 {
0207 GlobalVector direction(pftrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALEntrance).position().x(),
0208 pftrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALEntrance).position().y(),
0209 pftrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALEntrance).position().z());
0210 stip = IPTools::signedTransverseImpactParameter(thebuilder.build(*trackRef), direction, *pv)
0211 .second.significance();
0212 }
0213 pftrack.setSTIP(stip);
0214 PfTrColl->push_back(pftrack);
0215 }
0216 } else {
0217 reco::PFRecTrack pftrack(trackRef->charge(), reco::PFRecTrack::KF, i, trackRef);
0218 bool valid = false;
0219 if (trajinev_) {
0220 valid = pfTransformer_->addPoints(pftrack, *trackRef, (*Tj)[i]);
0221 } else {
0222 Trajectory FakeTraj;
0223 valid = pfTransformer_->addPoints(pftrack, *trackRef, FakeTraj);
0224 }
0225
0226 if (valid) {
0227 double stip = -999;
0228 const reco::PFTrajectoryPoint& atECAL = pftrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALEntrance);
0229 if (atECAL.isValid()) {
0230 GlobalVector direction(pftrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALEntrance).position().x(),
0231 pftrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALEntrance).position().y(),
0232 pftrack.extrapolatedPoint(reco::PFTrajectoryPoint::ECALEntrance).position().z());
0233 stip = IPTools::signedTransverseImpactParameter(thebuilder.build(*trackRef), direction, *pv)
0234 .second.significance();
0235 }
0236 pftrack.setSTIP(stip);
0237 PfTrColl->push_back(pftrack);
0238 }
0239 }
0240 }
0241 }
0242 iEvent.put(std::move(PfTrColl));
0243 }
0244
0245
0246 void PFTrackProducer::beginRun(const edm::Run& run, const EventSetup& iSetup) {
0247 auto const& magneticField = iSetup.getData(magneticFieldToken_);
0248 pfTransformer_ = std::make_unique<PFTrackTransformer>(math::XYZVector(magneticField.inTesla(GlobalPoint(0, 0, 0))));
0249 if (!trajinev_)
0250 pfTransformer_->OnlyProp();
0251 }
0252
0253
0254 void PFTrackProducer::endRun(const edm::Run& run, const EventSetup& iSetup) { pfTransformer_.reset(); }