CalibTracker/src/APVValidationPlots.cc

0001 // system include files
0002 #include <iostream>
0003 #include <memory>
0004 #include <sstream>
0005 #include <vector>
0006
0007 // user include files
0008 #include "FWCore/Framework/interface/Frameworkfwd.h"
0009 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/Framework/interface/MakerMacros.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include "FWCore/ServiceRegistry/interface/Service.h"
0014 #include "DataFormats/SiStripCluster/interface/SiStripClusterCollection.h"
0015 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0016 #include "DataFormats/Common/interface/DetSetVector.h"
0017 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0018 #include "DQMServices/Core/interface/DQMStore.h"
0019
0020 // ROOT includes
0021 #include "TCanvas.h"
0022 #include "TFile.h"
0023 #include "TH1.h"
0024 #include "TH1D.h"
0025 #include "TH1F.h"
0026 #include "TH2F.h"
0027 #include "TProfile.h"
0028 #include "TStyle.h"
0029 #include "TTree.h"
0030
0031 //
0032 // class decleration
0033 //
0034
0035 class APVValidationPlots : public edm::one::EDAnalyzer<> {
0036 public:
0037   typedef dqm::legacy::DQMStore DQMStore;
0038   typedef dqm::legacy::MonitorElement MonitorElement;
0039
0040   explicit APVValidationPlots(const edm::ParameterSet&);
0041
0042 private:
0043   void analyze(const edm::Event&, const edm::EventSetup&) override;
0044   void endJob() override;
0045
0046   const std::string infilename;
0047   const std::string outfilename;
0048
0049   // ----------member data ---------------------------
0050 };
0051
0052 //
0053 // constructors and destructor
0054 //
0055 APVValidationPlots::APVValidationPlots(const edm::ParameterSet& iConfig)
0056     : infilename(iConfig.getUntrackedParameter<std::string>("inputFilename", "in.root")),
0057       outfilename(iConfig.getUntrackedParameter<std::string>("outputFilename", "out.root"))
0058
0059 {
0060   //now do what ever initialization is needed
0061 }
0062
0063 //
0064 // member functions
0065 //
0066
0067 // ------------ method called to for each event  ------------
0068 void APVValidationPlots::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {}
0069
0070 // ------------ method called once each job just after ending the event loop  ------------
0071 void APVValidationPlots::endJob() {
0072   std::ostringstream oss;
0073   oss << 1;  //runNumber
0074
0075   DQMStore* dqmStore = edm::Service<DQMStore>().operator->();
0076   dqmStore->setCurrentFolder("ChannelStatusPlots");
0077
0078   // Initialize histograms
0079   std::vector<std::string> subDetName;
0080   std::vector<unsigned int> nLayers;
0081   std::vector<std::string> layerName;
0082   subDetName.push_back("");
0083   subDetName.push_back("TIB");
0084   subDetName.push_back("TID");
0085   subDetName.push_back("TOB");
0086   subDetName.push_back("TEC");
0087   nLayers.push_back(0);
0088   nLayers.push_back(4);
0089   nLayers.push_back(3);
0090   nLayers.push_back(6);
0091   nLayers.push_back(9);
0092   layerName.push_back("");
0093   layerName.push_back("Layer");
0094   layerName.push_back("Disk");
0095   layerName.push_back("Layer");
0096   layerName.push_back("Disk");
0097
0098   std::string histoName;
0099   std::string histoTitle;
0100   // Histograms
0101   // indexes in these arrays are [SubDetId-2][LayerN]
0102   // histograms for [SubDetId-2][0] are global for the subdetector
0103   // histogram for [0][0] is global for the tracker
0104   TH2F* medianVsAbsoluteOccupancy[5][10];
0105   TH1F* medianOccupancy[5][10];
0106   TH1F* absoluteOccupancy[5][10];
0107   for (unsigned int i = 0; i < subDetName.size(); i++) {
0108     for (unsigned int j = 0; j <= nLayers[i]; j++) {
0109       histoName = "medianVsAbsoluteOccupancy" + subDetName[i];
0110       if (j != 0) {
0111         oss.str("");
0112         oss << j;
0113         histoName += layerName[i] + oss.str();
0114       }
0115       histoTitle = "Median APV occupancy vs. absolute APV occupancy";
0116       if (i != 0)
0117         histoTitle += " in " + subDetName[i];
0118       if (j != 0) {
0119         histoTitle += " " + layerName[i] + " " + oss.str();
0120       }
0121       MonitorElement* tmp = dqmStore->book2D(histoName.c_str(), histoTitle.c_str(), 1000, 0., 6., 1000, -1., 3.);
0122       medianVsAbsoluteOccupancy[i][j] = tmp->getTH2F();
0123       medianVsAbsoluteOccupancy[i][j]->Rebin2D(10, 10);
0124       medianVsAbsoluteOccupancy[i][j]->GetXaxis()->SetTitle("log_{10}(Abs. Occupancy)");
0125       medianVsAbsoluteOccupancy[i][j]->GetYaxis()->SetTitle("log_{10}(Median Occupancy)");
0126       //
0127       histoName = "medianOccupancy" + subDetName[i];
0128       if (j != 0) {
0129         oss.str("");
0130         oss << j;
0131         histoName += layerName[i] + oss.str();
0132       }
0133       histoTitle = "Median APV occupancy";
0134       if (i != 0)
0135         histoTitle += " in " + subDetName[i];
0136       if (j != 0) {
0137         histoTitle += " " + layerName[i] + " " + oss.str();
0138       }
0139       tmp = dqmStore->book1D(histoName.c_str(), histoTitle.c_str(), 1000, -1., 3.);
0140       medianOccupancy[i][j] = tmp->getTH1F();
0141       medianOccupancy[i][j]->GetXaxis()->SetTitle("log_{10}(Occupancy)");
0142       medianOccupancy[i][j]->GetYaxis()->SetTitle("APVs");
0143       //
0144       histoName = "absoluteOccupancy" + subDetName[i];
0145       if (j != 0) {
0146         oss.str("");
0147         oss << j;
0148         histoName += layerName[i] + oss.str();
0149       }
0150       histoTitle = "Absolute APV occupancy";
0151       if (i != 0)
0152         histoTitle += " in " + subDetName[i];
0153       if (j != 0) {
0154         histoTitle += " " + layerName[i] + " " + oss.str();
0155       }
0156       tmp = dqmStore->book1D(histoName.c_str(), histoTitle.c_str(), 1000, 0., 6.);
0157       absoluteOccupancy[i][j] = tmp->getTH1F();
0158       absoluteOccupancy[i][j]->GetXaxis()->SetTitle("log_{10}(Occupancy)");
0159       absoluteOccupancy[i][j]->GetYaxis()->SetTitle("APVs");
0160     }
0161   }
0162
0163   TFile* infile = new TFile(infilename.c_str(), "READ");
0164   TTree* intree = (TTree*)infile->Get("moduleOccupancy");
0165
0166   // Declaration of leaf types
0167   Int_t DetRawId;
0168   Int_t SubDetId;
0169   Int_t Layer_Ring;
0170   Int_t Disc;
0171   Int_t IsBack;
0172   Int_t IsExternalString;
0173   Int_t IsZMinusSide;
0174   Int_t RodStringPetal;
0175   Int_t IsStereo;
0176   Int_t ModulePosition;
0177   Int_t NumberOfStrips;
0178   Float_t APVGlobalPositionX;
0179   Float_t APVGlobalPositionY;
0180   Float_t APVGlobalPositionZ;
0181   Int_t APVNumber;
0182   Int_t APVAbsoluteOccupancy;
0183   Double_t APVMedianOccupancy;
0184   intree->SetBranchAddress("DetRawId", &DetRawId);
0185   intree->SetBranchAddress("SubDetId", &SubDetId);
0186   intree->SetBranchAddress("Layer_Ring", &Layer_Ring);
0187   intree->SetBranchAddress("Disc", &Disc);
0188   intree->SetBranchAddress("IsBack", &IsBack);
0189   intree->SetBranchAddress("IsExternalString", &IsExternalString);
0190   intree->SetBranchAddress("IsZMinusSide", &IsZMinusSide);
0191   intree->SetBranchAddress("RodStringPetal", &RodStringPetal);
0192   intree->SetBranchAddress("IsStereo", &IsStereo);
0193   intree->SetBranchAddress("ModuleNumber", &ModulePosition);
0194   intree->SetBranchAddress("NumberOfStrips", &NumberOfStrips);
0195   intree->SetBranchAddress("APVGlobalPositionX", &APVGlobalPositionX);
0196   intree->SetBranchAddress("APVGlobalPositionY", &APVGlobalPositionY);
0197   intree->SetBranchAddress("APVGlobalPositionZ", &APVGlobalPositionZ);
0198   intree->SetBranchAddress("APVNumber", &APVNumber);
0199   intree->SetBranchAddress("APVAbsoluteOccupancy", &APVAbsoluteOccupancy);
0200   intree->SetBranchAddress("APVMedianOccupancy", &APVMedianOccupancy);
0201
0202   for (int i = 0; i < intree->GetEntries(); i++) {
0203     intree->GetEntry(i);
0204
0205     double logMedianOccupancy = -1;
0206     double logAbsoluteOccupancy = -1;
0207
0208     if (APVMedianOccupancy > 0)
0209       logMedianOccupancy = log10(APVMedianOccupancy);
0210     if (APVAbsoluteOccupancy > 0)
0211       logAbsoluteOccupancy = log10(APVAbsoluteOccupancy);
0212
0213     // The layer/disk information is stored in Layer_Ring for TIB/TOB and in Disc for TID/TEC
0214     unsigned int layer = 0;
0215     if (SubDetId == 3 || SubDetId == 5)
0216       layer = Layer_Ring;
0217     else
0218       layer = Disc;
0219
0220     // Fill histograms for all the tracker
0221     medianVsAbsoluteOccupancy[0][0]->Fill(logAbsoluteOccupancy, logMedianOccupancy);
0222     medianOccupancy[0][0]->Fill(logMedianOccupancy);
0223     absoluteOccupancy[0][0]->Fill(logAbsoluteOccupancy);
0224     // Fill summary histograms for each subdetector
0225     medianVsAbsoluteOccupancy[SubDetId - 2][0]->Fill(logAbsoluteOccupancy, logMedianOccupancy);
0226     medianOccupancy[SubDetId - 2][0]->Fill(logMedianOccupancy);
0227     absoluteOccupancy[SubDetId - 2][0]->Fill(logAbsoluteOccupancy);
0228     // Fill histograms for each layer/disk
0229     medianVsAbsoluteOccupancy[SubDetId - 2][layer]->Fill(logAbsoluteOccupancy, logMedianOccupancy);
0230     medianOccupancy[SubDetId - 2][layer]->Fill(logMedianOccupancy);
0231     absoluteOccupancy[SubDetId - 2][layer]->Fill(logAbsoluteOccupancy);
0232   }
0233
0234   dqmStore->cd();
0235   dqmStore->save(outfilename, "ChannelStatusPlots");
0236 }
0237
0238 //define this as a plug-in
0239 DEFINE_FWK_MODULE(APVValidationPlots);