Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2023-03-17 11:23:02 

0001 #include "RecoTracker/TkSeedingLayers/interface/SeedingLayerSetsBuilder.h"
0002 #include "FWCore/Framework/interface/Event.h"
0003 #include "FWCore/Framework/interface/EventSetup.h"
0004 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0005 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0006 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0007 #include "RecoTracker/TkDetLayers/interface/GeometricSearchTracker.h"
0008 #include "RecoTracker/Record/interface/TrackerRecoGeometryRecord.h"
0009 
0010 #include "TrackingTools/Records/interface/TransientRecHitRecord.h"
0011 #include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHit.h"
0012 #include "RecoTracker/TransientTrackingRecHit/interface/TkTransientTrackingRecHitBuilder.h"
0013 
0014 #include "TrackingTools/DetLayers/interface/BarrelDetLayer.h"
0015 #include "TrackingTools/DetLayers/interface/ForwardDetLayer.h"
0016 
0017 #include "RecoLocalTracker/SiStripClusterizer/interface/ClusterChargeCut.h"
0018 
0019 #include "HitExtractorPIX.h"
0020 #include "HitExtractorSTRP.h"
0021 
0022 #include <iostream>
0023 #include <sstream>
0024 #include <ostream>
0025 #include <fstream>
0026 #include <map>
0027 
0028 using namespace ctfseeding;
0029 using namespace std;
0030 
0031 SeedingLayerSetsBuilder::SeedingLayerId SeedingLayerSetsBuilder::nameToEnumId(const std::string& name) {
0032   GeomDetEnumerators::SubDetector subdet = GeomDetEnumerators::invalidDet;
0033   TrackerDetSide side = TrackerDetSide::Barrel;
0034   int idLayer = 0;
0035 
0036   size_t index;
0037   //
0038   // BPIX
0039   //
0040   if ((index = name.find("BPix")) != string::npos) {
0041     subdet = GeomDetEnumerators::PixelBarrel;
0042     side = TrackerDetSide::Barrel;
0043     idLayer = atoi(name.substr(index + 4, 1).c_str());
0044   }
0045   //
0046   // FPIX
0047   //
0048   else if ((index = name.find("FPix")) != string::npos) {
0049     subdet = GeomDetEnumerators::PixelEndcap;
0050     idLayer = atoi(name.substr(index + 4).c_str());
0051     if (name.find("pos") != string::npos) {
0052       side = TrackerDetSide::PosEndcap;
0053     } else {
0054       side = TrackerDetSide::NegEndcap;
0055     }
0056   }
0057   //
0058   // TIB
0059   //
0060   else if ((index = name.find("TIB")) != string::npos) {
0061     subdet = GeomDetEnumerators::TIB;
0062     side = TrackerDetSide::Barrel;
0063     idLayer = atoi(name.substr(index + 3, 1).c_str());
0064   }
0065   //
0066   // TID
0067   //
0068   else if ((index = name.find("TID")) != string::npos) {
0069     subdet = GeomDetEnumerators::TID;
0070     idLayer = atoi(name.substr(index + 3, 1).c_str());
0071     if (name.find("pos") != string::npos) {
0072       side = TrackerDetSide::PosEndcap;
0073     } else {
0074       side = TrackerDetSide::NegEndcap;
0075     }
0076   }
0077   //
0078   // TOB
0079   //
0080   else if ((index = name.find("TOB")) != string::npos) {
0081     subdet = GeomDetEnumerators::TOB;
0082     side = TrackerDetSide::Barrel;
0083     idLayer = atoi(name.substr(index + 3, 1).c_str());
0084   }
0085   //
0086   // TEC
0087   //
0088   else if ((index = name.find("TEC")) != string::npos) {
0089     subdet = GeomDetEnumerators::TEC;
0090     idLayer = atoi(name.substr(index + 3, 1).c_str());
0091     if (name.find("pos") != string::npos) {
0092       side = TrackerDetSide::PosEndcap;
0093     } else {
0094       side = TrackerDetSide::NegEndcap;
0095     }
0096   }
0097   return std::make_tuple(subdet, side, idLayer);
0098 }
0099 
0100 SeedingLayerSetsBuilder::LayerSpec::LayerSpec(unsigned short index,
0101                                               const std::string& layerName,
0102                                               const edm::ParameterSet& cfgLayer,
0103                                               edm::ConsumesCollector& iC)
0104     : nameIndex(index),
0105       hitBuilder(cfgLayer.getParameter<string>("TTRHBuilder")),
0106       hitBuilderToken(iC.esConsumes(edm::ESInputTag("", hitBuilder))) {
0107   usePixelHitProducer = false;
0108   if (cfgLayer.exists("HitProducer")) {
0109     pixelHitProducer = cfgLayer.getParameter<string>("HitProducer");
0110     usePixelHitProducer = true;
0111   }
0112 
0113   bool skipClusters = cfgLayer.exists("skipClusters");
0114   if (skipClusters) {
0115     LogDebug("SeedingLayerSetsBuilder") << layerName << " ready for skipping";
0116   } else {
0117     LogDebug("SeedingLayerSetsBuilder") << layerName << " not skipping ";
0118   }
0119 
0120   auto subdetData = nameToEnumId(layerName);
0121   subdet = std::get<0>(subdetData);
0122   side = std::get<1>(subdetData);
0123   idLayer = std::get<2>(subdetData);
0124   if (subdet == GeomDetEnumerators::PixelBarrel || subdet == GeomDetEnumerators::PixelEndcap) {
0125     extractor = std::make_unique<HitExtractorPIX>(side, idLayer, pixelHitProducer, iC);
0126   } else if (subdet != GeomDetEnumerators::invalidDet) {  // strip
0127     auto extr = std::make_unique<HitExtractorSTRP>(subdet, side, idLayer, clusterChargeCut(cfgLayer), iC);
0128     if (cfgLayer.exists("matchedRecHits")) {
0129       extr->useMatchedHits(cfgLayer.getParameter<edm::InputTag>("matchedRecHits"), iC);
0130     }
0131     if (cfgLayer.exists("rphiRecHits")) {
0132       extr->useRPhiHits(cfgLayer.getParameter<edm::InputTag>("rphiRecHits"), iC);
0133     }
0134     if (cfgLayer.exists("stereoRecHits")) {
0135       extr->useStereoHits(cfgLayer.getParameter<edm::InputTag>("stereoRecHits"), iC);
0136     }
0137     if (cfgLayer.exists("vectorRecHits")) {
0138       extr->useVectorHits(cfgLayer.getParameter<edm::InputTag>("vectorRecHits"), iC);
0139     }
0140     if (cfgLayer.exists("useRingSlector") && cfgLayer.getParameter<bool>("useRingSlector")) {
0141       extr->useRingSelector(cfgLayer.getParameter<int>("minRing"), cfgLayer.getParameter<int>("maxRing"));
0142     }
0143     bool useSimpleRphiHitsCleaner =
0144         cfgLayer.exists("useSimpleRphiHitsCleaner") ? cfgLayer.getParameter<bool>("useSimpleRphiHitsCleaner") : true;
0145     extr->useSimpleRphiHitsCleaner(useSimpleRphiHitsCleaner);
0146 
0147     double minAbsZ = cfgLayer.exists("MinAbsZ") ? cfgLayer.getParameter<double>("MinAbsZ") : 0.;
0148     if (minAbsZ > 0.) {
0149       extr->setMinAbsZ(minAbsZ);
0150     }
0151     if (skipClusters) {
0152       bool useProjection = cfgLayer.exists("useProjection") ? cfgLayer.getParameter<bool>("useProjection") : false;
0153       if (useProjection) {
0154         LogDebug("SeedingLayerSetsBuilder") << layerName << " will project partially masked matched rechit";
0155       } else {
0156         extr->setNoProjection();
0157       }
0158     }
0159     extractor = std::move(extr);
0160   }
0161   if (extractor && skipClusters) {
0162     extractor->useSkipClusters(cfgLayer.getParameter<edm::InputTag>("skipClusters"), iC);
0163   }
0164 }
0165 
0166 std::string SeedingLayerSetsBuilder::LayerSpec::print(const std::vector<std::string>& names) const {
0167   std::ostringstream str;
0168   str << "Layer=" << names[nameIndex] << ", hitBldr: " << hitBuilder;
0169 
0170   str << ", useRingSelector: ";
0171   HitExtractorSTRP* ext = nullptr;
0172   if ((ext = dynamic_cast<HitExtractorSTRP*>(extractor.get())) && ext->useRingSelector()) {
0173     auto minMaxRing = ext->getMinMaxRing();
0174     str << "true,"
0175         << " Rings: (" << std::get<0>(minMaxRing) << "," << std::get<1>(minMaxRing) << ")";
0176   } else
0177     str << "false";
0178 
0179   return str.str();
0180 }
0181 //FastSim specific constructor
0182 SeedingLayerSetsBuilder::SeedingLayerSetsBuilder(const edm::ParameterSet& cfg,
0183                                                  edm::ConsumesCollector& iC,
0184                                                  const edm::InputTag& fastsimHitTag)
0185     : SeedingLayerSetsBuilder(cfg, iC) {
0186   fastSimrecHitsToken_ = iC.consumes<FastTrackerRecHitCollection>(fastsimHitTag);
0187   trackerTopologyToken_ = iC.esConsumes();
0188 }
0189 SeedingLayerSetsBuilder::SeedingLayerSetsBuilder(const edm::ParameterSet& cfg, edm::ConsumesCollector&& iC)
0190     : SeedingLayerSetsBuilder(cfg, iC) {}
0191 SeedingLayerSetsBuilder::SeedingLayerSetsBuilder(const edm::ParameterSet& cfg, edm::ConsumesCollector& iC)
0192     : trackerToken_(iC.esConsumes()) {
0193   std::vector<std::string> namesPset = cfg.getParameter<std::vector<std::string> >("layerList");
0194   std::vector<std::vector<std::string> > layerNamesInSets = this->layerNamesInSets(namesPset);
0195   // debug printout of layers
0196   typedef std::vector<std::string>::const_iterator IS;
0197   typedef std::vector<std::vector<std::string> >::const_iterator IT;
0198   std::ostringstream str;
0199   // The following should not be set to cout
0200   //  for (IT it = layerNamesInSets.begin(); it != layerNamesInSets.end(); it++) {
0201   //    str << "SET: ";
0202   //    for (IS is = it->begin(); is != it->end(); is++)  str << *is <<" ";
0203   //    str << std::endl;
0204   //  }
0205   //  std::cout << str.str() << std::endl;
0206   if (layerNamesInSets.empty())
0207     theNumberOfLayersInSet = 0;
0208   else
0209     theNumberOfLayersInSet = layerNamesInSets[0].size();
0210 
0211   for (IT it = layerNamesInSets.begin(); it != layerNamesInSets.end(); it++) {
0212     if (it->size() != theNumberOfLayersInSet)
0213       throw cms::Exception("Configuration")
0214           << "Assuming all SeedingLayerSets to have same number of layers. LayerSet " << (it - layerNamesInSets.begin())
0215           << " has " << it->size() << " while 0th has " << theNumberOfLayersInSet;
0216     for (const std::string& layerName : *it) {
0217       auto found = std::find(theLayerNames.begin(), theLayerNames.end(), layerName);
0218       unsigned short layerIndex = 0;
0219       if (found != theLayerNames.end()) {
0220         layerIndex = found - theLayerNames.begin();
0221       } else {
0222         if (std::numeric_limits<unsigned short>::max() == theLayerNames.size()) {
0223           throw cms::Exception("Assert")
0224               << "Too many layers in " << __FILE__ << ":" << __LINE__
0225               << ", we may have to enlarge the index type from unsigned short to unsigned int";
0226         }
0227 
0228         layerIndex = theLayers.size();
0229         theLayers.emplace_back(theLayerNames.size(), layerName, layerConfig(layerName, cfg), iC);
0230         theLayerNames.push_back(layerName);
0231       }
0232       theLayerSetIndices.push_back(layerIndex);
0233     }
0234   }
0235   theLayerDets.resize(theLayers.size());
0236   theTTRHBuilders.resize(theLayers.size());
0237 
0238   // debug printout
0239   // The following should not be set to cout
0240   //for(const LayerSpec& layer: theLayers) {
0241   //  std::cout << layer.print(theLayerNames) << std::endl;
0242   //}
0243 }
0244 
0245 SeedingLayerSetsBuilder::~SeedingLayerSetsBuilder() {}
0246 
0247 void SeedingLayerSetsBuilder::fillDescriptions(edm::ParameterSetDescription& desc) {
0248   edm::ParameterSetDescription empty;
0249   empty.setAllowAnything();  // for now accept any parameter in the PSets, consider improving later
0250 
0251   desc.add<std::vector<std::string> >("layerList", {});
0252   desc.add<edm::ParameterSetDescription>("BPix", empty);
0253   desc.add<edm::ParameterSetDescription>("FPix", empty);
0254   desc.add<edm::ParameterSetDescription>("TIB", empty);
0255   desc.add<edm::ParameterSetDescription>("TID", empty);
0256   desc.add<edm::ParameterSetDescription>("TOB", empty);
0257   desc.add<edm::ParameterSetDescription>("TEC", empty);
0258   desc.add<edm::ParameterSetDescription>("MTIB", empty);
0259   desc.add<edm::ParameterSetDescription>("MTID", empty);
0260   desc.add<edm::ParameterSetDescription>("MTOB", empty);
0261   desc.add<edm::ParameterSetDescription>("MTEC", empty);
0262 }
0263 
0264 edm::ParameterSet SeedingLayerSetsBuilder::layerConfig(const std::string& nameLayer,
0265                                                        const edm::ParameterSet& cfg) const {
0266   edm::ParameterSet result;
0267 
0268   for (string::size_type iEnd = nameLayer.size(); iEnd > 0; --iEnd) {
0269     string name = nameLayer.substr(0, iEnd);
0270     if (cfg.exists(name))
0271       return cfg.getParameter<edm::ParameterSet>(name);
0272   }
0273   edm::LogError("SeedingLayerSetsBuilder")
0274       << "configuration for layer: " << nameLayer << " not found, job will probably crash!";
0275   return result;
0276 }
0277 
0278 vector<vector<string> > SeedingLayerSetsBuilder::layerNamesInSets(const vector<string>& namesPSet) {
0279   std::vector<std::vector<std::string> > result;
0280   for (std::vector<std::string>::const_iterator is = namesPSet.begin(); is < namesPSet.end(); ++is) {
0281     vector<std::string> layersInSet;
0282     string line = *is;
0283     string::size_type pos = 0;
0284     while (pos != string::npos) {
0285       pos = line.find('+');
0286       string layer = line.substr(0, pos);
0287       layersInSet.push_back(layer);
0288       line = line.substr(pos + 1, string::npos);
0289     }
0290     result.push_back(layersInSet);
0291   }
0292   return result;
0293 }
0294 
0295 void SeedingLayerSetsBuilder::updateEventSetup(const edm::EventSetup& es) {
0296   // We want to evaluate both in the first invocation (to properly
0297   // initialize ESWatcher), and this way we avoid one branch compared
0298   // to || (should be tiny effect)
0299   if (!(geometryWatcher_.check(es) || trhWatcher_.check(es)))
0300     return;
0301 
0302   const GeometricSearchTracker& tracker = es.getData(trackerToken_);
0303 
0304   const std::vector<BarrelDetLayer const*>& bpx = tracker.barrelLayers();
0305   const std::vector<BarrelDetLayer const*>& tib = tracker.tibLayers();
0306   const std::vector<BarrelDetLayer const*>& tob = tracker.tobLayers();
0307 
0308   const std::vector<ForwardDetLayer const*>& fpx_pos = tracker.posForwardLayers();
0309   const std::vector<ForwardDetLayer const*>& tid_pos = tracker.posTidLayers();
0310   const std::vector<ForwardDetLayer const*>& tec_pos = tracker.posTecLayers();
0311 
0312   const std::vector<ForwardDetLayer const*>& fpx_neg = tracker.negForwardLayers();
0313   const std::vector<ForwardDetLayer const*>& tid_neg = tracker.negTidLayers();
0314   const std::vector<ForwardDetLayer const*>& tec_neg = tracker.negTecLayers();
0315 
0316   for (const auto& layer : theLayers) {
0317     const DetLayer* detLayer = nullptr;
0318     int index = layer.idLayer - 1;
0319 
0320     if (layer.subdet == GeomDetEnumerators::PixelBarrel) {
0321       detLayer = bpx[index];
0322     } else if (layer.subdet == GeomDetEnumerators::PixelEndcap) {
0323       if (layer.side == TrackerDetSide::PosEndcap) {
0324         detLayer = fpx_pos[index];
0325       } else {
0326         detLayer = fpx_neg[index];
0327       }
0328     } else if (layer.subdet == GeomDetEnumerators::TIB) {
0329       detLayer = tib[index];
0330     } else if (layer.subdet == GeomDetEnumerators::TID) {
0331       if (layer.side == TrackerDetSide::PosEndcap) {
0332         detLayer = tid_pos[index];
0333       } else {
0334         detLayer = tid_neg[index];
0335       }
0336     } else if (layer.subdet == GeomDetEnumerators::TOB) {
0337       detLayer = tob[index];
0338     } else if (layer.subdet == GeomDetEnumerators::TEC) {
0339       if (layer.side == TrackerDetSide::PosEndcap) {
0340         detLayer = tec_pos[index];
0341       } else {
0342         detLayer = tec_neg[index];
0343       }
0344     } else {
0345       throw cms::Exception("Configuration") << "Did not find DetLayer for layer " << theLayerNames[layer.nameIndex];
0346     }
0347 
0348     theLayerDets[layer.nameIndex] = detLayer;
0349     theTTRHBuilders[layer.nameIndex] = &es.getData(layer.hitBuilderToken);
0350   }
0351 }
0352 
0353 std::vector<SeedingLayerSetsBuilder::SeedingLayerId> SeedingLayerSetsBuilder::layers() const {
0354   std::vector<SeedingLayerId> ret;
0355   ret.reserve(numberOfLayers());
0356   for (const auto& layer : theLayers) {
0357     ret.emplace_back(layer.subdet, layer.side, layer.idLayer);
0358   }
0359   return ret;
0360 }
0361 
0362 std::unique_ptr<SeedingLayerSetsHits> SeedingLayerSetsBuilder::hits(const edm::Event& ev, const edm::EventSetup& es) {
0363   updateEventSetup(es);
0364 
0365   auto ret = std::make_unique<SeedingLayerSetsHits>(
0366       theNumberOfLayersInSet, &theLayerSetIndices, &theLayerNames, &theLayerDets);
0367 
0368   for (auto& layer : theLayers) {
0369     ret->addHits(
0370         layer.nameIndex,
0371         layer.extractor->hits((const TkTransientTrackingRecHitBuilder&)(*theTTRHBuilders[layer.nameIndex]), ev, es));
0372   }
0373   ret->shrink_to_fit();
0374   return ret;
0375 }
0376 //new function for FastSim only
0377 std::unique_ptr<SeedingLayerSetsHits> SeedingLayerSetsBuilder::makeSeedingLayerSetsHitsforFastSim(
0378     const edm::Event& ev, const edm::EventSetup& es) {
0379   updateEventSetup(es);
0380 
0381   edm::Handle<FastTrackerRecHitCollection> fastSimrechits_;
0382   ev.getByToken(fastSimrecHitsToken_, fastSimrechits_);  //using FastSim RecHits
0383   const TrackerTopology* const tTopo = &es.getData(trackerTopologyToken_);
0384   SeedingLayerSetsHits::OwnedHits layerhits_;
0385 
0386   auto ret = std::make_unique<SeedingLayerSetsHits>(
0387       theNumberOfLayersInSet, &theLayerSetIndices, &theLayerNames, &theLayerDets);
0388 
0389   for (auto& layer : theLayers) {
0390     layerhits_.clear();
0391     for (auto& rh : *fastSimrechits_) {
0392       GeomDetEnumerators::SubDetector subdet = GeomDetEnumerators::invalidDet;
0393       TrackerDetSide side = TrackerDetSide::Barrel;
0394       int idLayer = 0;
0395       if ((rh.det()->geographicalId()).subdetId() == PixelSubdetector::PixelBarrel) {
0396         subdet = GeomDetEnumerators::PixelBarrel;
0397         side = TrackerDetSide::Barrel;
0398         idLayer = tTopo->pxbLayer(rh.det()->geographicalId());
0399       } else if ((rh.det()->geographicalId()).subdetId() == PixelSubdetector::PixelEndcap) {
0400         subdet = GeomDetEnumerators::PixelEndcap;
0401         idLayer = tTopo->pxfDisk(rh.det()->geographicalId());
0402         if (tTopo->pxfSide(rh.det()->geographicalId()) == 1)
0403           side = TrackerDetSide::NegEndcap;
0404         else
0405           side = TrackerDetSide::PosEndcap;
0406       }
0407 
0408       if (layer.subdet == subdet && layer.side == side && layer.idLayer == idLayer) {
0409         BaseTrackerRecHit const& b(rh);
0410         auto ptrHit = (BaseTrackerRecHit*)(b.clone());
0411         layerhits_.emplace_back(ptrHit);
0412       } else
0413         continue;
0414     }
0415     ret->addHits(layer.nameIndex, std::move(layerhits_));
0416   }
0417   ret->shrink_to_fit();
0418   return ret;
0419 }
[ 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!