Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2025-01-22 07:34:25 

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 "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0009 
0010 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
0011 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0012 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
0013 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
0014 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
0015 #include "DataFormats/TrackerRecHit2D/interface/ProjectedSiStripRecHit2D.h"
0016 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h"
0017 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHitCollection.h"
0018 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2DCollection.h"
0019 #include "DataFormats/Common/interface/Handle.h"
0020 #include "DataFormats/SiStripDetId/interface/SiStripDetId.h"
0021 #include "DataFormats/Common/interface/DetSetVector.h"
0022 #include "DataFormats/Common/interface/ValueMap.h"
0023 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0024 #include "DataFormats/Provenance/interface/ProductID.h"
0025 #include "DataFormats/Common/interface/ContainerMask.h"
0026 
0027 #include "DataFormats/DetId/interface/DetId.h"
0028 
0029 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0030 #include "DataFormats/TrackReco/interface/Track.h"
0031 #include "DataFormats/TrackerRecHit2D/interface/ClusterRemovalInfo.h"
0032 
0033 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0034 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0035 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0036 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0037 
0038 #include "RecoTracker/TransientTrackingRecHit/interface/Traj2TrackHits.h"
0039 
0040 //
0041 // class decleration
0042 //
0043 
0044 class HITrackClusterRemover : public edm::stream::EDProducer<> {
0045 public:
0046   HITrackClusterRemover(const edm::ParameterSet &iConfig);
0047   ~HITrackClusterRemover() override;
0048   void produce(edm::Event &iEvent, const edm::EventSetup &iSetup) override;
0049   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0050 
0051 private:
0052   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> const tTrackerGeom_;
0053   struct ParamBlock {
0054     ParamBlock() : isSet_(false), usesCharge_(false) {}
0055     ParamBlock(const edm::ParameterSet &iConfig)
0056         : isSet_(true),
0057           maxChi2_(iConfig.getParameter<double>("maxChi2")),
0058           maxCharge_(iConfig.getParameter<double>("maxCharge")),
0059           minGoodPixelCharge_(iConfig.getParameter<double>("minGoodPixelCharge")),
0060           minGoodStripCharge_(iConfig.getParameter<double>("minGoodStripCharge")),
0061           maxSize_(iConfig.getParameter<uint32_t>("maxSize")),
0062           usesCharge_(maxCharge_ > 0.),
0063           usesSize_(maxSize_ > 0.),
0064           cutOnPixelCharge_(minGoodPixelCharge_ > 0.),
0065           cutOnStripCharge_(minGoodStripCharge_ > 0.) {}
0066     bool isSet_;
0067     float maxChi2_, maxCharge_, minGoodPixelCharge_, minGoodStripCharge_;
0068     size_t maxSize_;
0069     bool usesCharge_, usesSize_, cutOnPixelCharge_, cutOnStripCharge_;
0070 
0071     static void fillPSetDescription(edm::ParameterSetDescription &desc) {
0072       desc.addOptional<double>("maxChi2");
0073       desc.add<double>("maxCharge", 0.);
0074       desc.add<double>("minGoodPixelCharge", 0.);
0075       desc.add<double>("minGoodStripCharge", 0.);
0076       desc.add<uint32_t>("maxSize", 0.);
0077     }
0078   };
0079   static const unsigned int NumberOfParamBlocks = 6;
0080 
0081   bool doTracks_;
0082   bool doStrip_, doPixel_;
0083   bool mergeOld_;
0084   typedef edm::ContainerMask<edmNew::DetSetVector<SiPixelCluster> > PixelMaskContainer;
0085   typedef edm::ContainerMask<edmNew::DetSetVector<SiStripCluster> > StripMaskContainer;
0086   edm::EDGetTokenT<edmNew::DetSetVector<SiPixelCluster> > pixelClusters_;
0087   edm::EDGetTokenT<edmNew::DetSetVector<SiStripCluster> > stripClusters_;
0088   edm::EDGetTokenT<reco::TrackCollection> tracks_;
0089   edm::EDGetTokenT<reco::ClusterRemovalInfo> oldRemovalInfo_;
0090   edm::EDGetTokenT<PixelMaskContainer> oldPxlMaskToken_;
0091   edm::EDGetTokenT<StripMaskContainer> oldStrMaskToken_;
0092   std::vector<edm::EDGetTokenT<edm::ValueMap<int> > > overrideTrkQuals_;
0093   edm::EDGetTokenT<SiStripRecHit2DCollection> rphiRecHitToken_, stereoRecHitToken_;
0094   //    edm::EDGetTokenT<SiPixelRecHitCollection> pixelRecHitsToken_;
0095 
0096   ParamBlock pblocks_[NumberOfParamBlocks];
0097   void readPSet(const edm::ParameterSet &iConfig,
0098                 const std::string &name,
0099                 int id1 = -1,
0100                 int id2 = -1,
0101                 int id3 = -1,
0102                 int id4 = -1,
0103                 int id5 = -1,
0104                 int id6 = -1);
0105 
0106   std::vector<uint8_t> pixels, strips;                  // avoid unneed alloc/dealloc of this
0107   edm::ProductID pixelSourceProdID, stripSourceProdID;  // ProdIDs refs must point to (for consistency tests)
0108 
0109   inline void process(const TrackingRecHit *hit, unsigned char chi2, const TrackerGeometry *tg);
0110   inline void process(const OmniClusterRef &cluRef, SiStripDetId &detid, bool fromTrack);
0111 
0112   template <typename T>
0113   std::unique_ptr<edmNew::DetSetVector<T> > cleanup(const edmNew::DetSetVector<T> &oldClusters,
0114                                                     const std::vector<uint8_t> &isGood,
0115                                                     reco::ClusterRemovalInfo::Indices &refs,
0116                                                     const reco::ClusterRemovalInfo::Indices *oldRefs);
0117 
0118   // Carries in full removal info about a given det from oldRefs
0119   void mergeOld(reco::ClusterRemovalInfo::Indices &refs, const reco::ClusterRemovalInfo::Indices &oldRefs);
0120 
0121   bool clusterWasteSolution_, doStripChargeCheck_, doPixelChargeCheck_;
0122   std::string stripRecHits_, pixelRecHits_;
0123   bool filterTracks_;
0124   int minNumberOfLayersWithMeasBeforeFiltering_;
0125   reco::TrackBase::TrackQuality trackQuality_;
0126   std::vector<bool> collectedStrips_;
0127   std::vector<bool> collectedPixels_;
0128 
0129   float sensorThickness(const SiStripDetId &detid) const;
0130 };
0131 
0132 using namespace std;
0133 using namespace edm;
0134 using namespace reco;
0135 
0136 void HITrackClusterRemover::readPSet(
0137     const edm::ParameterSet &iConfig, const std::string &name, int id1, int id2, int id3, int id4, int id5, int id6) {
0138   if (iConfig.exists(name)) {
0139     ParamBlock pblock(iConfig.getParameter<ParameterSet>(name));
0140     if (id1 == -1) {
0141       fill(pblocks_, pblocks_ + NumberOfParamBlocks, pblock);
0142     } else {
0143       pblocks_[id1] = pblock;
0144       if (id2 != -1)
0145         pblocks_[id2] = pblock;
0146       if (id3 != -1)
0147         pblocks_[id3] = pblock;
0148       if (id4 != -1)
0149         pblocks_[id4] = pblock;
0150       if (id5 != -1)
0151         pblocks_[id5] = pblock;
0152       if (id6 != -1)
0153         pblocks_[id6] = pblock;
0154     }
0155   }
0156 }
0157 
0158 HITrackClusterRemover::HITrackClusterRemover(const ParameterSet &iConfig)
0159     : tTrackerGeom_(esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>()),
0160       doTracks_(!iConfig.getParameter<edm::InputTag>("trajectories").label().empty()),
0161       doStrip_(iConfig.getParameter<bool>("doStrip")),
0162       doPixel_(iConfig.getParameter<bool>("doPixel")),
0163       mergeOld_(!iConfig.getParameter<edm::InputTag>("oldClusterRemovalInfo").label().empty()),
0164       clusterWasteSolution_(true),
0165       doStripChargeCheck_(iConfig.getParameter<bool>("doStripChargeCheck")),
0166       doPixelChargeCheck_(iConfig.getParameter<bool>("doPixelChargeCheck")),
0167       stripRecHits_(doStripChargeCheck_ ? iConfig.getParameter<std::string>("stripRecHits")
0168                                         : std::string("siStripMatchedRecHits")),
0169       pixelRecHits_(doPixelChargeCheck_ ? iConfig.getParameter<std::string>("pixelRecHits")
0170                                         : std::string("siPixelRecHits")) {
0171   mergeOld_ = mergeOld_ && !iConfig.getParameter<InputTag>("oldClusterRemovalInfo").label().empty();
0172   if (!iConfig.getParameter<edm::InputTag>("overrideTrkQuals").label().empty())
0173     overrideTrkQuals_.push_back(consumes<edm::ValueMap<int> >(iConfig.getParameter<InputTag>("overrideTrkQuals")));
0174   clusterWasteSolution_ = !iConfig.getParameter<bool>("clusterLessSolution");
0175   if ((doPixelChargeCheck_ && !doPixel_) || (doStripChargeCheck_ && !doStrip_))
0176     throw cms::Exception("Configuration Error")
0177         << "HITrackClusterRemover: Charge check asked without cluster collection ";
0178   if (doPixelChargeCheck_)
0179     throw cms::Exception("Configuration Error")
0180         << "HITrackClusterRemover: Pixel cluster charge check not yet implemented";
0181 
0182   if (doPixel_ && clusterWasteSolution_)
0183     produces<edmNew::DetSetVector<SiPixelCluster> >();
0184   if (doStrip_ && clusterWasteSolution_)
0185     produces<edmNew::DetSetVector<SiStripCluster> >();
0186   if (clusterWasteSolution_)
0187     produces<ClusterRemovalInfo>();
0188 
0189   assert(!clusterWasteSolution_);
0190 
0191   fill(pblocks_, pblocks_ + NumberOfParamBlocks, ParamBlock());
0192   readPSet(iConfig, "Common", -1);
0193   if (doPixel_) {
0194     readPSet(iConfig, "Pixel", 0, 1);
0195     readPSet(iConfig, "PXB", 0);
0196     readPSet(iConfig, "PXE", 1);
0197   }
0198   if (doStrip_) {
0199     readPSet(iConfig, "Strip", 2, 3, 4, 5);
0200     readPSet(iConfig, "StripInner", 2, 3);
0201     readPSet(iConfig, "StripOuter", 4, 5);
0202     readPSet(iConfig, "TIB", 2);
0203     readPSet(iConfig, "TID", 3);
0204     readPSet(iConfig, "TOB", 4);
0205     readPSet(iConfig, "TEC", 5);
0206   }
0207 
0208   bool usingCharge = false;
0209   for (size_t i = 0; i < NumberOfParamBlocks; ++i) {
0210     if (!pblocks_[i].isSet_)
0211       throw cms::Exception("Configuration Error")
0212           << "HITrackClusterRemover: Missing configuration for detector with subDetID = " << (i + 1);
0213     if (pblocks_[i].usesCharge_ && !usingCharge) {
0214       throw cms::Exception("Configuration Error") << "HITrackClusterRemover: Configuration for subDetID = " << (i + 1)
0215                                                   << " uses cluster charge, which is not enabled.";
0216     }
0217   }
0218 
0219   if (!clusterWasteSolution_) {
0220     produces<edm::ContainerMask<edmNew::DetSetVector<SiPixelCluster> > >();
0221     produces<edm::ContainerMask<edmNew::DetSetVector<SiStripCluster> > >();
0222   }
0223   trackQuality_ = reco::TrackBase::undefQuality;
0224   filterTracks_ = false;
0225   if (!iConfig.getParameter<std::string>("TrackQuality").empty()) {
0226     filterTracks_ = true;
0227     trackQuality_ = reco::TrackBase::qualityByName(iConfig.getParameter<std::string>("TrackQuality"));
0228     minNumberOfLayersWithMeasBeforeFiltering_ = iConfig.getParameter<int>("minNumberOfLayersWithMeasBeforeFiltering");
0229   }
0230 
0231   if (doTracks_)
0232     tracks_ = consumes<reco::TrackCollection>(iConfig.getParameter<InputTag>("trajectories"));
0233   if (doPixel_)
0234     pixelClusters_ = consumes<edmNew::DetSetVector<SiPixelCluster> >(iConfig.getParameter<InputTag>("pixelClusters"));
0235   if (doStrip_)
0236     stripClusters_ = consumes<edmNew::DetSetVector<SiStripCluster> >(iConfig.getParameter<InputTag>("stripClusters"));
0237   if (mergeOld_) {
0238     oldRemovalInfo_ = consumes<ClusterRemovalInfo>(iConfig.getParameter<InputTag>("oldClusterRemovalInfo"));
0239     oldPxlMaskToken_ = consumes<PixelMaskContainer>(iConfig.getParameter<InputTag>("oldClusterRemovalInfo"));
0240     oldStrMaskToken_ = consumes<StripMaskContainer>(iConfig.getParameter<InputTag>("oldClusterRemovalInfo"));
0241   }
0242 
0243   if (doStripChargeCheck_) {
0244     rphiRecHitToken_ = consumes<SiStripRecHit2DCollection>(InputTag(stripRecHits_, "rphiRecHit"));
0245     stereoRecHitToken_ = consumes<SiStripRecHit2DCollection>(InputTag(stripRecHits_, "stereoRecHit"));
0246   }
0247   //    if(doPixelChargeCheck_) pixelRecHitsToken_ = consumes<SiPixelRecHitCollection>(InputTag(pixelRecHits_));
0248 }
0249 
0250 HITrackClusterRemover::~HITrackClusterRemover() {}
0251 
0252 void HITrackClusterRemover::mergeOld(ClusterRemovalInfo::Indices &refs, const ClusterRemovalInfo::Indices &oldRefs) {
0253   for (size_t i = 0, n = refs.size(); i < n; ++i) {
0254     refs[i] = oldRefs[refs[i]];
0255   }
0256 }
0257 
0258 template <typename T>
0259 std::unique_ptr<edmNew::DetSetVector<T> > HITrackClusterRemover::cleanup(
0260     const edmNew::DetSetVector<T> &oldClusters,
0261     const std::vector<uint8_t> &isGood,
0262     reco::ClusterRemovalInfo::Indices &refs,
0263     const reco::ClusterRemovalInfo::Indices *oldRefs) {
0264   typedef typename edmNew::DetSetVector<T> DSV;
0265   typedef typename edmNew::DetSetVector<T>::FastFiller DSF;
0266   typedef typename edmNew::DetSet<T> DS;
0267   auto output = std::make_unique<DSV>();
0268   output->reserve(oldClusters.size(), oldClusters.dataSize());
0269 
0270   // cluster removal loop
0271   const T *firstOffset = &oldClusters.data().front();
0272   for (typename DSV::const_iterator itdet = oldClusters.begin(), enddet = oldClusters.end(); itdet != enddet; ++itdet) {
0273     DS oldDS = *itdet;
0274 
0275     if (oldDS.empty())
0276       continue;  // skip empty detsets
0277 
0278     uint32_t id = oldDS.detId();
0279     DSF outds(*output, id);
0280 
0281     for (typename DS::const_iterator it = oldDS.begin(), ed = oldDS.end(); it != ed; ++it) {
0282       uint32_t index = ((&*it) - firstOffset);
0283       if (isGood[index]) {
0284         outds.push_back(*it);
0285         refs.push_back(index);
0286         //std::cout << "HITrackClusterRemover::cleanup " << typeid(T).name() << " reference " << index << " to " << (refs.size() - 1) << std::endl;
0287       }
0288     }
0289     if (outds.empty())
0290       outds.abort();  // not write in an empty DSV
0291   }
0292 
0293   if (oldRefs != nullptr)
0294     mergeOld(refs, *oldRefs);
0295   return output;
0296 }
0297 
0298 float HITrackClusterRemover::sensorThickness(const SiStripDetId &detid) const {
0299   if (detid.subdetId() >= SiStripDetId::TIB) {
0300     if (detid.subdetId() == SiStripDetId::TOB)
0301       return 0.047;
0302     if (detid.moduleGeometry() == SiStripModuleGeometry::W5 || detid.moduleGeometry() == SiStripModuleGeometry::W6 ||
0303         detid.moduleGeometry() == SiStripModuleGeometry::W7)
0304       return 0.047;
0305     return 0.029;  // so it is TEC ring 1-4 or TIB or TOB;
0306   } else if (detid.subdetId() == PixelSubdetector::PixelBarrel)
0307     return 0.0285;
0308   else
0309     return 0.027;
0310 }
0311 
0312 void HITrackClusterRemover::process(OmniClusterRef const &ocluster, SiStripDetId &detid, bool fromTrack) {
0313   SiStripRecHit2D::ClusterRef cluster = ocluster.cluster_strip();
0314   if (cluster.id() != stripSourceProdID)
0315     throw cms::Exception("Inconsistent Data")
0316         << "HITrackClusterRemover: strip cluster ref from Product ID = " << cluster.id()
0317         << " does not match with source cluster collection (ID = " << stripSourceProdID << ")\n.";
0318 
0319   uint32_t subdet = detid.subdetId();
0320   assert(cluster.id() == stripSourceProdID);
0321   if (pblocks_[subdet - 1].usesSize_ && (cluster->amplitudes().size() > pblocks_[subdet - 1].maxSize_))
0322     return;
0323   if (!fromTrack) {
0324     if (pblocks_[subdet - 1].cutOnStripCharge_ &&
0325         (cluster->charge() > (pblocks_[subdet - 1].minGoodStripCharge_ * sensorThickness(detid))))
0326       return;
0327   }
0328 
0329   if (collectedStrips_.size() <= cluster.key())
0330     edm::LogError("BadCollectionSize") << collectedStrips_.size() << " is smaller than " << cluster.key();
0331 
0332   assert(collectedStrips_.size() > cluster.key());
0333   strips[cluster.key()] = false;
0334   if (!clusterWasteSolution_)
0335     collectedStrips_[cluster.key()] = true;
0336 }
0337 
0338 void HITrackClusterRemover::process(const TrackingRecHit *hit, unsigned char chi2, const TrackerGeometry *tg) {
0339   SiStripDetId detid = hit->geographicalId();
0340   uint32_t subdet = detid.subdetId();
0341 
0342   assert((subdet > 0) && (subdet <= NumberOfParamBlocks));
0343 
0344   // chi2 cut
0345   if (chi2 > Traj2TrackHits::toChi2x5(pblocks_[subdet - 1].maxChi2_))
0346     return;
0347 
0348   if (GeomDetEnumerators::isTrackerPixel(tg->geomDetSubDetector(subdet))) {
0349     if (!doPixel_)
0350       return;
0351     // this is a pixel, and i *know* it is
0352     const SiPixelRecHit *pixelHit = static_cast<const SiPixelRecHit *>(hit);
0353 
0354     SiPixelRecHit::ClusterRef cluster = pixelHit->cluster();
0355 
0356     if (cluster.id() != pixelSourceProdID)
0357       throw cms::Exception("Inconsistent Data")
0358           << "HITrackClusterRemover: pixel cluster ref from Product ID = " << cluster.id()
0359           << " does not match with source cluster collection (ID = " << pixelSourceProdID << ")\n.";
0360 
0361     assert(cluster.id() == pixelSourceProdID);
0362     //DBG// cout << "HIT NEW PIXEL DETID = " << detid.rawId() << ", Cluster [ " << cluster.key().first << " / " <<  cluster.key().second << " ] " << endl;
0363 
0364     // if requested, cut on cluster size
0365     if (pblocks_[subdet - 1].usesSize_ && (cluster->pixels().size() > pblocks_[subdet - 1].maxSize_))
0366       return;
0367 
0368     // mark as used
0369     pixels[cluster.key()] = false;
0370 
0371     //if(!clusterWasteSolution_) collectedPixel[detid.rawId()].insert(cluster);
0372     assert(collectedPixels_.size() > cluster.key());
0373     //assert(detid.rawId() == cluster->geographicalId()); //This condition fails
0374     if (!clusterWasteSolution_)
0375       collectedPixels_[cluster.key()] = true;
0376 
0377   } else {  // aka Strip
0378     if (!doStrip_)
0379       return;
0380     const type_info &hitType = typeid(*hit);
0381     if (hitType == typeid(SiStripRecHit2D)) {
0382       const SiStripRecHit2D *stripHit = static_cast<const SiStripRecHit2D *>(hit);
0383       //DBG//     cout << "Plain RecHit 2D: " << endl;
0384       process(stripHit->omniClusterRef(), detid, true);
0385     } else if (hitType == typeid(SiStripRecHit1D)) {
0386       const SiStripRecHit1D *hit1D = static_cast<const SiStripRecHit1D *>(hit);
0387       process(hit1D->omniClusterRef(), detid, true);
0388     } else if (hitType == typeid(SiStripMatchedRecHit2D)) {
0389       const SiStripMatchedRecHit2D *matchHit = static_cast<const SiStripMatchedRecHit2D *>(hit);
0390       //DBG//     cout << "Matched RecHit 2D: " << endl;
0391       process(matchHit->monoClusterRef(), detid, true);
0392       process(matchHit->stereoClusterRef(), detid, true);
0393     } else if (hitType == typeid(ProjectedSiStripRecHit2D)) {
0394       const ProjectedSiStripRecHit2D *projHit = static_cast<const ProjectedSiStripRecHit2D *>(hit);
0395       //DBG//     cout << "Projected RecHit 2D: " << endl;
0396       process(projHit->originalHit().omniClusterRef(), detid, true);
0397     } else
0398       throw cms::Exception("NOT IMPLEMENTED")
0399           << "Don't know how to handle " << hitType.name() << " on detid " << detid.rawId() << "\n";
0400   }
0401 }
0402 
0403 /*   Schematic picture of n-th step Iterative removal
0404  *   (that os removing clusters after n-th step tracking)
0405  *   clusters:    [ C1 ] -> [ C2 ] -> ... -> [ Cn ] -> [ Cn + 1 ]
0406  *                   ^                          ^           ^--- OUTPUT "new" ID
0407  *                   |-- before any removal     |----- Source clusters                                                                             
0408  *                   |-- OUTPUT "old" ID        |----- Hits in Traj. point here                       
0409  *                   |                          \----- Old ClusterRemovalInfo "new" ID                 
0410  *                   \-- Old ClusterRemovalInfo "old" ID                                                  
0411  */
0412 
0413 void HITrackClusterRemover::produce(Event &iEvent, const EventSetup &iSetup) {
0414   ProductID pixelOldProdID, stripOldProdID;
0415 
0416   const auto &tgh = &iSetup.getData(tTrackerGeom_);
0417 
0418   Handle<edmNew::DetSetVector<SiPixelCluster> > pixelClusters;
0419   if (doPixel_) {
0420     iEvent.getByToken(pixelClusters_, pixelClusters);
0421     pixelSourceProdID = pixelClusters.id();
0422   }
0423   //DBG// std::cout << "HITrackClusterRemover: Read pixel " << pixelClusters_.encode() << " = ID " << pixelSourceProdID << std::endl;
0424 
0425   Handle<edmNew::DetSetVector<SiStripCluster> > stripClusters;
0426   if (doStrip_) {
0427     iEvent.getByToken(stripClusters_, stripClusters);
0428     stripSourceProdID = stripClusters.id();
0429   }
0430   //DBG// std::cout << "HITrackClusterRemover: Read strip " << stripClusters_.encode() << " = ID " << stripSourceProdID << std::endl;
0431 
0432   std::unique_ptr<ClusterRemovalInfo> cri;
0433   if (clusterWasteSolution_) {
0434     if (doStrip_ && doPixel_)
0435       cri = std::make_unique<ClusterRemovalInfo>(pixelClusters, stripClusters);
0436     else if (doStrip_)
0437       cri = std::make_unique<ClusterRemovalInfo>(stripClusters);
0438     else if (doPixel_)
0439       cri = std::make_unique<ClusterRemovalInfo>(pixelClusters);
0440   }
0441 
0442   Handle<ClusterRemovalInfo> oldRemovalInfo;
0443   if (mergeOld_ && clusterWasteSolution_) {
0444     iEvent.getByToken(oldRemovalInfo_, oldRemovalInfo);
0445     // Check ProductIDs
0446     if ((oldRemovalInfo->stripNewRefProd().id() == stripClusters.id()) &&
0447         (oldRemovalInfo->pixelNewRefProd().id() == pixelClusters.id())) {
0448       cri->getOldClustersFrom(*oldRemovalInfo);
0449 
0450       pixelOldProdID = oldRemovalInfo->pixelRefProd().id();
0451       stripOldProdID = oldRemovalInfo->stripRefProd().id();
0452 
0453     } else {
0454       edm::EDConsumerBase::Labels labels;
0455       labelsForToken(oldRemovalInfo_, labels);
0456       throw cms::Exception("Inconsistent Data")
0457           << "HITrackClusterRemover: "
0458           << "Input collection product IDs are [pixel: " << pixelClusters.id() << ", strip: " << stripClusters.id()
0459           << "] \n"
0460           << "\t but the *old* ClusterRemovalInfo " << labels.productInstance << " refers as 'new product ids' to "
0461           << "[pixel: " << oldRemovalInfo->pixelNewRefProd().id()
0462           << ", strip: " << oldRemovalInfo->stripNewRefProd().id() << "]\n"
0463           << "NOTA BENE: when running HITrackClusterRemover with an old ClusterRemovalInfo the hits in the trajectory "
0464              "MUST be already re-keyed.\n";
0465     }
0466   } else {  // then Old == Source
0467     pixelOldProdID = pixelSourceProdID;
0468     stripOldProdID = stripSourceProdID;
0469   }
0470 
0471   if (doStrip_) {
0472     strips.resize(stripClusters->dataSize());
0473     fill(strips.begin(), strips.end(), true);
0474   }
0475   if (doPixel_) {
0476     pixels.resize(pixelClusters->dataSize());
0477     fill(pixels.begin(), pixels.end(), true);
0478   }
0479   if (mergeOld_) {
0480     edm::Handle<PixelMaskContainer> oldPxlMask;
0481     edm::Handle<StripMaskContainer> oldStrMask;
0482     iEvent.getByToken(oldPxlMaskToken_, oldPxlMask);
0483     iEvent.getByToken(oldStrMaskToken_, oldStrMask);
0484     LogDebug("HITrackClusterRemover") << "to merge in, " << oldStrMask->size() << " strp and " << oldPxlMask->size()
0485                                       << " pxl";
0486     oldStrMask->copyMaskTo(collectedStrips_);
0487     oldPxlMask->copyMaskTo(collectedPixels_);
0488     assert(stripClusters->dataSize() >= collectedStrips_.size());
0489     collectedStrips_.resize(stripClusters->dataSize(), false);
0490   } else {
0491     collectedStrips_.resize(stripClusters->dataSize(), false);
0492     collectedPixels_.resize(pixelClusters->dataSize(), false);
0493   }
0494 
0495   if (doTracks_) {
0496     Handle<reco::TrackCollection> tracks;
0497     iEvent.getByToken(tracks_, tracks);
0498 
0499     std::vector<Handle<edm::ValueMap<int> > > quals;
0500     if (!overrideTrkQuals_.empty()) {
0501       quals.resize(1);
0502       iEvent.getByToken(overrideTrkQuals_[0], quals[0]);
0503     }
0504     int it = 0;
0505     for (const auto &track : *tracks) {
0506       if (filterTracks_) {
0507         bool goodTk = true;
0508         if (!quals.empty()) {
0509           int qual = (*(quals[0])).get(it++);
0510           if (qual < 0) {
0511             goodTk = false;
0512           }
0513           //note that this does not work for some trackquals (goodIterative  or undefQuality)
0514           else
0515             goodTk = (qual & (1 << trackQuality_)) >> trackQuality_;
0516         } else
0517           goodTk = (track.quality(trackQuality_));
0518         if (!goodTk)
0519           continue;
0520         if (track.hitPattern().trackerLayersWithMeasurement() < minNumberOfLayersWithMeasBeforeFiltering_)
0521           continue;
0522       }
0523       auto const &chi2sX5 = track.extra()->chi2sX5();
0524       assert(chi2sX5.size() == track.recHitsSize());
0525       auto hb = track.recHitsBegin();
0526       for (unsigned int h = 0; h < track.recHitsSize(); h++) {
0527         auto hit = *(hb + h);
0528         if (!hit->isValid())
0529           continue;
0530         process(hit, chi2sX5[h], tgh);
0531       }
0532     }
0533   }
0534 
0535   if (doStripChargeCheck_) {
0536     edm::Handle<SiStripRecHit2DCollection> rechitsrphi;
0537     iEvent.getByToken(rphiRecHitToken_, rechitsrphi);
0538     const SiStripRecHit2DCollection::DataContainer *rphiRecHits = &(rechitsrphi).product()->data();
0539     for (SiStripRecHit2DCollection::DataContainer::const_iterator recHit = rphiRecHits->begin();
0540          recHit != rphiRecHits->end();
0541          recHit++) {
0542       SiStripDetId detid = recHit->geographicalId();
0543       process(recHit->omniClusterRef(), detid, false);
0544     }
0545     edm::Handle<SiStripRecHit2DCollection> rechitsstereo;
0546     iEvent.getByToken(stereoRecHitToken_, rechitsstereo);
0547     const SiStripRecHit2DCollection::DataContainer *stereoRecHits = &(rechitsstereo).product()->data();
0548     for (SiStripRecHit2DCollection::DataContainer::const_iterator recHit = stereoRecHits->begin();
0549          recHit != stereoRecHits->end();
0550          recHit++) {
0551       SiStripDetId detid = recHit->geographicalId();
0552       process(recHit->omniClusterRef(), detid, false);
0553     }
0554   }
0555   //    if(doPixelChargeCheck_) {
0556   //    edm::Handle<SiPixelRecHitCollection> pixelrechits;
0557   //    iEvent.getByToken(pixelRecHitsToken_,pixelrechits);
0558   //    }
0559 
0560   if (doPixel_ && clusterWasteSolution_) {
0561     OrphanHandle<edmNew::DetSetVector<SiPixelCluster> > newPixels = iEvent.put(
0562         cleanup(*pixelClusters, pixels, cri->pixelIndices(), mergeOld_ ? &oldRemovalInfo->pixelIndices() : nullptr));
0563     //DBG// std::cout << "HITrackClusterRemover: Wrote pixel " << newPixels.id() << " from " << pixelSourceProdID << std::endl;
0564     cri->setNewPixelClusters(newPixels);
0565   }
0566   if (doStrip_ && clusterWasteSolution_) {
0567     OrphanHandle<edmNew::DetSetVector<SiStripCluster> > newStrips = iEvent.put(
0568         cleanup(*stripClusters, strips, cri->stripIndices(), mergeOld_ ? &oldRemovalInfo->stripIndices() : nullptr));
0569     //DBG// std::cout << "HITrackClusterRemover: Wrote strip " << newStrips.id() << " from " << stripSourceProdID << std::endl;
0570     cri->setNewStripClusters(newStrips);
0571   }
0572 
0573   if (clusterWasteSolution_) {
0574     //      double fraction_pxl= cri->pixelIndices().size() / (double) pixels.size();
0575     //      double fraction_strp= cri->stripIndices().size() / (double) strips.size();
0576     //      edm::LogWarning("HITrackClusterRemover")<<" fraction: " << fraction_pxl <<" "<<fraction_strp;
0577     iEvent.put(std::move(cri));
0578   }
0579 
0580   pixels.clear();
0581   strips.clear();
0582 
0583   if (!clusterWasteSolution_) {
0584     //auto_ptr<edmNew::DetSetVector<SiPixelClusterRefNew> > removedPixelClsuterRefs(new edmNew::DetSetVector<SiPixelClusterRefNew>());
0585     //auto_ptr<edmNew::DetSetVector<SiStripRecHit1D::ClusterRef> > removedStripClsuterRefs(new edmNew::DetSetVector<SiStripRecHit1D::ClusterRef>());
0586 
0587     LogDebug("HITrackClusterRemover") << "total strip to skip: "
0588                                       << std::count(collectedStrips_.begin(), collectedStrips_.end(), true);
0589     // std::cout << "HITrackClusterRemover " <<"total strip to skip: "<<std::count(collectedStrips_.begin(),collectedStrips_.end(),true) <<std::endl;
0590     iEvent.put(std::make_unique<StripMaskContainer>(edm::RefProd<edmNew::DetSetVector<SiStripCluster> >(stripClusters),
0591                                                     collectedStrips_));
0592 
0593     LogDebug("HITrackClusterRemover") << "total pxl to skip: "
0594                                       << std::count(collectedPixels_.begin(), collectedPixels_.end(), true);
0595     iEvent.put(std::make_unique<PixelMaskContainer>(edm::RefProd<edmNew::DetSetVector<SiPixelCluster> >(pixelClusters),
0596                                                     collectedPixels_));
0597   }
0598   collectedStrips_.clear();
0599   collectedPixels_.clear();
0600 }
0601 
0602 void HITrackClusterRemover::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0603   edm::ParameterSetDescription desc;
0604 
0605   // Define the structure for optional ParameterSets like Common, Pixel, Strip, etc.
0606   edm::ParameterSetDescription paramBlockDesc;
0607   ParamBlock::fillPSetDescription(paramBlockDesc);
0608 
0609   // Add all possible ParameterSets (optional) using the above structure
0610   desc.addOptional<edm::ParameterSetDescription>("Common", paramBlockDesc);
0611   desc.addOptional<edm::ParameterSetDescription>("Pixel", paramBlockDesc);
0612   desc.addOptional<edm::ParameterSetDescription>("Strip", paramBlockDesc);
0613   desc.addOptional<edm::ParameterSetDescription>("PXB", paramBlockDesc);
0614   desc.addOptional<edm::ParameterSetDescription>("PXE", paramBlockDesc);
0615   desc.addOptional<edm::ParameterSetDescription>("StripInner", paramBlockDesc);
0616   desc.addOptional<edm::ParameterSetDescription>("StripOuter", paramBlockDesc);
0617   desc.addOptional<edm::ParameterSetDescription>("TIB", paramBlockDesc);
0618   desc.addOptional<edm::ParameterSetDescription>("TID", paramBlockDesc);
0619   desc.addOptional<edm::ParameterSetDescription>("TOB", paramBlockDesc);
0620   desc.addOptional<edm::ParameterSetDescription>("TEC", paramBlockDesc);
0621 
0622   desc.add<bool>("doStrip", true);
0623   desc.add<bool>("doPixel", true);
0624   desc.add<bool>("doStripChargeCheck", false);
0625   desc.add<bool>("doPixelChargeCheck", false);
0626   desc.add<std::string>("stripRecHits");
0627   desc.add<std::string>("pixelRecHits");
0628   desc.add<InputTag>("oldClusterRemovalInfo", edm::InputTag(""));
0629   desc.add<InputTag>("overrideTrkQuals", edm::InputTag(""));
0630   desc.add<bool>("clusterLessSolution", false);
0631   desc.add<std::string>("TrackQuality", "");
0632   desc.add<int>("minNumberOfLayersWithMeasBeforeFiltering", 0);
0633   desc.add<InputTag>("trajectories", edm::InputTag(""));
0634   desc.add<InputTag>("pixelClusters");
0635   desc.add<InputTag>("stripClusters");
0636   descriptions.addDefault(desc);
0637 }
0638 
0639 #include "FWCore/PluginManager/interface/ModuleDef.h"
0640 #include "FWCore/Framework/interface/MakerMacros.h"
0641 DEFINE_FWK_MODULE(HITrackClusterRemover);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!