CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
CosmicRegionalSeedGenerator.cc
Go to the documentation of this file.
1 #include <iostream>
2 #include <memory>
3 #include <string>
4 
8 
10 
11 #include <TMath.h>
12 
13 using namespace std;
14 using namespace trigger;
15 using namespace reco;
16 using namespace edm;
17 
19  : magfieldToken_(iC.esConsumes()) {
20  edm::LogInfo("CosmicRegionalSeedGenerator") << "Begin Run:: Constructing CosmicRegionalSeedGenerator";
21 
22  edm::ParameterSet regionPSet = conf.getParameter<edm::ParameterSet>("RegionPSet");
23  ptMin_ = regionPSet.getParameter<double>("ptMin");
24  rVertex_ = regionPSet.getParameter<double>("rVertex");
25  zVertex_ = regionPSet.getParameter<double>("zVertex");
26  deltaEta_ = regionPSet.getParameter<double>("deltaEtaRegion");
27  deltaPhi_ = regionPSet.getParameter<double>("deltaPhiRegion");
28 
29  edm::ParameterSet toolsPSet = conf.getParameter<edm::ParameterSet>("ToolsPSet");
30  propagatorToken_ = iC.esConsumes(edm::ESInputTag("", toolsPSet.getParameter<std::string>("thePropagatorName")));
31  regionBase_ = toolsPSet.getParameter<std::string>("regionBase");
32 
33  edm::ParameterSet collectionsPSet = conf.getParameter<edm::ParameterSet>("CollectionsPSet");
34  recoMuonsCollection_ = collectionsPSet.getParameter<edm::InputTag>("recoMuonsCollection");
35  recoTrackMuonsCollection_ = collectionsPSet.getParameter<edm::InputTag>("recoTrackMuonsCollection");
36  recoL2MuonsCollection_ = collectionsPSet.getParameter<edm::InputTag>("recoL2MuonsCollection");
37 
38  edm::ParameterSet regionInJetsCheckPSet = conf.getParameter<edm::ParameterSet>("RegionInJetsCheckPSet");
39  doJetsExclusionCheck_ = regionInJetsCheckPSet.getParameter<bool>("doJetsExclusionCheck");
40  deltaRExclusionSize_ = regionInJetsCheckPSet.getParameter<double>("deltaRExclusionSize");
41  jetsPtMin_ = regionInJetsCheckPSet.getParameter<double>("jetsPtMin");
42  recoCaloJetsCollection_ = regionInJetsCheckPSet.getParameter<edm::InputTag>("recoCaloJetsCollection");
47  measurementTrackerEventToken_ = iC.consumes<MeasurementTrackerEvent>(edm::InputTag("MeasurementTrackerEvent"));
48 
49  edm::LogInfo("CosmicRegionalSeedGenerator") << "Reco muons collection: " << recoMuonsCollection_ << "\n"
50  << "Reco tracks muons collection: " << recoTrackMuonsCollection_ << "\n"
51  << "Reco L2 muons collection: " << recoL2MuonsCollection_;
52 }
53 
54 std::vector<std::unique_ptr<TrackingRegion>> CosmicRegionalSeedGenerator::regions(const edm::Event& event,
55  const edm::EventSetup& es) const {
56  std::vector<std::unique_ptr<TrackingRegion>> result;
57 
61  event.getByToken(measurementTrackerEventToken_, hmte);
62  measurementTracker = hmte.product();
63  }
64 
65  //get the propagator
66  auto const& propagator = es.getData(propagatorToken_);
67  auto const& magfield = es.getData(magfieldToken_);
68 
69  //________________________________________
70  //
71  //Seeding on Sta muon (MC && Datas)
72  //________________________________________
73 
74  if (regionBase_ == "seedOnStaMuon" || regionBase_.empty()) {
75  LogDebug("CosmicRegionalSeedGenerator") << "Seeding on stand alone muons ";
76 
77  //get collections
78  //+++++++++++++++
79 
80  //get the muon collection
82  event.getByToken(recoMuonsToken_, muonsHandle);
83  if (!muonsHandle.isValid()) {
84  edm::LogError("CollectionNotFound") << "Error::No reco muons collection (" << recoMuonsCollection_
85  << ") in the event - Please verify the name of the muon collection";
86  return result;
87  }
88 
89  LogDebug("CosmicRegionalSeedGenerator") << "Muons collection size = " << muonsHandle->size();
90 
91  //get the jet collection
92  edm::Handle<CaloJetCollection> caloJetsHandle;
93  event.getByToken(recoCaloJetsToken_, caloJetsHandle);
94 
95  //definition of the region
96  //+++++++++++++++++++++++++
97 
98  int nmuons = 0;
99  for (reco::MuonCollection::const_iterator staMuon = muonsHandle->begin(); staMuon != muonsHandle->end();
100  ++staMuon) {
101  //select sta muons
102  if (!staMuon->isStandAloneMuon()) {
103  LogDebug("CosmicRegionalSeedGenerator") << "This muon is not a stand alone muon";
104  continue;
105  }
106 
107  //bit 25 as a coverage -1.4 < eta < 1.4
108  if (abs(staMuon->standAloneMuon()->eta()) > 1.5)
109  continue;
110 
111  //debug
112  nmuons++;
113  LogDebug("CosmicRegionalSeedGenerator") << "Muon stand alone found in the collection - in muons chambers: \n "
114  << "Position = " << staMuon->standAloneMuon()->outerPosition() << "\n "
115  << "Momentum = " << staMuon->standAloneMuon()->outerMomentum() << "\n "
116  << "Eta = " << staMuon->standAloneMuon()->eta() << "\n "
117  << "Phi = " << staMuon->standAloneMuon()->phi();
118 
119  //initial position, momentum, charge
120 
121  GlobalPoint initialRegionPosition(staMuon->standAloneMuon()->referencePoint().x(),
122  staMuon->standAloneMuon()->referencePoint().y(),
123  staMuon->standAloneMuon()->referencePoint().z());
124  GlobalVector initialRegionMomentum(staMuon->standAloneMuon()->momentum().x(),
125  staMuon->standAloneMuon()->momentum().y(),
126  staMuon->standAloneMuon()->momentum().z());
127  int charge = (int)staMuon->standAloneMuon()->charge();
128 
129  LogDebug("CosmicRegionalSeedGenerator") << "Initial region - Reference point of the sta muon: \n "
130  << "Position = " << initialRegionPosition << "\n "
131  << "Momentum = " << initialRegionMomentum << "\n "
132  << "Eta = " << initialRegionPosition.eta() << "\n "
133  << "Phi = " << initialRegionPosition.phi() << "\n "
134  << "Charge = " << charge;
135 
136  //propagation on the last layers of TOB
137  if (staMuon->standAloneMuon()->outerPosition().y() > 0)
138  initialRegionMomentum *= -1;
139  GlobalTrajectoryParameters glb_parameters(
140  initialRegionPosition, initialRegionMomentum, charge, propagator.magneticField());
141  FreeTrajectoryState fts(glb_parameters);
142  StateOnTrackerBound onBounds(&propagator);
143  TrajectoryStateOnSurface outer = onBounds(fts);
144 
145  if (!outer.isValid()) {
146  //edm::LogError("FailedPropagation") << "Trajectory state on surface not valid" ;
147  LogDebug("CosmicRegionalSeedGenerator") << "Trajectory state on surface not valid";
148  continue;
149  }
150 
151  //final position & momentum
152  GlobalPoint regionPosition = outer.globalPosition();
153  GlobalVector regionMom = outer.globalMomentum();
154 
155  LogDebug("CosmicRegionalSeedGenerator")
156  << "Region after propagation: \n "
157  << "Position = " << outer.globalPosition() << "\n "
158  << "Momentum = " << outer.globalMomentum() << "\n "
159  << "R = " << regionPosition.perp() << " ---- z = " << regionPosition.z() << "\n "
160  << "Eta = " << outer.globalPosition().eta() << "\n "
161  << "Phi = " << outer.globalPosition().phi();
162 
163  //step back
164  double stepBack = 1;
165  GlobalPoint center = regionPosition + stepBack * regionMom.unit();
166  GlobalVector v = stepBack * regionMom.unit();
167  LogDebug("CosmicRegionalSeedGenerator") << "Step back vector = " << v << "\n";
168 
169  //exclude region built in jets
170  if (doJetsExclusionCheck_) {
171  double delta_R_min = 1000.;
172  for (CaloJetCollection::const_iterator jet = caloJetsHandle->begin(); jet != caloJetsHandle->end(); jet++) {
173  if (jet->pt() < jetsPtMin_)
174  continue;
175 
176  double deta = center.eta() - jet->eta();
177  double dphi = fabs(center.phi() - jet->phi());
178  if (dphi > TMath::Pi())
179  dphi = 2 * TMath::Pi() - dphi;
180 
181  double delta_R = sqrt(deta * deta + dphi * dphi);
182  if (delta_R < delta_R_min)
183  delta_R_min = delta_R;
184 
185  } //end loop on jets
186 
187  if (delta_R_min < deltaRExclusionSize_) {
188  LogDebug("CosmicRegionalSeedGenerator") << "Region built too close from a jet";
189  continue;
190  }
191  } //end if doJetsExclusionCheck
192 
193  //definition of the region
194 
195  result.push_back(std::make_unique<CosmicTrackingRegion>(
196  (-1) * regionMom, center, ptMin_, rVertex_, zVertex_, deltaEta_, deltaPhi_, magfield, measurementTracker));
197 
198  LogDebug("CosmicRegionalSeedGenerator")
199  << "Final CosmicTrackingRegion \n "
200  << "Position = " << center << "\n "
201  << "Direction = " << result.back()->direction() << "\n "
202  << "Distance from the region on the layer = " << (regionPosition - center).mag() << "\n "
203  << "Eta = " << center.eta() << "\n "
204  << "Phi = " << center.phi();
205 
206  } //end loop on muons
207 
208  } //end if SeedOnStaMuon
209 
210  //________________________________________
211  //
212  //Seeding on cosmic muons (MC && Datas)
213  //________________________________________
214 
215  if (regionBase_ == "seedOnCosmicMuon") {
216  LogDebug("CosmicRegionalSeedGenerator") << "Seeding on cosmic muons tracks";
217 
218  //get collections
219  //+++++++++++++++
220 
221  //get the muon collection
222  edm::Handle<reco::TrackCollection> cosmicMuonsHandle;
223  event.getByToken(recoTrackMuonsToken_, cosmicMuonsHandle);
224  if (!cosmicMuonsHandle.isValid()) {
225  edm::LogError("CollectionNotFound")
226  << "Error::No cosmic muons collection (" << recoTrackMuonsCollection_
227  << ") in the event - Please verify the name of the muon reco track collection";
228  return result;
229  }
230 
231  LogDebug("CosmicRegionalSeedGenerator") << "Cosmic muons tracks collection size = " << cosmicMuonsHandle->size();
232 
233  //get the jet collection
234  edm::Handle<CaloJetCollection> caloJetsHandle;
235  event.getByToken(recoCaloJetsToken_, caloJetsHandle);
236 
237  //definition of the region
238  //+++++++++++++++++++++++++
239 
240  int nmuons = 0;
241  for (reco::TrackCollection::const_iterator cosmicMuon = cosmicMuonsHandle->begin();
242  cosmicMuon != cosmicMuonsHandle->end();
243  ++cosmicMuon) {
244  //bit 25 as a coverage -1.4 < eta < 1.4
245  if (abs(cosmicMuon->eta()) > 1.5)
246  continue;
247 
248  nmuons++;
249 
250  //initial position, momentum, charge
251  GlobalPoint initialRegionPosition(
252  cosmicMuon->referencePoint().x(), cosmicMuon->referencePoint().y(), cosmicMuon->referencePoint().z());
253  GlobalVector initialRegionMomentum(
254  cosmicMuon->momentum().x(), cosmicMuon->momentum().y(), cosmicMuon->momentum().z());
255  int charge = (int)cosmicMuon->charge();
256 
257  LogDebug("CosmicRegionalSeedGenerator") << "Position and momentum of the muon track in the muon chambers: \n "
258  << "x = " << cosmicMuon->outerPosition().x() << "\n "
259  << "y = " << cosmicMuon->outerPosition().y() << "\n "
260  << "y = " << cosmicMuon->pt() << "\n "
261  << "Initial region - Reference point of the cosmic muon track: \n "
262  << "Position = " << initialRegionPosition << "\n "
263  << "Momentum = " << initialRegionMomentum << "\n "
264  << "Eta = " << initialRegionPosition.eta() << "\n "
265  << "Phi = " << initialRegionPosition.phi() << "\n "
266  << "Charge = " << charge;
267 
268  //propagation on the last layers of TOB
269  if (cosmicMuon->outerPosition().y() > 0 && cosmicMuon->momentum().y() < 0)
270  initialRegionMomentum *= -1;
271  GlobalTrajectoryParameters glb_parameters(
272  initialRegionPosition, initialRegionMomentum, charge, propagator.magneticField());
273  FreeTrajectoryState fts(glb_parameters);
274  StateOnTrackerBound onBounds(&propagator);
275  TrajectoryStateOnSurface outer = onBounds(fts);
276 
277  if (!outer.isValid()) {
278  //edm::LogError("FailedPropagation") << "Trajectory state on surface not valid" ;
279  LogDebug("CosmicRegionalSeedGenerator") << "Trajectory state on surface not valid";
280  continue;
281  }
282 
283  //final position & momentum
284  GlobalPoint regionPosition = outer.globalPosition();
285  GlobalVector regionMom = outer.globalMomentum();
286 
287  LogDebug("CosmicRegionalSeedGenerator")
288  << "Region after propagation: \n "
289  << "Position = " << outer.globalPosition() << "\n "
290  << "Momentum = " << outer.globalMomentum() << "\n "
291  << "R = " << regionPosition.perp() << " ---- z = " << regionPosition.z() << "\n "
292  << "Eta = " << outer.globalPosition().eta() << "\n "
293  << "Phi = " << outer.globalPosition().phi();
294 
295  //step back
296  double stepBack = 1;
297  GlobalPoint center = regionPosition + stepBack * regionMom.unit();
298  GlobalVector v = stepBack * regionMom.unit();
299  LogDebug("CosmicRegionalSeedGenerator") << "Step back vector = " << v << "\n";
300 
301  //exclude region built in jets
302  if (doJetsExclusionCheck_) {
303  double delta_R_min = 1000.;
304  for (CaloJetCollection::const_iterator jet = caloJetsHandle->begin(); jet != caloJetsHandle->end(); jet++) {
305  if (jet->pt() < jetsPtMin_)
306  continue;
307 
308  double deta = center.eta() - jet->eta();
309  double dphi = fabs(center.phi() - jet->phi());
310  if (dphi > TMath::Pi())
311  dphi = 2 * TMath::Pi() - dphi;
312 
313  double delta_R = sqrt(deta * deta + dphi * dphi);
314  if (delta_R < delta_R_min)
315  delta_R_min = delta_R;
316 
317  } //end loop on jets
318 
319  if (delta_R_min < deltaRExclusionSize_) {
320  LogDebug("CosmicRegionalSeedGenerator") << "Region built too close from a jet";
321  continue;
322  }
323  } // end if doJetsExclusionCheck
324 
325  //definition of the region
326  result.push_back(std::make_unique<CosmicTrackingRegion>(
327  (-1) * regionMom, center, ptMin_, rVertex_, zVertex_, deltaEta_, deltaPhi_, magfield, measurementTracker));
328 
329  LogDebug("CosmicRegionalSeedGenerator")
330  << "Final CosmicTrackingRegion \n "
331  << "Position = " << center << "\n "
332  << "Direction = " << result.back()->direction() << "\n "
333  << "Distance from the region on the layer = " << (regionPosition - center).mag() << "\n "
334  << "Eta = " << center.eta() << "\n "
335  << "Phi = " << center.phi();
336 
337  } //end loop on muons
338 
339  } //end if SeedOnCosmicMuon
340 
341  //________________________________________
342  //
343  //Seeding on L2 muons (MC && Datas)
344  //________________________________________
345 
346  if (regionBase_ == "seedOnL2Muon") {
347  LogDebug("CosmicRegionalSeedGenerator") << "Seeding on L2 muons";
348 
349  //get collections
350  //+++++++++++++++
351 
352  //get the muon collection
354  event.getByToken(recoL2MuonsToken_, L2MuonsHandle);
355 
356  if (!L2MuonsHandle.isValid()) {
357  edm::LogError("CollectionNotFound") << "Error::No L2 muons collection (" << recoL2MuonsCollection_
358  << ") in the event - Please verify the name of the L2 muon collection";
359  return result;
360  }
361 
362  LogDebug("CosmicRegionalSeedGenerator") << "L2 muons collection size = " << L2MuonsHandle->size();
363 
364  //definition of the region
365  //+++++++++++++++++++++++++
366 
367  int nmuons = 0;
368  for (reco::RecoChargedCandidateCollection::const_iterator L2Muon = L2MuonsHandle->begin();
369  L2Muon != L2MuonsHandle->end();
370  ++L2Muon) {
371  reco::TrackRef tkL2Muon = L2Muon->get<reco::TrackRef>();
372 
373  //bit 25 as a coverage -1.4 < eta < 1.4
374  if (abs(tkL2Muon->eta()) > 1.5)
375  continue;
376 
377  nmuons++;
378 
379  //initial position, momentum, charge
380  GlobalPoint initialRegionPosition(
381  tkL2Muon->referencePoint().x(), tkL2Muon->referencePoint().y(), tkL2Muon->referencePoint().z());
382  GlobalVector initialRegionMomentum(tkL2Muon->momentum().x(), tkL2Muon->momentum().y(), tkL2Muon->momentum().z());
383  int charge = (int)tkL2Muon->charge();
384 
385  LogDebug("CosmicRegionalSeedGenerator") << "Position and momentum of the L2 muon track in the muon chambers: \n "
386  << "x = " << tkL2Muon->outerPosition().x() << "\n "
387  << "y = " << tkL2Muon->outerPosition().y() << "\n "
388  << "y = " << tkL2Muon->pt() << "\n "
389  << "Initial region - Reference point of the L2 muon track: \n "
390  << "Position = " << initialRegionPosition << "\n "
391  << "Momentum = " << initialRegionMomentum << "\n "
392  << "Eta = " << initialRegionPosition.eta() << "\n "
393  << "Phi = " << initialRegionPosition.phi() << "\n "
394  << "Charge = " << charge;
395 
396  //seeding only in the bottom
397  if (tkL2Muon->outerPosition().y() > 0) {
398  LogDebug("CosmicRegionalSeedGenerator") << "L2 muon in the TOP --- Region not created";
399  return result;
400  }
401 
402  GlobalTrajectoryParameters glb_parameters(
403  initialRegionPosition, initialRegionMomentum, charge, propagator.magneticField());
404  FreeTrajectoryState fts(glb_parameters);
405  StateOnTrackerBound onBounds(&propagator);
406  TrajectoryStateOnSurface outer = onBounds(fts);
407 
408  if (!outer.isValid()) {
409  //edm::LogError("FailedPropagation") << "Trajectory state on surface not valid" ;
410  LogDebug("CosmicRegionalSeedGenerator") << "Trajectory state on surface not valid";
411  continue;
412  }
413 
414  //final position & momentum
415  GlobalPoint regionPosition = outer.globalPosition();
416  GlobalVector regionMom = outer.globalMomentum();
417 
418  LogDebug("CosmicRegionalSeedGenerator")
419  << "Region after propagation: \n "
420  << "Position = " << outer.globalPosition() << "\n "
421  << "Momentum = " << outer.globalMomentum() << "\n "
422  << "R = " << regionPosition.perp() << " ---- z = " << regionPosition.z() << "\n "
423  << "Eta = " << outer.globalPosition().eta() << "\n "
424  << "Phi = " << outer.globalPosition().phi();
425 
426  //step back
427  double stepBack = 1;
428  GlobalPoint center = regionPosition + stepBack * regionMom.unit();
429  GlobalVector v = stepBack * regionMom.unit();
430  LogDebug("CosmicRegionalSeedGenerator") << "Step back vector = " << v << "\n";
431 
432  //definition of the region
433  result.push_back(std::make_unique<CosmicTrackingRegion>(
434  (-1) * regionMom, center, ptMin_, rVertex_, zVertex_, deltaEta_, deltaPhi_, magfield, measurementTracker));
435 
436  LogDebug("CosmicRegionalSeedGenerator")
437  << "Final L2TrackingRegion \n "
438  << "Position = " << center << "\n "
439  << "Direction = " << result.back()->direction() << "\n "
440  << "Distance from the region on the layer = " << (regionPosition - center).mag() << "\n "
441  << "Eta = " << center.eta() << "\n "
442  << "Phi = " << center.phi();
443 
444  } //end loop on muons
445 
446  } //end if SeedOnL2Muon
447 
448  return result;
449 }
edm::EDGetTokenT< reco::TrackCollection > recoTrackMuonsToken_
const double Pi
tuple propagator
edm::EDGetTokenT< reco::CaloJetCollection > recoCaloJetsToken_
T perp() const
Definition: PV3DBase.h:69
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
tuple measurementTracker
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
constexpr bool isUninitialized() const noexcept
Definition: EDGetToken.h:99
std::vector< Track > TrackCollection
collection of Tracks
Definition: TrackFwd.h:14
GlobalPoint globalPosition() const
edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > magfieldToken_
Log< level::Error, false > LogError
std::vector< Muon > MuonCollection
collection of Muon objects
Definition: MuonFwd.h:9
tuple result
Definition: mps_fire.py:311
bool getData(T &iHolder) const
Definition: EventSetup.h:122
std::vector< std::unique_ptr< TrackingRegion > > regions(const edm::Event &event, const edm::EventSetup &es) const override
edm::ESGetToken< Propagator, TrackingComponentsRecord > propagatorToken_
T sqrt(T t)
Definition: SSEVec.h:19
T z() const
Definition: PV3DBase.h:61
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
edm::EDGetTokenT< reco::MuonCollection > recoMuonsToken_
double delta_R(double eta1, double phi1, double eta2, double phi2)
Definition: AnglesUtil.h:106
bool isValid() const
Definition: HandleBase.h:70
edm::EDGetTokenT< MeasurementTrackerEvent > measurementTrackerEventToken_
Log< level::Info, false > LogInfo
Vector3DBase unit() const
Definition: Vector3DBase.h:54
std::vector< RecoChargedCandidate > RecoChargedCandidateCollection
collectin of RecoChargedCandidate objects
T const * product() const
Definition: Handle.h:70
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
T eta() const
Definition: PV3DBase.h:73
GlobalVector globalMomentum() const
tuple MeasurementTrackerEvent
edm::EDGetTokenT< reco::RecoChargedCandidateCollection > recoL2MuonsToken_
CosmicRegionalSeedGenerator(const edm::ParameterSet &conf, edm::ConsumesCollector &&iC)
ESGetTokenH3DDVariant esConsumes(std::string const &Reccord, edm::ConsumesCollector &)
Definition: DeDxTools.cc:283
std::vector< CaloJet > CaloJetCollection
collection of CaloJet objects
#define LogDebug(id)