Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 11:59:58

0001 // -*- C++ -*-
0002 //
0003 // Package:    CalibTracker/SiStripLorentzAnglePCLMonitor
0004 // Class:      SiStripLorentzAnglePCLMonitor
0005 //
0006 /**\class SiStripLorentzAnglePCLMonitor SiStripLorentzAnglePCLMonitor.cc CalibTracker/SiStripLorentzAnglePCLMonitor/plugins/SiStripLorentzAnglePCLMonitor.cc
0007 
0008  Description: class to book and fill histograms necessary for the online monitoring of the SiStripLorentzAngle
0009 
0010  Implementation:
0011      Largely taken from https://github.com/robervalwalsh/tracker-la/blob/master/SiStripLAMonitor.cc
0012 */
0013 //
0014 // Original Author:  musich
0015 //         Created:  Sun, 07 May 2023 16:57:10 GMT
0016 //
0017 //
0018 
0019 // system includes
0020 #include <string>
0021 #include <fmt/format.h>
0022 #include <fmt/printf.h>
0023 
0024 // user include files
0025 #include "CalibTracker/Records/interface/SiStripDependentRecords.h"
0026 #include "CalibTracker/SiStripCommon/interface/ShallowTools.h"
0027 #include "CalibTracker/SiStripLorentzAngle/interface/SiStripLorentzAngleCalibrationHelpers.h"
0028 #include "CalibTracker/SiStripLorentzAngle/interface/SiStripLorentzAngleCalibrationStruct.h"
0029 #include "CondFormats/SiStripObjects/interface/SiStripLorentzAngle.h"
0030 #include "CondFormats/SiStripObjects/interface/SiStripLatency.h"
0031 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0032 #include "DQMServices/Core/interface/MonitorElement.h"
0033 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
0034 #include "DataFormats/SiStripDetId/interface/SiStripDetId.h"
0035 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0036 #include "DataFormats/TrackReco/interface/Track.h"
0037 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
0038 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
0039 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
0040 #include "FWCore/Framework/interface/Event.h"
0041 #include "FWCore/Framework/interface/Frameworkfwd.h"
0042 #include "FWCore/Framework/interface/MakerMacros.h"
0043 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0044 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0045 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
0046 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0047 #include "MagneticField/Engine/interface/MagneticField.h"
0048 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0049 #include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterInfo.h"
0050 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
0051 
0052 // ROOT includes
0053 #include "TVector3.h"
0054 
0055 //
0056 // class declaration
0057 //
0058 
0059 class SiStripLorentzAnglePCLMonitor : public DQMEDAnalyzer {
0060 public:
0061   explicit SiStripLorentzAnglePCLMonitor(const edm::ParameterSet&);
0062   ~SiStripLorentzAnglePCLMonitor() override = default;
0063 
0064   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0065 
0066 private:
0067   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
0068   void analyze(const edm::Event&, const edm::EventSetup&) override;
0069   void dqmBeginRun(edm::Run const&, edm::EventSetup const&) override;
0070 
0071   std::string moduleLocationType(const uint32_t& mod, const TrackerTopology* tTopo);
0072 
0073   // ------------ member data ------------
0074   SiStripClusterInfo m_clusterInfo;
0075   SiStripLorentzAngleCalibrationHistograms iHists_;
0076 
0077   // for magnetic field conversion
0078   static constexpr float teslaToInverseGeV_ = 2.99792458e-3f;
0079 
0080   bool mismatchedBField_;
0081   bool mismatchedLatency_;
0082   const std::string folder_;
0083   const bool saveHistosMods_;
0084   const edm::EDGetTokenT<edm::View<reco::Track>> m_tracks_token;
0085   const edm::EDGetTokenT<TrajTrackAssociationCollection> m_association_token;
0086   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> m_tkGeomToken;
0087 
0088   const edm::ESGetToken<SiStripLatency, SiStripLatencyRcd> m_latencyTokenBR;
0089   const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> m_topoEsTokenBR;
0090   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> m_tkGeomTokenBR;
0091   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> m_magFieldTokenBR;
0092   const edm::ESGetToken<SiStripLorentzAngle, SiStripLorentzAngleDepRcd> m_lorentzAngleTokenBR;
0093 
0094   struct OnTrackCluster {
0095     uint32_t det;
0096     const SiStripCluster* cluster;
0097     const Trajectory* traj;
0098     const reco::Track* track;
0099     const TrajectoryMeasurement& measurement;
0100     OnTrackCluster(uint32_t detId,
0101                    const SiStripCluster* stripCluster,
0102                    const Trajectory* trajectory,
0103                    const reco::Track* track_,
0104                    const TrajectoryMeasurement& measurement_)
0105         : det{detId}, cluster{stripCluster}, traj{trajectory}, track{track_}, measurement{measurement_} {}
0106   };
0107 };
0108 
0109 SiStripLorentzAnglePCLMonitor::SiStripLorentzAnglePCLMonitor(const edm::ParameterSet& iConfig)
0110     : m_clusterInfo(consumesCollector()),
0111       mismatchedBField_{false},
0112       mismatchedLatency_{false},
0113       folder_(iConfig.getParameter<std::string>("folder")),
0114       saveHistosMods_(iConfig.getParameter<bool>("saveHistoMods")),
0115       m_tracks_token(consumes<edm::View<reco::Track>>(iConfig.getParameter<edm::InputTag>("Tracks"))),
0116       m_association_token(consumes<TrajTrackAssociationCollection>(iConfig.getParameter<edm::InputTag>("Tracks"))),
0117       m_tkGeomToken{esConsumes<>()},
0118       m_latencyTokenBR{esConsumes<edm::Transition::BeginRun>()},
0119       m_topoEsTokenBR{esConsumes<edm::Transition::BeginRun>()},
0120       m_tkGeomTokenBR{esConsumes<edm::Transition::BeginRun>()},
0121       m_magFieldTokenBR{esConsumes<edm::Transition::BeginRun>()},
0122       m_lorentzAngleTokenBR{esConsumes<edm::Transition::BeginRun>()} {}
0123 //
0124 // member functions
0125 //
0126 
0127 void SiStripLorentzAnglePCLMonitor::dqmBeginRun(edm::Run const& run, edm::EventSetup const& iSetup) {
0128   const auto& tkGeom = iSetup.getData(m_tkGeomTokenBR);
0129   const auto& magField = iSetup.getData(m_magFieldTokenBR);
0130   const auto& lorentzAngle = iSetup.getData(m_lorentzAngleTokenBR);
0131   const TrackerTopology* tTopo = &iSetup.getData(m_topoEsTokenBR);
0132 
0133   // fast cachecd access
0134   const auto& theMagField = 1.f / (magField.inverseBzAtOriginInGeV() * teslaToInverseGeV_);
0135 
0136   if (iHists_.bfield_.empty()) {
0137     iHists_.bfield_ = siStripLACalibration::fieldAsString(theMagField);
0138   } else {
0139     if (iHists_.bfield_ != siStripLACalibration::fieldAsString(theMagField)) {
0140       mismatchedBField_ = true;
0141     }
0142   }
0143 
0144   const SiStripLatency* apvlat = &iSetup.getData(m_latencyTokenBR);
0145   if (iHists_.apvmode_.empty()) {
0146     iHists_.apvmode_ = siStripLACalibration::apvModeAsString(apvlat);
0147   } else {
0148     if (iHists_.apvmode_ != siStripLACalibration::apvModeAsString(apvlat)) {
0149       mismatchedLatency_ = true;
0150     }
0151   }
0152 
0153   std::vector<uint32_t> c_rawid;
0154   std::vector<float> c_globalZofunitlocalY, c_localB, c_BdotY, c_driftx, c_drifty, c_driftz, c_lorentzAngle;
0155 
0156   auto dets = tkGeom.detsTIB();
0157   //dets.insert(dets.end(), tkGeom.detsTID().begin(), tkGeom.detsTID().end()); // no LA in endcaps
0158   dets.insert(dets.end(), tkGeom.detsTOB().begin(), tkGeom.detsTOB().end());
0159   //dets.insert(dets.end(), tkGeom.detsTEC().begin(), tkGeom.detsTEC().end()); // no LA in endcaps
0160 
0161   for (auto det : dets) {
0162     auto detid = det->geographicalId().rawId();
0163     const StripGeomDetUnit* stripDet = dynamic_cast<const StripGeomDetUnit*>(tkGeom.idToDet(det->geographicalId()));
0164     if (stripDet) {
0165       c_rawid.push_back(detid);
0166       c_globalZofunitlocalY.push_back(stripDet->toGlobal(LocalVector(0, 1, 0)).z());
0167       iHists_.orientation_[detid] = (stripDet->toGlobal(LocalVector(0, 1, 0)).z() < 0 ? -1 : 1);
0168       const auto locB = magField.inTesla(stripDet->surface().position());
0169       c_localB.push_back(locB.mag());
0170       c_BdotY.push_back(stripDet->surface().toLocal(locB).y());
0171       const auto drift = shallow::drift(stripDet, magField, lorentzAngle);
0172       c_driftx.push_back(drift.x());
0173       c_drifty.push_back(drift.y());
0174       c_driftz.push_back(drift.z());
0175       c_lorentzAngle.push_back(lorentzAngle.getLorentzAngle(detid));
0176       iHists_.la_db_[detid] = lorentzAngle.getLorentzAngle(detid);
0177       iHists_.moduleLocationType_[detid] = this->moduleLocationType(detid, tTopo);
0178     }
0179   }
0180 
0181   // Sorted DetId list gives max performance, anything else is worse
0182   std::sort(c_rawid.begin(), c_rawid.end());
0183 
0184   // initialize the hash map
0185   // in case it's not already initialized
0186   if (iHists_.hash_.size() == 0) {
0187     iHists_.hash_ = SiStripHashedDetId(c_rawid);
0188   }
0189 
0190   //reserve the size of the vector
0191   if (saveHistosMods_) {
0192     iHists_.h2_ct_w_m_.reserve(c_rawid.size());
0193     iHists_.h2_ct_var2_m_.reserve(c_rawid.size());
0194     iHists_.h2_ct_var3_m_.reserve(c_rawid.size());
0195 
0196     iHists_.h2_t_w_m_.reserve(c_rawid.size());
0197     iHists_.h2_t_var2_m_.reserve(c_rawid.size());
0198     iHists_.h2_t_var3_m_.reserve(c_rawid.size());
0199   }
0200 }
0201 
0202 std::string SiStripLorentzAnglePCLMonitor::moduleLocationType(const uint32_t& mod, const TrackerTopology* tTopo) {
0203   const SiStripDetId detid(mod);
0204   std::string subdet = "";
0205   unsigned int layer = 0;
0206   if (detid.subDetector() == SiStripDetId::TIB) {
0207     subdet = "TIB";
0208     layer = tTopo->layer(mod);
0209   }
0210 
0211   if (detid.subDetector() == SiStripDetId::TOB) {
0212     subdet = "TOB";
0213     layer = tTopo->layer(mod);
0214   }
0215 
0216   std::string type = (detid.stereo() ? "s" : "a");
0217   std::string d_l_t = Form("%s_L%d%s", subdet.c_str(), layer, type.c_str());
0218 
0219   if (layer == 0)
0220     return subdet;
0221   return d_l_t;
0222 }
0223 
0224 // ------------ method called for each event  ------------
0225 void SiStripLorentzAnglePCLMonitor::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0226   using namespace edm;
0227 
0228   // return immediately if the field is not consistent!
0229   if (mismatchedBField_)
0230     return;
0231 
0232   if (mismatchedLatency_)
0233     return;
0234 
0235   edm::Handle<edm::View<reco::Track>> tracks;
0236   iEvent.getByToken(m_tracks_token, tracks);
0237   edm::Handle<TrajTrackAssociationCollection> trajTrackAssociations;
0238   iEvent.getByToken(m_association_token, trajTrackAssociations);
0239 
0240   LogDebug(moduleDescription().moduleName()) << "I AM IN EVENT" << iEvent.id() << std::endl;
0241 
0242   std::vector<OnTrackCluster> clusters{};
0243 
0244   // first collect all the clusters
0245   for (const auto& assoc : *trajTrackAssociations) {
0246     const auto traj = assoc.key.get();
0247     const auto track = assoc.val.get();
0248 
0249     iHists_.h1_["track_pt"]->Fill(track->pt());
0250     iHists_.h1_["track_eta"]->Fill(track->eta());
0251     iHists_.h1_["track_phi"]->Fill(track->phi());
0252     iHists_.h1_["track_validhits"]->Fill(track->numberOfValidHits());
0253 
0254     const auto normChi2 = track->ndof() > 0 ? track->chi2() / track->ndof() : -1.;
0255     iHists_.h1_["track_chi2ndof"]->Fill(normChi2);
0256     iHists_.h2_["track_chi2xhits"]->Fill(normChi2, track->numberOfValidHits());
0257     iHists_.h2_["track_ptxhits"]->Fill(track->pt(), track->numberOfValidHits());
0258     iHists_.h2_["track_etaxhits"]->Fill(track->eta(), track->numberOfValidHits());
0259     iHists_.h2_["track_ptxchi2"]->Fill(track->pt(), normChi2);
0260     iHists_.h2_["track_ptxeta"]->Fill(track->pt(), track->eta());
0261     iHists_.h2_["track_etaxchi2"]->Fill(track->eta(), normChi2);
0262 
0263     edm::LogInfo("SiStripLorentzAnglePCLMonitor")
0264         << " track pT()" << track->pt() << " track eta()" << track->eta() << std::endl;
0265 
0266     for (const auto& meas : traj->measurements()) {
0267       const auto& trajState = meas.updatedState();
0268       if (!trajState.isValid())
0269         continue;
0270 
0271       // there can be 2 (stereo module), 1 (no stereo module), or 0 (no strip hit) clusters per measurement
0272       const auto trechit = meas.recHit()->hit();
0273       const auto simple1d = dynamic_cast<const SiStripRecHit1D*>(trechit);
0274       const auto simple = dynamic_cast<const SiStripRecHit2D*>(trechit);
0275       const auto matched = dynamic_cast<const SiStripMatchedRecHit2D*>(trechit);
0276       if (matched) {
0277         clusters.emplace_back(matched->monoId(), &matched->monoCluster(), traj, track, meas);
0278         clusters.emplace_back(matched->stereoId(), &matched->stereoCluster(), traj, track, meas);
0279       } else if (simple) {
0280         clusters.emplace_back(simple->geographicalId().rawId(), simple->cluster().get(), traj, track, meas);
0281       } else if (simple1d) {
0282         clusters.emplace_back(simple1d->geographicalId().rawId(), simple1d->cluster().get(), traj, track, meas);
0283       }
0284     }
0285   }
0286 
0287   for (const auto clus : clusters) {
0288     uint32_t c_nstrips = clus.cluster->amplitudes().size();
0289     m_clusterInfo.setCluster(*clus.cluster, clus.det);
0290     float c_variance = m_clusterInfo.variance();
0291     const auto& trajState = clus.measurement.updatedState();
0292     const auto trackDir = trajState.localDirection();
0293     float c_localdirx = trackDir.x();
0294     float c_localdiry = trackDir.y();
0295     float c_localdirz = trackDir.z();
0296     const auto hit = clus.measurement.recHit()->hit();
0297 
0298     // not yet needed (might be used for Backplane correction later on
0299     /*
0300       const auto& tkGeom = iSetup.getData(m_tkGeomToken);
0301       const auto stripDet = dynamic_cast<const StripGeomDetUnit*>(tkGeom.idToDet(hit->geographicalId()));
0302       float c_barycenter = stripDet->specificTopology().localPosition(clus.cluster->barycenter()).x();
0303       float c_localx = stripDet->toLocal(trajState.globalPosition()).x();
0304       float c_rhlocalx = hit->localPosition().x();
0305       float c_rhlocalxerr = hit->localPositionError().xx();
0306     */
0307 
0308     const uint32_t mod = hit->geographicalId().rawId();
0309 
0310     std::string locationtype = iHists_.moduleLocationType_[mod];
0311     if (locationtype.empty())
0312       return;
0313 
0314     // retrive the hashed index
0315     const auto& hashedIndex = iHists_.hash_.hashedIndex(mod);
0316 
0317     if (saveHistosMods_) {
0318       LogDebug("SiStripLorentzAnglePCLMonitor") << "module ID: " << mod << " hashedIndex: " << hashedIndex;
0319       iHists_.h1_["occupancyPerIndex"]->Fill(hashedIndex);
0320     }
0321 
0322     TVector3 localdir(c_localdirx, c_localdiry, c_localdirz);
0323     int sign = iHists_.orientation_[mod];
0324     float tantheta = TMath::Tan(localdir.Theta());
0325     float cosphi = TMath::Cos(localdir.Phi());
0326     float theta = localdir.Theta();
0327 
0328     iHists_.h1_[Form("%s_nstrips", locationtype.c_str())]->Fill(c_nstrips);
0329     iHists_.h1_[Form("%s_tanthetatrk", locationtype.c_str())]->Fill(sign * tantheta);
0330     iHists_.h1_[Form("%s_cosphitrk", locationtype.c_str())]->Fill(cosphi);
0331 
0332     // nstrips
0333     iHists_.h2_[Form("%s_tanthcosphtrk_nstrip", locationtype.c_str())]->Fill(sign * cosphi * tantheta, c_nstrips);
0334     iHists_.h2_[Form("%s_thetatrk_nstrip", locationtype.c_str())]->Fill(sign * theta * cosphi, c_nstrips);
0335 
0336     if (saveHistosMods_) {
0337       iHists_.h1_[Form("%s_%d_nstrips", locationtype.c_str(), mod)]->Fill(c_nstrips);
0338       iHists_.h1_[Form("%s_%d_tanthetatrk", locationtype.c_str(), mod)]->Fill(sign * tantheta);
0339       iHists_.h1_[Form("%s_%d_cosphitrk", locationtype.c_str(), mod)]->Fill(cosphi);
0340     }
0341 
0342     // variance for width == 2
0343     if (c_nstrips == 2) {
0344       iHists_.h1_[Form("%s_variance_w2", locationtype.c_str())]->Fill(c_variance);
0345       iHists_.h2_[Form("%s_tanthcosphtrk_var2", locationtype.c_str())]->Fill(sign * cosphi * tantheta, c_variance);
0346       iHists_.h2_[Form("%s_thcosphtrk_var2", locationtype.c_str())]->Fill(sign * cosphi * theta, c_variance);
0347 
0348       // not in PCL
0349       if (saveHistosMods_) {
0350         LogDebug("SiStripLorentzAnglePCLMonitor") << iHists_.h2_ct_var2_m_[hashedIndex]->getName();
0351         iHists_.h1_[Form("%s_%d_variance_w2", locationtype.c_str(), mod)]->Fill(c_variance);
0352         iHists_.h2_ct_var2_m_[hashedIndex]->Fill(sign * cosphi * tantheta, c_variance);
0353         iHists_.h2_t_var2_m_[hashedIndex]->Fill(sign * cosphi * theta, c_variance);
0354       }
0355     }
0356     // variance for width == 3
0357     if (c_nstrips == 3) {
0358       iHists_.h1_[Form("%s_variance_w3", locationtype.c_str())]->Fill(c_variance);
0359       iHists_.h2_[Form("%s_tanthcosphtrk_var3", locationtype.c_str())]->Fill(sign * cosphi * tantheta, c_variance);
0360       iHists_.h2_[Form("%s_thcosphtrk_var3", locationtype.c_str())]->Fill(sign * cosphi * theta, c_variance);
0361 
0362       // not in PCL
0363       if (saveHistosMods_) {
0364         iHists_.h1_[Form("%s_%d_variance_w3", locationtype.c_str(), mod)]->Fill(c_variance);
0365         iHists_.h2_ct_var3_m_[hashedIndex]->Fill(sign * cosphi * tantheta, c_variance);
0366         iHists_.h2_t_var3_m_[hashedIndex]->Fill(sign * cosphi * theta, c_variance);
0367       }
0368     }
0369 
0370     // not in PCL
0371     if (saveHistosMods_) {
0372       iHists_.h2_ct_w_m_[hashedIndex]->Fill(sign * cosphi * tantheta, c_nstrips);
0373       iHists_.h2_t_w_m_[hashedIndex]->Fill(sign * cosphi * theta, c_nstrips);
0374     }
0375   }
0376 }
0377 
0378 void SiStripLorentzAnglePCLMonitor::bookHistograms(DQMStore::IBooker& ibook,
0379                                                    edm::Run const& run,
0380                                                    edm::EventSetup const& iSetup) {
0381   std::string bvalue = (iHists_.bfield_ == "3.8") ? "B-ON" : "B-OFF";
0382   std::string folderToBook = fmt::format("{}/{}_{}", folder_, bvalue, iHists_.apvmode_);
0383 
0384   ibook.setCurrentFolder(folderToBook);
0385   edm::LogPrint(moduleDescription().moduleName()) << "booking in " << folderToBook;
0386 
0387   // prepare track histograms
0388   // clang-format off
0389   iHists_.h1_["track_pt"] = ibook.book1D("track_pt", "track p_{T};track p_{T} [GeV];n. tracks", 2000, 0, 1000);
0390   iHists_.h1_["track_eta"] = ibook.book1D("track_eta", "track #eta;track #eta;n. tracks", 100, -4, 4);
0391   iHists_.h1_["track_phi"] = ibook.book1D("track_phi", "track #phi;track #phi;n. tracks", 80, -3.2, 3.2);
0392   iHists_.h1_["track_validhits"] =
0393       ibook.book1D("track_validhits", "track n. valid hits;track n. valid hits;n. tracks", 50, 0, 50);
0394   iHists_.h1_["track_chi2ndof"] =
0395       ibook.book1D("track_chi2ndof", "track #chi^{2}/ndf;track #chi^{2}/ndf;n. tracks", 100, 0, 5);
0396   iHists_.h2_["track_chi2xhits"] =
0397       ibook.book2D("track_chi2xhits_2d",
0398                    "track track n. hits vs track #chi^{2}/ndf;track #chi^{2};track n. valid hits;tracks",
0399                    100, 0, 5, 50, 0, 50);
0400   iHists_.h2_["track_ptxhits"] = ibook.book2D(
0401       "track_ptxhits_2d", "track n. hits vs p_{T};track p_{T} [GeV];track n. valid hits;tracks", 200, 0, 100, 50, 0, 50);
0402   iHists_.h2_["track_etaxhits"] = ibook.book2D(
0403       "track_etaxhits_2d", "track n. hits vs track #eta;track #eta;track n. valid hits;tracks", 60, -3, 3, 50, 0, 50);
0404   iHists_.h2_["track_ptxchi2"] =
0405       ibook.book2D("track_ptxchi2_2d",
0406                    "track #chi^{2}/ndf vs track p_{T};track p_{T} [GeV]; track #chi^{2}/ndf;tracks",
0407                    200, 0, 100, 100, 0, 5);
0408   iHists_.h2_["track_ptxeta"] = ibook.book2D(
0409       "track_ptxeta_2d", "track #eta vs track p_{T};track p_{T} [GeV];track #eta;tracks", 200, 0, 100, 60, -3, 3);
0410   iHists_.h2_["track_etaxchi2"] = ibook.book2D(
0411       "track_etaxchi2_2d", "track #chi^{2}/ndf vs track #eta;track #eta;track #chi^{2};tracks", 60, -3, 3, 100, 0, 5);
0412   // clang-format on
0413 
0414   if (saveHistosMods_) {
0415     iHists_.h1_["occupancyPerIndex"] = ibook.book1D("ClusterOccupancyPerHashedIndex",
0416                                                     "cluster occupancy;hashed index;# clusters per module",
0417                                                     iHists_.hash_.size(),
0418                                                     -0.5,
0419                                                     iHists_.hash_.size() - 0.5);
0420   }
0421 
0422   // fill in the module types
0423   iHists_.nlayers_["TIB"] = 4;
0424   iHists_.nlayers_["TOB"] = 6;
0425   iHists_.modtypes_.push_back("s");
0426   iHists_.modtypes_.push_back("a");
0427 
0428   // prepare type histograms
0429   for (auto& layers : iHists_.nlayers_) {
0430     std::string subdet = layers.first;
0431     for (int l = 1; l <= layers.second; ++l) {
0432       ibook.setCurrentFolder(folderToBook + Form("/%s/L%d", subdet.c_str(), l));
0433       for (auto& t : iHists_.modtypes_) {
0434         // do not fill stereo where there aren't
0435         if (l > 2 && t == "s")
0436           continue;
0437         std::string locType = Form("%s_L%d%s", subdet.c_str(), l, t.c_str());
0438 
0439         // clang-format off
0440     const char* titles = Form("n.strips in %s;n.strips;n. clusters", locType.c_str());
0441         iHists_.h1_[Form("%s_nstrips", locType.c_str())] = ibook.book1D(Form("%s_nstrips", locType.c_str()), titles, 20, 0, 20);
0442 
0443     titles =  Form("tan(#theta_{trk}) in %s;tan(#theta_{trk});n. clusters", locType.c_str());
0444         iHists_.h1_[Form("%s_tanthetatrk", locType.c_str())] = ibook.book1D(Form("%s_tanthetatrk", locType.c_str()), titles, 300, -1.5, 1.5);
0445 
0446     titles =  Form("cos(#phi_{trk}) in %s;cos(#phi_{trk});n. clusters", locType.c_str());
0447         iHists_.h1_[Form("%s_cosphitrk", locType.c_str())] = ibook.book1D(Form("%s_cosphitrk", locType.c_str()), titles, 40, -1, 1);
0448 
0449     titles = Form("Cluster variance (w=2) in %s;cluster variance (w=2);n. clusters", locType.c_str());
0450         iHists_.h1_[Form("%s_variance_w2", locType.c_str())] = ibook.book1D(Form("%s_variance_w2", locType.c_str()), titles, 100, 0, 1);
0451 
0452     titles = Form("Cluster variance (w=3) in %s;cluster variance (w=3);n. clusters", locType.c_str());
0453         iHists_.h1_[Form("%s_variance_w3", locType.c_str())] = ibook.book1D(Form("%s_variance_w3", locType.c_str()), titles, 100, 0, 1);
0454 
0455     titles = Form("n. strips in %s vs tan(#theta_{trk})cos(#phi_{trk});tan(#theta_{trk})cos(#phi_{trk});n. strips;n. clusters", locType.c_str());
0456         iHists_.h2_[Form("%s_tanthcosphtrk_nstrip", locType.c_str())] = ibook.book2D(Form("%s_tanthcosphtrk_nstrip", locType.c_str()), titles, 360, -0.9, 0.9, 20, 0, 20);
0457 
0458     titles = Form("n. strips in %s vs #theta_{trk};#theta_{trk} [rad];n. strips;n. clusters", locType.c_str());
0459         iHists_.h2_[Form("%s_thetatrk_nstrip", locType.c_str())] = ibook.book2D(Form("%s_thetatrk_nstrip", locType.c_str()), titles, 360, -0.9, 0.9, 20, 0, 20);
0460 
0461     titles = Form("cluster variance (w=2) in %s vs tan(#theta_{trk})cos(#phi_{trk});tan(#theta_{trk})cos(#phi_{trk});cluster variance (w=2);n. clusters", locType.c_str());
0462         iHists_.h2_[Form("%s_tanthcosphtrk_var2", locType.c_str())] = ibook.book2D(Form("%s_tanthcosphtrk_var2", locType.c_str()), titles, 360, -0.9, 0.9, 50, 0, 1);
0463 
0464     titles =  Form("cluster variance (w=3) in %s vs tan(#theta_{trk})cos(#phi_{trk});tan(#theta_{trk})cos(#phi_{trk});cluster variance (w=3);n. clusters", locType.c_str());
0465         iHists_.h2_[Form("%s_tanthcosphtrk_var3", locType.c_str())] = ibook.book2D(Form("%s_tanthcosphtrk_var3", locType.c_str()), titles, 360, -0.9, 0.9, 50, 0, 1);
0466 
0467     titles =  Form("cluster variance (w=2) in %s vs #theta_{trk}cos(#phi_{trk});#theta_{trk}cos(#phi_{trk});cluster variance (w=2);n. clusters", locType.c_str());
0468         iHists_.h2_[Form("%s_thcosphtrk_var2", locType.c_str())] = ibook.book2D(Form("%s_thcosphtrk_var2", locType.c_str()), titles, 360, -0.9, 0.9, 50, 0, 1);
0469 
0470     titles = Form("cluster variance (w=3) in %s vs #theta_{trk}cos(#phi_{trk});#theta_{trk}cos(#phi_{trk});cluster variance (w=3);n. clusters", locType.c_str());
0471         iHists_.h2_[Form("%s_thcosphtrk_var3", locType.c_str())] = ibook.book2D(Form("%s_thcosphtrk_var3", locType.c_str()), titles, 360, -0.9, 0.9, 50, 0, 1);
0472         // clang-format on
0473       }
0474     }
0475   }
0476 
0477   // prepare module histograms
0478   if (saveHistosMods_) {
0479     ibook.setCurrentFolder(folderToBook + "/modules");
0480     for (const auto& [mod, locationType] : iHists_.moduleLocationType_) {
0481       ibook.setCurrentFolder(folderToBook + "/modules" + Form("/%s", locationType.c_str()));
0482       // histograms for each module
0483       iHists_.h1_[Form("%s_%d_nstrips", locationType.c_str(), mod)] =
0484           ibook.book1D(Form("%s_%d_nstrips", locationType.c_str(), mod), "", 10, 0, 10);
0485       iHists_.h1_[Form("%s_%d_tanthetatrk", locationType.c_str(), mod)] =
0486           ibook.book1D(Form("%s_%d_tanthetatrk", locationType.c_str(), mod), "", 40, -1., 1.);
0487       iHists_.h1_[Form("%s_%d_cosphitrk", locationType.c_str(), mod)] =
0488           ibook.book1D(Form("%s_%d_cosphitrk", locationType.c_str(), mod), "", 40, -1, 1);
0489       iHists_.h1_[Form("%s_%d_variance_w2", locationType.c_str(), mod)] =
0490           ibook.book1D(Form("%s_%d_variance_w2", locationType.c_str(), mod), "", 20, 0, 1);
0491       iHists_.h1_[Form("%s_%d_variance_w3", locationType.c_str(), mod)] =
0492           ibook.book1D(Form("%s_%d_variance_w3", locationType.c_str(), mod), "", 20, 0, 1);
0493     }
0494 
0495     int counter{0};
0496     SiStripHashedDetId::const_iterator iter = iHists_.hash_.begin();
0497     for (; iter != iHists_.hash_.end(); ++iter) {
0498       LogDebug("SiStripLorentzAnglePCLMonitor")
0499           << "detId: " << (*iter) << " hashed index: " << iHists_.hash_.hashedIndex((*iter));
0500       const auto& locationType = iHists_.moduleLocationType_[(*iter)];
0501       ibook.setCurrentFolder(folderToBook + "/modules" + Form("/%s", locationType.c_str()));
0502       iHists_.h2_ct_w_m_.push_back(
0503           ibook.book2D(Form("ct_w_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 10, 0, 10));
0504       iHists_.h2_t_w_m_.push_back(
0505           ibook.book2D(Form("t_w_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 10, 0, 10));
0506       iHists_.h2_ct_var2_m_.push_back(
0507           ibook.book2D(Form("ct_var2_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 20, 0, 1));
0508       iHists_.h2_ct_var3_m_.push_back(
0509           ibook.book2D(Form("ct_var3_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 20, 0, 1));
0510       iHists_.h2_t_var2_m_.push_back(
0511           ibook.book2D(Form("t_var2_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 20, 0, 1));
0512       iHists_.h2_t_var3_m_.push_back(
0513           ibook.book2D(Form("t_var3_m_%s_%d", locationType.c_str(), *iter), "", 90, -0.9, 0.9, 20, 0, 1));
0514       counter++;
0515     }
0516     edm::LogPrint(moduleDescription().moduleName())
0517         << __PRETTY_FUNCTION__ << " Booked " << counter << " module level histograms!";
0518   }  // if saveHistoMods
0519 }
0520 
0521 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0522 void SiStripLorentzAnglePCLMonitor::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0523   edm::ParameterSetDescription desc;
0524   desc.add<std::string>("folder", "AlCaReco/SiStripLorentzAngle")->setComment("DQM folder to write into");
0525   desc.add<bool>("saveHistoMods", false)->setComment("save module level hisotgrams. Warning! takes a lot of space!");
0526   desc.add<edm::InputTag>("Tracks", edm::InputTag("SiStripCalCosmics"))->setComment("input track collection");
0527   descriptions.addWithDefaultLabel(desc);
0528 }
0529 
0530 // define this as a plug-in
0531 DEFINE_FWK_MODULE(SiStripLorentzAnglePCLMonitor);