Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-10-26 02:21:17

0001 /******************************************
0002  *
0003  * This is a part of CTPPSDQM software.
0004  * Authors:
0005  *   F.Ferro INFN Genova
0006  *   Vladimir Popov (vladimir.popov@cern.ch)
0007  *
0008  *******************************************/
0009 
0010 #include "FWCore/Framework/interface/ESHandle.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/EventSetup.h"
0013 #include "FWCore/Framework/interface/MakerMacros.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/Utilities/interface/InputTag.h"
0016 
0017 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0018 #include "DQMServices/Core/interface/DQMStore.h"
0019 
0020 #include "DataFormats/Common/interface/DetSetVector.h"
0021 
0022 #include "CondFormats/PPSObjects/interface/CTPPSPixelIndices.h"
0023 #include "DataFormats/CTPPSDetId/interface/CTPPSDetId.h"
0024 #include "DataFormats/CTPPSDigi/interface/CTPPSPixelDigi.h"
0025 #include "DataFormats/CTPPSDigi/interface/CTPPSPixelDataError.h"
0026 #include "DataFormats/CTPPSReco/interface/CTPPSPixelCluster.h"
0027 #include "DataFormats/CTPPSReco/interface/CTPPSPixelLocalTrack.h"
0028 
0029 #include <string>
0030 
0031 //-----------------------------------------------------------------------------
0032 
0033 class CTPPSPixelDQMSource : public DQMEDAnalyzer {
0034 public:
0035   CTPPSPixelDQMSource(const edm::ParameterSet &ps);
0036   ~CTPPSPixelDQMSource() override;
0037 
0038 protected:
0039   void dqmBeginRun(edm::Run const &, edm::EventSetup const &) override;
0040   void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override;
0041   void analyze(edm::Event const &e, edm::EventSetup const &eSetup) override;
0042 
0043 private:
0044   unsigned int verbosity;
0045   long int nEvents = 0;
0046   edm::EDGetTokenT<edm::DetSetVector<CTPPSPixelDigi>> tokenDigi;
0047   edm::EDGetTokenT<edm::DetSetVector<CTPPSPixelDataError>> tokenError;
0048   edm::EDGetTokenT<edm::DetSetVector<CTPPSPixelCluster>> tokenCluster;
0049   edm::EDGetTokenT<edm::DetSetVector<CTPPSPixelLocalTrack>> tokenTrack;
0050 
0051   static constexpr int NArms = 2;
0052   static constexpr int NStationMAX = 3;  // in an arm
0053   static constexpr int NRPotsMAX = 6;    // per station
0054   static constexpr int NplaneMAX = 6;    // per RPot
0055   static constexpr int NROCsMAX = 6;     // per plane
0056   static constexpr int RPn_first = 3, RPn_last = 4;
0057   static constexpr int ADCMax = 256;
0058   static constexpr int StationIDMAX = 4;  // possible range of ID
0059   static constexpr int RPotsIDMAX = 8;    // possible range of ID
0060   static constexpr int NLocalTracksMAX = 20;
0061   static constexpr int hitMultMAX = 50;   // tuned
0062   static constexpr int ClusMultMAX = 10;  // tuned
0063   static constexpr int ClusterSizeMax = 9;
0064   static constexpr int errCodeSize = 15;
0065 
0066   static constexpr int mapXbins = 200;
0067   static constexpr int mapYbins = 240;
0068   static constexpr float mapYmin = -16.;
0069   static constexpr float mapYmax = 8.;
0070   const float mapXmin = 0. * TMath::Cos(18.4 / 180. * TMath::Pi());
0071   const float mapXmax = 30. * TMath::Cos(18.4 / 180. * TMath::Pi());
0072 
0073   CTPPSPixelIndices thePixIndices;
0074 
0075   int TrackFitDimension = 4;
0076 
0077   static constexpr int NRPotBinsInStation = RPn_last - RPn_first;
0078   static constexpr int NPlaneBins = NplaneMAX * NRPotBinsInStation;
0079 
0080   MonitorElement *hBX, *hBXshort, *h2AllPlanesActive, *hpixLTrack;
0081   MonitorElement *hpRPactive;
0082 
0083   MonitorElement *h2HitsMultipl[NArms][NStationMAX];
0084   MonitorElement *h2CluSize[NArms][NStationMAX];
0085 
0086   static constexpr int RPotsTotalNumber = NArms * NStationMAX * NRPotsMAX;
0087 
0088   int RPindexValid[RPotsTotalNumber];
0089   MonitorElement *h2trackXY0[RPotsTotalNumber];
0090   MonitorElement *h2ErrorCodeRP[RPotsTotalNumber];
0091   MonitorElement *h2ErrorCode;
0092 
0093   MonitorElement *htrackMult[RPotsTotalNumber];
0094   MonitorElement *htrackHits[RPotsTotalNumber];
0095   MonitorElement *hRPotActivPlanes[RPotsTotalNumber];
0096   MonitorElement *hRPotActivBX[RPotsTotalNumber];
0097   MonitorElement *hRPotActivBXroc[RPotsTotalNumber];
0098   MonitorElement *hRPotActivBXroc_3[RPotsTotalNumber];
0099   MonitorElement *hRPotActivBXroc_2[RPotsTotalNumber];
0100   MonitorElement *h2HitsMultROC[RPotsTotalNumber];
0101   MonitorElement *hp2HitsMultROC_LS[RPotsTotalNumber];
0102   MonitorElement *hHitsMult[RPotsTotalNumber][NplaneMAX];
0103   MonitorElement *h2xyHits[RPotsTotalNumber][NplaneMAX];
0104   MonitorElement *hp2xyADC[RPotsTotalNumber][NplaneMAX];
0105   MonitorElement *h2Efficiency[RPotsTotalNumber][NplaneMAX];
0106   MonitorElement *h2xyROCHits[RPotsTotalNumber * NplaneMAX][NROCsMAX];
0107   MonitorElement *hROCadc[RPotsTotalNumber * NplaneMAX][NROCsMAX];
0108   MonitorElement *hRPotActivBXall[RPotsTotalNumber];
0109   int HitsMultROC[RPotsTotalNumber * NplaneMAX][NROCsMAX];
0110   int HitsMultPlane[RPotsTotalNumber][NplaneMAX];
0111 
0112   // Flags for disabling set of plots
0113   bool offlinePlots = true;
0114   bool onlinePlots = true;
0115 
0116   // Flags for disabling plots of a plane
0117   bool isPlanePlotsTurnedOff[NArms][NStationMAX][NRPotsMAX][NplaneMAX] = {};
0118 
0119   unsigned int rpStatusWord = 0x8008;      // 220_fr_hr(stn2rp3)+ 210_fr_hr
0120   int RPstatus[StationIDMAX][RPotsIDMAX];  // symmetric in both arms
0121   int StationStatus[StationIDMAX];         // symmetric in both arms
0122   const int IndexNotValid = 0;
0123 
0124   int getRPindex(int arm, int station, int rp) {
0125     if (arm < 0 || station < 0 || rp < 0)
0126       return (IndexNotValid);
0127     if (arm > 1 || station >= NStationMAX || rp >= NRPotsMAX)
0128       return (IndexNotValid);
0129     int rc = (arm * NStationMAX + station) * NRPotsMAX + rp;
0130     return (rc);
0131   }
0132 
0133   int getPlaneIndex(int arm, int station, int rp, int plane) {
0134     if (plane < 0 || plane >= NplaneMAX)
0135       return (IndexNotValid);
0136     int rc = getRPindex(arm, station, rp);
0137     if (rc == IndexNotValid)
0138       return (IndexNotValid);
0139     return (rc * NplaneMAX + plane);
0140   }
0141 
0142   int getRPInStationBin(int rp) { return (rp - RPn_first + 1); }
0143 
0144   static constexpr int NRPglobalBins = 4;  // 2 arms w. 2 stations w. 1 RP
0145 
0146   int getRPglobalBin(int arm, int stn) {
0147     static constexpr int stationBinOrder[NStationMAX] = {0, 4, 1};
0148     return (arm * 2 + stationBinOrder[stn] + 1);
0149   }
0150 
0151   int prIndex(int rp, int plane)  // plane index in station
0152 
0153   {
0154     return ((rp - RPn_first) * NplaneMAX + plane);
0155   }
0156   int getDet(int id) { return (id >> DetId::kDetOffset) & 0xF; }
0157   int getPixPlane(int id) { return ((id >> 16) & 0x7); }
0158   //  int getSubdet(int id) { return ((id>>kSubdetOffset)&0x7); }
0159 
0160   float x0_MIN, x0_MAX, y0_MIN, y0_MAX;
0161 };
0162 
0163 //----------------------------------------------------------------------------------
0164 
0165 using namespace std;
0166 using namespace edm;
0167 
0168 //-------------------------------------------------------------------------------
0169 
0170 CTPPSPixelDQMSource::CTPPSPixelDQMSource(const edm::ParameterSet &ps)
0171     : verbosity(ps.getUntrackedParameter<unsigned int>("verbosity", 0)),
0172       rpStatusWord(ps.getUntrackedParameter<unsigned int>("RPStatusWord", 0x8008)) {
0173   tokenDigi = consumes<DetSetVector<CTPPSPixelDigi>>(ps.getUntrackedParameter<edm::InputTag>("tagRPixDigi"));
0174   tokenError = consumes<DetSetVector<CTPPSPixelDataError>>(ps.getUntrackedParameter<edm::InputTag>("tagRPixError"));
0175   tokenCluster = consumes<DetSetVector<CTPPSPixelCluster>>(ps.getUntrackedParameter<edm::InputTag>("tagRPixCluster"));
0176   tokenTrack = consumes<DetSetVector<CTPPSPixelLocalTrack>>(ps.getUntrackedParameter<edm::InputTag>("tagRPixLTrack"));
0177   offlinePlots = ps.getUntrackedParameter<bool>("offlinePlots", true);
0178   onlinePlots = ps.getUntrackedParameter<bool>("onlinePlots", true);
0179 
0180   vector<string> disabledPlanePlotsVec =
0181       ps.getUntrackedParameter<vector<string>>("turnOffPlanePlots", vector<string>());
0182 
0183   // Parse the strings in disabledPlanePlotsVec and set the flags in
0184   // isPlanePlotsTurnedOff
0185   for (auto s : disabledPlanePlotsVec) {
0186     // Check that the format is <arm>_<station>_<RP>_<Plane>
0187     if (count(s.begin(), s.end(), '_') != 3)
0188       throw cms::Exception("RPixPlaneCombinatoryTracking") << "Invalid string in turnOffPlanePlots: " << s;
0189     else {
0190       vector<string> armStationRpPlane;
0191       size_t pos = 0;
0192       while ((pos = s.find('_')) != string::npos) {
0193         armStationRpPlane.push_back(s.substr(0, pos));
0194         s.erase(0, pos + 1);
0195       }
0196       armStationRpPlane.push_back(s);
0197 
0198       int arm = stoi(armStationRpPlane.at(0));
0199       int station = stoi(armStationRpPlane.at(1));
0200       int rp = stoi(armStationRpPlane.at(2));
0201       int plane = stoi(armStationRpPlane.at(3));
0202 
0203       if (arm < NArms && station < NStationMAX && rp < NRPotsMAX && plane < NplaneMAX) {
0204         if (verbosity)
0205           LogPrint("CTPPSPixelDQMSource")
0206               << "Shutting off plots for: Arm " << arm << " Station " << station << " Rp " << rp << " Plane " << plane;
0207         isPlanePlotsTurnedOff[arm][station][rp][plane] = true;
0208       } else {
0209         throw cms::Exception("RPixPlaneCombinatoryTracking") << "Invalid string in turnOffPlanePlots: " << s;
0210       }
0211     }
0212   }
0213 }
0214 
0215 //----------------------------------------------------------------------------------
0216 
0217 CTPPSPixelDQMSource::~CTPPSPixelDQMSource() {}
0218 
0219 //--------------------------------------------------------------------------
0220 
0221 void CTPPSPixelDQMSource::dqmBeginRun(edm::Run const &run, edm::EventSetup const &) {
0222   if (verbosity)
0223     LogPrint("CTPPSPixelDQMSource") << "RPstatusWord= " << rpStatusWord;
0224   nEvents = 0;
0225 
0226   CTPPSPixelLocalTrack thePixelLocalTrack;
0227   TrackFitDimension = thePixelLocalTrack.dimension;
0228 
0229   for (int stn = 0; stn < StationIDMAX; stn++) {
0230     StationStatus[stn] = 0;
0231     for (int rp = 0; rp < RPotsIDMAX; rp++)
0232       RPstatus[stn][rp] = 0;
0233   }
0234 
0235   unsigned int rpSts = rpStatusWord << 1;
0236   for (int stn = 0; stn < NStationMAX; stn++) {
0237     int stns = 0;
0238     for (int rp = 0; rp < NRPotsMAX; rp++) {
0239       rpSts = (rpSts >> 1);
0240       RPstatus[stn][rp] = rpSts & 1;
0241       if (RPstatus[stn][rp] > 0)
0242         stns = 1;
0243     }
0244     StationStatus[stn] = stns;
0245   }
0246 
0247   for (int ind = 0; ind < 2 * 3 * NRPotsMAX; ind++)
0248     RPindexValid[ind] = 0;
0249 
0250   x0_MIN = y0_MIN = 1.0e06;
0251   x0_MAX = y0_MAX = -1.0e06;
0252 }
0253 
0254 //-------------------------------------------------------------------------------------
0255 
0256 void CTPPSPixelDQMSource::bookHistograms(DQMStore::IBooker &ibooker, edm::Run const &, edm::EventSetup const &) {
0257   ibooker.cd();
0258   ibooker.setCurrentFolder("CTPPS/TrackingPixel");
0259   char s[50];
0260   string armTitleShort, stnTitleShort;
0261 
0262   TAxis *yah1st = nullptr;
0263   TAxis *xaRPact = nullptr;
0264   TAxis *xah1trk = nullptr;
0265   if (onlinePlots) {
0266     hBX = ibooker.book1D("events per BX", "ctpps_pixel;Event.BX", 4002, -1.5, 4000. + 0.5);
0267     hBXshort = ibooker.book1D("events per BX(short)", "ctpps_pixel;Event.BX", 102, -1.5, 100. + 0.5);
0268 
0269     string str1st = "Pixel planes activity";
0270     h2AllPlanesActive = ibooker.book2DD(
0271         str1st, str1st + "(digi task);Plane #", NplaneMAX, 0, NplaneMAX, NRPglobalBins, 0.5, NRPglobalBins + 0.5);
0272     TH2D *h1st = h2AllPlanesActive->getTH2D();
0273     h1st->SetOption("colz");
0274     yah1st = h1st->GetYaxis();
0275 
0276     string str2 = "Pixel RP active";
0277     hpRPactive = ibooker.bookProfile(
0278         str2, str2 + " per event(digi task)", NRPglobalBins, 0.5, NRPglobalBins + 0.5, -0.1, 1.1, "");
0279     xaRPact = hpRPactive->getTProfile()->GetXaxis();
0280     hpRPactive->getTProfile()->SetOption("hist");
0281     hpRPactive->getTProfile()->SetMinimum(0.);
0282     hpRPactive->getTProfile()->SetMaximum(1.1);
0283 
0284     str2 = "Pixel Local Tracks";
0285     hpixLTrack = ibooker.bookProfile(
0286         str2, str2 + " per event", NRPglobalBins, 0.5, NRPglobalBins + 0.5, -0.1, NLocalTracksMAX, "");
0287 
0288     xah1trk = hpixLTrack->getTProfile()->GetXaxis();
0289     hpixLTrack->getTProfile()->GetYaxis()->SetTitle("average number of tracks per event");
0290     hpixLTrack->getTProfile()->SetOption("hist");
0291   }
0292   const float minErrCode = 25.;
0293   const string errCode[errCodeSize] = {
0294       "Invalid ROC                        ",  // error  25
0295       "Gap word",                             // error 26
0296       "Dummy word",                           // error 27
0297       "FIFO nearly full",                     // error 28
0298       "Channel Timeout",                      // error 29
0299       "TBM Trailer",                          // error 30
0300       "Event number mismatch",                // error 31
0301       "Invalid/no FED header",                // error 32
0302       "Invalid/no FED trailer",               // error 33
0303       "Size mismatch",                        // error 34
0304       "Conversion: inv. channel",             // error 35
0305       "Conversion: inv. ROC number",          // error 36
0306       "Conversion: inv. pixel address",       // error 37
0307       "-",
0308       "CRC"  //error 39
0309   };
0310   h2ErrorCode = ibooker.book2D("Errors in Unidentified Det",
0311                                "Errors in Unidentified Det;;fed",
0312                                errCodeSize,
0313                                minErrCode - 0.5,
0314                                minErrCode + float(errCodeSize) - 0.5,
0315                                2,
0316                                1461.5,
0317                                1463.5);
0318   for (unsigned int iBin = 1; iBin <= errCodeSize; iBin++)
0319     h2ErrorCode->setBinLabel(iBin, errCode[iBin - 1]);
0320   h2ErrorCode->setBinLabel(1, "1462", 2);
0321   h2ErrorCode->setBinLabel(2, "1463", 2);
0322   h2ErrorCode->getTH2F()->SetOption("colz");
0323 
0324   for (int arm = 0; arm < 2; arm++) {
0325     CTPPSDetId ID(CTPPSDetId::sdTrackingPixel, arm, 0);
0326     string sd, armTitle;
0327     ID.armName(sd, CTPPSDetId::nPath);
0328     ID.armName(armTitle, CTPPSDetId::nFull);
0329     ID.armName(armTitleShort, CTPPSDetId::nShort);
0330 
0331     ibooker.setCurrentFolder(sd);
0332 
0333     for (int stn = 0; stn < NStationMAX; stn++) {
0334       if (StationStatus[stn] == 0)
0335         continue;
0336       ID.setStation(stn);
0337       string stnd, stnTitle;
0338 
0339       CTPPSDetId(ID.stationId()).stationName(stnd, CTPPSDetId::nPath);
0340       CTPPSDetId(ID.stationId()).stationName(stnTitle, CTPPSDetId::nFull);
0341       CTPPSDetId(ID.stationId()).stationName(stnTitleShort, CTPPSDetId::nShort);
0342 
0343       ibooker.setCurrentFolder(stnd);
0344       //--------- RPots ---
0345       int pixBinW = 4;
0346       for (int rp = RPn_first; rp < RPn_last; rp++) {  // only installed pixel pots
0347         ID.setRP(rp);
0348         string rpd, rpTitle;
0349         CTPPSDetId(ID.rpId()).rpName(rpTitle, CTPPSDetId::nShort);
0350         string rpBinName = armTitleShort + "_" + stnTitleShort + "_" + rpTitle;
0351         if (onlinePlots) {
0352           yah1st->SetBinLabel(getRPglobalBin(arm, stn), rpBinName.c_str());
0353           xah1trk->SetBinLabel(getRPglobalBin(arm, stn), rpBinName.c_str());
0354           xaRPact->SetBinLabel(getRPglobalBin(arm, stn), rpBinName.c_str());
0355         }
0356         if (RPstatus[stn][rp] == 0)
0357           continue;
0358         int indexP = getRPindex(arm, stn, rp);
0359         RPindexValid[indexP] = 1;
0360 
0361         CTPPSDetId(ID.rpId()).rpName(rpTitle, CTPPSDetId::nFull);
0362         CTPPSDetId(ID.rpId()).rpName(rpd, CTPPSDetId::nPath);
0363 
0364         ibooker.setCurrentFolder(rpd);
0365 
0366         const float x0Maximum = 70.;
0367         const float y0Maximum = 15.;
0368         string st = "track intercept point";
0369         string st2 = ": " + stnTitle;
0370         h2trackXY0[indexP] = ibooker.book2D(
0371             st, st + st2 + ";x0;y0", int(x0Maximum) * 2, 0., x0Maximum, int(y0Maximum) * 4, -y0Maximum, y0Maximum);
0372         h2trackXY0[indexP]->getTH2F()->SetOption("colz");
0373         st = "Error Code";
0374         h2ErrorCodeRP[indexP] = ibooker.book2D(st,
0375                                                st + st2 + ";;plane",
0376                                                errCodeSize,
0377                                                minErrCode - 0.5,
0378                                                minErrCode + float(errCodeSize) - 0.5,
0379                                                6,
0380                                                -0.5,
0381                                                5.5);
0382         for (unsigned int iBin = 1; iBin <= errCodeSize; iBin++)
0383           h2ErrorCodeRP[indexP]->setBinLabel(iBin, errCode[iBin - 1]);
0384         h2ErrorCodeRP[indexP]->getTH2F()->SetOption("colz");
0385 
0386         st = "number of tracks per event";
0387         htrackMult[indexP] = ibooker.bookProfile(st,
0388                                                  rpTitle + ";number of tracks",
0389                                                  NLocalTracksMAX + 1,
0390                                                  -0.5,
0391                                                  NLocalTracksMAX + 0.5,
0392                                                  -0.5,
0393                                                  NLocalTracksMAX + 0.5,
0394                                                  "");
0395         htrackMult[indexP]->getTProfile()->SetOption("hist");
0396 
0397         hRPotActivPlanes[indexP] = ibooker.bookProfile("number of fired planes per event",
0398                                                        rpTitle + ";nPlanes;Probability",
0399                                                        NplaneMAX + 1,
0400                                                        -0.5,
0401                                                        NplaneMAX + 0.5,
0402                                                        -0.5,
0403                                                        NplaneMAX + 0.5,
0404                                                        "");
0405         hRPotActivPlanes[indexP]->getTProfile()->SetOption("hist");
0406 
0407         hp2HitsMultROC_LS[indexP] = ibooker.bookProfile2D("ROCs hits multiplicity per event vs LS",
0408                                                           rpTitle + ";LumiSection;Plane#___ROC#",
0409                                                           1000,
0410                                                           0.,
0411                                                           1000.,
0412                                                           NplaneMAX * NROCsMAX,
0413                                                           0.,
0414                                                           double(NplaneMAX * NROCsMAX),
0415                                                           0.,
0416                                                           rpixValues::ROCSizeInX *rpixValues::ROCSizeInY,
0417                                                           "");
0418         hp2HitsMultROC_LS[indexP]->getTProfile2D()->SetOption("colz");
0419         hp2HitsMultROC_LS[indexP]->getTProfile2D()->SetMinimum(1.0e-10);
0420         hp2HitsMultROC_LS[indexP]->getTProfile2D()->SetCanExtend(TProfile2D::kXaxis);
0421         TAxis *yahp2 = hp2HitsMultROC_LS[indexP]->getTProfile2D()->GetYaxis();
0422         for (int p = 0; p < NplaneMAX; p++) {
0423           sprintf(s, "plane%d_0", p);
0424           yahp2->SetBinLabel(p * NplaneMAX + 1, s);
0425           for (int r = 1; r < NROCsMAX; r++) {
0426             sprintf(s, "   %d_%d", p, r);
0427             yahp2->SetBinLabel(p * NplaneMAX + r + 1, s);
0428           }
0429         }
0430 
0431         if (onlinePlots) {
0432           string st3 = ";PlaneIndex(=pixelPot*PlaneMAX + plane)";
0433 
0434           st = "hit multiplicity in planes";
0435           h2HitsMultipl[arm][stn] = ibooker.book2DD(
0436               st, st + st2 + st3 + ";multiplicity", NPlaneBins, 0, NPlaneBins, hitMultMAX, 0, hitMultMAX);
0437           h2HitsMultipl[arm][stn]->getTH2D()->SetOption("colz");
0438 
0439           st = "cluster size in planes";
0440           h2CluSize[arm][stn] = ibooker.book2D(st,
0441                                                st + st2 + st3 + ";Cluster size",
0442                                                NPlaneBins,
0443                                                0,
0444                                                NPlaneBins,
0445                                                ClusterSizeMax + 1,
0446                                                0,
0447                                                ClusterSizeMax + 1);
0448           h2CluSize[arm][stn]->getTH2F()->SetOption("colz");
0449 
0450           st = "number of hits per track";
0451           htrackHits[indexP] = ibooker.bookProfile(st, rpTitle + ";number of hits", 5, 1.5, 6.5, -0.1, 1.1, "");
0452           htrackHits[indexP]->getTProfile()->SetOption("hist");
0453 
0454           h2HitsMultROC[indexP] = ibooker.bookProfile2D("ROCs hits multiplicity per event",
0455                                                         rpTitle + ";plane # ;ROC #",
0456                                                         NplaneMAX,
0457                                                         -0.5,
0458                                                         NplaneMAX - 0.5,
0459                                                         NROCsMAX,
0460                                                         -0.5,
0461                                                         NROCsMAX - 0.5,
0462                                                         0.,
0463                                                         rpixValues::ROCSizeInX * rpixValues::ROCSizeInY,
0464                                                         "");
0465           h2HitsMultROC[indexP]->getTProfile2D()->SetOption("colztext");
0466           h2HitsMultROC[indexP]->getTProfile2D()->SetMinimum(1.e-10);
0467 
0468           ibooker.setCurrentFolder(rpd + "/latency");
0469           hRPotActivBX[indexP] =
0470               ibooker.book1D("5 fired planes per BX", rpTitle + ";Event.BX", 4002, -1.5, 4000. + 0.5);
0471 
0472           hRPotActivBXroc[indexP] =
0473               ibooker.book1D("4 fired ROCs per BX", rpTitle + ";Event.BX", 4002, -1.5, 4000. + 0.5);
0474           hRPotActivBXroc_3[indexP] =
0475               ibooker.book1D("3 fired ROCs per BX", rpTitle + ";Event.BX", 4002, -1.5, 4000. + 0.5);
0476           hRPotActivBXroc_2[indexP] =
0477               ibooker.book1D("2 fired ROCs per BX", rpTitle + ";Event.BX", 4002, -1.5, 4000. + 0.5);
0478 
0479           hRPotActivBXall[indexP] = ibooker.book1D("hits per BX", rpTitle + ";Event.BX", 4002, -1.5, 4000. + 0.5);
0480         }
0481         int nbins = rpixValues::defaultDetSizeInX / pixBinW;
0482 
0483         for (int p = 0; p < NplaneMAX; p++) {
0484           if (isPlanePlotsTurnedOff[arm][stn][rp][p])
0485             continue;
0486           sprintf(s, "plane_%d", p);
0487           string pd = rpd + "/" + string(s);
0488           ibooker.setCurrentFolder(pd);
0489           string st1 = ": " + rpTitle + "_" + string(s);
0490 
0491           st = "adc average value";
0492           hp2xyADC[indexP][p] = ibooker.bookProfile2D(st,
0493                                                       st1 + ";pix col;pix row",
0494                                                       nbins,
0495                                                       0,
0496                                                       rpixValues::defaultDetSizeInX,
0497                                                       nbins,
0498                                                       0,
0499                                                       rpixValues::defaultDetSizeInX,
0500                                                       0.,
0501                                                       512.,
0502                                                       "");
0503           hp2xyADC[indexP][p]->getTProfile2D()->SetOption("colz");
0504 
0505           if (onlinePlots) {
0506             st = "hits position";
0507             h2xyHits[indexP][p] = ibooker.book2DD(st,
0508                                                   st1 + ";pix col;pix row",
0509                                                   rpixValues::defaultDetSizeInX,
0510                                                   0,
0511                                                   rpixValues::defaultDetSizeInX,
0512                                                   rpixValues::defaultDetSizeInX,
0513                                                   0,
0514                                                   rpixValues::defaultDetSizeInX);
0515             h2xyHits[indexP][p]->getTH2D()->SetOption("colz");
0516 
0517             st = "hits multiplicity";
0518             hHitsMult[indexP][p] =
0519                 ibooker.book1DD(st, st1 + ";number of hits;N / 1 hit", hitMultMAX + 1, -0.5, hitMultMAX + 0.5);
0520           }
0521 
0522           if (offlinePlots) {
0523             st = "plane efficiency";
0524             h2Efficiency[indexP][p] = ibooker.bookProfile2D(
0525                 st, st1 + ";x0;y0", mapXbins, mapXmin, mapXmax, mapYbins, mapYmin, mapYmax, 0, 1, "");
0526             h2Efficiency[indexP][p]->getTProfile2D()->SetOption("colz");
0527           }
0528         }  // end of for(int p=0; p<NplaneMAX;..
0529 
0530       }  // end for(int rp=0; rp<NRPotsMAX;...
0531     }    // end of for(int stn=0; stn<
0532   }      // end of for(int arm=0; arm<2;...
0533 
0534   return;
0535 }
0536 
0537 //-------------------------------------------------------------------------------
0538 
0539 void CTPPSPixelDQMSource::analyze(edm::Event const &event, edm::EventSetup const &eventSetup) {
0540   ++nEvents;
0541   int lumiId = event.getLuminosityBlock().id().luminosityBlock();
0542   if (lumiId < 0)
0543     lumiId = 0;
0544 
0545   int RPactivity[RPotsTotalNumber], RPdigiSize[RPotsTotalNumber];
0546   int pixRPTracks[RPotsTotalNumber];
0547 
0548   for (int rp = 0; rp < RPotsTotalNumber; rp++) {
0549     RPactivity[rp] = RPdigiSize[rp] = pixRPTracks[rp] = 0;
0550   }
0551 
0552   for (int ind = 0; ind < RPotsTotalNumber; ind++) {
0553     for (int p = 0; p < NplaneMAX; p++) {
0554       HitsMultPlane[ind][p] = 0;
0555     }
0556   }
0557   for (int ind = 0; ind < RPotsTotalNumber * NplaneMAX; ind++) {
0558     for (int roc = 0; roc < NROCsMAX; roc++) {
0559       HitsMultROC[ind][roc] = 0;
0560     }
0561   }
0562   Handle<DetSetVector<CTPPSPixelDigi>> pixDigi;
0563   event.getByToken(tokenDigi, pixDigi);
0564 
0565   Handle<DetSetVector<CTPPSPixelDataError>> pixError;
0566   event.getByToken(tokenError, pixError);
0567 
0568   Handle<DetSetVector<CTPPSPixelCluster>> pixClus;
0569   event.getByToken(tokenCluster, pixClus);
0570 
0571   Handle<DetSetVector<CTPPSPixelLocalTrack>> pixTrack;
0572   event.getByToken(tokenTrack, pixTrack);
0573 
0574   if (onlinePlots) {
0575     hBX->Fill(event.bunchCrossing());
0576     hBXshort->Fill(event.bunchCrossing());
0577   }
0578 
0579   if (pixTrack.isValid()) {
0580     for (const auto &ds_tr : *pixTrack) {
0581       int idet = getDet(ds_tr.id);
0582       if (idet != DetId::VeryForward) {
0583         if (verbosity > 1)
0584           LogPrint("CTPPSPixelDQMSource") << "not CTPPS: ds_tr.id" << ds_tr.id;
0585         continue;
0586       }
0587       CTPPSDetId theId(ds_tr.id);
0588       int arm = theId.arm() & 0x1;
0589       int station = theId.station() & 0x3;
0590       int rpot = theId.rp() & 0x7;
0591       int rpInd = getRPindex(arm, station, rpot);
0592 
0593       for (DetSet<CTPPSPixelLocalTrack>::const_iterator dit = ds_tr.begin(); dit != ds_tr.end(); ++dit) {
0594         ++pixRPTracks[rpInd];
0595         int nh_tr = (dit->ndf() + TrackFitDimension) / 2;
0596         if (onlinePlots) {
0597           for (int i = 0; i <= NplaneMAX; i++) {
0598             if (i == nh_tr)
0599               htrackHits[rpInd]->Fill(nh_tr, 1.);
0600             else
0601               htrackHits[rpInd]->Fill(i, 0.);
0602           }
0603         }
0604         float x0 = dit->x0();
0605         float y0 = dit->y0();
0606         h2trackXY0[rpInd]->Fill(x0, y0);
0607 
0608         if (x0_MAX < x0)
0609           x0_MAX = x0;
0610         if (y0_MAX < y0)
0611           y0_MAX = y0;
0612         if (x0_MIN > x0)
0613           x0_MIN = x0;
0614         if (y0_MIN > y0)
0615           y0_MIN = y0;
0616 
0617         if (offlinePlots) {
0618           edm::DetSetVector<CTPPSPixelFittedRecHit> fittedHits = dit->hits();
0619 
0620           std::map<int, int> numberOfPointPerPlaneEff;
0621           for (const auto &ds_frh : fittedHits) {
0622             int plane = getPixPlane(ds_frh.id);
0623             for (DetSet<CTPPSPixelFittedRecHit>::const_iterator frh_it = ds_frh.begin(); frh_it != ds_frh.end();
0624                  ++frh_it) {  // there should always be only one hit in each
0625                               // vector
0626               if (frh_it != ds_frh.begin())
0627                 if (verbosity > 1)
0628                   LogPrint("CTPPSPixelDQMSource") << "More than one FittedRecHit found in plane " << plane;
0629               if (frh_it->isRealHit())
0630                 for (int p = 0; p < NplaneMAX; p++) {
0631                   if (p != plane)
0632                     numberOfPointPerPlaneEff[p]++;
0633                 }
0634             }
0635           }
0636 
0637           if (verbosity > 1)
0638             for (auto planeAndHitsOnOthers : numberOfPointPerPlaneEff) {
0639               LogPrint("CTPPSPixelDQMSource")
0640                   << "For plane " << planeAndHitsOnOthers.first << ", " << planeAndHitsOnOthers.second
0641                   << " hits on other planes were found" << endl;
0642             }
0643 
0644           for (const auto &ds_frh : fittedHits) {
0645             int plane = getPixPlane(ds_frh.id);
0646             if (isPlanePlotsTurnedOff[arm][station][rpot][plane])
0647               continue;
0648             for (DetSet<CTPPSPixelFittedRecHit>::const_iterator frh_it = ds_frh.begin(); frh_it != ds_frh.end();
0649                  ++frh_it) {
0650               float frhX0 = frh_it->globalCoordinates().x() + frh_it->xResidual();
0651               float frhY0 = frh_it->globalCoordinates().y() + frh_it->yResidual();
0652               if (numberOfPointPerPlaneEff[plane] >= 3) {
0653                 if (frh_it->isRealHit())
0654                   h2Efficiency[rpInd][plane]->Fill(frhX0, frhY0, 1);
0655                 else
0656                   h2Efficiency[rpInd][plane]->Fill(frhX0, frhY0, 0);
0657               }
0658             }
0659           }
0660         }
0661       }
0662     }
0663   }  // end  if(pixTrack.isValid())
0664 
0665   bool valid = false;
0666   valid |= pixDigi.isValid();
0667   //  valid |= Clus.isValid();
0668 
0669   if (!valid && verbosity)
0670     LogPrint("CTPPSPixelDQMSource") << "No valid data in Event " << nEvents;
0671 
0672   if (pixDigi.isValid()) {
0673     for (const auto &ds_digi : *pixDigi) {
0674       int idet = getDet(ds_digi.id);
0675       if (idet != DetId::VeryForward) {
0676         if (verbosity > 1)
0677           LogPrint("CTPPSPixelDQMSource") << "not CTPPS: ds_digi.id" << ds_digi.id;
0678         continue;
0679       }
0680       //   int subdet = getSubdet(ds_digi.id);
0681 
0682       int plane = getPixPlane(ds_digi.id);
0683 
0684       CTPPSDetId theId(ds_digi.id);
0685       int arm = theId.arm() & 0x1;
0686       int station = theId.station() & 0x3;
0687       int rpot = theId.rp() & 0x7;
0688       int rpInd = getRPindex(arm, station, rpot);
0689       RPactivity[rpInd] = 1;
0690       ++RPdigiSize[rpInd];
0691 
0692       if (StationStatus[station] && RPstatus[station][rpot]) {
0693         if (onlinePlots) {
0694           h2HitsMultipl[arm][station]->Fill(prIndex(rpot, plane), ds_digi.data.size());
0695           h2AllPlanesActive->Fill(plane, getRPglobalBin(arm, station));
0696         }
0697         int index = getRPindex(arm, station, rpot);
0698         HitsMultPlane[index][plane] += ds_digi.data.size();
0699         if (RPindexValid[index]) {
0700           int nh = ds_digi.data.size();
0701           if (nh > hitMultMAX)
0702             nh = hitMultMAX;
0703           if (!isPlanePlotsTurnedOff[arm][station][rpot][plane])
0704             if (onlinePlots)
0705               hHitsMult[index][plane]->Fill(nh);
0706         }
0707         int rocHistIndex = getPlaneIndex(arm, station, rpot, plane);
0708 
0709         for (DetSet<CTPPSPixelDigi>::const_iterator dit = ds_digi.begin(); dit != ds_digi.end(); ++dit) {
0710           int row = dit->row();
0711           int col = dit->column();
0712           int adc = dit->adc();
0713 
0714           if (RPindexValid[index]) {
0715             if (!isPlanePlotsTurnedOff[arm][station][rpot][plane]) {
0716               if (onlinePlots)
0717                 h2xyHits[index][plane]->Fill(col, row);
0718               hp2xyADC[index][plane]->Fill(col, row, adc);
0719             }
0720             int colROC, rowROC;
0721             int trocId;
0722             if (!thePixIndices.transformToROC(col, row, trocId, colROC, rowROC)) {
0723               if (trocId >= 0 && trocId < NROCsMAX) {
0724                 ++HitsMultROC[rocHistIndex][trocId];
0725               }
0726             }
0727           }  // end if(RPindexValid[index]) {
0728         }
0729       }  // end  if(StationStatus[station]) {
0730     }    // end for(const auto &ds_digi : *pixDigi)
0731   }      // if(pixDigi.isValid()) {
0732 
0733   if (pixError.isValid()) {
0734     for (const auto &ds_error : *pixError) {
0735       int idet = getDet(ds_error.id);
0736       if (idet != DetId::VeryForward) {
0737         if (idet == 15) {  //dummy det id: store in a plot with fed info
0738 
0739           for (DetSet<CTPPSPixelDataError>::const_iterator dit = ds_error.begin(); dit != ds_error.end(); ++dit) {
0740             h2ErrorCode->Fill(dit->errorType(), dit->fedId());
0741           }
0742           continue;
0743         }
0744         if (verbosity > 1)
0745           LogPrint("CTPPSPixelDQMSource") << "not CTPPS: ds_error.id" << ds_error.id;
0746         continue;
0747       }
0748 
0749       int plane = getPixPlane(ds_error.id);
0750       CTPPSDetId theId(ds_error.id);
0751       int arm = theId.arm() & 0x1;
0752       int station = theId.station() & 0x3;
0753       int rpot = theId.rp() & 0x7;
0754       int rpInd = getRPindex(arm, station, rpot);
0755       RPactivity[rpInd] = 1;
0756 
0757       if (StationStatus[station] && RPstatus[station][rpot]) {
0758         int index = getRPindex(arm, station, rpot);
0759         for (DetSet<CTPPSPixelDataError>::const_iterator dit = ds_error.begin(); dit != ds_error.end(); ++dit) {
0760           if (RPindexValid[index]) {
0761             if (!isPlanePlotsTurnedOff[arm][station][rpot][plane]) {
0762               h2ErrorCodeRP[index]->Fill(dit->errorType(), plane);
0763             }
0764           }  // end if(RPindexValid[index]) {
0765         }
0766       }  // end  if(StationStatus[station]) {
0767     }    // end for(const auto &ds_error : *pixDigi)
0768   }      // if(pixError.isValid())
0769 
0770   if (pixClus.isValid() && onlinePlots)
0771     for (const auto &ds : *pixClus) {
0772       int idet = getDet(ds.id);
0773       if (idet != DetId::VeryForward && verbosity > 1) {
0774         LogPrint("CTPPSPixelDQMSource") << "not CTPPS: cluster.id" << ds.id;
0775         continue;
0776       }
0777 
0778       CTPPSDetId theId(ds.id);
0779       int plane = getPixPlane(ds.id);
0780       int arm = theId.arm() & 0x1;
0781       int station = theId.station() & 0x3;
0782       int rpot = theId.rp() & 0x7;
0783 
0784       if ((StationStatus[station] == 0) || (RPstatus[station][rpot] == 0))
0785         continue;
0786 
0787       for (const auto &p : ds) {
0788         int clusize = p.size();
0789 
0790         if (clusize > ClusterSizeMax)
0791           clusize = ClusterSizeMax;
0792 
0793         h2CluSize[arm][station]->Fill(prIndex(rpot, plane), clusize);
0794       }
0795     }  // end if(pixClus.isValid()) for(const auto &ds : *pixClus)
0796 
0797   bool allRPactivity = false;
0798   for (int rp = 0; rp < RPotsTotalNumber; rp++)
0799     if (RPactivity[rp] > 0)
0800       allRPactivity = true;
0801   for (int arm = 0; arm < 2; arm++) {
0802     for (int stn = 0; stn < NStationMAX; stn++) {
0803       for (int rp = 0; rp < NRPotsMAX; rp++) {
0804         int index = getRPindex(arm, stn, rp);
0805         if (RPindexValid[index] == 0)
0806           continue;
0807 
0808         if (onlinePlots)
0809           hpRPactive->Fill(getRPglobalBin(arm, stn), RPactivity[index]);
0810         //        if(RPactivity[index]==0) continue;
0811         if (!allRPactivity)
0812           continue;
0813         if (onlinePlots)
0814           hpixLTrack->Fill(getRPglobalBin(arm, stn), pixRPTracks[index]);
0815         int ntr = pixRPTracks[index];
0816         if (ntr > NLocalTracksMAX)
0817           ntr = NLocalTracksMAX;
0818         for (int i = 0; i <= NLocalTracksMAX; i++) {
0819           if (i == ntr)
0820             htrackMult[index]->Fill(ntr, 1.);
0821           else
0822             htrackMult[index]->Fill(i, 0.);
0823         }
0824 
0825         int np = 0;
0826         for (int p = 0; p < NplaneMAX; p++)
0827           if (HitsMultPlane[index][p] > 0)
0828             np++;
0829         for (int p = 0; p <= NplaneMAX; p++) {
0830           if (p == np)
0831             hRPotActivPlanes[index]->Fill(p, 1.);
0832           else
0833             hRPotActivPlanes[index]->Fill(p, 0.);
0834         }
0835         if (onlinePlots) {
0836           if (np >= 5)
0837             hRPotActivBX[index]->Fill(event.bunchCrossing());
0838           hRPotActivBXall[index]->Fill(event.bunchCrossing(), float(RPdigiSize[index]));
0839         }
0840 
0841         int planesFiredAtROC[NROCsMAX];  // how many planes registered hits on ROC r
0842         for (int r = 0; r < NROCsMAX; r++)
0843           planesFiredAtROC[r] = 0;
0844         for (int p = 0; p < NplaneMAX; p++) {
0845           int indp = getPlaneIndex(arm, stn, rp, p);
0846           for (int r = 0; r < NROCsMAX; r++)
0847             if (HitsMultROC[indp][r] > 0)
0848               ++planesFiredAtROC[r];
0849           for (int r = 0; r < NROCsMAX; r++) {
0850             if (onlinePlots)
0851               h2HitsMultROC[index]->Fill(p, r, HitsMultROC[indp][r]);
0852             hp2HitsMultROC_LS[index]->Fill(lumiId, p * NROCsMAX + r, HitsMultROC[indp][r]);
0853           }
0854         }
0855         int max = 0;
0856         for (int r = 0; r < NROCsMAX; r++)
0857           if (max < planesFiredAtROC[r])
0858             max = planesFiredAtROC[r];
0859         if (max >= 4 && onlinePlots)  // fill only if there are at least 4 aligned ROCs firing
0860           hRPotActivBXroc[index]->Fill(event.bunchCrossing());
0861         if (max >= 3 && onlinePlots)  // fill only if there are at least 3 aligned ROCs firing
0862           hRPotActivBXroc_3[index]->Fill(event.bunchCrossing());
0863         if (max >= 2 && onlinePlots)  // fill only if there are at least 2 aligned ROCs firing
0864           hRPotActivBXroc_2[index]->Fill(event.bunchCrossing());
0865       }  // end for(int rp=0; rp<NRPotsMAX; rp++) {
0866     }
0867   }  // end for(int arm=0; arm<2; arm++) {
0868 
0869   if ((nEvents % 100))
0870     return;
0871   if (verbosity)
0872     LogPrint("CTPPSPixelDQMSource") << "analyze event " << nEvents;
0873 }
0874 
0875 //---------------------------------------------------------------------------
0876 DEFINE_FWK_MODULE(CTPPSPixelDQMSource);