Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-06-04 22:36:04

0001 /****************************************************************************
0002 *
0003 * This is a part of CTPPSDQM software.
0004 * Authors:
0005 *   Jan Kašpar (jan.kaspar@gmail.com)
0006 *   Nicola Minafra
0007 *   Laurent Forthomme
0008 *   Christopher Misan
0009 *
0010 ****************************************************************************/
0011 
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 #include "FWCore/Framework/interface/Event.h"
0014 #include "FWCore/Framework/interface/EventSetup.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/Utilities/interface/ESGetToken.h"
0017 #include "FWCore/Utilities/interface/InputTag.h"
0018 #include "FWCore/Utilities/interface/Transition.h"
0019 #include "FWCore/Framework/interface/Run.h"
0020 
0021 #include "DQMServices/Core/interface/DQMOneEDAnalyzer.h"
0022 #include "DQMServices/Core/interface/DQMStore.h"
0023 
0024 #include "DataFormats/Provenance/interface/EventRange.h"
0025 #include "DataFormats/CTPPSDigi/interface/TotemVFATStatus.h"
0026 #include "DataFormats/CTPPSDigi/interface/TotemFEDInfo.h"
0027 #include "DataFormats/Common/interface/DetSetVector.h"
0028 
0029 #include "DataFormats/CTPPSDetId/interface/CTPPSDiamondDetId.h"
0030 #include "DataFormats/CTPPSDigi/interface/CTPPSDiamondDigi.h"
0031 
0032 #include "DataFormats/CTPPSReco/interface/CTPPSPixelLocalTrack.h"
0033 #include "DataFormats/CTPPSDetId/interface/CTPPSPixelDetId.h"
0034 
0035 #include "DataFormats/CTPPSReco/interface/CTPPSDiamondRecHit.h"
0036 #include "DataFormats/CTPPSReco/interface/CTPPSDiamondLocalTrack.h"
0037 #include "Geometry/VeryForwardGeometryBuilder/interface/CTPPSGeometry.h"
0038 #include "Geometry/Records/interface/VeryForwardRealGeometryRecord.h"
0039 
0040 #include <string>
0041 
0042 //----------------------------------------------------------------------------------------------------
0043 
0044 // Utility for efficiency computations
0045 bool channelAlignedWithTrack(const CTPPSGeometry* geom,
0046                              const CTPPSDiamondDetId& detid,
0047                              const CTPPSDiamondLocalTrack& localTrack,
0048                              const float tolerance = 1) {
0049   const DetGeomDesc* det = geom->sensor(detid);
0050   const float x_pos = det->translation().x(),
0051               x_width = 2.0 * det->getDiamondDimensions().xHalfWidth;  // parameters stand for half the size
0052   return ((x_pos + 0.5 * x_width > localTrack.x0() - localTrack.x0Sigma() - tolerance &&
0053            x_pos + 0.5 * x_width < localTrack.x0() + localTrack.x0Sigma() + tolerance) ||
0054           (x_pos - 0.5 * x_width > localTrack.x0() - localTrack.x0Sigma() - tolerance &&
0055            x_pos - 0.5 * x_width < localTrack.x0() + localTrack.x0Sigma() + tolerance) ||
0056           (x_pos - 0.5 * x_width < localTrack.x0() - localTrack.x0Sigma() - tolerance &&
0057            x_pos + 0.5 * x_width > localTrack.x0() + localTrack.x0Sigma() + tolerance));
0058 }
0059 
0060 namespace dds {
0061   struct Cache {
0062     std::unordered_map<unsigned int, std::unique_ptr<TH2F>> hitDistribution2dMap;
0063 
0064     std::unordered_map<unsigned int, unsigned long> hitsCounterMap;
0065   };
0066 }  // namespace dds
0067 
0068 class CTPPSDiamondDQMSource : public DQMOneEDAnalyzer<edm::LuminosityBlockCache<dds::Cache>> {
0069 public:
0070   CTPPSDiamondDQMSource(const edm::ParameterSet&);
0071 
0072 protected:
0073   void dqmBeginRun(const edm::Run&, const edm::EventSetup&) override;
0074   void dqmEndRun(edm::Run const&, edm::EventSetup const&) override;
0075   void bookHistograms(DQMStore::IBooker&, const edm::Run&, const edm::EventSetup&) override;
0076   void analyze(const edm::Event&, const edm::EventSetup&) override;
0077   std::shared_ptr<dds::Cache> globalBeginLuminosityBlock(const edm::LuminosityBlock&,
0078                                                          const edm::EventSetup&) const override;
0079   void globalEndLuminosityBlock(const edm::LuminosityBlock&, const edm::EventSetup&) override;
0080 
0081 private:
0082   // Helper method to transform coordinates to local coordinate system
0083   CTPPSGeometry::Vector transformToLocalCoordinates(const CTPPSDetId& detId_pot,
0084                                                     double x,
0085                                                     double y,
0086                                                     double z = 0) const {
0087     auto localVector = CTPPSGeometry::Vector(x, y, z);
0088     localVector -= diamTranslations_.at(detId_pot);
0089     localVector = diamRotations_.at(detId_pot).Inverse() * localVector;
0090     return localVector;
0091   }
0092 
0093   // Constants
0094   /// Number of seconds per lumisection: used to compute hit rates in Hz
0095   static constexpr double SEC_PER_LUMI_SECTION = 23.31;
0096   /// Channel ID of the VFAT that contains clock data
0097   static constexpr unsigned short CHANNEL_OF_VFAT_CLOCK = 30;
0098   /// Bin width of histograms showing hits and tracks (in mm)
0099   static constexpr double DISPLAY_RESOLUTION_FOR_HITS_MM = 0.1;
0100   static constexpr double INV_DISPLAY_RESOLUTION_FOR_HITS_MM = 1. / DISPLAY_RESOLUTION_FOR_HITS_MM;
0101   /// ns per HPTDC bin
0102   static constexpr double HPTDC_BIN_WIDTH_NS = 25. / 1024;
0103   static constexpr unsigned short CTPPS_PIXEL_STATION_ID = 2;
0104   static constexpr unsigned short CTPPS_PIXEL_FAR_RP_ID = 3;
0105   static constexpr unsigned short CTPPS_DIAMOND_CYL_STATION_ID = 1;
0106   static constexpr unsigned short CTPPS_DIAMOND_CYL_RP_ID = 6;
0107   static constexpr unsigned short CTPPS_DIAMOND_NUM_OF_CHANNELS = 12;
0108   static constexpr unsigned short CTPPS_FED_ID_45 = 583;
0109   static constexpr unsigned short CTPPS_FED_ID_56 = 582;
0110   static constexpr unsigned short HPTDC_0_CHANNEL = 6;
0111   static constexpr unsigned short HPTDC_1_CHANNEL = 7;
0112   /// Number of OOT indices monitored
0113   static constexpr unsigned int FIRST_RUN_W_PIXELS = 300000;
0114 
0115   bool perLSsaving_;  //to avoid nanoDQMIO crashing, driven by  DQMServices/Core/python/DQMStore_cfi.py
0116 
0117   /// plots related to the whole system
0118   struct GlobalPlots {
0119     GlobalPlots() = default;
0120     GlobalPlots(DQMStore::IBooker& ibooker);
0121   };
0122   /// plots related to one sector
0123   struct SectorPlots {
0124     // Tracks
0125     MonitorElement* trackCorrelation = nullptr;
0126     MonitorElement* trackCorrelationLowMultiplicity = nullptr;
0127     MonitorElement* digiSentPercentage = nullptr;
0128     SectorPlots() {}
0129     SectorPlots(DQMStore::IBooker& ibooker, unsigned int id, bool plotOnline);
0130   };
0131   /// plots related to one Diamond detector package
0132   struct PotPlots {
0133     std::unordered_map<unsigned int, MonitorElement*> activity_per_bx;
0134 
0135     MonitorElement* hitDistribution2d = nullptr;
0136     MonitorElement* hitDistribution2d_lumisection = nullptr;
0137     MonitorElement* hitDistribution2dOOT = nullptr;
0138     MonitorElement* hitDistribution2dOOT_le = nullptr;
0139     MonitorElement *activePlanes = nullptr, *activePlanesInclusive = nullptr;
0140 
0141     MonitorElement* trackDistribution = nullptr;
0142     MonitorElement* trackDistributionOOT = nullptr;
0143 
0144     std::unordered_map<unsigned int, MonitorElement*> pixelTomographyAll;
0145 
0146     MonitorElement *leadingEdgeCumulative_both = nullptr, *leadingEdgeCumulative_all = nullptr,
0147                    *leadingEdgeCumulative_le = nullptr, *trailingEdgeCumulative_te = nullptr;
0148     MonitorElement* timeOverThresholdCumulativePot = nullptr;  //, *leadingTrailingCorrelationPot = nullptr;
0149     MonitorElement* leadingWithoutTrailingCumulativePot = nullptr;
0150 
0151     MonitorElement* ECCheck = nullptr;
0152 
0153     MonitorElement* HPTDCErrorFlags_2D = nullptr;
0154     MonitorElement* MHComprensive = nullptr;
0155 
0156     MonitorElement* recHitTime = nullptr;
0157 
0158     // MonitorElement* clock_Digi1_le = nullptr;
0159     // MonitorElement* clock_Digi1_te = nullptr;
0160     // MonitorElement* clock_Digi3_le = nullptr;
0161     // MonitorElement* clock_Digi3_te = nullptr;
0162 
0163     unsigned int HitCounter, MHCounter, LeadingOnlyCounter, TrailingOnlyCounter, CompleteCounter;
0164 
0165     std::map<int, int> effTriplecountingChMap;
0166     std::map<int, int> effDoublecountingChMap;
0167     MonitorElement* EfficiencyOfChannelsInPot = nullptr;
0168     TH2F pixelTracksMap;
0169 
0170     // MonitorElement* TOTVsLS = nullptr;
0171     // MonitorElement* trackTimeVsLS = nullptr;
0172     MonitorElement* trackTimeVsBX = nullptr;
0173     // MonitorElement* trackTimeVsXAngle = nullptr;
0174 
0175     // MonitorElement* TOTVsLSProfile = nullptr;
0176     // MonitorElement* trackTimeVsLSProfile = nullptr;
0177     MonitorElement* trackTimeVsBXProfile = nullptr;
0178     // MonitorElement* trackTimeVsXAngleProfile = nullptr;
0179 
0180     PotPlots() = default;
0181     PotPlots(DQMStore::IBooker& ibooker,
0182              unsigned int id,
0183              unsigned int windowsNum,
0184              bool plotOnline,
0185              bool plotOffline,
0186              bool perLSsaving);
0187   };
0188   /// plots related to one Diamond plane
0189   struct PlanePlots {
0190     MonitorElement* digiProfileCumulativePerPlane = nullptr;
0191     MonitorElement* hitProfile = nullptr;
0192     MonitorElement* hit_multiplicity = nullptr;
0193 
0194     MonitorElement *leadingEdgeCumulativePerPlane_both = nullptr, *leadingEdgeCumulativePerPlane_le = nullptr,
0195                    *trailingEdgeCumulativePerPlane_te = nullptr;
0196     MonitorElement* TimeOverThresholdCumulativePerPlane = nullptr;
0197     MonitorElement* pixelTomography_far = nullptr;
0198     MonitorElement* EfficiencyWRTPixelsInPlane = nullptr;
0199 
0200     TH2F pixelTracksMapWithDiamonds;
0201 
0202     PlanePlots() = default;
0203     PlanePlots(DQMStore::IBooker& ibooker, unsigned int id, unsigned int windowsNum);
0204   };
0205   /// plots related to one Diamond channel
0206   struct ChannelPlots {
0207     std::unordered_map<unsigned int, MonitorElement*> activity_per_bx;
0208 
0209     MonitorElement* HPTDCErrorFlags = nullptr;
0210     MonitorElement *leadingEdgeCumulative_both = nullptr, *leadingEdgeCumulative_le = nullptr,
0211                    *trailingEdgeCumulative_te = nullptr;
0212     MonitorElement* TimeOverThresholdCumulativePerChannel = nullptr;
0213     //MonitorElement* LeadingTrailingCorrelationPerChannel = nullptr;
0214     MonitorElement* leadingWithoutTrailing = nullptr;
0215     MonitorElement* pixelTomography_far = nullptr;
0216     MonitorElement* hit_rate = nullptr;
0217     MonitorElement* recHitTime = nullptr;
0218 
0219     unsigned int HitCounter, MHCounter, LeadingOnlyCounter, TrailingOnlyCounter, CompleteCounter;
0220 
0221     ChannelPlots() = default;
0222     ChannelPlots(DQMStore::IBooker& ibooker, unsigned int id, unsigned int windowsNum);
0223   };
0224 
0225   void checkEventNumber(const CTPPSDiamondDetId&, const TotemFEDInfo&, const TotemVFATStatus&, PotPlots&, int&) const;
0226 
0227   edm::EDGetTokenT<edm::DetSetVector<TotemVFATStatus>> tokenStatus_;
0228   edm::EDGetTokenT<edm::DetSetVector<CTPPSPixelLocalTrack>> tokenPixelTrack_;
0229   edm::EDGetTokenT<edm::DetSetVector<CTPPSDiamondDigi>> tokenDigi_;
0230   edm::EDGetTokenT<edm::DetSetVector<CTPPSDiamondRecHit>> tokenDiamondHit_;
0231   edm::EDGetTokenT<edm::DetSetVector<CTPPSDiamondLocalTrack>> tokenDiamondTrack_;
0232   edm::EDGetTokenT<std::vector<TotemFEDInfo>> tokenFEDInfo_;
0233 
0234   edm::ESGetToken<CTPPSGeometry, VeryForwardRealGeometryRecord> ctppsGeometryRunToken_;
0235   edm::ESGetToken<CTPPSGeometry, VeryForwardRealGeometryRecord> ctppsGeometryEventToken_;
0236 
0237   bool excludeMultipleHits_;
0238   const bool extract_digi_info_;
0239   struct DiamondShifts {
0240     double global, withPixels;
0241   };
0242   std::unordered_map<CTPPSDetId, DiamondShifts> diamShifts_;
0243   std::unordered_map<CTPPSDetId, double> diamHalfWidths_;
0244   std::unordered_map<CTPPSDetId, ROOT::Math::Rotation3D> diamRotations_;
0245   std::unordered_map<CTPPSDetId, ROOT::Math::DisplacementVector3D<ROOT::Math::Cartesian3D<double>>> diamTranslations_;
0246   std::vector<std::pair<edm::EventRange, int>> runParameters_;
0247   int centralOOT_;
0248   unsigned int verbosity_;
0249   const bool plotOnline_;
0250   const bool plotOffline_;
0251   unsigned int windowsNum_;
0252   unsigned int trackCorrelationThreshold_;
0253 
0254   GlobalPlots globalPlot_;
0255   std::unordered_map<unsigned int, PotPlots> potPlots_;
0256   std::unordered_map<unsigned int, SectorPlots> sectorPlots_;
0257   std::unordered_map<unsigned int, PlanePlots> planePlots_;
0258   std::unordered_map<unsigned int, ChannelPlots> channelPlots_;
0259 
0260   int EC_difference_56_, EC_difference_45_;
0261 };
0262 
0263 //----------------------------------------------------------------------------------------------------
0264 
0265 CTPPSDiamondDQMSource::GlobalPlots::GlobalPlots(DQMStore::IBooker& ibooker) { ibooker.setCurrentFolder("CTPPS"); }
0266 
0267 //----------------------------------------------------------------------------------------------------
0268 
0269 CTPPSDiamondDQMSource::SectorPlots::SectorPlots(DQMStore::IBooker& ibooker, unsigned int id, bool plotOnline) {
0270   std::string path, title;
0271   CTPPSDiamondDetId(id).armName(path, CTPPSDiamondDetId::nPath);
0272   ibooker.setCurrentFolder(path);
0273 
0274   CTPPSDiamondDetId(id).armName(title, CTPPSDiamondDetId::nFull);
0275 
0276   trackCorrelation = ibooker.book2D("tracks correlation near-far",
0277                                     title + " tracks correlation near-far;track x 220nr_hr (mm);track x 220cyl (mm)",
0278                                     19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM,
0279                                     -1,
0280                                     18,
0281                                     19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM,
0282                                     -1,
0283                                     18);
0284   trackCorrelationLowMultiplicity = ibooker.book2D(
0285       "tracks correlation with low multiplicity near-far",
0286       title + " tracks correlation with low multiplicity near-far;track x 220nr_hr (mm);track x 220cyl(mm)",
0287       19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM,
0288       -1,
0289       18,
0290       19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM,
0291       -1,
0292       18);
0293 }
0294 
0295 //----------------------------------------------------------------------------------------------------
0296 CTPPSDiamondDQMSource::PotPlots::PotPlots(DQMStore::IBooker& ibooker,
0297                                           unsigned int id,
0298                                           unsigned int windowsNum,
0299                                           bool plotOnline,
0300                                           bool plotOffline,
0301                                           bool perLSsaving)
0302     : HitCounter(0),
0303       MHCounter(0),
0304       LeadingOnlyCounter(0),
0305       TrailingOnlyCounter(0),
0306       CompleteCounter(0),
0307       pixelTracksMap("Pixel track maps for efficiency", "Pixel track maps for efficiency", 25, 0, 25, 24, -6, 18) {
0308   std::string path, title;
0309   CTPPSDiamondDetId(id).rpName(path, CTPPSDiamondDetId::nPath);
0310   ibooker.setCurrentFolder(path);
0311 
0312   CTPPSDiamondDetId(id).rpName(title, CTPPSDiamondDetId::nFull);
0313 
0314   if (plotOnline) {
0315     hitDistribution2d_lumisection =
0316         ibooker.book2D("hits in planes lumisection",
0317                        title + " hits in planes in the last lumisection;plane number;x (mm)",
0318                        10,
0319                        -0.5,
0320                        4.5,
0321                        19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM,
0322                        -0.5,
0323                        18.5);
0324 
0325     hitDistribution2dOOT_le =
0326         ibooker.book2D("hits with OOT in planes (le only)",
0327                        title + " hits with OOT in planes (le only);plane number, OOT index;x (mm)",
0328                        1 + windowsNum * 4,
0329                        -1. / windowsNum,
0330                        4,
0331                        19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM,
0332                        -0.5,
0333                        18.5);
0334 
0335     activePlanesInclusive =
0336         ibooker.book1D("active planes inclusive",
0337                        title + " active planes, MH and le only included (per event);number of active planes",
0338                        6,
0339                        -0.5,
0340                        5.5);
0341 
0342     ECCheck = ibooker.book1D("optorxEC(8bit) - vfatEC", title + " EC Error;optorxEC-vfatEC", 50, -25, 25);
0343 
0344     EfficiencyOfChannelsInPot =
0345         ibooker.book2D("Efficiency in channels",
0346                        title + " Efficiency (%) in channels (diamonds only);plane number;ch number",
0347                        10,
0348                        -0.5,
0349                        4.5,
0350                        14,
0351                        -1,
0352                        13);
0353   }
0354 
0355   if (plotOffline && !perLSsaving) {
0356     ibooker.setCurrentFolder(path + "/timing_profiles");
0357     // TOTVsLS=ibooker.book2D("ToT vs LS",title +" ToT vs LS;LS;ToT(ns)",4000,0,4000, 200,5,25);
0358     // trackTimeVsLS=ibooker.book2D("track time vs LS",title+" track time vs LS;LS;track_time(ns)",4000,0,4000, 500, -25, 25);
0359     trackTimeVsBX =
0360         ibooker.book2D("track time vs BX", title + " track time vs BX;BX;track_time(ns)", 4000, 0, 4000, 500, -25, 25);
0361     // trackTimeVsXAngle = ibooker.book2D(
0362     //     "track time vs xangle", title + " track time vs xangle;xangle;track_time(ns)", 60, 120, 180, 500, -25, 25);
0363 
0364     // TOTVsLSProfile=ibooker.bookProfile("ToT vs LS profile",title+" ToT vs LS profile;LS;track_time(ns)", 500, -25, 25,4000,0,4000);
0365     // trackTimeVsLSProfile=ibooker.bookProfile("track time vs LS profile",title+" track time vs LS profile;LS;track_time(ns)", 500, -25, 25,4000,0,4000);
0366     trackTimeVsBXProfile = ibooker.bookProfile(
0367         "track time vs BX profile", title + " track time vs BX profile;BX;track_time(ns)", 500, -25, 25, 4000, 0, 4000);
0368     // trackTimeVsXAngleProfile = ibooker.bookProfile("track time vs xangle profile",
0369     //                                                title + " track time vs xangle profile;xangle;track_time(ns)",
0370     //                                                500,
0371     //                                                -25,
0372     //                                                25,
0373     //                                                60,
0374     //                                                120,
0375     //                                                180);
0376     ibooker.setCurrentFolder(path);
0377   }
0378 
0379   for (unsigned int i = 0; i < windowsNum; i++) {
0380     std::string window = std::to_string(i * 25) + "-" + std::to_string((i + 1) * 25);
0381     activity_per_bx[i] = ibooker.book1D(
0382         "activity per BX " + window, title + " Activity per BX " + window + " ns;Event.BX", 3600, -1.5, 3598. + 0.5);
0383     pixelTomographyAll[i] =
0384         ibooker.book2D("tomography pixel " + window,
0385                        title + " tomography with pixel " + window + " ns (all planes);x + 25*plane(mm);y (mm)",
0386                        100,
0387                        0,
0388                        100,
0389                        10,
0390                        -5,
0391                        5);
0392   }
0393 
0394   hitDistribution2d = ibooker.book2D("hits in planes",
0395                                      title + " hits in planes;plane number;x (mm)",
0396                                      10,
0397                                      -0.5,
0398                                      4.5,
0399                                      19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM,
0400                                      -0.5,
0401                                      18.5);
0402 
0403   hitDistribution2dOOT = ibooker.book2D("hits with OOT in planes",
0404                                         title + " hits with OOT in planes;plane number, OOT index;x (mm)",
0405                                         1 + windowsNum * 4,
0406                                         -1. / windowsNum,
0407                                         4,
0408                                         19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM,
0409                                         -0.5,
0410                                         18.5);
0411 
0412   {               // bin labelling (for clarity)
0413     int idx = 2;  // start counting at 1, first bin is empty
0414     for (int pl = 0; pl < 4; ++pl)
0415       for (unsigned int oot = 0; oot < windowsNum; ++oot) {
0416         const std::string bin_label =
0417             (oot == 0 ? "Plane " + std::to_string(pl) + ", " : "") + "OOT" + std::to_string(oot);
0418         hitDistribution2dOOT->setBinLabel(idx, bin_label);
0419         if (plotOnline)
0420           hitDistribution2dOOT_le->setBinLabel(idx, bin_label);
0421         ++idx;
0422       }
0423   }
0424 
0425   recHitTime = ibooker.book1D("recHit time", title + " recHit time; t (ns)", 500, -25, 25);
0426 
0427   activePlanes =
0428       ibooker.book1D("active planes", title + " active planes (per event);number of active planes", 6, -0.5, 5.5);
0429 
0430   trackDistribution =
0431       ibooker.book1D("tracks", title + " tracks;x (mm)", 19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM, -0.5, 18.5);
0432   trackDistributionOOT = ibooker.book2D("tracks with OOT",
0433                                         title + " tracks with OOT;plane number;x (mm)",
0434                                         9,
0435                                         -0.5,
0436                                         4,
0437                                         19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM,
0438                                         -0.5,
0439                                         18.5);
0440 
0441   leadingEdgeCumulative_both = ibooker.book1D("leading edge (le and te)",
0442                                               title + " leading edge (le and te) (recHits); leading edge (ns)",
0443                                               25 * windowsNum,
0444                                               0,
0445                                               25 * windowsNum);
0446   leadingEdgeCumulative_all = ibooker.book1D("leading edge (all)",
0447                                              title + " leading edge (with or without te) (DIGIs); leading edge (ns)",
0448                                              25 * windowsNum,
0449                                              0,
0450                                              25 * windowsNum);
0451   leadingEdgeCumulative_le = ibooker.book1D("leading edge (le only)",
0452                                             title + " leading edge (le only) (DIGIs); leading edge (ns)",
0453                                             25 * windowsNum,
0454                                             0,
0455                                             25 * windowsNum);
0456   trailingEdgeCumulative_te = ibooker.book1D("trailing edge (te only)",
0457                                              title + " trailing edge (te only) (DIGIs); trailing edge (ns)",
0458                                              25 * windowsNum,
0459                                              0,
0460                                              25 * windowsNum);
0461   timeOverThresholdCumulativePot =
0462       ibooker.book1D("time over threshold", title + " time over threshold;time over threshold (ns)", 250, -25, 100);
0463   // leadingTrailingCorrelationPot =
0464   //     ibooker.book2D("leading trailing correlation",
0465   //                    title + " leading trailing correlation;leading edge (ns);trailing edge (ns)",
0466   //                    75,
0467   //                    0,
0468   //                    75,
0469   //                    75,
0470   //                    0,
0471   //                    75);
0472 
0473   leadingWithoutTrailingCumulativePot =
0474       ibooker.book1D("event category", title + " leading edges without trailing;;%", 3, 0.5, 3.5);
0475   leadingWithoutTrailingCumulativePot->setBinLabel(1, "Leading only");
0476   leadingWithoutTrailingCumulativePot->setBinLabel(2, "Trailing only");
0477   leadingWithoutTrailingCumulativePot->setBinLabel(3, "Both");
0478 
0479   HPTDCErrorFlags_2D = ibooker.book2D("HPTDC Errors", title + " HPTDC Errors", 16, -0.5, 16.5, 9, -0.5, 8.5);
0480   for (unsigned short error_index = 1; error_index < 16; ++error_index)
0481     HPTDCErrorFlags_2D->setBinLabel(error_index, HPTDCErrorFlags::hptdcErrorName(error_index - 1));
0482   HPTDCErrorFlags_2D->setBinLabel(16, "Wrong EC");
0483 
0484   int tmpIndex = 0;
0485   HPTDCErrorFlags_2D->setBinLabel(++tmpIndex, "DB 0 TDC 18", /* axis */ 2);
0486   HPTDCErrorFlags_2D->setBinLabel(++tmpIndex, "DB 0 TDC 17", /* axis */ 2);
0487   HPTDCErrorFlags_2D->setBinLabel(++tmpIndex, "DB 0 TDC 16", /* axis */ 2);
0488   HPTDCErrorFlags_2D->setBinLabel(++tmpIndex, "DB 0 TDC 15", /* axis */ 2);
0489   HPTDCErrorFlags_2D->setBinLabel(++tmpIndex, "DB 1 TDC 18", /* axis */ 2);
0490   HPTDCErrorFlags_2D->setBinLabel(++tmpIndex, "DB 1 TDC 17", /* axis */ 2);
0491   HPTDCErrorFlags_2D->setBinLabel(++tmpIndex, "DB 1 TDC 16", /* axis */ 2);
0492   HPTDCErrorFlags_2D->setBinLabel(++tmpIndex, "DB 1 TDC 15", /* axis */ 2);
0493 
0494   MHComprensive =
0495       ibooker.book2D("MH in channels", title + " MH (%) in channels;plane number;ch number", 10, -0.5, 4.5, 14, -1, 13);
0496 
0497   // ibooker.setCurrentFolder( path+"/clock/" );
0498   // clock_Digi1_le = ibooker.book1D( "clock1 leading edge", title+" clock1;leading edge (ns)", 250, 0, 25 );
0499   // clock_Digi1_te = ibooker.book1D( "clock1 trailing edge", title+" clock1;trailing edge (ns)", 75, 0, 75 );
0500   // clock_Digi3_le = ibooker.book1D( "clock3 leading edge", title+" clock3;leading edge (ns)", 250, 0, 25 );
0501   // clock_Digi3_te = ibooker.book1D( "clock3 trailing edge", title+" clock3;trailing edge (ns)", 75, 0, 75 );
0502 }
0503 
0504 //----------------------------------------------------------------------------------------------------
0505 
0506 CTPPSDiamondDQMSource::PlanePlots::PlanePlots(DQMStore::IBooker& ibooker, unsigned int id, unsigned int windowsNum)
0507     : pixelTracksMapWithDiamonds("Pixel track maps for efficiency with coincidence",
0508                                  "Pixel track maps for efficiency with coincidence",
0509                                  25,
0510                                  0,
0511                                  25,
0512                                  24,
0513                                  -6,
0514                                  18) {
0515   std::string path, title;
0516   CTPPSDiamondDetId(id).planeName(path, CTPPSDiamondDetId::nPath);
0517   ibooker.setCurrentFolder(path);
0518 
0519   CTPPSDiamondDetId(id).planeName(title, CTPPSDiamondDetId::nFull);
0520 
0521   digiProfileCumulativePerPlane = ibooker.book1D("digi profile",
0522                                                  title + " digi profile; ch number",
0523                                                  CTPPS_DIAMOND_NUM_OF_CHANNELS,
0524                                                  -0.5,
0525                                                  CTPPS_DIAMOND_NUM_OF_CHANNELS - 0.5);
0526   hitProfile = ibooker.book1D(
0527       "hit profile", title + " hit profile;x (mm)", 19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM, -0.5, 18.5);
0528   hit_multiplicity = ibooker.book1D("channels per plane", title + " channels per plane; ch per plane", 13, -0.5, 12.5);
0529 
0530   leadingEdgeCumulativePerPlane_both = ibooker.book1D("leading edge (le and te)",
0531                                                       title + " leading edge (le and te) (recHits); leading edge (ns)",
0532                                                       25 * windowsNum,
0533                                                       0,
0534                                                       25 * windowsNum);
0535   leadingEdgeCumulativePerPlane_le = ibooker.book1D("leading edge (le only)",
0536                                                     title + " leading edge (le only) (DIGIs); leading edge (ns)",
0537                                                     25 * windowsNum,
0538                                                     0,
0539                                                     25 * windowsNum);
0540   trailingEdgeCumulativePerPlane_te = ibooker.book1D("trailing edge (te only)",
0541                                                      title + " trailing edge (te only) (DIGIs); trailing edge (ns)",
0542                                                      25 * windowsNum,
0543                                                      0,
0544                                                      25 * windowsNum);
0545   TimeOverThresholdCumulativePerPlane =
0546       ibooker.book1D("time over threshold", title + " time over threshold;time over threshold (ns)", 75, -25, 50);
0547 
0548   pixelTomography_far = ibooker.book2D("tomography pixel",
0549                                        title + " tomography with pixel;x + 25 OOT (mm);y (mm)",
0550                                        25 * windowsNum,
0551                                        0,
0552                                        25 * windowsNum,
0553                                        16,
0554                                        -8,
0555                                        8);
0556   EfficiencyWRTPixelsInPlane =
0557       ibooker.book2D("Efficiency wrt pixels", title + " Efficiency wrt pixels;x (mm);y (mm)", 25, 0, 25, 24, -6, 18);
0558 }
0559 
0560 //----------------------------------------------------------------------------------------------------
0561 
0562 CTPPSDiamondDQMSource::ChannelPlots::ChannelPlots(DQMStore::IBooker& ibooker, unsigned int id, unsigned int windowsNum)
0563     : HitCounter(0), MHCounter(0), LeadingOnlyCounter(0), TrailingOnlyCounter(0), CompleteCounter(0) {
0564   std::string path, title;
0565   CTPPSDiamondDetId(id).channelName(path, CTPPSDiamondDetId::nPath);
0566   ibooker.setCurrentFolder(path);
0567 
0568   CTPPSDiamondDetId(id).channelName(title, CTPPSDiamondDetId::nFull);
0569 
0570   leadingWithoutTrailing = ibooker.book1D("event category", title + " Event Category;;%", 3, 0.5, 3.5);
0571   leadingWithoutTrailing->setBinLabel(1, "Leading only");
0572   leadingWithoutTrailing->setBinLabel(2, "Trailing only");
0573   leadingWithoutTrailing->setBinLabel(3, "Full");
0574 
0575   for (unsigned int i = 0; i < windowsNum; i++) {
0576     std::string window = std::to_string(i * 25) + "-" + std::to_string((i + 1) * 25);
0577     activity_per_bx[i] = ibooker.book1D(
0578         "activity per BX " + window, title + " Activity per BX " + window + " ns;Event.BX", 3600, -1.5, 3598. + 0.5);
0579   }
0580 
0581   HPTDCErrorFlags = ibooker.book1D("hptdc_Errors", title + " HPTDC Errors", 16, -0.5, 16.5);
0582   for (unsigned short error_index = 1; error_index < 16; ++error_index)
0583     HPTDCErrorFlags->setBinLabel(error_index, HPTDCErrorFlags::hptdcErrorName(error_index - 1));
0584   HPTDCErrorFlags->setBinLabel(16, "MH  (%)");
0585 
0586   leadingEdgeCumulative_both = ibooker.book1D("leading edge (le and te)",
0587                                               title + " leading edge (le and te) (recHits); leading edge (ns)",
0588                                               25 * windowsNum,
0589                                               0,
0590                                               25 * windowsNum);
0591   leadingEdgeCumulative_le = ibooker.book1D("leading edge (le only)",
0592                                             title + " leading edge (le only) (DIGIs); leading edge (ns)",
0593                                             25 * windowsNum,
0594                                             0,
0595                                             25 * windowsNum);
0596   trailingEdgeCumulative_te = ibooker.book1D("trailing edge (te only)",
0597                                              title + " trailing edge (te only) (DIGIs); trailing edge (ns)",
0598                                              25 * windowsNum,
0599                                              0,
0600                                              25 * windowsNum);
0601   TimeOverThresholdCumulativePerChannel =
0602       ibooker.book1D("time over threshold", title + " time over threshold;time over threshold (ns)", 75, -25, 50);
0603   // LeadingTrailingCorrelationPerChannel =
0604   //     ibooker.book2D("leading trailing correlation",
0605   //                    title + " leading trailing correlation;leading edge (ns);trailing edge (ns)",
0606   //                    75,
0607   //                    0,
0608   //                    75,
0609   //                    75,
0610   //                    0,
0611   //                    75);
0612 
0613   pixelTomography_far = ibooker.book2D("tomography pixel",
0614                                        "tomography with pixel;x + 25 OOT (mm);y (mm)",
0615                                        25 * windowsNum,
0616                                        0,
0617                                        25 * windowsNum,
0618                                        16,
0619                                        -8,
0620                                        8);
0621 
0622   hit_rate = ibooker.book1D("hit rate", title + "hit rate;rate (Hz)", 40, 0, 20);
0623 
0624   recHitTime = ibooker.book1D("recHit Time", title + " recHit Time; t (ns)", 500, -25, 25);
0625 }
0626 
0627 //----------------------------------------------------------------------------------------------------
0628 
0629 CTPPSDiamondDQMSource::CTPPSDiamondDQMSource(const edm::ParameterSet& ps)
0630     : perLSsaving_(ps.getUntrackedParameter<bool>("perLSsaving", false)),
0631       tokenPixelTrack_(consumes<edm::DetSetVector<CTPPSPixelLocalTrack>>(
0632           ps.getUntrackedParameter<edm::InputTag>("tagPixelLocalTracks"))),
0633       tokenDiamondHit_(consumes<edm::DetSetVector<CTPPSDiamondRecHit>>(
0634           ps.getUntrackedParameter<edm::InputTag>("tagDiamondRecHits"))),
0635       tokenDiamondTrack_(consumes<edm::DetSetVector<CTPPSDiamondLocalTrack>>(
0636           ps.getUntrackedParameter<edm::InputTag>("tagDiamondLocalTracks"))),
0637       ctppsGeometryRunToken_(esConsumes<CTPPSGeometry, VeryForwardRealGeometryRecord, edm::Transition::BeginRun>()),
0638       ctppsGeometryEventToken_(esConsumes<CTPPSGeometry, VeryForwardRealGeometryRecord>()),
0639       excludeMultipleHits_(ps.getParameter<bool>("excludeMultipleHits")),
0640       extract_digi_info_(ps.getParameter<bool>("extractDigiInfo")),
0641       centralOOT_(-999),
0642       verbosity_(ps.getUntrackedParameter<unsigned int>("verbosity", 0)),
0643       plotOnline_(ps.getUntrackedParameter<bool>("plotOnline", true)),
0644       plotOffline_(ps.getUntrackedParameter<bool>("plotOffline", false)),
0645       windowsNum_(ps.getUntrackedParameter<unsigned int>("windowsNum", 3)),
0646       trackCorrelationThreshold_(ps.getUntrackedParameter<unsigned int>("trackCorrelationThreshold", 3)),
0647       EC_difference_56_(-500),
0648       EC_difference_45_(-500) {
0649   if (extract_digi_info_) {
0650     tokenStatus_ = consumes<edm::DetSetVector<TotemVFATStatus>>(ps.getUntrackedParameter<edm::InputTag>("tagStatus"));
0651     tokenFEDInfo_ = consumes<std::vector<TotemFEDInfo>>(ps.getUntrackedParameter<edm::InputTag>("tagFEDInfo"));
0652     tokenDigi_ = consumes<edm::DetSetVector<CTPPSDiamondDigi>>(ps.getUntrackedParameter<edm::InputTag>("tagDigi"));
0653   }
0654   for (const auto& pset : ps.getParameter<std::vector<edm::ParameterSet>>("offsetsOOT")) {
0655     runParameters_.emplace_back(
0656         std::make_pair(pset.getParameter<edm::EventRange>("validityRange"), pset.getParameter<int>("centralOOT")));
0657   }
0658 }
0659 
0660 //----------------------------------------------------------------------------------------------------
0661 
0662 void CTPPSDiamondDQMSource::dqmBeginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
0663   centralOOT_ = -999;
0664   for (const auto& oot : runParameters_) {
0665     if (edm::contains(oot.first, edm::EventID(iRun.run(), 0, 1))) {
0666       centralOOT_ = oot.second;
0667       break;
0668     }
0669   }
0670 
0671   // Get detector shifts from the geometry
0672   const CTPPSGeometry& geom = iSetup.getData(ctppsGeometryRunToken_);
0673   for (auto it = geom.beginRP(); it != geom.endRP(); ++it)
0674     if (CTPPSDiamondDetId::check(it->first)) {
0675       const CTPPSDiamondDetId diam_id(it->first);
0676       const auto diam = geom.sensor(it->first);
0677       diamShifts_[diam_id].global = diam->translation().x() - diam->getDiamondDimensions().xHalfWidth;
0678       diamRotations_[diam_id] = diam->rotation();
0679       diamTranslations_[diam_id] = diam->translation();
0680       diamHalfWidths_[diam_id] = diam->getDiamondDimensions().xHalfWidth;
0681 
0682       if (iRun.run() > FIRST_RUN_W_PIXELS) {  // pixel installed
0683         const CTPPSPixelDetId pixid(diam_id.arm(), CTPPS_PIXEL_STATION_ID, CTPPS_PIXEL_FAR_RP_ID);
0684         auto pix = geom.sensor(pixid);
0685         // Rough alignement of pixel detector for diamond tomography
0686         diamShifts_[diam_id].withPixels =
0687             pix->translation().x() - pix->getDiamondDimensions().xHalfWidth - diamShifts_[diam_id].global - 1.;
0688       }
0689     }
0690 }
0691 
0692 //----------------------------------------------------------------------------------------------------
0693 
0694 void CTPPSDiamondDQMSource::bookHistograms(DQMStore::IBooker& ibooker, const edm::Run&, const edm::EventSetup& iSetup) {
0695   ibooker.cd();
0696   ibooker.setCurrentFolder("CTPPS");
0697 
0698   globalPlot_ = GlobalPlots(ibooker);
0699 
0700   // book plots from the geometry
0701   const CTPPSGeometry& geom = iSetup.getData(ctppsGeometryRunToken_);
0702   for (auto it = geom.beginSensor(); it != geom.endSensor(); ++it) {
0703     if (!CTPPSDiamondDetId::check(it->first))
0704       continue;
0705     // per-channel plots
0706     const CTPPSDiamondDetId chId(it->first);
0707     if (plotOnline_ && channelPlots_.count(chId) == 0)
0708       channelPlots_[chId] = ChannelPlots(ibooker, chId, windowsNum_);
0709 
0710     // per-plane plots
0711     const CTPPSDiamondDetId plId(chId.planeId());
0712     if (planePlots_.count(plId) == 0)
0713       planePlots_[plId] = PlanePlots(ibooker, plId, windowsNum_);
0714     // per-pot plots
0715     const CTPPSDiamondDetId rpId(chId.rpId());
0716     if (potPlots_.count(rpId) == 0)
0717       potPlots_[rpId] = PotPlots(ibooker, rpId, windowsNum_, plotOnline_, plotOffline_, perLSsaving_);
0718 
0719     // per-sector plots
0720     const CTPPSDiamondDetId secId(chId.armId());
0721     if (plotOffline_ && sectorPlots_.count(secId) == 0)
0722       sectorPlots_[secId] = SectorPlots(ibooker, secId, plotOnline_);
0723   }
0724 }
0725 
0726 //----------------------------------------------------------------------------------------------------
0727 
0728 std::shared_ptr<dds::Cache> CTPPSDiamondDQMSource::globalBeginLuminosityBlock(const edm::LuminosityBlock&,
0729                                                                               const edm::EventSetup&) const {
0730   auto d = std::make_shared<dds::Cache>();
0731   d->hitDistribution2dMap.reserve(potPlots_.size());
0732   if (!perLSsaving_ && plotOnline_) {
0733     for (auto& plot : potPlots_) {
0734       d->hitDistribution2dMap[plot.first] = std::make_unique<TH2F>(
0735           "hits in planes lumisection",
0736           (std::string(plot.second.hitDistribution2d_lumisection->getTH2F()->GetTitle()) + ";plane number;x (mm)")
0737               .c_str(),
0738           10,
0739           -0.5,
0740           4.5,
0741           19. * INV_DISPLAY_RESOLUTION_FOR_HITS_MM,
0742           -0.5,
0743           18.5);
0744     }
0745   }
0746   return d;
0747 }
0748 
0749 //----------------------------------------------------------------------------------------------------
0750 
0751 void CTPPSDiamondDQMSource::analyze(const edm::Event& event, const edm::EventSetup& iSetup) {
0752   // get event data
0753 
0754   edm::Handle<edm::DetSetVector<TotemVFATStatus>> diamondVFATStatus;
0755   edm::Handle<edm::DetSetVector<CTPPSDiamondDigi>> diamondDigis;
0756   edm::Handle<std::vector<TotemFEDInfo>> fedInfo;
0757   if (extract_digi_info_) {
0758     event.getByToken(tokenStatus_, diamondVFATStatus);
0759     event.getByToken(tokenDigi_, diamondDigis);
0760     event.getByToken(tokenFEDInfo_, fedInfo);
0761   }
0762 
0763   edm::Handle<edm::DetSetVector<CTPPSPixelLocalTrack>> pixelTracks;
0764   event.getByToken(tokenPixelTrack_, pixelTracks);
0765 
0766   edm::Handle<edm::DetSetVector<CTPPSDiamondRecHit>> diamondRecHits;
0767   event.getByToken(tokenDiamondHit_, diamondRecHits);
0768 
0769   edm::Handle<edm::DetSetVector<CTPPSDiamondLocalTrack>> diamondLocalTracks;
0770   event.getByToken(tokenDiamondTrack_, diamondLocalTracks);
0771 
0772   const CTPPSGeometry* ctppsGeometry = &iSetup.getData(ctppsGeometryEventToken_);
0773 
0774   // check validity
0775   bool valid = true;
0776   if (extract_digi_info_) {  // drop DIGI-level validity checks if not monitored
0777     valid &= diamondVFATStatus.isValid();
0778     valid &= diamondDigis.isValid();
0779     valid &= fedInfo.isValid();
0780   }
0781   valid &= pixelTracks.isValid();
0782   valid &= diamondRecHits.isValid();
0783   valid &= diamondLocalTracks.isValid();
0784 
0785   if (!valid) {
0786     if (verbosity_)
0787       edm::LogProblem("CTPPSDiamondDQMSource")
0788           << "ERROR in CTPPSDiamondDQMSource::analyze > some of the required inputs are not valid. Skipping this "
0789              "event.\n"
0790           << "  DIGI-level: (checked? " << std::boolalpha << extract_digi_info_ << ")\n"
0791           << "    diamondVFATStatus.isValid = " << diamondVFATStatus.isValid() << "\n"
0792           << "    diamondDigis.isValid = " << diamondDigis.isValid() << "\n"
0793           << "    fedInfo.isValid = " << fedInfo.isValid() << "\n"
0794           << "  RECO-level:\n"
0795           << "    pixelTracks.isValid = " << pixelTracks.isValid() << "\n"
0796           << "    diamondRecHits.isValid = " << diamondRecHits.isValid() << "\n"
0797           << "    diamondLocalTracks.isValid = " << diamondLocalTracks.isValid();
0798     return;
0799   }
0800 
0801   //------------------------------
0802   // Sector Plots
0803   //------------------------------
0804   // Using CTPPSDiamondLocalTrack
0805   if (plotOffline_)
0806 
0807     // diamond timing detectors are located in:
0808     //   - 220cyl: "cylindrical pot" station  (id=1), in horizontal Roman Pot with id=6 ("far")
0809     //   - 220nr_hr: "220m"  station (id=2), in horizontal Roman Pot with id=2 ("near horizontal")
0810 
0811     for (const auto& tracks_220nr_hr : *diamondLocalTracks) {
0812       // to preprare correlation plot, we need to select tracks from nr_hr pot in 220m station
0813       const CTPPSDiamondDetId detId_220nr_hr(tracks_220nr_hr.detId());
0814 
0815       // selecting only tracks from 220nr station, realised as skipping tracks from 220cyl station
0816       if ((detId_220nr_hr.rp() == CTPPS_DIAMOND_CYL_RP_ID) &&
0817           (detId_220nr_hr.station() == CTPPS_DIAMOND_CYL_STATION_ID))
0818         continue;
0819 
0820       if (potPlots_.count(detId_220nr_hr.rpId()) == 0)
0821         continue;
0822       TH1F* trackHistoInTimeTmp = potPlots_[detId_220nr_hr.rpId()].trackDistribution->getTH1F();
0823 
0824       for (const auto& track_220nr_hr : tracks_220nr_hr) {
0825         if (!track_220nr_hr.isValid())
0826           continue;
0827 
0828         // get the bins from per-pot plots
0829         int startBin_220nr_hr = trackHistoInTimeTmp->FindBin(
0830             track_220nr_hr.x0() - diamShifts_[detId_220nr_hr.rpId()].global - track_220nr_hr.x0Sigma());
0831         int numOfBins_220nr_hr = 2 * track_220nr_hr.x0Sigma() * INV_DISPLAY_RESOLUTION_FOR_HITS_MM;
0832 
0833         for (const auto& tracks_220cyl : *diamondLocalTracks) {
0834           CTPPSDiamondDetId detId_220cyl(tracks_220cyl.detId());
0835 
0836           // select tracks in the same arm, but belonging to the cylindrical pot
0837           // that means skipping tracks from the opposite arm (skip if 220nr_hr.arm != 220cyl.arm)
0838           // and skipping tracks from the 220nr_hr pot  (skip if 220nr_hr == 220cyl.station)
0839           if (detId_220nr_hr.arm() != detId_220cyl.arm() || detId_220nr_hr.station() == detId_220cyl.station())
0840             continue;
0841 
0842           if (sectorPlots_.count(detId_220cyl.armId()) == 0)
0843             continue;
0844 
0845           TH2F* trackHistoTmp = sectorPlots_[detId_220cyl.armId()].trackCorrelation->getTH2F();
0846           TAxis* trackHistoTmpXAxis = trackHistoTmp->GetXaxis();
0847           TAxis* trackHistoTmpYAxis = trackHistoTmp->GetYaxis();
0848 
0849           for (const auto& track_220cyl : tracks_220cyl) {
0850             if (!track_220cyl.isValid())
0851               continue;
0852             int startBin_220cyl = trackHistoTmpYAxis->FindBin(
0853                 track_220cyl.x0() - diamShifts_[detId_220cyl.rpId()].global - track_220cyl.x0Sigma());
0854             int numOfBins_220cyl = 2 * track_220cyl.x0Sigma() * INV_DISPLAY_RESOLUTION_FOR_HITS_MM;
0855 
0856             // fill the correlation plot
0857             for (int i = 0; i < numOfBins_220nr_hr; ++i)
0858               for (int y = 0; y < numOfBins_220cyl; ++y) {
0859                 float track_220nr_hr_x = trackHistoTmpXAxis->GetBinCenter(startBin_220nr_hr + i);
0860                 float track_220cyl_x = trackHistoTmpYAxis->GetBinCenter(startBin_220cyl + y);
0861                 trackHistoTmp->Fill(track_220nr_hr_x, track_220cyl_x);
0862                 // fill low multiplicity histogram
0863                 if (tracks_220nr_hr.size() < 3 && tracks_220cyl.size() < trackCorrelationThreshold_)
0864                   sectorPlots_[detId_220cyl.armId()].trackCorrelationLowMultiplicity->Fill(track_220nr_hr_x,
0865                                                                                            track_220cyl_x);
0866               }
0867           }
0868         }
0869       }
0870     }
0871 
0872   //------------------------------
0873   // RP Plots
0874   //------------------------------
0875 
0876   //------------------------------
0877   // Correlation Plots
0878   //------------------------------
0879 
0880   if (extract_digi_info_) {
0881     // Using CTPPSDiamondDigi
0882     for (const auto& digis : *diamondDigis) {
0883       const CTPPSDiamondDetId detId(digis.detId()), detId_pot(detId.rpId());
0884       for (const auto& digi : digis) {
0885         if (detId.channel() == CHANNEL_OF_VFAT_CLOCK)
0886           continue;
0887         if (potPlots_.count(detId_pot) == 0)
0888           continue;
0889         //Leading without trailing investigation
0890         if (digi.leadingEdge() != 0 || digi.trailingEdge() != 0) {
0891           ++potPlots_[detId_pot].HitCounter;
0892           if (digi.leadingEdge() != 0) {
0893             potPlots_[detId_pot].leadingEdgeCumulative_all->Fill(HPTDC_BIN_WIDTH_NS * digi.leadingEdge());
0894           }
0895           if (digi.leadingEdge() != 0 && digi.trailingEdge() == 0) {
0896             ++potPlots_[detId_pot].LeadingOnlyCounter;
0897             potPlots_[detId_pot].leadingEdgeCumulative_le->Fill(HPTDC_BIN_WIDTH_NS * digi.leadingEdge());
0898           }
0899           if (digi.leadingEdge() == 0 && digi.trailingEdge() != 0) {
0900             ++potPlots_[detId_pot].TrailingOnlyCounter;
0901             potPlots_[detId_pot].trailingEdgeCumulative_te->Fill(HPTDC_BIN_WIDTH_NS * digi.trailingEdge());
0902           }
0903           if (digi.leadingEdge() != 0 && digi.trailingEdge() != 0) {
0904             ++potPlots_[detId_pot].CompleteCounter;
0905             // potPlots_[detId_pot].leadingTrailingCorrelationPot->Fill(HPTDC_BIN_WIDTH_NS * digi.leadingEdge(),
0906             //                                                          HPTDC_BIN_WIDTH_NS * digi.trailingEdge());
0907           }
0908         }
0909 
0910         // HPTDC Errors
0911         const HPTDCErrorFlags hptdcErrors = digi.hptdcErrorFlags();
0912         if (detId.channel() == HPTDC_0_CHANNEL ||
0913             detId.channel() == HPTDC_1_CHANNEL) {  // ch6 for HPTDC 0 and ch7 for HPTDC 1
0914           int verticalIndex = 2 * detId.plane() + (detId.channel() - HPTDC_0_CHANNEL);
0915           for (unsigned short hptdcErrorIndex = 1; hptdcErrorIndex < 16; ++hptdcErrorIndex)
0916             if (hptdcErrors.errorId(hptdcErrorIndex - 1))
0917               potPlots_[detId_pot].HPTDCErrorFlags_2D->Fill(hptdcErrorIndex, verticalIndex);
0918         }
0919         if (digi.multipleHit())
0920           ++potPlots_[detId_pot].MHCounter;
0921       }
0922     }
0923   }
0924 
0925   // EC Errors
0926   if (extract_digi_info_) {
0927     for (const auto& vfat_status : *diamondVFATStatus) {
0928       const CTPPSDiamondDetId detId(vfat_status.detId());
0929       for (const auto& status : vfat_status) {
0930         if (!status.isOK())
0931           continue;
0932         if (potPlots_.count(detId.rpId()) == 0)
0933           continue;
0934         if (channelPlots_.count(detId) == 0)
0935           continue;
0936 
0937         // Check Event Number
0938         for (const auto& optorx : *fedInfo) {
0939           if (detId.arm() == 1 && optorx.fedId() == CTPPS_FED_ID_56)
0940             checkEventNumber(detId, optorx, status, potPlots_[detId.rpId()], EC_difference_56_);
0941           else if (detId.arm() == 0 && optorx.fedId() == CTPPS_FED_ID_45)
0942             checkEventNumber(detId, optorx, status, potPlots_[detId.rpId()], EC_difference_45_);
0943         }
0944       }
0945     }
0946   }
0947 
0948   // Using CTPPSDiamondRecHit
0949   std::unordered_map<unsigned int, std::set<unsigned int>> planes;
0950   std::unordered_map<unsigned int, std::set<unsigned int>> planes_inclusive;
0951 
0952   auto lumiCache = luminosityBlockCache(event.getLuminosityBlock().index());
0953   for (const auto& rechits : *diamondRecHits) {
0954     const CTPPSDiamondDetId detId(rechits.detId()), detId_pot(detId.rpId());
0955 
0956     for (const auto& rechit : rechits) {
0957       planes_inclusive[detId_pot].insert(detId.plane());
0958       if (excludeMultipleHits_ && rechit.multipleHits() > 0)
0959         continue;
0960       if (rechit.toT() != 0 && centralOOT_ != -999 && rechit.ootIndex() == centralOOT_)
0961         planes[detId_pot].insert(detId.plane());
0962 
0963       if (potPlots_.count(detId_pot) == 0)
0964         continue;
0965 
0966       potPlots_[detId_pot].recHitTime->Fill(rechit.time());
0967       // coordinates of the rechit in the local coordinate system are needed for the rotated pot
0968       auto localRecHit = transformToLocalCoordinates(detId_pot, rechit.x(), rechit.y(), rechit.z());
0969       auto localRecHitLeftX = localRecHit.x() + diamHalfWidths_.at(detId_pot) - 0.5 * rechit.xWidth();
0970 
0971       float UFSDShift = 0.0;
0972       if (rechit.yWidth() < 3)
0973         UFSDShift = 0.5;  // Display trick for UFSD that have 2 pixels with same X
0974 
0975       if (rechit.toT() != 0 && centralOOT_ != -999 && rechit.ootIndex() == centralOOT_) {
0976         TH2F* hitHistoTmp = potPlots_[detId_pot].hitDistribution2d->getTH2F();
0977         TAxis* hitHistoTmpYAxis = hitHistoTmp->GetYaxis();
0978         int startBin = hitHistoTmpYAxis->FindBin(localRecHitLeftX);
0979         int numOfBins = rechit.xWidth() * INV_DISPLAY_RESOLUTION_FOR_HITS_MM;
0980         for (int i = 0; i < numOfBins; ++i)
0981           hitHistoTmp->Fill(detId.plane() + UFSDShift, hitHistoTmpYAxis->GetBinCenter(startBin + i));
0982 
0983         if (!perLSsaving_ && plotOnline_) {
0984           hitHistoTmp = lumiCache->hitDistribution2dMap[detId_pot].get();
0985           hitHistoTmpYAxis = hitHistoTmp->GetYaxis();
0986           startBin = hitHistoTmpYAxis->FindBin(localRecHitLeftX);
0987           numOfBins = rechit.xWidth() * INV_DISPLAY_RESOLUTION_FOR_HITS_MM;
0988           for (int i = 0; i < numOfBins; ++i)
0989             hitHistoTmp->Fill(detId.plane() + UFSDShift, hitHistoTmpYAxis->GetBinCenter(startBin + i));
0990         }
0991       }
0992 
0993       if (rechit.toT() > 0) {
0994         // Both
0995         potPlots_[detId_pot].leadingEdgeCumulative_both->Fill(rechit.time() + 25 * rechit.ootIndex());
0996         potPlots_[detId_pot].timeOverThresholdCumulativePot->Fill(rechit.toT());
0997 
0998         TH2F* hitHistoOOTTmp = potPlots_[detId_pot].hitDistribution2dOOT->getTH2F();
0999         TAxis* hitHistoOOTTmpYAxis = hitHistoOOTTmp->GetYaxis();
1000         int startBin = hitHistoOOTTmpYAxis->FindBin(localRecHitLeftX);
1001         int numOfBins = rechit.xWidth() * INV_DISPLAY_RESOLUTION_FOR_HITS_MM;
1002         for (int i = 0; i < numOfBins; ++i)
1003           hitHistoOOTTmp->Fill(detId.plane() + 1. / windowsNum_ * rechit.ootIndex(),
1004                                hitHistoOOTTmpYAxis->GetBinCenter(startBin + i));
1005 
1006       } else if (rechit.ootIndex() != CTPPSDiamondRecHit::TIMESLICE_WITHOUT_LEADING && plotOnline_) {
1007         // Only leading
1008         TH2F* hitHistoOOTTmp = potPlots_[detId_pot].hitDistribution2dOOT_le->getTH2F();
1009         TAxis* hitHistoOOTTmpYAxis = hitHistoOOTTmp->GetYaxis();
1010         int startBin = hitHistoOOTTmpYAxis->FindBin(localRecHitLeftX);
1011         int numOfBins = rechit.xWidth() * INV_DISPLAY_RESOLUTION_FOR_HITS_MM;
1012         for (int i = 0; i < numOfBins; ++i)
1013           hitHistoOOTTmp->Fill(detId.plane() + 1. / windowsNum_ * rechit.ootIndex(),
1014                                hitHistoOOTTmpYAxis->GetBinCenter(startBin + i));
1015       }
1016       if (rechit.ootIndex() != CTPPSDiamondRecHit::TIMESLICE_WITHOUT_LEADING &&
1017           potPlots_[detId_pot].activity_per_bx.count(rechit.ootIndex()) > 0)
1018         potPlots_[detId_pot].activity_per_bx.at(rechit.ootIndex())->Fill(event.bunchCrossing());
1019     }
1020   }
1021 
1022   for (const auto& plt : potPlots_) {
1023     plt.second.activePlanes->Fill(planes[plt.first].size());
1024     if (plotOnline_)
1025       plt.second.activePlanesInclusive->Fill(planes_inclusive[plt.first].size());
1026   }
1027 
1028   // Using CTPPSDiamondLocalTrack
1029   for (const auto& tracks : *diamondLocalTracks) {
1030     const CTPPSDiamondDetId detId(tracks.detId()), detId_pot(detId.rpId());
1031     const auto& x_shift = diamShifts_.at(detId_pot);
1032     // extract the rotation angle for the diamond pot
1033     auto cosRotAngle = diamRotations_.at(detId_pot) * CTPPSGeometry::Vector(1, 0, 0);
1034 
1035     for (const auto& track : tracks) {
1036       if (!track.isValid())
1037         continue;
1038       if (track.ootIndex() == CTPPSDiamondRecHit::TIMESLICE_WITHOUT_LEADING)
1039         continue;
1040       if (excludeMultipleHits_ && track.multipleHits() > 0)
1041         continue;
1042       if (potPlots_.count(detId_pot) == 0)
1043         continue;
1044 
1045       TH2F* trackHistoOOTTmp = potPlots_[detId_pot].trackDistributionOOT->getTH2F();
1046       TAxis* trackHistoOOTTmpYAxis = trackHistoOOTTmp->GetYaxis();
1047       int startBin = trackHistoOOTTmpYAxis->FindBin(track.x0() - x_shift.global - track.x0Sigma());
1048       int numOfBins = 2 * track.x0Sigma() * INV_DISPLAY_RESOLUTION_FOR_HITS_MM;
1049       for (int i = 0; i < numOfBins; ++i)
1050         trackHistoOOTTmp->Fill(track.ootIndex(), trackHistoOOTTmpYAxis->GetBinCenter(startBin + i));
1051 
1052       if (centralOOT_ != -999 && track.ootIndex() == centralOOT_) {
1053         TH1F* trackHistoInTimeTmp = potPlots_[detId_pot].trackDistribution->getTH1F();
1054 
1055         // X coordinate of the left edge of the track in the local coordinate system
1056         auto localTrackX =
1057             (track.x0() - diamTranslations_.at(detId_pot).x() + diamHalfWidths_.at(detId_pot) - track.x0Sigma()) /
1058             cosRotAngle.x();
1059 
1060         int startBin = trackHistoInTimeTmp->FindBin((localTrackX));
1061         int numOfBins = 2 * track.x0Sigma() * INV_DISPLAY_RESOLUTION_FOR_HITS_MM;
1062         for (int i = 0; i < numOfBins; ++i)
1063           trackHistoInTimeTmp->Fill(trackHistoInTimeTmp->GetBinCenter(startBin + i));
1064       }
1065       if (plotOffline_ && !perLSsaving_) {
1066         // potPlots_[detId_pot].trackTimeVsLS->Fill(event.luminosityBlock(),track.time());
1067         potPlots_[detId_pot].trackTimeVsBX->Fill(event.bunchCrossing(), track.time());
1068       }
1069     }
1070   }
1071 
1072   // Channel efficiency using CTPPSDiamondLocalTrack
1073   for (const auto& tracks : *diamondLocalTracks) {
1074     const CTPPSDiamondDetId detId(tracks.detId()), detId_pot(detId.rpId());
1075     for (const auto& track : tracks) {
1076       // Find hits and planes in the track
1077       int numOfHits = 0;
1078       std::set<int> planesInTrackSet;
1079       for (const auto& vec : *diamondRecHits) {
1080         const CTPPSDiamondDetId detid(vec.detId());
1081         if (detid.arm() != detId_pot.arm())
1082           continue;
1083 
1084         for (const auto& hit : vec) {
1085           // first check if the hit contributes to the track
1086           if (track.containsHit(hit)) {
1087             ++numOfHits;
1088             planesInTrackSet.insert(detid.plane());
1089           }
1090         }
1091       }
1092 
1093       if (numOfHits > 0 && numOfHits <= 10 && planesInTrackSet.size() > 2) {
1094         for (const auto& plt_vs_ch : channelPlots_) {  // loop over all channels registered in geometry
1095           const CTPPSDiamondDetId detId(plt_vs_ch.first);
1096           if (detId.rpId() != detId_pot)
1097             continue;
1098           const unsigned short map_index = detId.plane() * 100 + detId.channel();
1099           if (potPlots_[detId_pot].effDoublecountingChMap.count(map_index) == 0) {
1100             potPlots_[detId_pot].effTriplecountingChMap[map_index] = 0;
1101             potPlots_[detId_pot].effDoublecountingChMap[map_index] = 0;
1102           }
1103           if (channelAlignedWithTrack(ctppsGeometry, detId, track, 0.2)) {
1104             // Channel should fire
1105             ++potPlots_[detId_pot].effDoublecountingChMap[map_index];
1106             for (const auto& rechits : *diamondRecHits) {
1107               const CTPPSDiamondDetId detId_hit(rechits.detId());
1108               if (detId_hit == detId)
1109                 for (const auto& rechit : rechits)
1110                   if (track.containsHit(rechit, 1))  // Channel fired
1111                     ++potPlots_[detId_pot].effTriplecountingChMap[map_index];
1112             }
1113           }
1114         }
1115       }
1116     }
1117   }
1118 
1119   // Tomography of diamonds using pixel
1120   for (const auto& rechits : *diamondRecHits) {
1121     const CTPPSDiamondDetId detId(rechits.detId()), detId_pot(detId.rpId());
1122     const auto pix_shift = diamShifts_.at(detId_pot).withPixels;
1123     for (const auto& rechit : rechits) {
1124       if (excludeMultipleHits_ && rechit.multipleHits() > 0)
1125         continue;
1126       if (rechit.toT() == 0)
1127         continue;
1128       if (!pixelTracks.isValid())
1129         continue;
1130       if (potPlots_.count(detId_pot) == 0)
1131         continue;
1132 
1133       for (const auto& ds : *pixelTracks) {
1134         if (ds.size() > 1)
1135           continue;
1136         const CTPPSPixelDetId pixId(ds.detId());
1137         if (pixId.station() != CTPPS_PIXEL_STATION_ID || pixId.rp() != CTPPS_PIXEL_FAR_RP_ID)
1138           continue;
1139         for (const auto& lt : ds) {
1140           if (lt.isValid() && pixId.arm() == detId_pot.arm()) {
1141             if (rechit.ootIndex() != CTPPSDiamondRecHit::TIMESLICE_WITHOUT_LEADING && rechit.ootIndex() >= 0 &&
1142                 potPlots_[detId_pot].pixelTomographyAll.count(rechit.ootIndex()) > 0 && lt.x0() - pix_shift < 24)
1143               potPlots_[detId_pot]
1144                   .pixelTomographyAll.at(rechit.ootIndex())
1145                   ->Fill(lt.x0() - pix_shift + 25 * detId.plane(), lt.y0());
1146           }
1147         }
1148       }
1149     }
1150   }
1151 
1152   //------------------------------
1153   // Clock Plots
1154   //------------------------------
1155   // Commented out to save space in the DQM files, but code should be kept
1156   // for ( const auto& digis : *diamondDigis ) {
1157   //   const CTPPSDiamondDetId detId(digis.detId()), detId_pot(detId.rpId());
1158   //   if ( detId.channel() == CHANNEL_OF_VFAT_CLOCK ) {
1159   //     for ( const auto& digi : digis ) {
1160   //       if ( digi.leadingEdge() != 0 )  {
1161   //         if ( detId.plane() == 1 ) {
1162   //           potPlots_[detId_pot].clock_Digi1_le->Fill( HPTDC_BIN_WIDTH_NS * digi.leadingEdge() );
1163   //           potPlots_[detId_pot].clock_Digi1_te->Fill( HPTDC_BIN_WIDTH_NS * digi.trailingEdge() );
1164   //         }
1165   //         if ( detId.plane() == 3 ) {
1166   //           potPlots_[detId_pot].clock_Digi3_le->Fill( HPTDC_BIN_WIDTH_NS * digi.leadingEdge() );
1167   //           potPlots_[detId_pot].clock_Digi3_te->Fill( HPTDC_BIN_WIDTH_NS * digi.trailingEdge() );
1168   //         }
1169   //       }
1170   //     }
1171   //   }
1172   // }
1173 
1174   //------------------------------
1175   // Plane Plots
1176   //------------------------------
1177 
1178   // Using CTPPSDiamondDigi
1179   std::unordered_map<unsigned int, unsigned int> channelsPerPlane;
1180   if (extract_digi_info_) {
1181     for (const auto& digis : *diamondDigis) {
1182       const CTPPSDiamondDetId detId(digis.detId()), detId_plane(detId.planeId());
1183       for (const auto& digi : digis) {
1184         if (detId.channel() == CHANNEL_OF_VFAT_CLOCK)
1185           continue;
1186         if (planePlots_.count(detId_plane) == 0)
1187           continue;
1188 
1189         if (digi.leadingEdge() != 0) {
1190           planePlots_[detId_plane].digiProfileCumulativePerPlane->Fill(detId.channel());
1191           channelsPerPlane[detId_plane]++;
1192         }
1193 
1194         // Check dropped trailing edges
1195         if ((digi.trailingEdge() == 0) && (digi.leadingEdge() != 0)) {
1196           planePlots_[detId_plane].leadingEdgeCumulativePerPlane_le->Fill(HPTDC_BIN_WIDTH_NS * digi.leadingEdge());
1197         } else if ((digi.leadingEdge() == 0 && (digi.trailingEdge() != 0))) {  // check dropped leading edges
1198           planePlots_[detId_plane].trailingEdgeCumulativePerPlane_te->Fill(HPTDC_BIN_WIDTH_NS * digi.trailingEdge());
1199         }
1200       }
1201     }
1202   }
1203 
1204   for (const auto& plt : channelsPerPlane)
1205     planePlots_[plt.first].hit_multiplicity->Fill(plt.second);
1206 
1207   // Using CTPPSDiamondRecHit
1208   for (const auto& rechits : *diamondRecHits) {
1209     const CTPPSDiamondDetId detId(rechits.detId()), detId_plane(detId.planeId()), detId_pot(detId.rpId());
1210     for (const auto& rechit : rechits) {
1211       if (excludeMultipleHits_ && rechit.multipleHits() > 0)
1212         continue;
1213       if (rechit.toT() == 0)
1214         continue;
1215 
1216       if (rechit.ootIndex() != CTPPSDiamondRecHit::TIMESLICE_WITHOUT_LEADING) {
1217         planePlots_[detId_plane].leadingEdgeCumulativePerPlane_both->Fill(rechit.time() + 25 * rechit.ootIndex());
1218         planePlots_[detId_plane].TimeOverThresholdCumulativePerPlane->Fill(rechit.toT());
1219       }
1220       if (planePlots_.count(detId_plane) != 0) {
1221         if (centralOOT_ != -999 && rechit.ootIndex() == centralOOT_) {
1222           TH1F* hitHistoTmp = planePlots_[detId_plane].hitProfile->getTH1F();
1223           // coordinates of the rechit in the local coordinate system are needed for the rotated pot
1224           auto localRecHit = transformToLocalCoordinates(detId_pot, rechit.x(), rechit.y(), rechit.z());
1225           int startBin = hitHistoTmp->FindBin(localRecHit.x() + diamHalfWidths_.at(detId_pot) - 0.5 * rechit.xWidth());
1226           int numOfBins = rechit.xWidth() * INV_DISPLAY_RESOLUTION_FOR_HITS_MM;
1227           for (int i = 0; i < numOfBins; ++i)
1228             hitHistoTmp->Fill(hitHistoTmp->GetBinCenter(startBin + i));
1229         }
1230       }
1231     }
1232   }
1233 
1234   //Tomography of diamonds using pixel and Efficiency WRT Pixels
1235   for (const auto& ds : *pixelTracks) {
1236     const CTPPSPixelDetId pixId(ds.detId());
1237     if (pixId.station() != CTPPS_PIXEL_STATION_ID || pixId.rp() != CTPPS_PIXEL_FAR_RP_ID)
1238       continue;
1239     if (ds.size() > 1)
1240       continue;
1241     for (const auto& lt : ds) {
1242       if (lt.isValid()) {
1243         for (const auto& sh_vs_id : diamShifts_) {
1244           const CTPPSDiamondDetId& detId_pot = sh_vs_id.first;
1245           const auto& shift = sh_vs_id.second;
1246           if (detId_pot.arm() == pixId.arm())
1247             // For efficieny
1248             potPlots_[detId_pot].pixelTracksMap.Fill(lt.x0() - shift.withPixels, lt.y0());
1249         }
1250 
1251         std::set<CTPPSDiamondDetId> planesWitHits_set;
1252         for (const auto& rechits : *diamondRecHits) {
1253           const CTPPSDiamondDetId detId(rechits.detId()), detId_plane(detId.planeId()), detId_pot(detId.rpId());
1254           const auto pix_shift = diamShifts_.at(detId_pot).withPixels;
1255           for (const auto& rechit : rechits) {
1256             if (excludeMultipleHits_ && rechit.multipleHits() > 0)
1257               continue;
1258             if (rechit.ootIndex() == CTPPSDiamondRecHit::TIMESLICE_WITHOUT_LEADING || rechit.toT() == 0)
1259               continue;
1260             if (planePlots_.count(detId_plane) == 0)
1261               continue;
1262             if (pixId.arm() == detId_plane.arm() && lt.x0() - pix_shift < 24) {
1263               planePlots_[detId_plane].pixelTomography_far->Fill(lt.x0() - pix_shift + 25 * rechit.ootIndex(), lt.y0());
1264               if (centralOOT_ != -999 && rechit.ootIndex() == centralOOT_)
1265                 planesWitHits_set.insert(detId_plane);
1266             }
1267           }
1268         }
1269 
1270         for (auto& planeId : planesWitHits_set)
1271           planePlots_[planeId].pixelTracksMapWithDiamonds.Fill(lt.x0() - diamShifts_.at(planeId.rpId()).withPixels,
1272                                                                lt.y0());
1273       }
1274     }
1275   }
1276 
1277   //------------------------------
1278   // Channel Plots
1279   //------------------------------
1280 
1281   // digi profile cumulative
1282   if (extract_digi_info_) {
1283     for (const auto& digis : *diamondDigis) {
1284       const CTPPSDiamondDetId detId(digis.detId());
1285       for (const auto& digi : digis) {
1286         if (detId.channel() == CHANNEL_OF_VFAT_CLOCK)
1287           continue;
1288         if (channelPlots_.count(detId) != 0) {
1289           // HPTDC Errors
1290           const HPTDCErrorFlags hptdcErrors = digi.hptdcErrorFlags();
1291           for (unsigned short hptdcErrorIndex = 1; hptdcErrorIndex < 16; ++hptdcErrorIndex)
1292             if (hptdcErrors.errorId(hptdcErrorIndex - 1))
1293               channelPlots_[detId].HPTDCErrorFlags->Fill(hptdcErrorIndex);
1294           if (digi.multipleHit())
1295             ++channelPlots_[detId].MHCounter;
1296 
1297           // Check dropped trailing edges
1298           if (digi.leadingEdge() != 0 || digi.trailingEdge() != 0) {
1299             ++channelPlots_[detId].HitCounter;
1300             if (digi.trailingEdge() == 0) {
1301               ++channelPlots_[detId].LeadingOnlyCounter;
1302               channelPlots_[detId].leadingEdgeCumulative_le->Fill(HPTDC_BIN_WIDTH_NS * digi.leadingEdge());
1303             } else if (digi.leadingEdge() == 0) {
1304               ++channelPlots_[detId].TrailingOnlyCounter;
1305               channelPlots_[detId].trailingEdgeCumulative_te->Fill(HPTDC_BIN_WIDTH_NS * digi.trailingEdge());
1306             } else {
1307               ++channelPlots_[detId].CompleteCounter;
1308               // channelPlots_[detId].LeadingTrailingCorrelationPerChannel->Fill(HPTDC_BIN_WIDTH_NS * digi.leadingEdge(),
1309               //                                                                 HPTDC_BIN_WIDTH_NS * digi.trailingEdge());
1310             }
1311           }
1312         }
1313       }
1314     }
1315   }
1316 
1317   // Using CTPPSDiamondRecHit
1318 
1319   for (const auto& rechits : *diamondRecHits) {
1320     CTPPSDiamondDetId detId(rechits.detId());
1321     for (const auto& rechit : rechits) {
1322       if (excludeMultipleHits_ && rechit.multipleHits() > 0)
1323         continue;
1324       if (channelPlots_.count(detId) != 0) {
1325         channelPlots_[detId].recHitTime->Fill(rechit.time());
1326         if (rechit.ootIndex() != CTPPSDiamondRecHit::TIMESLICE_WITHOUT_LEADING && rechit.toT() != 0) {
1327           channelPlots_[detId].leadingEdgeCumulative_both->Fill(rechit.time() + 25 * rechit.ootIndex());
1328           channelPlots_[detId].TimeOverThresholdCumulativePerChannel->Fill(rechit.toT());
1329         }
1330         ++lumiCache->hitsCounterMap[detId];
1331 
1332         if (rechit.ootIndex() != CTPPSDiamondRecHit::TIMESLICE_WITHOUT_LEADING &&
1333             channelPlots_[detId].activity_per_bx.count(rechit.ootIndex()) > 0)
1334           channelPlots_[detId].activity_per_bx.at(rechit.ootIndex())->Fill(event.bunchCrossing());
1335       }
1336     }
1337   }
1338 
1339   // Tomography of diamonds using pixel
1340   for (const auto& rechits : *diamondRecHits) {
1341     const CTPPSDiamondDetId detId(rechits.detId()), detId_pot(detId.rpId());
1342     const auto shift_pix = diamShifts_.at(detId_pot).withPixels;
1343     for (const auto& rechit : rechits) {
1344       if (excludeMultipleHits_ && rechit.multipleHits() > 0)
1345         continue;
1346       if (rechit.ootIndex() == CTPPSDiamondRecHit::TIMESLICE_WITHOUT_LEADING || rechit.toT() == 0)
1347         continue;
1348       if (!pixelTracks.isValid())
1349         continue;
1350       if (channelPlots_.count(detId) == 0)
1351         continue;
1352 
1353       for (const auto& ds : *pixelTracks) {
1354         const CTPPSPixelDetId pixId(ds.detId());
1355         if (pixId.station() != CTPPS_PIXEL_STATION_ID || pixId.rp() != CTPPS_PIXEL_FAR_RP_ID)
1356           continue;
1357         if (ds.size() > 1)
1358           continue;
1359         for (const auto& lt : ds) {
1360           if (lt.isValid() && pixId.arm() == detId.arm() && lt.x0() - shift_pix < 24)
1361             channelPlots_[detId].pixelTomography_far->Fill(lt.x0() - shift_pix + 25 * rechit.ootIndex(), lt.y0());
1362         }
1363       }
1364     }
1365   }
1366 }
1367 
1368 //----------------------------------------------------------------------------------------------------
1369 
1370 void CTPPSDiamondDQMSource::globalEndLuminosityBlock(const edm::LuminosityBlock& iLumi, const edm::EventSetup&) {
1371   auto lumiCache = luminosityBlockCache(iLumi.index());
1372   if (!perLSsaving_) {
1373     if (plotOnline_)
1374       for (auto& plot : potPlots_)
1375         *(plot.second.hitDistribution2d_lumisection->getTH2F()) = *(lumiCache->hitDistribution2dMap[plot.first]);
1376 
1377     for (auto& plot : channelPlots_) {
1378       auto hitsCounterPerLumisection = lumiCache->hitsCounterMap[plot.first];
1379       if (hitsCounterPerLumisection != 0) {
1380         plot.second.hit_rate->Fill((double)hitsCounterPerLumisection / SEC_PER_LUMI_SECTION);
1381       }
1382 
1383       double HundredOverHitCounter = .0;
1384       if (plot.second.HitCounter != 0)
1385         HundredOverHitCounter = 100. / plot.second.HitCounter;
1386       plot.second.HPTDCErrorFlags->setBinContent(16, HundredOverHitCounter * plot.second.MHCounter);
1387       plot.second.leadingWithoutTrailing->setBinContent(1, HundredOverHitCounter * plot.second.LeadingOnlyCounter);
1388       plot.second.leadingWithoutTrailing->setBinContent(2, HundredOverHitCounter * plot.second.TrailingOnlyCounter);
1389       plot.second.leadingWithoutTrailing->setBinContent(3, HundredOverHitCounter * plot.second.CompleteCounter);
1390     }
1391 
1392     for (auto& plot : potPlots_) {
1393       double HundredOverHitCounterPot = 0.;
1394       if (plot.second.HitCounter != 0)
1395         HundredOverHitCounterPot = 100. / plot.second.HitCounter;
1396       plot.second.leadingWithoutTrailingCumulativePot->setBinContent(
1397           1, HundredOverHitCounterPot * plot.second.LeadingOnlyCounter);
1398       plot.second.leadingWithoutTrailingCumulativePot->setBinContent(
1399           2, HundredOverHitCounterPot * plot.second.TrailingOnlyCounter);
1400       plot.second.leadingWithoutTrailingCumulativePot->setBinContent(
1401           3, HundredOverHitCounterPot * plot.second.CompleteCounter);
1402 
1403       plot.second.MHComprensive->Reset();
1404       const CTPPSDiamondDetId rpId(plot.first);
1405       for (auto& chPlot : channelPlots_) {
1406         const CTPPSDiamondDetId chId(chPlot.first);
1407         if (chId.arm() == rpId.arm() && chId.rp() == rpId.rp()) {
1408           plot.second.MHComprensive->Fill(
1409               chId.plane(), chId.channel(), chPlot.second.HPTDCErrorFlags->getBinContent(16));
1410         }
1411       }
1412     }
1413     // Efficiencies of single channels
1414     if (plotOnline_)
1415       for (auto& plot : potPlots_) {
1416         plot.second.EfficiencyOfChannelsInPot->Reset();
1417         for (auto& element : plot.second.effTriplecountingChMap) {
1418           if (plot.second.effDoublecountingChMap[element.first] > 0) {
1419             int plane = element.first / 100;
1420             int channel = element.first % 100;
1421             double counted = element.second;
1422             double total = plot.second.effDoublecountingChMap[element.first];
1423             double efficiency = counted / total;
1424             //         double error = std::sqrt( efficiency * ( 1 - efficiency ) / total );
1425 
1426             plot.second.EfficiencyOfChannelsInPot->Fill(plane, channel, 100 * efficiency);
1427           }
1428         }
1429       }
1430 
1431     // Efficeincy wrt pixels  //TODO
1432     for (auto& plot : planePlots_) {
1433       TH2F* hitHistoTmp = plot.second.EfficiencyWRTPixelsInPlane->getTH2F();
1434 
1435       const CTPPSDiamondDetId detId(plot.first), detId_pot(detId.rpId());
1436       hitHistoTmp->Divide(&(plot.second.pixelTracksMapWithDiamonds), &(potPlots_[detId_pot].pixelTracksMap));
1437     }
1438   }  //perLSsaving
1439 }
1440 
1441 //----------------------------------------------------------------------------------------------------
1442 
1443 void CTPPSDiamondDQMSource::checkEventNumber(const CTPPSDiamondDetId& detId,
1444                                              const TotemFEDInfo& optorx,
1445                                              const TotemVFATStatus& status,
1446                                              CTPPSDiamondDQMSource::PotPlots& plots,
1447                                              int& EC_difference) const {
1448   const CTPPSDiamondDetId detId_pot(detId.rpId());
1449   if (plotOnline_)
1450     plots.ECCheck->Fill((int)((optorx.lv1() & 0xFF) - ((unsigned int)status.ec() & 0xFF)) & 0xFF);
1451   if ((static_cast<int>((optorx.lv1() & 0xFF) - status.ec()) != EC_difference) &&
1452       (static_cast<uint8_t>((optorx.lv1() & 0xFF) - status.ec()) < 128))
1453     EC_difference = static_cast<int>(optorx.lv1() & 0xFF) - (static_cast<unsigned int>(status.ec()) & 0xFF);
1454   if (EC_difference != 1 && EC_difference != -500 && std::abs(EC_difference) < 127) {
1455     if (detId.channel() == HPTDC_0_CHANNEL || detId.channel() == HPTDC_1_CHANNEL)
1456       plots.HPTDCErrorFlags_2D->Fill(16, 2 * detId.plane() + (detId.channel() - HPTDC_0_CHANNEL));
1457     if (verbosity_)
1458       edm::LogProblem("CTPPSDiamondDQMSource")
1459           << "FED " << optorx.fedId() << ": ECError at EV: 0x" << std::hex << optorx.lv1() << "\t\tVFAT EC: 0x"
1460           << static_cast<unsigned int>(status.ec()) << "\twith ID: " << std::dec << detId
1461           << "\tdiff: " << EC_difference;
1462   }
1463 }
1464 
1465 //----------------------------------------------------------------------------------------------------
1466 
1467 void CTPPSDiamondDQMSource::dqmEndRun(edm::Run const&, edm::EventSetup const&) {
1468   if (plotOffline_ && !perLSsaving_)
1469     for (const auto& rpPlots : potPlots_) {
1470       auto plots = rpPlots.second;
1471       // *(plots.TOTVsLSProfile->getTProfile())=*plots.TOTVsLS->getTH2F()->ProfileX();
1472       // *(plots.trackTimeVsLSProfile->getTProfile())=*plots.trackTimeVsLS->getTH2F()->ProfileX();
1473       *(plots.trackTimeVsBXProfile->getTProfile()) = *plots.trackTimeVsBX->getTH2F()->ProfileX();
1474       // *(plots.trackTimeVsXAngleProfile->getTProfile()) = *plots.trackTimeVsXAngle->getTH2F()->ProfileX();
1475     }
1476 }
1477 
1478 //----------------------------------------------------------------------------------------------------
1479 DEFINE_FWK_MODULE(CTPPSDiamondDQMSource);