Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:33:10

0001 #include <memory>
0002 
0003 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0004 
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/Framework/interface/MakerMacros.h"
0007 
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 
0010 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0011 #include "AnalysisDataFormats/TrackInfo/interface/TrackInfo.h"
0012 #include "AnalysisDataFormats/TrackInfo/interface/TrackInfoTrackAssociation.h"
0013 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
0014 
0015 #include <iostream>
0016 #include <string>
0017 
0018 using namespace edm;
0019 
0020 class TrackInfoAnalyzerExample : public edm::one::EDAnalyzer<> {
0021   edm::EDGetTokenT<reco::TrackInfoTrackAssociationCollection> TItkassociatorCollectionToken;
0022   edm::EDGetTokenT<reco::TrackCollection> tkCollectionToken;
0023 
0024 public:
0025   TrackInfoAnalyzerExample(const edm::ParameterSet& pset);
0026 
0027   ~TrackInfoAnalyzerExample() {}
0028 
0029   virtual void analyze(const edm::Event& event, const edm::EventSetup& setup) {
0030     using namespace reco;
0031 
0032     //get TrackInfoTrackAssociationCollection from the event
0033     edm::Handle<reco::TrackInfoTrackAssociationCollection> TItkassociatorCollection;
0034     event.getByToken(TItkassociatorCollectionToken, TItkassociatorCollection);
0035 
0036     // get track collection from the event
0037     edm::Handle<reco::TrackCollection> tkCollection;
0038     event.getByToken(tkCollectionToken, tkCollection);
0039     edm::LogInfo("TrackInfoAnalyzerExample") << "track info collection size " << TItkassociatorCollection->size();
0040     edm::LogInfo("TrackInfoAnalyzerExample") << "track collection size " << tkCollection->size();
0041 
0042     // loop on the tracks
0043     for (unsigned int track = 0; track < tkCollection->size(); ++track) {
0044       //build the ref to the track
0045       reco::TrackRef trackref = reco::TrackRef(tkCollection, track);
0046       edm::LogInfo("TrackInfoAnalyzerExample") << "Track pt" << trackref->pt();
0047 
0048       //get the ref to the trackinfo
0049       reco::TrackInfoRef trackinforef = (*TItkassociatorCollection.product())[trackref];
0050 
0051       // get additional track information from trackinfo:
0052 
0053       //the seed:
0054       const TrajectorySeed seed = trackinforef->seed();
0055       edm::LogInfo("TrackInfoAnalyzerExample") << "N hits in the seed: " << seed.nHits();
0056       edm::LogInfo("TrackInfoAnalyzerExample")
0057           << "Starting state position" << seed.startingState().parameters().position();
0058       edm::LogInfo("TrackInfoAnalyzerExample")
0059           << "Starting state direction" << seed.startingState().parameters().momentum();
0060 
0061       //local angle for a specific hit
0062       TrackingRecHitRef rechitref = trackref->recHit(2);
0063       if (rechitref->isValid()) {
0064         const LocalVector localdir = trackinforef->localTrackMomentum(Combined, rechitref);
0065         edm::LogInfo("TrackInfoAnalyzerExample")
0066             << "Local x-z plane angle of 3rd hit:" << atan2(localdir.x(), localdir.z());
0067       }
0068 
0069       // loop on all the track hits
0070       reco::TrackInfo::TrajectoryInfo::const_iterator iter;
0071       for (iter = trackinforef->trajStateMap().begin(); iter != trackinforef->trajStateMap().end(); ++iter) {
0072         //trajectory local direction and position on detector
0073         LocalVector statedirection = (trackinforef->stateOnDet(Combined, (*iter).first)->parameters()).momentum();
0074         LocalPoint stateposition = (trackinforef->stateOnDet(Combined, (*iter).first)->parameters()).position();
0075         edm::LogInfo("TrackInfoAnalyzerExample") << "LocalMomentum: " << statedirection;
0076         edm::LogInfo("TrackInfoAnalyzerExample") << "LocalPosition: " << stateposition;
0077         edm::LogInfo("TrackInfoAnalyzerExample")
0078             << "Local x-z plane angle: " << atan2(statedirection.x(), statedirection.z());
0079         if (trackinforef->type((*iter).first) == Matched) {  // get the direction for the components
0080           edm::LogInfo("TrackInfoAnalyzerExample")
0081               << "LocalMomentum (mono): " << trackinforef->localTrackMomentumOnMono(Combined, (*iter).first);
0082           edm::LogInfo("TrackInfoAnalyzerExample")
0083               << "LocalMomentum (stereo): " << trackinforef->localTrackMomentumOnStereo(Combined, (*iter).first);
0084         } else if (trackinforef->type((*iter).first) == Projected) {  //one should be 0
0085           edm::LogInfo("TrackInfoAnalyzerExample")
0086               << "LocalMomentum (mono): " << trackinforef->localTrackMomentumOnMono(Combined, (*iter).first);
0087           edm::LogInfo("TrackInfoAnalyzerExample")
0088               << "LocalMomentum (stereo): " << trackinforef->localTrackMomentumOnStereo(Combined, (*iter).first);
0089         }
0090         //hit position on detector
0091         if (((*iter).first)->isValid())
0092           edm::LogInfo("TrackInfoAnalyzerExample") << "LocalPosition (rechit): " << ((*iter).first)->localPosition();
0093       }
0094     }
0095   }
0096 };
0097 
0098 TrackInfoAnalyzerExample::TrackInfoAnalyzerExample(const edm::ParameterSet& pset)
0099     : TItkassociatorCollectionToken(
0100           consumes<reco::TrackInfoTrackAssociationCollection>(pset.getParameter<edm::InputTag>("TrackInfo"))),
0101       tkCollectionToken(consumes<reco::TrackCollection>(pset.getParameter<edm::InputTag>("Tracks"))) {}
0102 
0103 DEFINE_FWK_MODULE(TrackInfoAnalyzerExample);