Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-02-20 03:46:00

0001 // -*- C++ -*-
0002 //
0003 // Package:    Validation/SiTrackerPhase2V
0004 // Class:      Phase2OTValidateTTStub
0005 
0006 /**
0007  * This class is part of the Phase 2 Tracker validation framework. It validates the
0008  * association of tracking particles to stubs and evaluates stub reconstruction performance
0009  * by generating detailed histograms, including residuals for key parameters.
0010  * 
0011  * Usage:
0012  * To generate histograms from this code, run the test configuration files provided
0013  * in the DQM/SiTrackerPhase2/test directory. The generated histograms can then be
0014  * analyzed or visualized to study stub performance.
0015  */
0016 
0017 // Original Author:
0018 
0019 // Updated by: Brandi Skipworth, 2025
0020 
0021 // system include files
0022 #include <memory>
0023 #include <numeric>
0024 #include <vector>
0025 
0026 // user include files
0027 #include "SimDataFormats/Associations/interface/TTStubAssociationMap.h"
0028 
0029 #include "L1Trigger/TrackFindingTracklet/interface/Settings.h"
0030 
0031 #include "FWCore/Framework/interface/ESHandle.h"
0032 #include "FWCore/Framework/interface/Event.h"
0033 #include "FWCore/Framework/interface/Frameworkfwd.h"
0034 #include "FWCore/Framework/interface/MakerMacros.h"
0035 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0036 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0037 #include "FWCore/ServiceRegistry/interface/Service.h"
0038 #include "FWCore/Utilities/interface/EDGetToken.h"
0039 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0040 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0041 
0042 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0043 #include "DataFormats/L1TrackTrigger/interface/TTTypes.h"
0044 #include "DataFormats/L1TrackTrigger/interface/TTStub.h"
0045 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0046 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0047 
0048 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0049 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0050 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0051 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
0052 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0053 #include "Geometry/CommonTopologies/interface/PixelGeomDetUnit.h"
0054 #include "Geometry/CommonTopologies/interface/PixelTopology.h"
0055 
0056 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0057 #include "DQMServices/Core/interface/DQMStore.h"
0058 
0059 class Phase2OTValidateTTStub : public DQMEDAnalyzer {
0060 public:
0061   explicit Phase2OTValidateTTStub(const edm::ParameterSet&);
0062   ~Phase2OTValidateTTStub() override;
0063   void analyze(const edm::Event&, const edm::EventSetup&) override;
0064   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
0065   void dqmBeginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) override;
0066   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0067   float phiOverBendCorrection(bool isBarrel,
0068                               float stub_z,
0069                               float stub_r,
0070                               const TrackerTopology* tTopo,
0071                               uint32_t detid,
0072                               const GeomDetUnit* det0,
0073                               const GeomDetUnit* det1);
0074   std::vector<double> getTPDerivedCoords(edm::Ptr<TrackingParticle> associatedTP,
0075                                          bool isBarrel,
0076                                          float stub_z,
0077                                          float stub_r) const;
0078   // TTStub stacks
0079   // Global position of the stubs
0080   MonitorElement* Stub_RZ = nullptr;  // TTStub #rho vs. z
0081 
0082   // delta_z hists
0083   MonitorElement* z_res_isPS_barrel = nullptr;
0084   MonitorElement* z_res_is2S_barrel = nullptr;
0085 
0086   // delta_r hists
0087   MonitorElement* r_res_isPS_fw_endcap = nullptr;
0088   MonitorElement* r_res_is2S_fw_endcap = nullptr;
0089   MonitorElement* r_res_isPS_bw_endcap = nullptr;
0090   MonitorElement* r_res_is2S_bw_endcap = nullptr;
0091 
0092   // delta_phi hists
0093   MonitorElement* phi_res_isPS_barrel = nullptr;
0094   MonitorElement* phi_res_is2S_barrel = nullptr;
0095   MonitorElement* phi_res_fw_endcap = nullptr;
0096   MonitorElement* phi_res_bw_endcap = nullptr;
0097   std::vector<MonitorElement*> phi_res_barrel_layers;
0098   std::vector<MonitorElement*> phi_res_fw_endcap_discs;
0099   std::vector<MonitorElement*> phi_res_bw_endcap_discs;
0100 
0101   // delta_bend hists
0102   MonitorElement* bend_res_fw_endcap = nullptr;
0103   MonitorElement* bend_res_bw_endcap = nullptr;
0104   std::vector<MonitorElement*> bend_res_barrel_layers;
0105   std::vector<MonitorElement*> bend_res_fw_endcap_discs;
0106   std::vector<MonitorElement*> bend_res_bw_endcap_discs;
0107 
0108   std::vector<MonitorElement*>* phi_res_vec = nullptr;
0109   std::vector<MonitorElement*>* bend_res_vec = nullptr;
0110 
0111 private:
0112   edm::ParameterSet conf_;
0113   edm::EDGetTokenT<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>> tagTTStubsToken_;
0114   edm::EDGetTokenT<TTStubAssociationMap<Ref_Phase2TrackerDigi_>> ttStubMCTruthToken_;
0115   std::string topFolderName_;
0116   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geomToken_;
0117   const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> topoToken_;
0118   const TrackerGeometry* tkGeom_ = nullptr;
0119   const TrackerTopology* tTopo_ = nullptr;
0120   double TP_minPt;
0121   double TP_maxEta;
0122   double TP_maxVtxZ;
0123   double TP_maxD0;
0124   double TP_maxDxy;
0125 };
0126 
0127 // constructors and destructor
0128 Phase2OTValidateTTStub::Phase2OTValidateTTStub(const edm::ParameterSet& iConfig)
0129     : conf_(iConfig),
0130       geomToken_(esConsumes<TrackerGeometry, TrackerDigiGeometryRecord, edm::Transition::BeginRun>()),
0131       topoToken_(esConsumes<TrackerTopology, TrackerTopologyRcd, edm::Transition::BeginRun>()) {
0132   // now do what ever initialization is needed
0133   topFolderName_ = conf_.getParameter<std::string>("TopFolderName");
0134   tagTTStubsToken_ =
0135       consumes<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>>(conf_.getParameter<edm::InputTag>("TTStubs"));
0136   ttStubMCTruthToken_ =
0137       consumes<TTStubAssociationMap<Ref_Phase2TrackerDigi_>>(conf_.getParameter<edm::InputTag>("MCTruthStubInputTag"));
0138   TP_minPt = conf_.getParameter<double>("TP_minPt");
0139   TP_maxEta = conf_.getParameter<double>("TP_maxEta");
0140   TP_maxVtxZ = conf_.getParameter<double>("TP_maxVtxZ");
0141   TP_maxD0 = conf_.getParameter<double>("TP_maxD0");
0142   TP_maxDxy = conf_.getParameter<double>("TP_maxDxy");
0143 }
0144 
0145 Phase2OTValidateTTStub::~Phase2OTValidateTTStub() {
0146   // do anything here that needs to be done at desctruction time
0147   // (e.g. close files, deallocate resources etc.)
0148 }
0149 
0150 void Phase2OTValidateTTStub::dqmBeginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
0151   tkGeom_ = &(iSetup.getData(geomToken_));
0152   tTopo_ = &(iSetup.getData(topoToken_));
0153 
0154   // Clear existing histograms
0155   phi_res_barrel_layers.clear();
0156   bend_res_barrel_layers.clear();
0157   phi_res_fw_endcap_discs.clear();
0158   bend_res_fw_endcap_discs.clear();
0159   phi_res_bw_endcap_discs.clear();
0160   bend_res_bw_endcap_discs.clear();
0161 
0162   // Resize vectors and set elements to nullptr
0163   phi_res_barrel_layers.resize(trklet::N_LAYER, nullptr);
0164   bend_res_barrel_layers.resize(trklet::N_LAYER, nullptr);
0165   phi_res_fw_endcap_discs.resize(trklet::N_DISK, nullptr);
0166   bend_res_fw_endcap_discs.resize(trklet::N_DISK, nullptr);
0167   phi_res_bw_endcap_discs.resize(trklet::N_DISK, nullptr);
0168   bend_res_bw_endcap_discs.resize(trklet::N_DISK, nullptr);
0169 }
0170 // member functions
0171 
0172 // Calculate the correction factor for tilted modules in barrel
0173 float Phase2OTValidateTTStub::phiOverBendCorrection(bool isBarrel,
0174                                                     float stub_z,
0175                                                     float stub_r,
0176                                                     const TrackerTopology* tTopo,
0177                                                     uint32_t detid,
0178                                                     const GeomDetUnit* det0,
0179                                                     const GeomDetUnit* det1) {
0180   // Get R0, R1, Z0, Z1 values
0181   float R0 = det0->position().perp();
0182   float R1 = det1->position().perp();
0183   float Z0 = det0->position().z();
0184   float Z1 = det1->position().z();
0185 
0186   bool isTiltedBarrel = (isBarrel && tTopo->tobSide(detid) != 3);
0187 
0188   float tiltAngle = 0;  // Initialize to 0 (meaning no tilt, in the endcaps)
0189   if (isTiltedBarrel) {
0190     float deltaR = std::abs(R1 - R0);
0191     float deltaZ = (R1 - R0 > 0)
0192                        ? (Z1 - Z0)
0193                        : -(Z1 - Z0);  // if module parallel, tilt angle should be π/2 and deltaZ would approach zero
0194     // fill histograms here
0195     tiltAngle = atan(deltaR / std::abs(deltaZ));
0196   }
0197 
0198   float correction;
0199   if (isBarrel && tTopo->tobSide(detid) != 3) {  // Assuming this condition represents tiltedBarrel
0200     correction = cos(tiltAngle) * std::abs(stub_z) / stub_r + sin(tiltAngle);
0201   } else if (isBarrel) {
0202     correction = 1;
0203   } else {
0204     correction =
0205         std::abs(stub_z) /
0206         stub_r;  // if tiltAngle = 0, stub (not module) is parallel to the beam line, if tiltAngle = 90, stub is perpendicular to beamline
0207   }
0208 
0209   return correction;
0210 }
0211 
0212 // Compute derived coordinates (z, phi, r) for tracking particle (TP)
0213 std::vector<double> Phase2OTValidateTTStub::getTPDerivedCoords(edm::Ptr<TrackingParticle> associatedTP,
0214                                                                bool isBarrel,
0215                                                                float stub_z,
0216                                                                float stub_r) const {
0217   double tp_phi = -99;
0218   double tp_r = -99;
0219   double tp_z = -99;
0220 
0221   trklet::Settings settings;
0222   double bfield_ = settings.bfield();
0223   double c_ = settings.c();
0224 
0225   // Get values from the tracking particle associatedTP
0226   double tp_pt = associatedTP->pt();
0227   double tp_charge = associatedTP->charge();
0228   float tp_z0 = associatedTP->vertex().z();
0229   double tp_t = associatedTP->tanl();
0230   double tp_rinv = (tp_charge * bfield_) / (tp_pt);
0231 
0232   if (isBarrel) {
0233     tp_r = stub_r;
0234     tp_phi = associatedTP->p4().phi() - std::asin(tp_r * tp_rinv * c_ / 2.0E2);
0235     tp_phi = reco::reducePhiRange(tp_phi);
0236     tp_z = tp_z0 + (2.0E2 / c_) * tp_t * (1 / tp_rinv) * std::asin(tp_r * tp_rinv * c_ / 2.0E2);
0237   } else {
0238     tp_z = stub_z;
0239     tp_phi = associatedTP->p4().phi() - (tp_z - tp_z0) * tp_rinv * c_ / 2.0E2 / tp_t;
0240     tp_phi = reco::reducePhiRange(tp_phi);
0241     tp_r = 2.0E2 / tp_rinv / c_ * std::sin((tp_z - tp_z0) * tp_rinv * c_ / 2.0E2 / tp_t);
0242   }
0243 
0244   std::vector<double> tpDerived_coords{tp_z, tp_phi, tp_r};
0245   return tpDerived_coords;
0246 }
0247 
0248 // ------------ method called for each event  ------------
0249 void Phase2OTValidateTTStub::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0250   /// Track Trigger Stubs
0251   edm::Handle<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>> Phase2TrackerDigiTTStubHandle;
0252   iEvent.getByToken(tagTTStubsToken_, Phase2TrackerDigiTTStubHandle);
0253   edm::Handle<TTStubAssociationMap<Ref_Phase2TrackerDigi_>> MCTruthTTStubHandle;
0254   iEvent.getByToken(ttStubMCTruthToken_, MCTruthTTStubHandle);
0255 
0256   trklet::Settings settings;
0257   double bfield_ = settings.bfield();
0258   double c_ = settings.c();
0259 
0260   /// Loop over input Stubs for basic histogram filling (e.g., Stub_RZ)
0261   typename edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>::const_iterator inputIter;
0262   typename edmNew::DetSet<TTStub<Ref_Phase2TrackerDigi_>>::const_iterator contentIter;
0263   // Adding protection
0264   if (!Phase2TrackerDigiTTStubHandle.isValid() || !MCTruthTTStubHandle.isValid()) {
0265     edm::LogError("Phase2OTValidateTTStub") << "Invalid handle(s) detected.";
0266     return;
0267   }
0268 
0269   for (inputIter = Phase2TrackerDigiTTStubHandle->begin(); inputIter != Phase2TrackerDigiTTStubHandle->end();
0270        ++inputIter) {
0271     for (contentIter = inputIter->begin(); contentIter != inputIter->end(); ++contentIter) {
0272       /// Make reference stub
0273       edm::Ref<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>, TTStub<Ref_Phase2TrackerDigi_>> tempStubRef =
0274           edmNew::makeRefTo(Phase2TrackerDigiTTStubHandle, contentIter);
0275 
0276       /// Get det ID (place of the stub)
0277       //  tempStubRef->getDetId() gives the stackDetId, not rawId
0278       DetId detIdStub = tkGeom_->idToDet((tempStubRef->clusterRef(0))->getDetId())->geographicalId();
0279 
0280       /// Get trigger displacement/offset
0281       //double rawBend = tempStubRef->rawBend();
0282       //double bendOffset = tempStubRef->bendOffset();
0283 
0284       /// Define position stub by position inner cluster
0285       MeasurementPoint mp = (tempStubRef->clusterRef(0))->findAverageLocalCoordinates();
0286       const GeomDet* theGeomDet = tkGeom_->idToDet(detIdStub);
0287       Global3DPoint posStub = theGeomDet->surface().toGlobal(theGeomDet->topology().localPosition(mp));
0288 
0289       Stub_RZ->Fill(posStub.z(), posStub.perp());
0290     }
0291   }
0292 
0293   // Loop over geometric detectors
0294   for (auto gd = tkGeom_->dets().begin(); gd != tkGeom_->dets().end(); gd++) {
0295     DetId detid = (*gd)->geographicalId();
0296 
0297     // Check if detid belongs to TOB or TID subdetectors
0298     if (detid.subdetId() != StripSubdetector::TOB && detid.subdetId() != StripSubdetector::TID)
0299       continue;
0300 
0301     // Process only the lower part of the stack
0302     if (!tTopo_->isLower(detid))
0303       continue;
0304 
0305     // Get the stack DetId
0306     DetId stackDetid = tTopo_->stack(detid);
0307 
0308     // Check if the stackDetid exists in TTStubHandle
0309     if (Phase2TrackerDigiTTStubHandle->find(stackDetid) == Phase2TrackerDigiTTStubHandle->end())
0310       continue;
0311 
0312     // Get the DetSets of the Clusters
0313     edmNew::DetSet<TTStub<Ref_Phase2TrackerDigi_>> stubs = (*Phase2TrackerDigiTTStubHandle)[stackDetid];
0314 
0315     // Calculate detector module positions
0316     const GeomDetUnit* det0 = tkGeom_->idToDetUnit(detid);
0317     const GeomDetUnit* det1 = tkGeom_->idToDetUnit(tTopo_->partnerDetId(detid));
0318     if (!det0 || !det1) {
0319       edm::LogError("Phase2OTValidateTTStub") << "Error: det0 or det1 is null";
0320       continue;
0321     }
0322     float modMinR = std::min(det0->position().perp(), det1->position().perp());
0323     float modMaxR = std::max(det0->position().perp(), det1->position().perp());
0324     float modMinZ = std::min(det0->position().z(), det1->position().z());
0325     float modMaxZ = std::max(det0->position().z(), det1->position().z());
0326     float sensorSpacing = sqrt((modMaxR - modMinR) * (modMaxR - modMinR) + (modMaxZ - modMinZ) * (modMaxZ - modMinZ));
0327 
0328     // Calculate strip pitch
0329     const PixelGeomDetUnit* theGeomDetUnit = dynamic_cast<const PixelGeomDetUnit*>(det0);
0330     if (!theGeomDetUnit) {
0331       edm::LogError("Phase2OTValidateTTStub") << "Error: theGeomDetUnit is null";
0332       continue;
0333     }
0334     const PixelTopology& topo = theGeomDetUnit->specificTopology();
0335     float stripPitch = topo.pitch().first;
0336 
0337     // Loop over input stubs
0338     for (auto stubIter = stubs.begin(); stubIter != stubs.end(); ++stubIter) {
0339       // Create reference to the stub
0340       edm::Ref<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>, TTStub<Ref_Phase2TrackerDigi_>> tempStubPtr =
0341           edmNew::makeRefTo(Phase2TrackerDigiTTStubHandle, stubIter);
0342 
0343       // Check if the stub is genuine
0344       if (!MCTruthTTStubHandle->isGenuine(tempStubPtr))
0345         continue;
0346 
0347       // Get det ID from the stub
0348       DetId detIdStub = tkGeom_->idToDet((tempStubPtr->clusterRef(0))->getDetId())->geographicalId();
0349 
0350       // Retrieve geometrical detector
0351       const GeomDet* theGeomDet = tkGeom_->idToDet(detIdStub);
0352       if (!theGeomDet) {
0353         edm::LogError("Phase2OTValidateTTStub") << "Error: theGeomDet is null";
0354         continue;
0355       }
0356 
0357       // Retrieve tracking particle associated with TTStub
0358       edm::Ptr<TrackingParticle> associatedTP = MCTruthTTStubHandle->findTrackingParticlePtr(tempStubPtr);
0359       if (associatedTP.isNull())
0360         continue;
0361 
0362       // Determine layer and subdetector information
0363       int isBarrel = 0;
0364       int layer = -999999;
0365       if (detid.subdetId() == StripSubdetector::TOB) {
0366         isBarrel = 1;
0367         layer = static_cast<int>(tTopo_->layer(detid));
0368       } else if (detid.subdetId() == StripSubdetector::TID) {
0369         isBarrel = 0;
0370         layer = static_cast<int>(tTopo_->layer(detid));
0371       } else {
0372         edm::LogVerbatim("Tracklet") << "WARNING -- neither TOB nor TID stub, shouldn't happen...";
0373         layer = -1;
0374       }
0375 
0376       bool isPSmodule = tkGeom_->getDetectorType(detid) == TrackerGeometry::ModuleType::Ph2PSP;
0377 
0378       // Calculate local coordinates of clusters
0379       MeasurementPoint innerClusterCoords = tempStubPtr->clusterRef(0)->findAverageLocalCoordinatesCentered();
0380       MeasurementPoint outerClusterCoords = tempStubPtr->clusterRef(1)->findAverageLocalCoordinatesCentered();
0381 
0382       // Convert local coordinates to global positions
0383       Global3DPoint innerClusterGlobalPos =
0384           theGeomDet->surface().toGlobal(theGeomDet->topology().localPosition(innerClusterCoords));
0385       Global3DPoint outerClusterGlobalPos =
0386           theGeomDet->surface().toGlobal(theGeomDet->topology().localPosition(outerClusterCoords));
0387 
0388       // Determine maximum Z positions of the stubs
0389       float stub_maxZ = std::max(innerClusterGlobalPos.z(), outerClusterGlobalPos.z());
0390 
0391       // Stub parameters
0392       float stub_phi = innerClusterGlobalPos.phi();
0393       float stub_r = innerClusterGlobalPos.perp();
0394       float stub_z = innerClusterGlobalPos.z();
0395 
0396       // Tracking particle parameters
0397       int tp_charge = associatedTP->charge();
0398       float tp_pt = associatedTP->p4().pt();
0399       float tp_eta = associatedTP->p4().eta();
0400       float tp_d0 = associatedTP->d0();
0401       float tp_vx = associatedTP->vx();
0402       float tp_vy = associatedTP->vy();
0403       float tp_vz = associatedTP->vz();
0404       float tp_dxy = std::sqrt(tp_vx * tp_vx + tp_vy * tp_vy);
0405 
0406       if (tp_charge == 0)
0407         continue;
0408       if (tp_pt < TP_minPt)
0409         continue;
0410       if (std::abs(tp_eta) > TP_maxEta)
0411         continue;
0412       if (std::abs(tp_vz) > TP_maxVtxZ)
0413         continue;
0414       if (std::abs(tp_d0) > TP_maxD0)
0415         continue;
0416       if (std::abs(tp_dxy) > TP_maxDxy)
0417         continue;
0418 
0419       // Derived coordinates
0420       std::vector<double> tpDerivedCoords = getTPDerivedCoords(associatedTP, isBarrel, stub_z, stub_r);
0421       float tp_z = tpDerivedCoords[0];
0422       float tp_phi = tpDerivedCoords[1];
0423       float tp_r = tpDerivedCoords[2];
0424 
0425       // Trigger information
0426       float trigBend = tempStubPtr->bendFE();
0427       if (!isBarrel && stub_maxZ < 0.0) {
0428         trigBend = -trigBend;
0429       }
0430 
0431       float correctionValue = phiOverBendCorrection(isBarrel, stub_z, stub_r, tTopo_, detid, det0, det1);
0432       float trackBend =
0433           -(sensorSpacing * stub_r * bfield_ * c_ * tp_charge) / (stripPitch * 2.0E2 * tp_pt * correctionValue);
0434 
0435       float bendRes = trackBend - trigBend;
0436       float zRes = tp_z - stub_z;
0437       float phiRes = tp_phi - stub_phi;
0438       float rRes = tp_r - stub_r;
0439 
0440       // Histograms for z_res, phi_res, and r_res based on module type and location
0441       if (isBarrel == 1) {
0442         if (isPSmodule) {
0443           z_res_isPS_barrel->Fill(zRes);
0444           phi_res_isPS_barrel->Fill(phiRes);
0445         } else {
0446           z_res_is2S_barrel->Fill(zRes);
0447           phi_res_is2S_barrel->Fill(phiRes);
0448         }
0449       } else {
0450         if (stub_maxZ > 0) {
0451           if (isPSmodule) {
0452             r_res_isPS_fw_endcap->Fill(rRes);
0453           } else {
0454             r_res_is2S_fw_endcap->Fill(rRes);
0455           }
0456         } else {
0457           if (isPSmodule) {
0458             r_res_isPS_bw_endcap->Fill(rRes);
0459           } else {
0460             r_res_is2S_bw_endcap->Fill(rRes);
0461           }
0462         }
0463       }
0464 
0465       // Ensure that the vectors are correctly assigned before use
0466       if (isBarrel == 1) {
0467         phi_res_vec = &phi_res_barrel_layers;
0468         bend_res_vec = &bend_res_barrel_layers;
0469       } else {
0470         if (stub_maxZ > 0) {
0471           // Forward endcap
0472           bend_res_fw_endcap->Fill(bendRes);
0473           phi_res_fw_endcap->Fill(phiRes);
0474           phi_res_vec = &phi_res_fw_endcap_discs;
0475           bend_res_vec = &bend_res_fw_endcap_discs;
0476         } else {
0477           // Backward endcap
0478           bend_res_bw_endcap->Fill(bendRes);
0479           phi_res_bw_endcap->Fill(phiRes);
0480           phi_res_vec = &phi_res_bw_endcap_discs;
0481           bend_res_vec = &bend_res_bw_endcap_discs;
0482         }
0483       }
0484       // Fill the appropriate histogram based on layer/disc
0485       if (layer >= 1 && layer <= trklet::N_LAYER) {
0486         (*bend_res_vec)[layer - 1]->Fill(bendRes);
0487         (*phi_res_vec)[layer - 1]->Fill(phiRes);
0488       }
0489     }
0490   }
0491 }  // end of method
0492 
0493 // ------------ method called when starting to processes a run  ------------
0494 void Phase2OTValidateTTStub::bookHistograms(DQMStore::IBooker& iBooker,
0495                                             edm::Run const& run,
0496                                             edm::EventSetup const& es) {
0497   edm::ParameterSet psTTStub_RZ = conf_.getParameter<edm::ParameterSet>("TH2TTStub_RZ");
0498   edm::ParameterSet ps_2S_Res = conf_.getParameter<edm::ParameterSet>("TH1_2S_Res");
0499   edm::ParameterSet ps_PS_Res = conf_.getParameter<edm::ParameterSet>("TH1_PS_Res");
0500   edm::ParameterSet psPhi_Res = conf_.getParameter<edm::ParameterSet>("TH1Phi_Res");
0501   edm::ParameterSet psBend_Res = conf_.getParameter<edm::ParameterSet>("TH1Bend_Res");
0502   std::string HistoName;
0503   iBooker.setCurrentFolder(topFolderName_);
0504   // 2D histogram for stub_RZ
0505   HistoName = "Stub_RZ";
0506   Stub_RZ = iBooker.book2D(HistoName,
0507                            HistoName,
0508                            psTTStub_RZ.getParameter<int32_t>("Nbinsx"),
0509                            psTTStub_RZ.getParameter<double>("xmin"),
0510                            psTTStub_RZ.getParameter<double>("xmax"),
0511                            psTTStub_RZ.getParameter<int32_t>("Nbinsy"),
0512                            psTTStub_RZ.getParameter<double>("ymin"),
0513                            psTTStub_RZ.getParameter<double>("ymax"));
0514 
0515   iBooker.setCurrentFolder(topFolderName_ + "/Residual");
0516   // z-res for PS modules
0517   HistoName = "#Delta z Barrel PS modules";
0518   z_res_isPS_barrel = iBooker.book1D(HistoName,
0519                                      HistoName,
0520                                      ps_PS_Res.getParameter<int32_t>("Nbinsx"),
0521                                      ps_PS_Res.getParameter<double>("xmin"),
0522                                      ps_PS_Res.getParameter<double>("xmax"));
0523   z_res_isPS_barrel->setAxisTitle("tp_z - stub_z", 1);
0524   z_res_isPS_barrel->setAxisTitle("events ", 2);
0525 
0526   // z-res for 2S modules
0527   HistoName = "#Delta z Barrel 2S modules";
0528   z_res_is2S_barrel = iBooker.book1D(HistoName,
0529                                      HistoName,
0530                                      ps_2S_Res.getParameter<int32_t>("Nbinsx"),
0531                                      ps_2S_Res.getParameter<double>("xmin"),
0532                                      ps_2S_Res.getParameter<double>("xmax"));
0533   z_res_is2S_barrel->setAxisTitle("tp_z - stub_z [cm]", 1);
0534   z_res_is2S_barrel->setAxisTitle("events ", 2);
0535 
0536   // r-res for fw endcap PS modules
0537   HistoName = "#Delta r FW Endcap PS modules";
0538   r_res_isPS_fw_endcap = iBooker.book1D(HistoName,
0539                                         HistoName,
0540                                         ps_PS_Res.getParameter<int32_t>("Nbinsx"),
0541                                         ps_PS_Res.getParameter<double>("xmin"),
0542                                         ps_PS_Res.getParameter<double>("xmax"));
0543   r_res_isPS_fw_endcap->setAxisTitle("tp_r - stub_r [cm]", 1);
0544   r_res_isPS_fw_endcap->setAxisTitle("events ", 2);
0545 
0546   // r-res for fw endcap 2S modules
0547   HistoName = "#Delta r FW Endcap 2S modules";
0548   r_res_is2S_fw_endcap = iBooker.book1D(HistoName,
0549                                         HistoName,
0550                                         ps_2S_Res.getParameter<int32_t>("Nbinsx"),
0551                                         ps_2S_Res.getParameter<double>("xmin"),
0552                                         ps_2S_Res.getParameter<double>("xmax"));
0553   r_res_is2S_fw_endcap->setAxisTitle("tp_r - stub_r [cm]", 1);
0554   r_res_is2S_fw_endcap->setAxisTitle("events ", 2);
0555 
0556   // r-res for bw endcap PS modules
0557   HistoName = "#Delta r BW Endcap PS modules";
0558   r_res_isPS_bw_endcap = iBooker.book1D(HistoName,
0559                                         HistoName,
0560                                         ps_PS_Res.getParameter<int32_t>("Nbinsx"),
0561                                         ps_PS_Res.getParameter<double>("xmin"),
0562                                         ps_PS_Res.getParameter<double>("xmax"));
0563   r_res_isPS_bw_endcap->setAxisTitle("tp_r - stub_r [cm]", 1);
0564   r_res_isPS_bw_endcap->setAxisTitle("events ", 2);
0565 
0566   // r-res for bw endcap 2S modules
0567   HistoName = "#Delta r BW Endcap 2S modules";
0568   r_res_is2S_bw_endcap = iBooker.book1D(HistoName,
0569                                         HistoName,
0570                                         ps_2S_Res.getParameter<int32_t>("Nbinsx"),
0571                                         ps_2S_Res.getParameter<double>("xmin"),
0572                                         ps_2S_Res.getParameter<double>("xmax"));
0573   r_res_is2S_bw_endcap->setAxisTitle("tp_r - stub_r [cm]", 1);
0574   r_res_is2S_bw_endcap->setAxisTitle("events ", 2);
0575 
0576   // histograms for phi_res and bend_res
0577   HistoName = "#Delta #phi Barrel PS modules";
0578   phi_res_isPS_barrel = iBooker.book1D(HistoName,
0579                                        HistoName,
0580                                        psPhi_Res.getParameter<int32_t>("Nbinsx"),
0581                                        psPhi_Res.getParameter<double>("xmin"),
0582                                        psPhi_Res.getParameter<double>("xmax"));
0583   phi_res_isPS_barrel->setAxisTitle("tp_phi - stub_phi", 1);
0584   phi_res_isPS_barrel->setAxisTitle("events", 2);
0585 
0586   HistoName = "#Delta #phi Barrel 2S modules";
0587   phi_res_is2S_barrel = iBooker.book1D(HistoName,
0588                                        HistoName,
0589                                        psPhi_Res.getParameter<int32_t>("Nbinsx"),
0590                                        psPhi_Res.getParameter<double>("xmin"),
0591                                        psPhi_Res.getParameter<double>("xmax"));
0592 
0593   HistoName = "#Delta #phi FW Endcap";
0594   phi_res_fw_endcap = iBooker.book1D(HistoName,
0595                                      HistoName,
0596                                      psPhi_Res.getParameter<int32_t>("Nbinsx"),
0597                                      psPhi_Res.getParameter<double>("xmin"),
0598                                      psPhi_Res.getParameter<double>("xmax"));
0599 
0600   HistoName = "#Delta #phi BW Endcap";
0601   phi_res_bw_endcap = iBooker.book1D(HistoName,
0602                                      HistoName,
0603                                      psPhi_Res.getParameter<int32_t>("Nbinsx"),
0604                                      psPhi_Res.getParameter<double>("xmin"),
0605                                      psPhi_Res.getParameter<double>("xmax"));
0606 
0607   HistoName = "#Delta bend FW Endcap";
0608   bend_res_fw_endcap = iBooker.book1D(HistoName,
0609                                       HistoName,
0610                                       psBend_Res.getParameter<int32_t>("Nbinsx"),
0611                                       psBend_Res.getParameter<double>("xmin"),
0612                                       psBend_Res.getParameter<double>("xmax"));
0613 
0614   HistoName = "#Delta bend BW Endcap";
0615   bend_res_bw_endcap = iBooker.book1D(HistoName,
0616                                       HistoName,
0617                                       psBend_Res.getParameter<int32_t>("Nbinsx"),
0618                                       psBend_Res.getParameter<double>("xmin"),
0619                                       psBend_Res.getParameter<double>("xmax"));
0620 
0621   // barrel layers
0622   for (int i = 0; i < trklet::N_LAYER; ++i) {
0623     std::string HistoName = "#Delta #phi Barrel L" + std::to_string(i + 1);
0624     phi_res_barrel_layers[i] = iBooker.book1D(HistoName,
0625                                               HistoName,
0626                                               psPhi_Res.getParameter<int32_t>("Nbinsx"),
0627                                               psPhi_Res.getParameter<double>("xmin"),
0628                                               psPhi_Res.getParameter<double>("xmax"));
0629     phi_res_barrel_layers[i]->setAxisTitle("tp_phi - stub_phi", 1);
0630     phi_res_barrel_layers[i]->setAxisTitle("events", 2);
0631 
0632     HistoName = "#Delta bend Barrel L" + std::to_string(i + 1);
0633     bend_res_barrel_layers[i] = iBooker.book1D(HistoName,
0634                                                HistoName,
0635                                                psBend_Res.getParameter<int32_t>("Nbinsx"),
0636                                                psBend_Res.getParameter<double>("xmin"),
0637                                                psBend_Res.getParameter<double>("xmax"));
0638     bend_res_barrel_layers[i]->setAxisTitle("tp_bend - stub_bend", 1);
0639     bend_res_barrel_layers[i]->setAxisTitle("events", 2);
0640   }
0641 
0642   // endcap discs
0643   for (int i = 0; i < trklet::N_DISK; ++i) {
0644     std::string HistoName = "#Delta #phi FW Endcap D" + std::to_string(i + 1);
0645     phi_res_fw_endcap_discs[i] = iBooker.book1D(HistoName,
0646                                                 HistoName,
0647                                                 psPhi_Res.getParameter<int32_t>("Nbinsx"),
0648                                                 psPhi_Res.getParameter<double>("xmin"),
0649                                                 psPhi_Res.getParameter<double>("xmax"));
0650     phi_res_fw_endcap_discs[i]->setAxisTitle("tp_phi - stub_phi", 1);
0651     phi_res_fw_endcap_discs[i]->setAxisTitle("events", 2);
0652 
0653     HistoName = "#Delta bend FW Endcap D" + std::to_string(i + 1);
0654     bend_res_fw_endcap_discs[i] = iBooker.book1D(HistoName,
0655                                                  HistoName,
0656                                                  psBend_Res.getParameter<int32_t>("Nbinsx"),
0657                                                  psBend_Res.getParameter<double>("xmin"),
0658                                                  psBend_Res.getParameter<double>("xmax"));
0659     bend_res_fw_endcap_discs[i]->setAxisTitle("tp_bend - stub_bend", 1);
0660     bend_res_fw_endcap_discs[i]->setAxisTitle("events", 2);
0661 
0662     HistoName = "#Delta #phi BW Endcap D" + std::to_string(i + 1);
0663     phi_res_bw_endcap_discs[i] = iBooker.book1D(HistoName,
0664                                                 HistoName,
0665                                                 psPhi_Res.getParameter<int32_t>("Nbinsx"),
0666                                                 psPhi_Res.getParameter<double>("xmin"),
0667                                                 psPhi_Res.getParameter<double>("xmax"));
0668     phi_res_bw_endcap_discs[i]->setAxisTitle("tp_phi - stub_phi", 1);
0669     phi_res_bw_endcap_discs[i]->setAxisTitle("events", 2);
0670 
0671     HistoName = "#Delta bend BW Endcap D" + std::to_string(i + 1);
0672     bend_res_bw_endcap_discs[i] = iBooker.book1D(HistoName,
0673                                                  HistoName,
0674                                                  psBend_Res.getParameter<int32_t>("Nbinsx"),
0675                                                  psBend_Res.getParameter<double>("xmin"),
0676                                                  psBend_Res.getParameter<double>("xmax"));
0677     bend_res_bw_endcap_discs[i]->setAxisTitle("tp_bend - stub_bend", 1);
0678     bend_res_bw_endcap_discs[i]->setAxisTitle("events", 2);
0679   }
0680 }
0681 
0682 void Phase2OTValidateTTStub::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0683   // Phase2OTValidateTTStub
0684   edm::ParameterSetDescription desc;
0685   {
0686     edm::ParameterSetDescription psd0;
0687     psd0.add<int>("Nbinsx", 900);
0688     psd0.add<double>("xmax", 300);
0689     psd0.add<double>("xmin", -300);
0690     psd0.add<int>("Nbinsy", 900);
0691     psd0.add<double>("ymax", 120);
0692     psd0.add<double>("ymin", 0);
0693     desc.add<edm::ParameterSetDescription>("TH2TTStub_RZ", psd0);
0694   }
0695   {
0696     edm::ParameterSetDescription psd0;
0697     psd0.add<int>("Nbinsx", 99);
0698     psd0.add<double>("xmax", 5.5);
0699     psd0.add<double>("xmin", -5.5);
0700     desc.add<edm::ParameterSetDescription>("TH1_2S_Res", psd0);
0701   }
0702   {
0703     edm::ParameterSetDescription psd0;
0704     psd0.add<int>("Nbinsx", 99);
0705     psd0.add<double>("xmax", 2.0);
0706     psd0.add<double>("xmin", -2.0);
0707     desc.add<edm::ParameterSetDescription>("TH1_PS_Res", psd0);
0708   }
0709   {
0710     edm::ParameterSetDescription psd0;
0711     psd0.add<int>("Nbinsx", 599);
0712     psd0.add<double>("xmax", 0.1);
0713     psd0.add<double>("xmin", -0.1);
0714     desc.add<edm::ParameterSetDescription>("TH1Phi_Res", psd0);
0715   }
0716   {
0717     edm::ParameterSetDescription psd0;
0718     psd0.add<int>("Nbinsx", 59);
0719     psd0.add<double>("xmax", 5.0);
0720     psd0.add<double>("xmin", -5.5);
0721     desc.add<edm::ParameterSetDescription>("TH1Bend_Res", psd0);
0722   }
0723 
0724   desc.add<std::string>("TopFolderName", "TrackerPhase2OTStubV");
0725   desc.add<edm::InputTag>("TTStubs", edm::InputTag("TTStubsFromPhase2TrackerDigis", "StubAccepted"));
0726   desc.add<edm::InputTag>("trackingParticleToken", edm::InputTag("mix", "MergedTrackTruth"));
0727   desc.add<edm::InputTag>("MCTruthStubInputTag", edm::InputTag("TTStubAssociatorFromPixelDigis", "StubAccepted"));
0728   desc.add<edm::InputTag>("MCTruthClusterInputTag",
0729                           edm::InputTag("TTClusterAssociatorFromPixelDigis", "ClusterInclusive"));
0730   desc.add<int>("TP_minNStub", 4);
0731   desc.add<int>("TP_minNLayersStub", 4);
0732   desc.add<double>("TP_minPt", 2.0);
0733   desc.add<double>("TP_maxEta", 2.4);
0734   desc.add<double>("TP_maxVtxZ", 15.0);
0735   desc.add<double>("TP_maxD0", 1.0);
0736   desc.add<double>("TP_maxDxy", 1.0);
0737   descriptions.add("Phase2OTValidateTTStub", desc);
0738   // or use the following to generate the label from the module's C++ type
0739   //descriptions.addWithDefaultLabel(desc);
0740 }
0741 DEFINE_FWK_MODULE(Phase2OTValidateTTStub);