CMS 3D CMS Logo

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

#include <RecoMuon/MuonIdentification/src/CSCTimingExtractor.cc>

Classes

class  TimeMeasurement
 

Public Member Functions

 CSCTimingExtractor (const edm::ParameterSet &, MuonSegmentMatcher *segMatcher)
 Constructor. More...
 
void fillTiming (TimeMeasurementSequence &tmSequence, reco::TrackRef muonTrack, const edm::Event &iEvent, const edm::EventSetup &iSetup)
 
 ~CSCTimingExtractor ()
 Destructor. More...
 

Private Attributes

edm::InputTag CSCSegmentTags_
 
bool debug
 
unsigned int theHitsMin_
 
MuonSegmentMatchertheMatcher
 
double thePruneCut_
 
std::unique_ptr< MuonServiceProxytheService
 
double theStripError_
 
double theStripTimeOffset_
 
double theWireError_
 
double theWireTimeOffset_
 
bool UseStripTime
 
bool UseWireTime
 

Detailed Description

Extracts timing information associated to a muon track

Description: Produce timing information for a muon track using CSC hits from segments used to build the track

Definition at line 54 of file CSCTimingExtractor.h.

Constructor & Destructor Documentation

CSCTimingExtractor::CSCTimingExtractor ( const edm::ParameterSet iConfig,
MuonSegmentMatcher segMatcher 
)

Constructor.

Definition at line 69 of file CSCTimingExtractor.cc.

References edm::ParameterSet::getParameter(), theMatcher, and theService.

70  :
71  thePruneCut_(iConfig.getParameter<double>("PruneCut")),
72  theStripTimeOffset_(iConfig.getParameter<double>("CSCStripTimeOffset")),
73  theWireTimeOffset_(iConfig.getParameter<double>("CSCWireTimeOffset")),
74  theStripError_(iConfig.getParameter<double>("CSCStripError")),
75  theWireError_(iConfig.getParameter<double>("CSCWireError")),
76  UseWireTime(iConfig.getParameter<bool>("UseWireTime")),
77  UseStripTime(iConfig.getParameter<bool>("UseStripTime")),
78  debug(iConfig.getParameter<bool>("debug"))
79 {
80  edm::ParameterSet serviceParameters = iConfig.getParameter<edm::ParameterSet>("ServiceParameters");
81  theService = std::make_unique<MuonServiceProxy>(serviceParameters);
82  theMatcher = segMatcher;
83 }
T getParameter(std::string const &) const
std::unique_ptr< MuonServiceProxy > theService
MuonSegmentMatcher * theMatcher
CSCTimingExtractor::~CSCTimingExtractor ( )

Destructor.

Definition at line 86 of file CSCTimingExtractor.cc.

87 {
88 }

Member Function Documentation

void CSCTimingExtractor::fillTiming ( TimeMeasurementSequence tmSequence,
reco::TrackRef  muonTrack,
const edm::Event iEvent,
const edm::EventSetup iSetup 
)

Definition at line 97 of file CSCTimingExtractor.cc.

References gather_cfg::cout, debug, CSCTimingExtractor::TimeMeasurement::distIP, TimeMeasurementSequence::dstnc, edm::EventSetup::get(), mps_fire::i, GlobalTrackingGeometry::idToDet(), TimeMeasurementSequence::local_t0, mag(), MuonSegmentMatcher::matchCSC(), edm::ESHandle< T >::product(), Propagator::propagateWithPath(), PhotonConversionTrajectorySeedProducerFromQuadruplets_cfi::propagator, mathSSE::sqrt(), GeomDet::surface(), theMatcher, thePruneCut_, theService, theStripError_, theStripTimeOffset_, MuonServiceProxy::theTrackingGeometry, theWireError_, theWireTimeOffset_, CSCTimingExtractor::TimeMeasurement::timeCorr, GeomDet::toGlobal(), TimeMeasurementSequence::totalWeightInvbeta, TimeMeasurementSequence::totalWeightTimeVtx, UseStripTime, UseWireTime, TimeMeasurementSequence::weightInvbeta, CSCTimingExtractor::TimeMeasurement::weightInvbeta, TimeMeasurementSequence::weightTimeVtx, and CSCTimingExtractor::TimeMeasurement::weightTimeVtx.

