Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-05-21 22:39:40

0001 
0002 #ifndef DTCalibValidation_H
0003 #define DTCalibValidation_H
0004 
0005 /** \class DTCalibValidation
0006  *  Analysis on DT residuals to validate the kFactor
0007  *
0008  *
0009  *  \author G. Mila - INFN Torino
0010  */
0011 
0012 #include "FWCore/Framework/interface/Frameworkfwd.h"
0013 #include "DataFormats/MuonDetId/interface/DTSuperLayerId.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/Framework/interface/EDAnalyzer.h"
0016 #include "FWCore/Framework/interface/LuminosityBlock.h"
0017 
0018 #include "DQMServices/Core/interface/DQMStore.h"
0019 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0020 
0021 #include "DataFormats/MuonDetId/interface/DTWireId.h"
0022 #include "DataFormats/DTRecHit/interface/DTRecHitCollection.h"
0023 #include "DataFormats/DTRecHit/interface/DTRecSegment2DCollection.h"
0024 #include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h"
0025 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0026 
0027 #include "FWCore/Framework/interface/ESHandle.h"
0028 
0029 #include <string>
0030 #include <map>
0031 #include <vector>
0032 
0033 // To remove into CMSSW versions before 20X
0034 // To add into CMSSW versions before 20X
0035 //class DaqMonitorBEInterface;
0036 
0037 class DTGeometry;
0038 class DTChamber;
0039 
0040 // FR class DTCalibValidation: public edm::EDAnalyzer{
0041 class DTCalibValidation : public DQMEDAnalyzer {
0042 public:
0043   /// Constructor
0044   DTCalibValidation(const edm::ParameterSet& pset);
0045 
0046   /// Destructor
0047   ~DTCalibValidation() override;
0048 
0049   /// BeginRun
0050   void dqmBeginRun(const edm::Run& r, const edm::EventSetup& c) override;
0051 
0052   // Operations
0053   void analyze(const edm::Event& event, const edm::EventSetup& setup) override;
0054 
0055 protected:
0056   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
0057 
0058 private:
0059   // Switch for verbosity
0060   //bool debug;
0061   edm::ParameterSet parameters;
0062   int wrongSegment;
0063   int rightSegment;
0064   int nevent;
0065   // the analysis type
0066   bool detailedAnalysis;
0067   // the geometry
0068   edm::ESGetToken<DTGeometry, MuonGeometryRecord> muonGeomToken_;
0069   const DTGeometry* dtGeom;
0070 
0071   // Lable of 1D rechits in the event
0072   edm::EDGetTokenT<DTRecHitCollection> recHits1DToken_;
0073   // Lable of 2D segments in the event
0074   edm::EDGetTokenT<DTRecSegment2DCollection> segment2DToken_;
0075   // Lable of 4D segments in the event
0076   edm::EDGetTokenT<DTRecSegment4DCollection> segment4DToken_;
0077 
0078   // Return a map between DTRecHit1DPair and wireId
0079   std::map<DTWireId, std::vector<DTRecHit1DPair> > map1DRecHitsPerWire(const DTRecHitCollection* dt1DRecHitPairs);
0080 
0081   // Return a map between DTRecHit1D and wireId
0082   std::map<DTWireId, std::vector<DTRecHit1D> > map1DRecHitsPerWire(const DTRecSegment2DCollection* segment2Ds);
0083 
0084   // Return a map between DTRecHit1D and wireId
0085   std::map<DTWireId, std::vector<DTRecHit1D> > map1DRecHitsPerWire(const DTRecSegment4DCollection* segment4Ds);
0086 
0087   template <typename type>
0088   const type* findBestRecHit(const DTLayer* layer,
0089                              DTWireId wireId,
0090                              const std::vector<type>& recHits,
0091                              const float simHitDist);
0092 
0093   // Compute the distance from wire (cm) of a hits in a DTRecHit1DPair
0094   float recHitDistFromWire(const DTRecHit1DPair& hitPair, const DTLayer* layer);
0095   // Compute the distance from wire (cm) of a hits in a DTRecHit1D
0096   float recHitDistFromWire(const DTRecHit1D& recHit, const DTLayer* layer);
0097   // Compute the position with respect to the wire (cm) of a hits in a DTRecHit1DPair
0098   float recHitPosition(
0099       const DTRecHit1DPair& hitPair, const DTLayer* layer, const DTChamber* chamber, float segmPos, int sl);
0100   // Compute the position with respect to the wire (cm) of a hits in a DTRecHit1D
0101   float recHitPosition(const DTRecHit1D& recHit, const DTLayer* layer, const DTChamber* chamber, float segmPos, int sl);
0102 
0103   // Does the real job
0104   template <typename type>
0105   void compute(const DTGeometry* dtGeom,
0106                const DTRecSegment4D& segment,
0107                const std::map<DTWireId, std::vector<type> >& recHitsPerWire,
0108                int step);
0109 
0110   // Book a set of histograms for a give chamber
0111   void bookHistos(DTSuperLayerId slId, int step);
0112   // Fill a set of histograms for a give chamber
0113   void fillHistos(DTSuperLayerId slId,
0114                   float distance,
0115                   float residualOnDistance,
0116                   float position,
0117                   float residualOnPosition,
0118                   int step);
0119 
0120   std::map<std::pair<DTSuperLayerId, int>, std::vector<MonitorElement*> > histosPerSL;
0121 };
0122 #endif
0123 
0124 /* Local Variables: */
0125 /* show-trailing-whitespace: t */
0126 /* truncate-lines: t */
0127 /* End: */