Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-09-07 04:34:59

0001 // system include files
0002 #include <map>
0003 #include <memory>
0004 #include <string>
0005 #include <iostream>
0006 #include <fstream>
0007 #include <sstream>
0008 #include <vector>
0009 
0010 // user include files
0011 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0012 #include "CondFormats/DataRecord/interface/HcalRespCorrsRcd.h"
0013 #include "CondFormats/DataRecord/interface/L1GtTriggerMenuRcd.h"
0014 #include "CondFormats/HcalObjects/interface/HcalRespCorrs.h"
0015 #include "CondFormats/L1TObjects/interface/L1GtTriggerMenu.h"
0016 
0017 #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
0018 #include "DataFormats/CaloTowers/interface/CaloTowerDetId.h"
0019 #include "DataFormats/Candidate/interface/Candidate.h"
0020 #include "DataFormats/Common/interface/Ref.h"
0021 #include "DataFormats/DetId/interface/DetId.h"
0022 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0023 #include "DataFormats/FEDRawData/interface/FEDRawData.h"
0024 #include "DataFormats/FEDRawData/interface/FEDRawDataCollection.h"
0025 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
0026 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0027 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0028 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
0029 #include "DataFormats/HcalDigi/interface/HcalCalibrationEventTypes.h"
0030 #include "DataFormats/JetReco/interface/Jet.h"
0031 #include "DataFormats/JetReco/interface/CaloJet.h"
0032 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
0033 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
0034 #include "DataFormats/L1GlobalTrigger/interface/L1GtfeWord.h"
0035 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
0036 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerRecord.h"
0037 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerObjectMap.h"
0038 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerObjectMapRecord.h"
0039 #include "DataFormats/Provenance/interface/StableProvenance.h"
0040 #include "DataFormats/TrackReco/interface/Track.h"
0041 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
0042 
0043 #include "EventFilter/HcalRawToDigi/interface/HcalDCCHeader.h"
0044 
0045 #include "FWCore/Common/interface/TriggerNames.h"
0046 #include "FWCore/Framework/interface/Frameworkfwd.h"
0047 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0048 #include "FWCore/Framework/interface/Event.h"
0049 #include "FWCore/Framework/interface/EventSetup.h"
0050 #include "FWCore/Framework/interface/Run.h"
0051 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0052 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0053 #include "FWCore/ServiceRegistry/interface/Service.h"
0054 
0055 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0056 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0057 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0058 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0059 
0060 #include "TFile.h"
0061 #include "TH1.h"
0062 #include "TH2.h"
0063 #include "TTree.h"
0064 
0065 //
0066 // class declaration
0067 //
0068 namespace cms {
0069   class Analyzer_minbias : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::SharedResources> {
0070   public:
0071     explicit Analyzer_minbias(const edm::ParameterSet&);
0072     ~Analyzer_minbias() override;
0073 
0074     void beginJob() override;
0075     void analyze(edm::Event const&, edm::EventSetup const&) override;
0076     void beginRun(edm::Run const&, edm::EventSetup const&) override;
0077     void endRun(edm::Run const&, edm::EventSetup const&) override;
0078     void endJob() override;
0079 
0080   private:
0081     // ----------member data ---------------------------
0082     std::string fOutputFileName_;
0083     bool theRecalib_;
0084     std::string hcalfile_;
0085     std::ofstream* myout_hcal;
0086 
0087     // names of modules, producing object collections
0088     TFile* hOutputFile;
0089     TTree* myTree;
0090     TH1F* hCalo1[73][43];
0091     TH1F* hCalo2[73][43];
0092     TH1F* hCalo1mom2[73][43];
0093     TH1F* hCalo2mom2[73][43];
0094     TH1F* hbheNoiseE;
0095     TH1F* hbheSignalE;
0096     TH1F* hfNoiseE;
0097     TH1F* hfSignalE;
0098 
0099     TH2F* hHBHEsize_vs_run;
0100     TH2F* hHFsize_vs_run;
0101     // Root tree members
0102     int nevent_run;
0103     int mydet, mysubd, depth, iphi, ieta;
0104     float phi, eta;
0105     float mom0_MB, mom1_MB, mom2_MB, mom3_MB, mom4_MB, occup;
0106     float mom0_Noise, mom1_Noise, mom2_Noise, mom3_Noise, mom4_Noise;
0107     float mom0_Diff, mom1_Diff, mom2_Diff, mom3_Diff, mom4_Diff;
0108 
0109     // Noise subtraction
0110 
0111     double meannoise_pl[73][43], meannoise_min[73][43];
0112     double noise_pl[73][43], noise_min[73][43];
0113 
0114     // counters
0115 
0116     double nevent;
0117     double theMBFillDetMapPl0[5][5][73][43];
0118     double theMBFillDetMapPl1[5][5][73][43];
0119     double theMBFillDetMapPl2[5][5][73][43];
0120     double theMBFillDetMapPl4[5][5][73][43];
0121 
0122     double theMBFillDetMapMin0[5][5][73][43];
0123     double theMBFillDetMapMin1[5][5][73][43];
0124     double theMBFillDetMapMin2[5][5][73][43];
0125     double theMBFillDetMapMin4[5][5][73][43];
0126 
0127     double theNSFillDetMapPl0[5][5][73][43];
0128     double theNSFillDetMapPl1[5][5][73][43];
0129     double theNSFillDetMapPl2[5][5][73][43];
0130     double theNSFillDetMapPl4[5][5][73][43];
0131 
0132     double theNSFillDetMapMin0[5][5][73][43];
0133     double theNSFillDetMapMin1[5][5][73][43];
0134     double theNSFillDetMapMin2[5][5][73][43];
0135     double theNSFillDetMapMin4[5][5][73][43];
0136 
0137     double theDFFillDetMapPl0[5][5][73][43];
0138     double theDFFillDetMapPl1[5][5][73][43];
0139     double theDFFillDetMapPl2[5][5][73][43];
0140     double theDFFillDetMapMin0[5][5][73][43];
0141     double theDFFillDetMapMin1[5][5][73][43];
0142     double theDFFillDetMapMin2[5][5][73][43];
0143 
0144     const edm::EDGetTokenT<HBHERecHitCollection> tok_hbhe_;
0145     const edm::EDGetTokenT<HORecHitCollection> tok_ho_;
0146     const edm::EDGetTokenT<HFRecHitCollection> tok_hf_;
0147 
0148     const edm::EDGetTokenT<FEDRawDataCollection> tok_data_;
0149 
0150     const edm::EDGetTokenT<HBHERecHitCollection> tok_hbheNoise_;
0151     const edm::EDGetTokenT<HORecHitCollection> tok_hoNoise_;
0152     const edm::EDGetTokenT<HFRecHitCollection> tok_hfNoise_;
0153 
0154     //
0155     const edm::EDGetTokenT<L1GlobalTriggerReadoutRecord> tok_gtRec_;
0156     const edm::EDGetTokenT<HBHERecHitCollection> tok_hbheNorm_;
0157 
0158     const edm::ESGetToken<HcalRespCorrs, HcalRespCorrsRcd> tok_respCorr_;
0159     const edm::ESGetToken<L1GtTriggerMenu, L1GtTriggerMenuRcd> tok_l1gt_;
0160   };
0161 }  // namespace cms
0162 
0163 //
0164 // constructors and destructor
0165 //
0166 namespace cms {
0167   Analyzer_minbias::Analyzer_minbias(const edm::ParameterSet& iConfig)
0168       : fOutputFileName_(iConfig.getUntrackedParameter<std::string>("HistOutFile")),
0169         theRecalib_(iConfig.getParameter<bool>("Recalib")),
0170         tok_hbhe_(consumes<HBHERecHitCollection>(iConfig.getParameter<edm::InputTag>("hbheInputMB"))),
0171         tok_ho_(consumes<HORecHitCollection>(iConfig.getParameter<edm::InputTag>("hoInputMB"))),
0172         tok_hf_(consumes<HFRecHitCollection>(iConfig.getParameter<edm::InputTag>("hfInputMB"))),
0173         tok_data_(consumes<FEDRawDataCollection>(edm::InputTag(iConfig.getParameter<std::string>("InputLabel")))),
0174         tok_hbheNoise_(consumes<HBHERecHitCollection>(iConfig.getParameter<edm::InputTag>("hbheInputNoise"))),
0175         tok_hoNoise_(consumes<HORecHitCollection>(iConfig.getParameter<edm::InputTag>("hoInputNoise"))),
0176         tok_hfNoise_(consumes<HFRecHitCollection>(iConfig.getParameter<edm::InputTag>("hfInputNoise"))),
0177         tok_gtRec_(consumes<L1GlobalTriggerReadoutRecord>(edm::InputTag("gtDigisAlCaMB"))),
0178         tok_hbheNorm_(consumes<HBHERecHitCollection>(edm::InputTag("hbhereco"))),
0179         tok_respCorr_(esConsumes<HcalRespCorrs, HcalRespCorrsRcd>()),
0180         tok_l1gt_(esConsumes<L1GtTriggerMenu, L1GtTriggerMenuRcd>()) {
0181     usesResource(TFileService::kSharedResource);
0182     // get name of output file with histogramms
0183     // get names of modules, producing object collections
0184     // some of the label names are hardcodded..
0185     //
0186     for (int i = 0; i < 73; i++) {
0187       for (int j = 0; j < 43; j++) {
0188         noise_min[i][j] = 0.;
0189         noise_pl[i][j] = 0.;
0190       }
0191     }
0192   }
0193 
0194   Analyzer_minbias::~Analyzer_minbias() {
0195     // do anything here that needs to be done at desctruction time
0196     // (e.g. close files, deallocate resources etc.)
0197   }
0198 
0199   void Analyzer_minbias::beginRun(const edm::Run&, const edm::EventSetup&) { nevent_run = 0; }
0200   void Analyzer_minbias::endRun(const edm::Run& r, const edm::EventSetup&) {
0201     edm::LogVerbatim("AnalyzerMB") << " Runnumber " << r.run() << " Nevents  " << nevent_run;
0202   }
0203 
0204   void Analyzer_minbias::beginJob() {
0205     edm::Service<TFileService> fs;
0206     myTree = fs->make<TTree>("RecJet", "RecJet Tree");
0207     myTree->Branch("mydet", &mydet, "mydet/I");
0208     myTree->Branch("mysubd", &mysubd, "mysubd/I");
0209     myTree->Branch("depth", &depth, "depth/I");
0210     myTree->Branch("ieta", &ieta, "ieta/I");
0211     myTree->Branch("iphi", &iphi, "iphi/I");
0212     myTree->Branch("eta", &eta, "eta/F");
0213     myTree->Branch("phi", &phi, "phi/F");
0214 
0215     myTree->Branch("mom0_MB", &mom0_MB, "mom0_MB/F");
0216     myTree->Branch("mom1_MB", &mom1_MB, "mom1_MB/F");
0217     myTree->Branch("mom2_MB", &mom2_MB, "mom2_MB/F");
0218     myTree->Branch("mom4_MB", &mom4_MB, "mom4_MB/F");
0219 
0220     myTree->Branch("mom0_Noise", &mom0_Noise, "mom0_Noise/F");
0221     myTree->Branch("mom1_Noise", &mom1_Noise, "mom1_Noise/F");
0222     myTree->Branch("mom2_Noise", &mom2_Noise, "mom2_Noise/F");
0223     myTree->Branch("mom4_Noise", &mom4_Noise, "mom4_Noise/F");
0224 
0225     myTree->Branch("mom0_Diff", &mom0_Diff, "mom0_Diff/F");
0226     myTree->Branch("mom1_Diff", &mom1_Diff, "mom1_Diff/F");
0227     myTree->Branch("mom2_Diff", &mom2_Diff, "mom2_Diff/F");
0228 
0229     myTree->Branch("occup", &occup, "occup/F");
0230 
0231     edm::LogVerbatim("AnalyzerMB") << " Before ordering Histos ";
0232 
0233     char str0[32];
0234     char str1[32];
0235 
0236     char str10[32];
0237     char str11[32];
0238 
0239     int k = 0;
0240     nevent = 0;
0241     // Size of collections
0242 
0243     hHBHEsize_vs_run =
0244         fs->make<TH2F>("hHBHEsize_vs_run", "hHBHEsize_vs_run", 500, 111500., 112000., 6101, -100.5, 6000.5);
0245     hHFsize_vs_run = fs->make<TH2F>("hHFsize_vs_run", "hHFsize_vs_run", 500, 111500., 112000., 6101, -100.5, 6000.5);
0246 
0247     for (int i = 1; i < 73; i++) {
0248       for (int j = 1; j < 43; j++) {
0249         meannoise_pl[i][j] = 0.;
0250         meannoise_min[i][j] = 0.;
0251 
0252         k = i * 1000 + j;
0253         sprintf(str0, "mpl%d", k);
0254         sprintf(str1, "mmin%d", k);
0255 
0256         sprintf(str10, "vpl%d", k);
0257         sprintf(str11, "vmin%d", k);
0258         if (j < 30) {
0259           // first order moment
0260           hCalo1[i][j] = fs->make<TH1F>(str0, "h0", 320, -10., 10.);
0261           hCalo2[i][j] = fs->make<TH1F>(str1, "h1", 320, -10., 10.);
0262 
0263           // second order moment
0264           hCalo1mom2[i][j] = fs->make<TH1F>(str10, "h10", 320, 0., 20.);
0265           hCalo2mom2[i][j] = fs->make<TH1F>(str11, "h11", 320, 0., 20.);
0266         } else {
0267           // HF
0268           // first order moment
0269           if (j < 40) {
0270             hCalo1[i][j] = fs->make<TH1F>(str0, "h0", 320, -10., 10.);
0271             hCalo2[i][j] = fs->make<TH1F>(str1, "h1", 320, -10., 10.);
0272             //
0273             // second order moment
0274             hCalo1mom2[i][j] = fs->make<TH1F>(str10, "h10", 320, 0., 40.);
0275             hCalo2mom2[i][j] = fs->make<TH1F>(str11, "h11", 320, 0., 40.);
0276           } else {
0277             hCalo1[i][j] = fs->make<TH1F>(str0, "h0", 320, -10., 10.);
0278             hCalo2[i][j] = fs->make<TH1F>(str1, "h1", 320, -10., 10.);
0279 
0280             // second order moment
0281             hCalo1mom2[i][j] = fs->make<TH1F>(str10, "h10", 320, 0., 120.);
0282             hCalo2mom2[i][j] = fs->make<TH1F>(str11, "h11", 320, 0., 120.);
0283           }
0284         }  // HE/HF boundary
0285       }  // j
0286     }  // i
0287 
0288     hbheNoiseE = fs->make<TH1F>("hbheNoiseE", "hbheNoiseE", 320, -10., 10.);
0289     hfNoiseE = fs->make<TH1F>("hfNoiseE", "hfNoiseE", 320, -10., 10.);
0290     hbheSignalE = fs->make<TH1F>("hbheSignalE", "hbheSignalE", 320, -10., 10.);
0291     hfSignalE = fs->make<TH1F>("hfSignalE", "hfSignalE", 320, -10., 10.);
0292 
0293     edm::LogVerbatim("AnalyzerMB") << " After ordering Histos ";
0294 
0295     std::string ccc = "noise_0.dat";
0296 
0297     myout_hcal = new std::ofstream(ccc.c_str());
0298     if (!myout_hcal)
0299       edm::LogVerbatim("AnalyzerMB") << " Output file not open!!! ";
0300 
0301     //
0302     for (int i = 0; i < 5; i++) {
0303       for (int j = 0; j < 5; j++) {
0304         for (int k = 0; k < 73; k++) {
0305           for (int l = 0; l < 43; l++) {
0306             theMBFillDetMapPl0[i][j][k][l] = 0.;
0307             theMBFillDetMapPl1[i][j][k][l] = 0.;
0308             theMBFillDetMapPl2[i][j][k][l] = 0.;
0309             theMBFillDetMapPl4[i][j][k][l] = 0.;
0310 
0311             theMBFillDetMapMin0[i][j][k][l] = 0.;
0312             theMBFillDetMapMin1[i][j][k][l] = 0.;
0313             theMBFillDetMapMin2[i][j][k][l] = 0.;
0314             theMBFillDetMapMin4[i][j][k][l] = 0.;
0315 
0316             theNSFillDetMapPl0[i][j][k][l] = 0.;
0317             theNSFillDetMapPl1[i][j][k][l] = 0.;
0318             theNSFillDetMapPl2[i][j][k][l] = 0.;
0319             theNSFillDetMapPl4[i][j][k][l] = 0.;
0320 
0321             theNSFillDetMapMin0[i][j][k][l] = 0.;
0322             theNSFillDetMapMin1[i][j][k][l] = 0.;
0323             theNSFillDetMapMin2[i][j][k][l] = 0.;
0324             theNSFillDetMapMin4[i][j][k][l] = 0.;
0325 
0326             theDFFillDetMapPl0[i][j][k][l] = 0.;
0327             theDFFillDetMapPl1[i][j][k][l] = 0.;
0328             theDFFillDetMapPl2[i][j][k][l] = 0.;
0329             theDFFillDetMapMin0[i][j][k][l] = 0.;
0330             theDFFillDetMapMin1[i][j][k][l] = 0.;
0331             theDFFillDetMapMin2[i][j][k][l] = 0.;
0332           }
0333         }
0334       }
0335     }
0336 
0337     return;
0338   }
0339   //
0340   //  EndJob
0341   //
0342   void Analyzer_minbias::endJob() {
0343     int ii = 0;
0344 
0345     for (int i = 1; i < 5; i++) {
0346       for (int j = 1; j < 5; j++) {
0347         for (int k = 1; k < 73; k++) {
0348           for (int l = 1; l < 43; l++) {
0349             if (theMBFillDetMapPl0[i][j][k][l] > 0) {
0350               mom0_MB = theMBFillDetMapPl0[i][j][k][l];
0351               mom1_MB = theMBFillDetMapPl1[i][j][k][l];
0352               mom2_MB = theMBFillDetMapPl2[i][j][k][l];
0353               mom4_MB = theMBFillDetMapPl4[i][j][k][l];
0354               mom0_Noise = theNSFillDetMapPl0[i][j][k][l];
0355               mom1_Noise = theNSFillDetMapPl1[i][j][k][l];
0356               mom2_Noise = theNSFillDetMapPl2[i][j][k][l];
0357               mom4_Noise = theNSFillDetMapPl4[i][j][k][l];
0358               mom0_Diff = theDFFillDetMapPl0[i][j][k][l];
0359               mom1_Diff = theDFFillDetMapPl1[i][j][k][l];
0360               mom2_Diff = theDFFillDetMapPl2[i][j][k][l];
0361 
0362               mysubd = i;
0363               depth = j;
0364               ieta = l;
0365               iphi = k;
0366               edm::LogVerbatim("AnalyzerMB") << " Result Plus= " << mysubd << " " << ieta << " " << iphi << " mom0  "
0367                                              << mom0_MB << " mom1 " << mom1_MB << " mom2 " << mom2_MB;
0368               myTree->Fill();
0369               ii++;
0370             }  // Pl > 0
0371 
0372             if (theMBFillDetMapMin0[i][j][k][l] > 0) {
0373               mom0_MB = theMBFillDetMapMin0[i][j][k][l];
0374               mom1_MB = theMBFillDetMapMin1[i][j][k][l];
0375               mom2_MB = theMBFillDetMapMin2[i][j][k][l];
0376               mom4_MB = theMBFillDetMapMin4[i][j][k][l];
0377               mom0_Noise = theNSFillDetMapMin0[i][j][k][l];
0378               mom1_Noise = theNSFillDetMapMin1[i][j][k][l];
0379               mom2_Noise = theNSFillDetMapMin2[i][j][k][l];
0380               mom4_Noise = theNSFillDetMapMin4[i][j][k][l];
0381               mom0_Diff = theDFFillDetMapMin0[i][j][k][l];
0382               mom1_Diff = theDFFillDetMapMin1[i][j][k][l];
0383               mom2_Diff = theDFFillDetMapMin2[i][j][k][l];
0384 
0385               mysubd = i;
0386               depth = j;
0387               ieta = -1 * l;
0388               iphi = k;
0389               edm::LogVerbatim("AnalyzerMB") << " Result Minus= " << mysubd << " " << ieta << " " << iphi << " mom0  "
0390                                              << mom0_MB << " mom1 " << mom1_MB << " mom2 " << mom2_MB;
0391               myTree->Fill();
0392               ii++;
0393 
0394             }  // Min>0
0395           }  // ieta
0396         }  // iphi
0397       }  // depth
0398     }  //subd
0399 
0400     edm::LogVerbatim("AnalyzerMB") << " Number of cells " << ii;
0401 
0402     for (int i = 1; i < 73; i++) {
0403       for (int j = 1; j < 43; j++) {
0404         hCalo1[i][j]->Write();
0405         hCalo2[i][j]->Write();
0406         hCalo1mom2[i][j]->Write();
0407         hCalo2mom2[i][j]->Write();
0408       }
0409     }
0410 
0411     edm::LogVerbatim("AnalyzerMB") << " File is closed ";
0412 
0413     return;
0414   }
0415 
0416   //
0417   // member functions
0418   //
0419 
0420   // ------------ method called to produce the data  ------------
0421   void Analyzer_minbias::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0422     edm::LogVerbatim("AnalyzerMB") << " Start Analyzer_minbias::analyze " << nevent;
0423     nevent++;
0424     nevent_run++;
0425 
0426     float rnnum = (float)iEvent.run();
0427 
0428     std::vector<edm::StableProvenance const*> theProvenance;
0429     iEvent.getAllStableProvenance(theProvenance);
0430 
0431     for (auto const& provenance : theProvenance) {
0432       edm::LogVerbatim("AnalyzerMB") << " Print all process/modulelabel/product names " << provenance->processName()
0433                                      << " , " << provenance->moduleLabel() << " , "
0434                                      << provenance->productInstanceName();
0435     }
0436     const HcalRespCorrs* myRecalib = nullptr;
0437     if (theRecalib_) {
0438       myRecalib = &iSetup.getData(tok_respCorr_);
0439     }  // theRecalib_
0440 
0441     // Noise part for HB HE
0442 
0443     double tmpNSFillDetMapPl1[5][5][73][43];
0444     double tmpNSFillDetMapMin1[5][5][73][43];
0445 
0446     for (int i = 0; i < 5; i++) {
0447       for (int j = 0; j < 5; j++) {
0448         for (int k = 0; k < 73; k++) {
0449           for (int l = 0; l < 43; l++) {
0450             tmpNSFillDetMapPl1[i][j][k][l] = 0.;
0451             tmpNSFillDetMapMin1[i][j][k][l] = 0.;
0452           }
0453         }
0454       }
0455     }
0456 
0457     const edm::Handle<HBHERecHitCollection> hbheNormal = iEvent.getHandle(tok_hbheNorm_);
0458     if (!hbheNormal.isValid()) {
0459       edm::LogWarning("AnalyzerMB") << " hbheNormal failed ";
0460     } else {
0461       edm::LogVerbatim("AnalyzerMB") << " The size of the normal collection " << hbheNormal->size();
0462     }
0463 
0464     const edm::Handle<HBHERecHitCollection> hbheNS = iEvent.getHandle(tok_hbheNoise_);
0465 
0466     if (!hbheNS.isValid()) {
0467       edm::LogWarning("AnalyzerMB") << "HcalCalibAlgos: Error! can't get hbhe"
0468                                     << " product! No HBHE MS ";
0469       return;
0470     }
0471 
0472     const HBHERecHitCollection HithbheNS = *(hbheNS.product());
0473     edm::LogVerbatim("AnalyzerMB") << " HBHE NS size of collection " << HithbheNS.size();
0474     hHBHEsize_vs_run->Fill(rnnum, (float)HithbheNS.size());
0475 
0476     if (HithbheNS.size() != 5184) {
0477       edm::LogWarning("AnalyzerMB") << " HBHE problem " << rnnum << " " << HithbheNS.size();
0478     }
0479     const edm::Handle<HBHERecHitCollection> hbheMB = iEvent.getHandle(tok_hbhe_);
0480 
0481     if (!hbheMB.isValid()) {
0482       edm::LogWarning("AnalyzerMB") << "HcalCalibAlgos: Error! can't get hbhe"
0483                                     << " product! No HBHE MB";
0484     }
0485 
0486     const HBHERecHitCollection HithbheMB = *(hbheMB.product());
0487     edm::LogVerbatim("AnalyzerMB") << " HBHE MB size of collection " << HithbheMB.size();
0488     if (HithbheMB.size() != 5184) {
0489       edm::LogWarning("AnalyzerMB") << " HBHE problem " << rnnum << " " << HithbheMB.size();
0490     }
0491 
0492     const edm::Handle<HFRecHitCollection> hfNS = iEvent.getHandle(tok_hfNoise_);
0493 
0494     if (!hfNS.isValid()) {
0495       edm::LogWarning("AnalyzerMB") << "HcalCalibAlgos: Error! can't get hf"
0496                                     << " product! No HF NS ";
0497     }
0498 
0499     const HFRecHitCollection HithfNS = *(hfNS.product());
0500     edm::LogVerbatim("AnalyzerMB") << " HFE NS size of collection " << HithfNS.size();
0501     hHFsize_vs_run->Fill(rnnum, (float)HithfNS.size());
0502     if (HithfNS.size() != 1728) {
0503       edm::LogWarning("AnalyzerMB") << " HF problem " << rnnum << " " << HithfNS.size();
0504     }
0505 
0506     const edm::Handle<HFRecHitCollection> hfMB = iEvent.getHandle(tok_hf_);
0507 
0508     if (!hfMB.isValid()) {
0509       edm::LogWarning("AnalyzerMB") << "HcalCalibAlgos: Error! can't get hf"
0510                                     << " product! No HF MB";
0511     }
0512 
0513     const HFRecHitCollection HithfMB = *(hfMB.product());
0514     edm::LogVerbatim("AnalyzerMB") << " HF MB size of collection " << HithfMB.size();
0515     if (HithfMB.size() != 1728) {
0516       edm::LogWarning("AnalyzerMB") << " HF problem " << rnnum << " " << HithfMB.size();
0517     }
0518 
0519     for (HBHERecHitCollection::const_iterator hbheItr = HithbheNS.begin(); hbheItr != HithbheNS.end(); hbheItr++) {
0520       // Recalibration of energy
0521       float icalconst = 1.;
0522       DetId mydetid = hbheItr->id().rawId();
0523       if (theRecalib_)
0524         icalconst = myRecalib->getValues(mydetid)->getValue();
0525 
0526       HBHERecHit aHit(hbheItr->id(), hbheItr->energy() * icalconst, hbheItr->time());
0527 
0528       double energyhit = aHit.energy();
0529 
0530       DetId id = (*hbheItr).detid();
0531       HcalDetId hid = HcalDetId(id);
0532 
0533       int mysu = ((hid).rawId() >> 25) & 0x7;
0534       if (hid.ieta() > 0) {
0535         theNSFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0536             theNSFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] + 1.;
0537         theNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0538             theNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] + energyhit;
0539         theNSFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0540             theNSFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] + pow(energyhit, 2);
0541         theNSFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0542             theNSFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] + pow(energyhit, 4);
0543 
0544         tmpNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = energyhit;
0545 
0546       } else {
0547         theNSFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0548             theNSFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + 1.;
0549         theNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0550             theNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + energyhit;
0551         theNSFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0552             theNSFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + pow(energyhit, 2);
0553         theNSFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0554             theNSFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + pow(energyhit, 4);
0555 
0556         tmpNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = energyhit;
0557       }
0558 
0559       if (hid.depth() == 1) {
0560         hbheNoiseE->Fill(energyhit);
0561 
0562         if (energyhit < -2.)
0563           edm::LogVerbatim("AnalyzerMB") << " Run " << rnnum << " ieta,iphi " << hid.ieta() << " " << hid.iphi()
0564                                          << energyhit;
0565 
0566       }  // depth=1
0567 
0568     }  // HBHE_NS
0569 
0570     // Signal part for HB HE
0571 
0572     for (HBHERecHitCollection::const_iterator hbheItr = HithbheMB.begin(); hbheItr != HithbheMB.end(); hbheItr++) {
0573       // Recalibration of energy
0574       float icalconst = 1.;
0575       DetId mydetid = hbheItr->id().rawId();
0576       if (theRecalib_)
0577         icalconst = myRecalib->getValues(mydetid)->getValue();
0578 
0579       HBHERecHit aHit(hbheItr->id(), hbheItr->energy() * icalconst, hbheItr->time());
0580 
0581       double energyhit = aHit.energy();
0582 
0583       DetId id = (*hbheItr).detid();
0584       HcalDetId hid = HcalDetId(id);
0585 
0586       int mysu = ((hid).rawId() >> 25) & 0x7;
0587       if (hid.ieta() > 0) {
0588         theMBFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] += 1.;
0589         theMBFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] += energyhit;
0590         theMBFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] += pow(energyhit, 2);
0591         theMBFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] += pow(energyhit, 4);
0592         float mydiff = energyhit - tmpNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()];
0593 
0594         theDFFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0595             theDFFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] + mydiff;
0596         theDFFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0597             theDFFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] + pow(mydiff, 2);
0598       } else {
0599         theMBFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0600             theMBFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + 1.;
0601         theMBFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0602             theMBFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + energyhit;
0603         theMBFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0604             theMBFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + pow(energyhit, 2);
0605         theMBFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0606             theMBFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + pow(energyhit, 4);
0607 
0608         float mydiff = energyhit - tmpNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()];
0609         theDFFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0610             theDFFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()] + mydiff;
0611         theDFFillDetMapMin2[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0612             theDFFillDetMapMin2[mysu][hid.depth()][hid.iphi()][hid.ieta()] + pow(mydiff, 2);
0613       }
0614 
0615       if (hid.depth() == 1) {
0616         hbheSignalE->Fill(energyhit);
0617 
0618         if (hid.ieta() > 0) {
0619           hCalo1[hid.iphi()][hid.ieta()]->Fill(energyhit);
0620           hCalo1mom2[hid.iphi()][hid.ieta()]->Fill(pow(energyhit, 2));
0621         } else {
0622           hCalo2[hid.iphi()][abs(hid.ieta())]->Fill(energyhit);
0623           hCalo2mom2[hid.iphi()][abs(hid.ieta())]->Fill(pow(energyhit, 2));
0624         }  // eta><0
0625 
0626       }  // depth=1
0627 
0628     }  // HBHE_MB
0629 
0630     // HF
0631 
0632     for (HFRecHitCollection::const_iterator hbheItr = HithfNS.begin(); hbheItr != HithfNS.end(); hbheItr++) {
0633       // Recalibration of energy
0634       float icalconst = 1.;
0635       DetId mydetid = hbheItr->id().rawId();
0636       if (theRecalib_)
0637         icalconst = myRecalib->getValues(mydetid)->getValue();
0638 
0639       HFRecHit aHit(hbheItr->id(), hbheItr->energy() * icalconst, hbheItr->time());
0640 
0641       double energyhit = aHit.energy();
0642       //
0643       // Remove PMT hits
0644       //
0645       DetId id = (*hbheItr).detid();
0646       HcalDetId hid = HcalDetId(id);
0647 
0648       if (fabs(energyhit) > 40.)
0649         continue;
0650 
0651       int mysu = hid.subdetId();
0652       if (hid.ieta() > 0) {
0653         theNSFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0654             theNSFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] + 1.;
0655         theNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0656             theNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] + energyhit;
0657         theNSFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0658             theNSFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] + pow(energyhit, 2);
0659         theNSFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0660             theNSFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] + pow(energyhit, 4);
0661 
0662         tmpNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = energyhit;
0663 
0664       } else {
0665         theNSFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0666             theNSFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + 1.;
0667         theNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0668             theNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + energyhit;
0669         theNSFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0670             theNSFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + pow(energyhit, 2);
0671         theNSFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0672             theNSFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + pow(energyhit, 4);
0673 
0674         tmpNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()] = energyhit;
0675       }
0676 
0677       if (hid.depth() == 1) {
0678         hfNoiseE->Fill(energyhit);
0679 
0680       }  // depth=1
0681 
0682     }  // HBHE_NS
0683 
0684     // Signal part for HB HE
0685 
0686     for (HFRecHitCollection::const_iterator hbheItr = HithfMB.begin(); hbheItr != HithfMB.end(); hbheItr++) {
0687       // Recalibration of energy
0688       float icalconst = 1.;
0689       DetId mydetid = hbheItr->id().rawId();
0690       if (theRecalib_)
0691         icalconst = myRecalib->getValues(mydetid)->getValue();
0692 
0693       HFRecHit aHit(hbheItr->id(), hbheItr->energy() * icalconst, hbheItr->time());
0694 
0695       double energyhit = aHit.energy();
0696       //
0697       // Remove PMT hits
0698       //
0699       if (fabs(energyhit) > 40.)
0700         continue;
0701 
0702       DetId id = (*hbheItr).detid();
0703       HcalDetId hid = HcalDetId(id);
0704 
0705       int mysu = ((hid).rawId() >> 25) & 0x7;
0706       if (hid.ieta() > 0) {
0707         theMBFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0708             theMBFillDetMapPl0[mysu][hid.depth()][hid.iphi()][hid.ieta()] + 1.;
0709         theMBFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0710             theMBFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] + energyhit;
0711         theMBFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0712             theMBFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] + pow(energyhit, 2);
0713         theMBFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0714             theMBFillDetMapPl4[mysu][hid.depth()][hid.iphi()][hid.ieta()] + pow(energyhit, 4);
0715 
0716         theDFFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0717             theDFFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()] + energyhit -
0718             tmpNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()];
0719         theDFFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0720             theDFFillDetMapPl2[mysu][hid.depth()][hid.iphi()][hid.ieta()] +
0721             pow((energyhit - tmpNSFillDetMapPl1[mysu][hid.depth()][hid.iphi()][hid.ieta()]), 2);
0722       } else {
0723         theMBFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0724             theMBFillDetMapMin0[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + 1.;
0725         theMBFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0726             theMBFillDetMapMin1[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + energyhit;
0727         theMBFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0728             theMBFillDetMapMin2[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + pow(energyhit, 2);
0729         theMBFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] =
0730             theMBFillDetMapMin4[mysu][hid.depth()][hid.iphi()][abs(hid.ieta())] + pow(energyhit, 4);
0731 
0732         theDFFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0733             theDFFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()] + energyhit -
0734             tmpNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()];
0735         theDFFillDetMapMin2[mysu][hid.depth()][hid.iphi()][hid.ieta()] =
0736             theDFFillDetMapMin2[mysu][hid.depth()][hid.iphi()][hid.ieta()] +
0737             pow((energyhit - tmpNSFillDetMapMin1[mysu][hid.depth()][hid.iphi()][hid.ieta()]), 2);
0738       }
0739 
0740       if (hid.depth() == 1) {
0741         hfSignalE->Fill(energyhit);
0742 
0743         if (hid.ieta() > 0) {
0744           hCalo1[hid.iphi()][hid.ieta()]->Fill(energyhit);
0745           hCalo1mom2[hid.iphi()][hid.ieta()]->Fill(pow(energyhit, 2));
0746         } else {
0747           hCalo2[hid.iphi()][abs(hid.ieta())]->Fill(energyhit);
0748           hCalo2mom2[hid.iphi()][abs(hid.ieta())]->Fill(pow(energyhit, 2));
0749         }  // eta><0
0750 
0751       }  // depth=1
0752 
0753     }  // HF_MB
0754 
0755     edm::LogVerbatim("AnalyzerMB") << " Event is finished ";
0756   }
0757 }  // namespace cms
0758 
0759 #include "FWCore/Framework/interface/MakerMacros.h"
0760 
0761 using cms::Analyzer_minbias;
0762 
0763 DEFINE_FWK_MODULE(Analyzer_minbias);