CMS 3D CMS Logo

List of all members | Public Member Functions | Private Attributes
SETPatternRecognition Class Reference

#include <SETPatternRecognition.h>

Inheritance diagram for SETPatternRecognition:
MuonSeedVPatternRecognition

Public Member Functions

void produce (const edm::Event &event, const edm::EventSetup &eSetup, std::vector< MuonRecHitContainer > &result) override
 
bool segmentCleaning (const DetId &detId, const LocalPoint &localPosition, const LocalError &localError, const LocalVector &localDirection, const LocalError &localDirectionError, const double &chi2, const int &ndf)
 
 SETPatternRecognition (const edm::ParameterSet &pset, edm::ConsumesCollector &iC)
 
void setServiceProxy (MuonServiceProxy *service)
 
 ~SETPatternRecognition () override
 
- Public Member Functions inherited from MuonSeedVPatternRecognition
 MuonSeedVPatternRecognition (const edm::ParameterSet &pset)
 
virtual ~MuonSeedVPatternRecognition ()
 

Private Attributes

edm::InputTag CSCRecSegmentLabel
 
edm::EDGetTokenT< CSCSegmentCollectioncscToken
 
edm::InputTag DTRecSegmentLabel
 
edm::EDGetTokenT< DTRecSegment4DCollectiondtToken
 
int maxActiveChambers
 
double minLocalSegmentAngle
 
double outsideChamberErrorScale
 
edm::InputTag RPCRecSegmentLabel
 
edm::EDGetTokenT< RPCRecHitCollectionrpcToken
 
MuonServiceProxytheService
 
bool useRPCs
 

Additional Inherited Members

- Public Types inherited from MuonSeedVPatternRecognition
typedef MuonTransientTrackingRecHit::ConstMuonRecHitPointer ConstMuonRecHitPointer
 
typedef MuonTransientTrackingRecHit::MuonRecHitContainer MuonRecHitContainer
 
typedef MuonTransientTrackingRecHit::MuonRecHitPointer MuonRecHitPointer
 
- Protected Attributes inherited from MuonSeedVPatternRecognition
bool enableCSCMeasurement
 Enable the CSC measurement. More...
 
bool enableDTMeasurement
 Enable the DT measurement. More...
 
bool enableME0Measurement
 Enable the ME0 measurement. More...
 
edm::InputTag theCSCRecSegmentLabel
 the name of the CSC rec hits collection More...
 
edm::InputTag theDTRecSegmentLabel
 the name of the DT rec hits collection More...
 
edm::InputTag theME0RecSegmentLabel
 the name of the ME0 rec hits collection More...
 

Detailed Description

I. Bloch, E. James, S. Stoynev

Definition at line 13 of file SETPatternRecognition.h.

Constructor & Destructor Documentation

◆ SETPatternRecognition()

SETPatternRecognition::SETPatternRecognition ( const edm::ParameterSet pset,
edm::ConsumesCollector iC 
)
explicit

Definition at line 17 of file SETPatternRecognition.cc.

References edm::ConsumesCollector::consumes(), CSCRecSegmentLabel, cscToken, DTRecSegmentLabel, dtToken, edm::ParameterSet::getParameter(), maxActiveChambers, metname, minLocalSegmentAngle, outsideChamberErrorScale, edm::parameterSet(), RPCRecSegmentLabel, rpcToken, and useRPCs.

