Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-03-14 23:24:27

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/DetSetVectorNew.h"
0019 #include "DataFormats/Common/interface/ContainerMask.h"
0020 #include "DataFormats/Provenance/interface/ProductID.h"
0021 
0022 #include "DataFormats/TrackReco/interface/Track.h"
0023 #include "DataFormats/TrackerRecHit2D/interface/ClusterRemovalInfo.h"
0024 
0025 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0026 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
0027 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0028 
0029 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0030 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0031 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0032 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0033 
0034 //
0035 // class decleration
0036 //
0037 
0038 class HLTTrackClusterRemoverNew final : public edm::stream::EDProducer<> {
0039 public:
0040   HLTTrackClusterRemoverNew(const edm::ParameterSet &iConfig);
0041   ~HLTTrackClusterRemoverNew() override;
0042   void produce(edm::Event &iEvent, const edm::EventSetup &iSetup) override;
0043 
0044 private:
0045   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> const tTrackerGeom_;
0046   struct ParamBlock {
0047     ParamBlock() : isSet_(false), usesCharge_(false) {}
0048     ParamBlock(const edm::ParameterSet &iConfig)
0049         : isSet_(true),
0050           usesCharge_(iConfig.exists("maxCharge")),
0051           usesSize_(iConfig.exists("maxSize")),
0052           cutOnPixelCharge_(iConfig.exists("minGoodPixelCharge")),
0053           cutOnStripCharge_(iConfig.exists("minGoodStripCharge")),
0054           maxChi2_(iConfig.getParameter<double>("maxChi2")),
0055           maxCharge_(usesCharge_ ? iConfig.getParameter<double>("maxCharge") : 0),
0056           minGoodPixelCharge_(cutOnPixelCharge_ ? iConfig.getParameter<double>("minGoodPixelCharge") : 0),
0057           minGoodStripCharge_(cutOnStripCharge_ ? iConfig.getParameter<double>("minGoodStripCharge") : 0),
0058           maxSize_(usesSize_ ? iConfig.getParameter<uint32_t>("maxSize") : 0) {}
0059     bool isSet_, usesCharge_, usesSize_, cutOnPixelCharge_, cutOnStripCharge_;
0060     float maxChi2_, maxCharge_, minGoodPixelCharge_, minGoodStripCharge_;
0061     size_t maxSize_;
0062   };
0063   static const unsigned int NumberOfParamBlocks = 6;
0064 
0065   bool doTracks_;
0066   bool doStrip_, doPixel_;
0067   bool mergeOld_;
0068 
0069   typedef edm::ContainerMask<edmNew::DetSetVector<SiPixelCluster> > PixelMaskContainer;
0070   typedef edm::ContainerMask<edmNew::DetSetVector<SiStripCluster> > StripMaskContainer;
0071   edm::EDGetTokenT<edmNew::DetSetVector<SiPixelCluster> > pixelClusters_;
0072   edm::EDGetTokenT<edmNew::DetSetVector<SiStripCluster> > stripClusters_;
0073   edm::EDGetTokenT<PixelMaskContainer> oldPxlMaskToken_;
0074   edm::EDGetTokenT<StripMaskContainer> oldStrMaskToken_;
0075   edm::EDGetTokenT<std::vector<Trajectory> > trajectories_;
0076 
0077   ParamBlock pblocks_[NumberOfParamBlocks];
0078   void readPSet(const edm::ParameterSet &iConfig,
0079                 const std::string &name,
0080                 int id1 = -1,
0081                 int id2 = -1,
0082                 int id3 = -1,
0083                 int id4 = -1,
0084                 int id5 = -1,
0085                 int id6 = -1);
0086 
0087   std::vector<uint8_t> pixels, strips;                  // avoid unneed alloc/dealloc of this
0088   edm::ProductID pixelSourceProdID, stripSourceProdID;  // ProdIDs refs must point to (for consistency tests)
0089 
0090   inline void process(const TrackingRecHit *hit, float chi2, const TrackerGeometry *tg);
0091   inline void process(const OmniClusterRef &cluRef, uint32_t subdet);
0092 
0093   template <typename T>
0094   std::unique_ptr<edmNew::DetSetVector<T> > cleanup(const edmNew::DetSetVector<T> &oldClusters,
0095                                                     const std::vector<uint8_t> &isGood,
0096                                                     reco::ClusterRemovalInfo::Indices &refs,
0097                                                     const reco::ClusterRemovalInfo::Indices *oldRefs);
0098 
0099   // Carries in full removal info about a given det from oldRefs
0100   void mergeOld(reco::ClusterRemovalInfo::Indices &refs, const reco::ClusterRemovalInfo::Indices &oldRefs);
0101 
0102   bool makeProducts_;
0103   bool doStripChargeCheck_, doPixelChargeCheck_;
0104   std::vector<bool> collectedRegStrips_;
0105   std::vector<bool> collectedPixels_;
0106 };
0107 
0108 using namespace std;
0109 using namespace edm;
0110 using namespace reco;
0111 
0112 void HLTTrackClusterRemoverNew::readPSet(
0113     const edm::ParameterSet &iConfig, const std::string &name, int id1, int id2, int id3, int id4, int id5, int id6) {
0114   if (iConfig.exists(name)) {
0115     ParamBlock pblock(iConfig.getParameter<ParameterSet>(name));
0116     if (id1 == -1) {
0117       fill(pblocks_, pblocks_ + NumberOfParamBlocks, pblock);
0118     } else {
0119       pblocks_[id1] = pblock;
0120       if (id2 != -1)
0121         pblocks_[id2] = pblock;
0122       if (id3 != -1)
0123         pblocks_[id3] = pblock;
0124       if (id4 != -1)
0125         pblocks_[id4] = pblock;
0126       if (id5 != -1)
0127         pblocks_[id5] = pblock;
0128       if (id6 != -1)
0129         pblocks_[id6] = pblock;
0130     }
0131   }
0132 }
0133 
0134 HLTTrackClusterRemoverNew::HLTTrackClusterRemoverNew(const ParameterSet &iConfig)
0135     : tTrackerGeom_(esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>()),
0136       doTracks_(iConfig.exists("trajectories")),
0137       doStrip_(iConfig.existsAs<bool>("doStrip") ? iConfig.getParameter<bool>("doStrip") : true),
0138       doPixel_(iConfig.existsAs<bool>("doPixel") ? iConfig.getParameter<bool>("doPixel") : true),
0139       mergeOld_(false),
0140       makeProducts_(true),
0141       doStripChargeCheck_(
0142           iConfig.existsAs<bool>("doStripChargeCheck") ? iConfig.getParameter<bool>("doStripChargeCheck") : false),
0143       doPixelChargeCheck_(
0144           iConfig.existsAs<bool>("doPixelChargeCheck") ? iConfig.getParameter<bool>("doPixelChargeCheck") : false)
0145 
0146 {
0147   if (iConfig.exists("oldClusterRemovalInfo")) {
0148     oldPxlMaskToken_ = consumes<PixelMaskContainer>(iConfig.getParameter<InputTag>("oldClusterRemovalInfo"));
0149     oldStrMaskToken_ = consumes<StripMaskContainer>(iConfig.getParameter<InputTag>("oldClusterRemovalInfo"));
0150     if (not(iConfig.getParameter<InputTag>("oldClusterRemovalInfo") == edm::InputTag()))
0151       mergeOld_ = true;
0152   }
0153 
0154   if ((doPixelChargeCheck_ && !doPixel_) || (doStripChargeCheck_ && !doStrip_))
0155     throw cms::Exception("Configuration Error")
0156         << "HLTTrackClusterRemoverNew: Charge check asked without cluster collection ";
0157   if (doPixelChargeCheck_)
0158     throw cms::Exception("Configuration Error")
0159         << "HLTTrackClusterRemoverNew: Pixel cluster charge check not yet implemented";
0160 
0161   fill(pblocks_, pblocks_ + NumberOfParamBlocks, ParamBlock());
0162   readPSet(iConfig, "Common", -1);
0163   if (doPixel_) {
0164     readPSet(iConfig, "Pixel", 0, 1);
0165     readPSet(iConfig, "PXB", 0);
0166     readPSet(iConfig, "PXE", 1);
0167   }
0168   if (doStrip_) {
0169     readPSet(iConfig, "Strip", 2, 3, 4, 5);
0170     readPSet(iConfig, "StripInner", 2, 3);
0171     readPSet(iConfig, "StripOuter", 4, 5);
0172     readPSet(iConfig, "TIB", 2);
0173     readPSet(iConfig, "TID", 3);
0174     readPSet(iConfig, "TOB", 4);
0175     readPSet(iConfig, "TEC", 5);
0176   }
0177 
0178   bool usingCharge = false;
0179   for (size_t i = 0; i < NumberOfParamBlocks; ++i) {
0180     if (!pblocks_[i].isSet_)
0181       throw cms::Exception("Configuration Error")
0182           << "HLTTrackClusterRemoverNew: Missing configuration for detector with subDetID = " << (i + 1);
0183     if (pblocks_[i].usesCharge_ && !usingCharge) {
0184       throw cms::Exception("Configuration Error")
0185           << "HLTTrackClusterRemoverNew: Configuration for subDetID = " << (i + 1)
0186           << " uses cluster charge, which is not enabled.";
0187     }
0188   }
0189 
0190   //    trajectories_ = consumes<vector<Trajectory> >(iConfig.getParameter<InputTag>("trajectories"));
0191   if (doTracks_)
0192     trajectories_ = consumes<vector<Trajectory> >(iConfig.getParameter<InputTag>("trajectories"));
0193   if (doPixel_)
0194     pixelClusters_ = consumes<edmNew::DetSetVector<SiPixelCluster> >(iConfig.getParameter<InputTag>("pixelClusters"));
0195   if (doStrip_)
0196     stripClusters_ = consumes<edmNew::DetSetVector<SiStripCluster> >(iConfig.getParameter<InputTag>("stripClusters"));
0197   if (mergeOld_) {
0198     oldPxlMaskToken_ = consumes<PixelMaskContainer>(iConfig.getParameter<InputTag>("oldClusterRemovalInfo"));
0199     oldStrMaskToken_ = consumes<StripMaskContainer>(iConfig.getParameter<InputTag>("oldClusterRemovalInfo"));
0200   }
0201 
0202   //produces<edmNew::DetSetVector<SiPixelClusterRefNew> >();
0203   //produces<edmNew::DetSetVector<SiStripRecHit1D::ClusterRegionalRef> >();
0204 
0205   produces<edm::ContainerMask<edmNew::DetSetVector<SiPixelCluster> > >();
0206   produces<edm::ContainerMask<edmNew::DetSetVector<SiStripCluster> > >();
0207 }
0208 
0209 HLTTrackClusterRemoverNew::~HLTTrackClusterRemoverNew() {}
0210 
0211 void HLTTrackClusterRemoverNew::mergeOld(ClusterRemovalInfo::Indices &refs,
0212                                          const ClusterRemovalInfo::Indices &oldRefs) {
0213   for (size_t i = 0, n = refs.size(); i < n; ++i) {
0214     refs[i] = oldRefs[refs[i]];
0215   }
0216 }
0217 
0218 template <typename T>
0219 std::unique_ptr<edmNew::DetSetVector<T> > HLTTrackClusterRemoverNew::cleanup(
0220     const edmNew::DetSetVector<T> &oldClusters,
0221     const std::vector<uint8_t> &isGood,
0222     reco::ClusterRemovalInfo::Indices &refs,
0223     const reco::ClusterRemovalInfo::Indices *oldRefs) {
0224   typedef typename edmNew::DetSetVector<T> DSV;
0225   typedef typename edmNew::DetSetVector<T>::FastFiller DSF;
0226   typedef typename edmNew::DetSet<T> DS;
0227   auto output = std::make_unique<DSV>();
0228   output->reserve(oldClusters.size(), oldClusters.dataSize());
0229 
0230   unsigned int countOld = 0;
0231   unsigned int countNew = 0;
0232 
0233   // cluster removal loop
0234   const T *firstOffset = &oldClusters.data().front();
0235   for (typename DSV::const_iterator itdet = oldClusters.begin(), enddet = oldClusters.end(); itdet != enddet; ++itdet) {
0236     DS oldDS = *itdet;
0237 
0238     if (oldDS.empty())
0239       continue;  // skip empty detsets
0240 
0241     uint32_t id = oldDS.detId();
0242     DSF outds(*output, id);
0243 
0244     for (typename DS::const_iterator it = oldDS.begin(), ed = oldDS.end(); it != ed; ++it) {
0245       uint32_t index = ((&*it) - firstOffset);
0246       countOld++;
0247       if (isGood[index]) {
0248         outds.push_back(*it);
0249         countNew++;
0250         refs.push_back(index);
0251         //std::cout << "HLTTrackClusterRemoverNew::cleanup " << typeid(T).name() << " reference " << index << " to " << (refs.size() - 1) << std::endl;
0252       }
0253     }
0254     if (outds.empty())
0255       outds.abort();  // not write in an empty DSV
0256   }
0257   //    double fraction = countNew  / (double) countOld;
0258   //    std::cout<<"fraction: "<<fraction<<std::endl;
0259   if (oldRefs != nullptr)
0260     mergeOld(refs, *oldRefs);
0261   return output;
0262 }
0263 
0264 void HLTTrackClusterRemoverNew::process(OmniClusterRef const &clusterReg, uint32_t subdet) {
0265   if (clusterReg.id() != stripSourceProdID)
0266     throw cms::Exception("Inconsistent Data")
0267         << "HLTTrackClusterRemoverNew: strip cluster ref from Product ID = " << clusterReg.id()
0268         << " does not match with source cluster collection (ID = " << stripSourceProdID << ")\n.";
0269 
0270   if (collectedRegStrips_.size() <= clusterReg.key()) {
0271     edm::LogError("BadCollectionSize") << collectedRegStrips_.size() << " is smaller than " << clusterReg.key();
0272 
0273     assert(collectedRegStrips_.size() > clusterReg.key());
0274   }
0275   collectedRegStrips_[clusterReg.key()] = true;
0276 }
0277 
0278 void HLTTrackClusterRemoverNew::process(const TrackingRecHit *hit, float chi2, const TrackerGeometry *tg) {
0279   DetId detid = hit->geographicalId();
0280   uint32_t subdet = detid.subdetId();
0281 
0282   assert((subdet > 0) && (subdet <= NumberOfParamBlocks));
0283 
0284   // chi2 cut
0285   if (chi2 > pblocks_[subdet - 1].maxChi2_)
0286     return;
0287 
0288   if (GeomDetEnumerators::isTrackerPixel(tg->geomDetSubDetector(subdet))) {
0289     //      std::cout<<"process pxl hit"<<std::endl;
0290     if (!doPixel_)
0291       return;
0292     // this is a pixel, and i *know* it is
0293     const SiPixelRecHit *pixelHit = static_cast<const SiPixelRecHit *>(hit);
0294 
0295     SiPixelRecHit::ClusterRef cluster = pixelHit->cluster();
0296     if (cluster.id() != pixelSourceProdID)
0297       throw cms::Exception("Inconsistent Data")
0298           << "HLTTrackClusterRemoverNew: pixel cluster ref from Product ID = " << cluster.id()
0299           << " does not match with source cluster collection (ID = " << pixelSourceProdID << ")\n.";
0300 
0301     assert(cluster.id() == pixelSourceProdID);
0302     //DBG// cout << "HIT NEW PIXEL DETID = " << detid.rawId() << ", Cluster [ " << cluster.key().first << " / " <<  cluster.key().second << " ] " << endl;
0303 
0304     // if requested, cut on cluster size
0305     if (pblocks_[subdet - 1].usesSize_ && (cluster->pixels().size() > pblocks_[subdet - 1].maxSize_))
0306       return;
0307 
0308     // mark as used
0309     //pixels[cluster.key()] = false;
0310     assert(collectedPixels_.size() > cluster.key());
0311     collectedPixels_[cluster.key()] = true;
0312   } else {  // aka Strip
0313     if (!doStrip_)
0314       return;
0315     const type_info &hitType = typeid(*hit);
0316     if (hitType == typeid(SiStripRecHit2D)) {
0317       const SiStripRecHit2D *stripHit = static_cast<const SiStripRecHit2D *>(hit);
0318       //DBG//     cout << "Plain RecHit 2D: " << endl;
0319       process(stripHit->omniClusterRef(), subdet);
0320       //      int clusCharge=0;
0321       //      for ( auto cAmp : stripHit->omniClusterRef().stripCluster().amplitudes() ) clusCharge+=cAmp;
0322       //      std::cout << "[HLTTrackClusterRemoverNew::process (SiStripRecHit2D) chi2: " << chi2 << " [" << subdet << " --> charge: " << clusCharge << "]" << std::endl;
0323     } else if (hitType == typeid(SiStripRecHit1D)) {
0324       const SiStripRecHit1D *hit1D = static_cast<const SiStripRecHit1D *>(hit);
0325       process(hit1D->omniClusterRef(), subdet);
0326       //      int clusCharge=0;
0327       //      for ( auto cAmp : hit1D->omniClusterRef().stripCluster().amplitudes() ) clusCharge+=cAmp;
0328       //      std::cout << "[HLTTrackClusterRemoverNew::process (SiStripRecHit1D) chi2: " << chi2 << " [" << subdet << " --> charge: " << clusCharge << "]" << std::endl;
0329     } else if (hitType == typeid(SiStripMatchedRecHit2D)) {
0330       const SiStripMatchedRecHit2D *matchHit = static_cast<const SiStripMatchedRecHit2D *>(hit);
0331       //DBG//     cout << "Matched RecHit 2D: " << endl;
0332       process(matchHit->monoClusterRef(), subdet);
0333       //        int clusCharge=0;
0334       //        for ( auto cAmp : matchHit->monoClusterRef().stripCluster().amplitudes() ) clusCharge+=cAmp;
0335       //        std::cout << "[HLTTrackClusterRemoverNew::process (SiStripMatchedRecHit2D:mono) chi2: " << chi2 << " [" << subdet << " --> charge: " << clusCharge << "]" << std::endl;
0336 
0337       process(matchHit->stereoClusterRef(), subdet);
0338       //        clusCharge=0;
0339       //        for ( auto cAmp : matchHit->stereoClusterRef().stripCluster().amplitudes() ) clusCharge+=cAmp;
0340       //        std::cout << "[HLTTrackClusterRemoverNew::process (SiStripMatchedRecHit2D:stereo) chi2: " << chi2 << " [" << subdet << " --> charge: " << clusCharge << "]" << std::endl;
0341 
0342     } else if (hitType == typeid(ProjectedSiStripRecHit2D)) {
0343       const ProjectedSiStripRecHit2D *projHit = static_cast<const ProjectedSiStripRecHit2D *>(hit);
0344       //DBG//     cout << "Projected RecHit 2D: " << endl;
0345       process(projHit->originalHit().omniClusterRef(), subdet);
0346       //        int clusCharge=0;
0347       //        for ( auto cAmp : projHit->originalHit().omniClusterRef().stripCluster().amplitudes() ) clusCharge+=cAmp;
0348       //        std::cout << "[HLTTrackClusterRemoverNew::process (ProjectedSiStripRecHit2D) chi2: " << chi2 << " [" << subdet << " --> charge: " << clusCharge << "]" << std::endl;
0349     } else
0350       throw cms::Exception("NOT IMPLEMENTED")
0351           << "Don't know how to handle " << hitType.name() << " on detid " << detid.rawId() << "\n";
0352   }
0353 }
0354 
0355 /*   Schematic picture of n-th step Iterative removal
0356  *   (that os removing clusters after n-th step tracking)
0357  *   clusters:    [ C1 ] -> [ C2 ] -> ... -> [ Cn ] -> [ Cn + 1 ]
0358  *                   ^                          ^           ^--- OUTPUT "new" ID
0359  *                   |-- before any removal     |----- Source clusters                                                                             
0360  *                   |-- OUTPUT "old" ID        |----- Hits in Traj. point here                       
0361  *                   |                          \----- Old ClusterRemovalInfo "new" ID                 
0362  *                   \-- Old ClusterRemovalInfo "old" ID                                                  
0363  */
0364 
0365 void HLTTrackClusterRemoverNew::produce(Event &iEvent, const EventSetup &iSetup) {
0366   ProductID pixelOldProdID, stripOldProdID;
0367 
0368   const auto &tgh = &iSetup.getData(tTrackerGeom_);
0369 
0370   edm::Handle<edmNew::DetSetVector<SiPixelCluster> > pixelClusters;
0371   if (doPixel_) {
0372     iEvent.getByToken(pixelClusters_, pixelClusters);
0373     pixelSourceProdID = pixelClusters.id();
0374   }
0375 
0376   edm::Handle<edmNew::DetSetVector<SiStripCluster> > stripClusters;
0377   if (doStrip_) {
0378     iEvent.getByToken(stripClusters_, stripClusters);
0379     stripSourceProdID = stripClusters.id();
0380   }
0381 
0382   //Handle<TrajTrackAssociationCollection> trajectories;
0383   edm::Handle<vector<Trajectory> > trajectories;
0384   iEvent.getByToken(trajectories_, trajectories);
0385 
0386   if (mergeOld_) {
0387     edm::Handle<PixelMaskContainer> oldPxlMask;
0388     edm::Handle<StripMaskContainer> oldStrMask;
0389     iEvent.getByToken(oldPxlMaskToken_, oldPxlMask);
0390     iEvent.getByToken(oldStrMaskToken_, oldStrMask);
0391     LogDebug("TrackClusterRemover") << "to merge in, " << oldStrMask->size() << " strp and " << oldPxlMask->size()
0392                                     << " pxl";
0393     oldStrMask->copyMaskTo(collectedRegStrips_);
0394     oldPxlMask->copyMaskTo(collectedPixels_);
0395     collectedRegStrips_.resize(stripClusters->dataSize(), false);
0396   } else {
0397     collectedRegStrips_.resize(stripClusters->dataSize(), false);
0398     collectedPixels_.resize(pixelClusters->dataSize(), false);
0399   }
0400 
0401   //for (TrajTrackAssociationCollection::const_iterator it = trajectories->begin(), ed = trajectories->end(); it != ed; ++it)  {
0402   //    const Trajectory &tj = * it->key;
0403 
0404   for (std::vector<Trajectory>::const_iterator it = trajectories->begin(), ed = trajectories->end(); it != ed; ++it) {
0405     const Trajectory &tj = *it;
0406     const std::vector<TrajectoryMeasurement> &tms = tj.measurements();
0407 
0408     std::vector<TrajectoryMeasurement>::const_iterator itm, endtm;
0409     for (itm = tms.begin(), endtm = tms.end(); itm != endtm; ++itm) {
0410       const TrackingRecHit *hit = itm->recHit()->hit();
0411       if (!hit->isValid())
0412         continue;
0413       //        std::cout<<"process hit"<<std::endl;
0414       process(hit, itm->estimate(), tgh);
0415     }
0416   }
0417 
0418   //    std::cout << " => collectedRegStrips_: " << collectedRegStrips_.size() << std::endl;
0419   //    std::cout << " total strip to skip (before charge check): "<<std::count(collectedRegStrips_.begin(),collectedRegStrips_.end(),true) << std::endl;
0420 
0421   // checks only cluster already found! creates regression
0422   if (doStripChargeCheck_) {
0423     //      std::cout << "[HLTTrackClusterRemoverNew::produce] doStripChargeCheck_: " << (doStripChargeCheck_ ? "true" : "false") << " stripClusters: " << stripClusters->size() << std::endl;
0424 
0425     auto const &clusters = stripClusters->data();
0426     for (auto const &item : stripClusters->ids()) {
0427       if (!item.isValid())
0428         continue;  // not umpacked
0429 
0430       DetId detid = item.id;
0431       uint32_t subdet = detid.subdetId();
0432       if (!pblocks_[subdet - 1].cutOnStripCharge_)
0433         continue;
0434 
0435       //    std::cout << " i: " << i << " --> detid: " << detid << " --> subdet: " << subdet << std::endl;
0436 
0437       for (int i = item.offset; i < item.offset + int(item.size); ++i) {
0438         int clusCharge = 0;
0439         for (auto cAmp : clusters[i].amplitudes())
0440           clusCharge += cAmp;
0441 
0442         //  if (clusCharge < pblocks_[subdet-1].minGoodStripCharge_) std::cout << " clusCharge: " << clusCharge << std::endl;
0443         if (clusCharge < pblocks_[subdet - 1].minGoodStripCharge_)
0444           collectedRegStrips_[i] = true;  // (|= does not work!)
0445       }
0446     }
0447   }
0448 
0449   //    std::cout << " => collectedRegStrips_: " << collectedRegStrips_.size() << std::endl;
0450   //    std::cout << " total strip to skip: "<<std::count(collectedRegStrips_.begin(),collectedRegStrips_.end(),true) << std::endl;
0451   LogDebug("TrackClusterRemover") << "total strip to skip: "
0452                                   << std::count(collectedRegStrips_.begin(), collectedRegStrips_.end(), true);
0453   // std::cout << "TrackClusterRemover" <<" total strip to skip: "<<std::count(collectedRegStrips_.begin(),collectedRegStrips_.end(),true)<<std::endl;
0454   iEvent.put(std::make_unique<StripMaskContainer>(edm::RefProd<edmNew::DetSetVector<SiStripCluster> >(stripClusters),
0455                                                   collectedRegStrips_));
0456 
0457   LogDebug("TrackClusterRemover") << "total pxl to skip: "
0458                                   << std::count(collectedPixels_.begin(), collectedPixels_.end(), true);
0459   iEvent.put(std::make_unique<PixelMaskContainer>(edm::RefProd<edmNew::DetSetVector<SiPixelCluster> >(pixelClusters),
0460                                                   collectedPixels_));
0461 
0462   collectedRegStrips_.clear();
0463   collectedPixels_.clear();
0464 }
0465 
0466 #include "FWCore/PluginManager/interface/ModuleDef.h"
0467 #include "FWCore/Framework/interface/MakerMacros.h"
0468 DEFINE_FWK_MODULE(HLTTrackClusterRemoverNew);