File indexing completed on 2024-04-06 12:28:07
0001 #include "RecoTracker/FinalTrackSelectors/interface/TrackCollectionCloner.h"
0002
0003 #include "FWCore/Framework/interface/global/EDProducer.h"
0004 #include "FWCore/Utilities/interface/InputTag.h"
0005
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0008 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0009
0010 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0011 #include "DataFormats/TrackReco/interface/Track.h"
0012
0013 #include "DataFormats/TrackCandidate/interface/TrackCandidate.h"
0014
0015 #include <vector>
0016 #include <algorithm>
0017 #include <string>
0018 #include <iostream>
0019 #include <memory>
0020
0021 using namespace reco;
0022
0023 namespace {
0024 class TrackCollectionFilterCloner final : public edm::global::EDProducer<> {
0025 public:
0026
0027 explicit TrackCollectionFilterCloner(const edm::ParameterSet& iConfig);
0028
0029 ~TrackCollectionFilterCloner() override;
0030
0031
0032 using CandidateToDuplicate = std::vector<std::pair<int, int>>;
0033
0034 using RecHitContainer = edm::OwnVector<TrackingRecHit>;
0035
0036 using MVACollection = std::vector<float>;
0037 using QualityMaskCollection = std::vector<unsigned char>;
0038
0039 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0040
0041 private:
0042
0043 void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0044
0045 private:
0046 TrackCollectionCloner collectionCloner;
0047 const TrackCollectionCloner::Tokens originalTrackSource_;
0048
0049 const edm::EDGetTokenT<MVACollection> originalMVAValsToken_;
0050 const edm::EDGetTokenT<QualityMaskCollection> originalQualValsToken_;
0051
0052 const reco::TrackBase::TrackQuality minQuality_;
0053 };
0054 }
0055
0056 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0057 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
0058 #include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h"
0059 #include "TrackingTools/PatternTools/interface/ClusterRemovalRefSetter.h"
0060
0061 #include "FWCore/Framework/interface/Event.h"
0062 #include "CommonTools/Statistics/interface/ChiSquaredProbability.h"
0063
0064 #include <tuple>
0065 #include <array>
0066 #include "CommonTools/Utils/interface/DynArray.h"
0067
0068 namespace {
0069 void TrackCollectionFilterCloner::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0070 edm::ParameterSetDescription desc;
0071 desc.add<edm::InputTag>("originalSource", edm::InputTag());
0072 desc.add<edm::InputTag>("originalMVAVals", edm::InputTag());
0073 desc.add<edm::InputTag>("originalQualVals", edm::InputTag());
0074 desc.add<std::string>("minQuality", "loose");
0075 TrackCollectionCloner::fill(desc);
0076 descriptions.add("TrackCollectionFilterCloner", desc);
0077 }
0078
0079 TrackCollectionFilterCloner::TrackCollectionFilterCloner(const edm::ParameterSet& iConfig)
0080 : collectionCloner(producesCollector(), iConfig, true),
0081 originalTrackSource_(iConfig.getParameter<edm::InputTag>("originalSource"), consumesCollector()),
0082 originalMVAValsToken_(consumes<MVACollection>(iConfig.getParameter<edm::InputTag>("originalMVAVals"))),
0083 originalQualValsToken_(
0084 consumes<QualityMaskCollection>(iConfig.getParameter<edm::InputTag>("originalQualVals"))),
0085 minQuality_(reco::TrackBase::qualityByName(iConfig.getParameter<std::string>("minQuality"))) {
0086 produces<MVACollection>("MVAValues");
0087 produces<QualityMaskCollection>("QualityMasks");
0088 }
0089
0090 TrackCollectionFilterCloner::~TrackCollectionFilterCloner() {}
0091
0092 void TrackCollectionFilterCloner::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0093 TrackCollectionCloner::Producer producer(iEvent, collectionCloner);
0094
0095
0096 auto const& originalsTracks = originalTrackSource_.tracks(iEvent);
0097
0098 auto nTracks = originalsTracks.size();
0099
0100 edm::Handle<MVACollection> originalMVAStore;
0101 iEvent.getByToken(originalMVAValsToken_, originalMVAStore);
0102 assert((*originalMVAStore).size() == nTracks);
0103
0104 edm::Handle<QualityMaskCollection> originalQualStore;
0105 iEvent.getByToken(originalQualValsToken_, originalQualStore);
0106 assert((*originalQualStore).size() == nTracks);
0107
0108
0109 unsigned char qualMask = ~0;
0110 if (minQuality_ != reco::TrackBase::undefQuality)
0111 qualMask = 1 << minQuality_;
0112
0113
0114 std::vector<unsigned int> selId;
0115 auto pmvas = std::make_unique<MVACollection>();
0116 auto pquals = std::make_unique<QualityMaskCollection>();
0117
0118 auto k = 0U;
0119 for (auto j = 0U; j < nTracks; ++j) {
0120 if (!(qualMask & (*originalQualStore)[j]))
0121 continue;
0122
0123 selId.push_back(j);
0124 pmvas->push_back((*originalMVAStore)[j]);
0125 pquals->push_back((*originalQualStore)[j]);
0126
0127 ++k;
0128 }
0129
0130
0131 auto nsel = k;
0132 auto isel = 0U;
0133 assert(producer.selTracks_->size() == isel);
0134 producer(originalTrackSource_, selId);
0135 assert(producer.selTracks_->size() == nsel);
0136
0137 for (; isel < nsel; ++isel) {
0138 auto& otk = (*producer.selTracks_)[isel];
0139 otk.setQualityMask((*pquals)[isel]);
0140 }
0141 assert(producer.selTracks_->size() == pmvas->size());
0142 assert(producer.selTracks_->size() == pquals->size());
0143
0144 iEvent.put(std::move(pmvas), "MVAValues");
0145 iEvent.put(std::move(pquals), "QualityMasks");
0146 }
0147
0148 }
0149
0150 #include "FWCore/PluginManager/interface/ModuleDef.h"
0151 #include "FWCore/Framework/interface/MakerMacros.h"
0152
0153 DEFINE_FWK_MODULE(TrackCollectionFilterCloner);