19  .getParameter<ParameterSet>("FilterParameters")) {
20  const string metname = "Muon|RecoMuon|SETPatternRecognition";
21  // Parameter set for the Builder
22  ParameterSet trajectoryBuilderParameters = parameterSet.getParameter<ParameterSet>("SETTrajBuilderParameters");
23  // The inward-outward fitter (starts from seed state)
24  ParameterSet filterPSet = trajectoryBuilderParameters.getParameter<ParameterSet>("FilterParameters");
25  maxActiveChambers = filterPSet.getParameter<int>("maxActiveChambers");
26  useRPCs = filterPSet.getParameter<bool>("EnableRPCMeasurement");
27  DTRecSegmentLabel = filterPSet.getParameter<edm::InputTag>("DTRecSegmentLabel");
28  CSCRecSegmentLabel = filterPSet.getParameter<edm::InputTag>("CSCRecSegmentLabel");
29  RPCRecSegmentLabel = filterPSet.getParameter<edm::InputTag>("RPCRecSegmentLabel");
30 
31  outsideChamberErrorScale = filterPSet.getParameter<double>("OutsideChamberErrorScale");
32  minLocalSegmentAngle = filterPSet.getParameter<double>("MinLocalSegmentAngle");
33  //----
37 }
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
edm::EDGetTokenT< DTRecSegment4DCollection > dtToken
edm::EDGetTokenT< RPCRecHitCollection > rpcToken
edm::EDGetTokenT< CSCSegmentCollection > cscToken
const std::string metname
ParameterSet const & parameterSet(StableProvenance const &provenance, ProcessHistory const &history)
Definition: Provenance.cc:11
MuonSeedVPatternRecognition(const edm::ParameterSet &pset)

◆ ~SETPatternRecognition()

SETPatternRecognition::~SETPatternRecognition ( )
inlineoverride

Definition at line 16 of file SETPatternRecognition.h.

16 {}

Member Function Documentation

◆ produce()

void SETPatternRecognition::produce ( const edm::Event event,
const edm::EventSetup eSetup,
std::vector< MuonRecHitContainer > &  result 
)
overridevirtual

Output is a cluster, with possibly more than one hit per layer

Implements MuonSeedVPatternRecognition.

Definition at line 40 of file SETPatternRecognition.cc.

References SplitLinear::begin, isoTrack_cff::chi2, dtChamberEfficiency_cfi::cscSegments, cscToken, dtToken, HLT_2024v14_cff::dXclusBoxMax, HLT_2024v14_cff::dYclusBoxMax, mps_fire::end, edm::eventsetup::heterocontainer::insert(), ALPAKA_ACCELERATOR_NAMESPACE::vertexFinder::it, maxActiveChambers, metname, Pi, dtTriggerPhase2PrimitiveDigis_cfi::rpcRecHits, rpcToken, HLT_2024v14_cff::seeds, segmentCleaning(), findQualityFiles::size, MuonTransientTrackingRecHit::specificBuild(), groupFilesInBlocks::temp, theService, and useRPCs.

Referenced by SETMuonSeedProducer::produce().

