test
CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MuonTrackProducer.cc
Go to the documentation of this file.
1 //
2 // modified & integrated by Giovanni Abbiendi
3 // from code by Arun Luthra: UserCode/luthra/MuonTrackSelector/src/MuonTrackSelector.cc
4 //
13 #include <sstream>
14 
18 
20  muonsToken(consumes<reco::MuonCollection>(parset.getParameter< edm::InputTag >("muonsTag"))),
21  inputDTRecSegment4DToken_(consumes<DTRecSegment4DCollection>(parset.getParameter<edm::InputTag>("inputDTRecSegment4DCollection"))),
22  inputCSCSegmentToken_(consumes<CSCSegmentCollection>(parset.getParameter<edm::InputTag>("inputCSCSegmentCollection"))),
23  selectionTags(parset.getParameter< std::vector<std::string> >("selectionTags")),
24  trackType(parset.getParameter< std::string >("trackType")),
25  parset_(parset)
26 {
27  edm::LogVerbatim("MuonTrackProducer") << "constructing MuonTrackProducer" << parset_.dump();
28  produces<reco::TrackCollection>();
29  produces<reco::TrackExtraCollection>();
30  produces<TrackingRecHitCollection>();
31 }
32 
34 }
35 
37 {
41 
43  iSetup.get<TrackerTopologyRcd>().get(httopo);
44  const TrackerTopology& ttopo = *httopo;
45 
46  std::auto_ptr<reco::TrackCollection> selectedTracks(new reco::TrackCollection);
47  std::auto_ptr<reco::TrackExtraCollection> selectedTrackExtras( new reco::TrackExtraCollection() );
48  std::auto_ptr<TrackingRecHitCollection> selectedTrackHits( new TrackingRecHitCollection() );
49 
53 
56 
57  edm::LogVerbatim("MuonTrackProducer") <<"\nThere are "<< dtSegmentCollectionH_->size()<<" DT segments.";
58  unsigned int index_dt_segment = 0;
60  segment != dtSegmentCollectionH_->end(); ++segment , index_dt_segment++) {
61  LocalPoint segmentLocalPosition = segment->localPosition();
62  LocalVector segmentLocalDirection = segment->localDirection();
63  LocalError segmentLocalPositionError = segment->localPositionError();
64  LocalError segmentLocalDirectionError = segment->localDirectionError();
65  DetId geoid = segment->geographicalId();
66  DTChamberId dtdetid = DTChamberId(geoid);
67  int wheel = dtdetid.wheel();
68  int station = dtdetid.station();
69  int sector = dtdetid.sector();
70 
71  float segmentX = segmentLocalPosition.x();
72  float segmentY = segmentLocalPosition.y();
73  float segmentdXdZ = segmentLocalDirection.x()/segmentLocalDirection.z();
74  float segmentdYdZ = segmentLocalDirection.y()/segmentLocalDirection.z();
75  float segmentXerr = sqrt(segmentLocalPositionError.xx());
76  float segmentYerr = sqrt(segmentLocalPositionError.yy());
77  float segmentdXdZerr = sqrt(segmentLocalDirectionError.xx());
78  float segmentdYdZerr = sqrt(segmentLocalDirectionError.yy());
79 
80  edm::LogVerbatim("MuonTrackProducer")
81  <<"\nDT segment index :"<<index_dt_segment
82  <<"\nchamber Wh:"<<wheel<<",St:"<<station<<",Se:"<<sector
83  <<"\nLocal Position (X,Y)=("<<segmentX<<","<<segmentY<<") +/- ("<<segmentXerr<<","<<segmentYerr<<"), "
84  <<"Local Direction (dXdZ,dYdZ)=("<<segmentdXdZ<<","<<segmentdYdZ<<") +/- ("<<segmentdXdZerr<<","<<segmentdYdZerr<<")";
85  }
86 
87  edm::LogVerbatim("MuonTrackProducer") <<"\nThere are "<< cscSegmentCollectionH_->size()<<" CSC segments.";
88  unsigned int index_csc_segment = 0;
90  segment != cscSegmentCollectionH_->end(); ++segment , index_csc_segment++) {
91  LocalPoint segmentLocalPosition = segment->localPosition();
92  LocalVector segmentLocalDirection = segment->localDirection();
93  LocalError segmentLocalPositionError = segment->localPositionError();
94  LocalError segmentLocalDirectionError = segment->localDirectionError();
95 
96  DetId geoid = segment->geographicalId();
97  CSCDetId cscdetid = CSCDetId(geoid);
98  int endcap = cscdetid.endcap();
99  int station = cscdetid.station();
100  int ring = cscdetid.ring();
101  int chamber = cscdetid.chamber();
102 
103  float segmentX = segmentLocalPosition.x();
104  float segmentY = segmentLocalPosition.y();
105  float segmentdXdZ = segmentLocalDirection.x()/segmentLocalDirection.z();
106  float segmentdYdZ = segmentLocalDirection.y()/segmentLocalDirection.z();
107  float segmentXerr = sqrt(segmentLocalPositionError.xx());
108  float segmentYerr = sqrt(segmentLocalPositionError.yy());
109  float segmentdXdZerr = sqrt(segmentLocalDirectionError.xx());
110  float segmentdYdZerr = sqrt(segmentLocalDirectionError.yy());
111 
112  edm::LogVerbatim("MuonTrackProducer")
113  <<"\nCSC segment index :"<<index_csc_segment
114  <<"\nchamber Endcap:"<<endcap<<",St:"<<station<<",Ri:"<<ring<<",Ch:"<<chamber
115  <<"\nLocal Position (X,Y)=("<<segmentX<<","<<segmentY<<") +/- ("<<segmentXerr<<","<<segmentYerr<<"), "
116  <<"Local Direction (dXdZ,dYdZ)=("<<segmentdXdZ<<","<<segmentdYdZ<<") +/- ("<<segmentdXdZerr<<","<<segmentdYdZerr<<")";
117  }
118 
119  edm::LogVerbatim("MuonTrackProducer") <<"\nThere are "<< muonCollectionH->size() <<" reco::Muons.";
120  unsigned int muon_index = 0;
121  for(reco::MuonCollection::const_iterator muon = muonCollectionH->begin();
122  muon != muonCollectionH->end(); ++muon, muon_index++) {
123  edm::LogVerbatim("MuonTrackProducer") <<"\n******* muon index : "<<muon_index;
124 
125  std::vector<bool> isGood;
126  for(unsigned int index=0; index<selectionTags.size(); ++index) {
127  isGood.push_back(false);
128 
130  isGood[index] = muon::isGoodMuon(*muon, muonType);
131  }
132 
133  bool isGoodResult=true;
134  for(unsigned int index=0; index<isGood.size(); ++index) {
135  edm::LogVerbatim("MuonTrackProducer") << "selectionTag = "<<selectionTags[index]<< ": "<<isGood[index]<<"\n";
136  isGoodResult *= isGood[index];
137  }
138 
139  if (isGoodResult) {
140  // new copy of Track
141  reco::TrackRef trackref;
142  if (trackType == "innerTrack") {
143  if (muon->innerTrack().isNonnull()) trackref = muon->innerTrack();
144  else continue;
145  }
146  else if (trackType == "outerTrack") {
147  if (muon->outerTrack().isNonnull()) trackref = muon->outerTrack();
148  else continue;
149  }
150  else if (trackType == "globalTrack") {
151  if (muon->globalTrack().isNonnull()) trackref = muon->globalTrack();
152  else continue;
153  }
154  else if (trackType == "innerTrackPlusSegments") {
155  if (muon->innerTrack().isNonnull()) trackref = muon->innerTrack();
156  else continue;
157  }
158 
159  const reco::Track* trk = &(*trackref);
160  // pointer to old track:
161  reco::Track* newTrk = new reco::Track(*trk);
162 
163  newTrk->setExtra( reco::TrackExtraRef( rTrackExtras, idx++ ) );
164  PropagationDirection seedDir = trk->seedDirection();
165  // new copy of track Extras
166  reco::TrackExtra * newExtra = new reco::TrackExtra( trk->outerPosition(), trk->outerMomentum(),
167  trk->outerOk(), trk->innerPosition(),
168  trk->innerMomentum(), trk->innerOk(),
169  trk->outerStateCovariance(), trk->outerDetId(),
170  trk->innerStateCovariance(), trk->innerDetId() , seedDir ) ;
171 
172  // new copy of the silicon hits; add hit refs to Extra and hits to hit collection
173 
174  // edm::LogVerbatim("MuonTrackProducer")<<"\n printing initial hit_pattern";
175  // trk->hitPattern().print();
176  unsigned int nHitsToAdd = 0;
177  for (trackingRecHit_iterator iHit = trk->recHitsBegin(); iHit != trk->recHitsEnd(); iHit++) {
178  TrackingRecHit* hit = (*iHit)->clone();
179  selectedTrackHits->push_back( hit );
180  ++nHitsToAdd;
181  }
182  newExtra->setHits( rHits, hidx, nHitsToAdd );
183  hidx += nHitsToAdd;
184  if (trackType == "innerTrackPlusSegments") {
185 
186  int wheel, station, sector;
187  int endcap, /*station, */ ring, chamber;
188 
189  edm::LogVerbatim("MuonTrackProducer") <<"Number of chambers: "<<muon->matches().size()
190  <<", arbitrated: "<<muon->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
191  unsigned int index_chamber = 0;
192 
193  for(std::vector<reco::MuonChamberMatch>::const_iterator chamberMatch = muon->matches().begin();
194  chamberMatch != muon->matches().end(); ++chamberMatch, index_chamber++) {
195  std::stringstream chamberStr;
196  chamberStr <<"\nchamber index: "<<index_chamber;
197 
198  int subdet = chamberMatch->detector();
199  DetId did = chamberMatch->id;
200 
201  if (subdet == MuonSubdetId::DT) {
202  DTChamberId dtdetid = DTChamberId(did);
203  wheel = dtdetid.wheel();
204  station = dtdetid.station();
205  sector = dtdetid.sector();
206  chamberStr << ", DT chamber Wh:"<<wheel<<",St:"<<station<<",Se:"<<sector;
207  }
208  else if (subdet == MuonSubdetId::CSC) {
209  CSCDetId cscdetid = CSCDetId(did);
210  endcap = cscdetid.endcap();
211  station = cscdetid.station();
212  ring = cscdetid.ring();
213  chamber = cscdetid.chamber();
214  chamberStr << ", CSC chamber End:"<<endcap<<",St:"<<station<<",Ri:"<<ring<<",Ch:"<<chamber;
215  }
216 
217  chamberStr << ", Number of segments: "<<chamberMatch->segmentMatches.size();
218  edm::LogVerbatim("MuonTrackProducer") << chamberStr.str();
219 
220  unsigned int index_segment = 0;
221 
222  for(std::vector<reco::MuonSegmentMatch>::const_iterator segmentMatch = chamberMatch->segmentMatches.begin();
223  segmentMatch != chamberMatch->segmentMatches.end(); ++segmentMatch, index_segment++) {
224 
225  float segmentX = segmentMatch->x;
226  float segmentY = segmentMatch->y ;
227  float segmentdXdZ = segmentMatch->dXdZ;
228  float segmentdYdZ = segmentMatch->dYdZ;
229  float segmentXerr = segmentMatch->xErr;
230  float segmentYerr = segmentMatch->yErr;
231  float segmentdXdZerr = segmentMatch->dXdZErr;
232  float segmentdYdZerr = segmentMatch->dYdZErr;
233 
234  CSCSegmentRef segmentCSC = segmentMatch->cscSegmentRef;
235  DTRecSegment4DRef segmentDT = segmentMatch->dtSegmentRef;
236 
237  bool segment_arbitrated_Ok = (segmentMatch->isMask(reco::MuonSegmentMatch::BestInChamberByDR) &&
238  segmentMatch->isMask(reco::MuonSegmentMatch::BelongsToTrackByDR));
239 
240  std::string ARBITRATED(" ***Arbitrated Off*** ");
241  if (segment_arbitrated_Ok) ARBITRATED = " ***ARBITRATED OK*** ";
242 
243  if (subdet == MuonSubdetId::DT) {
244  edm::LogVerbatim("MuonTrackProducer")
245  <<"\n\t segment index: "<<index_segment << ARBITRATED
246  <<"\n\t Local Position (X,Y)=("<<segmentX<<","<<segmentY<<") +/- ("<<segmentXerr<<","<<segmentYerr<<"), "
247  <<"\n\t Local Direction (dXdZ,dYdZ)=("<<segmentdXdZ<<","<<segmentdYdZ<<") +/- ("<<segmentdXdZerr<<","<<segmentdYdZerr<<")";
248 
249  if (!segment_arbitrated_Ok) continue;
250 
251  if (segmentDT.get() != 0) {
252  const DTRecSegment4D* segment = segmentDT.get();
253 
254  edm::LogVerbatim("MuonTrackProducer")<<"\t ===> MATCHING with DT segment with index = "<<segmentDT.key();
255 
256  if(segment->hasPhi()) {
257  const DTChamberRecSegment2D* phiSeg = segment->phiSegment();
258  std::vector<const TrackingRecHit*> phiHits = phiSeg->recHits();
259  unsigned int nHitsAdded = 0;
260  for(std::vector<const TrackingRecHit*>::const_iterator ihit = phiHits.begin();
261  ihit != phiHits.end(); ++ihit) {
262  TrackingRecHit* seghit = (*ihit)->clone();
263  newTrk->appendHitPattern(*seghit, ttopo);
264  // edm::LogVerbatim("MuonTrackProducer")<<"hit pattern for position "<<index_hit<<" set to:";
265  // newTrk->hitPattern().printHitPattern(index_hit, std::cout);
266  selectedTrackHits->push_back( seghit );
267  ++nHitsAdded;
268  }
269  newExtra->setHits( rHits, hidx, nHitsAdded );
270  hidx += nHitsAdded;
271  }
272 
273  if(segment->hasZed()) {
274  const DTSLRecSegment2D* zSeg = (*segment).zSegment();
275  std::vector<const TrackingRecHit*> zedHits = zSeg->recHits();
276  unsigned int nHitsAdded = 0;
277  for(std::vector<const TrackingRecHit*>::const_iterator ihit = zedHits.begin();
278  ihit != zedHits.end(); ++ihit) {
279  TrackingRecHit* seghit = (*ihit)->clone();
280  newTrk->appendHitPattern(*seghit, ttopo);
281  // edm::LogVerbatim("MuonTrackProducer")<<"hit pattern for position "<<index_hit<<" set to:";
282  // newTrk->hitPattern().printHitPattern(index_hit, std::cout);
283  selectedTrackHits->push_back( seghit );
284  ++nHitsAdded;
285  }
286  newExtra->setHits( rHits, hidx, nHitsAdded );
287  hidx += nHitsAdded;
288  }
289  } else edm::LogWarning("MuonTrackProducer")<<"\n***WARNING: UNMATCHED DT segment ! \n";
290  } // if (subdet == MuonSubdetId::DT)
291 
292  else if (subdet == MuonSubdetId::CSC) {
293  edm::LogVerbatim("MuonTrackProducer")
294  <<"\n\t segment index: "<<index_segment << ARBITRATED
295  <<"\n\t Local Position (X,Y)=("<<segmentX<<","<<segmentY<<") +/- ("<<segmentXerr<<","<<segmentYerr<<"), "
296  <<"\n\t Local Direction (dXdZ,dYdZ)=("<<segmentdXdZ<<","<<segmentdYdZ<<") +/- ("<<segmentdXdZerr<<","<<segmentdYdZerr<<")";
297 
298  if (!segment_arbitrated_Ok) continue;
299 
300  if (segmentCSC.get() != 0) {
301  const CSCSegment* segment = segmentCSC.get();
302 
303  edm::LogVerbatim("MuonTrackProducer")<<"\t ===> MATCHING with CSC segment with index = "<<segmentCSC.key();
304 
305  std::vector<const TrackingRecHit*> hits = segment->recHits();
306  unsigned int nHitsAdded = 0;
307  for(std::vector<const TrackingRecHit*>::const_iterator ihit = hits.begin();
308  ihit != hits.end(); ++ihit) {
309  TrackingRecHit* seghit = (*ihit)->clone();
310  newTrk->appendHitPattern(*seghit, ttopo);
311  // edm::LogVerbatim("MuonTrackProducer")<<"hit pattern for position "<<index_hit<<" set to:";
312  // newTrk->hitPattern().printHitPattern(index_hit, std::cout);
313  selectedTrackHits->push_back( seghit );
314  ++nHitsAdded;
315  }
316  newExtra->setHits( rHits, hidx, nHitsAdded );
317  hidx += nHitsAdded;
318  } else edm::LogWarning("MuonTrackProducer")<<"\n***WARNING: UNMATCHED CSC segment ! \n";
319  } // else if (subdet == MuonSubdetId::CSC)
320 
321  } // loop on vector<MuonSegmentMatch>
322  } // loop on vector<MuonChamberMatch>
323  } // if (trackType == "innerTrackPlusSegments")
324 
325  // edm::LogVerbatim("MuonTrackProducer")<<"\n printing final hit_pattern";
326  // newTrk->hitPattern().print();
327 
328  selectedTracks->push_back( *newTrk );
329  selectedTrackExtras->push_back( *newExtra );
330 
331  } // if (isGoodResult)
332  } // loop on reco::MuonCollection
333 
334  iEvent.put(selectedTracks);
335  iEvent.put(selectedTrackExtras);
336  iEvent.put(selectedTrackHits);
337 }
int chamber() const
Definition: CSCDetId.h:68
float xx() const
Definition: LocalError.h:24
edm::EDGetTokenT< DTRecSegment4DCollection > inputDTRecSegment4DToken_
std::vector< std::string > selectionTags
std::string dump(unsigned int indent=0) const
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:462
void setHits(TrackingRecHitRefProd const &prod, unsigned firstH, unsigned int nH)
const DTChamberRecSegment2D * phiSegment() const
The superPhi segment: 0 if no phi projection available.
edm::EDGetTokenT< reco::MuonCollection > muonsToken
std::vector< Track > TrackCollection
collection of Tracks
Definition: TrackFwd.h:14
T y() const
Definition: PV3DBase.h:63
bool innerOk() const
return true if the innermost hit is valid
Definition: Track.h:50
PropagationDirection
key_type key() const
Accessor for product key.
Definition: Ref.h:264
const math::XYZPoint & outerPosition() const
position of the outermost hit
Definition: Track.h:65
std::vector< Muon > MuonCollection
collection of Muon objects
Definition: MuonFwd.h:9
SelectionType
Selector type.
Definition: MuonSelectors.h:17
const math::XYZPoint & innerPosition() const
position of the innermost hit
Definition: Track.h:55
int endcap() const
Definition: CSCDetId.h:93
int iEvent
Definition: GenABIO.cc:230
static const int CSC
Definition: MuonSubdetId.h:13
edm::Handle< CSCSegmentCollection > cscSegmentCollectionH_
float yy() const
Definition: LocalError.h:26
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:121
virtual std::vector< const TrackingRecHit * > recHits() const
Access to component RecHits (if any)
T sqrt(T t)
Definition: SSEVec.h:18
CovarianceMatrix outerStateCovariance() const
outermost trajectory state curvilinear errors
Definition: Track.h:75
T z() const
Definition: PV3DBase.h:64
edm::Handle< DTRecSegment4DCollection > dtSegmentCollectionH_
unsigned int outerDetId() const
DetId of the detector on which surface the outermost state is located.
Definition: Track.h:94
static const unsigned int BestInChamberByDR
T const * get() const
Returns C++ pointer to the item.
Definition: Ref.h:244
bool hasPhi() const
Does it have the Phi projection?
virtual void produce(edm::Event &, const edm::EventSetup &)
trackingRecHit_iterator recHitsBegin() const
Iterator to first hit on the track.
Definition: Track.h:104
virtual std::vector< const TrackingRecHit * > recHits() const
Access to component RecHits (if any)
Definition: CSCSegment.cc:30
RefProd< PROD > getRefBeforePut()
Definition: Event.h:141
bool isGoodMuon(const reco::Muon &muon, SelectionType type, reco::Muon::ArbitrationType arbitrationType=reco::Muon::SegmentAndTrackArbitration)
main GoodMuon wrapper call
virtual TrackingRecHit * clone() const =0
MuonTrackProducer(const edm::ParameterSet &)
int ring() const
Definition: CSCDetId.h:75
bool hasZed() const
Does it have the Z projection?
Definition: DetId.h:18
std::vector< TrackExtra > TrackExtraCollection
collection of TrackExtra objects
Definition: TrackExtraFwd.h:11
edm::OwnVector< TrackingRecHit > TrackingRecHitCollection
collection of TrackingRecHits
void setExtra(const TrackExtraRef &ref)
set reference to &quot;extra&quot; object
Definition: Track.h:184
const math::XYZVector & outerMomentum() const
momentum vector at the outermost hit position
Definition: Track.h:70
bool outerOk() const
return true if the outermost hit is valid
Definition: Track.h:45
SelectionType selectionTypeFromString(const std::string &label)
Definition: MuonSelectors.cc:9
tuple idx
DEBUGGING if hasattr(process,&quot;trackMonIterativeTracking2012&quot;): print &quot;trackMonIterativeTracking2012 D...
edm::Handle< reco::MuonCollection > muonCollectionH
const T & get() const
Definition: EventSetup.h:56
CovarianceMatrix innerStateCovariance() const
innermost trajectory state curvilinear errors
Definition: Track.h:80
bool appendHitPattern(const TrackingRecHit &hit, const TrackerTopology &ttopo)
append a single hit to the HitPattern
Definition: TrackBase.h:455
static const unsigned int BelongsToTrackByDR
const edm::ParameterSet parset_
int sector() const
Definition: DTChamberId.h:61
const math::XYZVector & innerMomentum() const
momentum vector at the innermost hit position
Definition: Track.h:60
PropagationDirection seedDirection() const
direction of how the hits were sorted in the original seed
Definition: Track.h:204
int station() const
Definition: CSCDetId.h:86
static const int DT
Definition: MuonSubdetId.h:12
int station() const
Return the station number.
Definition: DTChamberId.h:51
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:45
T x() const
Definition: PV3DBase.h:62
unsigned int innerDetId() const
DetId of the detector on which surface the innermost state is located.
Definition: Track.h:99
edm::EDGetTokenT< CSCSegmentCollection > inputCSCSegmentToken_
boost::remove_cv< typename boost::remove_reference< argument_type >::type >::type key_type
Definition: Ref.h:168
TrackingRecHitCollection::base::const_iterator trackingRecHit_iterator
iterator over a vector of reference to TrackingRecHit in the same collection
trackingRecHit_iterator recHitsEnd() const
Iterator to last hit on the track.
Definition: Track.h:109