RecoMuon/src/MuonSeedTrack.cc

0001 /**
0002  *  Class: MuonSeedTrack
0003  *
0004  *
0005  *
0006  *  Authors :
0007  *  \author Adam Everett - Purdue University
0008  *
0009  */
0010
0011 #include "Validation/RecoMuon/src/MuonSeedTrack.h"
0012
0013 #include "FWCore/Framework/interface/ConsumesCollector.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015
0016 #include "DataFormats/TrackReco/interface/Track.h"
0017 #include "DataFormats/TrackReco/interface/TrackExtra.h"
0018
0019 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0020
0021 #include "TrackingTools/DetLayers/interface/DetLayer.h"
0022 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0023 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
0024 #include "TrackingTools/GeomPropagators/interface/TrackerBounds.h"
0025
0026 #include "RecoMuon/TrackingTools/interface/MuonServiceProxy.h"
0027 #include "RecoMuon/TrackingTools/interface/MuonUpdatorAtVertex.h"
0028 #include "RecoMuon/TrackingTools/interface/MuonPatternRecoDumper.h"
0029
0030 using namespace reco;
0031 using namespace edm;
0032 using namespace std;
0033
0034 typedef TrajectoryStateOnSurface TSOS;
0035
0036 //
0037 // constructors
0038 //
0039 MuonSeedTrack::MuonSeedTrack(const edm::ParameterSet& pset) {
0040   // service parameters
0041   ParameterSet serviceParameters = pset.getParameter<ParameterSet>("ServiceParameters");
0042   theService = new MuonServiceProxy(serviceParameters, consumesCollector());
0043
0044   ParameterSet updatorPar = pset.getParameter<ParameterSet>("MuonUpdatorAtVertexParameters");
0045   //theSeedPropagatorName = updatorPar.getParameter<string>("Propagator");
0046   theSeedsLabel = pset.getParameter<InputTag>("MuonSeed");
0047   theSeedsToken = consumes<TrajectorySeedCollection>(theSeedsLabel);
0048   theUpdatorAtVtx = new MuonUpdatorAtVertex(updatorPar, theService);
0049
0050   theAllowNoVtxFlag = pset.getUntrackedParameter<bool>("AllowNoVertex", false);
0051
0052   //register products
0053   setAlias(pset.getParameter<std::string>("@module_label"));
0054   produces<reco::TrackCollection>().setBranchAlias(theAlias + "Tracks");
0055 }
0056
0057 //
0058 // destructor
0059 //
0060 MuonSeedTrack::~MuonSeedTrack() {
0061   if (theService)
0062     delete theService;
0063   if (theUpdatorAtVtx)
0064     delete theUpdatorAtVtx;
0065 }
0066
0067 //
0068 // member functions
0069 //
0070
0071 /*!
0072  * For each seed, make a (fake) reco::Track
0073  */
0074 void MuonSeedTrack::produce(edm::Event& event, const edm::EventSetup& eventSetup) {
0075   using namespace edm;
0076
0077   // Update the services
0078   theService->update(eventSetup);
0079
0080   // the track collectios; they will be loaded in the event
0081   unique_ptr<reco::TrackCollection> trackCollection(new reco::TrackCollection());
0082   // ... and its reference into the event
0083   reco::TrackRefProd trackCollectionRefProd = event.getRefBeforePut<reco::TrackCollection>();
0084
0085   Handle<TrajectorySeedCollection> seeds;
0086   event.getByToken(theSeedsToken, seeds);
0087
0088   for (TrajectorySeedCollection::const_iterator iSeed = seeds->begin(); iSeed != seeds->end(); iSeed++) {
0089     pair<bool, reco::Track> resultOfTrackExtrapAtPCA = buildTrackAtPCA(*iSeed);
0090     if (!resultOfTrackExtrapAtPCA.first)
0091       continue;
0092     // take the "bare" track at PCA
0093     reco::Track& track = resultOfTrackExtrapAtPCA.second;
0094     // fill the TrackCollection
0095     trackCollection->push_back(track);
0096   }
0097
0098   event.put(std::move(trackCollection));
0099 }
0100
0101 /*!
0102  * empty method
0103  */
0104 void MuonSeedTrack::beginJob() {}
0105
0106 /*!
0107  * empty method
0108  */
0109 void MuonSeedTrack::endJob() {}
0110
0111 /*!
0112  * Get the TrajectoryStateOnSurface from the TrajectorySeed
0113  */
0114 TrajectoryStateOnSurface MuonSeedTrack::getSeedTSOS(const TrajectorySeed& seed) const {
0115   // Get the Trajectory State on Det (persistent version of a TSOS) from the seed
0116   PTrajectoryStateOnDet pTSOD = seed.startingState();
0117
0118   // Transform it in a TrajectoryStateOnSurface
0119
0120   DetId seedDetId(pTSOD.detId());
0121
0122   const GeomDet* gdet = theService->trackingGeometry()->idToDet(seedDetId);
0123
0124   TrajectoryStateOnSurface initialState =
0125       trajectoryStateTransform::transientState(pTSOD, &(gdet->surface()), &*theService->magneticField());
0126
0127   /*
0128   // Get the layer on which the seed relies
0129   const DetLayer *initialLayer = theService->detLayerGeometry()->idToLayer( seedDetId );
0130
0131   PropagationDirection detLayerOrder = oppositeToMomentum;
0132
0133   // ask for compatible layers
0134   vector<const DetLayer*> detLayers;
0135   detLayers = initialLayer->compatibleLayers( *initialState.freeState(),detLayerOrder);
0136
0137   TrajectoryStateOnSurface result = initialState;
0138   if(detLayers.size()){
0139     const DetLayer* finalLayer = detLayers.back();
0140     const TrajectoryStateOnSurface propagatedState = theService->propagator(theSeedPropagatorName)->propagate(initialState, finalLayer->surface());
0141     if(propagatedState.isValid())
0142       result = propagatedState;
0143   }
0144
0145   return result;
0146   */
0147
0148   return initialState;
0149 }
0150
0151 /*!
0152  * First calls getSeedTSOS, then propagates to the vertex.  Creates a
0153  * reco::Track from the propagated initial FreeTrajectoryState.
0154  */
0155 pair<bool, reco::Track> MuonSeedTrack::buildTrackAtPCA(const TrajectorySeed& seed) const {
0156   const string metname = "MuonSeedTrack";
0157
0158   MuonPatternRecoDumper debug;
0159
0160   TSOS seedTSOS = getSeedTSOS(seed);
0161   // This is needed to extrapolate the tsos at vertex
0162   LogTrace(metname) << "Propagate to PCA...";
0163   pair<bool, FreeTrajectoryState> extrapolationResult = theUpdatorAtVtx->propagateToNominalLine(seedTSOS);
0164   FreeTrajectoryState ftsAtVtx;
0165
0166   if (extrapolationResult.first) {
0167     ftsAtVtx = extrapolationResult.second;
0168   } else {
0169     if (TrackerBounds::isInside(seedTSOS.globalPosition())) {
0170       LogWarning(metname) << "Track in the Tracker: taking the innermost state instead of the state at PCA";
0171       ftsAtVtx = *seedTSOS.freeState();
0172     } else {
0173       if (theAllowNoVtxFlag) {
0174         LogWarning(metname) << "Propagation to PCA failed, taking the innermost state instead of the state at PCA";
0175         ftsAtVtx = *seedTSOS.freeState();
0176       } else {
0177         LogWarning(metname) << "Stand Alone track: this track will be rejected";
0178         return pair<bool, reco::Track>(false, reco::Track());
0179       }
0180     }
0181   }
0182
0183   LogTrace(metname) << "TSOS after the extrapolation at vtx";
0184   LogTrace(metname) << debug.dumpFTS(ftsAtVtx);
0185
0186   GlobalPoint pca = ftsAtVtx.position();
0187   math::XYZPoint persistentPCA(pca.x(), pca.y(), pca.z());
0188   GlobalVector p = ftsAtVtx.momentum();
0189   math::XYZVector persistentMomentum(p.x(), p.y(), p.z());
0190
0191   double dummyNDOF = 1.0;
0192   double dummyChi2 = 1.0;
0193
0194   reco::Track track(
0195       dummyChi2, dummyNDOF, persistentPCA, persistentMomentum, ftsAtVtx.charge(), ftsAtVtx.curvilinearError());
0196
0197   return pair<bool, reco::Track>(true, track);
0198 }