42  {
43  const string metname = "Muon|RecoMuon|SETMuonSeedSeed";
44 
45  //---- Build collection of all segments
46  MuonRecHitContainer muonRecHits;
47  MuonRecHitContainer muonRecHits_DT2D_hasPhi;
48  MuonRecHitContainer muonRecHits_DT2D_hasZed;
49  MuonRecHitContainer muonRecHits_RPC;
50 
51  // ********************************************;
52  // Get the DT-Segment collection from the Event
53  // ********************************************;
54 
56  event.getByToken(dtToken, dtRecHits);
57  std::vector<DTChamberId> chambers_DT;
58  std::vector<DTChamberId>::const_iterator chIt_DT;
59  for (DTRecSegment4DCollection::const_iterator rechit = dtRecHits->begin(); rechit != dtRecHits->end(); ++rechit) {
60  bool insert = true;
61  for (chIt_DT = chambers_DT.begin(); chIt_DT != chambers_DT.end(); ++chIt_DT) {
62  if (((*rechit).chamberId().wheel()) == ((*chIt_DT).wheel()) &&
63  ((*rechit).chamberId().station() == (*chIt_DT).station()) &&
64  ((*rechit).chamberId().sector() == (*chIt_DT).sector())) {
65  insert = false;
66  }
67  }
68  if (insert) {
69  chambers_DT.push_back((*rechit).chamberId());
70  }
71  if (segmentCleaning((*rechit).geographicalId(),
72  rechit->localPosition(),
73  rechit->localPositionError(),
74  rechit->localDirection(),
75  rechit->localDirectionError(),
76  rechit->chi2(),
77  rechit->degreesOfFreedom())) {
78  continue;
79  }
80  if ((rechit->hasZed() && rechit->hasPhi())) {
81  muonRecHits.push_back(MuonTransientTrackingRecHit::specificBuild(
82  theService->trackingGeometry()->idToDet((*rechit).geographicalId()), &*rechit));
83  } else if (rechit->hasZed()) {
84  muonRecHits_DT2D_hasZed.push_back(MuonTransientTrackingRecHit::specificBuild(
85  theService->trackingGeometry()->idToDet((*rechit).geographicalId()), &*rechit));
86  } else if (rechit->hasPhi()) { // safeguard
87  muonRecHits_DT2D_hasPhi.push_back(MuonTransientTrackingRecHit::specificBuild(
88  theService->trackingGeometry()->idToDet((*rechit).geographicalId()), &*rechit));
89  } else {
90  //std::cout<<"Warning in "<<metname<<": DT segment which claims to have neither phi nor Z."<<std::endl;
91  }
92  }
93  //std::cout<<"DT done"<<std::endl;
94 
95  // ********************************************;
96  // Get the CSC-Segment collection from the event
97  // ********************************************;
98 
100  event.getByToken(cscToken, cscSegments);
101  std::vector<CSCDetId> chambers_CSC;
102  std::vector<CSCDetId>::const_iterator chIt_CSC;
103  for (CSCSegmentCollection::const_iterator rechit = cscSegments->begin(); rechit != cscSegments->end(); ++rechit) {
104  bool insert = true;
105  for (chIt_CSC = chambers_CSC.begin(); chIt_CSC != chambers_CSC.end(); ++chIt_CSC) {
106  if (((*rechit).cscDetId().chamber() == (*chIt_CSC).chamber()) &&
107  ((*rechit).cscDetId().station() == (*chIt_CSC).station()) &&
108  ((*rechit).cscDetId().ring() == (*chIt_CSC).ring()) &&
109  ((*rechit).cscDetId().endcap() == (*chIt_CSC).endcap())) {
110  insert = false;
111  }
112  }
113  if (insert) {
114  chambers_CSC.push_back((*rechit).cscDetId().chamberId());
115  }
116  if (segmentCleaning((*rechit).geographicalId(),
117  rechit->localPosition(),
118  rechit->localPositionError(),
119  rechit->localDirection(),
120  rechit->localDirectionError(),
121  rechit->chi2(),
122  rechit->degreesOfFreedom())) {
123  continue;
124  }
125  muonRecHits.push_back(MuonTransientTrackingRecHit::specificBuild(
126  theService->trackingGeometry()->idToDet((*rechit).geographicalId()), &*rechit));
127  }
128  //std::cout<<"CSC done"<<std::endl;
129 
130  // ********************************************;
131  // Get the RPC-Hit collection from the event
132  // ********************************************;
133 
135  event.getByToken(rpcToken, rpcRecHits);
136  if (useRPCs) {
137  for (RPCRecHitCollection::const_iterator rechit = rpcRecHits->begin(); rechit != rpcRecHits->end(); ++rechit) {
138  // RPCs are special
139  const LocalVector localDirection(0., 0., 1.);
140  const LocalError localDirectionError(0., 0., 0.);
141  const double chi2 = 1.;
142  const int ndf = 1;
143  if (segmentCleaning((*rechit).geographicalId(),
144  rechit->localPosition(),
145  rechit->localPositionError(),
146  localDirection,
147  localDirectionError,
148  chi2,
149  ndf)) {
150  continue;
151  }
152  muonRecHits_RPC.push_back(MuonTransientTrackingRecHit::specificBuild(
153  theService->trackingGeometry()->idToDet((*rechit).geographicalId()), &*rechit));
154  }
155  }
156  //std::cout<<"RPC done"<<std::endl;
157  //
158  if (int(chambers_DT.size() + chambers_CSC.size()) > maxActiveChambers) {
159  // std::cout <<" Too many active chambers : nDT = "<<chambers_DT.size()<<
160  // " nCSC = "<<chambers_CSC.size()<<" Skip them all."<<std::endl;
161  edm::LogWarning("tooManyActiveChambers") << " Too many active chambers : nDT = " << chambers_DT.size()
162  << " nCSC = " << chambers_CSC.size() << " Skip them all.";
163  muonRecHits.clear();
164  muonRecHits_DT2D_hasPhi.clear();
165  muonRecHits_DT2D_hasZed.clear();
166  muonRecHits_RPC.clear();
167  }
168  //---- Find "pre-clusters" from all segments; these contain potential muon candidates
169 
170  //---- From all the hits (i.e. segments; sometimes "rechits" is also used with the same meaning;
171  //---- this convention has meaning in the global reconstruction though could be misleading
172  //---- from a local reconstruction point of view; "local rechits" are used in the backward fit only)
173  //---- make clusters of hits; a trajectory could contain hits from one cluster only
174 
175  // the clustering procedure is very similar to the one used in the segment reconstruction
176 
177  bool useDT2D_hasPhi = true;
178  bool useDT2D_hasZed = true;
179  double dXclusBoxMax = 0.60; // phi - can be as large as 15 - 20 degrees for 6 GeV muons
180  double dYclusBoxMax = 0.;
181 
182  // this is the main selection criteria; the value of 0.02 rad seems wide enough to
183  // contain any hit from a passing muon and still narrow enough to remove good part of
184  // possible "junk" hits
185  // (Comment: it may be better to allow maximum difference between any two hits in a trajectory
186  // to be 0.02 or 0.04 or ...; currently the requirement below is imposed on two consecutive hits)
187 
188  dYclusBoxMax = 0.02; // theta // hardoded - remove it!
189 
190  // X and Y are distance variables - we use eta and phi here
191 
192  float dXclus = 0.0;
193  float dXclus_box = 0.0;
194  float dYclus_box = 0.0;
195 
197 
198  std::vector<MuonRecHitContainer> seeds;
199 
200  std::vector<float> running_meanX;
201  std::vector<float> running_meanY;
202 
203  std::vector<float> seed_minX;
204  std::vector<float> seed_maxX;
205  std::vector<float> seed_minY;
206  std::vector<float> seed_maxY;
207 
208  // split rechits into subvectors and return vector of vectors:
209  // Loop over rechits
210  // Create one seed per hit
211  for (MuonRecHitContainer::const_iterator it = muonRecHits.begin(); it != muonRecHits.end(); ++it) {
212  // try to avoid using 2D DT segments. We will add them later to the
213  // clusters they are most likely to belong to. Might need to add them
214  // to more than just one cluster, if we find them to be consistent with
215  // more than one. This would lead to an implicit sharing of hits between
216  // SA muon candidates.
217 
218  temp.clear();
219 
220  temp.push_back((*it));
221 
222  seeds.push_back(temp);
223 
224  // First added hit in seed defines the mean to which the next hit is compared
225  // for this seed.
226 
227  running_meanX.push_back((*it)->globalPosition().phi());
228  running_meanY.push_back((*it)->globalPosition().theta());
229 
230  // set min/max X and Y for box containing the hits in the precluster:
231  seed_minX.push_back((*it)->globalPosition().phi());
232  seed_maxX.push_back((*it)->globalPosition().phi());
233  seed_minY.push_back((*it)->globalPosition().theta());
234  seed_maxY.push_back((*it)->globalPosition().theta());
235  }
236 
237  // merge clusters that are too close
238  // measure distance between final "running mean"
239  for (unsigned int NNN = 0; NNN < seeds.size(); ++NNN) {
240  for (unsigned int MMM = NNN + 1; MMM < seeds.size(); ++MMM) {
241  if (running_meanX[MMM] == 999999. || running_meanX[NNN] == 999999.) {
242  // LogDebug("CSC") << "CSCSegmentST::clusterHits: Warning: Skipping used seeds, this should happen - inform developers!\n";
243  //std::cout<<"We should never see this line now!!!"<<std::endl;
244  continue; //skip seeds that have been used
245  }
246 
247  // Some complications for using phi as a clustering variable due to wrap-around (-pi = pi)
248  // Define temporary mean, min, and max variables for the cluster which could be merged (NNN)
249  double temp_meanX = running_meanX[NNN];
250  double temp_minX = seed_minX[NNN];
251  double temp_maxX = seed_maxX[NNN];
252 
253  // check if the difference between the two phi values is greater than pi
254  // if so, need to shift temporary values by 2*pi to make a valid comparison
255  dXclus = running_meanX[NNN] - running_meanX[MMM];
256  if (dXclus > TMath::Pi()) {
257  temp_meanX = temp_meanX - 2. * TMath::Pi();
258  temp_minX = temp_minX - 2. * TMath::Pi();
259  temp_maxX = temp_maxX - 2. * TMath::Pi();
260  }
261  if (dXclus < -TMath::Pi()) {
262  temp_meanX = temp_meanX + 2. * TMath::Pi();
263  temp_minX = temp_minX + 2. * TMath::Pi();
264  temp_maxX = temp_maxX + 2. * TMath::Pi();
265  }
266 
267  // // calculate cut criteria for simple running mean distance cut:
268  // // not sure that these values are really used anywhere
269 
270  // calculate minmal distance between precluster boxes containing the hits:
271  // use the temp variables from above for phi of the NNN cluster
272  if (temp_meanX > running_meanX[MMM])
273  dXclus_box = temp_minX - seed_maxX[MMM];
274  else
275  dXclus_box = seed_minX[MMM] - temp_maxX;
276  if (running_meanY[NNN] > running_meanY[MMM])
277  dYclus_box = seed_minY[NNN] - seed_maxY[MMM];
278  else
279  dYclus_box = seed_minY[MMM] - seed_maxY[NNN];
280 
281  if (dXclus_box < dXclusBoxMax && dYclus_box < dYclusBoxMax) {
282  // merge clusters!
283  // merge by adding seed NNN to seed MMM and erasing seed NNN
284 
285  // calculate running mean for the merged seed:
286  // use the temp variables from above for phi of the NNN cluster
287  running_meanX[MMM] = (temp_meanX * seeds[NNN].size() + running_meanX[MMM] * seeds[MMM].size()) /
288  (seeds[NNN].size() + seeds[MMM].size());
289  running_meanY[MMM] = (running_meanY[NNN] * seeds[NNN].size() + running_meanY[MMM] * seeds[MMM].size()) /
290  (seeds[NNN].size() + seeds[MMM].size());
291 
292  // update min/max X and Y for box containing the hits in the merged cluster:
293  // use the temp variables from above for phi of the NNN cluster
294  if (temp_minX <= seed_minX[MMM])
295  seed_minX[MMM] = temp_minX;
296  if (temp_maxX > seed_maxX[MMM])
297  seed_maxX[MMM] = temp_maxX;
298  if (seed_minY[NNN] <= seed_minY[MMM])
299  seed_minY[MMM] = seed_minY[NNN];
300  if (seed_maxY[NNN] > seed_maxY[MMM])
301  seed_maxY[MMM] = seed_maxY[NNN];
302 
303  // now check to see if the running mean has moved outside of the allowed -pi to pi region
304  // if so, then adjust shift all values up or down by 2 * pi
305  if (running_meanX[MMM] > TMath::Pi()) {
306  running_meanX[MMM] = running_meanX[MMM] - 2. * TMath::Pi();
307  seed_minX[MMM] = seed_minX[MMM] - 2. * TMath::Pi();
308  seed_maxX[MMM] = seed_maxX[MMM] - 2. * TMath::Pi();
309  }
310  if (running_meanX[MMM] < -TMath::Pi()) {
311  running_meanX[MMM] = running_meanX[MMM] + 2. * TMath::Pi();
312  seed_minX[MMM] = seed_minX[MMM] + 2. * TMath::Pi();
313  seed_maxX[MMM] = seed_maxX[MMM] + 2. * TMath::Pi();
314  }
315 
316  // add seed NNN to MMM (lower to larger number)
317  seeds[MMM].insert(seeds[MMM].end(), seeds[NNN].begin(), seeds[NNN].end());
318 
319  // mark seed NNN as used (at the moment just set running mean to 999999.)
320  running_meanX[NNN] = 999999.;
321  running_meanY[NNN] = 999999.;
322  // we have merged a seed (NNN) to the highter seed (MMM) - need to contimue to
323  // next seed (NNN+1)
324  break;
325  }
326  }
327  }
328  bool tooCloseClusters = false;
329  if (seeds.size() > 1) {
330  std::vector<double> seedTheta(seeds.size());
331  for (unsigned int iSeed = 0; iSeed < seeds.size(); ++iSeed) {
332  seedTheta[iSeed] = seeds[iSeed][0]->globalPosition().theta();
333  if (iSeed) {
334  double dTheta = fabs(seedTheta[iSeed] - seedTheta[iSeed - 1]);
335  if (dTheta < 0.5) { //? should be something more clever
336  tooCloseClusters = true;
337  break;
338  }
339  }
340  }
341  }
342 
343  // have formed clusters from all hits except for 2D DT segments. Now add the 2D segments to the
344  // compatible clusters. For this we compare the mean cluster postition with the
345  // 2D segment position. We should use the valid coordinate only and use the bad coordinate
346  // as a cross check.
347  for (unsigned int NNN = 0; NNN < seeds.size(); ++NNN) {
348  if (running_meanX[NNN] == 999999.)
349  continue; //skip seeds that have been marked as used up in merging
350 
351  // We have a valid cluster - loop over all 2D segments.
352  if (useDT2D_hasZed) {
353  for (MuonRecHitContainer::const_iterator it2 = muonRecHits_DT2D_hasZed.begin();
354  it2 != muonRecHits_DT2D_hasZed.end();
355  ++it2) {
356  // check that global theta of 2-D segment lies within cluster box plus or minus allowed slop
357  if (((*it2)->globalPosition().theta() < seed_maxY[NNN] + dYclusBoxMax) &&
358  ((*it2)->globalPosition().theta() > seed_minY[NNN] - dYclusBoxMax)) {
359  // check that global phi of 2-D segment (assumed to be center of chamber since no phi hit info)
360  // matches with cluster box plus or minus allowed slop given that the true phi value could be
361  // anywhere within a given chamber (+/- 5 degrees ~ 0.09 radians from center)
362  if (!((((*it2)->globalPosition().phi() + 0.09) < (seed_minX[NNN] - dXclusBoxMax) &&
363  ((*it2)->globalPosition().phi() - 0.09) < (seed_minX[NNN] - dXclusBoxMax)) ||
364  (((*it2)->globalPosition().phi() + 0.09) > (seed_maxX[NNN] + dXclusBoxMax) &&
365  // we have checked that the 2Dsegment is within tight theta boundaries and loose phi boundaries of the current cluster -> add it
366  ((*it2)->globalPosition().phi() - 0.09) > (seed_maxX[NNN] + dXclusBoxMax)))) {
367  seeds[NNN].push_back((*it2));
368  }
369  }
370  }
371  }
372 
373  // put DT hasphi loop here
374  if (useDT2D_hasPhi) {
375  for (MuonRecHitContainer::const_iterator it2 = muonRecHits_DT2D_hasPhi.begin();
376  it2 != muonRecHits_DT2D_hasPhi.end();
377  ++it2) {
378  if (((*it2)->globalPosition().phi() < seed_maxX[NNN] + dXclusBoxMax) &&
379  ((*it2)->globalPosition().phi() > seed_minX[NNN] - dXclusBoxMax)) {
380  if (!((((*it2)->globalPosition().theta() + 0.3) < (seed_minY[NNN] - dYclusBoxMax) &&
381  ((*it2)->globalPosition().theta() - 0.3) < (seed_minY[NNN] - dYclusBoxMax)) ||
382  (((*it2)->globalPosition().theta() + 0.3) > (seed_maxY[NNN] + dYclusBoxMax) &&
383  // we have checked that the 2Dsegment is within tight phi boundaries and loose theta boundaries of the current cluster -> add it
384  ((*it2)->globalPosition().theta() - 0.3) > (seed_maxY[NNN] + dYclusBoxMax)))) {
385  seeds[NNN].push_back((*it2)); // warning - neeed eta/theta switch here
386  }
387  }
388  }
389  } // DT2D_hastPhi loop
390 
391  // put RPC loop here
392  int secondCh = 0;
393  DetId detId_prev;
394  if (seeds[NNN].size() > 1) { // actually we should check how many chambers with measurements are present
395  for (unsigned int iRH = 0; iRH < seeds[NNN].size(); ++iRH) {
396  if (iRH && detId_prev != seeds[NNN][iRH]->hit()->geographicalId()) {
397  ++secondCh;
398  break;
399  }
400  detId_prev = seeds[NNN][iRH]->hit()->geographicalId();
401  }
402  }
403 
404  if (useRPCs && !secondCh && !tooCloseClusters) {
405  for (MuonRecHitContainer::const_iterator it2 = muonRecHits_RPC.begin(); it2 != muonRecHits_RPC.end(); ++it2) {
406  if (((*it2)->globalPosition().phi() < seed_maxX[NNN] + dXclusBoxMax) &&
407  ((*it2)->globalPosition().phi() > seed_minX[NNN] - dXclusBoxMax)) {
408  if (!((((*it2)->globalPosition().theta() + 0.3) < (seed_minY[NNN] - dYclusBoxMax) &&
409  ((*it2)->globalPosition().theta() - 0.3) < (seed_minY[NNN] - dYclusBoxMax)) ||
410  (((*it2)->globalPosition().theta() + 0.3) > (seed_maxY[NNN] + dYclusBoxMax) &&
411  // we have checked that the 2Dsegment is within tight phi boundaries and loose theta boundaries of the current cluster -> add it
412  ((*it2)->globalPosition().theta() - 0.3) > (seed_maxY[NNN] + dYclusBoxMax)))) {
413  seeds[NNN].push_back((*it2)); // warning - neeed eta/theta switch here
414  }
415  }
416  }
417  } // RPC loop
418  }
419 
420  // hand over the final seeds to the output
421  // would be more elegant if we could do the above step with
422  // erasing the merged ones, rather than the
423  for (unsigned int NNN = 0; NNN < seeds.size(); ++NNN) {
424  if (running_meanX[NNN] == 999999.)
425  continue; //skip seeds that have been marked as used up in merging
426  //std::cout<<"Next Cluster..."<<std::endl;
427  segments_clusters.push_back(seeds[NNN]);
428  }
429 }
size
Write out results.
const double Pi
edm::EDGetTokenT< DTRecSegment4DCollection > dtToken
edm::EDGetTokenT< RPCRecHitCollection > rpcToken
edm::EDGetTokenT< CSCSegmentCollection > cscToken
const std::string metname
C::const_iterator const_iterator
constant access iterator type
Definition: RangeMap.h:43
bool segmentCleaning(const DetId &detId, const LocalPoint &localPosition, const LocalError &localError, const LocalVector &localDirection, const LocalError &localDirectionError, const double &chi2, const int &ndf)
bool insert(Storage &iStorage, ItemType *iItem, const IdTag &iIdTag)
Definition: HCMethods.h:50
MuonServiceProxy * theService
Definition: DetId.h:17
MuonTransientTrackingRecHit::MuonRecHitContainer MuonRecHitContainer
Log< level::Warning, false > LogWarning
static MuonRecHitPointer specificBuild(const GeomDet *geom, const TrackingRecHit *rh)

