TrackInfoProducer/test/TrackInfoAnalyzer.cc

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 //#include "DataFormats/DetId/interface/DetId.h"
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     //       cout << "beginJob" <<endl;
0029     //     edm::ESHandle<TrackerGeometry> tkgeom;
0030     //    c.get<TrackerDigiGeometryRecord>().get( tkgeom );
0031     //     tracker=&(* tkgeom);
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) {  //analyze
0046
0047     using namespace reco;
0048
0049     //Retrieve tracker topology from geometry
0050     const TrackerTopology* const tTopo = &setup.getData(tkTopoToken);
0051
0052     //std::cout << "\nEvent ID = "<< event.id() << std::endl ;
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       //const reco::TrackInfo::TrajectoryInfo tinfo=track->trajstate();
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);