99 {
100 
101  if (debug)
102  std::cout << " *** CSC Timimng Extractor ***" << std::endl;
103 
104  theService->update(iSetup);
105 
106  const GlobalTrackingGeometry *theTrackingGeometry = &*theService->trackingGeometry();
107 
108  // get the propagator
110  iSetup.get<TrackingComponentsRecord>().get("SteppingHelixPropagatorAny", propagator);
111  const Propagator *propag = propagator.product();
112 
113  math::XYZPoint pos=muonTrack->innerPosition();
114  math::XYZVector mom=muonTrack->innerMomentum();
115  if (sqrt(muonTrack->innerPosition().mag2()) > sqrt(muonTrack->outerPosition().mag2())){
116  pos=muonTrack->outerPosition();
117  mom=-1*muonTrack->outerMomentum();
118  }
119  GlobalPoint posp(pos.x(), pos.y(), pos.z());
120  GlobalVector momv(mom.x(), mom.y(), mom.z());
121  FreeTrajectoryState muonFTS(posp, momv, (TrackCharge)muonTrack->charge(), theService->magneticField().product());
122 
123  // get the CSC segments that were used to construct the muon
124  std::vector<const CSCSegment*> range = theMatcher->matchCSC(*muonTrack,iEvent);
125 
126  // create a collection on TimeMeasurements for the track
127  std::vector<TimeMeasurement> tms;
128  for (std::vector<const CSCSegment*>::iterator rechit = range.begin(); rechit!=range.end();++rechit) {
129 
130  // Create the ChamberId
131  DetId id = (*rechit)->geographicalId();
132  CSCDetId chamberId(id.rawId());
133  // int station = chamberId.station();
134 
135  if (!(*rechit)->specificRecHits().size()) continue;
136 
137  const std::vector<CSCRecHit2D> hits2d = (*rechit)->specificRecHits();
138 
139  // store all the hits from the segment
140  for (std::vector<CSCRecHit2D>::const_iterator hiti=hits2d.begin(); hiti!=hits2d.end(); hiti++) {
141 
142  const GeomDet* cscDet = theTrackingGeometry->idToDet(hiti->geographicalId());
143  TimeMeasurement thisHit;
144 
145  std::pair< TrajectoryStateOnSurface, double> tsos;
146  tsos=propag->propagateWithPath(muonFTS,cscDet->surface());
147 
148  double dist;
149  if (tsos.first.isValid()) dist = tsos.second+posp.mag();
150  else dist = cscDet->toGlobal(hiti->localPosition()).mag();
151 
152  thisHit.distIP = dist;
153  if (UseStripTime) {
154  thisHit.weightInvbeta = dist*dist/(theStripError_*theStripError_*30.*30.);
155  thisHit.weightTimeVtx = 1./(theStripError_*theStripError_);
156  thisHit.timeCorr = hiti->tpeak()-theStripTimeOffset_;
157  tms.push_back(thisHit);
158  }
159 
160  if (UseWireTime) {
161  thisHit.weightInvbeta = dist*dist/(theWireError_*theWireError_*30.*30.);
162  thisHit.weightTimeVtx = 1./(theWireError_*theWireError_);
163  thisHit.timeCorr = hiti->wireTime()-theWireTimeOffset_;
164  tms.push_back(thisHit);
165  }
166 
167 // std::cout << " CSC Hit. Dist= " << dist << " Time= " << thisHit.timeCorr
168 // << " invBeta= " << (1.+thisHit.timeCorr/dist*30.) << std::endl;
169  }
170 
171  } // rechit
172 
173  bool modified = false;
174  std::vector <double> dstnc, local_t0, hitWeightInvbeta, hitWeightTimeVtx;
175  double totalWeightInvbeta=0;
176  double totalWeightTimeVtx=0;
177 
178  // Now loop over the measurements, calculate 1/beta and cut away outliers
179  do {
180 
181  modified = false;
182  dstnc.clear();
183  local_t0.clear();
184  hitWeightInvbeta.clear();
185  hitWeightTimeVtx.clear();
186 
187  totalWeightInvbeta=0;
188  totalWeightTimeVtx=0;
189 
190  for (std::vector<TimeMeasurement>::iterator tm=tms.begin(); tm!=tms.end(); ++tm) {
191  dstnc.push_back(tm->distIP);
192  local_t0.push_back(tm->timeCorr);
193  hitWeightInvbeta.push_back(tm->weightInvbeta);
194  hitWeightTimeVtx.push_back(tm->weightTimeVtx);
195  totalWeightInvbeta+=tm->weightInvbeta;
196  totalWeightTimeVtx+=tm->weightTimeVtx;
197  }
198 
199  if (totalWeightInvbeta==0) break;
200 
201  // calculate the value and error of 1/beta and timeVtx from the complete set of 1D hits
202  if (debug)
203  std::cout << " Points for global fit: " << dstnc.size() << std::endl;
204 
205  double invbeta=0,invbetaErr=0;
206  double timeVtx=0,timeVtxErr=0;
207 
208  for (unsigned int i=0;i<dstnc.size();i++) {
209  invbeta+=(1.+local_t0.at(i)/dstnc.at(i)*30.)*hitWeightInvbeta.at(i)/totalWeightInvbeta;
210  timeVtx+=local_t0.at(i)*hitWeightTimeVtx.at(i)/totalWeightTimeVtx;
211  }
212 
213  double chimax=0.;
214  std::vector<TimeMeasurement>::iterator tmmax;
215 
216  // Calculate the inv beta and time at vertex dispersion
217  double diff_ibeta,diff_tvtx;
218  for (unsigned int i=0;i<dstnc.size();i++) {
219  diff_ibeta=(1.+local_t0.at(i)/dstnc.at(i)*30.)-invbeta;
220  diff_ibeta=diff_ibeta*diff_ibeta*hitWeightInvbeta.at(i);
221  diff_tvtx=local_t0.at(i)-timeVtx;
222  diff_tvtx=diff_tvtx*diff_tvtx*hitWeightTimeVtx.at(i);
223  invbetaErr+=diff_ibeta;
224  timeVtxErr+=diff_tvtx;
225 
226  // decide if we cut away time at vertex outliers or inverse beta outliers
227  // currently not configurable.
228  if (diff_tvtx>chimax) {
229  tmmax=tms.begin()+i;
230  chimax=diff_tvtx;
231  }
232  }
233 
234  // cut away the outliers
235  if (chimax>thePruneCut_) {
236  tms.erase(tmmax);
237  modified=true;
238  }
239 
240  if (debug) {
241  double cf = 1./(dstnc.size()-1);
242  invbetaErr=sqrt(invbetaErr/totalWeightInvbeta*cf);
243  timeVtxErr=sqrt(timeVtxErr/totalWeightTimeVtx*cf);
244  std::cout << " Measured 1/beta: " << invbeta << " +/- " << invbetaErr << std::endl;
245  std::cout << " Measured time: " << timeVtx << " +/- " << timeVtxErr << std::endl;
246  }
247 
248  } while (modified);
249 
250  for (unsigned int i=0;i<dstnc.size();i++) {
251  tmSequence.dstnc.push_back(dstnc.at(i));
252  tmSequence.local_t0.push_back(local_t0.at(i));
253  tmSequence.weightInvbeta.push_back(hitWeightInvbeta.at(i));
254  tmSequence.weightTimeVtx.push_back(hitWeightTimeVtx.at(i));
255  }
256 
257  tmSequence.totalWeightInvbeta=totalWeightInvbeta;
258  tmSequence.totalWeightTimeVtx=totalWeightTimeVtx;
259 
260 }
std::vector< double > local_t0
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:54
std::vector< const CSCSegment * > matchCSC(const reco::Track &muon, const edm::Event &event)
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:42
std::unique_ptr< MuonServiceProxy > theService
int TrackCharge
Definition: TrackCharge.h:4
MuonSegmentMatcher * theMatcher
T sqrt(T t)
Definition: SSEVec.h:18
virtual std::pair< TrajectoryStateOnSurface, double > propagateWithPath(const FreeTrajectoryState &, const Surface &) const final
Definition: Propagator.cc:15
std::vector< double > weightInvbeta
Definition: DetId.h:18
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:30
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
std::vector< double > weightTimeVtx
const T & get() const
Definition: EventSetup.h:55
const GeomDet * idToDet(DetId) const override
std::vector< double > dstnc
T const * product() const
Definition: ESHandle.h:86