◆ segmentCleaning()

bool SETPatternRecognition::segmentCleaning ( const DetId detId,
const LocalPoint localPosition,
const LocalError localError,
const LocalVector localDirection,
const LocalError localDirectionError,
const double &  chi2,
const int &  ndf 
)

Definition at line 431 of file SETPatternRecognition.cc.

References Surface::bounds(), hcalRecHitTable_cff::detId, Bounds::inside(), minLocalSegmentAngle, outsideChamberErrorScale, GeomDet::surface(), theService, and PV3DBase< T, PVType, FrameType >::z().

Referenced by produce().

437  {
438  // drop segments which are "bad"
439  bool dropTheSegment = true;
440  const GeomDet* geomDet = theService->trackingGeometry()->idToDet(detId);
441  // only segments whithin the boundaries of the chamber
442  bool insideCh = geomDet->surface().bounds().inside(localPosition, localError, outsideChamberErrorScale);
443 
444  // Don't use segments (nearly) parallel to the chamberi;
445  // the direction vector is normalized (R=1)
446  bool parallelSegment = fabs(localDirection.z()) > minLocalSegmentAngle ? false : true;
447 
448  if (insideCh && !parallelSegment) {
449  dropTheSegment = false;
450  }
451  // use chi2 too? (DT, CSCs, RPCs; 2D, 4D;...)
452 
453  return dropTheSegment;
454 }
virtual bool inside(const Local3DPoint &) const =0
Determine if the point is inside the bounds.
T z() const
Definition: PV3DBase.h:61
MuonServiceProxy * theService
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
const Bounds & bounds() const
Definition: Surface.h:87

