CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
SeedFinderSelector.cc
Go to the documentation of this file.
2 
3 // framework
6 
7 // track reco
21 // data formats
23 
25  : trackingRegion_(nullptr)
26  , eventSetup_(nullptr)
27  , measurementTracker_(nullptr)
28  , measurementTrackerLabel_(cfg.getParameter<std::string>("measurementTracker"))
29 {
30  if(cfg.exists("pixelTripletGeneratorFactory"))
31  {
32  const edm::ParameterSet & tripletConfig = cfg.getParameter<edm::ParameterSet>("pixelTripletGeneratorFactory");
33  pixelTripletGenerator_.reset(HitTripletGeneratorFromPairAndLayersFactory::get()->create(tripletConfig.getParameter<std::string>("ComponentName"),tripletConfig,consumesCollector));
34  }
35 
36  if(cfg.exists("MultiHitGeneratorFactory"))
37  {
38  const edm::ParameterSet & tripletConfig = cfg.getParameter<edm::ParameterSet>("MultiHitGeneratorFactory");
39  multiHitGenerator_.reset(MultiHitGeneratorFromPairAndLayersFactory::get()->create(tripletConfig.getParameter<std::string>("ComponentName"),tripletConfig));
40  }
41 
42  if(cfg.exists("CAHitTripletGeneratorFactory"))
43  {
44  const edm::ParameterSet & tripletConfig = cfg.getParameter<edm::ParameterSet>("CAHitTripletGeneratorFactory");
45  CAHitTriplGenerator_ = std::make_unique<CAHitTripletGenerator>(tripletConfig,consumesCollector);
46  seedingLayers_ = std::make_unique<SeedingLayerSetsBuilder>(cfg, consumesCollector, edm::InputTag("fastTrackerRecHits")); //calling the new FastSim specific constructor to make SeedingLayerSetsHits pointer for triplet iterations
47  layerPairs_ = cfg.getParameter<std::vector<unsigned>>("layerPairs"); //allowed layer pairs for CA triplets
48  }
49 
50  if(cfg.exists("CAHitQuadrupletGeneratorFactory"))
51  {
52  const edm::ParameterSet & quadrupletConfig = cfg.getParameter<edm::ParameterSet>("CAHitQuadrupletGeneratorFactory");
53  CAHitQuadGenerator_ = std::make_unique<CAHitQuadrupletGenerator>(quadrupletConfig, consumesCollector);
54  seedingLayers_ = std::make_unique<SeedingLayerSetsBuilder>(cfg, consumesCollector, edm::InputTag("fastTrackerRecHits")); //calling the new FastSim specific constructor to make SeedingLayerSetsHits pointer for quadruplet iterations
55  layerPairs_ = cfg.getParameter<std::vector<unsigned>>("layerPairs"); //allowed layer pairs for CA quadruplets
56  }
57 
59  {
60  throw cms::Exception("FastSimTracking") << "It is forbidden to specify together 'pixelTripletGeneratorFactory', 'CAHitTripletGeneratorFactory' and 'MultiHitGeneratorFactory' in configuration of SeedFinderSelection";
61  }
63  {
64  throw cms::Exception("FastSimTracking") << "It is forbidden to specify 'CAHitQuadrupletGeneratorFactory' together with 'pixelTripletGeneratorFactory', 'CAHitTripletGeneratorFactory' or 'MultiHitGeneratorFactory' in configuration of SeedFinderSelection";
65  }
66 }
67 
68 
70 
72 {
73  eventSetup_ = &es;
74 
75  edm::ESHandle<MeasurementTracker> measurementTrackerHandle;
76  es.get<CkfComponentsRecord>().get(measurementTrackerLabel_, measurementTrackerHandle);
78  measurementTracker_ = &(*measurementTrackerHandle);
79 
81  {
82  multiHitGenerator_->initES(es);
83  }
84 
85  //for CA triplet iterations
87  seedingLayer = seedingLayers_->makeSeedingLayerSetsHitsforFastSim(ev, es);
88  CAHitTriplGenerator_->initEvent(ev,es);
89  }
90  //for CA quadruplet iterations
92  seedingLayer = seedingLayers_->makeSeedingLayerSetsHitsforFastSim(ev, es);
93  CAHitQuadGenerator_->initEvent(ev,es);
94  }
95 }
96 
97 
98 bool SeedFinderSelector::pass(const std::vector<const FastTrackerRecHit *>& hits) const
99 {
101  {
102  throw cms::Exception("FastSimTracking") << "ERROR: event not initialized";
103  }
104  if(!trackingRegion_)
105  {
106  throw cms::Exception("FastSimTracking") << "ERROR: trackingRegion not set";
107  }
108 
109 
110  // check the inner 2 hits
111  if(hits.size() < 2)
112  {
113  throw cms::Exception("FastSimTracking") << "SeedFinderSelector::pass requires at least 2 hits";
114  }
115  const DetLayer * firstLayer = measurementTracker_->geometricSearchTracker()->detLayer(hits[0]->det()->geographicalId());
116  const DetLayer * secondLayer = measurementTracker_->geometricSearchTracker()->detLayer(hits[1]->det()->geographicalId());
117 
118  std::vector<BaseTrackerRecHit const *> firstHits{hits[0]};
119  std::vector<BaseTrackerRecHit const *> secondHits{hits[1]};
120 
121  const RecHitsSortedInPhi fhm(firstHits, trackingRegion_->origin(), firstLayer);
122  const RecHitsSortedInPhi shm(secondHits, trackingRegion_->origin(), secondLayer);
123 
124  HitDoublets result(fhm,shm);
125  HitPairGeneratorFromLayerPair::doublets(*trackingRegion_,*firstLayer,*secondLayer,fhm,shm,*eventSetup_,0,result);
126 
127  if(result.empty())
128  {
129  return false;
130  }
131 
132  // check the inner 3 hits
134  {
135  if(hits.size() < 3)
136  {
137  throw cms::Exception("FastSimTracking") << "For the given configuration, SeedFinderSelector::pass requires at least 3 hits";
138  }
139  const DetLayer * thirdLayer = measurementTracker_->geometricSearchTracker()->detLayer(hits[2]->det()->geographicalId());
140  std::vector<const DetLayer *> thirdLayerDetLayer(1,thirdLayer);
141  std::vector<BaseTrackerRecHit const *> thirdHits{hits[2]};
142  const RecHitsSortedInPhi thm(thirdHits,trackingRegion_->origin(), thirdLayer);
143  const RecHitsSortedInPhi * thmp =&thm;
144 
146  {
147  OrderedHitTriplets tripletresult;
148  pixelTripletGenerator_->hitTriplets(*trackingRegion_,tripletresult,*eventSetup_,result,&thmp,thirdLayerDetLayer,1);
149  return !tripletresult.empty();
150  }
151  else if(multiHitGenerator_)
152  {
153  OrderedMultiHits tripletresult;
154  multiHitGenerator_->hitTriplets(*trackingRegion_,tripletresult,*eventSetup_,result,&thmp,thirdLayerDetLayer,1);
155  return !tripletresult.empty();
156  }
157  //new for Phase1
158  else if(CAHitTriplGenerator_)
159  {
160  if(!seedingLayer)
161  throw cms::Exception("FastSimTracking") << "ERROR: SeedingLayers pointer not set for CATripletGenerator";
162 
164  //constructing IntermediateHitDoublets to be passed onto CAHitTripletGenerator::hitNtuplets()
165  IntermediateHitDoublets ihd(&layers);
166  const TrackingRegion& tr_ = *trackingRegion_;
167  auto filler = ihd.beginRegion(&tr_);
168 
169  //Forming doublets for CA triplets from the allowed layer pair combinations:(0,1),(1,2)
170  std::string hitPair[2];
171  for(int i=0; i<2; i++){
173  hitPair[0] = Layer_name(hits[i]);
174  hitPair[1] = Layer_name(hits[i+1]);
175 
176  bool found;
178  found = false;
179  for(const auto p : layerPairs_){
180  pairCandidate = ls.slice(p,p+2);
181  if(hitPair[0] == pairCandidate[0].name() && hitPair[1] == pairCandidate[1].name()){
182  found = true;
183  break;
184  }
185  }
186  if(found)
187  break;
188  }
189  assert(found == true);
190  const DetLayer * fLayer = measurementTracker_->geometricSearchTracker()->detLayer(hits[i]->det()->geographicalId());
191  const DetLayer * sLayer = measurementTracker_->geometricSearchTracker()->detLayer(hits[i+1]->det()->geographicalId());
192  std::vector<BaseTrackerRecHit const *> fHits{hits[i]};
193  std::vector<BaseTrackerRecHit const *> sHits{hits[i+1]};
194 
195  //Important: doublets to be added to the cache
196  auto& layerCache = filler.layerHitMapCache();
197  const RecHitsSortedInPhi& firsthm = *layerCache.add(pairCandidate[0], std::make_unique<RecHitsSortedInPhi>(fHits, trackingRegion_->origin(),fLayer));
198  const RecHitsSortedInPhi& secondhm = *layerCache.add(pairCandidate[1], std::make_unique<RecHitsSortedInPhi>(sHits, trackingRegion_->origin(),sLayer));
199  HitDoublets res(firsthm,secondhm);
200  HitPairGeneratorFromLayerPair::doublets(*trackingRegion_,*fLayer,*sLayer,firsthm,secondhm,*eventSetup_,0,res);
201  filler.addDoublets(pairCandidate, std::move(res)); //fill the formed doublet
202  }
203  std::vector<OrderedHitSeeds> tripletresult;
204  tripletresult.resize(ihd.regionSize());
205  for(auto& ntuplet : tripletresult)
206  ntuplet.reserve(3);
207  CAHitTriplGenerator_->hitNtuplets(ihd,tripletresult,*eventSetup_,*seedingLayer); //calling the function from the class, modifies tripletresult
208  return !tripletresult.empty();
209  }
210  }
211  //new for Phase1
213  {
214  if(hits.size() < 4)
215  {
216  throw cms::Exception("FastSimTracking") << "For the given configuration, SeedFinderSelector::pass requires at least 4 hits";
217  }
218 
219  if(!seedingLayer)
220  throw cms::Exception("FastSimTracking") << "ERROR: SeedingLayers pointer not set for CAHitQuadrupletGenerator";
221 
223  //constructing IntermediateHitDoublets to be passed onto CAHitQuadrupletGenerator::hitNtuplets()
224  IntermediateHitDoublets ihd(&layers);
225  const TrackingRegion& tr_ = *trackingRegion_;
226  auto filler = ihd.beginRegion(&tr_);
227 
228  //Forming doublets for CA quadruplets from the allowed layer pair combinations:(0,1),(1,2),(2,3)
229  std::string hitPair[2];
230  for(int i=0; i<3; i++){
232  hitPair[0] = Layer_name(hits[i]);
233  hitPair[1] = Layer_name(hits[i+1]);
234 
235  bool found;
237  found = false;
238  for(const auto p : layerPairs_){
239  pairCandidate = ls.slice(p,p+2);
240  if(hitPair[0] == pairCandidate[0].name() && hitPair[1] == pairCandidate[1].name()){
241  found = true;
242  break;
243  }
244  }
245  if(found)
246  break;
247  }
248  assert(found == true);
249  const DetLayer * fLayer = measurementTracker_->geometricSearchTracker()->detLayer(hits[i]->det()->geographicalId());
250  const DetLayer * sLayer = measurementTracker_->geometricSearchTracker()->detLayer(hits[i+1]->det()->geographicalId());
251  std::vector<BaseTrackerRecHit const *> fHits{hits[i]};
252  std::vector<BaseTrackerRecHit const *> sHits{hits[i+1]};
253 
254  //Important: doublets to be added to the cache
255  auto& layerCache = filler.layerHitMapCache();
256  const RecHitsSortedInPhi& firsthm = *layerCache.add(pairCandidate[0], std::make_unique<RecHitsSortedInPhi>(fHits, trackingRegion_->origin(), fLayer));
257  const RecHitsSortedInPhi& secondhm = *layerCache.add(pairCandidate[1], std::make_unique<RecHitsSortedInPhi>(sHits, trackingRegion_->origin(), sLayer));
258  HitDoublets res(firsthm,secondhm);
259  HitPairGeneratorFromLayerPair::doublets(*trackingRegion_,*fLayer,*sLayer,firsthm,secondhm,*eventSetup_,0,res);
260  filler.addDoublets(pairCandidate, std::move(res)); //fill the formed doublet
261  }
262 
263  std::vector<OrderedHitSeeds> quadrupletresult;
264  quadrupletresult.resize(ihd.regionSize());
265  for(auto& ntuplet : quadrupletresult)
266  ntuplet.reserve(4);
267  CAHitQuadGenerator_->hitNtuplets(ihd,quadrupletresult,*eventSetup_,*seedingLayer); //calling the function from the class, modifies quadrupletresult
268  return !quadrupletresult.empty();
269  }
270 
271  return true;
272 
273 }
274 
275 //new for Phase1
277 {
278  std::string name="UNKWN";
279  const TrackerTopology* const tTopo = trackerTopology.product();
280  int idLayer = 0;
281  if( (hit->det()->geographicalId()).subdetId() == PixelSubdetector::PixelBarrel){
282  idLayer = tTopo->pxbLayer(hit->det()->geographicalId());
283  name = "BPix"+std::to_string(idLayer);
284  }
285  else if ((hit->det()->geographicalId()).subdetId() == PixelSubdetector::PixelEndcap){
286  idLayer = tTopo->pxfDisk(hit->det()->geographicalId());
287  if(tTopo->pxfSide(hit->det()->geographicalId())==1)
288  name = "FPix"+std::to_string(idLayer)+"_neg";
289  else
290  name = "FPix"+std::to_string(idLayer)+"_pos";
291  }
292  return name;
293 }
T getParameter(std::string const &) const
bool pass(const std::vector< const FastTrackerRecHit * > &hits) const
std::vector< LayerSetAndLayers > layers(const SeedingLayerSetsHits &sets)
Definition: LayerTriplets.cc:4
std::unique_ptr< CAHitQuadrupletGenerator > CAHitQuadGenerator_
const edm::EventSetup * eventSetup_
GlobalPoint const & origin() const
def create(alignables, pedeDump, additionalData, outputFile, config)
std::unique_ptr< SeedingLayerSetsHits > seedingLayer
unsigned int pxfDisk(const DetId &id) const
bool exists(std::string const &parameterName) const
checks if a parameter exists
bool ev
SeedFinderSelector(const edm::ParameterSet &cfg, edm::ConsumesCollector &&consumesCollector)
RegionFiller beginRegion(const TrackingRegion *region)
#define nullptr
std::unique_ptr< HitTripletGeneratorFromPairAndLayers > pixelTripletGenerator_
Definition: Electron.h:4
std::string Layer_name(const FastTrackerRecHit *hit) const
const GeomDet * det() const
std::vector< unsigned > layerPairs_
const TrackingRegion * trackingRegion_
const DetLayer * detLayer(const DetId &id) const
obsolete method. Use idToLayer() instead.
DetId geographicalId() const
The label of this GeomDet.
Definition: GeomDet.h:79
HitDoublets doublets(const TrackingRegion &reg, const edm::Event &ev, const edm::EventSetup &es, Layers layers)
edm::ESHandle< TrackerTopology > trackerTopology
unsigned int pxbLayer(const DetId &id) const
def ls(path, rec=False)
Definition: eostools.py:348
const T & get() const
Definition: EventSetup.h:55
void initEvent(const edm::Event &ev, const edm::EventSetup &es)
std::unique_ptr< SeedingLayerSetsBuilder > seedingLayers_
const MeasurementTracker * measurementTracker_
unsigned int pxfSide(const DetId &id) const
const GeometricSearchTracker * geometricSearchTracker() const
const std::string measurementTrackerLabel_
T const * product() const
Definition: ESHandle.h:86
std::unique_ptr< CAHitTripletGenerator > CAHitTriplGenerator_
def move(src, dest)
Definition: eostools.py:510
T get(const Candidate &c)
Definition: component.h:55
std::unique_ptr< MultiHitGeneratorFromPairAndLayers > multiHitGenerator_