CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
PFClusterTimeSelector.cc
Go to the documentation of this file.
2 #include <memory>
3 
7 
8 using namespace std;
9 using namespace edm;
10 
12  : clusters_(consumes<reco::PFClusterCollection>(iConfig.getParameter<edm::InputTag>("src"))) {
13  std::vector<edm::ParameterSet> cuts = iConfig.getParameter<std::vector<edm::ParameterSet> >("cuts");
14  for (const auto& cut : cuts) {
15  CutInfo info;
16  info.depth = cut.getParameter<double>("depth");
17  info.minE = cut.getParameter<double>("minEnergy");
18  info.maxE = cut.getParameter<double>("maxEnergy");
19  info.minTime = cut.getParameter<double>("minTime");
20  info.maxTime = cut.getParameter<double>("maxTime");
21  info.endcap = cut.getParameter<bool>("endcap");
22  cutInfo_.push_back(info);
23  }
24 
25  produces<reco::PFClusterCollection>();
26  produces<reco::PFClusterCollection>("OOT");
27 }
28 
31  iEvent.getByToken(clusters_, clusters);
32  auto out = std::make_unique<reco::PFClusterCollection>();
33  auto outOOT = std::make_unique<reco::PFClusterCollection>();
34 
35  for (const auto& cluster : *clusters) {
36  const double energy = cluster.energy();
37  const double time = cluster.time();
38  const double depth = cluster.depth();
39  const PFLayer::Layer layer = cluster.layer();
40  for (const auto& info : cutInfo_) {
41  if (energy < info.minE || energy > info.maxE)
42  continue;
43  if (depth < 0.9 * info.depth || depth > 1.1 * info.depth)
44  continue;
45  if ((info.endcap &&
46  (layer == PFLayer::ECAL_BARREL || layer == PFLayer::HCAL_BARREL1 || layer == PFLayer::HCAL_BARREL2)) ||
47  (((!info.endcap) && (layer == PFLayer::ECAL_ENDCAP || layer == PFLayer::HCAL_ENDCAP))))
48  continue;
49 
50  if (time > info.minTime && time < info.maxTime)
51  out->push_back(cluster);
52  else
53  outOOT->push_back(cluster);
54 
55  break;
56  }
57  }
58 
59  iEvent.put(std::move(out));
60  iEvent.put(std::move(outOOT), "OOT");
61 }
62 
64 
65 // ------------ method called once each job just before starting event loop ------------
T getParameter(std::string const &) const
static const TGPicture * info(bool iBackgroundIsBlack)
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:131
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:525
std::vector< CutInfo > cutInfo_
PFClusterTimeSelector(const edm::ParameterSet &)
~PFClusterTimeSelector() override
void beginRun(const edm::Run &run, const edm::EventSetup &es) override
int iEvent
Definition: GenABIO.cc:224
void produce(edm::Event &iEvent, const edm::EventSetup &iSetup) override
Layer
layer definition
Definition: PFLayer.h:29
fixed size matrix
HLT enums.
std::vector< PFCluster > PFClusterCollection
collection of PFCluster objects
Definition: PFClusterFwd.h:9
edm::EDGetTokenT< reco::PFClusterCollection > clusters_
def move(src, dest)
Definition: eostools.py:511
Definition: Run.h:45