◆ setServiceProxy()

void SETPatternRecognition::setServiceProxy ( MuonServiceProxy service)
inline

Definition at line 22 of file SETPatternRecognition.h.

References theService.

Referenced by SETMuonSeedProducer::SETMuonSeedProducer().

22 { theService = service; }
MuonServiceProxy * theService

Member Data Documentation

◆ CSCRecSegmentLabel

edm::InputTag SETPatternRecognition::CSCRecSegmentLabel
private

Definition at line 37 of file SETPatternRecognition.h.

Referenced by SETPatternRecognition().

◆ cscToken

edm::EDGetTokenT<CSCSegmentCollection> SETPatternRecognition::cscToken
private

Definition at line 41 of file SETPatternRecognition.h.

Referenced by produce(), and SETPatternRecognition().

◆ DTRecSegmentLabel

edm::InputTag SETPatternRecognition::DTRecSegmentLabel
private

Definition at line 36 of file SETPatternRecognition.h.

Referenced by SETPatternRecognition().

◆ dtToken

edm::EDGetTokenT<DTRecSegment4DCollection> SETPatternRecognition::dtToken
private

Definition at line 40 of file SETPatternRecognition.h.

Referenced by produce(), and SETPatternRecognition().

◆ maxActiveChambers

int SETPatternRecognition::maxActiveChambers
private

