PhysicsObjectsMonitoring/src/PhysicsObjectsMonitor.cc

0001 /** \class PhysicsObjectsMonitor
0002  *  Analyzer of the StandAlone muon tracks
0003  *
0004  *  \author M. Mulders - CERN <martijn.mulders@cern.ch>
0005  *  Based on STAMuonAnalyzer by R. Bellan - INFN Torino
0006  *<riccardo.bellan@cern.ch>
0007  */
0008
0009 #include "DQM/PhysicsObjectsMonitoring/interface/PhysicsObjectsMonitor.h"
0010
0011 // Collaborating Class Header
0012 #include "FWCore/Framework/interface/ESHandle.h"
0013 #include "FWCore/Framework/interface/Event.h"
0014 #include "FWCore/Framework/interface/MakerMacros.h"
0015
0016 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0017
0018 #include "DataFormats/TrackReco/interface/Track.h"
0019 #include "RecoMuon/TrackingTools/interface/MuonPatternRecoDumper.h"
0020 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0021
0022 #include "DataFormats/MuonDetId/interface/MuonSubdetId.h"
0023
0024 #include "DQMServices/Core/interface/DQMStore.h"
0025 #include <FWCore/MessageLogger/interface/MessageLogger.h>
0026
0027 using namespace std;
0028 using namespace edm;
0029
0030 /// Constructor
0031 PhysicsObjectsMonitor::PhysicsObjectsMonitor(const ParameterSet &pset) {
0032   theSTAMuonLabel = pset.getUntrackedParameter<string>("StandAloneTrackCollectionLabel");
0033   theSeedCollectionLabel = pset.getUntrackedParameter<string>("MuonSeedCollectionLabel");
0034   theDataType = pset.getUntrackedParameter<string>("DataType");
0035   magFiledToken_ = esConsumes();
0036   if (theDataType != "RealData" && theDataType != "SimData")
0037     edm::LogInfo("PhysicsObjectsMonitor") << "Error in Data Type!!" << endl;
0038
0039   if (theDataType == "SimData") {
0040     edm::LogInfo("PhysicsObjectsMonitor") << "Sorry! Running this package on simulation is no longer supported! ";
0041   }
0042
0043   // set Token(-s)
0044   theSTAMuonToken_ =
0045       consumes<reco::TrackCollection>(pset.getUntrackedParameter<string>("StandAloneTrackCollectionLabel"));
0046 }
0047
0048 /// Destructor
0049 PhysicsObjectsMonitor::~PhysicsObjectsMonitor() {}
0050
0051 void PhysicsObjectsMonitor::bookHistograms(DQMStore::IBooker &iBooker, edm::Run const &, edm::EventSetup const &) {
0052   iBooker.setCurrentFolder("PhysicsObjects/MuonReconstruction");
0053
0054   hPres = iBooker.book1D("pTRes", "pT Resolution", 100, -2, 2);
0055   h1_Pres = iBooker.book1D("invPTRes", "1/pT Resolution", 100, -2, 2);
0056
0057   charge = iBooker.book1D("charge", "track charge", 5, -2.5, 2.5);
0058   ptot = iBooker.book1D("ptot", "track momentum", 50, 0, 50);
0059   pt = iBooker.book1D("pt", "track pT", 100, 0, 50);
0060   px = iBooker.book1D("px ", "track px", 100, -50, 50);
0061   py = iBooker.book1D("py", "track py", 100, -50, 50);
0062   pz = iBooker.book1D("pz", "track pz", 100, -50, 50);
0063   Nmuon = iBooker.book1D("Nmuon", "Number of muon tracks", 11, -.5, 10.5);
0064   Nrechits = iBooker.book1D("Nrechits", "Number of RecHits/Segments on track", 21, -.5, 21.5);
0065   NDThits = iBooker.book1D("NDThits", "Number of DT Hits/Segments on track", 31, -.5, 31.5);
0066   NCSChits = iBooker.book1D("NCSChits", "Number of CSC Hits/Segments on track", 31, -.5, 31.5);
0067   NRPChits = iBooker.book1D("NRPChits", "Number of RPC hits on track", 11, -.5, 11.5);
0068
0069   DTvsCSC = iBooker.book2D("DTvsCSC", "Number of DT vs CSC hits on track", 29, -.5, 28.5, 29, -.5, 28.5);
0070   TH2F *root_ob = DTvsCSC->getTH2F();
0071   root_ob->SetXTitle("Number of DT hits");
0072   root_ob->SetYTitle("Number of CSC hits");
0073 }
0074
0075 void PhysicsObjectsMonitor::analyze(const Event &event, const EventSetup &eventSetup) {
0076   edm::LogInfo("PhysicsObjectsMonitor") << "Run: " << event.id().run() << " Event: " << event.id().event();
0077   MuonPatternRecoDumper debug;
0078
0079   // Get the RecTrack collection from the event
0080   Handle<reco::TrackCollection> staTracks;
0081   event.getByToken(theSTAMuonToken_, staTracks);
0082
0083   const auto &theMGField = eventSetup.getHandle(magFiledToken_);
0084
0085   double recPt = 0.;
0086   double simPt = 0.;
0087
0088   reco::TrackCollection::const_iterator staTrack;
0089
0090   edm::LogInfo("PhysicsObjectsMonitor") << "Reconstructed tracks: " << staTracks->size() << endl;
0091   Nmuon->Fill(staTracks->size());
0092   for (staTrack = staTracks->begin(); staTrack != staTracks->end(); ++staTrack) {
0093     reco::TransientTrack track(*staTrack, &*theMGField);
0094
0095     int nrechits = 0;
0096     int nDThits = 0;
0097     int nCSChits = 0;
0098     int nRPChits = 0;
0099
0100     for (trackingRecHit_iterator it = track.recHitsBegin(); it != track.recHitsEnd(); it++) {
0101       if ((*it)->isValid()) {
0102         edm::LogInfo("PhysicsObjectsMonitor") << "Analyzer:  Aha this looks like a Rechit!" << std::endl;
0103         if ((*it)->geographicalId().subdetId() == MuonSubdetId::DT) {
0104           nDThits++;
0105         } else if ((*it)->geographicalId().subdetId() == MuonSubdetId::CSC) {
0106           nCSChits++;
0107         } else if ((*it)->geographicalId().subdetId() == MuonSubdetId::RPC) {
0108           nRPChits++;
0109         } else {
0110           edm::LogInfo("PhysicsObjectsMonitor") << "This is an UNKNOWN hit !! " << std::endl;
0111         }
0112         nrechits++;
0113       }
0114     }
0115
0116     Nrechits->Fill(nrechits);
0117     NDThits->Fill(nDThits);
0118     NCSChits->Fill(nCSChits);
0119     DTvsCSC->Fill(nDThits, nCSChits);
0120     NRPChits->Fill(nRPChits);
0121
0122     debug.dumpFTS(track.impactPointTSCP().theState());
0123
0124     recPt = track.impactPointTSCP().momentum().perp();
0125     edm::LogInfo("PhysicsObjectsMonitor")
0126         << " p: " << track.impactPointTSCP().momentum().mag() << " pT: " << recPt << endl;
0127     pt->Fill(recPt);
0128     ptot->Fill(track.impactPointTSCP().momentum().mag());
0129     charge->Fill(track.impactPointTSCP().charge());
0130     px->Fill(track.impactPointTSCP().momentum().x());
0131     py->Fill(track.impactPointTSCP().momentum().y());
0132     pz->Fill(track.impactPointTSCP().momentum().z());
0133   }
0134   edm::LogInfo("PhysicsObjectsMonitor") << "---" << endl;
0135   if (recPt && theDataType == "SimData") {
0136     hPres->Fill((recPt - simPt) / simPt);
0137     h1_Pres->Fill((1 / recPt - 1 / simPt) / (1 / simPt));
0138   }
0139 }
0140
0141 DEFINE_FWK_MODULE(PhysicsObjectsMonitor);