File indexing completed on 2022-10-14 01:44:38
0001 #include "FWCore/Framework/interface/stream/EDProducer.h"
0002 #include "FWCore/Framework/interface/Run.h"
0003 #include "FWCore/Framework/interface/EventSetup.h"
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "FWCore/Utilities/interface/InputTag.h"
0008 #include "FWCore/Utilities/interface/EDGetToken.h"
0009 #include "FWCore/Utilities/interface/EDPutToken.h"
0010 #include "DataFormats/Common/interface/Handle.h"
0011
0012 #include "SimTracker/TrackTriggerAssociation/interface/TTTypes.h"
0013 #include "SimTracker/TrackTriggerAssociation/interface/StubAssociation.h"
0014 #include "L1Trigger/TrackTrigger/interface/Setup.h"
0015
0016 #include <vector>
0017 #include <map>
0018 #include <set>
0019 #include <algorithm>
0020 #include <iterator>
0021 #include <cmath>
0022
0023 using namespace std;
0024 using namespace edm;
0025
0026 namespace tt {
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037 class StubAssociator : public stream::EDProducer<> {
0038 public:
0039 explicit StubAssociator(const ParameterSet&);
0040 ~StubAssociator() override {}
0041
0042 private:
0043 void beginRun(const Run&, const EventSetup&) override;
0044 void produce(Event&, const EventSetup&) override;
0045 void endJob() {}
0046
0047
0048 const Setup* setup_ = nullptr;
0049
0050 EDGetTokenT<TTStubDetSetVec> getTokenTTStubDetSetVec_;
0051
0052 EDGetTokenT<TTClusterAssMap> getTokenTTClusterAssMap_;
0053
0054 EDPutTokenT<StubAssociation> putTokenReconstructable_;
0055
0056 EDPutTokenT<StubAssociation> putTokenSelection_;
0057
0058 ESGetToken<Setup, SetupRcd> esGetTokenSetup_;
0059 };
0060
0061 StubAssociator::StubAssociator(const ParameterSet& iConfig) {
0062
0063 getTokenTTStubDetSetVec_ = consumes<TTStubDetSetVec>(iConfig.getParameter<InputTag>("InputTagTTStubDetSetVec"));
0064 getTokenTTClusterAssMap_ = consumes<TTClusterAssMap>(iConfig.getParameter<InputTag>("InputTagTTClusterAssMap"));
0065 putTokenReconstructable_ = produces<StubAssociation>(iConfig.getParameter<string>("BranchReconstructable"));
0066 putTokenSelection_ = produces<StubAssociation>(iConfig.getParameter<string>("BranchSelection"));
0067
0068 esGetTokenSetup_ = esConsumes<Setup, SetupRcd, Transition::BeginRun>();
0069 }
0070
0071 void StubAssociator::beginRun(const Run& iRun, const EventSetup& iSetup) {
0072 setup_ = &iSetup.getData(esGetTokenSetup_);
0073 }
0074
0075 void StubAssociator::produce(Event& iEvent, const EventSetup& iSetup) {
0076
0077 Handle<TTStubDetSetVec> handleTTStubDetSetVec;
0078 iEvent.getByToken<TTStubDetSetVec>(getTokenTTStubDetSetVec_, handleTTStubDetSetVec);
0079 Handle<TTClusterAssMap> handleTTClusterAssMap;
0080 Handle<TTStubAssMap> handleTTStubAssMap;
0081 iEvent.getByToken<TTClusterAssMap>(getTokenTTClusterAssMap_, handleTTClusterAssMap);
0082 map<TPPtr, vector<TTStubRef>> mapTPPtrsTTStubRefs;
0083 auto isNonnull = [](const TPPtr& tpPtr) { return tpPtr.isNonnull(); };
0084 for (TTStubDetSetVec::const_iterator ttModule = handleTTStubDetSetVec->begin();
0085 ttModule != handleTTStubDetSetVec->end();
0086 ttModule++) {
0087 for (TTStubDetSet::const_iterator ttStub = ttModule->begin(); ttStub != ttModule->end(); ttStub++) {
0088 const TTStubRef ttStubRef = makeRefTo(handleTTStubDetSetVec, ttStub);
0089 set<TPPtr> tpPtrs;
0090 for (unsigned int iClus = 0; iClus < 2; iClus++) {
0091 const vector<TPPtr>& assocPtrs =
0092 handleTTClusterAssMap->findTrackingParticlePtrs(ttStubRef->clusterRef(iClus));
0093 copy_if(assocPtrs.begin(), assocPtrs.end(), inserter(tpPtrs, tpPtrs.begin()), isNonnull);
0094 }
0095 for (const TPPtr& tpPtr : tpPtrs)
0096 mapTPPtrsTTStubRefs[tpPtr].push_back(ttStubRef);
0097 }
0098 }
0099
0100 StubAssociation reconstructable(setup_);
0101 StubAssociation selection(setup_);
0102 for (const auto& p : mapTPPtrsTTStubRefs) {
0103 if (!setup_->useForReconstructable(*p.first) || !setup_->reconstructable(p.second))
0104 continue;
0105 reconstructable.insert(p.first, p.second);
0106 if (setup_->useForAlgEff(*p.first))
0107 selection.insert(p.first, p.second);
0108 }
0109 iEvent.emplace(putTokenReconstructable_, move(reconstructable));
0110 iEvent.emplace(putTokenSelection_, move(selection));
0111 }
0112
0113 }
0114
0115 DEFINE_FWK_MODULE(tt::StubAssociator);