Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-06-22 22:49:02

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