SubCollectionProducers/src/SeedClusterRemover.cc

0001 #include "FWCore/Framework/interface/Frameworkfwd.h"
0002 #include "FWCore/Framework/interface/stream/EDProducer.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0005 #include "FWCore/Utilities/interface/InputTag.h"
0006
0007 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
0008 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0009 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
0010 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
0011 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
0012 #include "DataFormats/TrackerRecHit2D/interface/ProjectedSiStripRecHit2D.h"
0013 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h"
0014 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHitCollection.h"
0015 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2DCollection.h"
0016 #include "DataFormats/Common/interface/Handle.h"
0017 #include "DataFormats/Common/interface/DetSetVector.h"
0018 #include "DataFormats/Common/interface/ValueMap.h"
0019 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0020 #include "DataFormats/Provenance/interface/ProductID.h"
0021 #include "DataFormats/Common/interface/ContainerMask.h"
0022
0023 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0024 #include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h"
0025
0026 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0027 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0028 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0029 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0030
0031 //
0032 // class decleration
0033 //
0034
0035 class SeedClusterRemover : public edm::stream::EDProducer<> {
0036 public:
0037   SeedClusterRemover(const edm::ParameterSet &iConfig);
0038   void produce(edm::Event &iEvent, const edm::EventSetup &iSetup) override;
0039
0040 private:
0041   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> const tTrackerGeom_;
0042   struct ParamBlock {
0043     ParamBlock() : isSet_(false), usesCharge_(false) {}
0044     ParamBlock(const edm::ParameterSet &iConfig)
0045         : isSet_(true),
0046           usesCharge_(iConfig.exists("maxCharge")),
0047           usesSize_(iConfig.exists("maxSize")),
0048           maxChi2_(iConfig.getParameter<double>("maxChi2")),
0049           maxCharge_(usesCharge_ ? iConfig.getParameter<double>("maxCharge") : 0),
0050           maxSize_(usesSize_ ? iConfig.getParameter<uint32_t>("maxSize") : 0) {}
0051     bool isSet_, usesCharge_, usesSize_;
0052     float maxChi2_, maxCharge_;
0053     size_t maxSize_;
0054   };
0055   static const unsigned int NumberOfParamBlocks = 6;
0056
0057   bool doStrip_, doPixel_;
0058   bool mergeOld_;
0059
0060   typedef edm::ContainerMask<edmNew::DetSetVector<SiPixelCluster> > PixelMaskContainer;
0061   typedef edm::ContainerMask<edmNew::DetSetVector<SiStripCluster> > StripMaskContainer;
0062   edm::EDGetTokenT<edmNew::DetSetVector<SiPixelCluster> > pixelClusters_;
0063   edm::EDGetTokenT<edmNew::DetSetVector<SiStripCluster> > stripClusters_;
0064   edm::EDGetTokenT<PixelMaskContainer> oldPxlMaskToken_;
0065   edm::EDGetTokenT<StripMaskContainer> oldStrMaskToken_;
0066   edm::EDGetTokenT<TrajectorySeedCollection> trajectories_;
0067
0068   ParamBlock pblocks_[NumberOfParamBlocks];
0069   void readPSet(const edm::ParameterSet &iConfig,
0070                 const std::string &name,
0071                 int id1 = -1,
0072                 int id2 = -1,
0073                 int id3 = -1,
0074                 int id4 = -1,
0075                 int id5 = -1,
0076                 int id6 = -1);
0077
0078   std::vector<uint8_t> pixels, strips;                  // avoid unneed alloc/dealloc of this
0079   edm::ProductID pixelSourceProdID, stripSourceProdID;  // ProdIDs refs must point to (for consistency tests)
0080
0081   inline void process(const TrackingRecHit *hit, float chi2, const TrackerGeometry *tg);
0082   inline void process(const OmniClusterRef &cluRef, uint32_t subdet);
0083
0084   std::vector<bool> collectedStrips_;
0085   std::vector<bool> collectedPixels_;
0086 };
0087
0088 using namespace std;
0089 using namespace edm;
0090
0091 void SeedClusterRemover::readPSet(
0092     const edm::ParameterSet &iConfig, const std::string &name, int id1, int id2, int id3, int id4, int id5, int id6) {
0093   if (iConfig.exists(name)) {
0094     ParamBlock pblock(iConfig.getParameter<ParameterSet>(name));
0095     if (id1 == -1) {
0096       fill(pblocks_, pblocks_ + NumberOfParamBlocks, pblock);
0097     } else {
0098       pblocks_[id1] = pblock;
0099       if (id2 != -1)
0100         pblocks_[id2] = pblock;
0101       if (id3 != -1)
0102         pblocks_[id3] = pblock;
0103       if (id4 != -1)
0104         pblocks_[id4] = pblock;
0105       if (id5 != -1)
0106         pblocks_[id5] = pblock;
0107       if (id6 != -1)
0108         pblocks_[id6] = pblock;
0109     }
0110   }
0111 }
0112
0113 SeedClusterRemover::SeedClusterRemover(const ParameterSet &iConfig)
0114     : tTrackerGeom_(esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>()),
0115       doStrip_(iConfig.existsAs<bool>("doStrip") ? iConfig.getParameter<bool>("doStrip") : true),
0116       doPixel_(iConfig.existsAs<bool>("doPixel") ? iConfig.getParameter<bool>("doPixel") : true),
0117       mergeOld_(iConfig.exists("oldClusterRemovalInfo")) {
0118   fill(pblocks_, pblocks_ + NumberOfParamBlocks, ParamBlock());
0119   readPSet(iConfig, "Common", -1);
0120   if (doPixel_) {
0121     readPSet(iConfig, "Pixel", 0, 1);
0122     readPSet(iConfig, "PXB", 0);
0123     readPSet(iConfig, "PXE", 1);
0124   }
0125   if (doStrip_) {
0126     readPSet(iConfig, "Strip", 2, 3, 4, 5);
0127     readPSet(iConfig, "StripInner", 2, 3);
0128     readPSet(iConfig, "StripOuter", 4, 5);
0129     readPSet(iConfig, "TIB", 2);
0130     readPSet(iConfig, "TID", 3);
0131     readPSet(iConfig, "TOB", 4);
0132     readPSet(iConfig, "TEC", 5);
0133   }
0134
0135   bool usingCharge = false;
0136   for (size_t i = 0; i < NumberOfParamBlocks; ++i) {
0137     if (!pblocks_[i].isSet_)
0138       throw cms::Exception("Configuration Error")
0139           << "SeedClusterRemover: Missing configuration for detector with subDetID = " << (i + 1);
0140     if (pblocks_[i].usesCharge_ && !usingCharge) {
0141       throw cms::Exception("Configuration Error") << "SeedClusterRemover: Configuration for subDetID = " << (i + 1)
0142                                                   << " uses cluster charge, which is not enabled.";
0143     }
0144   }
0145
0146   trajectories_ = consumes<TrajectorySeedCollection>(iConfig.getParameter<InputTag>("trajectories"));
0147   if (doPixel_)
0148     pixelClusters_ = consumes<edmNew::DetSetVector<SiPixelCluster> >(iConfig.getParameter<InputTag>("pixelClusters"));
0149   if (doStrip_)
0150     stripClusters_ = consumes<edmNew::DetSetVector<SiStripCluster> >(iConfig.getParameter<InputTag>("stripClusters"));
0151   if (mergeOld_) {
0152     oldPxlMaskToken_ = consumes<PixelMaskContainer>(iConfig.getParameter<InputTag>("oldClusterRemovalInfo"));
0153     oldStrMaskToken_ = consumes<StripMaskContainer>(iConfig.getParameter<InputTag>("oldClusterRemovalInfo"));
0154   }
0155
0156   produces<edm::ContainerMask<edmNew::DetSetVector<SiPixelCluster> > >();
0157   produces<edm::ContainerMask<edmNew::DetSetVector<SiStripCluster> > >();
0158 }
0159
0160 void SeedClusterRemover::process(OmniClusterRef const &ocluster, uint32_t subdet) {
0161   SiStripRecHit2D::ClusterRef cluster = ocluster.cluster_strip();
0162   if (cluster.id() != stripSourceProdID)
0163     throw cms::Exception("Inconsistent Data")
0164         << "SeedClusterRemover: strip cluster ref from Product ID = " << cluster.id()
0165         << " does not match with source cluster collection (ID = " << stripSourceProdID << ")\n.";
0166
0167   assert(cluster.id() == stripSourceProdID);
0168   if (pblocks_[subdet - 1].usesSize_ && (cluster->amplitudes().size() > pblocks_[subdet - 1].maxSize_))
0169     return;
0170
0171   strips[cluster.key()] = false;
0172   assert(collectedStrips_.size() > cluster.key());
0173   collectedStrips_[cluster.key()] = true;
0174 }
0175
0176 void SeedClusterRemover::process(const TrackingRecHit *hit, float chi2, const TrackerGeometry *tg) {
0177   DetId detid = hit->geographicalId();
0178   uint32_t subdet = detid.subdetId();
0179
0180   assert((subdet > 0) && (subdet <= NumberOfParamBlocks));
0181
0182   // chi2 cut
0183   if (chi2 > pblocks_[subdet - 1].maxChi2_)
0184     return;
0185
0186   if (GeomDetEnumerators::isTrackerPixel(tg->geomDetSubDetector(subdet))) {
0187     if (!doPixel_)
0188       return;
0189     // this is a pixel, and i *know* it is
0190     const SiPixelRecHit *pixelHit = static_cast<const SiPixelRecHit *>(hit);
0191
0192     SiPixelRecHit::ClusterRef cluster = pixelHit->cluster();
0193
0194     if (cluster.id() != pixelSourceProdID)
0195       throw cms::Exception("Inconsistent Data")
0196           << "SeedClusterRemover: pixel cluster ref from Product ID = " << cluster.id()
0197           << " does not match with source cluster collection (ID = " << pixelSourceProdID << ")\n.";
0198
0199     assert(cluster.id() == pixelSourceProdID);
0200     //DBG// cout << "HIT NEW PIXEL DETID = " << detid.rawId() << ", Cluster [ " << cluster.key().first << " / " <<  cluster.key().second << " ] " << endl;
0201
0202     // if requested, cut on cluster size
0203     if (pblocks_[subdet - 1].usesSize_ && (cluster->pixels().size() > pblocks_[subdet - 1].maxSize_))
0204       return;
0205
0206     // mark as used
0207     pixels[cluster.key()] = false;
0208
0209     assert(collectedPixels_.size() > cluster.key());
0210     collectedPixels_[cluster.key()] = true;
0211
0212   } else {  // aka Strip
0213     if (!doStrip_)
0214       return;
0215     const type_info &hitType = typeid(*hit);
0216     if (hitType == typeid(SiStripRecHit2D)) {
0217       const SiStripRecHit2D *stripHit = static_cast<const SiStripRecHit2D *>(hit);
0218       //DBG//     cout << "Plain RecHit 2D: " << endl;
0219       process(stripHit->omniClusterRef(), subdet);
0220     } else if (hitType == typeid(SiStripRecHit1D)) {
0221       const SiStripRecHit1D *hit1D = static_cast<const SiStripRecHit1D *>(hit);
0222       process(hit1D->omniClusterRef(), subdet);
0223     } else if (hitType == typeid(SiStripMatchedRecHit2D)) {
0224       const SiStripMatchedRecHit2D *matchHit = static_cast<const SiStripMatchedRecHit2D *>(hit);
0225       //DBG//     cout << "Matched RecHit 2D: " << endl;
0226       process(matchHit->monoClusterRef(), subdet);
0227       process(matchHit->stereoClusterRef(), subdet);
0228     } else if (hitType == typeid(ProjectedSiStripRecHit2D)) {
0229       const ProjectedSiStripRecHit2D *projHit = static_cast<const ProjectedSiStripRecHit2D *>(hit);
0230       //DBG//     cout << "Projected RecHit 2D: " << endl;
0231       process(projHit->originalHit().omniClusterRef(), subdet);
0232     } else
0233       throw cms::Exception("NOT IMPLEMENTED")
0234           << "Don't know how to handle " << hitType.name() << " on detid " << detid.rawId() << "\n";
0235   }
0236 }
0237
0238 void SeedClusterRemover::produce(Event &iEvent, const EventSetup &iSetup) {
0239   const auto &tgh = &iSetup.getData(tTrackerGeom_);
0240
0241   Handle<edmNew::DetSetVector<SiPixelCluster> > pixelClusters;
0242   if (doPixel_) {
0243     iEvent.getByToken(pixelClusters_, pixelClusters);
0244     pixelSourceProdID = pixelClusters.id();
0245   }
0246   //DBG// std::cout << "SeedClusterRemover: Read pixel " << pixelClusters_.encode() << " = ID " << pixelSourceProdID << std::endl;
0247
0248   Handle<edmNew::DetSetVector<SiStripCluster> > stripClusters;
0249   if (doStrip_) {
0250     iEvent.getByToken(stripClusters_, stripClusters);
0251     stripSourceProdID = stripClusters.id();
0252   }
0253   //DBG// std::cout << "SeedClusterRemover: Read strip " << stripClusters_.encode() << " = ID " << stripSourceProdID << std::endl;
0254
0255   if (doStrip_) {
0256     strips.resize(stripClusters->dataSize());
0257     fill(strips.begin(), strips.end(), true);
0258   }
0259   if (doPixel_) {
0260     pixels.resize(pixelClusters->dataSize());
0261     fill(pixels.begin(), pixels.end(), true);
0262   }
0263   if (mergeOld_) {
0264     edm::Handle<PixelMaskContainer> oldPxlMask;
0265     edm::Handle<StripMaskContainer> oldStrMask;
0266     iEvent.getByToken(oldPxlMaskToken_, oldPxlMask);
0267     iEvent.getByToken(oldStrMaskToken_, oldStrMask);
0268     LogDebug("SeedClusterRemover") << "to merge in, " << oldStrMask->size() << " strp and " << oldPxlMask->size()
0269                                    << " pxl";
0270     oldStrMask->copyMaskTo(collectedStrips_);
0271     oldPxlMask->copyMaskTo(collectedPixels_);
0272     assert(stripClusters->dataSize() >= collectedStrips_.size());
0273     collectedStrips_.resize(stripClusters->dataSize(), false);  // for ondemand
0274   } else {
0275     collectedStrips_.resize(stripClusters->dataSize(), false);
0276     collectedPixels_.resize(pixelClusters->dataSize(), false);
0277   }
0278
0279   edm::Handle<TrajectorySeedCollection> seeds;
0280   iEvent.getByToken(trajectories_, seeds);
0281
0282   for (auto const &seed : (*seeds)) {
0283     for (auto const &hit : seed.recHits()) {
0284       if (!hit.isValid())
0285         continue;
0286       process(&hit, 0., tgh);
0287     }
0288   }
0289
0290   pixels.clear();
0291   strips.clear();
0292
0293   LogDebug("SeedClusterRemover") << "total strip to skip: "
0294                                  << std::count(collectedStrips_.begin(), collectedStrips_.end(), true);
0295   iEvent.put(std::make_unique<StripMaskContainer>(edm::RefProd<edmNew::DetSetVector<SiStripCluster> >(stripClusters),
0296                                                   collectedStrips_));
0297
0298   LogDebug("SeedClusterRemover") << "total pxl to skip: "
0299                                  << std::count(collectedPixels_.begin(), collectedPixels_.end(), true);
0300   iEvent.put(std::make_unique<PixelMaskContainer>(edm::RefProd<edmNew::DetSetVector<SiPixelCluster> >(pixelClusters),
0301                                                   collectedPixels_));
0302
0303   collectedStrips_.clear();
0304   collectedPixels_.clear();
0305 }
0306
0307 #include "FWCore/PluginManager/interface/ModuleDef.h"
0308 #include "FWCore/Framework/interface/MakerMacros.h"
0309 DEFINE_FWK_MODULE(SeedClusterRemover);