File indexing completed on 2023-03-17 10:40:52
0001 #include <memory>
0002
0003 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0004
0005 #include "FWCore/Framework/interface/ESHandle.h"
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "FWCore/Framework/interface/MakerMacros.h"
0008
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010 #include "AnalysisDataFormats/TrackInfo/interface/TrackInfo.h"
0011 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0012 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
0013
0014 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0015 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0016 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
0017 #include <iostream>
0018 #include <string>
0019 #include "TH1F.h"
0020 #include "TFile.h"
0021
0022 using namespace edm;
0023
0024 class TrackInfoAnalyzer : public edm::one::EDAnalyzer<> {
0025 public:
0026 TrackInfoAnalyzer(const edm::ParameterSet& pset);
0027 void beginJob() {
0028
0029
0030
0031
0032 output = new TFile(filename.c_str(), "RECREATE");
0033 tib1int = new TH1F("Tib1Int", "Tib Layer 1 Int", 100, -1, 1);
0034 tib1ext = new TH1F("Tib1Ext", "Tib Layer 1 Ext", 100, -1, 1);
0035 tib2int = new TH1F("Tib2Int", "Tib Layer 2 Int", 100, -1, 1);
0036 tib2ext = new TH1F("Tib2Ext", "Tib Layer 2 Ext", 100, -1, 1);
0037 tib3int = new TH1F("Tib3Int", "Tib Layer 3 Int", 100, -0.5, 0.5);
0038 tib3ext = new TH1F("Tib3Ext", "Tib Layer 3 Ext", 100, -0.5, 0.5);
0039 tib4int = new TH1F("Tib4Int", "Tib Layer 4 Int", 100, -0.5, 0.5);
0040 tib4ext = new TH1F("Tib4Ext", "Tib Layer 4 Ext", 100, -0.5, 0.5);
0041 }
0042
0043 ~TrackInfoAnalyzer() { delete output; }
0044
0045 virtual void analyze(const edm::Event& event, const edm::EventSetup& setup) {
0046
0047 using namespace reco;
0048
0049
0050 const TrackerTopology* const tTopo = &setup.getData(tkTopoToken);
0051
0052
0053 edm::Handle<reco::TrackInfoCollection> trackCollection;
0054 event.getByToken(trackCollectionToken, trackCollection);
0055
0056 reco::TrackInfoCollection tC = *(trackCollection.product());
0057
0058 edm::LogInfo("TrackInfoAnalyzer") << "number of infos " << tC.size();
0059 for (reco::TrackInfoCollection::iterator track = tC.begin(); track != tC.end(); ++track) {
0060
0061 reco::TrackInfo::TrajectoryInfo::const_iterator iter;
0062 edm::LogInfo("TrackInfoAnalyzer") << "N hits in the seed: " << track->seed().nHits();
0063 edm::LogInfo("TrackInfoAnalyzer") << "Starting state " << track->seed().startingState().parameters().position();
0064 if (track->trajStateMap().size() > 0) {
0065 for (iter = track->trajStateMap().begin(); iter != track->trajStateMap().end(); ++iter) {
0066 edm::LogInfo("TrackInfoAnalyzer")
0067 << "LocalMomentum: " << (track->stateOnDet(Combined, (*iter).first)->parameters()).momentum();
0068 edm::LogInfo("TrackInfoAnalyzer")
0069 << "LocalPosition: " << (track->stateOnDet(Combined, (*iter).first)->parameters()).position();
0070 edm::LogInfo("TrackInfoAnalyzer") << "LocalPosition (rechit): " << ((*iter).first)->localPosition();
0071 DetId id = (*iter).first->geographicalId();
0072 unsigned int iSubDet = StripSubdetector(id).subdetId();
0073 if (iSubDet == StripSubdetector::TIB) {
0074 int layer = tTopo->tibLayer(id);
0075 unsigned int order = tTopo->tibOrder(id);
0076 if (layer == 1) {
0077 if (order == 0)
0078 tib1int->Fill((track->stateOnDet(Combined, (*iter).first)->parameters()).position().x() -
0079 ((*iter).first)->localPosition().x());
0080 else if (order == 1)
0081 tib1ext->Fill((track->stateOnDet(Combined, (*iter).first)->parameters()).position().x() -
0082 ((*iter).first)->localPosition().x());
0083 } else if (layer == 2) {
0084 if (order == 0)
0085 tib2int->Fill((track->stateOnDet(Combined, (*iter).first)->parameters()).position().x() -
0086 ((*iter).first)->localPosition().x());
0087 else if (order == 1)
0088 tib2ext->Fill((track->stateOnDet(Combined, (*iter).first)->parameters()).position().x() -
0089 ((*iter).first)->localPosition().x());
0090 } else if (layer == 3) {
0091 if (order == 0)
0092 tib3int->Fill((track->stateOnDet(Combined, (*iter).first)->parameters()).position().x() -
0093 ((*iter).first)->localPosition().x());
0094 else if (order == 1)
0095 tib3ext->Fill((track->stateOnDet(Combined, (*iter).first)->parameters()).position().x() -
0096 ((*iter).first)->localPosition().x());
0097 } else if (layer == 4) {
0098 if (order == 0)
0099 tib4int->Fill((track->stateOnDet(Combined, (*iter).first)->parameters()).position().x() -
0100 ((*iter).first)->localPosition().x());
0101 else if (order == 1)
0102 tib4ext->Fill((track->stateOnDet(Combined, (*iter).first)->parameters()).position().x() -
0103 ((*iter).first)->localPosition().x());
0104 }
0105 }
0106 }
0107 }
0108 }
0109 }
0110 void endJob() {
0111 output->Write();
0112 output->Close();
0113 }
0114
0115 private:
0116 TFile* output;
0117 TH1F* tib1int;
0118 TH1F* tib1ext;
0119 TH1F* tib2int;
0120 TH1F* tib2ext;
0121 TH1F* tib3int;
0122 TH1F* tib3ext;
0123 TH1F* tib4int;
0124 TH1F* tib4ext;
0125 const std::string filename;
0126 const edm::EDGetTokenT<reco::TrackInfoCollection> trackCollectionToken;
0127 const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> tkTopoToken;
0128 };
0129
0130 TrackInfoAnalyzer::TrackInfoAnalyzer(const edm::ParameterSet& pset)
0131 : filename(pset.getParameter<std::string>("OutFileName")),
0132 trackCollectionToken(consumes<reco::TrackInfoCollection>(pset.getParameter<edm::InputTag>("TrackInfo"))),
0133 tkTopoToken(esConsumes()) {}
0134
0135 DEFINE_FWK_MODULE(TrackInfoAnalyzer);