Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:26:20

0001 // File: ReadPixClusters.cc
0002 // Description: To test the pixel clusters.
0003 // Author: Danek Kotlinski
0004 // Creation Date:  Initial version. 3/06
0005 // Modify to work with CMSSW620, 8/13, CMSSW700, 10/13 d.k.
0006 // Add ByToken data access.
0007 //--------------------------------------------
0008 #include <memory>
0009 #include <string>
0010 #include <iostream>
0011 
0012 #include "DataFormats/Common/interface/Handle.h"
0013 #include "FWCore/Framework/interface/ESHandle.h"
0014 #include "FWCore/Framework/interface/EventSetup.h"
0015 
0016 #include "FWCore/Framework/interface/Frameworkfwd.h"
0017 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0018 
0019 #include "FWCore/Framework/interface/Event.h"
0020 #include "FWCore/Framework/interface/MakerMacros.h"
0021 
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "FWCore/ServiceRegistry/interface/Service.h"
0024 #include "FWCore/Utilities/interface/InputTag.h"
0025 
0026 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0027 #include "DataFormats/Common/interface/DetSetVector.h"
0028 #include "DataFormats/Common/interface/Ref.h"
0029 #include "DataFormats/DetId/interface/DetId.h"
0030 
0031 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0032 #include "DataFormats/TrackerCommon/interface/PixelBarrelName.h"
0033 
0034 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0035 #include "Geometry/CommonDetUnit/interface/PixelGeomDetType.h"
0036 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0037 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0038 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0039 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0040 #include "Geometry/CommonTopologies/interface/PixelTopology.h"
0041 
0042 // For L1
0043 #include "L1Trigger/GlobalTriggerAnalyzer/interface/L1GtUtils.h"
0044 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetupFwd.h"
0045 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetup.h"
0046 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
0047 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerObjectMapRecord.h"
0048 
0049 // For HLT
0050 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0051 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0052 #include "DataFormats/Common/interface/TriggerResults.h"
0053 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0054 #include "FWCore/Common/interface/TriggerNames.h"
0055 
0056 // To use root histos
0057 #include "FWCore/ServiceRegistry/interface/Service.h"
0058 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0059 
0060 // For ROOT
0061 #include <TROOT.h>
0062 #include <TChain.h>
0063 #include <TFile.h>
0064 #include <TF1.h>
0065 #include <TH2F.h>
0066 #include <TH1F.h>
0067 
0068 #define NEW_ID
0069 #ifdef NEW_ID
0070 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0071 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0072 #else
0073 #include "DataFormats/SiPixelDetId/interface/PXBDetId.h"
0074 #include "DataFormats/SiPixelDetId/interface/PXFDetId.h"
0075 #endif
0076 
0077 #define HISTOS
0078 
0079 #include "DataFormats/Provenance/interface/RunLumiEventNumber.h"
0080 
0081 using namespace std;
0082 
0083 //=============================================================================
0084 
0085 class ReadPixClusters : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0086 public:
0087   explicit ReadPixClusters(const edm::ParameterSet &conf);
0088   virtual ~ReadPixClusters();
0089   virtual void analyze(const edm::Event &e, const edm::EventSetup &c) override;
0090   virtual void beginJob() override;
0091   virtual void endJob() override;
0092 
0093 private:
0094   edm::InputTag src_;
0095   edm::EDGetTokenT<edmNew::DetSetVector<SiPixelCluster>> tPixelCluster;
0096   edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> trackerTopoToken_;
0097   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> trackerGeomToken_;
0098   bool printLocal;
0099   int countEvents, countAllEvents;
0100   double sumClusters;
0101 
0102   //TFile* hFile;
0103   TH1F *hdetunit;
0104   TH1F *hpixid, *hpixsubid, *hlayerid, *hladder1id, *hladder2id, *hladder3id, *hz1id, *hz2id, *hz3id;
0105 
0106   TH1F *hcharge1, *hcharge2, *hcharge3, *hcharge4, *hcharge5;
0107   TH1F *hpixcharge1, *hpixcharge2, *hpixcharge3, *hpixcharge4, *hpixcharge5;
0108   TH1F *hcols1, *hcols2, *hcols3, *hrows1, *hrows2, *hrows3;
0109   TH1F *hpcols1, *hpcols2, *hpcols3, *hprows1, *hprows2, *hprows3;
0110   TH1F *hsize1, *hsize2, *hsize3, *hsize4, *hsize5, *hsizex1, *hsizex2, *hsizex3, *hsizex4, *hsizex5, *hsizey1,
0111       *hsizey2, *hsizey3, *hsizey4, *hsizey5;
0112 
0113   TH1F *hclusPerDet1, *hclusPerDet2, *hclusPerDet3;
0114   TH1F *hpixPerDet1, *hpixPerDet2, *hpixPerDet3;
0115   TH1F *hpixPerLink1, *hpixPerLink2, *hpixPerLink3;
0116   TH1F *hclusPerLay1, *hclusPerLay2, *hclusPerLay3;
0117   TH1F *hpixPerLay1, *hpixPerLay2, *hpixPerLay3;
0118   TH1F *hdetsPerLay1, *hdetsPerLay2, *hdetsPerLay3;
0119   TH1F *hclus, *hclusBPix, *hclusFPix, *hdigis, *hdigisB, *hdigisF;
0120 
0121   TH1F *hdetr, *hdetz;
0122   //   TH1F *hcolsB,  *hrowsB,  *hcolsF,  *hrowsF;
0123   //   TH2F *htest1, *htest2;
0124   TH2F *hDetMap1, *hDetMap2, *hDetMap3;
0125   TH2F *hpDetMap1, *hpDetMap2, *hpDetMap3;
0126   TH2F *hpixDetMap1, *hpixDetMap2, *hpixDetMap3;
0127   TH2F *hcluDetMap1, *hcluDetMap2, *hcluDetMap3;
0128 
0129   TH1F *hevent, *horbit, *hlumi, *hlumi0, *hlumi1;
0130   TH1F *hbx, *hbx0, *hbx1;
0131   TH1F *hdets, *hmbits1, *hmbits2, *hmbits3, *hmaxPixPerDet;
0132 
0133   TH1F *hclusPerDisk1, *hclusPerDisk2, *hclusPerDisk3, *hclusPerDisk4;
0134   TH1F *htest;
0135 };
0136 
0137 /////////////////////////////////////////////////////////////////
0138 // Contructor, empty.
0139 ReadPixClusters::ReadPixClusters(edm::ParameterSet const &conf) : src_(conf.getParameter<edm::InputTag>("src")) {
0140   usesResource(TFileService::kSharedResource);
0141   printLocal = conf.getUntrackedParameter<bool>("Verbosity", false);
0142   //src_ =  conf.getParameter<edm::InputTag>( "src" );
0143   edm::LogPrint("ReadPixClusters") << " Construct " << printLocal;
0144   tPixelCluster = consumes<edmNew::DetSetVector<SiPixelCluster>>(src_);
0145   trackerTopoToken_ = esConsumes<TrackerTopology, TrackerTopologyRcd>();
0146   trackerGeomToken_ = esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>();
0147 }
0148 // Virtual destructor needed.
0149 ReadPixClusters::~ReadPixClusters() = default;
0150 
0151 // ------------ method called at the begining   ------------
0152 void ReadPixClusters::beginJob() {
0153   edm::LogPrint("ReadPixClusters") << "Initialize PixelClusterTest " << printLocal;
0154 
0155 #ifdef HISTOS
0156 
0157   // NEW way to use root (from 2.0.0?)
0158   edm::Service<TFileService> fs;
0159 
0160   //=====================================================================
0161 
0162   hladder1id = fs->make<TH1F>("hladder1id", "Ladder L1 id", 23, -11.5, 11.5);
0163   hladder2id = fs->make<TH1F>("hladder2id", "Ladder L2 id", 35, -17.5, 17.5);
0164   hladder3id = fs->make<TH1F>("hladder3id", "Ladder L3 id", 47, -23.5, 23.5);
0165   hz1id = fs->make<TH1F>("hz1id", "Z-index id L1", 11, -5.5, 5.5);
0166   hz2id = fs->make<TH1F>("hz2id", "Z-index id L2", 11, -5.5, 5.5);
0167   hz3id = fs->make<TH1F>("hz3id", "Z-index id L3", 11, -5.5, 5.5);
0168 
0169   int sizeH = 200;
0170   float lowH = -0.5;
0171   float highH = 199.5;
0172 
0173   hclusPerDet1 = fs->make<TH1F>("hclusPerDet1", "Clus per det l1", sizeH, lowH, highH);
0174   hclusPerDet2 = fs->make<TH1F>("hclusPerDet2", "Clus per det l2", sizeH, lowH, highH);
0175   hclusPerDet3 = fs->make<TH1F>("hclusPerDet3", "Clus per det l3", sizeH, lowH, highH);
0176 
0177   sizeH = 1000;
0178   highH = 1999.5;
0179   hpixPerDet1 = fs->make<TH1F>("hpixPerDet1", "Pix per det l1", sizeH, lowH, highH);
0180   hpixPerDet2 = fs->make<TH1F>("hpixPerDet2", "Pix per det l2", sizeH, lowH, highH);
0181   hpixPerDet3 = fs->make<TH1F>("hpixPerDet3", "Pix per det l3", sizeH, lowH, highH);
0182 
0183   sizeH = 1000;
0184   highH = 999.5;
0185   hpixPerLink1 = fs->make<TH1F>("hpixPerLink1", "Pix per link l1", sizeH, lowH, highH);
0186   hpixPerLink2 = fs->make<TH1F>("hpixPerLink2", "Pix per link l2", sizeH, lowH, highH);
0187   hpixPerLink3 = fs->make<TH1F>("hpixPerLink3", "Pix per link l3", sizeH, lowH, highH);
0188 
0189   sizeH = 2000;
0190   highH = 1999.5;
0191   hclusPerLay1 = fs->make<TH1F>("hclusPerLay1", "Clus per layer l1", sizeH, lowH, highH);
0192   hclusPerLay2 = fs->make<TH1F>("hclusPerLay2", "Clus per layer l2", sizeH, lowH, highH);
0193   hclusPerLay3 = fs->make<TH1F>("hclusPerLay3", "Clus per layer l3", sizeH, lowH, highH);
0194 
0195   hclusPerDisk1 = fs->make<TH1F>("hclusPerDisk1", "Clus per disk1", sizeH, lowH, highH);
0196   hclusPerDisk2 = fs->make<TH1F>("hclusPerDisk2", "Clus per disk2", sizeH, lowH, highH);
0197   hclusPerDisk3 = fs->make<TH1F>("hclusPerDisk3", "Clus per disk3", sizeH, lowH, highH);
0198   hclusPerDisk4 = fs->make<TH1F>("hclusPerDisk4", "Clus per disk4", sizeH, lowH, highH);
0199 
0200   sizeH = 2000;
0201   highH = 9999.5;
0202   hpixPerLay1 = fs->make<TH1F>("hpixPerLay1", "Pix per layer l1", sizeH, lowH, highH);
0203   hpixPerLay2 = fs->make<TH1F>("hpixPerLay2", "Pix per layer l2", sizeH, lowH, highH);
0204   hpixPerLay3 = fs->make<TH1F>("hpixPerLay3", "Pix per layer l3", sizeH, lowH, highH);
0205 
0206   hclus = fs->make<TH1F>("hclus", "Clus per event", sizeH, lowH, highH);
0207   hdigis = fs->make<TH1F>("hdigis", "All Digis in clus per event", 2000, lowH, 19999.5);
0208   hdigisB = fs->make<TH1F>("hdigisB", "BPix Digis in clus per event", 8000, lowH, 19999.5);
0209   hdigisF = fs->make<TH1F>("hdigisF", "FPix Digis in clus per event", 2000, lowH, 7999.5);
0210   hclusBPix = fs->make<TH1F>("hclusBPix", "Bpix Clus per event", 2000, 0., 2000.);
0211   hclusFPix = fs->make<TH1F>("hclusFPix", "Fpix Clus per event", 2000, 0., 2000.);
0212   hdets = fs->make<TH1F>("hdets", "Dets per event", 2000, -0.5, 1999.5);
0213 
0214   hmaxPixPerDet = fs->make<TH1F>("hmaxPixPerDet", "Max pixels per det", 1000, -0.5, 999.5);
0215 
0216   sizeH = 1000;
0217   highH = 1999.5;
0218 
0219   hdetsPerLay1 = fs->make<TH1F>("hdetsPerLay1", "Full dets per layer l1", 161, -0.5, 160.5);
0220   hdetsPerLay3 = fs->make<TH1F>("hdetsPerLay3", "Full dets per layer l3", 353, -0.5, 352.5);
0221   hdetsPerLay2 = fs->make<TH1F>("hdetsPerLay2", "Full dets per layer l2", 257, -0.5, 256.5);
0222 
0223   sizeH = 120;
0224   lowH = 0.;
0225   highH = 121.0;                                                             // charge limit in kelec
0226   hcharge1 = fs->make<TH1F>("hcharge1", "Clu charge l1", sizeH, 0., highH);  //in ke
0227   hcharge2 = fs->make<TH1F>("hcharge2", "Clu charge l2", sizeH, 0., highH);
0228   hcharge3 = fs->make<TH1F>("hcharge3", "Clu charge l3", sizeH, 0., highH);
0229   hcharge4 = fs->make<TH1F>("hcharge4", "Clu charge d1", sizeH, 0., highH);
0230   hcharge5 = fs->make<TH1F>("hcharge5", "Clu charge d2", sizeH, 0., highH);
0231 
0232   sizeH = 90;
0233   highH = 61.0;                                                                    // charge limit in kelec
0234   hpixcharge1 = fs->make<TH1F>("hpixcharge1", "Pix charge l1", sizeH, 0., highH);  //in ke
0235   hpixcharge2 = fs->make<TH1F>("hpixcharge2", "Pix charge l2", sizeH, 0., highH);
0236   hpixcharge3 = fs->make<TH1F>("hpixcharge3", "Pix charge l3", sizeH, 0., highH);
0237   hpixcharge4 = fs->make<TH1F>("hpixcharge4", "Pix charge d1", sizeH, 0., highH);
0238   hpixcharge5 = fs->make<TH1F>("hpixcharge5", "Pix charge d2", sizeH, 0., highH);
0239 
0240   hcols1 = fs->make<TH1F>("hcols1", "Layer 1 cols", 500, -0.5, 499.5);
0241   hcols2 = fs->make<TH1F>("hcols2", "Layer 2 cols", 500, -0.5, 499.5);
0242   hcols3 = fs->make<TH1F>("hcols3", "Layer 3 cols", 500, -0.5, 499.5);
0243 
0244   hrows1 = fs->make<TH1F>("hrows1", "Layer 1 rows", 200, -0.5, 199.5);
0245   hrows2 = fs->make<TH1F>("hrows2", "Layer 2 rows", 200, -0.5, 199.5);
0246   hrows3 = fs->make<TH1F>("hrows3", "layer 3 rows", 200, -0.5, 199.5);
0247 
0248   hpcols1 = fs->make<TH1F>("hpcols1", "Layer 1 pix cols", 500, -0.5, 499.5);
0249   hpcols2 = fs->make<TH1F>("hpcols2", "Layer 2 pix cols", 500, -0.5, 499.5);
0250   hpcols3 = fs->make<TH1F>("hpcols3", "Layer 3 pix cols", 500, -0.5, 499.5);
0251 
0252   hprows1 = fs->make<TH1F>("hprows1", "Layer 1 pix rows", 200, -0.5, 199.5);
0253   hprows2 = fs->make<TH1F>("hprows2", "Layer 2 pix rows", 200, -0.5, 199.5);
0254   hprows3 = fs->make<TH1F>("hprows3", "layer 3 pix rows", 200, -0.5, 199.5);
0255 
0256   sizeH = 1000;
0257   highH = 999.5;  // charge limit in kelec
0258   hsize1 = fs->make<TH1F>("hsize1", "layer 1 clu size", sizeH, -0.5, highH);
0259   hsize2 = fs->make<TH1F>("hsize2", "layer 2 clu size", sizeH, -0.5, highH);
0260   hsize3 = fs->make<TH1F>("hsize3", "layer 3 clu size", sizeH, -0.5, highH);
0261   hsize4 = fs->make<TH1F>("hsize4", "disk 1 clu size", sizeH, -0.5, highH);
0262   hsize5 = fs->make<TH1F>("hsize5", "disk 2 clu size", sizeH, -0.5, highH);
0263 
0264   hsizex1 = fs->make<TH1F>("hsizex1", "lay1 clu size in x", 10, -0.5, 9.5);
0265   hsizex2 = fs->make<TH1F>("hsizex2", "lay2 clu size in x", 10, -0.5, 9.5);
0266   hsizex3 = fs->make<TH1F>("hsizex3", "lay3 clu size in x", 10, -0.5, 9.5);
0267   hsizex4 = fs->make<TH1F>("hsizex4", "d1 clu size in x", 10, -0.5, 9.5);
0268   hsizex5 = fs->make<TH1F>("hsizex5", "d2 clu size in x", 10, -0.5, 9.5);
0269 
0270   hsizey1 = fs->make<TH1F>("hsizey1", "lay1 clu size in y", 20, -0.5, 19.5);
0271   hsizey2 = fs->make<TH1F>("hsizey2", "lay2 clu size in y", 20, -0.5, 19.5);
0272   hsizey3 = fs->make<TH1F>("hsizey3", "lay3 clu size in y", 20, -0.5, 19.5);
0273   hsizey4 = fs->make<TH1F>("hsizey4", "d1 clu size in y", 20, -0.5, 19.5);
0274   hsizey5 = fs->make<TH1F>("hsizey5", "d2 clu size in y", 20, -0.5, 19.5);
0275 
0276   hevent = fs->make<TH1F>("hevent", "event", 1000, 0, 10000000.);
0277   horbit = fs->make<TH1F>("horbit", "orbit", 100, 0, 100000000.);
0278 
0279   hlumi1 = fs->make<TH1F>("hlumi1", "lumi", 2000, 0, 2000.);
0280   hlumi0 = fs->make<TH1F>("hlumi0", "lumi", 2000, 0, 2000.);
0281   hlumi = fs->make<TH1F>("hlumi", "lumi", 2000, 0, 2000.);
0282   hbx1 = fs->make<TH1F>("hbx1", "bx", 4000, 0, 4000.);
0283   hbx0 = fs->make<TH1F>("hbx0", "bx", 4000, 0, 4000.);
0284   hbx = fs->make<TH1F>("hbx", "bx", 4000, 0, 4000.);
0285 
0286   hDetMap1 = fs->make<TH2F>("hDetMap1", " ", 9, -4.5, 4.5, 21, -10.5, 10.5);
0287   hDetMap1->SetOption("colz");
0288   hDetMap2 = fs->make<TH2F>("hDetMap2", " ", 9, -4.5, 4.5, 33, -16.5, 16.5);
0289   hDetMap2->SetOption("colz");
0290   hDetMap3 = fs->make<TH2F>("hDetMap3", " ", 9, -4.5, 4.5, 45, -22.5, 22.5);
0291   hDetMap3->SetOption("colz");
0292 
0293   hpDetMap1 = fs->make<TH2F>("hpDetMap1", " ", 9, -4.5, 4.5, 21, -10.5, 10.5);
0294   hpDetMap1->SetOption("colz");
0295   hpDetMap2 = fs->make<TH2F>("hpDetMap2", " ", 9, -4.5, 4.5, 33, -16.5, 16.5);
0296   hpDetMap2->SetOption("colz");
0297   hpDetMap3 = fs->make<TH2F>("hpDetMap3", " ", 9, -4.5, 4.5, 45, -22.5, 22.5);
0298   hpDetMap3->SetOption("colz");
0299 
0300   hpixDetMap1 = fs->make<TH2F>("hpixDetMap1", "pix det layer 1", 416, 0., 416., 160, 0., 160.);
0301   hpixDetMap2 = fs->make<TH2F>("hpixDetMap2", "pix det layer 2", 416, 0., 416., 160, 0., 160.);
0302   hpixDetMap3 = fs->make<TH2F>("hpixDetMap3", "pix det layer 3", 416, 0., 416., 160, 0., 160.);
0303 
0304   hcluDetMap1 = fs->make<TH2F>("hcluDetMap1", "clu det layer 1", 416, 0., 416., 160, 0., 160.);
0305   hcluDetMap2 = fs->make<TH2F>("hcluDetMap2", "clu det layer 1", 416, 0., 416., 160, 0., 160.);
0306   hcluDetMap3 = fs->make<TH2F>("hcluDetMap3", "clu det layer 1", 416, 0., 416., 160, 0., 160.);
0307   htest = fs->make<TH1F>("htest", "FPix R", 300, -15., 15.);
0308 
0309 #endif
0310 
0311   countEvents = 0;
0312   countAllEvents = 0;
0313   sumClusters = 0.;
0314 }
0315 // ------------ method called to at the end of the job  ------------
0316 void ReadPixClusters::endJob() {
0317   sumClusters = sumClusters / float(countEvents);
0318   edm::LogPrint("ReadPixClusters") << " End PixelClusTest, events all/with hits=  " << countAllEvents << "/"
0319                                    << countEvents << " " << sumClusters << " " << printLocal;
0320 }
0321 //////////////////////////////////////////////////////////////////
0322 // Functions that gets called by framework every event
0323 void ReadPixClusters::analyze(const edm::Event &e, const edm::EventSetup &es) {
0324   using namespace edm;
0325 
0326   // Get event setup
0327   edm::ESHandle<TrackerGeometry> geom = es.getHandle(trackerGeomToken_);
0328   const TrackerGeometry &theTracker(*geom);
0329 
0330 #ifdef NEW_ID
0331   //Retrieve tracker topology from geometry
0332   edm::ESHandle<TrackerTopology> tTopo = es.getHandle(trackerTopoToken_);
0333 #endif
0334 
0335   countAllEvents++;
0336   RunNumber_t const run = e.id().run();
0337   EventNumber_t const event = e.id().event();
0338   LuminosityBlockNumber_t const lumiBlock = e.luminosityBlock();
0339 
0340   int bx = e.bunchCrossing();
0341   int orbit = e.orbitNumber();
0342 
0343   // Get Cluster Collection from InputTag
0344   edm::Handle<edmNew::DetSetVector<SiPixelCluster>> clusters;
0345   e.getByToken(tPixelCluster, clusters);
0346 
0347   const edmNew::DetSetVector<SiPixelCluster> &input = *clusters;
0348   int numOf = input.size();
0349 
0350   hbx0->Fill(float(bx));
0351   hlumi0->Fill(float(lumiBlock));
0352 
0353   if (printLocal)
0354     edm::LogPrint("ReadPixClusters") << "run " << run << " event " << event << " bx " << bx << " lumi " << lumiBlock
0355                                      << " orbit " << orbit << " " << numOf;
0356 
0357   hdets->Fill(float(numOf));  // number of modules with pix
0358 
0359   // Select events with pixels
0360   //if(numOf<1) return; // skip events with  pixel dets
0361   //if(numOf<4) return; // skip events with few pixel dets
0362 
0363   hevent->Fill(float(event));
0364   hlumi->Fill(float(lumiBlock));
0365   hbx->Fill(float(bx));
0366   horbit->Fill(float(orbit));
0367 
0368   countEvents++;
0369   int numberOfClusters = 0;
0370   int numberOfPixels = 0;
0371   int numberOfDetUnits1 = 0;
0372   int numOfClustersPerDet1 = 0;
0373   int numOfClustersPerLay1 = 0;
0374   int numberOfDetUnits2 = 0;
0375   int numOfClustersPerDet2 = 0;
0376   int numOfClustersPerLay2 = 0;
0377   int numberOfDetUnits3 = 0;
0378   int numOfClustersPerDet3 = 0;
0379   int numOfClustersPerLay3 = 0;
0380 
0381   int numOfPixPerLay1 = 0;
0382   int numOfPixPerLay2 = 0;
0383   int numOfPixPerLay3 = 0;
0384 
0385   int numOfPixPerDet1 = 0;
0386   int numOfPixPerDet2 = 0;
0387   int numOfPixPerDet3 = 0;
0388 
0389   int numOfPixPerLink11 = 0;
0390   int numOfPixPerLink12 = 0;
0391   int numOfPixPerLink21 = 0;
0392   int numOfPixPerLink22 = 0;
0393   //SK:unused  int numOfPixPerLink3=0;
0394 
0395   int maxClusPerDet = 0;
0396   int maxPixPerDet = 0;
0397   unsigned int maxPixPerClu = 0;
0398 
0399   int numOfClustersPerDisk1 = 0;
0400   int numOfClustersPerDisk2 = 0;
0401   int numOfClustersPerDisk3 = 0;
0402   int numOfClustersPerDisk4 = 0;
0403   int numOfPixPerDisk1 = 0;
0404   int numOfPixPerDisk2 = 0;
0405   int numOfPixPerDisk3 = 0;
0406   int numOfPixPerDisk4 = 0;
0407 
0408   float aveCharge1 = 0., aveCharge2 = 0., aveCharge3 = 0., aveCharge4 = 0., aveCharge5 = 0.;
0409 
0410   static int module1[416][160] = {{0}};
0411   static int module2[416][160] = {{0}};
0412   static int module3[416][160] = {{0}};
0413 
0414   // get vector of detunit ids
0415   //--- Loop over detunits.
0416   edmNew::DetSetVector<SiPixelCluster>::const_iterator DSViter = input.begin();
0417   for (; DSViter != input.end(); DSViter++) {
0418     bool valid = false;
0419     unsigned int detid = DSViter->detId();
0420     // Det id
0421     DetId detId = DetId(detid);             // Get the Detid object
0422     unsigned int detType = detId.det();     // det type, pixel=1
0423     unsigned int subid = detId.subdetId();  //subdetector type, barrel=1
0424 
0425     if (printLocal)
0426       edm::LogPrint("ReadPixClusters") << "Det: " << detId.rawId() << " " << detId.null() << " " << detType << " "
0427                                        << subid;
0428 
0429 #ifdef HISTOS
0430       //hdetunit->Fill(float(detid));
0431       //hpixid->Fill(float(detType));
0432       //hpixsubid->Fill(float(subid));
0433 #endif  // HISTOS
0434 
0435     if (detType != 1)
0436       continue;  // look only at pixels
0437 
0438     //const GeomDetUnit * genericDet = geom->idToDet(detId);
0439     //const PixelGeomDetUnit * pixDet =
0440     //dynamic_cast<const PixelGeomDetUnit*>(genericDet);
0441 
0442     // Get the geom-detector
0443     const PixelGeomDetUnit *theGeomDet = dynamic_cast<const PixelGeomDetUnit *>(theTracker.idToDet(detId));
0444     double detZ = theGeomDet->surface().position().z();
0445     double detR = theGeomDet->surface().position().perp();
0446 
0447     //const BoundPlane& plane = theGeomDet->surface(); //for transf.
0448 
0449     //double detThick = theGeomDet->specificSurface().bounds().thickness();
0450     //int cols = theGeomDet->specificTopology().ncolumns();
0451     //int rows = theGeomDet->specificTopology().nrows();
0452 
0453     const PixelTopology *topol = &(theGeomDet->specificTopology());
0454 
0455     // barrel ids
0456     unsigned int layerC = 0;
0457     unsigned int ladderC = 0;
0458     unsigned int zindex = 0;
0459     int shell = 0;      // shell id // Shell { mO = 1, mI = 2 , pO =3 , pI =4 };
0460     int sector = 0;     // 1-8
0461     int ladder = 0;     // 1-22
0462     int layer = 0;      // 1-3
0463     int module = 0;     // 1-4
0464     bool half = false;  //
0465 
0466     // Endcap ids
0467     unsigned int disk = 0;     //1,2,3
0468     unsigned int blade = 0;    //1-24
0469     unsigned int moduleF = 0;  // plaquette 1,2,3,4
0470     unsigned int side = 0;     //size=1 for -z, 2 for +z
0471     unsigned int panel = 0;    //panel=1
0472 
0473     edmNew::DetSet<SiPixelCluster>::const_iterator clustIt;
0474 
0475     // Subdet id, pix barrel=1, forward=2
0476     if (subid == 2) {  // forward
0477 
0478 #ifdef NEW_ID
0479       disk = tTopo->pxfDisk(detid);      //1,2,3
0480       blade = tTopo->pxfBlade(detid);    //1-24
0481       zindex = tTopo->pxfModule(detid);  //
0482       side = tTopo->pxfSide(detid);      //size=1 for -z, 2 for +z
0483       panel = tTopo->pxfPanel(detid);    //panel=1
0484 #else
0485       PXFDetId pdetId = PXFDetId(detid);
0486       disk = pdetId.disk();       //1,2,3
0487       blade = pdetId.blade();     //1-24
0488       moduleF = pdetId.module();  // plaquette
0489       side = pdetId.side();       //size=1 for -z, 2 for +z
0490       panel = pdetId.panel();     //panel=1
0491 #endif
0492 
0493       if (printLocal)
0494         edm::LogPrint("ReadPixClusters") << " forward det, disk " << disk << ", blade " << blade << ", module "
0495                                          << moduleF << ", side " << side << ", panel " << panel << " pos = " << detZ
0496                                          << " " << detR;
0497 
0498       bool fpixInner = (((panel == 1) && (moduleF <= 2)) || ((panel == 2) && (moduleF <= 1)));  // make split at 10cm
0499       if (fpixInner)
0500         htest->Fill(detR);
0501 
0502     } else if (subid == 1) {  // barrel
0503 
0504 #ifdef HISTOS
0505       //hdetr->Fill(detR);
0506       //hdetz->Fill(detZ);
0507 #endif  // HISTOS
0508 
0509       //hcolsB->Fill(float(cols));
0510       //hrowsB->Fill(float(rows));
0511 
0512 #ifdef NEW_ID
0513       layerC = tTopo->pxbLayer(detid);
0514       ladderC = tTopo->pxbLadder(detid);
0515       zindex = tTopo->pxbModule(detid);
0516       PixelBarrelName pbn(detid);
0517 #else
0518       PXBDetId pdetId = PXBDetId(detid);
0519       //unsigned int detTypeP=pdetId.det();
0520       //unsigned int subidP=pdetId.subdetId();
0521       // Barell layer = 1,2,3
0522       layerC = pdetId.layer();
0523       // Barrel ladder id 1-20,32,44.
0524       ladderC = pdetId.ladder();
0525       // Barrel Z-index=1,8
0526       zindex = pdetId.module();
0527       // Convert to online
0528       PixelBarrelName pbn(pdetId);
0529 #endif
0530 
0531       // Shell { mO = 1, mI = 2 , pO =3 , pI =4 };
0532       PixelBarrelName::Shell sh = pbn.shell();  //enum
0533       sector = pbn.sectorName();
0534       ladder = pbn.ladderName();
0535       layer = pbn.layerName();
0536       module = pbn.moduleName();
0537       half = pbn.isHalfModule();
0538       shell = int(sh);
0539       // change the module sign for z<0
0540       if (shell == 1 || shell == 2)
0541         module = -module;
0542       // change ladeer sign for Outer )x<0)
0543       if (shell == 1 || shell == 3)
0544         ladder = -ladder;
0545 
0546       if (printLocal) {
0547         edm::LogPrint("ReadPixClusters") << " Barrel layer, ladder, module " << layerC << " " << ladderC << " "
0548                                          << zindex << " " << sh << "(" << shell << ") " << sector << " " << layer << " "
0549                                          << ladder << " " << module << " " << half;
0550         //edm::LogPrint("ReadPixClusters")<<" Barrel det, thick "<<detThick<<" "
0551         //  <<" layer, ladder, module "
0552         //  <<layer<<" "<<ladder<<" "<<zindex;
0553         //edm::LogPrint("ReadPixClusters")<<" col/row, pitch "<<cols<<" "<<rows<<" "
0554         //  <<pitchX<<" "<<pitchY;
0555       }
0556 
0557     }  // if subid
0558 
0559     if (printLocal) {
0560       edm::LogPrint("ReadPixClusters") << "List clusters : ";
0561       edm::LogPrint("ReadPixClusters") << "Num Charge Size SizeX SizeY X Y Xmin Xmax Ymin Ymax Edge";
0562     }
0563 
0564     // Loop over clusters
0565     for (clustIt = DSViter->begin(); clustIt != DSViter->end(); clustIt++) {
0566       sumClusters++;
0567       numberOfClusters++;
0568       float ch = float(clustIt->charge()) / 1000.;  // convert ke to electrons
0569       int size = clustIt->size();
0570       int sizeX = clustIt->sizeX();  //x=row=rfi, starting from CMSSW6 this is limited to 64
0571       int sizeY = clustIt->sizeY();  //y=col=z_global, starting from CMSSW6 this is limited to 64
0572       float x = clustIt->x();        // cluster position as float (int+0.5)
0573       float y = clustIt->y();        // analog average
0574 
0575       // Returns int index of the cluster min/max
0576       // This parameters are affected by the 64*64 limit on the cluster grid, so max-min<=64.
0577       int minPixelRow = clustIt->minPixelRow();  //x
0578       int maxPixelRow = clustIt->maxPixelRow();
0579       int minPixelCol = clustIt->minPixelCol();  //y
0580       int maxPixelCol = clustIt->maxPixelCol();
0581 
0582       //unsigned int geoId = clustIt->geographicalId(); // always 0?!
0583 
0584       // edge method moved to topologu class
0585       bool edgeHitX = (topol->isItEdgePixelInX(minPixelRow)) || (topol->isItEdgePixelInX(maxPixelRow));
0586       bool edgeHitY = (topol->isItEdgePixelInY(minPixelCol)) || (topol->isItEdgePixelInY(maxPixelCol));
0587 
0588       if (printLocal)
0589         edm::LogPrint("ReadPixClusters") << numberOfClusters << " " << ch << " " << size << " " << sizeX << " " << sizeY
0590                                          << " " << x << " " << y << " " << minPixelRow << " " << maxPixelRow << " "
0591                                          << minPixelCol << " " << maxPixelCol << " " << edgeHitX << " " << edgeHitY;
0592 
0593       // Get the pixels in the Cluster
0594       const vector<SiPixelCluster::Pixel> &pixelsVec = clustIt->pixels();
0595       if (printLocal)
0596         edm::LogPrint("ReadPixClusters") << " Pixels in this cluster ";
0597 
0598       bool bigInX = false, bigInY = false;
0599       bool edgeHitX2 = false;  // edge method moved
0600       bool edgeHitY2 = false;  // to topologu class
0601       //bool edgeInX = false; // edge method moved
0602       //bool edgeInY = false; // to topology class
0603       //SK:unused      bool cluBigInX = false; // does this clu include a big pixel
0604       //SK:unused      bool cluBigInY = false; // does this clu include a big pixel
0605       //int noisy = 0;
0606 
0607       if (pixelsVec.size() > maxPixPerClu)
0608         maxPixPerClu = pixelsVec.size();
0609 
0610       // Look at pixels in this cluster. ADC is calibrated, in electrons
0611       for (unsigned int i = 0; i < pixelsVec.size(); ++i) {  // loop over pixels
0612         numberOfPixels++;
0613         float pixx = pixelsVec[i].x;  // index as float=iteger, row index, 0-159
0614         float pixy = pixelsVec[i].y;  // same, col index, 0-415
0615         float adc = (float(pixelsVec[i].adc) / 1000.);
0616 
0617         //int chan = PixelChannelIdentifier::pixelToChannel(int(pixx),int(pixy));
0618         //bool binInX = (RectangularPixelTopology::isItBigPixelInX(int(pixx)));
0619         //bool bigInY = (RectangularPixelTopology::isItBigPixelInY(int(pixy)));
0620         //int roc = rocId(int(pixy),int(pixx));  // column, row
0621 
0622 #ifdef HISTOS
0623         // Pixel histos
0624         if (subid == 1) {  // barrel
0625           if (layer == 1) {
0626             numOfPixPerDet1++;
0627             numOfPixPerLay1++;
0628             valid = valid || true;
0629             hpixcharge1->Fill(adc);
0630             hpixDetMap1->Fill(pixy, pixx);
0631             hpDetMap1->Fill(float(module), float(ladder));
0632             module1[int(pixx)][int(pixy)]++;
0633 
0634             hpcols1->Fill(pixy);
0635             hprows1->Fill(pixx);
0636 
0637             if (pixx < 80.)
0638               numOfPixPerLink11++;
0639             else
0640               numOfPixPerLink12++;
0641 
0642           } else if (layer == 2) {
0643             numOfPixPerDet2++;
0644             numOfPixPerLay2++;
0645 
0646             hpcols2->Fill(pixy);
0647             hprows2->Fill(pixx);
0648 
0649             if (pixx < 80.)
0650               numOfPixPerLink21++;
0651             else
0652               numOfPixPerLink22++;
0653 
0654             hpixcharge2->Fill(adc);
0655             hpixDetMap2->Fill(pixy, pixx);
0656             hpDetMap2->Fill(float(module), float(ladder));
0657             module2[int(pixx)][int(pixy)]++;
0658 
0659           } else if (layer == 3) {
0660             numOfPixPerDet3++;
0661             numOfPixPerLay3++;
0662             valid = valid || true;
0663             hpixcharge3->Fill(adc);
0664             hpixDetMap3->Fill(pixy, pixx);
0665             hpDetMap3->Fill(float(module), float(ladder));
0666             module3[int(pixx)][int(pixy)]++;
0667 
0668             hpcols3->Fill(pixy);
0669             hprows3->Fill(pixx);
0670 
0671           }  // if layer
0672 
0673         } else if (subid == 2) {  // endcap
0674           // pixels
0675 
0676           if (disk == 1) {  // disk1 -+z
0677             if (side == 1)
0678               numOfPixPerDisk2++;  // d1,-z
0679             else if (side == 2)
0680               numOfPixPerDisk3++;  // d1, +z
0681             else
0682               edm::LogPrint("ReadPixClusters") << " unknown side " << side;
0683 
0684             hpixcharge4->Fill(adc);
0685 
0686           } else if (disk == 2) {  // disk2 -+z
0687 
0688             if (side == 1)
0689               numOfPixPerDisk1++;  // d2, -z
0690             else if (side == 2)
0691               numOfPixPerDisk4++;  // d2, +z
0692             else
0693               edm::LogPrint("ReadPixClusters") << " unknown side " << side;
0694 
0695             hpixcharge5->Fill(adc);
0696 
0697           } else
0698             edm::LogPrint("ReadPixClusters") << " unknown disk " << disk;
0699 
0700         }  // end if subdet (pixel loop)
0701 
0702 #endif  // HISTOS
0703 
0704         // This looks if there is an edge pixel in the cluster
0705         bool edgeInX = topol->isItEdgePixelInX(int(pixx));
0706         bool edgeInY = topol->isItEdgePixelInY(int(pixy));
0707         if (edgeInX)
0708           edgeHitX2 = true;
0709         if (edgeInY)
0710           edgeHitY2 = true;
0711 
0712         if (printLocal)
0713           edm::LogPrint("ReadPixClusters") << i << " " << pixx << " " << pixy << " " << adc << " " << bigInX << " "
0714                                            << bigInY << " " << edgeInX << " " << edgeInY;
0715 
0716         //SK:unused if(bigInX) cluBigInX=true;
0717         //SK:unused if(bigInY) cluBigInY=true;
0718 
0719       }  // pixel loop
0720 
0721       // if(edgeHitX ||  edgeHitX2 || edgeHitY ||  edgeHitY2)
0722       //    edm::LogPrint("ReadPixClusters")<<" egde "<<edgeHitX<<" "<<edgeHitX2<<" "<<edgeHitY<<" "<<edgeHitY2<<" "
0723       //        <<size<<" "<<sizeX<<" "<<sizeY<<" "
0724       //        <<x<<" "<<y<<" "<<minPixelRow<<" "<<maxPixelRow<<" "<<minPixelCol<<" "
0725       //        <<maxPixelCol<<" "<<edgeHitX<<" "<<edgeHitY;
0726 
0727       // This will happen  for clusters wider than 64
0728       if ((edgeHitX != edgeHitX2) && sizeX < 64)
0729         edm::LogPrint("ReadPixClusters") << " wrong egdeX " << edgeHitX << " "
0730                                          << edgeHitX2
0731                                          //<<" "<<event<<" "<<detid<<" "<<numberOfClusters
0732                                          << " size = " << size << " " << sizeX << " " << sizeY << " " << x << " " << y
0733                                          << " " << minPixelRow << " " << maxPixelRow << " " << minPixelCol << " "
0734                                          << maxPixelCol << " " << edgeHitX << " " << edgeHitY;
0735 
0736       if ((edgeHitY != edgeHitY2) && sizeY < 64)
0737         edm::LogPrint("ReadPixClusters") << " wrong egdeY " << edgeHitY << " "
0738                                          << edgeHitY2
0739                                          //<<" "<<event<<" "<<detid<<" "<<numberOfClusters
0740                                          << " size = " << size << " " << sizeX << " " << sizeY << " " << x << " " << y
0741                                          << " " << minPixelRow << " " << maxPixelRow << " " << minPixelCol << " "
0742                                          << maxPixelCol << " " << edgeHitX << " " << edgeHitY;
0743 
0744 #ifdef HISTOS
0745 
0746       // Cluster histos
0747       if (subid == 1) {  // barrel
0748         //if (subid==1) {  // barrel
0749         if (layer == 1) {  // layer
0750 
0751           hDetMap1->Fill(float(module), float(ladder));
0752           hcluDetMap1->Fill(y, x);
0753           hcharge1->Fill(ch);
0754           hcols1->Fill(y);
0755           hrows1->Fill(x);
0756           hsize1->Fill(float(size));
0757           hsizex1->Fill(float(sizeX));
0758           hsizey1->Fill(float(sizeY));
0759           numOfClustersPerDet1++;
0760           numOfClustersPerLay1++;
0761 
0762           aveCharge1 += ch;
0763 
0764         } else if (layer == 2) {
0765           hDetMap2->Fill(float(module), float(ladder));
0766           hcluDetMap2->Fill(y, x);
0767           hcharge2->Fill(ch);
0768           hcols2->Fill(y);
0769           hrows2->Fill(x);
0770           hsize2->Fill(float(size));
0771           hsizex2->Fill(float(sizeX));
0772           hsizey2->Fill(float(sizeY));
0773           numOfClustersPerDet2++;
0774           numOfClustersPerLay2++;
0775 
0776           aveCharge2 += ch;
0777 
0778         } else if (layer == 3) {
0779           hDetMap3->Fill(float(module), float(ladder));
0780           hcluDetMap3->Fill(y, x);
0781           hcharge3->Fill(ch);
0782           hcols3->Fill(y);
0783           hrows3->Fill(x);
0784           hsize3->Fill(float(size));
0785           hsizex3->Fill(float(sizeX));
0786           hsizey3->Fill(float(sizeY));
0787           numOfClustersPerDet3++;
0788           numOfClustersPerLay3++;
0789 
0790           aveCharge3 += ch;
0791 
0792         }  // end if bpix layer
0793 
0794       } else if (subid == 2) {  // endcap
0795 
0796         //edm::LogPrint("ReadPixClusters")<<disk<<" "<<side;
0797         if (disk == 1) {  // disk1 -+z
0798           if (side == 1)
0799             numOfClustersPerDisk2++;  // d1,-z
0800           else if (side == 2)
0801             numOfClustersPerDisk3++;  // d1, +z
0802           else
0803             edm::LogPrint("ReadPixClusters") << " unknown side " << side;
0804 
0805           hcharge4->Fill(ch);
0806           aveCharge4 += ch;
0807           hsize4->Fill(float(size));
0808           hsizex4->Fill(float(sizeX));
0809           hsizey4->Fill(float(sizeY));
0810 
0811         } else if (disk == 2) {  // disk2 -+z
0812 
0813           if (side == 1)
0814             numOfClustersPerDisk1++;  // d2, -z
0815           else if (side == 2)
0816             numOfClustersPerDisk4++;  // d2, +z
0817           else
0818             edm::LogPrint("ReadPixClusters") << " unknown side " << side;
0819 
0820           hcharge5->Fill(ch);
0821           aveCharge5 += ch;
0822           hsize5->Fill(float(size));
0823           hsizex5->Fill(float(sizeX));
0824           hsizey5->Fill(float(sizeY));
0825 
0826         } else
0827           edm::LogPrint("ReadPixClusters") << " unknown disk " << disk;  // end fpix disk
0828 
0829       }  // end if barrel/forward
0830 
0831 #endif  // HISTOS
0832 
0833     }  // clusters
0834 
0835     if (numOfClustersPerDet1 > maxClusPerDet)
0836       maxClusPerDet = numOfClustersPerDet1;
0837     if (numOfClustersPerDet2 > maxClusPerDet)
0838       maxClusPerDet = numOfClustersPerDet2;
0839     if (numOfClustersPerDet3 > maxClusPerDet)
0840       maxClusPerDet = numOfClustersPerDet3;
0841 
0842     if (printLocal) {
0843       if (layer == 1)
0844         edm::LogPrint("ReadPixClusters") << "Lay1: number of clusters per det = " << numOfClustersPerDet1;
0845       else if (layer == 2)
0846         edm::LogPrint("ReadPixClusters") << "Lay2: number of clusters per det = " << numOfClustersPerDet1;
0847       else if (layer == 3)
0848         edm::LogPrint("ReadPixClusters") << "Lay3: number of clusters per det = " << numOfClustersPerDet1;
0849     }  // end if printLocal
0850 
0851 #ifdef HISTOS
0852     if (subid == 1) {  // barrel
0853 
0854       // Det histos
0855       if (layer == 1) {
0856         hladder1id->Fill(float(ladder));
0857         hz1id->Fill(float(module));
0858         ++numberOfDetUnits1;
0859         hclusPerDet1->Fill(float(numOfClustersPerDet1));
0860         hpixPerDet1->Fill(float(numOfPixPerDet1));
0861         if (numOfPixPerDet1 > maxPixPerDet)
0862           maxPixPerDet = numOfPixPerDet1;
0863         numOfClustersPerDet1 = 0;
0864         numOfPixPerDet1 = 0;
0865         //if(numOfPixPerLink11>798 || numOfPixPerLink12>798) select=true;
0866         hpixPerLink1->Fill(float(numOfPixPerLink11));
0867         hpixPerLink1->Fill(float(numOfPixPerLink12));
0868         numOfPixPerLink11 = 0;
0869         numOfPixPerLink12 = 0;
0870 
0871       } else if (layer == 2) {
0872         hladder2id->Fill(float(ladder));
0873         hz2id->Fill(float(module));
0874         ++numberOfDetUnits2;
0875         hclusPerDet2->Fill(float(numOfClustersPerDet2));
0876         hpixPerDet2->Fill(float(numOfPixPerDet2));
0877         if (numOfPixPerDet2 > maxPixPerDet)
0878           maxPixPerDet = numOfPixPerDet2;
0879         numOfClustersPerDet2 = 0;
0880         numOfPixPerDet2 = 0;
0881         hpixPerLink2->Fill(float(numOfPixPerLink21));
0882         hpixPerLink2->Fill(float(numOfPixPerLink22));
0883         numOfPixPerLink21 = 0;
0884         numOfPixPerLink22 = 0;
0885 
0886       } else if (layer == 3) {
0887         hladder3id->Fill(float(ladder));
0888         hz3id->Fill(float(module));
0889         ++numberOfDetUnits3;
0890         hclusPerDet3->Fill(float(numOfClustersPerDet3));
0891         hpixPerDet3->Fill(float(numOfPixPerDet3));
0892         if (numOfPixPerDet3 > maxPixPerDet)
0893           maxPixPerDet = numOfPixPerDet3;
0894         numOfClustersPerDet3 = 0;
0895         numOfPixPerDet3 = 0;
0896         //SK:unused numOfPixPerLink3=0;
0897 
0898       }  // layer
0899 
0900     }  // end barrel/forward
0901 
0902 #endif  // HISTOS
0903 
0904   }  // detunits loop
0905 
0906   if (0) {
0907     edm::LogPrint("ReadPixClusters") << "run " << run << " event " << event << " bx " << bx << " lumi " << lumiBlock
0908                                      << " orbit " << orbit << " num " << countEvents;
0909     edm::LogPrint("ReadPixClusters") << "Num of pix " << numberOfPixels << " num of clus " << numberOfClusters
0910                                      << " num of dets " << numOf << " max clus per det " << maxClusPerDet
0911                                      << " max pix per clu " << maxPixPerClu << " count " << countEvents;
0912     edm::LogPrint("ReadPixClusters") << "Number of clusters per      Lay1,2,3: " << numOfClustersPerLay1 << " "
0913                                      << numOfClustersPerLay2 << " " << numOfClustersPerLay3;
0914     edm::LogPrint("ReadPixClusters") << "Number of pixels per        Lay1,2,3: " << numOfPixPerLay1 << " "
0915                                      << numOfPixPerLay2 << " " << numOfPixPerLay3;
0916     edm::LogPrint("ReadPixClusters") << "Number of dets with clus in Lay1,2,3: " << numberOfDetUnits1 << " "
0917                                      << numberOfDetUnits2 << " " << numberOfDetUnits3;
0918     edm::LogPrint("ReadPixClusters") << "Number of clus in disks 1,2,3,4     : " << numOfClustersPerDisk1 << " "
0919                                      << numOfClustersPerDisk2 << " " << numOfClustersPerDisk3 << " "
0920                                      << numOfClustersPerDisk4;
0921     aveCharge1 /= float(numOfClustersPerLay1);
0922     aveCharge2 /= float(numOfClustersPerLay2);
0923     aveCharge3 /= float(numOfClustersPerLay3);
0924     aveCharge4 /= float(numOfClustersPerDisk2 + numOfClustersPerDisk3);
0925     aveCharge5 /= float(numOfClustersPerDisk1 + numOfClustersPerDisk4);
0926     edm::LogPrint("ReadPixClusters") << " Average charge " << aveCharge1 << " " << aveCharge2 << " " << aveCharge3
0927                                      << " " << aveCharge4 << " " << aveCharge5;
0928   }
0929 
0930 #ifdef HISTOS
0931 
0932   hdigis->Fill(float(numberOfPixels));  // all pix
0933   int tmp1 = numOfPixPerLay1 + numOfPixPerLay2 + numOfPixPerLay3;
0934   hdigisB->Fill(float(tmp1));  // pix in bpix
0935   tmp1 = numOfPixPerDisk1 + numOfPixPerDisk2 + numOfPixPerDisk3 + numOfPixPerDisk4;
0936   hdigisF->Fill(float(tmp1));  // pix in fpix
0937 
0938   hclus->Fill(float(numberOfClusters));  // clusters fpix+bpix
0939   tmp1 = numOfClustersPerLay1 + numOfClustersPerLay2 + numOfClustersPerLay3;
0940   hclusBPix->Fill(float(tmp1));  // clusters in bpix
0941   int tmp2 = numOfClustersPerDisk1 + numOfClustersPerDisk2 + numOfClustersPerDisk3 + numOfClustersPerDisk4;
0942   hclusFPix->Fill(float(tmp2));  // clusters in fpix
0943 
0944   hclusPerLay1->Fill(float(numOfClustersPerLay1));
0945   hclusPerLay2->Fill(float(numOfClustersPerLay2));
0946   hclusPerLay3->Fill(float(numOfClustersPerLay3));
0947   hpixPerLay1->Fill(float(numOfPixPerLay1));
0948   hpixPerLay2->Fill(float(numOfPixPerLay2));
0949   hpixPerLay3->Fill(float(numOfPixPerLay3));
0950   if (numOfClustersPerLay1 > 0)
0951     hdetsPerLay1->Fill(float(numberOfDetUnits1));
0952   if (numOfClustersPerLay2 > 0)
0953     hdetsPerLay2->Fill(float(numberOfDetUnits2));
0954   if (numOfClustersPerLay3 > 0)
0955     hdetsPerLay3->Fill(float(numberOfDetUnits3));
0956 
0957   hclusPerDisk1->Fill(float(numOfClustersPerDisk1));
0958   hclusPerDisk2->Fill(float(numOfClustersPerDisk2));
0959   hclusPerDisk3->Fill(float(numOfClustersPerDisk3));
0960   hclusPerDisk4->Fill(float(numOfClustersPerDisk4));
0961 
0962   hmaxPixPerDet->Fill(float(maxPixPerDet));
0963 
0964   // Select only events with more tha 3 clusters
0965   if (numberOfClusters > 3) {
0966     hbx1->Fill(float(bx));
0967     hlumi1->Fill(float(lumiBlock));
0968   }  // if num of clusters
0969 
0970 #endif  // HISTOS
0971 
0972 }  // end
0973 
0974 //define this as a plug-in
0975 DEFINE_FWK_MODULE(ReadPixClusters);