Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2023-10-25 10:06:44 

0001 // -*- C++ -*-
0002 //
0003 // Package:    HGCalValidation/HGCalMTRecoStudy
0004 // Class:      HGCalMTRecoStudy
0005 //
0006 /**\class HGCalMTRecoStudy HGCalMTRecoStudy.cc Validation/HGCalValidation/test/HGCalMTRecoStudy.cc
0007 // Derived from : HGCalRecHitStudy.cc
0008 
0009  Description: [one line class summary]
0010 
0011  Implementation:
0012      [Notes on implementation]
0013 */
0014 //
0015 // Original Author:  Indranil Das
0016 //         Created:  Thu, 06 Oct 2022 06:18:11 GMT
0017 //
0018 //
0019 
0020 // system include files
0021 #include <cmath>
0022 #include <iostream>
0023 #include <fstream>
0024 #include <map>
0025 #include <vector>
0026 #include <string>
0027 
0028 // user include files
0029 #include "FWCore/Framework/interface/Frameworkfwd.h"
0030 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0031 #include "FWCore/Framework/interface/Event.h"
0032 #include "FWCore/Framework/interface/EventSetup.h"
0033 #include "FWCore/Framework/interface/ESHandle.h"
0034 #include "FWCore/Utilities/interface/InputTag.h"
0035 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0036 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0037 #include "FWCore/ServiceRegistry/interface/Service.h"
0038 
0039 #include "FWCore/Framework/interface/MakerMacros.h"
0040 #include "FWCore/Framework/interface/ESTransientHandle.h"
0041 #include "FWCore/Framework/interface/ModuleFactory.h"
0042 #include "FWCore/Framework/interface/ESProducer.h"
0043 
0044 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0045 
0046 #include "DataFormats/DetId/interface/DetId.h"
0047 #include "DataFormats/ForwardDetId/interface/HGCSiliconDetId.h"
0048 #include "DataFormats/ForwardDetId/interface/HGCScintillatorDetId.h"
0049 #include "DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"
0050 
0051 #include "RecoLocalCalo/HGCalRecAlgos/interface/RecHitTools.h"
0052 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0053 #include "Geometry/HGCalGeometry/interface/HGCalGeometry.h"
0054 #include "Geometry/HGCalCommonData/interface/HGCalDDDConstants.h"
0055 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0056 
0057 #include "TH1D.h"
0058 #include "TH2D.h"
0059 #include "TGraph.h"
0060 
0061 //
0062 // class declaration
0063 //
0064 
0065 class HGCalMTRecoStudy : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0066 public:
0067   //Implemented following Validation/HGCalValidation/test/HGCalRecHitStudy.cc
0068   explicit HGCalMTRecoStudy(const edm::ParameterSet &);
0069   ~HGCalMTRecoStudy() override = default;
0070 
0071   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0072 
0073 private:
0074   void analyze(const edm::Event &, const edm::EventSetup &) override;
0075 
0076   // ----------member data ---------------------------
0077   edm::ESGetToken<CaloGeometry, CaloGeometryRecord> caloGeomToken_;
0078   hgcal::RecHitTools rhtools_;
0079   const std::string nameDetector_;
0080   const edm::EDGetTokenT<HGCRecHitCollection> recHitSource_;
0081   const edm::ESGetToken<HGCalDDDConstants, IdealGeometryRecord> tok_hgcaldd_;
0082   const edm::ESGetToken<HGCalGeometry, IdealGeometryRecord> tok_hgcGeom_;
0083 
0084   const std::string layers_;
0085   std::vector<Int_t> layerList;
0086 
0087   // For rechittool z positions. The 0 and 1 are for -ve and +ve, respectively.
0088   TH1D *hEF;
0089   TH1D *hECN;
0090   TH1D *hECK;
0091   TH1D *hESc;
0092 
0093   std::vector<TH1D *> hELossLayer0;
0094   std::vector<TH1D *> hELossLayer1;
0095   std::vector<TH1D *> hNCellsLayer0;
0096   std::vector<TH1D *> hNCellsLayer1;
0097 
0098   std::vector<TH2D *> hXYhitsF0;
0099   std::vector<TH2D *> hXYhitsCN0;
0100   std::vector<TH2D *> hXYhitsCK0;
0101   std::vector<TH2D *> hXYhitsB0;
0102 
0103   std::vector<TH2D *> hXYhitsF1;
0104   std::vector<TH2D *> hXYhitsCN1;
0105   std::vector<TH2D *> hXYhitsCK1;
0106   std::vector<TH2D *> hXYhitsB1;
0107 
0108   std::vector<TH2D *> hEPhitsF0;
0109   std::vector<TH2D *> hEPhitsCN0;
0110   std::vector<TH2D *> hEPhitsCK0;
0111   std::vector<TH2D *> hEPhitsB0;
0112 
0113   std::vector<TH2D *> hEPhitsF1;
0114   std::vector<TH2D *> hEPhitsCN1;
0115   std::vector<TH2D *> hEPhitsCK1;
0116   std::vector<TH2D *> hEPhitsB1;
0117 
0118   std::vector<TH2D *> hXYFailhitsF0;
0119   std::vector<TH2D *> hXYFailhitsCN0;
0120   std::vector<TH2D *> hXYFailhitsCK0;
0121   std::vector<TH2D *> hXYFailhitsB0;
0122 
0123   std::vector<TH2D *> hXYFailhitsF1;
0124   std::vector<TH2D *> hXYFailhitsCN1;
0125   std::vector<TH2D *> hXYFailhitsCK1;
0126   std::vector<TH2D *> hXYFailhitsB1;
0127 
0128   std::vector<TH2D *> hEPFailhitsF0;
0129   std::vector<TH2D *> hEPFailhitsCN0;
0130   std::vector<TH2D *> hEPFailhitsCK0;
0131   std::vector<TH2D *> hEPFailhitsB0;
0132 
0133   std::vector<TH2D *> hEPFailhitsF1;
0134   std::vector<TH2D *> hEPFailhitsCN1;
0135   std::vector<TH2D *> hEPFailhitsCK1;
0136   std::vector<TH2D *> hEPFailhitsB1;
0137 
0138   std::vector<TH1D *> hELossLayerF0;
0139   std::vector<TH1D *> hELossLayerCN0;
0140   std::vector<TH1D *> hELossLayerCK0;
0141   std::vector<TH1D *> hELossLayerB0;
0142 
0143   std::vector<TH1D *> hELossLayerF1;
0144   std::vector<TH1D *> hELossLayerCN1;
0145   std::vector<TH1D *> hELossLayerCK1;
0146   std::vector<TH1D *> hELossLayerB1;
0147 
0148   // For rechittool z positions. The 0 and 1 are for -ve and +ve, respectively.
0149   std::vector<TGraph *> grXYhitsF0;
0150   std::vector<TGraph *> grXYhitsCN0;
0151   std::vector<TGraph *> grXYhitsCK0;
0152   std::vector<TGraph *> grXYhitsAR0;
0153   std::vector<TGraph *> grXYhitsB0;
0154   int ixyF0[50], ixyCN0[50], ixyCK0[50], ixyAR0[50], ixyB0[50];
0155 
0156   std::vector<TGraph *> grXYhitsF1;
0157   std::vector<TGraph *> grXYhitsCN1;
0158   std::vector<TGraph *> grXYhitsCK1;
0159   std::vector<TGraph *> grXYhitsAR1;
0160   std::vector<TGraph *> grXYhitsB1;
0161   int ixyF1[50], ixyCN1[50], ixyCK1[50], ixyAR1[50], ixyB1[50];
0162   /////////////////////////////////
0163 
0164   // For rechittool z positions. The 0 and 1 are for -ve and +ve, respectively.
0165   std::vector<TGraph *> grEtaPhihitsF0;
0166   std::vector<TGraph *> grEtaPhihitsCN0;
0167   std::vector<TGraph *> grEtaPhihitsCK0;
0168   std::vector<TGraph *> grEtaPhihitsB0;
0169   int iepF0[50], iepCN0[50], iepCK0[50], iepB0[50];
0170 
0171   std::vector<TGraph *> grEtaPhihitsF1;
0172   std::vector<TGraph *> grEtaPhihitsCN1;
0173   std::vector<TGraph *> grEtaPhihitsCK1;
0174   std::vector<TGraph *> grEtaPhihitsB1;
0175   int iepF1[50], iepCN1[50], iepCK1[50], iepB1[50];
0176 };
0177 
0178 //
0179 // constructors and destructor
0180 //
0181 HGCalMTRecoStudy::HGCalMTRecoStudy(const edm::ParameterSet &iConfig)
0182     : nameDetector_(iConfig.getParameter<std::string>("detectorName")),
0183       recHitSource_(consumes<HGCRecHitCollection>(iConfig.getParameter<edm::InputTag>("source"))),
0184       tok_hgcaldd_(esConsumes<HGCalDDDConstants, IdealGeometryRecord, edm::Transition::BeginRun>(
0185           edm::ESInputTag{"", nameDetector_})),
0186       tok_hgcGeom_(esConsumes<HGCalGeometry, IdealGeometryRecord>(edm::ESInputTag{"", nameDetector_})),
0187       layers_(iConfig.getParameter<std::string>("layerList")) {
0188   layerList.clear();
0189 
0190   if (layers_.find("-") != std::string::npos) {
0191     std::vector<std::string> tokens;
0192     std::stringstream check1(layers_);
0193     std::string intermediate;
0194     while (getline(check1, intermediate, '-'))
0195       tokens.push_back(intermediate);
0196     int minLayer = (stoi(tokens[0]) < 1) ? 1 : stoi(tokens[0]);
0197     int maxLayer = (stoi(tokens[1]) > 47) ? 47 : stoi(tokens[1]);
0198     for (int i = minLayer; i <= maxLayer; i++) {
0199       layerList.push_back(i);
0200       //std::cout << tokens[i] << '\n';
0201     }
0202     tokens.clear();
0203 
0204   } else if (layers_.find(",") != std::string::npos) {
0205     std::vector<std::string> tokens;
0206     std::stringstream check1(layers_);
0207     std::string intermediate;
0208     while (getline(check1, intermediate, ','))
0209       tokens.push_back(intermediate);
0210     for (unsigned int i = 0; i < tokens.size(); i++) {
0211       if (stoi(tokens[i]) >= 1 and stoi(tokens[i]) <= 47)
0212         layerList.push_back(stoi(tokens[i]));
0213       //std::cout << tokens[i] << '\n';
0214     }
0215     tokens.clear();
0216 
0217   } else {
0218     if (stoi(layers_) >= 1 and stoi(layers_) <= 47)
0219       layerList.push_back(stoi(layers_));
0220   }
0221 
0222   usesResource(TFileService::kSharedResource);
0223   edm::Service<TFileService> fs;
0224 
0225   hEF = fs->make<TH1D>("hEF", "hEF", 1000, 0., 50.);
0226   hECN = fs->make<TH1D>("hECN", "hECN", 1000, 0., 50.);
0227   hECK = fs->make<TH1D>("hECK", "hECK", 1000, 0., 50.);
0228   hESc = fs->make<TH1D>("hESc", "hESc", 1000, 0., 50.);
0229 
0230   for (unsigned int i = 0; i < 47; i++) {
0231     hELossLayer0.emplace_back(
0232         fs->make<TH1D>(Form("hELossLayer0_%02d", i + 1), Form("Rec:ELoss0 for layer %d", i + 1), 500, 0., 5000.));
0233     hELossLayer1.emplace_back(
0234         fs->make<TH1D>(Form("hELossLayer1_%02d", i + 1), Form("Rec:ELoss1 for layer %d", i + 1), 500, 0., 5000.));
0235     hNCellsLayer0.emplace_back(
0236         fs->make<TH1D>(Form("hNCellsLayer0_%02d", i + 1), Form("Rec:NCells0 for layer %d", i + 1), 200, -0.5, 199.5));
0237     hNCellsLayer1.emplace_back(
0238         fs->make<TH1D>(Form("hNCellsLayer1_%02d", i + 1), Form("Rec:NCells1 for layer %d", i + 1), 200, -0.5, 199.5));
0239   }
0240 
0241   for (unsigned int i = 0; i < layerList.size(); i++) {
0242     hXYhitsF0.emplace_back(fs->make<TH2D>(Form("hXYhitsF0_layer_%02d", layerList[i]),
0243                                           Form("Rec:HitsF0 in XY for layer %d", layerList[i]),
0244                                           600,
0245                                           -300.,
0246                                           300.,
0247                                           600,
0248                                           -300.,
0249                                           300.));
0250     hXYhitsCN0.emplace_back(fs->make<TH2D>(Form("hXYhitsCN0_layer_%02d", layerList[i]),
0251                                            Form("Rec:HitsCN0 in XY for layer %d", layerList[i]),
0252                                            600,
0253                                            -300.,
0254                                            300.,
0255                                            600,
0256                                            -300.,
0257                                            300.));
0258     hXYhitsCK0.emplace_back(fs->make<TH2D>(Form("hXYhitsCK0_layer_%02d", layerList[i]),
0259                                            Form("Rec:HitsCK0 in XY for layer %d", layerList[i]),
0260                                            600,
0261                                            -300.,
0262                                            300.,
0263                                            600,
0264                                            -300.,
0265                                            300.));
0266     hXYhitsB0.emplace_back(fs->make<TH2D>(Form("hXYhitsB0_layer_%02d", layerList[i]),
0267                                           Form("Rec:HitsB0 in XY for layer %d", layerList[i]),
0268                                           600,
0269                                           -300.,
0270                                           300.,
0271                                           600,
0272                                           -300.,
0273                                           300.));
0274     hXYhitsF1.emplace_back(fs->make<TH2D>(Form("hXYhitsF1_layer_%02d", layerList[i]),
0275                                           Form("Rec:HitsF1 in XY for layer %d", layerList[i]),
0276                                           600,
0277                                           -300.,
0278                                           300.,
0279                                           600,
0280                                           -300.,
0281                                           300.));
0282     hXYhitsCN1.emplace_back(fs->make<TH2D>(Form("hXYhitsCN1_layer_%02d", layerList[i]),
0283                                            Form("Rec:HitsCN1 in XY for layer %d", layerList[i]),
0284                                            600,
0285                                            -300.,
0286                                            300.,
0287                                            600,
0288                                            -300.,
0289                                            300.));
0290     hXYhitsCK1.emplace_back(fs->make<TH2D>(Form("hXYhitsCK1_layer_%02d", layerList[i]),
0291                                            Form("Rec:HitsCK1 in XY for layer %d", layerList[i]),
0292                                            600,
0293                                            -300.,
0294                                            300.,
0295                                            600,
0296                                            -300.,
0297                                            300.));
0298     hXYhitsB1.emplace_back(fs->make<TH2D>(Form("hXYhitsB1_layer_%02d", layerList[i]),
0299                                           Form("Rec:HitsB1 in XY for layer %d", layerList[i]),
0300                                           600,
0301                                           -300.,
0302                                           300.,
0303                                           600,
0304                                           -300.,
0305                                           300.));
0306 
0307     hEPhitsF0.emplace_back(fs->make<TH2D>(Form("hEPhitsF0_layer_%02d", layerList[i]),
0308                                           Form("Rec:HitsF0 in EP for layer %d", layerList[i]),
0309                                           640,
0310                                           -3.2,
0311                                           3.2,
0312                                           640,
0313                                           -3.2,
0314                                           3.2));
0315     hEPhitsCN0.emplace_back(fs->make<TH2D>(Form("hEPhitsCN0_layer_%02d", layerList[i]),
0316                                            Form("Rec:HitsCN0 in EP for layer %d", layerList[i]),
0317                                            640,
0318                                            -3.2,
0319                                            3.2,
0320                                            640,
0321                                            -3.2,
0322                                            3.2));
0323     hEPhitsCK0.emplace_back(fs->make<TH2D>(Form("hEPhitsCK0_layer_%02d", layerList[i]),
0324                                            Form("Rec:HitsCK0 in EP for layer %d", layerList[i]),
0325                                            640,
0326                                            -3.2,
0327                                            3.2,
0328                                            640,
0329                                            -3.2,
0330                                            3.2));
0331     hEPhitsB0.emplace_back(fs->make<TH2D>(Form("hEPhitsB0_layer_%02d", layerList[i]),
0332                                           Form("Rec:HitsB0 in EP for layer %d", layerList[i]),
0333                                           640,
0334                                           -3.2,
0335                                           3.2,
0336                                           640,
0337                                           -3.2,
0338                                           3.2));
0339     hEPhitsF1.emplace_back(fs->make<TH2D>(Form("hEPhitsF1_layer_%02d", layerList[i]),
0340                                           Form("Rec:HitsF1 in EP for layer %d", layerList[i]),
0341                                           640,
0342                                           -3.2,
0343                                           3.2,
0344                                           640,
0345                                           -3.2,
0346                                           3.2));
0347     hEPhitsCN1.emplace_back(fs->make<TH2D>(Form("hEPhitsCN1_layer_%02d", layerList[i]),
0348                                            Form("Rec:HitsCN1 in EP for layer %d", layerList[i]),
0349                                            640,
0350                                            -3.2,
0351                                            3.2,
0352                                            640,
0353                                            -3.2,
0354                                            3.2));
0355     hEPhitsCK1.emplace_back(fs->make<TH2D>(Form("hEPhitsCK1_layer_%02d", layerList[i]),
0356                                            Form("Rec:HitsCK1 in EP for layer %d", layerList[i]),
0357                                            640,
0358                                            -3.2,
0359                                            3.2,
0360                                            640,
0361                                            -3.2,
0362                                            3.2));
0363     hEPhitsB1.emplace_back(fs->make<TH2D>(Form("hEPhitsB1_layer_%02d", layerList[i]),
0364                                           Form("Rec:HitsB1 in EP for layer %d", layerList[i]),
0365                                           640,
0366                                           -3.2,
0367                                           3.2,
0368                                           640,
0369                                           -3.2,
0370                                           3.2));
0371 
0372     hXYFailhitsF0.emplace_back(fs->make<TH2D>(Form("hXYFailhitsF0_layer_%02d", layerList[i]),
0373                                               Form("Rec:FailhitsF0 in XY for layer %d", layerList[i]),
0374                                               600,
0375                                               -300.,
0376                                               300.,
0377                                               600,
0378                                               -300.,
0379                                               300.));
0380     hXYFailhitsCN0.emplace_back(fs->make<TH2D>(Form("hXYFailhitsCN0_layer_%02d", layerList[i]),
0381                                                Form("Rec:FailhitsCN0 in XY for layer %d", layerList[i]),
0382                                                600,
0383                                                -300.,
0384                                                300.,
0385                                                600,
0386                                                -300.,
0387                                                300.));
0388     hXYFailhitsCK0.emplace_back(fs->make<TH2D>(Form("hXYFailhitsCK0_layer_%02d", layerList[i]),
0389                                                Form("Rec:FailhitsCK0 in XY for layer %d", layerList[i]),
0390                                                600,
0391                                                -300.,
0392                                                300.,
0393                                                600,
0394                                                -300.,
0395                                                300.));
0396     hXYFailhitsB0.emplace_back(fs->make<TH2D>(Form("hXYFailhitsB0_layer_%02d", layerList[i]),
0397                                               Form("Rec:FailhitsB0 in XY for layer %d", layerList[i]),
0398                                               600,
0399                                               -300.,
0400                                               300.,
0401                                               600,
0402                                               -300.,
0403                                               300.));
0404     hXYFailhitsF1.emplace_back(fs->make<TH2D>(Form("hXYFailhitsF1_layer_%02d", layerList[i]),
0405                                               Form("Rec:FailhitsF1 in XY for layer %d", layerList[i]),
0406                                               600,
0407                                               -300.,
0408                                               300.,
0409                                               600,
0410                                               -300.,
0411                                               300.));
0412     hXYFailhitsCN1.emplace_back(fs->make<TH2D>(Form("hXYFailhitsCN1_layer_%02d", layerList[i]),
0413                                                Form("Rec:FailhitsCN1 in XY for layer %d", layerList[i]),
0414                                                600,
0415                                                -300.,
0416                                                300.,
0417                                                600,
0418                                                -300.,
0419                                                300.));
0420     hXYFailhitsCK1.emplace_back(fs->make<TH2D>(Form("hXYFailhitsCK1_layer_%02d", layerList[i]),
0421                                                Form("Rec:FailhitsCK1 in XY for layer %d", layerList[i]),
0422                                                600,
0423                                                -300.,
0424                                                300.,
0425                                                600,
0426                                                -300.,
0427                                                300.));
0428     hXYFailhitsB1.emplace_back(fs->make<TH2D>(Form("hXYFailhitsB1_layer_%02d", layerList[i]),
0429                                               Form("Rec:FailhitsB1 in XY for layer %d", layerList[i]),
0430                                               600,
0431                                               -300.,
0432                                               300.,
0433                                               600,
0434                                               -300.,
0435                                               300.));
0436 
0437     hEPFailhitsF0.emplace_back(fs->make<TH2D>(Form("hEPFailhitsF0_layer_%02d", layerList[i]),
0438                                               Form("Rec:FailhitsF0 in EP for layer %d", layerList[i]),
0439                                               640,
0440                                               -3.2,
0441                                               3.2,
0442                                               640,
0443                                               -3.2,
0444                                               3.2));
0445     hEPFailhitsCN0.emplace_back(fs->make<TH2D>(Form("hEPFailhitsCN0_layer_%02d", layerList[i]),
0446                                                Form("Rec:FailhitsCN0 in EP for layer %d", layerList[i]),
0447                                                640,
0448                                                -3.2,
0449                                                3.2,
0450                                                640,
0451                                                -3.2,
0452                                                3.2));
0453     hEPFailhitsCK0.emplace_back(fs->make<TH2D>(Form("hEPFailhitsCK0_layer_%02d", layerList[i]),
0454                                                Form("Rec:FailhitsCK0 in EP for layer %d", layerList[i]),
0455                                                640,
0456                                                -3.2,
0457                                                3.2,
0458                                                640,
0459                                                -3.2,
0460                                                3.2));
0461     hEPFailhitsB0.emplace_back(fs->make<TH2D>(Form("hEPFailhitsB0_layer_%02d", layerList[i]),
0462                                               Form("Rec:FailhitsB0 in EP for layer %d", layerList[i]),
0463                                               640,
0464                                               -3.2,
0465                                               3.2,
0466                                               640,
0467                                               -3.2,
0468                                               3.2));
0469     hEPFailhitsF1.emplace_back(fs->make<TH2D>(Form("hEPFailhitsF1_layer_%02d", layerList[i]),
0470                                               Form("Rec:FailhitsF1 in EP for layer %d", layerList[i]),
0471                                               640,
0472                                               -3.2,
0473                                               3.2,
0474                                               640,
0475                                               -3.2,
0476                                               3.2));
0477     hEPFailhitsCN1.emplace_back(fs->make<TH2D>(Form("hEPFailhitsCN1_layer_%02d", layerList[i]),
0478                                                Form("Rec:FailhitsCN1 in EP for layer %d", layerList[i]),
0479                                                640,
0480                                                -3.2,
0481                                                3.2,
0482                                                640,
0483                                                -3.2,
0484                                                3.2));
0485     hEPFailhitsCK1.emplace_back(fs->make<TH2D>(Form("hEPFailhitsCK1_layer_%02d", layerList[i]),
0486                                                Form("Rec:FailhitsCK1 in EP for layer %d", layerList[i]),
0487                                                640,
0488                                                -3.2,
0489                                                3.2,
0490                                                640,
0491                                                -3.2,
0492                                                3.2));
0493     hEPFailhitsB1.emplace_back(fs->make<TH2D>(Form("hEPFailhitsB1_layer_%02d", layerList[i]),
0494                                               Form("Rec:FailhitsB1 in EP for layer %d", layerList[i]),
0495                                               640,
0496                                               -3.2,
0497                                               3.2,
0498                                               640,
0499                                               -3.2,
0500                                               3.2));
0501 
0502     hELossLayerF0.emplace_back(fs->make<TH1D>(Form("hELossLayerF0_layer_%02d", layerList[i]),
0503                                               Form("Rec:ELossF0 in XY for layer %d", layerList[i]),
0504                                               1000,
0505                                               0.,
0506                                               1000.));
0507     hELossLayerCN0.emplace_back(fs->make<TH1D>(Form("hELossLayerCN0_layer_%02d", layerList[i]),
0508                                                Form("Rec:ELossCN0 in XY for layer %d", layerList[i]),
0509                                                1000,
0510                                                0.,
0511                                                1000.));
0512     hELossLayerCK0.emplace_back(fs->make<TH1D>(Form("hELossLayerCK0_layer_%02d", layerList[i]),
0513                                                Form("Rec:ELossCK0 in XY for layer %d", layerList[i]),
0514                                                1000,
0515                                                0.,
0516                                                1000.));
0517     hELossLayerB0.emplace_back(fs->make<TH1D>(Form("hELossLayerB0_layer_%02d", layerList[i]),
0518                                               Form("Rec:ELossB0 in XY for layer %d", layerList[i]),
0519                                               1000,
0520                                               0.,
0521                                               1000.));
0522     hELossLayerF1.emplace_back(fs->make<TH1D>(Form("hELossLayerF1_layer_%02d", layerList[i]),
0523                                               Form("Rec:ELossF1 in XY for layer %d", layerList[i]),
0524                                               1000,
0525                                               0.,
0526                                               1000.));
0527     hELossLayerCN1.emplace_back(fs->make<TH1D>(Form("hELossLayerCN1_layer_%02d", layerList[i]),
0528                                                Form("Rec:ELossCN1 in XY for layer %d", layerList[i]),
0529                                                1000,
0530                                                0.,
0531                                                1000.));
0532     hELossLayerCK1.emplace_back(fs->make<TH1D>(Form("hELossLayerCK1_layer_%02d", layerList[i]),
0533                                                Form("ELossCK1 in XY for layer %d", layerList[i]),
0534                                                1000,
0535                                                0.,
0536                                                1000.));
0537     hELossLayerB1.emplace_back(fs->make<TH1D>(Form("hELossLayerB1_layer_%02d", layerList[i]),
0538                                               Form("Rec:ELossB1 in XY for layer %d", layerList[i]),
0539                                               1000,
0540                                               0.,
0541                                               1000.));
0542   }
0543 
0544   for (unsigned int i = 0; i < layerList.size(); i++) {
0545     grXYhitsF0.emplace_back(fs->make<TGraph>(0));
0546     grXYhitsF0[i]->SetNameTitle(Form("grXYhitsF0_layer_%02d", layerList[i]),
0547                                 Form("Rec:HitsF0 in XY for layer %d", layerList[i]));
0548     grXYhitsCN0.emplace_back(fs->make<TGraph>(0));
0549     grXYhitsCN0[i]->SetNameTitle(Form("grXYhitsCN0_layer_%02d", layerList[i]),
0550                                  Form("Rec:HitsCN0 in XY for layer %d", layerList[i]));
0551     grXYhitsCK0.emplace_back(fs->make<TGraph>(0));
0552     grXYhitsCK0[i]->SetNameTitle(Form("grXYhitsCK0_layer_%02d", layerList[i]),
0553                                  Form("Rec:HitsCK0 in XY for layer %d", layerList[i]));
0554     grXYhitsB0.emplace_back(fs->make<TGraph>(0));
0555     grXYhitsB0[i]->SetNameTitle(Form("grXYhitsB0_layer_%02d", layerList[i]),
0556                                 Form("Rec:HitsB0 in XY for layer %d", layerList[i]));
0557     grXYhitsAR0.emplace_back(fs->make<TGraph>(0));
0558     grXYhitsAR0[i]->SetNameTitle(Form("grXYhitsAR0_layer_%02d", layerList[i]),
0559                                  Form("Rec:HitsAR0 in XY for layer %d", layerList[i]));
0560     ixyF0[i] = 0;
0561     ixyCN0[i] = 0;
0562     ixyCK0[i] = 0;
0563     ixyB0[i] = 0;
0564     ixyAR0[i] = 0;
0565 
0566     grXYhitsF1.emplace_back(fs->make<TGraph>(0));
0567     grXYhitsF1[i]->SetNameTitle(Form("grXYhitsF1_layer_%02d", layerList[i]),
0568                                 Form("Rec:HitsF1 in XY for layer %d", layerList[i]));
0569     grXYhitsCN1.emplace_back(fs->make<TGraph>(0));
0570     grXYhitsCN1[i]->SetNameTitle(Form("grXYhitsCN1_layer_%02d", layerList[i]),
0571                                  Form("Rec:HitsCN1 in XY for layer %d", layerList[i]));
0572     grXYhitsCK1.emplace_back(fs->make<TGraph>(0));
0573     grXYhitsCK1[i]->SetNameTitle(Form("grXYhitsCK1_layer_%02d", layerList[i]),
0574                                  Form("Rec:HitsCK1 in XY for layer %d", layerList[i]));
0575     grXYhitsB1.emplace_back(fs->make<TGraph>(0));
0576     grXYhitsB1[i]->SetNameTitle(Form("grXYhitsB1_layer_%02d", layerList[i]),
0577                                 Form("Rec:HitsB1 in XY for layer %d", layerList[i]));
0578     grXYhitsAR1.emplace_back(fs->make<TGraph>(0));
0579     grXYhitsAR1[i]->SetNameTitle(Form("grXYhitsAR1_layer_%02d", layerList[i]),
0580                                  Form("Rec:HitsAR1 in XY for layer %d", layerList[i]));
0581     ixyF1[i] = 0;
0582     ixyCN1[i] = 0;
0583     ixyCK1[i] = 0;
0584     ixyB1[i] = 0;
0585     ixyAR1[i] = 0;
0586 
0587     grEtaPhihitsF0.emplace_back(fs->make<TGraph>(0));
0588     grEtaPhihitsF0[i]->SetNameTitle(Form("grEtaPhihitsF0_layer_%02d", layerList[i]),
0589                                     Form("Rec:HitsF0 in XY for layer %d", layerList[i]));
0590     grEtaPhihitsCN0.emplace_back(fs->make<TGraph>(0));
0591     grEtaPhihitsCN0[i]->SetNameTitle(Form("grEtaPhihitsCN0_layer_%02d", layerList[i]),
0592                                      Form("Rec:HitsCN0 in XY for layer %d", layerList[i]));
0593     grEtaPhihitsCK0.emplace_back(fs->make<TGraph>(0));
0594     grEtaPhihitsCK0[i]->SetNameTitle(Form("grEtaPhihitsCK0_layer_%02d", layerList[i]),
0595                                      Form("Rec:HitsCK0 in XY for layer %d", layerList[i]));
0596     grEtaPhihitsB0.emplace_back(fs->make<TGraph>(0));
0597     grEtaPhihitsB0[i]->SetNameTitle(Form("grEtaPhihitsB0_layer_%02d", layerList[i]),
0598                                     Form("Rec:HitsB0 in XY for layer %d", layerList[i]));
0599     iepF0[i] = 0;
0600     iepCN0[i] = 0;
0601     iepCK0[i] = 0;
0602     iepB0[i] = 0;
0603 
0604     grEtaPhihitsF1.emplace_back(fs->make<TGraph>(0));
0605     grEtaPhihitsF1[i]->SetNameTitle(Form("grEtaPhihitsF1_layer_%02d", layerList[i]),
0606                                     Form("Rec:HitsF1 in XY for layer %d", layerList[i]));
0607     grEtaPhihitsCN1.emplace_back(fs->make<TGraph>(0));
0608     grEtaPhihitsCN1[i]->SetNameTitle(Form("grEtaPhihitsCN1_layer_%02d", layerList[i]),
0609                                      Form("Rec:HitsCN1 in XY for layer %d", layerList[i]));
0610     grEtaPhihitsCK1.emplace_back(fs->make<TGraph>(0));
0611     grEtaPhihitsCK1[i]->SetNameTitle(Form("grEtaPhihitsCK1_layer_%02d", layerList[i]),
0612                                      Form("Rec:HitsCK1 in XY for layer %d", layerList[i]));
0613     grEtaPhihitsB1.emplace_back(fs->make<TGraph>(0));
0614     grEtaPhihitsB1[i]->SetNameTitle(Form("grEtaPhihitsB1_layer_%02d", layerList[i]),
0615                                     Form("Rec:HitsB1 in XY for layer %d", layerList[i]));
0616     iepF1[i] = 0;
0617     iepCN1[i] = 0;
0618     iepCK1[i] = 0;
0619     iepB1[i] = 0;
0620   }
0621 
0622   caloGeomToken_ = esConsumes<CaloGeometry, CaloGeometryRecord>();
0623 }
0624 
0625 //
0626 // member functions
0627 //
0628 
0629 // ------------ method called for each event  ------------
0630 void HGCalMTRecoStudy::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0631   using namespace edm;
0632   const CaloGeometry &geomCalo = iSetup.getData(caloGeomToken_);
0633   rhtools_.setGeometry(geomCalo);
0634 
0635   int verbosity_ = 0;
0636   double ElossLayer0[47], ElossLayer1[47];
0637   int nCellsLayer0[47], nCellsLayer1[47];
0638   for (int ilayer = 0; ilayer < 47; ilayer++) {
0639     ElossLayer0[ilayer] = ElossLayer1[ilayer] = 0.0;
0640     nCellsLayer0[ilayer] = nCellsLayer1[ilayer] = 0;
0641   }
0642 
0643   const edm::ESHandle<HGCalGeometry> &geom = iSetup.getHandle(tok_hgcGeom_);
0644   if (!geom.isValid())
0645     edm::LogWarning("HGCalValidation") << "Cannot get valid HGCalGeometry Object for " << nameDetector_;
0646   const HGCalGeometry *geom0 = geom.product();
0647 
0648   const edm::Handle<HGCRecHitCollection> &theRecHitContainers = iEvent.getHandle(recHitSource_);
0649   if (theRecHitContainers.isValid()) {
0650     if (verbosity_ > 0)
0651       edm::LogVerbatim("HGCalValidation") << nameDetector_ << " with " << theRecHitContainers->size() << " element(s)";
0652 
0653     for (const auto &it : *(theRecHitContainers.product())) {
0654       DetId detId = it.id();
0655       GlobalPoint global1 = rhtools_.getPosition(detId);
0656       double energy = it.energy() * 1.e3;
0657 
0658       std::vector<int>::iterator ilyr =
0659           std::find(layerList.begin(), layerList.end(), rhtools_.getLayerWithOffset(detId));
0660 
0661       if (geom0->topology().valid(detId)) {
0662         if (rhtools_.isSilicon(detId)) {
0663           HGCSiliconDetId id(it.id());
0664 
0665           if (id.type() == HGCSiliconDetId::HGCalFine) {
0666             hEF->Fill(energy);
0667 
0668             if (ilyr != layerList.cend()) {
0669               int il = std::distance(layerList.begin(), ilyr);
0670               if (global1.z() < 0.0) {
0671                 grXYhitsF0[il]->SetPoint(ixyF0[il]++, global1.x(), global1.y());
0672                 grEtaPhihitsF0[il]->SetPoint(iepF0[il]++, global1.eta(), global1.phi());
0673                 hXYhitsF0[il]->Fill(global1.x(), global1.y());
0674                 hEPhitsF0[il]->Fill(global1.eta(), global1.phi());
0675                 hELossLayerF0[il]->Fill(energy);
0676               } else {
0677                 grXYhitsF1[il]->SetPoint(ixyF1[il]++, global1.x(), global1.y());
0678                 grEtaPhihitsF1[il]->SetPoint(iepF1[il]++, global1.eta(), global1.phi());
0679                 hXYhitsF1[il]->Fill(global1.x(), global1.y());
0680                 hEPhitsF1[il]->Fill(global1.eta(), global1.phi());
0681                 hELossLayerF1[il]->Fill(energy);
0682               }
0683             }  //isvalid layer
0684           }
0685           if (id.type() == HGCSiliconDetId::HGCalCoarseThin) {
0686             hECN->Fill(energy);
0687 
0688             if (ilyr != layerList.cend()) {
0689               int il = std::distance(layerList.begin(), ilyr);
0690               if (global1.z() < 0.0) {
0691                 grXYhitsCN0[il]->SetPoint(ixyCN0[il]++, global1.x(), global1.y());
0692                 grEtaPhihitsCN0[il]->SetPoint(iepCN0[il]++, global1.eta(), global1.phi());
0693                 hXYhitsCN0[il]->Fill(global1.x(), global1.y());
0694                 hEPhitsCN0[il]->Fill(global1.eta(), global1.phi());
0695                 hELossLayerCN0[il]->Fill(energy);
0696               } else {
0697                 grXYhitsCN1[il]->SetPoint(ixyCN1[il]++, global1.x(), global1.y());
0698                 grEtaPhihitsCN1[il]->SetPoint(iepCN1[il]++, global1.eta(), global1.phi());
0699                 hXYhitsCN1[il]->Fill(global1.x(), global1.y());
0700                 hEPhitsCN1[il]->Fill(global1.eta(), global1.phi());
0701                 hELossLayerCN1[il]->Fill(energy);
0702               }
0703             }  //isvalid layer
0704           }
0705           if (id.type() == HGCSiliconDetId::HGCalCoarseThick) {  //case 2 :
0706             hECK->Fill(energy);
0707 
0708             if (ilyr != layerList.cend()) {
0709               int il = std::distance(layerList.begin(), ilyr);
0710               if (global1.z() < 0.0) {
0711                 grXYhitsCK0[il]->SetPoint(ixyCK0[il]++, global1.x(), global1.y());
0712                 grEtaPhihitsCK0[il]->SetPoint(iepCK0[il]++, global1.eta(), global1.phi());
0713                 hXYhitsCK0[il]->Fill(global1.x(), global1.y());
0714                 hEPhitsCK0[il]->Fill(global1.eta(), global1.phi());
0715                 hELossLayerCK0[il]->Fill(energy);
0716               } else {
0717                 grXYhitsCK1[il]->SetPoint(ixyCK1[il]++, global1.x(), global1.y());
0718                 grEtaPhihitsCK1[il]->SetPoint(iepCK1[il]++, global1.eta(), global1.phi());
0719                 hXYhitsCK1[il]->Fill(global1.x(), global1.y());
0720                 hEPhitsCK1[il]->Fill(global1.eta(), global1.phi());
0721                 hELossLayerCK1[il]->Fill(energy);
0722               }
0723             }  //isvalid layer
0724           }
0725           //The following line by Pruthvi to number the cells with id0 and detId
0726           if (rhtools_.getCell(detId).first + rhtools_.getCell(detId).second <= 2) {
0727             if (ilyr != layerList.cend()) {
0728               int il = std::distance(layerList.begin(), ilyr);
0729               if (global1.z() < 0.0)
0730                 grXYhitsAR0[il]->SetPoint(ixyAR0[il]++, global1.x(), global1.y());
0731               else
0732                 grXYhitsAR1[il]->SetPoint(ixyAR1[il]++, global1.x(), global1.y());
0733             }  //isvalid layer
0734           }
0735         } else if (rhtools_.isScintillator(detId)) {
0736           //HGCScintillatorDetId id(itHit->id());
0737           //int il = rhtools_.getLayerWithOffset(detId);
0738           hESc->Fill(energy);
0739 
0740           if (ilyr != layerList.cend()) {
0741             int il = std::distance(layerList.begin(), ilyr);
0742             if (global1.z() < 0.0) {
0743               grXYhitsB0[il]->SetPoint(ixyB0[il]++, global1.x(), global1.y());
0744               grEtaPhihitsB0[il]->SetPoint(iepB0[il]++, global1.eta(), global1.phi());
0745               hXYhitsB0[il]->Fill(global1.x(), global1.y());
0746               hEPhitsB0[il]->Fill(global1.eta(), global1.phi());
0747               hELossLayerB0[il]->Fill(energy);
0748             } else {
0749               grXYhitsB1[il]->SetPoint(ixyB1[il]++, global1.x(), global1.y());
0750               grEtaPhihitsB1[il]->SetPoint(iepB1[il]++, global1.eta(), global1.phi());
0751               hXYhitsB1[il]->Fill(global1.x(), global1.y());
0752               hEPhitsB1[il]->Fill(global1.eta(), global1.phi());
0753               hELossLayerB1[il]->Fill(energy);
0754             }
0755           }  //isvalid layer
0756 
0757         }  //Silicon or scintillator
0758 
0759         int ilayer = rhtools_.getLayerWithOffset(detId) - 1;
0760         if (global1.z() < 0.0) {
0761           ElossLayer0[ilayer] += energy;
0762           nCellsLayer0[ilayer]++;
0763         } else {
0764           ElossLayer1[ilayer] += energy;
0765           nCellsLayer1[ilayer]++;
0766         }
0767 
0768         ///.................
0769       } else {  //valid topology else invalid
0770 
0771         if (rhtools_.isSilicon(detId)) {
0772           HGCSiliconDetId id(it.id());
0773           //int il = rhtools_.getLayerWithOffset(detId);
0774           if (id.type() == HGCSiliconDetId::HGCalFine) {
0775             if (ilyr != layerList.cend()) {
0776               int il = std::distance(layerList.begin(), ilyr);
0777               if (global1.z() < 0.0) {
0778                 hXYFailhitsF0[il]->Fill(global1.x(), global1.y());
0779                 hEPFailhitsF0[il]->Fill(global1.eta(), global1.phi());
0780               } else {
0781                 hXYFailhitsF1[il]->Fill(global1.x(), global1.y());
0782                 hEPFailhitsF1[il]->Fill(global1.eta(), global1.phi());
0783               }
0784             }  //isvalid layer
0785           }
0786           if (id.type() == HGCSiliconDetId::HGCalCoarseThin) {
0787             if (ilyr != layerList.cend()) {
0788               int il = std::distance(layerList.begin(), ilyr);
0789               if (global1.z() < 0.0) {
0790                 hXYFailhitsCN0[il]->Fill(global1.x(), global1.y());
0791                 hEPFailhitsCN0[il]->Fill(global1.eta(), global1.phi());
0792               } else {
0793                 hXYFailhitsCN1[il]->Fill(global1.x(), global1.y());
0794                 hEPFailhitsCN1[il]->Fill(global1.eta(), global1.phi());
0795               }
0796             }  //isvalid layer
0797           }
0798           if (id.type() == HGCSiliconDetId::HGCalCoarseThick) {  //case 2 :
0799 
0800             if (ilyr != layerList.cend()) {
0801               int il = std::distance(layerList.begin(), ilyr);
0802               if (global1.z() < 0.0) {
0803                 hXYFailhitsCK0[il]->Fill(global1.x(), global1.y());
0804                 hEPFailhitsCK0[il]->Fill(global1.eta(), global1.phi());
0805               } else {
0806                 hXYFailhitsCK1[il]->Fill(global1.x(), global1.y());
0807                 hEPFailhitsCK1[il]->Fill(global1.eta(), global1.phi());
0808               }
0809             }  //isvalid layer
0810           }
0811         } else if (rhtools_.isScintillator(detId)) {
0812           //int il = rhtools_.getLayerWithOffset(detId);
0813           if (ilyr != layerList.cend()) {
0814             int il = std::distance(layerList.begin(), ilyr);
0815             if (global1.z() < 0.0) {
0816               hXYFailhitsB0[il]->Fill(global1.x(), global1.y());
0817               hEPFailhitsB0[il]->Fill(global1.eta(), global1.phi());
0818             } else {
0819               hXYFailhitsB1[il]->Fill(global1.x(), global1.y());
0820               hEPFailhitsB1[il]->Fill(global1.eta(), global1.phi());
0821             }
0822           }  //isvalid layer
0823 
0824         }  //Silicon or scintillator
0825 
0826       }  //invalid topology
0827     }    //loop over iterator
0828   }      //is Valid container
0829 
0830   for (int i = 0; i < 47; i++) {
0831     if (ElossLayer0[i] > 0.0)
0832       hELossLayer0[i]->Fill(ElossLayer0[i]);
0833     if (nCellsLayer0[i] > 0)
0834       hNCellsLayer0[i]->Fill(nCellsLayer0[i]);
0835 
0836     if (ElossLayer1[i] > 0.0)
0837       hELossLayer1[i]->Fill(ElossLayer1[i]);
0838     if (nCellsLayer1[i] > 0)
0839       hNCellsLayer1[i]->Fill(nCellsLayer1[i]);
0840   }
0841 }
0842 
0843 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0844 void HGCalMTRecoStudy::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0845   edm::ParameterSetDescription desc;
0846   desc.add<std::string>("detectorName", "HGCalEESensitive");
0847   desc.add<edm::InputTag>("source", edm::InputTag("HGCalRecHit", "HGCEERecHits"));
0848   desc.add<std::string>("layerList", "1");
0849   descriptions.add("hgcalMTRecoStudy", desc);
0850 }
0851 
0852 //define this as a plug-in
0853 DEFINE_FWK_MODULE(HGCalMTRecoStudy);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!