PixelLowPtUtilities/src/PixelTripletLowPtGenerator.cc

0001 #include "RecoTracker/PixelLowPtUtilities/interface/PixelTripletLowPtGenerator.h"
0002 #include "RecoTracker/PixelLowPtUtilities/interface/ThirdHitPrediction.h"
0003 #include "RecoTracker/PixelLowPtUtilities/interface/TripletFilter.h"
0004 #include "RecoTracker/PixelLowPtUtilities/interface/HitInfo.h"
0005
0006 #include "RecoTracker/TkMSParametrization/interface/PixelRecoPointRZ.h"
0007 #include "RecoTracker/TkHitPairs/interface/HitPairGeneratorFromLayerPair.h"
0008
0009 #include "FWCore/Framework/interface/ConsumesCollector.h"
0010
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "DataFormats/SiPixelCluster/interface/SiPixelClusterShapeCache.h"
0013
0014 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0015 #include "RecoTracker/TkMSParametrization/interface/MultipleScatteringParametrisationMaker.h"
0016 #include "TrackingTools/Records/interface/TransientRecHitRecord.h"
0017 #include "RecoTracker/Record/interface/TrackerMultipleScatteringRecord.h"
0018 #include "RecoTracker/Record/interface/CkfComponentsRecord.h"
0019
0020 #undef Debug
0021
0022 using namespace std;
0023
0024 /*****************************************************************************/
0025 PixelTripletLowPtGenerator::PixelTripletLowPtGenerator(const edm::ParameterSet& cfg, edm::ConsumesCollector& iC)
0026     : HitTripletGeneratorFromPairAndLayers(),  // no theMaxElement used in this class
0027       m_geomToken(iC.esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>()),
0028       m_topoToken(iC.esConsumes<TrackerTopology, TrackerTopologyRcd>()),
0029       m_magfieldToken(iC.esConsumes()),
0030       m_ttrhBuilderToken(iC.esConsumes(edm::ESInputTag("", cfg.getParameter<string>("TTRHBuilder")))),
0031       m_msmakerToken(iC.esConsumes()),
0032       m_clusterFilterToken(iC.esConsumes(edm::ESInputTag("", "ClusterShapeHitFilter"))),
0033       theTracker(nullptr),
0034       theClusterShapeCacheToken(
0035           iC.consumes<SiPixelClusterShapeCache>(cfg.getParameter<edm::InputTag>("clusterShapeCacheSrc"))) {
0036   checkMultipleScattering = cfg.getParameter<bool>("checkMultipleScattering");
0037   nSigMultipleScattering = cfg.getParameter<double>("nSigMultipleScattering");
0038   checkClusterShape = cfg.getParameter<bool>("checkClusterShape");
0039   rzTolerance = cfg.getParameter<double>("rzTolerance");
0040   maxAngleRatio = cfg.getParameter<double>("maxAngleRatio");
0041 }
0042
0043 /*****************************************************************************/
0044 PixelTripletLowPtGenerator::~PixelTripletLowPtGenerator() {}
0045
0046 /*****************************************************************************/
0047
0048 /*****************************************************************************/
0049 void PixelTripletLowPtGenerator::getTracker(const edm::EventSetup& es) {
0050   if (theTracker == nullptr) {
0051     // Get tracker geometry
0052     theTracker = &es.getData(m_geomToken);
0053   }
0054
0055   if (!theFilter) {
0056     theFilter = std::make_unique<TripletFilter>(&es.getData(m_clusterFilterToken));
0057   }
0058 }
0059
0060 /*****************************************************************************/
0061 GlobalPoint PixelTripletLowPtGenerator::getGlobalPosition(const TrackingRecHit* recHit) {
0062   DetId detId = recHit->geographicalId();
0063
0064   return theTracker->idToDet(detId)->toGlobal(recHit->localPosition());
0065 }
0066
0067 /*****************************************************************************/
0068 void PixelTripletLowPtGenerator::hitTriplets(const TrackingRegion& region,
0069                                              OrderedHitTriplets& result,
0070                                              const edm::Event& ev,
0071                                              const edm::EventSetup& es,
0072                                              const SeedingLayerSetsHits::SeedingLayerSet& pairLayers,
0073                                              const std::vector<SeedingLayerSetsHits::SeedingLayer>& thirdLayers) {
0074   //Retrieve tracker topology from geometry
0075   const TrackerTopology* tTopo = &es.getData(m_topoToken);
0076   const auto& magfield = es.getData(m_magfieldToken);
0077   const auto& ttrhBuilder = es.getData(m_ttrhBuilderToken);
0078   const auto& msmaker = es.getData(m_msmakerToken);
0079
0080   edm::Handle<SiPixelClusterShapeCache> clusterShapeCache;
0081   ev.getByToken(theClusterShapeCacheToken, clusterShapeCache);
0082
0083   // Generate pairs
0084   OrderedHitPairs pairs;
0085   pairs.reserve(30000);
0086   thePairGenerator->hitPairs(region, pairs, ev, es, pairLayers);
0087
0088   if (pairs.empty())
0089     return;
0090
0091   int size = thirdLayers.size();
0092
0093   // Set aliases
0094   const RecHitsSortedInPhi** thirdHitMap = new const RecHitsSortedInPhi*[size];
0095   for (int il = 0; il < size; il++)
0096     thirdHitMap[il] = &(*theLayerCache)(thirdLayers[il], region);
0097
0098   // Get tracker
0099   getTracker(es);
0100
0101   // Look at all generated pairs
0102   for (OrderedHitPairs::const_iterator ip = pairs.begin(); ip != pairs.end(); ip++) {
0103     // Fill rechits and points
0104     vector<const TrackingRecHit*> recHits(3);
0105     vector<GlobalPoint> points(3);
0106
0107     recHits[0] = (*ip).inner()->hit();
0108     recHits[1] = (*ip).outer()->hit();
0109
0110 #ifdef Debug
0111     cerr << " RecHits " + HitInfo::getInfo(*recHits[0]) + HitInfo::getInfo(*recHits[1]) << endl;
0112 #endif
0113
0114     for (int i = 0; i < 2; i++)
0115       points[i] = getGlobalPosition(recHits[i]);
0116
0117     // Initialize helix prediction
0118     ThirdHitPrediction thePrediction(
0119         region, points[0], points[1], magfield, ttrhBuilder, nSigMultipleScattering, maxAngleRatio);
0120
0121     // Look at all layers
0122     for (int il = 0; il < size; il++) {
0123       const DetLayer* layer = thirdLayers[il].detLayer();
0124
0125 #ifdef Debug
0126       cerr << "  check layer " << layer->subDetector() << " " << layer->location() << endl;
0127 #endif
0128
0129       // Get ranges for the third hit
0130       float phi[2], rz[2];
0131       thePrediction.getRanges(layer, phi, rz);
0132
0133       PixelRecoRange<float> phiRange(phi[0], phi[1]);
0134       PixelRecoRange<float> rzRange(rz[0] - rzTolerance, rz[1] + rzTolerance);
0135
0136       // Get third hit candidates from cache
0137       typedef RecHitsSortedInPhi::Hit Hit;
0138       vector<Hit> thirdHits = thirdHitMap[il]->hits(phiRange.min(), phiRange.max());
0139       typedef vector<Hit>::const_iterator IH;
0140
0141       for (IH th = thirdHits.begin(), eh = thirdHits.end(); th < eh; ++th) {
0142         // Fill rechit and point
0143         recHits[2] = (*th)->hit();
0144         points[2] = getGlobalPosition(recHits[2]);
0145
0146 #ifdef Debug
0147         cerr << "  third hit " + HitInfo::getInfo(*recHits[2]) << endl;
0148 #endif
0149
0150         // Check if third hit is compatible with multiple scattering
0151         vector<GlobalVector> globalDirs;
0152         if (thePrediction.isCompatibleWithMultipleScattering(points[2], recHits, globalDirs, msmaker) == false) {
0153 #ifdef Debug
0154           cerr << "  not compatible: multiple scattering" << endl;
0155 #endif
0156           if (checkMultipleScattering)
0157             continue;
0158         }
0159
0160         // Convert to localDirs
0161         /*
0162         vector<LocalVector> localDirs;
0163         vector<GlobalVector>::const_iterator globalDir = globalDirs.begin();
0164         for(vector<const TrackingRecHit *>::const_iterator
0165                                             recHit  = recHits.begin();
0166                                             recHit != recHits.end(); recHit++)
0167         {
0168           localDirs.push_back(theTracker->idToDet(
0169                              (*recHit)->geographicalId())->toLocal(*globalDir));
0170           globalDir++;
0171         }
0172 */
0173
0174         // Check if the cluster shapes are compatible with thrusts
0175         if (checkClusterShape) {
0176           if (!theFilter->checkTrack(recHits, globalDirs, tTopo, *clusterShapeCache)) {
0177 #ifdef Debug
0178             cerr << "  not compatible: cluster shape" << endl;
0179 #endif
0180             continue;
0181           }
0182         }
0183
0184         // All checks passed, put triplet back
0185         result.push_back(OrderedHitTriplet((*ip).inner(), (*ip).outer(), *th));
0186       }
0187     }
0188   }
0189   delete[] thirdHitMap;
0190
0191   return;
0192 }
0193 void PixelTripletLowPtGenerator::hitTriplets(const TrackingRegion& region,
0194                                              OrderedHitTriplets& result,
0195                                              const edm::EventSetup& es,
0196                                              const HitDoublets& doublets,
0197                                              const RecHitsSortedInPhi** thirdHitMap,
0198                                              const std::vector<const DetLayer*>& thirdLayerDetLayer,
0199                                              const int nThirdLayers) {
0200   throw cms::Exception("Error") << "PixelTripletLowPtGenerator::hitTriplets is not implemented \n";
0201 }