Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:26:56 

0001 /**
0002  *  Class: GlobalMuonTrackMatcher
0003  *
0004  * 
0005  *  
0006  *  \author Chang Liu - Purdue University
0007  *  \author Norbert Neumeister - Purdue University
0008  *  \author Adam Everett - Purdue University
0009  *  \author Edwin Antillon - Purdue University
0010  *
0011  */
0012 
0013 #include "RecoMuon/GlobalTrackingTools/interface/GlobalMuonTrackMatcher.h"
0014 
0015 //---------------
0016 // C++ Headers --
0017 //---------------
0018 
0019 //-------------------------------
0020 // Collaborating Class Headers --
0021 //-------------------------------
0022 
0023 #include "FWCore/ServiceRegistry/interface/Service.h"
0024 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0025 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0026 
0027 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0028 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0029 #include "TrackingTools/GeomPropagators/interface/StateOnTrackerBound.h"
0030 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0031 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0032 
0033 #include "RecoMuon/TrackingTools/interface/MuonServiceProxy.h"
0034 
0035 #include "DataFormats/TrackReco/interface/Track.h"
0036 #include "DataFormats/TrajectoryState/interface/LocalTrajectoryParameters.h"
0037 
0038 #include "DataFormats/GeometrySurface/interface/TangentPlane.h"
0039 #include "DataFormats/Math/interface/deltaR.h"
0040 
0041 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0042 
0043 using namespace std;
0044 using namespace reco;
0045 
0046 //
0047 // constructor
0048 //
0049 GlobalMuonTrackMatcher::GlobalMuonTrackMatcher(const edm::ParameterSet& par, const MuonServiceProxy* service)
0050     : theService(service) {
0051   theMinP = par.getParameter<double>("MinP");
0052   theMinPt = par.getParameter<double>("MinPt");
0053   thePt_threshold1 = par.getParameter<double>("Pt_threshold1");
0054   thePt_threshold2 = par.getParameter<double>("Pt_threshold2");
0055   theEta_threshold = par.getParameter<double>("Eta_threshold");
0056   theChi2_1 = par.getParameter<double>("Chi2Cut_1");
0057   theChi2_2 = par.getParameter<double>("Chi2Cut_2");
0058   theChi2_3 = par.getParameter<double>("Chi2Cut_3");
0059   theLocChi2 = par.getParameter<double>("LocChi2Cut");
0060   theDeltaD_1 = par.getParameter<double>("DeltaDCut_1");
0061   theDeltaD_2 = par.getParameter<double>("DeltaDCut_2");
0062   theDeltaD_3 = par.getParameter<double>("DeltaDCut_3");
0063   theDeltaR_1 = par.getParameter<double>("DeltaRCut_1");
0064   theDeltaR_2 = par.getParameter<double>("DeltaRCut_2");
0065   theDeltaR_3 = par.getParameter<double>("DeltaRCut_3");
0066   theQual_1 = par.getParameter<double>("Quality_1");
0067   theQual_2 = par.getParameter<double>("Quality_2");
0068   theQual_3 = par.getParameter<double>("Quality_3");
0069   theOutPropagatorName = par.getParameter<string>("Propagator");
0070 }
0071 
0072 //
0073 // destructor
0074 //
0075 GlobalMuonTrackMatcher::~GlobalMuonTrackMatcher() {}
0076 
0077 //
0078 // check if two tracks are compatible with the *tight* criteria
0079 //
0080 bool GlobalMuonTrackMatcher::matchTight(const TrackCand& sta, const TrackCand& track) const {
0081   std::pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface> tsosPair = convertToTSOSMuHit(sta, track);
0082 
0083   double chi2 = match_Chi2(tsosPair.first, tsosPair.second);
0084   if (chi2 > 0. && chi2 < theChi2_2)
0085     return true;
0086 
0087   double distance = match_d(tsosPair.first, tsosPair.second);
0088   if (distance > 0. && distance < theDeltaD_2)
0089     return true;
0090 
0091   //double deltaR = match_Rpos(tsosPair.first,tsosPair.second);
0092   //if ( deltaR > 0. && deltaR < theDeltaR_3 ) return true;
0093 
0094   return false;
0095 }
0096 
0097 //
0098 // check if two tracks are compatible
0099 //
0100 double GlobalMuonTrackMatcher::match(const TrackCand& sta,
0101                                      const TrackCand& track,
0102                                      int matchOption,
0103                                      int surfaceOption) const {
0104   std::pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface> tsosPair;
0105   if (surfaceOption == 0)
0106     tsosPair = convertToTSOSMuHit(sta, track);
0107   if (surfaceOption == 1)
0108     tsosPair = convertToTSOSTkHit(sta, track);
0109   if (surfaceOption != 0 && surfaceOption != 1)
0110     return -1.0;
0111 
0112   if (matchOption == 0) {
0113     // chi^2
0114     return match_Chi2(tsosPair.first, tsosPair.second);
0115   } else if (matchOption == 1) {
0116     // distance
0117     return match_d(tsosPair.first, tsosPair.second);
0118   } else if (matchOption == 2) {
0119     // deltaR
0120     return match_Rpos(tsosPair.first, tsosPair.second);
0121   } else if (matchOption == 3) {
0122     return match_dist(tsosPair.first, tsosPair.second);
0123   } else {
0124     return -1.0;
0125   }
0126 }
0127 
0128 //
0129 // choose the track from a TrackCollection which best
0130 // matches a given standalone muon track
0131 //
0132 std::vector<GlobalMuonTrackMatcher::TrackCand>::const_iterator GlobalMuonTrackMatcher::matchOne(
0133     const TrackCand& sta, const std::vector<TrackCand>& tracks) const {
0134   if (tracks.empty())
0135     return tracks.end();
0136 
0137   double minChi2 = 1000.0;
0138   vector<TrackCand>::const_iterator result = tracks.end();
0139   for (vector<TrackCand>::const_iterator is = tracks.begin(); is != tracks.end(); ++is) {
0140     // propagate to common surface
0141     std::pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface> tsosPair = convertToTSOSMuHit(sta, *is);
0142 
0143     // calculate chi^2 of local track parameters
0144     double chi2 = match_Chi2(tsosPair.first, tsosPair.second);
0145     if (chi2 > 0. && chi2 <= minChi2) {
0146       minChi2 = chi2;
0147       result = is;
0148     }
0149   }
0150 
0151   return result;
0152 }
0153 
0154 //
0155 // choose a vector of tracks from a TrackCollection that are compatible
0156 // with a given standalone track. The order of checks for compatability
0157 // * for low momentum: use chi2 selection
0158 // * high momentum: use direction or local position
0159 //
0160 vector<GlobalMuonTrackMatcher::TrackCand> GlobalMuonTrackMatcher::match(const TrackCand& sta,
0161                                                                         const vector<TrackCand>& tracks) const {
0162   const string category = "GlobalMuonTrackMatcher";
0163 
0164   vector<TrackCand> result;
0165 
0166   if (tracks.empty())
0167     return result;
0168 
0169   typedef std::pair<TrackCand, TrajectoryStateOnSurface> TrackCandWithTSOS;
0170   vector<TrackCandWithTSOS> cands;
0171   int iiTk = 1;
0172   TrajectoryStateOnSurface muonTSOS;
0173 
0174   LogTrace(category) << "   ***" << endl << "STA Muon pT " << sta.second->pt();
0175   LogTrace(category) << "   Tk in Region " << tracks.size() << endl;
0176 
0177   for (vector<TrackCand>::const_iterator is = tracks.begin(); is != tracks.end(); ++is, iiTk++) {
0178     // propagate to a common surface
0179     std::pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface> tsosPair = convertToTSOSMuHit(sta, *is);
0180     LogTrace(category) << "    Tk " << iiTk << " of " << tracks.size() << "  ConvertToMuHitSurface muon isValid "
0181                        << tsosPair.first.isValid() << " tk isValid " << tsosPair.second.isValid() << endl;
0182     if (tsosPair.first.isValid())
0183       muonTSOS = tsosPair.first;
0184     cands.push_back(TrackCandWithTSOS(*is, tsosPair.second));
0185   }
0186 
0187   // initialize variables
0188   double min_chisq = 999999;
0189   double min_d = 999999;
0190   double min_de = 999999;
0191   double min_r_pos = 999999;
0192   std::vector<bool> passes(cands.size(), false);
0193   int jj = 0;
0194 
0195   int iTkCand = 1;
0196   for (vector<TrackCandWithTSOS>::const_iterator ii = cands.begin(); ii != cands.end(); ++ii, jj++, iTkCand++) {
0197     // tracks that are able not able propagate to a common surface
0198     if (!muonTSOS.isValid() || !(*ii).second.isValid())
0199       continue;
0200 
0201     // calculate matching variables
0202     double distance = match_d(muonTSOS, (*ii).second);
0203     double chi2 = match_Chi2(muonTSOS, (*ii).second);
0204     double loc_chi2 = match_dist(muonTSOS, (*ii).second);
0205     double deltaR = match_Rpos(muonTSOS, (*ii).second);
0206 
0207     LogTrace(category) << "   iTk " << iTkCand << " of " << cands.size() << " eta "
0208                        << (*ii).second.globalPosition().eta() << " phi " << (*ii).second.globalPosition().phi() << endl;
0209     LogTrace(category) << "    distance " << distance << " distance cut "
0210                        << " " << endl;
0211     LogTrace(category) << "    chi2 " << chi2 << " chi2 cut "
0212                        << " " << endl;
0213     LogTrace(category) << "    loc_chi2 " << loc_chi2 << " locChi2 cut "
0214                        << " " << endl;
0215     LogTrace(category) << "    deltaR " << deltaR << " deltaR cut "
0216                        << " " << endl;
0217 
0218     if ((*ii).second.globalMomentum().perp() < thePt_threshold1) {
0219       LogTrace(category) << "    Enters  a1" << endl;
0220 
0221       if ((chi2 > 0 && fabs((*ii).second.globalMomentum().eta()) < theEta_threshold && chi2 < theChi2_1) ||
0222           (distance > 0 && distance / (*ii).first.second->pt() < theDeltaD_1 && loc_chi2 > 0 &&
0223            loc_chi2 < theLocChi2)) {
0224         LogTrace(category) << "    Passes a1" << endl;
0225         result.push_back((*ii).first);
0226         passes[jj] = true;
0227       }
0228     }
0229     if ((passes[jj] == false) && (*ii).second.globalMomentum().perp() < thePt_threshold2) {
0230       LogTrace(category) << "    Enters a2" << endl;
0231       if ((chi2 > 0 && chi2 < theChi2_2) || (distance > 0 && distance < theDeltaD_2)) {
0232         LogTrace(category) << "    Passes a2" << endl;
0233         result.push_back((*ii).first);
0234         passes[jj] = true;
0235       }
0236     } else {
0237       LogTrace(category) << "    Enters a3" << endl;
0238       if (distance > 0 && distance < theDeltaD_3 && deltaR > 0 && deltaR < theDeltaR_1) {
0239         LogTrace(category) << "    Passes a3" << endl;
0240         result.push_back((*ii).first);
0241         passes[jj] = true;
0242       }
0243     }
0244 
0245     if (passes[jj]) {
0246       if (distance < min_d)
0247         min_d = distance;
0248       if (loc_chi2 < min_de)
0249         min_de = loc_chi2;
0250       if (deltaR < min_r_pos)
0251         min_r_pos = deltaR;
0252       if (chi2 < min_chisq)
0253         min_chisq = chi2;
0254     }
0255   }
0256 
0257   // re-initialize mask counter
0258   jj = 0;
0259 
0260   if (result.empty()) {
0261     LogTrace(category) << "   Stage 1 returned 0 results";
0262     for (vector<TrackCandWithTSOS>::const_iterator is = cands.begin(); is != cands.end(); ++is, jj++) {
0263       double deltaR = match_Rpos(muonTSOS, (*is).second);
0264 
0265       if (muonTSOS.isValid() && (*is).second.isValid()) {
0266         // check matching between tracker and muon tracks using dEta cut looser then dPhi cut
0267         LogTrace(category) << "    Stage 2 deltaR " << deltaR << " deltaEta "
0268                            << fabs((*is).second.globalPosition().eta() - muonTSOS.globalPosition().eta() <
0269                                    1.5 * theDeltaR_2)
0270                            << " deltaPhi "
0271                            << (fabs(deltaPhi((*is).second.globalPosition().barePhi(),
0272                                              muonTSOS.globalPosition().barePhi())) < theDeltaR_2)
0273                            << endl;
0274 
0275         if (fabs((*is).second.globalPosition().eta() - muonTSOS.globalPosition().eta()) < 1.5 * theDeltaR_2 &&
0276             fabs(deltaPhi((*is).second.globalPosition().barePhi(), muonTSOS.globalPosition().barePhi())) <
0277                 theDeltaR_2) {
0278           result.push_back((*is).first);
0279           passes[jj] = true;
0280         }
0281       }
0282 
0283       if (passes[jj]) {
0284         double distance = match_d(muonTSOS, (*is).second);
0285         double chi2 = match_Chi2(muonTSOS, (*is).second);
0286         double loc_chi2 = match_dist(muonTSOS, (*is).second);
0287         if (distance < min_d)
0288           min_d = distance;
0289         if (loc_chi2 < min_de)
0290           min_de = loc_chi2;
0291         if (deltaR < min_r_pos)
0292           min_r_pos = deltaR;
0293         if (chi2 < min_chisq)
0294           min_chisq = chi2;
0295       }
0296     }
0297   }
0298 
0299   for (vector<TrackCand>::const_iterator iTk = result.begin(); iTk != result.end(); ++iTk) {
0300     LogTrace(category) << "   -----" << endl
0301                        << "selected pt " << iTk->second->pt() << " eta " << iTk->second->eta() << " phi "
0302                        << iTk->second->phi() << endl;
0303   }
0304 
0305   if (result.size() < 2)
0306     return result;
0307   else
0308     result.clear();
0309 
0310   LogTrace(category) << "   Cleaning matched candiates" << endl;
0311 
0312   // re-initialize mask counter
0313   jj = 0;
0314 
0315   for (vector<TrackCandWithTSOS>::const_iterator is = cands.begin(); is != cands.end(); ++is, jj++) {
0316     if (!passes[jj])
0317       continue;
0318 
0319     double distance = match_d(muonTSOS, (*is).second);
0320     double chi2 = match_Chi2(muonTSOS, (*is).second);
0321     //unused    double loc_chi2 = match_dist(muonTSOS,(*is).second);
0322     double deltaR = match_Rpos(muonTSOS, (*is).second);
0323 
0324     // compute quality as the relative ratio to the minimum found for each variable
0325 
0326     int qual = (int)(chi2 / min_chisq + distance / min_d + deltaR / min_r_pos);
0327     int n_min =
0328         ((chi2 / min_chisq == 1) ? 1 : 0) + ((distance / min_d == 1) ? 1 : 0) + ((deltaR / min_r_pos == 1) ? 1 : 0);
0329 
0330     if (n_min == 3) {
0331       result.push_back((*is).first);
0332     }
0333 
0334     if (n_min == 2 && qual < theQual_1) {
0335       result.push_back((*is).first);
0336     }
0337 
0338     if (n_min == 1 && qual < theQual_2) {
0339       result.push_back((*is).first);
0340     }
0341 
0342     if (n_min == 0 && qual < theQual_3) {
0343       result.push_back((*is).first);
0344     }
0345   }
0346 
0347   for (vector<TrackCand>::const_iterator iTk = result.begin(); iTk != result.end(); ++iTk) {
0348     LogTrace(category) << "   -----" << endl
0349                        << "selected pt " << iTk->second->pt() << " eta " << iTk->second->eta() << " phi "
0350                        << iTk->second->phi() << endl;
0351   }
0352 
0353   return result;
0354 }
0355 
0356 //
0357 // propagate the two track candidates to the tracker bound surface
0358 //
0359 std::pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface> GlobalMuonTrackMatcher::convertToTSOSTk(
0360     const TrackCand& staCand, const TrackCand& tkCand) const {
0361   const string category = "GlobalMuonTrackMatcher";
0362 
0363   TrajectoryStateOnSurface empty;
0364 
0365   TransientTrack muTT(*staCand.second, &*theService->magneticField(), theService->trackingGeometry());
0366   TrajectoryStateOnSurface impactMuTSOS = muTT.impactPointState();
0367 
0368   TrajectoryStateOnSurface outerTkTsos;
0369   if (tkCand.second.isNonnull()) {
0370     // make sure the tracker track has enough momentum to reach the muon chambers
0371     if (!(tkCand.second->p() < theMinP || tkCand.second->pt() < theMinPt)) {
0372       outerTkTsos = trajectoryStateTransform::outerStateOnSurface(
0373           *tkCand.second, *theService->trackingGeometry(), &*theService->magneticField());
0374     }
0375   }
0376 
0377   if (!impactMuTSOS.isValid() || !outerTkTsos.isValid())
0378     return pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface>(empty, empty);
0379 
0380   // define StateOnTrackerBound object
0381   StateOnTrackerBound fromInside(&*theService->propagator(theOutPropagatorName));
0382 
0383   // extrapolate to outer tracker surface
0384   TrajectoryStateOnSurface tkTsosFromMu = fromInside(impactMuTSOS);
0385   TrajectoryStateOnSurface tkTsosFromTk = fromInside(outerTkTsos);
0386 
0387   if (!samePlane(tkTsosFromMu, tkTsosFromTk)) {
0388     // propagate tracker track to same surface as muon
0389     bool same1, same2;
0390     TrajectoryStateOnSurface newTkTsosFromTk, newTkTsosFromMu;
0391     if (tkTsosFromMu.isValid())
0392       newTkTsosFromTk = theService->propagator(theOutPropagatorName)->propagate(outerTkTsos, tkTsosFromMu.surface());
0393     same1 = samePlane(newTkTsosFromTk, tkTsosFromMu);
0394     LogTrace(category) << "Propagating to same tracker surface (Mu):" << same1;
0395     if (!same1) {
0396       if (tkTsosFromTk.isValid())
0397         newTkTsosFromMu = theService->propagator(theOutPropagatorName)->propagate(impactMuTSOS, tkTsosFromTk.surface());
0398       same2 = samePlane(newTkTsosFromMu, tkTsosFromTk);
0399       LogTrace(category) << "Propagating to same tracker surface (Tk):" << same2;
0400     }
0401     if (same1)
0402       tkTsosFromTk = newTkTsosFromTk;
0403     else if (same2)
0404       tkTsosFromMu = newTkTsosFromMu;
0405     else {
0406       LogTrace(category) << "Could not propagate Muon and Tracker track to the same tracker bound!";
0407       return pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface>(empty, empty);
0408     }
0409   }
0410 
0411   return pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface>(tkTsosFromMu, tkTsosFromTk);
0412 }
0413 
0414 //
0415 // propagate the two track candidates to the surface of the innermost muon hit
0416 //
0417 std::pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface> GlobalMuonTrackMatcher::convertToTSOSMuHit(
0418     const TrackCand& staCand, const TrackCand& tkCand) const {
0419   const string category = "GlobalMuonTrackMatcher";
0420   TrajectoryStateOnSurface empty;
0421   TransientTrack muTT(*staCand.second, &*theService->magneticField(), theService->trackingGeometry());
0422   TrajectoryStateOnSurface innerMuTSOS = muTT.innermostMeasurementState();
0423   TrajectoryStateOnSurface outerTkTsos, innerTkTsos;
0424   if (tkCand.second.isNonnull()) {
0425     // make sure the tracker track has enough momentum to reach the muon chambers
0426     if (!(tkCand.second->p() < theMinP || tkCand.second->pt() < theMinPt)) {
0427       TrajectoryStateOnSurface innerTkTsos;
0428 
0429       outerTkTsos = trajectoryStateTransform::outerStateOnSurface(
0430           *tkCand.second, *theService->trackingGeometry(), &*theService->magneticField());
0431       innerTkTsos = trajectoryStateTransform::innerStateOnSurface(
0432           *tkCand.second, *theService->trackingGeometry(), &*theService->magneticField());
0433       // for cosmics, outer-most referst to last traversed layer
0434       if ((innerMuTSOS.globalPosition() - outerTkTsos.globalPosition()).mag() >
0435           (innerMuTSOS.globalPosition() - innerTkTsos.globalPosition()).mag())
0436         outerTkTsos = innerTkTsos;
0437     }
0438   }
0439 
0440   if (!innerMuTSOS.isValid() || !outerTkTsos.isValid()) {
0441     LogTrace(category) << "A TSOS validity problem! MuTSOS " << innerMuTSOS.isValid() << " TkTSOS "
0442                        << outerTkTsos.isValid();
0443     return pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface>(empty, empty);
0444   }
0445 
0446   const Surface& refSurface = innerMuTSOS.surface();
0447   TrajectoryStateOnSurface tkAtMu =
0448       theService->propagator(theOutPropagatorName)->propagate(*outerTkTsos.freeState(), refSurface);
0449 
0450   if (!tkAtMu.isValid()) {
0451     LogTrace(category) << "Could not propagate Muon and Tracker track to the same muon hit surface!";
0452     return pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface>(empty, empty);
0453   }
0454 
0455   return pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface>(innerMuTSOS, tkAtMu);
0456 }
0457 
0458 //
0459 // propagate the two track candidates to the surface of the outermost tracker hit
0460 //
0461 std::pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface> GlobalMuonTrackMatcher::convertToTSOSTkHit(
0462     const TrackCand& staCand, const TrackCand& tkCand) const {
0463   const string category = "GlobalMuonTrackMatcher";
0464 
0465   TrajectoryStateOnSurface empty;
0466 
0467   TransientTrack muTT(*staCand.second, &*theService->magneticField(), theService->trackingGeometry());
0468   TrajectoryStateOnSurface impactMuTSOS = muTT.impactPointState();
0469   TrajectoryStateOnSurface innerMuTSOS = muTT.innermostMeasurementState();
0470 
0471   TrajectoryStateOnSurface outerTkTsos, innerTkTsos;
0472   if (tkCand.second.isNonnull()) {
0473     // make sure the tracker track has enough momentum to reach the muon chambers
0474     if (!(tkCand.second->p() < theMinP || tkCand.second->pt() < theMinPt)) {
0475       outerTkTsos = trajectoryStateTransform::outerStateOnSurface(
0476           *tkCand.second, *theService->trackingGeometry(), &*theService->magneticField());
0477       innerTkTsos = trajectoryStateTransform::innerStateOnSurface(
0478           *tkCand.second, *theService->trackingGeometry(), &*theService->magneticField());
0479 
0480       // for cosmics, outer-most referst to last traversed layer
0481       if ((innerMuTSOS.globalPosition() - outerTkTsos.globalPosition()).mag() >
0482           (innerMuTSOS.globalPosition() - innerTkTsos.globalPosition()).mag())
0483         outerTkTsos = innerTkTsos;
0484     }
0485   }
0486 
0487   if (!impactMuTSOS.isValid() || !outerTkTsos.isValid()) {
0488     LogTrace(category) << "A TSOS validity problem! MuTSOS " << impactMuTSOS.isValid() << " TkTSOS "
0489                        << outerTkTsos.isValid();
0490     return pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface>(empty, empty);
0491   }
0492 
0493   const Surface& refSurface = outerTkTsos.surface();
0494   TrajectoryStateOnSurface muAtTk =
0495       theService->propagator(theOutPropagatorName)->propagate(*impactMuTSOS.freeState(), refSurface);
0496 
0497   if (!muAtTk.isValid()) {
0498     LogTrace(category) << "Could not propagate Muon and Tracker track to the same tracker hit surface!";
0499     return pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface>(empty, empty);
0500   }
0501 
0502   return pair<TrajectoryStateOnSurface, TrajectoryStateOnSurface>(muAtTk, outerTkTsos);
0503 }
0504 
0505 //
0506 //
0507 //
0508 bool GlobalMuonTrackMatcher::samePlane(const TrajectoryStateOnSurface& tsos1,
0509                                        const TrajectoryStateOnSurface& tsos2) const {
0510   if (!tsos1.isValid() || !tsos2.isValid())
0511     return false;
0512 
0513   if (fabs(match_D(tsos1, tsos2) - match_d(tsos1, tsos2)) > 0.1)
0514     return false;
0515 
0516   const float maxtilt = 0.999;
0517   const float maxdist = 0.01;  // in cm
0518 
0519   auto p1(tsos1.surface().tangentPlane(tsos1.localPosition()));
0520   auto p2(tsos2.surface().tangentPlane(tsos2.localPosition()));
0521 
0522   bool returnValue = ((fabs(p1->normalVector().dot(p2->normalVector())) > maxtilt) ||
0523                       (fabs((p1->toLocal(p2->position())).z()) < maxdist))
0524                          ? true
0525                          : false;
0526 
0527   return returnValue;
0528 }
0529 
0530 //
0531 // calculate Chi^2 of two trajectory states
0532 //
0533 double GlobalMuonTrackMatcher::match_Chi2(const TrajectoryStateOnSurface& tsos1,
0534                                           const TrajectoryStateOnSurface& tsos2) const {
0535   const string category = "GlobalMuonTrackMatcher";
0536   //LogTrace(category) << "match_Chi2 sanity check: " << tsos1.isValid() << " " << tsos2.isValid();
0537   if (!tsos1.isValid() || !tsos2.isValid())
0538     return -1.;
0539 
0540   AlgebraicVector5 v(tsos1.localParameters().vector() - tsos2.localParameters().vector());
0541   AlgebraicSymMatrix55 m(tsos1.localError().matrix() + tsos2.localError().matrix());
0542 
0543   //LogTrace(category) << "match_Chi2 vector v " << v;
0544 
0545   bool ierr = !m.Invert();
0546 
0547   if (ierr) {
0548     edm::LogInfo(category) << "Error inverting covariance matrix";
0549     return -1;
0550   }
0551 
0552   double est = ROOT::Math::Similarity(v, m);
0553 
0554   //LogTrace(category) << "Chi2 " << est;
0555 
0556   return est;
0557 }
0558 
0559 //
0560 // calculate Delta_R of two track candidates at the IP
0561 //
0562 double GlobalMuonTrackMatcher::match_R_IP(const TrackCand& staCand, const TrackCand& tkCand) const {
0563   double dR = 99.0;
0564   if (tkCand.second.isNonnull()) {
0565     dR = (deltaR<double>(staCand.second->eta(), staCand.second->phi(), tkCand.second->eta(), tkCand.second->phi()));
0566   }
0567 
0568   return dR;
0569 }
0570 
0571 //
0572 // calculate Delta_R of two trajectory states
0573 //
0574 double GlobalMuonTrackMatcher::match_Rmom(const TrajectoryStateOnSurface& sta,
0575                                           const TrajectoryStateOnSurface& tk) const {
0576   if (!sta.isValid() || !tk.isValid())
0577     return -1;
0578   return (deltaR<double>(
0579       sta.globalMomentum().eta(), sta.globalMomentum().phi(), tk.globalMomentum().eta(), tk.globalMomentum().phi()));
0580 }
0581 
0582 //
0583 // calculate Delta_R of two trajectory states
0584 //
0585 double GlobalMuonTrackMatcher::match_Rpos(const TrajectoryStateOnSurface& sta,
0586                                           const TrajectoryStateOnSurface& tk) const {
0587   if (!sta.isValid() || !tk.isValid())
0588     return -1;
0589   return (deltaR<double>(
0590       sta.globalPosition().eta(), sta.globalPosition().phi(), tk.globalPosition().eta(), tk.globalPosition().phi()));
0591 }
0592 
0593 //
0594 // calculate the distance in global position of two trajectory states
0595 //
0596 double GlobalMuonTrackMatcher::match_D(const TrajectoryStateOnSurface& sta, const TrajectoryStateOnSurface& tk) const {
0597   if (!sta.isValid() || !tk.isValid())
0598     return -1;
0599   return (sta.globalPosition() - tk.globalPosition()).mag();
0600 }
0601 
0602 //
0603 // calculate the distance in local position of two trajectory states
0604 //
0605 double GlobalMuonTrackMatcher::match_d(const TrajectoryStateOnSurface& sta, const TrajectoryStateOnSurface& tk) const {
0606   if (!sta.isValid() || !tk.isValid())
0607     return -1;
0608   return (sta.localPosition() - tk.localPosition()).mag();
0609 }
0610 
0611 //
0612 // calculate the chi2 of the distance in local position of two
0613 // trajectory states including local errors
0614 //
0615 double GlobalMuonTrackMatcher::match_dist(const TrajectoryStateOnSurface& sta,
0616                                           const TrajectoryStateOnSurface& tk) const {
0617   const string category = "GlobalMuonTrackMatcher";
0618 
0619   if (!sta.isValid() || !tk.isValid())
0620     return -1;
0621 
0622   AlgebraicMatrix22 m;
0623   m(0, 0) = tk.localError().positionError().xx() + sta.localError().positionError().xx();
0624   m(1, 0) = m(0, 1) = tk.localError().positionError().xy() + sta.localError().positionError().xy();
0625   m(1, 1) = tk.localError().positionError().yy() + sta.localError().positionError().yy();
0626 
0627   AlgebraicVector2 v;
0628   v[0] = tk.localPosition().x() - sta.localPosition().x();
0629   v[1] = tk.localPosition().y() - sta.localPosition().y();
0630 
0631   if (!m.Invert()) {
0632     LogTrace(category) << "Error inverting local matrix ";
0633     return -1;
0634   }
0635 
0636   return ROOT::Math::Similarity(v, m);
0637 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!