SiPixelDigitizer/test/PixelDigisTest.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:    PixelDigisTest
0004 // Class:      PixelDigisTest
0005 //
0006 /**\class PixelDigisTest PixelDigisTest.cc
0007
0008  Description: Test pixel digis.
0009  Barrel & Forward digis. Uses root histos.
0010  Adopted for the new simLinks.
0011  Added the online detector index. d.k. 11/09
0012  Works with CMSSW_7
0013  New detector ID.
0014  Modified to use "byToken"
0015
0016 */
0017 //
0018 // Original Author:  d.k.
0019 //         Created:  Jan CET 2006
0020 //
0021 //
0022 // system include files
0023 #include <memory>
0024
0025 // user include files
0026 #include "FWCore/Framework/interface/Frameworkfwd.h"
0027 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0028
0029 #include "FWCore/Framework/interface/Event.h"
0030 #include "FWCore/Framework/interface/MakerMacros.h"
0031 #include "DataFormats/Common/interface/Handle.h"
0032
0033 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0034 #include "FWCore/ServiceRegistry/interface/Service.h"
0035 #include "FWCore/Framework/interface/ESHandle.h"
0036 #include "FWCore/Utilities/interface/InputTag.h"
0037
0038 // my includes
0039 //#include "Geometry/CommonDetUnit/interface/GeomDet.h"
0040
0041 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0042
0043 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0044 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0045 #include "Geometry/CommonDetUnit/interface/PixelGeomDetType.h"
0046
0047 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0048 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0049 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0050 #include "DataFormats/TrackerCommon/interface/PixelBarrelName.h"
0051
0052 // data formats
0053 #include "DataFormats/Common/interface/DetSetVector.h"
0054 #include "DataFormats/SiPixelDigi/interface/PixelDigi.h"
0055 //#include "DataFormats/SiPixelDigi/interface/PixelDigiCollection.h"
0056 #include "DataFormats/DetId/interface/DetId.h"
0057
0058 #include "SimDataFormats/TrackerDigiSimLink/interface/PixelDigiSimLink.h"
0059
0060 // For the big pixel recongnition
0061 #include "Geometry/TrackerGeometryBuilder/interface/RectangularPixelTopology.h"
0062
0063 // for simulated Tracker hits
0064 //#include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0065
0066 // For L1
0067 #include "L1Trigger/GlobalTriggerAnalyzer/interface/L1GtUtils.h"
0068 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetupFwd.h"
0069 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetup.h"
0070 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
0071 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerObjectMapRecord.h"
0072
0073 // For HLT
0074 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0075 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0076 #include "DataFormats/Common/interface/TriggerResults.h"
0077 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0078 #include "FWCore/Common/interface/TriggerNames.h"
0079
0080 // To use root histos
0081 #include "FWCore/ServiceRegistry/interface/Service.h"
0082 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0083
0084 // For ROOT
0085 #include <TROOT.h>
0086 #include <TChain.h>
0087 #include <TFile.h>
0088 #include <TF1.h>
0089 #include <TH2F.h>
0090 #include <TH1F.h>
0091
0092 #define HISTOS
0093 #define L1
0094 #define HLT
0095
0096 using namespace std;
0097
0098 // Enable this to look at simlinks (link simhit->digis)
0099 // Can be used only with simulated data.
0100 //#define USE_SIM_LINKS
0101
0102 //
0103 // class decleration
0104 //
0105
0106 class PixelDigisTest : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0107 public:
0108   explicit PixelDigisTest(const edm::ParameterSet &);
0109   ~PixelDigisTest() override;
0110   void beginJob() override;
0111   void analyze(const edm::Event &, const edm::EventSetup &) override;
0112   void endJob() override;
0113
0114 private:
0115   // ----------member data ---------------------------
0116   bool PRINT;
0117   edm::EDGetTokenT<edm::DetSetVector<PixelDigi>> tPixelDigi;
0118
0119   edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> trackerTopoToken;
0120   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> trackerGeomToken;
0121
0122 #ifdef USE_SIM_LINKS
0123   edm::EDGetTokenT<edm::DetSetVector<PixelDigiSimLink>> tPixelDigiSimLink;
0124 #endif
0125
0126 #ifdef HISTOS
0127
0128   //TFile* hFile;
0129   TH1F *hdetunit;
0130   TH1F *heloss1, *heloss2, *heloss3;
0131   TH1F *hneloss1, *hneloss2, *hneloss3;
0132   TH1F *helossF1, *helossF2;
0133   TH1F *hpixid, *hpixsubid, *hlayerid, *hshellid, *hsectorid, *hladder1id, *hladder2id, *hladder3id, *hz1id, *hz2id,
0134       *hz3id;
0135   TH1F *hcols1, *hcols2, *hcols3, *hrows1, *hrows2, *hrows3;
0136   TH1F *hcolsF1, *hcolsF2, *hcolsF3, *hrowsF1, *hrowsF2, *hrowsF3;
0137   TH1F *hdigisPerDet1, *hdigisPerDet2, *hdigisPerDet3;
0138   TH1F *hdigisPerLay1, *hdigisPerLay2, *hdigisPerLay3;
0139   TH1F *hdetsPerLay1, *hdetsPerLay2, *hdetsPerLay3;
0140   TH1F *hdigisPerDetF1, *hdigisPerDetF2, *hdigisPerDetF3;
0141   TH1F *hdigisPerLayF1, *hdigisPerLayF2, *hdigisPerLayF3;
0142   TH1F *hdetsPerLayF1, *hdetsPerLayF2, *hdetsPerLayF3;
0143   TH1F *hdetr, *hdetz, *hdetrF, *hdetzF;
0144   TH1F *hcolsB, *hrowsB, *hcolsF, *hrowsF;
0145   TH1F *hcols1big, *hrows1big, *heloss1bigx, *heloss1bigy;
0146   TH1F *hsimlinks, *hfract;
0147   TH1F *hblade1, *hblade2;
0148
0149   //TH2F *htest, *htest2;
0150   TH2F *hdetMap3, *hdetMap2, *hdetMap1, *hpixMap1, *hpixMap2, *hpixMap3, *hpixMapNoise;
0151
0152   TH1F *hevent, *hlumi, *horbit, *hbx0, *hlumi0, *hlumi1, *hbx1, *hbx2, *hbx3, *hbx4, *hbx5, *hbx6;
0153   TH1F *hdets, *hdigis, *hdigis1, *hdigis2, *hdigis3, *hdigis4, *hdigis5;
0154
0155 #endif
0156
0157   edm::InputTag src_;
0158 };
0159
0160 //
0161 // constants, enums and typedefs
0162 //
0163
0164 //
0165 // static data member definitions
0166 //
0167
0168 //
0169 // constructors and destructor
0170 //
0171 PixelDigisTest::PixelDigisTest(const edm::ParameterSet &iConfig) {
0172   usesResource(TFileService::kSharedResource);
0173
0174   //We put this here for the moment since there is no better place
0175   //edm::Service<MonitorDaemon> daemon;
0176   //daemon.operator->();
0177
0178   PRINT = iConfig.getUntrackedParameter<bool>("Verbosity", false);
0179   src_ = iConfig.getParameter<edm::InputTag>("src");
0180   tPixelDigi = consumes<edm::DetSetVector<PixelDigi>>(src_);
0181
0182   trackerTopoToken = esConsumes<TrackerTopology, TrackerTopologyRcd>();
0183   trackerGeomToken = esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>();
0184
0185 #ifdef USE_SIM_LINKS
0186   tPixelDigiSimLink = consumes<edm::DetSetVector<PixelDigiSimLink>>(src_);
0187 #endif
0188
0189   cout << " Construct PixelDigisTest " << endl;
0190 }
0191
0192 PixelDigisTest::~PixelDigisTest() {
0193   // do anything here that needs to be done at desctruction time
0194   // (e.g. close files, deallocate resources etc.)
0195   cout << " Destroy PixelDigisTest " << endl;
0196 }
0197
0198 //
0199 // member functions
0200 //
0201 // ------------ method called at the begining   ------------
0202 void PixelDigisTest::beginJob() {
0203   using namespace edm;
0204   cout << "Initialize PixelDigisTest " << endl;
0205
0206 #ifdef HISTOS
0207   edm::Service<TFileService> fs;
0208
0209   // Histos go to a subdirectory "PixRecHits")
0210   //TFileDirectory subDir = fs->mkdir( "mySubDirectory" );
0211   //TFileDirectory subSubDir = subDir.mkdir( "mySubSubDirectory" );
0212
0213   // put here whatever you want to do at the beginning of the job
0214   //hFile = new TFile ( "histo.root", "RECREATE" );
0215
0216   hdetunit = fs->make<TH1F>("hdetunit", "Det unit", 1000, 302000000., 302300000.);
0217   hpixid = fs->make<TH1F>("hpixid", "Pix det id", 10, 0., 10.);
0218   hpixsubid = fs->make<TH1F>("hpixsubid", "Pix Barrel id", 10, 0., 10.);
0219   hlayerid = fs->make<TH1F>("hlayerid", "Pix layer id", 10, 0., 10.);
0220   hsectorid = fs->make<TH1F>("hsectorid", "Pix sector ", 10, 0., 10.);
0221   hshellid = fs->make<TH1F>("hshellid", "Shell", 5, 0., 5.);
0222   hladder1id = fs->make<TH1F>("hladder1id", "Ladder L1 id", 23, -11.5, 11.5);
0223   hladder2id = fs->make<TH1F>("hladder2id", "Ladder L2 id", 35, -17.5, 17.5);
0224   hladder3id = fs->make<TH1F>("hladder3id", "Ladder L3 id", 47, -23.5, 23.5);
0225   hz1id = fs->make<TH1F>("hz1id", "Z-index id L1", 11, -5.5, 5.5);
0226   hz2id = fs->make<TH1F>("hz2id", "Z-index id L2", 11, -5.5, 5.5);
0227   hz3id = fs->make<TH1F>("hz3id", "Z-index id L3", 11, -5.5, 5.5);
0228
0229   hdigisPerDet1 = fs->make<TH1F>("hdigisPerDet1", "Digis per det l1", 200, -0.5, 199.5);
0230   hdigisPerDet2 = fs->make<TH1F>("hdigisPerDet2", "Digis per det l2", 200, -0.5, 199.5);
0231   hdigisPerDet3 = fs->make<TH1F>("hdigisPerDet3", "Digis per det l3", 200, -0.5, 199.5);
0232   hdigisPerLay1 = fs->make<TH1F>("hdigisPerLay1", "Digis per layer l1", 200, -0.5, 199.5);
0233   hdigisPerLay2 = fs->make<TH1F>("hdigisPerLay2", "Digis per layer l2", 200, -0.5, 199.5);
0234   hdigisPerLay3 = fs->make<TH1F>("hdigisPerLay3", "Digis per layer l3", 200, -0.5, 199.5);
0235   hdetsPerLay1 = fs->make<TH1F>("hdetsPerLay1", "Full dets per layer l1", 161, -0.5, 160.5);
0236   hdetsPerLay3 = fs->make<TH1F>("hdetsPerLay3", "Full dets per layer l3", 353, -0.5, 352.5);
0237   hdetsPerLay2 = fs->make<TH1F>("hdetsPerLay2", "Full dets per layer l2", 257, -0.5, 256.5);
0238
0239   hdigisPerDetF1 = fs->make<TH1F>("hdigisPerDetF1", "Digis per det d1", 200, -0.5, 199.5);
0240   hdigisPerDetF2 = fs->make<TH1F>("hdigisPerDetF2", "Digis per det d2", 200, -0.5, 199.5);
0241   hdigisPerLayF1 = fs->make<TH1F>("hdigisPerLayF1", "Digis per layer d1", 2000, -0.5, 1999.5);
0242   hdigisPerLayF2 = fs->make<TH1F>("hdigisPerLayF2", "Digis per layer d2", 2000, -0.5, 1999.5);
0243   hdetsPerLayF1 = fs->make<TH1F>("hdetsPerLayF1", "Full dets per layer d1", 161, -0.5, 160.5);
0244   hdetsPerLayF2 = fs->make<TH1F>("hdetsPerLayF2", "Full dets per layer d2", 257, -0.5, 256.5);
0245
0246   heloss1 = fs->make<TH1F>("heloss1", "Pix charge l1", 256, 0., 256.);
0247   heloss2 = fs->make<TH1F>("heloss2", "Pix charge l2", 256, 0., 256.);
0248   heloss3 = fs->make<TH1F>("heloss3", "Pix charge l3", 256, 0., 256.);
0249   hneloss1 = fs->make<TH1F>("hneloss1", "Pix noise charge l1", 256, 0., 256.);
0250   hneloss2 = fs->make<TH1F>("hneloss2", "Pix noise charge l2", 256, 0., 256.);
0251   hneloss3 = fs->make<TH1F>("hneloss3", "Pix noise charge l3", 256, 0., 256.);
0252   heloss1bigx = fs->make<TH1F>("heloss1bigx", "L1 big-x pix", 256, 0., 256.);
0253   heloss1bigy = fs->make<TH1F>("heloss1bigy", "L1 big-y pix", 256, 0., 256.);
0254
0255   hcols1 = fs->make<TH1F>("hcols1", "Layer 1 cols", 500, -1.5, 498.5);
0256   hcols2 = fs->make<TH1F>("hcols2", "Layer 2 cols", 500, -1.5, 498.5);
0257   hcols3 = fs->make<TH1F>("hcols3", "Layer 3 cols", 500, -1.5, 498.5);
0258   hcols1big = fs->make<TH1F>("hcols1big", "Layer 1 big cols", 500, -1.5, 498.5);
0259
0260   hrows1 = fs->make<TH1F>("hrows1", "Layer 1 rows", 200, -1.5, 198.5);
0261   hrows2 = fs->make<TH1F>("hrows2", "Layer 2 rows", 200, -1.5, 198.5);
0262   hrows3 = fs->make<TH1F>("hrows3", "layer 3 rows", 200, -1.5, 198.5);
0263   hrows1big = fs->make<TH1F>("hrows1big", "Layer 1 big rows", 200, -1.5, 198.5);
0264
0265   hblade1 = fs->make<TH1F>("hblade1", "blade num, disk1", 24, 0., 24.);
0266   hblade2 = fs->make<TH1F>("hblade2", "blade num, disk2", 24, 0., 24.);
0267
0268   helossF1 = fs->make<TH1F>("helossF1", "Pix charge d1", 100, 0., 300.);
0269   helossF2 = fs->make<TH1F>("helossF2", "Pix charge d2", 100, 0., 300.);
0270   hcolsF1 = fs->make<TH1F>("hcolsF1", "Disk 1 cols", 500, -1.5, 498.5);
0271   hcolsF2 = fs->make<TH1F>("hcolsF2", "Disk 2 cols", 500, -1.5, 498.5);
0272   hrowsF1 = fs->make<TH1F>("hrowsF1", "Disk 1 rows", 200, -1.5, 198.5);
0273   hrowsF2 = fs->make<TH1F>("hrowsF2", "Disk 2 rows", 200, -1.5, 198.5);
0274
0275   hdetr = fs->make<TH1F>("hdetr", "det r", 150, 0., 15.);
0276   hdetz = fs->make<TH1F>("hdetz", "det z", 520, -26., 26.);
0277   hdetrF = fs->make<TH1F>("hdetrF", "det r", 150, 0., 15.);
0278   hdetzF = fs->make<TH1F>("hdetzF", "det z", 700, -70., 70.);
0279
0280   hcolsB = fs->make<TH1F>("hcolsB", "cols per bar det", 450, 0., 450.);
0281   hrowsB = fs->make<TH1F>("hrowsB", "rows per bar det", 200, 0., 200.);
0282   hcolsF = fs->make<TH1F>("hcolsF", "cols per for det", 300, 0., 300.);
0283   hrowsF = fs->make<TH1F>("hrowsF", "rows per for det", 200, 0., 200.);
0284
0285   hsimlinks = fs->make<TH1F>("hsimlinks", " track ids", 200, 0., 200.);
0286   hfract = fs->make<TH1F>("hfract", " track rractions", 100, 0., 1.);
0287   //                                             mod      ladder
0288   hdetMap1 = fs->make<TH2F>("hdetMap1", " ", 9, -4.5, 4.5, 21, -10.5, 10.5);
0289   hdetMap1->SetOption("colz");
0290   hdetMap2 = fs->make<TH2F>("hdetMap2", " ", 9, -4.5, 4.5, 33, -16.5, 16.5);
0291   hdetMap2->SetOption("colz");
0292   hdetMap3 = fs->make<TH2F>("hdetMap3", " ", 9, -4.5, 4.5, 45, -22.5, 22.5);
0293   hdetMap3->SetOption("colz");
0294   hpixMap1 = fs->make<TH2F>("hpixMap1", " ", 416, 0., 416., 160, 0., 160.);
0295   hpixMap1->SetOption("colz");
0296   hpixMap2 = fs->make<TH2F>("hpixMap2", " ", 416, 0., 416., 160, 0., 160.);
0297   hpixMap2->SetOption("colz");
0298   hpixMap3 = fs->make<TH2F>("hpixMap3", " ", 416, 0., 416., 160, 0., 160.);
0299   hpixMap3->SetOption("colz");
0300   hpixMapNoise = fs->make<TH2F>("hpixMapNoise", " ", 416, 0., 416., 160, 0., 160.);
0301   hpixMapNoise->SetOption("colz");
0302
0303   //htest = fs->make<TH2F>("htest"," ",10,0.,10.,20,0.,20.);
0304   //htest2 = fs->make<TH2F>("htest2"," ",10,0.,10.,300,0.,300.);
0305   //htest->SetOption("colz");
0306   //htest2->SetOption("colz");
0307
0308   hevent = fs->make<TH1F>("hevent", "event", 1000, 0, 10000000.);
0309   horbit = fs->make<TH1F>("horbit", "orbit", 100, 0, 100000000.);
0310
0311   hlumi1 = fs->make<TH1F>("hlumi1", "lumi", 2000, 0, 2000.);
0312   hlumi0 = fs->make<TH1F>("hlumi0", "lumi", 2000, 0, 2000.);
0313
0314   hbx6 = fs->make<TH1F>("hbx6", "bx", 4000, 0, 4000.);
0315   hbx5 = fs->make<TH1F>("hbx5", "bx", 4000, 0, 4000.);
0316   hbx4 = fs->make<TH1F>("hbx4", "bx", 4000, 0, 4000.);
0317   hbx3 = fs->make<TH1F>("hbx3", "bx", 4000, 0, 4000.);
0318   hbx2 = fs->make<TH1F>("hbx2", "bx", 4000, 0, 4000.);
0319   hbx1 = fs->make<TH1F>("hbx1", "bx", 4000, 0, 4000.);
0320   hbx0 = fs->make<TH1F>("hbx0", "bx", 4000, 0, 4000.);
0321
0322   hdets = fs->make<TH1F>("hdets", "Dets with hits", 2000, -0.5, 1999.5);
0323   hdigis = fs->make<TH1F>("hdigis", "All Digis", 2000, -0.5, 1999.5);
0324   hdigis1 = fs->make<TH1F>("hdigis1", "All Digis for full events", 2000, -0.5, 1999.5);
0325   hdigis2 = fs->make<TH1F>("hdigis2", "BPix Digis", 2000, -0.5, 1999.5);
0326   hdigis3 = fs->make<TH1F>("hdigis3", "Fpix Digis", 2000, -0.5, 1999.5);
0327   hdigis4 = fs->make<TH1F>("hdigis4", "All Digis - on bunch", 2000, -0.5, 1999.5);
0328   hdigis5 = fs->make<TH1F>("hdigis5", "All Digis - off bunch ", 2000, -0.5, 1999.5);
0329 #endif
0330 }
0331
0332 // ------------ method called to produce the data  ------------
0333 void PixelDigisTest::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0334   //Retrieve tracker topology from geometry
0335   edm::ESHandle<TrackerTopology> tTopo = iSetup.getHandle(trackerTopoToken);
0336
0337   using namespace edm;
0338   if (PRINT)
0339     cout << " Analyze PixelDigisTest " << endl;
0340
0341   //  int run       = iEvent.id().run();
0342   int event = iEvent.id().event();
0343   int lumiBlock = iEvent.luminosityBlock();
0344   int bx = iEvent.bunchCrossing();
0345   int orbit = iEvent.orbitNumber();
0346
0347   hbx0->Fill(float(bx));
0348   hlumi0->Fill(float(lumiBlock));
0349
0350   // eliminate bunches with beam
0351   bool bunch = false;
0352   if (bx == 410 || bx == 460 || bx == 510)
0353     bunch = true;
0354   else if (bx == 560 || bx == 610 || bx == 660)
0355     bunch = true;
0356   else if (bx == 1292 || bx == 1342 || bx == 1392)
0357     bunch = true;
0358   else if (bx == 1454 || bx == 1504 || bx == 1554)
0359     bunch = true;
0360   else if (bx == 2501 || bx == 2601)
0361     bunch = true;
0362   else if (bx == 3080 || bx == 3030 || bx == 3180)
0363     bunch = true;
0364
0365   if (bx >= 1 && bx <= 351) {
0366     if ((bx % 50) == 1)
0367       bunch = true;
0368   } else if (bx >= 892 && bx <= 1245) {
0369     if (((bx - 892) % 50) == 0)
0370       bunch = true;
0371     else if (((bx - 895) % 50) == 0)
0372       bunch = true;
0373   } else if (bx >= 1786 && bx <= 2286) {
0374     if (((bx - 1786) % 50) == 0)
0375       bunch = true;
0376   } else if (bx >= 2671 && bx <= 3021) {
0377     if (((bx - 2671) % 50) == 0)
0378       bunch = true;
0379   }
0380
0381   if (bunch) {
0382     //cout<<" reject "<<bx<<endl;
0383     hbx2->Fill(float(bx));
0384   } else {
0385     if (bx == 892)
0386       cout << " something wrong" << endl;
0387     if (bx == 1245)
0388       cout << " something wrong" << endl;
0389     if (bx == 3021)
0390       cout << " something wrong" << endl;
0391     if (bx == 2286)
0392       cout << " something wrong" << endl;
0393   }
0394
0395   // Get digis
0396   edm::Handle<edm::DetSetVector<PixelDigi>> pixelDigis;
0397   iEvent.getByToken(tPixelDigi, pixelDigis);
0398
0399 #ifdef USE_SIM_LINKS
0400   // Get simlink data
0401   edm::Handle<edm::DetSetVector<PixelDigiSimLink>> pixelSimLinks;
0402   iEvent.getByToken(tPixelDigiSimLink, pixelSimLinks);
0403 #endif
0404
0405   // Get event setup (to get global transformation)
0406   edm::ESHandle<TrackerGeometry> geom = iSetup.getHandle(trackerGeomToken);
0407   const TrackerGeometry &theTracker(*geom);
0408
0409   int numberOfDetUnits = 0;
0410   int totalNumOfDigis = 0;
0411
0412   int numberOfDetUnits1 = 0;
0413   int totalNumOfDigis1 = 0;
0414   int numberOfDetUnits2 = 0;
0415   int totalNumOfDigis2 = 0;
0416   int numberOfDetUnits3 = 0;
0417   int totalNumOfDigis3 = 0;
0418   int numOfDigisPerDet1 = 0;
0419   int numOfDigisPerDet2 = 0;
0420   int numOfDigisPerDet3 = 0;
0421
0422   int numberOfDetUnitsF1 = 0;
0423   int totalNumOfDigisF1 = 0;
0424   int numberOfDetUnitsF2 = 0;
0425   int totalNumOfDigisF2 = 0;
0426   int numOfDigisPerDetF1 = 0;
0427   int numOfDigisPerDetF2 = 0;
0428
0429   // Iterate on detector units
0430   edm::DetSetVector<PixelDigi>::const_iterator DSViter;
0431   for (DSViter = pixelDigis->begin(); DSViter != pixelDigis->end(); DSViter++) {
0432     bool valid = false;
0433     unsigned int detid = DSViter->id;  // = rawid
0434     DetId detId(detid);
0435     //const GeomDetUnit      * geoUnit = geom->idToDetUnit( detId );
0436     //const PixelGeomDetUnit * pixDet  = dynamic_cast<const PixelGeomDetUnit*>(geoUnit);
0437     unsigned int detType = detId.det();     // det type, tracker=1
0438     unsigned int subid = detId.subdetId();  //subdetector type, barrel=1
0439
0440     if (PRINT)
0441       cout << "Det: " << detId.rawId() << " " << detId.null() << " " << detType << " " << subid << endl;
0442
0443 #ifdef HISTOS
0444     hdetunit->Fill(float(detid));
0445     hpixid->Fill(float(detType));
0446     hpixsubid->Fill(float(subid));
0447 #endif
0448
0449     if (detType != 1)
0450       continue;  // look only at tracker
0451     ++numberOfDetUnits;
0452
0453     // Get the geom-detector
0454     const PixelGeomDetUnit *theGeomDet = dynamic_cast<const PixelGeomDetUnit *>(theTracker.idToDet(detId));
0455     double detZ = theGeomDet->surface().position().z();
0456     double detR = theGeomDet->surface().position().perp();
0457     //const BoundPlane plane = theGeomDet->surface(); // does not work
0458
0459     //     int cols = theGeomDet->specificTopology().ncolumns();
0460     //     int rows = theGeomDet->specificTopology().nrows();
0461     //     float pitchX = theGeomDet->specificTopology().pitch().first;
0462     //     float pitchY = theGeomDet->specificTopology().pitch().second;
0463
0464     const PixelTopology &topology = theGeomDet->specificTopology();
0465     int cols = topology.ncolumns();
0466     int rows = topology.nrows();
0467     float pitchX = topology.pitch().first;
0468     float pitchY = topology.pitch().second;
0469
0470     unsigned int layerC = 0;
0471     unsigned int ladderC = 0;
0472     unsigned int zindex = 0;
0473
0474     // Shell { mO = 1, mI = 2 , pO =3 , pI =4 };
0475     int shell = 0;      // shell id
0476     int sector = 0;     // 1-8
0477     int ladder = 0;     // 1-22
0478     int layer = 0;      // 1-3
0479     int module = 0;     // 1-4
0480     bool half = false;  //
0481
0482     unsigned int disk = 0;   //1,2,3
0483     unsigned int blade = 0;  //1-24
0484     unsigned int side = 0;   //size=1 for -z, 2 for +z
0485     unsigned int panel = 0;  //panel=1
0486
0487     // Subdet it, pix barrel=1, forward=2
0488     if (subid == 2) {  // forward
0489
0490 #ifdef HISTOS
0491       hdetrF->Fill(detR);
0492       hdetzF->Fill(detZ);
0493       hcolsF->Fill(float(cols));
0494       hrowsF->Fill(float(rows));
0495 #endif
0496
0497       disk = tTopo->pxfDisk(detid);      //1,2,3
0498       blade = tTopo->pxfBlade(detid);    //1-24
0499       zindex = tTopo->pxfModule(detid);  //
0500       side = tTopo->pxfSide(detid);      //size=1 for -z, 2 for +z
0501       panel = tTopo->pxfPanel(detid);    //panel=1
0502
0503       if (PRINT) {
0504         cout << "Forward det " << subid << ", disk " << disk << ", blade " << blade << ", module " << zindex
0505              << ", side " << side << ", panel " << panel << endl;
0506         cout << " col/row, pitch " << cols << " " << rows << " " << pitchX << " " << pitchY << endl;
0507       }
0508
0509     } else if (subid == 1) {  // Barrel
0510
0511       // Barell layer = 1,2,3
0512       layerC = tTopo->pxbLayer(detid);
0513       // Barrel ladder id 1-20,32,44.
0514       ladderC = tTopo->pxbLadder(detid);
0515       // Barrel Z-index=1,8
0516       zindex = tTopo->pxbModule(detid);
0517       // Convert to online
0518       PixelBarrelName pbn(detid);
0519       // Shell { mO = 1, mI = 2 , pO =3 , pI =4 };
0520       PixelBarrelName::Shell sh = pbn.shell();  //enum
0521       sector = pbn.sectorName();
0522       ladder = pbn.ladderName();
0523       layer = pbn.layerName();
0524       module = pbn.moduleName();
0525       half = pbn.isHalfModule();
0526       shell = int(sh);
0527       // change the module sign for z<0
0528       if (shell == 1 || shell == 2)
0529         module = -module;
0530       // change ladeer sign for Outer )x<0)
0531       if (shell == 1 || shell == 3)
0532         ladder = -ladder;
0533
0534       if (PRINT) {
0535         cout << " Barrel layer, ladder, module " << layerC << " " << ladderC << " " << zindex << " " << sh << "("
0536              << shell << ") " << sector << " " << layer << " " << ladder << " " << module << " " << half << endl;
0537         //cout<<" Barrel det, thick "<<detThick<<" "
0538         //  <<" layer, ladder, module "
0539         //  <<layer<<" "<<ladder<<" "<<zindex<<endl;
0540         //cout<<" col/row, pitch "<<cols<<" "<<rows<<" "
0541         //  <<pitchX<<" "<<pitchY<<endl;
0542       }
0543
0544     }  // end fb-bar
0545
0546 #ifdef USE_SIM_LINKS
0547     // Look at simlink information (simulated data only)
0548
0549     int numberOfSimLinks = 0;
0550     edm::DetSetVector<PixelDigiSimLink>::const_iterator isearch = pixelSimLinks->find(detid);
0551
0552     if (isearch != pixelSimLinks->end()) {  //if it is not empty
0553       edm::DetSet<PixelDigiSimLink> link_detset = (*pixelSimLinks)[detid];
0554       edm::DetSet<PixelDigiSimLink>::const_iterator it;
0555       // Loop over DigisSimLink in this det unit
0556       for (it = link_detset.data.begin(); it != link_detset.data.end(); it++) {
0557         numberOfSimLinks++;
0558         // these methods should be declared const, fixed by M.P.
0559         // wait for next releasse and then uncomment
0560         unsigned int chan = it->channel();
0561         unsigned int simTrack = it->SimTrackId();
0562         float frac = it->fraction();
0563 #ifdef HISTOS
0564         hsimlinks->Fill(float(simTrack));
0565         hfract->Fill(frac);
0566 #endif
0567
0568         if (PRINT)
0569           cout << " Sim link " << numberOfSimLinks << " " << chan << " " << simTrack << " " << frac << endl;
0570       }  // end simlink det loop
0571
0572     }  // end simlink if
0573
0574 #endif  // USE_SIM_LINKS
0575
0576     unsigned int numberOfDigis = 0;
0577
0578     // Look at digis now
0579     edm::DetSet<PixelDigi>::const_iterator di;
0580     for (di = DSViter->data.begin(); di != DSViter->data.end(); di++) {
0581       //for(di = begin; di != end; di++) {
0582
0583       numberOfDigis++;
0584       totalNumOfDigis++;
0585       int adc = di->adc();     // charge, modifued to unsiged short
0586       int col = di->column();  // column
0587       int row = di->row();     // row
0588       //int tof = di->time();    // tof always 0, method deleted
0589
0590       // channel index needed to look for the simlink to simtracks
0591       int channel = PixelChannelIdentifier::pixelToChannel(row, col);
0592       if (PRINT)
0593         cout << numberOfDigis << " Col: " << col << " Row: " << row << " ADC: " << adc << " channel = " << channel
0594              << endl;
0595
0596       if (col > 415)
0597         cout << " Error: column index too large " << col << " Barrel layer, ladder, module " << layer << " " << ladder
0598              << " " << zindex << endl;
0599       if (row > 159)
0600         cout << " Error: row index too large " << row << endl;
0601
0602 #ifdef HISTOS
0603       if (layer == 1) {
0604         bool noise = false;  //(ladder==6) || (module==-2) || (col==364) || (row==1);
0605         if (!noise) {
0606           valid = valid || true;
0607           heloss1->Fill(float(adc));
0608           hcols1->Fill(float(col));
0609           hrows1->Fill(float(row));
0610           hpixMap1->Fill(float(col), float(row));
0611           totalNumOfDigis1++;
0612           //htest2->Fill(float(module),float(adc));
0613           numOfDigisPerDet1++;
0614
0615           //old        if(RectangularPixelTopology::isItBigPixelInX(row)) {
0616           //new    if(topology.isItBigPixelInX(row)) {
0617           //         //cout<<" big in X "<<row<<endl;
0618           //         heloss1bigx->Fill(float(adc));
0619           //         hrows1big->Fill(float(row));
0620           //       }
0621           //old    if(RectangularPixelTopology::isItBigPixelInY(col)) {
0622           //new    if(topology.isItBigPixelInY(col)) {
0623           //         //cout<<" big in Y "<<col<<endl;
0624           //         heloss1bigy->Fill(float(adc));
0625           //         hcols1big->Fill(float(col));
0626           //       }
0627
0628         }  // noise
0629       } else if (layer == 2) {
0630         // look for the noisy pixel
0631         bool noise = false;  // (ladder==6) && (module==-2) && (col==364) && (row==1);
0632         if (noise) {
0633           //cout<<" noise pixel "<<layer<<" "<<sector<<" "<<shell<<endl;
0634           hpixMapNoise->Fill(float(col), float(row));
0635           hneloss2->Fill(float(adc));
0636         } else {
0637           valid = valid || true;
0638           heloss2->Fill(float(adc));
0639           hcols2->Fill(float(col));
0640           hrows2->Fill(float(row));
0641           hpixMap2->Fill(float(col), float(row));
0642           totalNumOfDigis2++;
0643           numOfDigisPerDet2++;
0644         }  // noise
0645       } else if (layer == 3) {
0646         bool noise = false;  //(ladder==6) || (module==-2) || (col==364) || (row==1);
0647         if (!noise) {
0648           valid = valid || true;
0649           heloss3->Fill(float(adc));
0650           hcols3->Fill(float(col));
0651           hrows3->Fill(float(row));
0652           hpixMap3->Fill(float(col), float(row));
0653           totalNumOfDigis3++;
0654           numOfDigisPerDet3++;
0655         }  // noise
0656       } else if (disk == 1) {
0657         bool noise = false;  //(ladder==6) || (module==-2) || (col==364) || (row==1);
0658         if (!noise) {
0659           valid = valid || true;
0660           helossF1->Fill(float(adc));
0661           hcolsF1->Fill(float(col));
0662           hrowsF1->Fill(float(row));
0663           totalNumOfDigisF1++;
0664           numOfDigisPerDetF1++;
0665         }  // noise
0666
0667       } else if (disk == 2) {
0668         bool noise = false;  //(ladder==6) || (module==-2) || (col==364) || (row==1);
0669         if (!noise) {
0670           valid = valid || true;
0671           helossF2->Fill(float(adc));
0672           hcolsF2->Fill(float(col));
0673           hrowsF2->Fill(float(row));
0674           totalNumOfDigisF2++;
0675           numOfDigisPerDetF2++;
0676         }  // noise
0677       }    // end if layer
0678 #endif
0679
0680     }  // end for digis in detunit
0681        //if(PRINT)
0682        //cout<<" for det "<<detid<<" digis = "<<numberOfDigis<<endl;
0683
0684 #ifdef HISTOS
0685     // Some histos
0686     if (valid) {         // to igore noise pixels
0687       if (subid == 2) {  // forward
0688
0689         hdetrF->Fill(detR);
0690         hdetzF->Fill(detZ);
0691         hcolsF->Fill(float(cols));
0692         hrowsF->Fill(float(rows));
0693
0694         if (disk == 1) {
0695           hblade1->Fill(float(blade));
0696           ++numberOfDetUnitsF1;
0697           hdigisPerDetF1->Fill(float(numOfDigisPerDetF1));
0698           numOfDigisPerDetF1 = 0;
0699
0700         } else if (disk == 2) {
0701           hblade2->Fill(float(blade));
0702           ++numberOfDetUnitsF2;
0703           hdigisPerDetF2->Fill(float(numOfDigisPerDetF2));
0704           numOfDigisPerDetF2 = 0;
0705         }  // if disk
0706
0707       } else if (subid == 1) {  // barrel
0708
0709         hdetr->Fill(detR);
0710         hdetz->Fill(detZ);
0711         hcolsB->Fill(float(cols));
0712         hrowsB->Fill(float(rows));
0713
0714         hlayerid->Fill(float(layer));
0715         hshellid->Fill(float(shell));
0716         hsectorid->Fill(float(sector));
0717
0718         if (layer == 1) {
0719           hladder1id->Fill(float(ladder));
0720           hz1id->Fill(float(module));
0721           hdetMap1->Fill(float(module), float(ladder));
0722           ++numberOfDetUnits1;
0723           hdigisPerDet1->Fill(float(numOfDigisPerDet1));
0724           numOfDigisPerDet1 = 0;
0725
0726         } else if (layer == 2) {
0727           hladder2id->Fill(float(ladder));
0728           hz2id->Fill(float(module));
0729           hdetMap2->Fill(float(module), float(ladder));
0730           ++numberOfDetUnits2;
0731           hdigisPerDet2->Fill(float(numOfDigisPerDet2));
0732           numOfDigisPerDet2 = 0;
0733
0734         } else if (layer == 3) {
0735           hladder3id->Fill(float(ladder));
0736           hz3id->Fill(float(module));
0737           hdetMap3->Fill(float(module), float(ladder));
0738           ++numberOfDetUnits3;
0739           hdigisPerDet3->Fill(float(numOfDigisPerDet3));
0740           numOfDigisPerDet3 = 0;
0741
0742         }  // layer
0743       }    // if bpix
0744     }      // if valid
0745 #endif
0746
0747   }  // end for det-units
0748
0749   if (PRINT)
0750     cout << " Number of full det-units = " << numberOfDetUnits << " total digis = " << totalNumOfDigis << endl;
0751   hdets->Fill(float(numberOfDetUnits));
0752   hdigis->Fill(float(totalNumOfDigis));
0753
0754   if (numberOfDetUnits > 0) {
0755     hevent->Fill(float(event));
0756     hlumi1->Fill(float(lumiBlock));
0757     hbx1->Fill(float(bx));
0758     horbit->Fill(float(orbit));
0759
0760     hdigis1->Fill(float(totalNumOfDigis));
0761     float tmp = float(totalNumOfDigis1) + float(totalNumOfDigis2) + float(totalNumOfDigis3);
0762     hdigis2->Fill(tmp);
0763     tmp = float(totalNumOfDigisF1) + float(totalNumOfDigisF2);
0764     hdigis3->Fill(tmp);
0765
0766     if (bunch)
0767       hdigis4->Fill(float(totalNumOfDigis));
0768     else
0769       hdigis5->Fill(float(totalNumOfDigis));
0770
0771     if (numberOfDetUnits > 20)
0772       hbx3->Fill(float(bx));
0773
0774     if (totalNumOfDigis > 100)
0775       hbx4->Fill(float(bx));
0776     else if (totalNumOfDigis > 4)
0777       hbx5->Fill(float(bx));
0778     else
0779       hbx6->Fill(float(bx));
0780   }
0781
0782 #ifdef HISTOS
0783   hdigisPerLay1->Fill(float(totalNumOfDigis1));
0784   hdigisPerLay2->Fill(float(totalNumOfDigis2));
0785   hdigisPerLay3->Fill(float(totalNumOfDigis3));
0786   if (totalNumOfDigis1 > 0)
0787     hdetsPerLay1->Fill(float(numberOfDetUnits1));
0788   if (totalNumOfDigis2 > 0)
0789     hdetsPerLay2->Fill(float(numberOfDetUnits2));
0790   if (totalNumOfDigis3 > 0)
0791     hdetsPerLay3->Fill(float(numberOfDetUnits3));
0792
0793   hdigisPerLayF1->Fill(float(totalNumOfDigisF1));
0794   hdigisPerLayF2->Fill(float(totalNumOfDigisF2));
0795   hdetsPerLayF1->Fill(float(numberOfDetUnitsF1));
0796   hdetsPerLayF2->Fill(float(numberOfDetUnitsF2));
0797 #endif
0798 }
0799 // ------------ method called to at the end of the job  ------------
0800 void PixelDigisTest::endJob() {
0801   cout << " End PixelDigisTest " << endl;
0802   //hFile->Write();
0803   //hFile->Close();
0804 }
0805
0806 //define this as a plug-in
0807 DEFINE_FWK_MODULE(PixelDigisTest);