File indexing completed on 2023-10-25 10:00:48
0001 #include "FWCore/Framework/interface/MakerMacros.h"
0002 #include "RecoLocalTracker/SubCollectionProducers/interface/PixelClusterSelectorTopBottom.h"
0003
0004 void PixelClusterSelectorTopBottom::produce(edm::StreamID, edm::Event& event, const edm::EventSetup& setup) const {
0005 edm::Handle<SiPixelClusterCollectionNew> input;
0006 event.getByToken(token_, input);
0007
0008 const TrackerGeometry& theTracker = setup.getData(tTrackerGeom_);
0009
0010 auto output = std::make_unique<SiPixelClusterCollectionNew>();
0011
0012 for (SiPixelClusterCollectionNew::const_iterator clustSet = input->begin(); clustSet != input->end(); ++clustSet) {
0013 edmNew::DetSet<SiPixelCluster>::const_iterator clustIt = clustSet->begin();
0014 edmNew::DetSet<SiPixelCluster>::const_iterator end = clustSet->end();
0015
0016 DetId detIdObject(clustSet->detId());
0017 edmNew::DetSetVector<SiPixelCluster>::FastFiller spc(*output, detIdObject);
0018 const PixelGeomDetUnit* theGeomDet = dynamic_cast<const PixelGeomDetUnit*>(theTracker.idToDet(detIdObject));
0019 const PixelTopology* topol = (&(theGeomDet->specificTopology()));
0020
0021 for (; clustIt != end; ++clustIt) {
0022 LocalPoint lpclust = topol->localPosition(MeasurementPoint((*clustIt).x(), (*clustIt).y()));
0023 GlobalPoint GPclust = theGeomDet->surface().toGlobal(Local3DPoint(lpclust.x(), lpclust.y(), lpclust.z()));
0024 double clustY = GPclust.y();
0025 if ((clustY * y_) > 0) {
0026 spc.push_back(*clustIt);
0027 }
0028 }
0029 }
0030 event.put(std::move(output));
0031 }
0032
0033 DEFINE_FWK_MODULE(PixelClusterSelectorTopBottom);