Definition at line 33 of file SETPatternRecognition.h.

Referenced by produce(), and SETPatternRecognition().

◆ minLocalSegmentAngle

double SETPatternRecognition::minLocalSegmentAngle
private

Definition at line 45 of file SETPatternRecognition.h.

Referenced by segmentCleaning(), and SETPatternRecognition().

◆ outsideChamberErrorScale

double SETPatternRecognition::outsideChamberErrorScale
private

Definition at line 44 of file SETPatternRecognition.h.

Referenced by segmentCleaning(), and SETPatternRecognition().

◆ RPCRecSegmentLabel

edm::InputTag SETPatternRecognition::RPCRecSegmentLabel
private

Definition at line 38 of file SETPatternRecognition.h.

Referenced by SETPatternRecognition().

◆ rpcToken

edm::EDGetTokenT<RPCRecHitCollection> SETPatternRecognition::rpcToken
private

Definition at line 42 of file SETPatternRecognition.h.

Referenced by produce(), and SETPatternRecognition().

◆ theService

MuonServiceProxy* SETPatternRecognition::theService
private

Definition at line 48 of file SETPatternRecognition.h.

Referenced by produce(), segmentCleaning(), and setServiceProxy().

◆ useRPCs

bool SETPatternRecognition::useRPCs
private

Definition at line 34 of file SETPatternRecognition.h.

Referenced by produce(), and SETPatternRecognition().