Member Data Documentation

edm::InputTag CSCTimingExtractor::CSCSegmentTags_
private

Definition at line 78 of file CSCTimingExtractor.h.

bool CSCTimingExtractor::debug
private

Definition at line 87 of file CSCTimingExtractor.h.

Referenced by fillTiming().

unsigned int CSCTimingExtractor::theHitsMin_
private

Definition at line 79 of file CSCTimingExtractor.h.

MuonSegmentMatcher* CSCTimingExtractor::theMatcher
private

Definition at line 90 of file CSCTimingExtractor.h.

Referenced by CSCTimingExtractor(), and fillTiming().

double CSCTimingExtractor::thePruneCut_
private

Definition at line 80 of file CSCTimingExtractor.h.

Referenced by fillTiming().

std::unique_ptr<MuonServiceProxy> CSCTimingExtractor::theService
private

Definition at line 89 of file CSCTimingExtractor.h.

Referenced by CSCTimingExtractor(), and fillTiming().

double CSCTimingExtractor::theStripError_
private

Definition at line 83 of file CSCTimingExtractor.h.

Referenced by fillTiming().

double CSCTimingExtractor::theStripTimeOffset_
private

Definition at line 81 of file CSCTimingExtractor.h.

Referenced by fillTiming().

double CSCTimingExtractor::theWireError_
private

Definition at line 84 of file CSCTimingExtractor.h.

Referenced by fillTiming().

double CSCTimingExtractor::theWireTimeOffset_
private

Definition at line 82 of file CSCTimingExtractor.h.

Referenced by fillTiming().

bool CSCTimingExtractor::UseStripTime
private

Definition at line 86 of file CSCTimingExtractor.h.

Referenced by fillTiming().

bool CSCTimingExtractor::UseWireTime
private

Definition at line 85 of file CSCTimingExtractor.h.

Referenced by fillTiming().