forked from cms-sw/cmssw
/
MTDRecHitProducer.cc
96 lines (71 loc) · 3.61 KB
/
MTDRecHitProducer.cc
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
#include "FWCore/Framework/interface/stream/EDProducer.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/EventSetup.h"
#include "DataFormats/Common/interface/Handle.h"
#include "DataFormats/FTLDigi/interface/FTLDigiCollections.h"
#include "DataFormats/FTLRecHit/interface/FTLRecHitCollections.h"
#include "RecoLocalFastTime/FTLCommonAlgos/interface/MTDRecHitAlgoBase.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
class MTDRecHitProducer : public edm::stream::EDProducer<> {
public:
explicit MTDRecHitProducer(const edm::ParameterSet& ps);
~MTDRecHitProducer() override;
void produce(edm::Event& evt, const edm::EventSetup& es) override;
private:
const edm::EDGetTokenT<FTLUncalibratedRecHitCollection> ftlbURecHits_; // collection of barrel digis
const edm::EDGetTokenT<FTLUncalibratedRecHitCollection> ftleURecHits_; // collection of endcap digis
const std::string ftlbInstance_; // instance name of barrel hits
const std::string ftleInstance_; // instance name of endcap hits
std::unique_ptr<MTDRecHitAlgoBase> barrel_,endcap_;
};
MTDRecHitProducer::MTDRecHitProducer(const edm::ParameterSet& ps) :
ftlbURecHits_( consumes<FTLUncalibratedRecHitCollection>( ps.getParameter<edm::InputTag>("barrelUncalibratedRecHits") ) ),
ftleURecHits_( consumes<FTLUncalibratedRecHitCollection>( ps.getParameter<edm::InputTag>("endcapUncalibratedRecHits") ) ),
ftlbInstance_( ps.getParameter<std::string>("BarrelHitsName") ),
ftleInstance_( ps.getParameter<std::string>("EndcapHitsName") )
{
produces< FTLRecHitCollection >(ftlbInstance_);
produces< FTLRecHitCollection >(ftleInstance_);
auto sumes = consumesCollector();
const edm::ParameterSet& barrel = ps.getParameterSet("barrel");
const std::string& barrelAlgo = barrel.getParameter<std::string>("algoName");
barrel_.reset( MTDRecHitAlgoFactory::get()->create(barrelAlgo, barrel, sumes) );
const edm::ParameterSet& endcap = ps.getParameterSet("endcap");
const std::string& endcapAlgo = endcap.getParameter<std::string>("algoName");
endcap_.reset( MTDRecHitAlgoFactory::get()->create(endcapAlgo, endcap, sumes) );
}
MTDRecHitProducer::~MTDRecHitProducer() {
}
void
MTDRecHitProducer::produce(edm::Event& evt, const edm::EventSetup& es) {
// tranparently get things from event setup
barrel_->getEventSetup(es);
endcap_->getEventSetup(es);
barrel_->getEvent(evt);
endcap_->getEvent(evt);
// prepare output
auto barrelRechits = std::make_unique<FTLRecHitCollection>();
auto endcapRechits = std::make_unique<FTLRecHitCollection>();
edm::Handle< FTLUncalibratedRecHitCollection > hBarrel;
evt.getByToken( ftlbURecHits_, hBarrel );
barrelRechits->reserve(hBarrel->size()/2);
for(const auto& uhit : *hBarrel) {
uint32_t flags = FTLRecHit::kGood;
auto rechit = std::move( barrel_->makeRecHit(uhit, flags) );
if( flags == FTLRecHit::kGood ) barrelRechits->push_back( std::move(rechit) );
}
edm::Handle< FTLUncalibratedRecHitCollection > hEndcap;
evt.getByToken( ftleURecHits_, hEndcap );
endcapRechits->reserve(hEndcap->size()/2);
for(const auto& uhit : *hEndcap) {
uint32_t flags = FTLRecHit::kGood;
auto rechit = std::move( endcap_->makeRecHit(uhit, flags) );
if( flags == FTLRecHit::kGood ) endcapRechits->push_back( std::move(rechit) );
}
// put the collection of recunstructed hits in the event
evt.put(std::move(barrelRechits), ftlbInstance_);
evt.put(std::move(endcapRechits), ftleInstance_);
}
#include "FWCore/Framework/interface/MakerMacros.h"
DEFINE_FWK_MODULE( MTDRecHitProducer );