SpecialSeedGenerators/src/BeamHaloPairGenerator.cc

0001 #include "RecoTracker/SpecialSeedGenerators/interface/BeamHaloPairGenerator.h"
0002 typedef SeedingHitSet::ConstRecHitPointer SeedingHit;
0003
0004 #include "FWCore/Framework/interface/ConsumesCollector.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "TrackingTools/TransientTrackingRecHit/interface/SeedingLayerSetsHits.h"
0007
0008 BeamHaloPairGenerator::BeamHaloPairGenerator(const edm::ParameterSet& conf, edm::ConsumesCollector& iC)
0009     : theSeedingLayerToken(iC.consumes<SeedingLayerSetsHits>(conf.getParameter<edm::InputTag>("LayerSrc"))) {
0010   edm::LogInfo("CtfSpecialSeedGenerator|BeamHaloPairGenerator") << "Constructing BeamHaloPairGenerator";
0011   theMaxTheta = conf.getParameter<double>("maxTheta");
0012   theMaxTheta = fabs(sin(theMaxTheta));
0013 }
0014
0015 const OrderedSeedingHits& BeamHaloPairGenerator::run(const TrackingRegion& region,
0016                                                      const edm::Event& e,
0017                                                      const edm::EventSetup& es) {
0018   hitPairs.clear();
0019   edm::Handle<SeedingLayerSetsHits> hlayers;
0020   e.getByToken(theSeedingLayerToken, hlayers);
0021   const SeedingLayerSetsHits& layers = *hlayers;
0022   if (layers.numberOfLayersInSet() != 2)
0023     throw cms::Exception("CtfSpecialSeedGenerator")
0024         << "You are using " << layers.numberOfLayersInSet() << " layers in set instead of 2 ";
0025   for (SeedingLayerSetsHits::SeedingLayerSet ls : layers) {
0026     auto innerHits = region.hits(ls[0]);
0027     auto outerHits = region.hits(ls[1]);
0028
0029     for (auto iOuterHit = outerHits.begin(); iOuterHit != outerHits.end(); iOuterHit++) {
0030       for (auto iInnerHit = innerHits.begin(); iInnerHit != innerHits.end(); iInnerHit++) {
0031         //do something in there... if necessary
0032         SeedingHitSet::ConstRecHitPointer crhpi = &(**iInnerHit);
0033         SeedingHitSet::ConstRecHitPointer crhpo = &(**iOuterHit);
0034         GlobalVector d = crhpo->globalPosition() - crhpi->globalPosition();
0035         double ABSsinDtheta = fabs(sin(d.theta()));
0036         LogDebug("BeamHaloPairGenerator")
0037             << "position1: " << crhpo->globalPosition() << " position2: " << crhpi->globalPosition()
0038             << " |sin(Dtheta)|: " << ABSsinDtheta << ((ABSsinDtheta > theMaxTheta) ? " skip" : " keep");
0039
0040         if (ABSsinDtheta > theMaxTheta) {
0041           continue;
0042         }
0043
0044         hitPairs.push_back(OrderedHitPair(crhpi, crhpo));
0045       }
0046     }
0047   }
0048   return hitPairs;
0049 }