00001 #include "RecoTracker/TrackProducer/plugins/GsfTrackProducer.h" 00002 // system include files 00003 #include <memory> 00004 // user include files 00005 #include "FWCore/Framework/interface/Frameworkfwd.h" 00006 #include "FWCore/Framework/interface/MakerMacros.h" 00007 #include "FWCore/MessageLogger/interface/MessageLogger.h" 00008 #include "TrackingTools/GeomPropagators/interface/Propagator.h" 00009 #include "TrackingTools/PatternTools/interface/Trajectory.h" 00010 00011 #include "TrackingTools/GsfTracking/interface/TrajGsfTrackAssociation.h" 00012 00013 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" 00014 #include "DataFormats/GsfTrackReco/interface/GsfTrackFwd.h" 00015 #include "DataFormats/GsfTrackReco/interface/GsfTrackExtra.h" 00016 #include "DataFormats/GsfTrackReco/interface/GsfComponent5D.h" 00017 00018 GsfTrackProducer::GsfTrackProducer(const edm::ParameterSet& iConfig): 00019 GsfTrackProducerBase(iConfig.getParameter<bool>("TrajectoryInEvent"), 00020 iConfig.getParameter<bool>("useHitsSplitting")), 00021 theAlgo(iConfig) 00022 { 00023 setConf(iConfig); 00024 setSrc( iConfig.getParameter<edm::InputTag>( "src" ), iConfig.getParameter<edm::InputTag>( "beamSpot" )); 00025 setAlias( iConfig.getParameter<std::string>( "@module_label" ) ); 00026 // string a = alias_; 00027 // a.erase(a.size()-6,a.size()); 00028 //register your products 00029 produces<reco::GsfTrackCollection>().setBranchAlias( alias_ + "GsfTracks" ); 00030 produces<reco::TrackExtraCollection>().setBranchAlias( alias_ + "TrackExtras" ); 00031 produces<reco::GsfTrackExtraCollection>().setBranchAlias( alias_ + "GsfTrackExtras" ); 00032 produces<TrackingRecHitCollection>().setBranchAlias( alias_ + "RecHits" ); 00033 produces<std::vector<Trajectory> >() ; 00034 produces<TrajGsfTrackAssociationCollection>(); 00035 00036 } 00037 00038 00039 void GsfTrackProducer::produce(edm::Event& theEvent, const edm::EventSetup& setup) 00040 { 00041 edm::LogInfo("GsfTrackProducer") << "Analyzing event number: " << theEvent.id() << "\n"; 00042 // 00043 // create empty output collections 00044 // 00045 std::auto_ptr<TrackingRecHitCollection> outputRHColl (new TrackingRecHitCollection); 00046 std::auto_ptr<reco::GsfTrackCollection> outputTColl(new reco::GsfTrackCollection); 00047 std::auto_ptr<reco::TrackExtraCollection> outputTEColl(new reco::TrackExtraCollection); 00048 std::auto_ptr<reco::GsfTrackExtraCollection> outputGsfTEColl(new reco::GsfTrackExtraCollection); 00049 std::auto_ptr<std::vector<Trajectory> > outputTrajectoryColl(new std::vector<Trajectory>); 00050 00051 // 00052 //declare and get stuff to be retrieved from ES 00053 // 00054 edm::ESHandle<TrackerGeometry> theG; 00055 edm::ESHandle<MagneticField> theMF; 00056 edm::ESHandle<TrajectoryFitter> theFitter; 00057 edm::ESHandle<Propagator> thePropagator; 00058 edm::ESHandle<MeasurementTracker> theMeasTk; 00059 edm::ESHandle<TransientTrackingRecHitBuilder> theBuilder; 00060 getFromES(setup,theG,theMF,theFitter,thePropagator,theMeasTk,theBuilder); 00061 00062 // 00063 //declare and get TrackColection to be retrieved from the event 00064 // 00065 AlgoProductCollection algoResults; 00066 reco::BeamSpot bs; 00067 try{ 00068 edm::Handle<TrackCandidateCollection> theTCCollection; 00069 getFromEvt(theEvent,theTCCollection,bs); 00070 00071 // 00072 //run the algorithm 00073 // 00074 LogDebug("GsfTrackProducer") << "run the algorithm" << "\n"; 00075 theAlgo.runWithCandidate(theG.product(), theMF.product(), *theTCCollection, 00076 theFitter.product(), thePropagator.product(), theBuilder.product(), bs, algoResults); 00077 } catch (cms::Exception &e){ edm::LogInfo("GsfTrackProducer") << "cms::Exception caught!!!" << "\n" << e << "\n"; throw; } 00078 // 00079 //put everything in the event 00080 putInEvt(theEvent, thePropagator.product(), theMeasTk.product(), outputRHColl, outputTColl, outputTEColl, outputGsfTEColl, 00081 outputTrajectoryColl, algoResults, bs); 00082 LogDebug("GsfTrackProducer") << "end" << "\n"; 00083 } 00084 00085 00086 // std::vector<reco::TransientTrack> GsfTrackProducer::getTransient(edm::Event& theEvent, const edm::EventSetup& setup) 00087 // { 00088 // edm::LogInfo("GsfTrackProducer") << "Analyzing event number: " << theEvent.id() << "\n"; 00089 // // 00090 // // create empty output collections 00091 // // 00092 // std::vector<reco::TransientTrack> ttks; 00093 00094 // // 00095 // //declare and get stuff to be retrieved from ES 00096 // // 00097 // edm::ESHandle<TrackerGeometry> theG; 00098 // edm::ESHandle<MagneticField> theMF; 00099 // edm::ESHandle<TrajectoryFitter> theFitter; 00100 // edm::ESHandle<Propagator> thePropagator; 00101 // edm::ESHandle<TransientTrackingRecHitBuilder> theBuilder; 00102 // getFromES(setup,theG,theMF,theFitter,thePropagator,theBuilder); 00103 00104 // // 00105 // //declare and get TrackColection to be retrieved from the event 00106 // // 00107 // AlgoProductCollection algoResults; 00108 // try{ 00109 // edm::Handle<TrackCandidateCollection> theTCCollection; 00110 // getFromEvt(theEvent,theTCCollection); 00111 00112 // // 00113 // //run the algorithm 00114 // // 00115 // LogDebug("GsfTrackProducer") << "run the algorithm" << "\n"; 00116 // theAlgo.runWithCandidate(theG.product(), theMF.product(), *theTCCollection, 00117 // theFitter.product(), thePropagator.product(), theBuilder.product(), algoResults); 00118 // } catch (cms::Exception &e){ edm::LogInfo("GsfTrackProducer") << "cms::Exception caught!!!" << "\n" << e << "\n"; throw;} 00119 00120 00121 // for (AlgoProductCollection::iterator prod=algoResults.begin();prod!=algoResults.end(); prod++){ 00122 // ttks.push_back( reco::TransientTrack(*((*prod).second),thePropagator.product()->magneticField() )); 00123 // } 00124 00125 // LogDebug("GsfTrackProducer") << "end" << "\n"; 00126 00127 // return ttks; 00128 // } 00129 00130 00131 // void 00132 // GsfTrackProducer::putInEvt(edm::Event& evt, 00133 // std::auto_ptr<TrackingRecHitCollection>& selHits, 00134 // std::auto_ptr<reco::GsfTrackCollection>& selTracks, 00135 // std::auto_ptr<reco::TrackExtraCollection>& selTrackExtras, 00136 // std::auto_ptr<reco::GsfTrackExtraCollection>& selGsfTrackExtras, 00137 // std::auto_ptr<std::vector<Trajectory> >& selTrajectories, 00138 // AlgoProductCollection& algoResults) 00139 // { 00140 00141 // TrackingRecHitRefProd rHits = evt.getRefBeforePut<TrackingRecHitCollection>(); 00142 // reco::TrackExtraRefProd rTrackExtras = evt.getRefBeforePut<reco::TrackExtraCollection>(); 00143 // reco::GsfTrackExtraRefProd rGsfTrackExtras = evt.getRefBeforePut<reco::GsfTrackExtraCollection>(); 00144 // reco::GsfTrackRefProd rTracks = evt.getRefBeforePut<reco::GsfTrackCollection>(); 00145 00146 // edm::Ref<reco::TrackExtraCollection>::key_type idx = 0; 00147 // edm::Ref<reco::TrackExtraCollection>::key_type hidx = 0; 00148 // edm::Ref<reco::GsfTrackExtraCollection>::key_type idxGsf = 0; 00149 // for(AlgoProductCollection::iterator i=algoResults.begin(); i!=algoResults.end();i++){ 00150 // Trajectory * theTraj = (*i).first; 00151 // if(trajectoryInEvent_) selTrajectories->push_back(*theTraj); 00152 // const TrajectoryFitter::RecHitContainer& transHits = theTraj->recHits(); 00153 00154 // reco::GsfTrack * theTrack = (*i).second.first; 00155 // PropagationDirection seedDir = (*i).second.second; 00156 // // if( ) { 00157 // reco::GsfTrack t = * theTrack; 00158 // selTracks->push_back( t ); 00159 00160 // //sets the outermost and innermost TSOSs 00161 // TrajectoryStateOnSurface outertsos; 00162 // TrajectoryStateOnSurface innertsos; 00163 // unsigned int innerId, outerId; 00164 // if (theTraj->direction() == alongMomentum) { 00165 // outertsos = theTraj->lastMeasurement().updatedState(); 00166 // innertsos = theTraj->firstMeasurement().updatedState(); 00167 // outerId = theTraj->lastMeasurement().recHit()->geographicalId().rawId(); 00168 // innerId = theTraj->firstMeasurement().recHit()->geographicalId().rawId(); 00169 // } else { 00170 // outertsos = theTraj->firstMeasurement().updatedState(); 00171 // innertsos = theTraj->lastMeasurement().updatedState(); 00172 // outerId = theTraj->firstMeasurement().recHit()->geographicalId().rawId(); 00173 // innerId = theTraj->lastMeasurement().recHit()->geographicalId().rawId(); 00174 // } 00175 00176 // GlobalPoint v = outertsos.globalParameters().position(); 00177 // GlobalVector p = outertsos.globalParameters().momentum(); 00178 // math::XYZVector outmom( p.x(), p.y(), p.z() ); 00179 // math::XYZPoint outpos( v.x(), v.y(), v.z() ); 00180 // v = innertsos.globalParameters().position(); 00181 // p = innertsos.globalParameters().momentum(); 00182 // math::XYZVector inmom( p.x(), p.y(), p.z() ); 00183 // math::XYZPoint inpos( v.x(), v.y(), v.z() ); 00184 00185 // reco::TrackExtraRef teref= reco::TrackExtraRef ( rTrackExtras, idx ++ ); 00186 // reco::GsfTrack & track = selTracks->back(); 00187 // track.setExtra( teref ); 00188 // selTrackExtras->push_back( reco::TrackExtra (outpos, outmom, true, 00189 // inpos, inmom, true, 00190 // outertsos.curvilinearError(), outerId, 00191 // innertsos.curvilinearError(), innerId, 00192 // seedDir,theTraj->seedRef())); 00193 00194 // reco::TrackExtra & tx = selTrackExtras->back(); 00195 // size_t i = 0; 00196 // for( TrajectoryFitter::RecHitContainer::const_iterator j = transHits.begin(); 00197 // j != transHits.end(); j ++ ) { 00198 // TrackingRecHit * hit = (**j).hit()->clone(); 00199 // track.setHitPattern( * hit, i ++ ); 00200 // selHits->push_back( hit ); 00201 // tx.add( TrackingRecHitRef( rHits, hidx ++ ) ); 00202 // } 00203 00204 // //build the GsfTrackExtra 00205 // std::vector<reco::GsfComponent5D> outerStates; 00206 // outerStates.reserve(outertsos.components().size()); 00207 // fillStates(outertsos,outerStates); 00208 // std::vector<reco::GsfComponent5D> innerStates; 00209 // innerStates.reserve(innertsos.components().size()); 00210 // fillStates(innertsos,innerStates); 00211 00212 // reco::GsfTrackExtraRef terefGsf = reco::GsfTrackExtraRef ( rGsfTrackExtras, idxGsf ++ ); 00213 // track.setGsfExtra( terefGsf ); 00214 // selGsfTrackExtras->push_back( reco::GsfTrackExtra (outerStates, outertsos.localParameters().pzSign(), 00215 // innerStates, innertsos.localParameters().pzSign())); 00216 00217 // delete theTrack; 00218 // delete theTraj; 00219 // } 00220 00221 // evt.put( selTracks ); 00222 // evt.put( selTrackExtras ); 00223 // evt.put( selGsfTrackExtras ); 00224 // evt.put( selHits ); 00225 // if(trajectoryInEvent_) evt.put(selTrajectories); 00226 // } 00227 00228 00229 // void 00230 // GsfTrackProducer::fillStates (TrajectoryStateOnSurface tsos, 00231 // std::vector<reco::GsfComponent5D>& states) const 00232 // { 00233 // // std::cout << "in fill states" << std::endl; 00234 // // if ( !tsos.isValid() ) { 00235 // // std::cout << std::endl << std::endl << "invalid tsos" << std::endl; 00236 // // return; 00237 // // } 00238 // reco::GsfComponent5D::ParameterVector pLocS; 00239 // reco::GsfComponent5D::CovarianceMatrix cLocS; 00240 // std::vector<TrajectoryStateOnSurface> components(tsos.components()); 00241 // for ( std::vector<TrajectoryStateOnSurface>::const_iterator i=components.begin(); 00242 // i!=components.end(); ++i ) { 00243 // // if ( !(*i).isValid() ) { 00244 // // std::cout << std::endl << "invalid component" << std::endl; 00245 // // continue; 00246 // // } 00247 // // Unneeded hack ... now we have SMatrix in tracking too 00248 // // const AlgebraicVector& pLoc = i->localParameters().vector(); 00249 // // for ( int j=0; j<reco::GsfTrackExtra::dimension; ++j ) pLocS(j) = pLoc[j]; 00250 // // const AlgebraicSymMatrix& cLoc = i->localError().matrix(); 00251 // // for ( int j1=0; j1<reco::GsfTrack::dimension; ++j1 ) 00252 // // for ( int j2=0; j2<=j1; ++j2 ) cLocS(j1,j2) = cLoc[j1][j2]; 00253 // // states.push_back(reco::GsfComponent5D(i->weight(),pLocS,cLocS)); 00254 00255 // states.push_back(reco::GsfComponent5D(i->weight(),i->localParameters().vector(),i->localError().matrix())); 00256 // } 00257 // // std::cout << "end fill states" << std::endl; 00258 // }