File indexing completed on 2023-03-17 11:13:43
0001 #include "L1Trigger/TrackFindingTMTT/interface/StubFEWindows.h"
0002 #include "L1Trigger/TrackFindingTMTT/interface/TrackerModule.h"
0003 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0004 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0005
0006 #include <algorithm>
0007 #include <utility>
0008
0009 using namespace std;
0010
0011 namespace tmtt {
0012
0013
0014
0015 StubFEWindows::StubFEWindows(const edm::ParameterSet& pSetStubAlgo) {
0016 numTiltedLayerRings_ = pSetStubAlgo.getParameter<vector<double>>("NTiltedRings");
0017 windowSizeBarrelLayers_ = pSetStubAlgo.getParameter<vector<double>>("BarrelCut");
0018 const auto& pSetTiltedLayer = pSetStubAlgo.getParameter<vector<edm::ParameterSet>>("TiltedBarrelCutSet");
0019 const auto& pSetEncapDisks = pSetStubAlgo.getParameter<vector<edm::ParameterSet>>("EndcapCutSet");
0020 windowSizeTiltedLayersRings_.reserve(pSetTiltedLayer.size());
0021 for (const auto& pSet : pSetTiltedLayer) {
0022 windowSizeTiltedLayersRings_.emplace_back(pSet.getParameter<vector<double>>("TiltedCut"));
0023 }
0024 windowSizeEndcapDisksRings_.reserve(pSetEncapDisks.size());
0025 for (const auto& pSet : pSetEncapDisks) {
0026 windowSizeEndcapDisksRings_.emplace_back(pSet.getParameter<vector<double>>("EndcapCut"));
0027 }
0028 }
0029
0030
0031
0032 void StubFEWindows::setZero() {
0033 std::fill(windowSizeBarrelLayers_.begin(), windowSizeBarrelLayers_.end(), 0.);
0034 for (auto& x : windowSizeEndcapDisksRings_)
0035 std::fill(x.begin(), x.end(), 0.);
0036 for (auto& y : windowSizeTiltedLayersRings_)
0037 std::fill(y.begin(), y.end(), 0.);
0038 }
0039
0040
0041
0042 const double* StubFEWindows::storedWindowSize(const TrackerTopology* trackerTopo, const DetId& detId) const {
0043
0044
0045 const double* storedHalfWindow = nullptr;
0046 if (detId.subdetId() == StripSubdetector::TOB) {
0047 unsigned int layer = trackerTopo->layer(detId);
0048 unsigned int ladder = trackerTopo->tobRod(detId);
0049 int type = 2 * trackerTopo->tobSide(detId) - 3;
0050 double corr = 0;
0051
0052 if (type != TrackerModule::BarrelModuleType::flat) {
0053
0054 corr = (numTiltedLayerRings_.at(layer) + 1) / 2.;
0055
0056 ladder = corr - (corr - ladder) * type;
0057 storedHalfWindow = &(windowSizeTiltedLayersRings_.at(layer).at(ladder));
0058 } else {
0059
0060 storedHalfWindow = &(windowSizeBarrelLayers_.at(layer));
0061 }
0062
0063 } else if (detId.subdetId() == StripSubdetector::TID) {
0064
0065 unsigned int wheel = trackerTopo->tidWheel(detId);
0066 unsigned int ring = trackerTopo->tidRing(detId);
0067 storedHalfWindow = &(windowSizeEndcapDisksRings_.at(wheel).at(ring));
0068 }
0069 return storedHalfWindow;
0070 }
0071
0072 double* StubFEWindows::storedWindowSize(const TrackerTopology* trackerTopo, const DetId& detId) {
0073
0074
0075
0076 return const_cast<double*>(std::as_const(*this).storedWindowSize(trackerTopo, detId));
0077 }
0078 }