Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2025-01-31 02:20:04 

0001 /** \file GlobalHitsProdHist.cc
0002  *
0003  *  See header file for description of class
0004  *
0005  *  \author M. Strang SUNY-Buffalo
0006  */
0007 
0008 #include "FWCore/Utilities/interface/Exception.h"
0009 #include "Geometry/HcalTowerAlgo/interface/HcalGeometry.h"
0010 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0011 #include "Validation/GlobalHits/interface/GlobalHitsProdHist.h"
0012 using CLHEP::micrometer;
0013 using CLHEP::millimeter;
0014 
0015 GlobalHitsProdHist::GlobalHitsProdHist(const edm::ParameterSet &iPSet)
0016     : fName(""),
0017       verbosity(0),
0018       frequency(0),
0019       vtxunit(0),
0020       getAllProvenances(false),
0021       printProvenanceInfo(false),
0022       G4VtxSrc_Token_(consumes<edm::SimVertexContainer>((iPSet.getParameter<edm::InputTag>("G4VtxSrc")))),
0023       G4TrkSrc_Token_(consumes<edm::SimTrackContainer>(iPSet.getParameter<edm::InputTag>("G4TrkSrc"))),
0024       tGeomToken_(esConsumes()),
0025       cscGeomToken_(esConsumes()),
0026       dtGeomToken_(esConsumes()),
0027       rpcGeomToken_(esConsumes()),
0028       caloGeomToken_(esConsumes()),
0029       count(0) {
0030   std::string MsgLoggerCat = "GlobalHitsProdHist_GlobalHitsProdHist";
0031 
0032   // get information from parameter set
0033   fName = iPSet.getUntrackedParameter<std::string>("Name");
0034   verbosity = iPSet.getUntrackedParameter<int>("Verbosity");
0035   frequency = iPSet.getUntrackedParameter<int>("Frequency");
0036   vtxunit = iPSet.getUntrackedParameter<int>("VtxUnit");
0037   edm::ParameterSet m_Prov = iPSet.getParameter<edm::ParameterSet>("ProvenanceLookup");
0038   getAllProvenances = m_Prov.getUntrackedParameter<bool>("GetAllProvenances");
0039   printProvenanceInfo = m_Prov.getUntrackedParameter<bool>("PrintProvenanceInfo");
0040 
0041   // get Labels to use to extract information
0042   PxlBrlLowSrc_ = iPSet.getParameter<edm::InputTag>("PxlBrlLowSrc");
0043   PxlBrlHighSrc_ = iPSet.getParameter<edm::InputTag>("PxlBrlHighSrc");
0044   PxlFwdLowSrc_ = iPSet.getParameter<edm::InputTag>("PxlFwdLowSrc");
0045   PxlFwdHighSrc_ = iPSet.getParameter<edm::InputTag>("PxlFwdHighSrc");
0046 
0047   SiTIBLowSrc_ = iPSet.getParameter<edm::InputTag>("SiTIBLowSrc");
0048   SiTIBHighSrc_ = iPSet.getParameter<edm::InputTag>("SiTIBHighSrc");
0049   SiTOBLowSrc_ = iPSet.getParameter<edm::InputTag>("SiTOBLowSrc");
0050   SiTOBHighSrc_ = iPSet.getParameter<edm::InputTag>("SiTOBHighSrc");
0051   SiTIDLowSrc_ = iPSet.getParameter<edm::InputTag>("SiTIDLowSrc");
0052   SiTIDHighSrc_ = iPSet.getParameter<edm::InputTag>("SiTIDHighSrc");
0053   SiTECLowSrc_ = iPSet.getParameter<edm::InputTag>("SiTECLowSrc");
0054   SiTECHighSrc_ = iPSet.getParameter<edm::InputTag>("SiTECHighSrc");
0055 
0056   MuonCscSrc_ = iPSet.getParameter<edm::InputTag>("MuonCscSrc");
0057   MuonDtSrc_ = iPSet.getParameter<edm::InputTag>("MuonDtSrc");
0058   MuonRpcSrc_ = iPSet.getParameter<edm::InputTag>("MuonRpcSrc");
0059 
0060   ECalEBSrc_ = iPSet.getParameter<edm::InputTag>("ECalEBSrc");
0061   ECalEESrc_ = iPSet.getParameter<edm::InputTag>("ECalEESrc");
0062   ECalESSrc_ = iPSet.getParameter<edm::InputTag>("ECalESSrc");
0063 
0064   HCalSrc_ = iPSet.getParameter<edm::InputTag>("HCalSrc");
0065 
0066   // use value of first digit to determine default output level (inclusive)
0067   // 0 is none, 1 is basic, 2 is fill output, 3 is gather output
0068   verbosity %= 10;
0069 
0070   // print out Parameter Set information being used
0071   if (verbosity >= 0) {
0072     edm::LogInfo(MsgLoggerCat)
0073         << "\n===============================\n"
0074         << "Initialized as EDProducer with parameter values:\n"
0075         << "    Name          = " << fName << "\n"
0076         << "    Verbosity     = " << verbosity << "\n"
0077         << "    Frequency     = " << frequency << "\n"
0078         << "    VtxUnit       = " << vtxunit << "\n"
0079         << "    GetProv       = " << getAllProvenances << "\n"
0080         << "    PrintProv     = " << printProvenanceInfo << "\n"
0081         << "    PxlBrlLowSrc  = " << PxlBrlLowSrc_.label() << ":" << PxlBrlLowSrc_.instance() << "\n"
0082         << "    PxlBrlHighSrc = " << PxlBrlHighSrc_.label() << ":" << PxlBrlHighSrc_.instance() << "\n"
0083         << "    PxlFwdLowSrc  = " << PxlFwdLowSrc_.label() << ":" << PxlBrlLowSrc_.instance() << "\n"
0084         << "    PxlFwdHighSrc = " << PxlFwdHighSrc_.label() << ":" << PxlBrlHighSrc_.instance() << "\n"
0085         << "    SiTIBLowSrc   = " << SiTIBLowSrc_.label() << ":" << SiTIBLowSrc_.instance() << "\n"
0086         << "    SiTIBHighSrc  = " << SiTIBHighSrc_.label() << ":" << SiTIBHighSrc_.instance() << "\n"
0087         << "    SiTOBLowSrc   = " << SiTOBLowSrc_.label() << ":" << SiTOBLowSrc_.instance() << "\n"
0088         << "    SiTOBHighSrc  = " << SiTOBHighSrc_.label() << ":" << SiTOBHighSrc_.instance() << "\n"
0089         << "    SiTIDLowSrc   = " << SiTIDLowSrc_.label() << ":" << SiTIDLowSrc_.instance() << "\n"
0090         << "    SiTIDHighSrc  = " << SiTIDHighSrc_.label() << ":" << SiTIDHighSrc_.instance() << "\n"
0091         << "    SiTECLowSrc   = " << SiTECLowSrc_.label() << ":" << SiTECLowSrc_.instance() << "\n"
0092         << "    SiTECHighSrc  = " << SiTECHighSrc_.label() << ":" << SiTECHighSrc_.instance() << "\n"
0093         << "    MuonCscSrc    = " << MuonCscSrc_.label() << ":" << MuonCscSrc_.instance() << "\n"
0094         << "    MuonDtSrc     = " << MuonDtSrc_.label() << ":" << MuonDtSrc_.instance() << "\n"
0095         << "    MuonRpcSrc    = " << MuonRpcSrc_.label() << ":" << MuonRpcSrc_.instance() << "\n"
0096         << "    ECalEBSrc     = " << ECalEBSrc_.label() << ":" << ECalEBSrc_.instance() << "\n"
0097         << "    ECalEESrc     = " << ECalEESrc_.label() << ":" << ECalEESrc_.instance() << "\n"
0098         << "    ECalESSrc     = " << ECalESSrc_.label() << ":" << ECalESSrc_.instance() << "\n"
0099         << "    HCalSrc       = " << HCalSrc_.label() << ":" << HCalSrc_.instance() << "\n"
0100         << "===============================\n";
0101   }
0102 
0103   // create histograms
0104   Char_t hname[200];
0105   Char_t htitle[200];
0106 
0107   // MCGeant
0108   sprintf(hname, "hMCRGP1");
0109   histName_.push_back(hname);
0110   sprintf(htitle, "RawGenParticles");
0111   hMCRGP[0] = new TH1F(hname, htitle, 100, 0., 5000.);
0112   sprintf(hname, "hMCRGP2");
0113   histName_.push_back(hname);
0114   hMCRGP[1] = new TH1F(hname, htitle, 100, 0., 500.);
0115   for (Int_t i = 0; i < 2; ++i) {
0116     hMCRGP[i]->GetXaxis()->SetTitle("Number of Raw Generated Particles");
0117     hMCRGP[i]->GetYaxis()->SetTitle("Count");
0118     histMap_[hMCRGP[i]->GetName()] = hMCRGP[i];
0119   }
0120 
0121   sprintf(hname, "hMCG4Vtx1");
0122   histName_.push_back(hname);
0123   sprintf(htitle, "G4 Vertices");
0124   hMCG4Vtx[0] = new TH1F(hname, htitle, 100, 0., 50000.);
0125   sprintf(hname, "hMCG4Vtx2");
0126   histName_.push_back(hname);
0127   hMCG4Vtx[1] = new TH1F(hname, htitle, 100, -0.5, 99.5);
0128   for (Int_t i = 0; i < 2; ++i) {
0129     hMCG4Vtx[i]->GetXaxis()->SetTitle("Number of Vertices");
0130     hMCG4Vtx[i]->GetYaxis()->SetTitle("Count");
0131     histMap_[hMCG4Vtx[i]->GetName()] = hMCG4Vtx[i];
0132   }
0133 
0134   sprintf(hname, "hMCG4Trk1");
0135   histName_.push_back(hname);
0136   sprintf(htitle, "G4 Tracks");
0137   hMCG4Trk[0] = new TH1F(hname, htitle, 150, 0., 15000.);
0138   sprintf(hname, "hMCG4Trk2");
0139   histName_.push_back(hname);
0140   hMCG4Trk[1] = new TH1F(hname, htitle, 150, -0.5, 99.5);
0141   for (Int_t i = 0; i < 2; ++i) {
0142     hMCG4Trk[i]->GetXaxis()->SetTitle("Number of Tracks");
0143     hMCG4Trk[i]->GetYaxis()->SetTitle("Count");
0144     histMap_[hMCG4Trk[i]->GetName()] = hMCG4Trk[i];
0145   }
0146 
0147   sprintf(hname, "hGeantVtxX1");
0148   histName_.push_back(hname);
0149   sprintf(htitle, "Geant vertex x/micrometer");
0150   hGeantVtxX[0] = new TH1F(hname, htitle, 100, -8000000., 8000000.);
0151   sprintf(hname, "hGeantVtxX2");
0152   histName_.push_back(hname);
0153   hGeantVtxX[1] = new TH1F(hname, htitle, 100, -50., 50.);
0154   for (Int_t i = 0; i < 2; ++i) {
0155     hGeantVtxX[i]->GetXaxis()->SetTitle("x of Vertex (um)");
0156     hGeantVtxX[i]->GetYaxis()->SetTitle("Count");
0157     histMap_[hGeantVtxX[i]->GetName()] = hGeantVtxX[i];
0158   }
0159 
0160   sprintf(hname, "hGeantVtxY1");
0161   histName_.push_back(hname);
0162   sprintf(htitle, "Geant vertex y/micrometer");
0163   hGeantVtxY[0] = new TH1F(hname, htitle, 100, -8000000, 8000000.);
0164   sprintf(hname, "hGeantVtxY2");
0165   histName_.push_back(hname);
0166   hGeantVtxY[1] = new TH1F(hname, htitle, 100, -50., 50.);
0167   for (Int_t i = 0; i < 2; ++i) {
0168     hGeantVtxY[i]->GetXaxis()->SetTitle("y of Vertex (um)");
0169     hGeantVtxY[i]->GetYaxis()->SetTitle("Count");
0170     histMap_[hGeantVtxY[i]->GetName()] = hGeantVtxY[i];
0171   }
0172 
0173   sprintf(hname, "hGeantVtxZ1");
0174   histName_.push_back(hname);
0175   sprintf(htitle, "Geant vertex z/millimeter");
0176   hGeantVtxZ[0] = new TH1F(hname, htitle, 100, -11000., 11000.);
0177   sprintf(hname, "hGeantVtxZ2");
0178   histName_.push_back(hname);
0179   hGeantVtxZ[1] = new TH1F(hname, htitle, 100, -250., 250.);
0180   for (Int_t i = 0; i < 2; ++i) {
0181     hGeantVtxZ[i]->GetXaxis()->SetTitle("z of Vertex (mm)");
0182     hGeantVtxZ[i]->GetYaxis()->SetTitle("Count");
0183     histMap_[hGeantVtxZ[i]->GetName()] = hGeantVtxZ[i];
0184   }
0185 
0186   sprintf(hname, "hGeantTrkPt");
0187   histName_.push_back(hname);
0188   sprintf(htitle, "Geant track pt/GeV");
0189   hGeantTrkPt = new TH1F(hname, htitle, 100, 0., 200.);
0190   hGeantTrkPt->GetXaxis()->SetTitle("pT of Track (GeV)");
0191   hGeantTrkPt->GetYaxis()->SetTitle("Count");
0192   histMap_[hGeantTrkPt->GetName()] = hGeantTrkPt;
0193 
0194   sprintf(hname, "hGeantTrkE");
0195   histName_.push_back(hname);
0196   sprintf(htitle, "Geant track E/GeV");
0197   hGeantTrkE = new TH1F(hname, htitle, 100, 0., 5000.);
0198   hGeantTrkE->GetXaxis()->SetTitle("E of Track (GeV)");
0199   hGeantTrkE->GetYaxis()->SetTitle("Count");
0200   histMap_[hGeantTrkE->GetName()] = hGeantTrkE;
0201 
0202   // ECal
0203   sprintf(hname, "hCaloEcal1");
0204   histName_.push_back(hname);
0205   sprintf(htitle, "Ecal hits");
0206   hCaloEcal[0] = new TH1F(hname, htitle, 100, 0., 10000.);
0207   sprintf(hname, "hCaloEcal2");
0208   histName_.push_back(hname);
0209   hCaloEcal[1] = new TH1F(hname, htitle, 100, -0.5, 99.5);
0210 
0211   sprintf(hname, "hCaloEcalE1");
0212   histName_.push_back(hname);
0213   sprintf(htitle, "Ecal hits, energy/GeV");
0214   hCaloEcalE[0] = new TH1F(hname, htitle, 100, 0., 10.);
0215   sprintf(hname, "hCaloEcalE2");
0216   histName_.push_back(hname);
0217   hCaloEcalE[1] = new TH1F(hname, htitle, 100, 0., 0.1);
0218 
0219   sprintf(hname, "hCaloEcalToF1");
0220   histName_.push_back(hname);
0221   sprintf(htitle, "Ecal hits, ToF/ns");
0222   hCaloEcalToF[0] = new TH1F(hname, htitle, 100, 0., 1000.);
0223   sprintf(hname, "hCaloEcalToF2");
0224   histName_.push_back(hname);
0225   hCaloEcalToF[1] = new TH1F(hname, htitle, 100, 0., 100.);
0226 
0227   for (Int_t i = 0; i < 2; ++i) {
0228     hCaloEcal[i]->GetXaxis()->SetTitle("Number of Hits");
0229     hCaloEcal[i]->GetYaxis()->SetTitle("Count");
0230     histMap_[hCaloEcal[i]->GetName()] = hCaloEcal[i];
0231     hCaloEcalE[i]->GetXaxis()->SetTitle("Energy of Hits (GeV)");
0232     hCaloEcalE[i]->GetYaxis()->SetTitle("Count");
0233     histMap_[hCaloEcalE[i]->GetName()] = hCaloEcalE[i];
0234     hCaloEcalToF[i]->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0235     hCaloEcalToF[i]->GetYaxis()->SetTitle("Count");
0236     histMap_[hCaloEcalToF[i]->GetName()] = hCaloEcalToF[i];
0237   }
0238 
0239   sprintf(hname, "hCaloEcalPhi");
0240   histName_.push_back(hname);
0241   sprintf(htitle, "Ecal hits, phi/rad");
0242   hCaloEcalPhi = new TH1F(hname, htitle, 100, -3.2, 3.2);
0243   hCaloEcalPhi->GetXaxis()->SetTitle("Phi of Hits (rad)");
0244   hCaloEcalPhi->GetYaxis()->SetTitle("Count");
0245   histMap_[hCaloEcalPhi->GetName()] = hCaloEcalPhi;
0246 
0247   sprintf(hname, "hCaloEcalEta");
0248   histName_.push_back(hname);
0249   sprintf(htitle, "Ecal hits, eta");
0250   hCaloEcalEta = new TH1F(hname, htitle, 100, -5.5, 5.5);
0251   hCaloEcalEta->GetXaxis()->SetTitle("Eta of Hits");
0252   hCaloEcalEta->GetYaxis()->SetTitle("Count");
0253   histMap_[hCaloEcalEta->GetName()] = hCaloEcalEta;
0254 
0255   sprintf(hname, "hCaloPreSh1");
0256   histName_.push_back(hname);
0257   sprintf(htitle, "PreSh hits");
0258   hCaloPreSh[0] = new TH1F(hname, htitle, 100, 0., 10000.);
0259   sprintf(hname, "hCaloPreSh2");
0260   histName_.push_back(hname);
0261   hCaloPreSh[1] = new TH1F(hname, htitle, 100, -0.5, 99.5);
0262 
0263   sprintf(hname, "hCaloPreShE1");
0264   histName_.push_back(hname);
0265   sprintf(htitle, "PreSh hits, energy/GeV");
0266   hCaloPreShE[0] = new TH1F(hname, htitle, 100, 0., 10.);
0267   sprintf(hname, "hCaloPreShE2");
0268   histName_.push_back(hname);
0269   hCaloPreShE[1] = new TH1F(hname, htitle, 100, 0., 0.1);
0270 
0271   sprintf(hname, "hCaloPreShToF1");
0272   histName_.push_back(hname);
0273   sprintf(htitle, "PreSh hits, ToF/ns");
0274   hCaloPreShToF[0] = new TH1F(hname, htitle, 100, 0., 1000.);
0275   sprintf(hname, "hCaloPreShToF2");
0276   histName_.push_back(hname);
0277   hCaloPreShToF[1] = new TH1F(hname, htitle, 100, 0., 100.);
0278 
0279   for (Int_t i = 0; i < 2; ++i) {
0280     hCaloPreSh[i]->GetXaxis()->SetTitle("Number of Hits");
0281     hCaloPreSh[i]->GetYaxis()->SetTitle("Count");
0282     histMap_[hCaloPreSh[i]->GetName()] = hCaloPreSh[i];
0283     hCaloPreShE[i]->GetXaxis()->SetTitle("Energy of Hits (GeV)");
0284     hCaloPreShE[i]->GetYaxis()->SetTitle("Count");
0285     histMap_[hCaloPreShE[i]->GetName()] = hCaloPreShE[i];
0286     hCaloPreShToF[i]->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0287     hCaloPreShToF[i]->GetYaxis()->SetTitle("Count");
0288     histMap_[hCaloPreShToF[i]->GetName()] = hCaloPreShToF[i];
0289   }
0290 
0291   sprintf(hname, "hCaloPreShPhi");
0292   histName_.push_back(hname);
0293   sprintf(htitle, "PreSh hits, phi/rad");
0294   hCaloPreShPhi = new TH1F(hname, htitle, 100, -3.2, 3.2);
0295   hCaloPreShPhi->GetXaxis()->SetTitle("Phi of Hits (rad)");
0296   hCaloPreShPhi->GetYaxis()->SetTitle("Count");
0297   histMap_[hCaloPreShPhi->GetName()] = hCaloPreShPhi;
0298 
0299   sprintf(hname, "hCaloPreShEta");
0300   histName_.push_back(hname);
0301   sprintf(htitle, "PreSh hits, eta");
0302   hCaloPreShEta = new TH1F(hname, htitle, 100, -5.5, 5.5);
0303   hCaloPreShEta->GetXaxis()->SetTitle("Eta of Hits");
0304   hCaloPreShEta->GetYaxis()->SetTitle("Count");
0305   histMap_[hCaloPreShEta->GetName()] = hCaloPreShEta;
0306 
0307   // Hcal
0308   sprintf(hname, "hCaloHcal1");
0309   histName_.push_back(hname);
0310   sprintf(htitle, "Hcal hits");
0311   hCaloHcal[0] = new TH1F(hname, htitle, 100, 0., 10000.);
0312   sprintf(hname, "hCaloHcal2");
0313   histName_.push_back(hname);
0314   hCaloHcal[1] = new TH1F(hname, htitle, 100, -0.5, 99.5);
0315 
0316   sprintf(hname, "hCaloHcalE1");
0317   histName_.push_back(hname);
0318   sprintf(htitle, "Hcal hits, energy/GeV");
0319   hCaloHcalE[0] = new TH1F(hname, htitle, 100, 0., 10.);
0320   sprintf(hname, "hCaloHcalE2");
0321   histName_.push_back(hname);
0322   hCaloHcalE[1] = new TH1F(hname, htitle, 100, 0., 0.1);
0323 
0324   sprintf(hname, "hCaloHcalToF1");
0325   histName_.push_back(hname);
0326   sprintf(htitle, "Hcal hits, ToF/ns");
0327   hCaloHcalToF[0] = new TH1F(hname, htitle, 100, 0., 1000.);
0328   sprintf(hname, "hCaloHcalToF2");
0329   histName_.push_back(hname);
0330   hCaloHcalToF[1] = new TH1F(hname, htitle, 100, 0., 100.);
0331 
0332   for (Int_t i = 0; i < 2; ++i) {
0333     hCaloHcal[i]->GetXaxis()->SetTitle("Number of Hits");
0334     hCaloHcal[i]->GetYaxis()->SetTitle("Count");
0335     histMap_[hCaloHcal[i]->GetName()] = hCaloHcal[i];
0336     hCaloHcalE[i]->GetXaxis()->SetTitle("Energy of Hits (GeV)");
0337     hCaloHcalE[i]->GetYaxis()->SetTitle("Count");
0338     histMap_[hCaloHcalE[i]->GetName()] = hCaloHcalE[i];
0339     hCaloHcalToF[i]->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0340     hCaloHcalToF[i]->GetYaxis()->SetTitle("Count");
0341     histMap_[hCaloHcalToF[i]->GetName()] = hCaloHcalToF[i];
0342   }
0343 
0344   sprintf(hname, "hCaloHcalPhi");
0345   histName_.push_back(hname);
0346   sprintf(htitle, "Hcal hits, phi/rad");
0347   hCaloHcalPhi = new TH1F(hname, htitle, 100, -3.2, 3.2);
0348   hCaloHcalPhi->GetXaxis()->SetTitle("Phi of Hits (rad)");
0349   hCaloHcalPhi->GetYaxis()->SetTitle("Count");
0350   histMap_[hCaloHcalPhi->GetName()] = hCaloHcalPhi;
0351 
0352   sprintf(hname, "hCaloHcalEta");
0353   histName_.push_back(hname);
0354   sprintf(htitle, "Hcal hits, eta");
0355   hCaloHcalEta = new TH1F(hname, htitle, 100, -5.5, 5.5);
0356   hCaloHcalEta->GetXaxis()->SetTitle("Eta of Hits");
0357   hCaloHcalEta->GetYaxis()->SetTitle("Count");
0358   histMap_[hCaloHcalEta->GetName()] = hCaloHcalEta;
0359 
0360   // tracker
0361   sprintf(hname, "hTrackerPx1");
0362   histName_.push_back(hname);
0363   sprintf(htitle, "Pixel hits");
0364   hTrackerPx[0] = new TH1F(hname, htitle, 100, 0., 10000.);
0365   sprintf(hname, "hTrackerPx2");
0366   histName_.push_back(hname);
0367   hTrackerPx[1] = new TH1F(hname, htitle, 100, -0.5, 99.5);
0368   for (Int_t i = 0; i < 2; ++i) {
0369     hTrackerPx[i]->GetXaxis()->SetTitle("Number of Pixel Hits");
0370     hTrackerPx[i]->GetYaxis()->SetTitle("Count");
0371     histMap_[hTrackerPx[i]->GetName()] = hTrackerPx[i];
0372   }
0373 
0374   sprintf(hname, "hTrackerPxPhi");
0375   histName_.push_back(hname);
0376   sprintf(htitle, "Pixel hits phi/rad");
0377   hTrackerPxPhi = new TH1F(hname, htitle, 100, -3.2, 3.2);
0378   hTrackerPxPhi->GetXaxis()->SetTitle("Phi of Hits (rad)");
0379   hTrackerPxPhi->GetYaxis()->SetTitle("Count");
0380   histMap_[hTrackerPxPhi->GetName()] = hTrackerPxPhi;
0381 
0382   sprintf(hname, "hTrackerPxEta");
0383   histName_.push_back(hname);
0384   sprintf(htitle, "Pixel hits eta");
0385   hTrackerPxEta = new TH1F(hname, htitle, 100, -3.5, 3.5);
0386   hTrackerPxEta->GetXaxis()->SetTitle("Eta of Hits");
0387   hTrackerPxEta->GetYaxis()->SetTitle("Count");
0388   histMap_[hTrackerPxEta->GetName()] = hTrackerPxEta;
0389 
0390   sprintf(hname, "hTrackerPxBToF");
0391   histName_.push_back(hname);
0392   sprintf(htitle, "Pixel barrel hits, ToF/ns");
0393   hTrackerPxBToF = new TH1F(hname, htitle, 100, 0., 40.);
0394   hTrackerPxBToF->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0395   hTrackerPxBToF->GetYaxis()->SetTitle("Count");
0396   histMap_[hTrackerPxBToF->GetName()] = hTrackerPxBToF;
0397 
0398   sprintf(hname, "hTrackerPxBR");
0399   histName_.push_back(hname);
0400   sprintf(htitle, "Pixel barrel hits, R/cm");
0401   hTrackerPxBR = new TH1F(hname, htitle, 100, 0., 50.);
0402   hTrackerPxBR->GetXaxis()->SetTitle("R of Hits (cm)");
0403   hTrackerPxBR->GetYaxis()->SetTitle("Count");
0404   histMap_[hTrackerPxBR->GetName()] = hTrackerPxBR;
0405 
0406   sprintf(hname, "hTrackerPxFToF");
0407   histName_.push_back(hname);
0408   sprintf(htitle, "Pixel forward hits, ToF/ns");
0409   hTrackerPxFToF = new TH1F(hname, htitle, 100, 0., 50.);
0410   hTrackerPxFToF->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0411   hTrackerPxFToF->GetYaxis()->SetTitle("Count");
0412   histMap_[hTrackerPxFToF->GetName()] = hTrackerPxFToF;
0413 
0414   sprintf(hname, "hTrackerPxFZ");
0415   histName_.push_back(hname);
0416   sprintf(htitle, "Pixel forward hits, Z/cm");
0417   hTrackerPxFZ = new TH1F(hname, htitle, 200, -100., 100.);
0418   hTrackerPxFZ->GetXaxis()->SetTitle("Z of Hits (cm)");
0419   hTrackerPxFZ->GetYaxis()->SetTitle("Count");
0420   histMap_[hTrackerPxFZ->GetName()] = hTrackerPxFZ;
0421 
0422   sprintf(hname, "hTrackerSi1");
0423   histName_.push_back(hname);
0424   sprintf(htitle, "Silicon hits");
0425   hTrackerSi[0] = new TH1F(hname, htitle, 100, 0., 10000.);
0426   sprintf(hname, "hTrackerSi2");
0427   histName_.push_back(hname);
0428   hTrackerSi[1] = new TH1F(hname, htitle, 100, -0.5, 99.5);
0429   for (Int_t i = 0; i < 2; ++i) {
0430     hTrackerSi[i]->GetXaxis()->SetTitle("Number of Silicon Hits");
0431     hTrackerSi[i]->GetYaxis()->SetTitle("Count");
0432     histMap_[hTrackerSi[i]->GetName()] = hTrackerSi[i];
0433   }
0434 
0435   sprintf(hname, "hTrackerSiPhi");
0436   histName_.push_back(hname);
0437   sprintf(htitle, "Silicon hits phi/rad");
0438   hTrackerSiPhi = new TH1F(hname, htitle, 100, -3.2, 3.2);
0439   hTrackerSiPhi->GetXaxis()->SetTitle("Phi of Hits (rad)");
0440   hTrackerSiPhi->GetYaxis()->SetTitle("Count");
0441   histMap_[hTrackerSiPhi->GetName()] = hTrackerSiPhi;
0442 
0443   sprintf(hname, "hTrackerSiEta");
0444   histName_.push_back(hname);
0445   sprintf(htitle, "Silicon hits eta");
0446   hTrackerSiEta = new TH1F(hname, htitle, 100, -3.5, 3.5);
0447   hTrackerSiEta->GetXaxis()->SetTitle("Eta of Hits");
0448   hTrackerSiEta->GetYaxis()->SetTitle("Count");
0449   histMap_[hTrackerSiEta->GetName()] = hTrackerSiEta;
0450 
0451   sprintf(hname, "hTrackerSiBToF");
0452   histName_.push_back(hname);
0453   sprintf(htitle, "Silicon barrel hits, ToF/ns");
0454   hTrackerSiBToF = new TH1F(hname, htitle, 100, 0., 50.);
0455   hTrackerSiBToF->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0456   hTrackerSiBToF->GetYaxis()->SetTitle("Count");
0457   histMap_[hTrackerSiBToF->GetName()] = hTrackerSiBToF;
0458 
0459   sprintf(hname, "hTrackerSiBR");
0460   histName_.push_back(hname);
0461   sprintf(htitle, "Silicon barrel hits, R/cm");
0462   hTrackerSiBR = new TH1F(hname, htitle, 100, 0., 200.);
0463   hTrackerSiBR->GetXaxis()->SetTitle("R of Hits (cm)");
0464   hTrackerSiBR->GetYaxis()->SetTitle("Count");
0465   histMap_[hTrackerSiBR->GetName()] = hTrackerSiBR;
0466 
0467   sprintf(hname, "hTrackerSiFToF");
0468   histName_.push_back(hname);
0469   sprintf(htitle, "Silicon forward hits, ToF/ns");
0470   hTrackerSiFToF = new TH1F(hname, htitle, 100, 0., 75.);
0471   hTrackerSiFToF->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0472   hTrackerSiFToF->GetYaxis()->SetTitle("Count");
0473   histMap_[hTrackerSiFToF->GetName()] = hTrackerSiFToF;
0474 
0475   sprintf(hname, "hTrackerSiFZ");
0476   histName_.push_back(hname);
0477   sprintf(htitle, "Silicon forward hits, Z/cm");
0478   hTrackerSiFZ = new TH1F(hname, htitle, 200, -300., 300.);
0479   hTrackerSiFZ->GetXaxis()->SetTitle("Z of Hits (cm)");
0480   hTrackerSiFZ->GetYaxis()->SetTitle("Count");
0481   histMap_[hTrackerSiFZ->GetName()] = hTrackerSiFZ;
0482 
0483   // muon
0484   sprintf(hname, "hMuon1");
0485   histName_.push_back(hname);
0486   sprintf(htitle, "Muon hits");
0487   hMuon[0] = new TH1F(hname, htitle, 100, 0., 10000.);
0488   sprintf(hname, "hMuon2");
0489   histName_.push_back(hname);
0490   hMuon[1] = new TH1F(hname, htitle, 100, -0.5, 99.5);
0491   for (Int_t i = 0; i < 2; ++i) {
0492     hMuon[i]->GetXaxis()->SetTitle("Number of Muon Hits");
0493     hMuon[i]->GetYaxis()->SetTitle("Count");
0494     histMap_[hMuon[i]->GetName()] = hMuon[i];
0495   }
0496 
0497   sprintf(hname, "hMuonPhi");
0498   histName_.push_back(hname);
0499   sprintf(htitle, "Muon hits phi/rad");
0500   hMuonPhi = new TH1F(hname, htitle, 100, -3.2, 3.2);
0501   hMuonPhi->GetXaxis()->SetTitle("Phi of Hits (rad)");
0502   hMuonPhi->GetYaxis()->SetTitle("Count");
0503   histMap_[hMuonPhi->GetName()] = hMuonPhi;
0504 
0505   sprintf(hname, "hMuonEta");
0506   histName_.push_back(hname);
0507   sprintf(htitle, "Muon hits eta");
0508   hMuonEta = new TH1F(hname, htitle, 100, -3.5, 3.5);
0509   hMuonEta->GetXaxis()->SetTitle("Eta of Hits");
0510   hMuonEta->GetYaxis()->SetTitle("Count");
0511   histMap_[hMuonEta->GetName()] = hMuonEta;
0512 
0513   sprintf(hname, "hMuonCscToF1");
0514   histName_.push_back(hname);
0515   sprintf(htitle, "Muon CSC hits, ToF/ns");
0516   hMuonCscToF[0] = new TH1F(hname, htitle, 100, 0., 250.);
0517   sprintf(hname, "hMuonCscToF2");
0518   histName_.push_back(hname);
0519   hMuonCscToF[1] = new TH1F(hname, htitle, 100, 0., 50.);
0520   for (Int_t i = 0; i < 2; ++i) {
0521     hMuonCscToF[i]->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0522     hMuonCscToF[i]->GetYaxis()->SetTitle("Count");
0523     histMap_[hMuonCscToF[i]->GetName()] = hMuonCscToF[i];
0524   }
0525 
0526   sprintf(hname, "hMuonCscZ");
0527   histName_.push_back(hname);
0528   sprintf(htitle, "Muon CSC hits, Z/cm");
0529   hMuonCscZ = new TH1F(hname, htitle, 200, -1500., 1500.);
0530   hMuonCscZ->GetXaxis()->SetTitle("Z of Hits (cm)");
0531   hMuonCscZ->GetYaxis()->SetTitle("Count");
0532   histMap_[hMuonCscZ->GetName()] = hMuonCscZ;
0533 
0534   sprintf(hname, "hMuonDtToF1");
0535   histName_.push_back(hname);
0536   sprintf(htitle, "Muon DT hits, ToF/ns");
0537   hMuonDtToF[0] = new TH1F(hname, htitle, 100, 0., 250.);
0538   sprintf(hname, "hMuonDtToF2");
0539   histName_.push_back(hname);
0540   hMuonDtToF[1] = new TH1F(hname, htitle, 100, 0., 50.);
0541   for (Int_t i = 0; i < 2; ++i) {
0542     hMuonDtToF[i]->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0543     hMuonDtToF[i]->GetYaxis()->SetTitle("Count");
0544     histMap_[hMuonDtToF[i]->GetName()] = hMuonDtToF[i];
0545   }
0546 
0547   sprintf(hname, "hMuonDtR");
0548   histName_.push_back(hname);
0549   sprintf(htitle, "Muon DT hits, R/cm");
0550   hMuonDtR = new TH1F(hname, htitle, 100, 0., 1500.);
0551   hMuonDtR->GetXaxis()->SetTitle("R of Hits (cm)");
0552   hMuonDtR->GetYaxis()->SetTitle("Count");
0553   histMap_[hMuonDtR->GetName()] = hMuonDtR;
0554 
0555   sprintf(hname, "hMuonRpcFToF1");
0556   histName_.push_back(hname);
0557   sprintf(htitle, "Muon RPC forward hits, ToF/ns");
0558   hMuonRpcFToF[0] = new TH1F(hname, htitle, 100, 0., 250.);
0559   sprintf(hname, "hMuonRpcFToF2");
0560   histName_.push_back(hname);
0561   hMuonRpcFToF[1] = new TH1F(hname, htitle, 100, 0., 50.);
0562   for (Int_t i = 0; i < 2; ++i) {
0563     hMuonRpcFToF[i]->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0564     hMuonRpcFToF[i]->GetYaxis()->SetTitle("Count");
0565     histMap_[hMuonRpcFToF[i]->GetName()] = hMuonRpcFToF[i];
0566   }
0567 
0568   sprintf(hname, "hMuonRpcFZ");
0569   histName_.push_back(hname);
0570   sprintf(htitle, "Muon RPC forward hits, Z/cm");
0571   hMuonRpcFZ = new TH1F(hname, htitle, 201, -1500., 1500.);
0572   hMuonRpcFZ->GetXaxis()->SetTitle("Z of Hits (cm)");
0573   hMuonRpcFZ->GetYaxis()->SetTitle("Count");
0574   histMap_[hMuonRpcFZ->GetName()] = hMuonRpcFZ;
0575 
0576   sprintf(hname, "hMuonRpcBToF1");
0577   histName_.push_back(hname);
0578   sprintf(htitle, "Muon RPC barrel hits, ToF/ns");
0579   hMuonRpcBToF[0] = new TH1F(hname, htitle, 100, 0., 250.);
0580   sprintf(hname, "hMuonRpcBToF2");
0581   histName_.push_back(hname);
0582   hMuonRpcBToF[1] = new TH1F(hname, htitle, 100, 0., 50.);
0583   for (Int_t i = 0; i < 2; ++i) {
0584     hMuonRpcBToF[i]->GetXaxis()->SetTitle("Time of Flight of Hits (ns)");
0585     hMuonRpcBToF[i]->GetYaxis()->SetTitle("Count");
0586     histMap_[hMuonRpcBToF[i]->GetName()] = hMuonRpcBToF[i];
0587   }
0588 
0589   sprintf(hname, "hMuonRpcBR");
0590   histName_.push_back(hname);
0591   sprintf(htitle, "Muon RPC barrel hits, R/cm");
0592   hMuonRpcBR = new TH1F(hname, htitle, 100, 0., 1500.);
0593   hMuonRpcBR->GetXaxis()->SetTitle("R of Hits (cm)");
0594   hMuonRpcBR->GetYaxis()->SetTitle("Count");
0595   histMap_[hMuonRpcBR->GetName()] = hMuonRpcBR;
0596 
0597   // create persistent objects
0598   for (std::size_t i = 0; i < histName_.size(); ++i) {
0599     produces<TH1F, edm::Transition::EndRun>(histName_[i]).setBranchAlias(histName_[i]);
0600   }
0601 }
0602 
0603 GlobalHitsProdHist::~GlobalHitsProdHist() {}
0604 
0605 void GlobalHitsProdHist::beginJob() { return; }
0606 
0607 void GlobalHitsProdHist::endJob() {
0608   std::string MsgLoggerCat = "GlobalHitsProdHist_endJob";
0609   if (verbosity >= 0)
0610     edm::LogInfo(MsgLoggerCat) << "Terminating having processed " << count << " events.";
0611   return;
0612 }
0613 
0614 void GlobalHitsProdHist::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0615   std::string MsgLoggerCat = "GlobalHitsProdHist_produce";
0616 
0617   // keep track of number of events processed
0618   ++count;
0619 
0620   // get event id information
0621   edm::RunNumber_t nrun = iEvent.id().run();
0622   edm::EventNumber_t nevt = iEvent.id().event();
0623 
0624   if (verbosity > 0) {
0625     edm::LogInfo(MsgLoggerCat) << "Processing run " << nrun << ", event " << nevt << " (" << count << " events total)";
0626   } else if (verbosity == 0) {
0627     if (nevt % frequency == 0 || nevt == 1) {
0628       edm::LogInfo(MsgLoggerCat) << "Processing run " << nrun << ", event " << nevt << " (" << count
0629                                  << " events total)";
0630     }
0631   }
0632 
0633   // look at information available in the event
0634   if (getAllProvenances) {
0635     std::vector<const edm::StableProvenance *> AllProv;
0636     iEvent.getAllStableProvenance(AllProv);
0637 
0638     if (verbosity >= 0)
0639       edm::LogInfo(MsgLoggerCat) << "Number of Provenances = " << AllProv.size();
0640 
0641     if (printProvenanceInfo && (verbosity >= 0)) {
0642       TString eventout("\nProvenance info:\n");
0643 
0644       for (unsigned int i = 0; i < AllProv.size(); ++i) {
0645         eventout += "\n       ******************************";
0646         eventout += "\n       Module       : ";
0647         eventout += AllProv[i]->moduleLabel();
0648         eventout += "\n       ProductID    : ";
0649         eventout += AllProv[i]->productID().id();
0650         eventout += "\n       ClassName    : ";
0651         eventout += AllProv[i]->className();
0652         eventout += "\n       InstanceName : ";
0653         eventout += AllProv[i]->productInstanceName();
0654         eventout += "\n       BranchName   : ";
0655         eventout += AllProv[i]->branchName();
0656       }
0657       eventout += "\n       ******************************\n";
0658       edm::LogInfo(MsgLoggerCat) << eventout << "\n";
0659       printProvenanceInfo = false;
0660     }
0661     getAllProvenances = false;
0662   }
0663 
0664   // call fill functions
0665   // gather G4MC information from event
0666   fillG4MC(iEvent);
0667   // gather Tracker information from event
0668   fillTrk(iEvent, iSetup);
0669   // gather muon information from event
0670   fillMuon(iEvent, iSetup);
0671   // gather Ecal information from event
0672   fillECal(iEvent, iSetup);
0673   // gather Hcal information from event
0674   fillHCal(iEvent, iSetup);
0675 
0676   if (verbosity > 0)
0677     edm::LogInfo(MsgLoggerCat) << "Done gathering data from event.";
0678 
0679   return;
0680 }
0681 
0682 void GlobalHitsProdHist::endRunProduce(edm::Run &iRun, const edm::EventSetup &iSetup) {
0683   std::string MsgLoggerCat = "GlobalHitsProdHist_endRun";
0684 
0685   TString eventout;
0686   TString eventoutw;
0687   bool warning = false;
0688 
0689   if (verbosity > 0)
0690     edm::LogInfo(MsgLoggerCat) << "\nStoring histograms.";
0691 
0692   // store persistent objects
0693   std::map<std::string, TH1F *>::iterator iter;
0694   for (std::size_t i = 0; i < histName_.size(); ++i) {
0695     iter = histMap_.find(histName_[i]);
0696     if (iter != histMap_.end()) {
0697       std::unique_ptr<TH1F> hist1D(iter->second);
0698       eventout += "\n Storing histogram " + histName_[i];
0699       iRun.put(std::move(hist1D), histName_[i]);
0700     } else {
0701       warning = true;
0702       eventoutw += "\n Unable to find histogram with name " + histName_[i];
0703     }
0704   }
0705 
0706   if (verbosity > 0) {
0707     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
0708     if (warning)
0709       edm::LogWarning(MsgLoggerCat) << eventoutw << "\n";
0710   }
0711   return;
0712 }
0713 
0714 //==================fill and store functions================================
0715 void GlobalHitsProdHist::fillG4MC(edm::Event &iEvent) {
0716   std::string MsgLoggerCat = "GlobalHitsProdHist_fillG4MC";
0717 
0718   TString eventout;
0719   if (verbosity > 0)
0720     eventout = "\nGathering info:";
0721 
0722   //////////////////////
0723   // get MC information
0724   /////////////////////
0725   edm::Handle<edm::HepMCProduct> HepMCEvt;
0726   std::vector<edm::Handle<edm::HepMCProduct>> AllHepMCEvt;
0727 
0728   //iEvent.getManyByType(AllHepMCEvt);
0729   throw cms::Exception("UnsupportedFunction") << "GlobalHitsProdHist::fillG4MC: "
0730                                               << "getManyByType has not been supported by the Framework since 2015. "
0731                                               << "This module has been broken since then. Maybe it should be deleted. "
0732                                               << "Another possibility is to upgrade to use GetterOfProducts instead.";
0733 
0734   // loop through products and extract VtxSmearing if available. Any of them
0735   // should have the information needed
0736   for (unsigned int i = 0; i < AllHepMCEvt.size(); ++i) {
0737     HepMCEvt = AllHepMCEvt[i];
0738     if ((HepMCEvt.provenance()->productDescription()).moduleLabel() == "generatorSmeared")
0739       break;
0740   }
0741 
0742   if (!HepMCEvt.isValid()) {
0743     edm::LogWarning(MsgLoggerCat) << "Unable to find HepMCProduct in event!";
0744     return;
0745   } else {
0746     eventout += "\n          Using HepMCProduct: ";
0747     eventout += (HepMCEvt.provenance()->productDescription()).moduleLabel();
0748   }
0749   const HepMC::GenEvent *MCEvt = HepMCEvt->GetEvent();
0750   nRawGenPart = MCEvt->particles_size();
0751 
0752   if (verbosity > 1) {
0753     eventout += "\n          Number of Raw Particles collected:......... ";
0754     eventout += nRawGenPart;
0755   }
0756 
0757   if (hMCRGP[0])
0758     hMCRGP[0]->Fill((float)nRawGenPart);
0759   if (hMCRGP[1])
0760     hMCRGP[1]->Fill((float)nRawGenPart);
0761 
0762   ////////////////////////////
0763   // get G4Vertex information
0764   ////////////////////////////
0765   // convert unit stored in SimVertex to mm
0766   float unit = 0.;
0767   if (vtxunit == 0)
0768     unit = 1.;  // already in mm
0769   if (vtxunit == 1)
0770     unit = 10.;  // stored in cm, convert to mm
0771 
0772   edm::Handle<edm::SimVertexContainer> G4VtxContainer;
0773   iEvent.getByToken(G4VtxSrc_Token_, G4VtxContainer);
0774   if (!G4VtxContainer.isValid()) {
0775     edm::LogWarning(MsgLoggerCat) << "Unable to find SimVertex in event!";
0776     return;
0777   }
0778   int i = 0;
0779   edm::SimVertexContainer::const_iterator itVtx;
0780   for (itVtx = G4VtxContainer->begin(); itVtx != G4VtxContainer->end(); ++itVtx) {
0781     ++i;
0782 
0783     const math::XYZTLorentzVector G4Vtx1(
0784         itVtx->position().x(), itVtx->position().y(), itVtx->position().z(), itVtx->position().e());
0785 
0786     double G4Vtx[4];
0787     G4Vtx1.GetCoordinates(G4Vtx);
0788 
0789     if (hGeantVtxX[0])
0790       hGeantVtxX[0]->Fill((G4Vtx[0] * unit) / micrometer);
0791     if (hGeantVtxX[1])
0792       hGeantVtxX[1]->Fill((G4Vtx[0] * unit) / micrometer);
0793 
0794     if (hGeantVtxY[0])
0795       hGeantVtxY[0]->Fill((G4Vtx[1] * unit) / micrometer);
0796     if (hGeantVtxY[1])
0797       hGeantVtxY[1]->Fill((G4Vtx[1] * unit) / micrometer);
0798 
0799     if (hGeantVtxZ[0])
0800       hGeantVtxZ[0]->Fill((G4Vtx[2] * unit) / millimeter);
0801     if (hGeantVtxZ[1])
0802       hGeantVtxZ[1]->Fill((G4Vtx[2] * unit) / millimeter);
0803   }
0804 
0805   if (verbosity > 1) {
0806     eventout += "\n          Number of G4Vertices collected:............ ";
0807     eventout += i;
0808   }
0809 
0810   if (hMCG4Vtx[0])
0811     hMCG4Vtx[0]->Fill((float)i);
0812   if (hMCG4Vtx[1])
0813     hMCG4Vtx[1]->Fill((float)i);
0814 
0815   ///////////////////////////
0816   // get G4Track information
0817   ///////////////////////////
0818   edm::Handle<edm::SimTrackContainer> G4TrkContainer;
0819   iEvent.getByToken(G4TrkSrc_Token_, G4TrkContainer);
0820   if (!G4TrkContainer.isValid()) {
0821     edm::LogWarning(MsgLoggerCat) << "Unable to find SimTrack in event!";
0822     return;
0823   }
0824   i = 0;
0825   edm::SimTrackContainer::const_iterator itTrk;
0826   for (itTrk = G4TrkContainer->begin(); itTrk != G4TrkContainer->end(); ++itTrk) {
0827     ++i;
0828 
0829     const math::XYZTLorentzVector G4Trk1(
0830         itTrk->momentum().x(), itTrk->momentum().y(), itTrk->momentum().z(), itTrk->momentum().e());
0831     double G4Trk[4];
0832     G4Trk1.GetCoordinates(G4Trk);
0833 
0834     if (hGeantTrkPt)
0835       hGeantTrkPt->Fill(sqrt(G4Trk[0] * G4Trk[0] + G4Trk[1] * G4Trk[1]));
0836     if (hGeantTrkE)
0837       hGeantTrkE->Fill(G4Trk[3]);
0838   }
0839 
0840   if (verbosity > 1) {
0841     eventout += "\n          Number of G4Tracks collected:.............. ";
0842     eventout += i;
0843   }
0844 
0845   if (hMCG4Trk[0])
0846     hMCG4Trk[0]->Fill((float)i);
0847   if (hMCG4Trk[1])
0848     hMCG4Trk[1]->Fill((float)i);
0849 
0850   if (verbosity > 0)
0851     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
0852 
0853   return;
0854 }
0855 
0856 void GlobalHitsProdHist::fillTrk(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0857   nPxlHits = 0;
0858   std::string MsgLoggerCat = "GlobalHitsProdHist_fillTrk";
0859 
0860   TString eventout;
0861   if (verbosity > 0)
0862     eventout = "\nGathering info:";
0863 
0864   // access the tracker geometry
0865   const auto &theTrackerGeometry = iSetup.getHandle(tGeomToken_);
0866   if (!theTrackerGeometry.isValid()) {
0867     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerDigiGeometryRecord in event!";
0868     return;
0869   }
0870   const TrackerGeometry &theTracker(*theTrackerGeometry);
0871 
0872   // iterator to access containers
0873   edm::PSimHitContainer::const_iterator itHit;
0874 
0875   ///////////////////////////////
0876   // get Pixel Barrel information
0877   ///////////////////////////////
0878   edm::PSimHitContainer thePxlBrlHits;
0879   // extract low container
0880   edm::Handle<edm::PSimHitContainer> PxlBrlLowContainer;
0881   iEvent.getByToken(PxlBrlLowSrc_Token_, PxlBrlLowContainer);
0882   if (!PxlBrlLowContainer.isValid()) {
0883     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsPixelBarrelLowTof in event!";
0884     return;
0885   }
0886   // extract high container
0887   edm::Handle<edm::PSimHitContainer> PxlBrlHighContainer;
0888   iEvent.getByToken(PxlBrlHighSrc_Token_, PxlBrlHighContainer);
0889   if (!PxlBrlHighContainer.isValid()) {
0890     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsPixelBarrelHighTof in event!";
0891     return;
0892   }
0893   // place both containers into new container
0894   thePxlBrlHits.insert(thePxlBrlHits.end(), PxlBrlLowContainer->begin(), PxlBrlLowContainer->end());
0895   thePxlBrlHits.insert(thePxlBrlHits.end(), PxlBrlHighContainer->begin(), PxlBrlHighContainer->end());
0896 
0897   // cycle through new container
0898   int i = 0, j = 0;
0899   for (itHit = thePxlBrlHits.begin(); itHit != thePxlBrlHits.end(); ++itHit) {
0900     ++i;
0901 
0902     // create a DetId from the detUnitId
0903     DetId theDetUnitId(itHit->detUnitId());
0904     int detector = theDetUnitId.det();
0905     int subdetector = theDetUnitId.subdetId();
0906 
0907     // check that expected detector is returned
0908     if ((detector == dTrk) && (subdetector == sdPxlBrl)) {
0909       // get the GeomDetUnit from the geometry using theDetUnitID
0910       const GeomDetUnit *theDet = theTracker.idToDetUnit(theDetUnitId);
0911 
0912       if (!theDet) {
0913         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from PxlBrlHits for Hit " << i;
0914         continue;
0915       }
0916 
0917       ++j;
0918 
0919       // get the Surface of the hit (knows how to go from local <-> global)
0920       const BoundPlane &bSurface = theDet->surface();
0921 
0922       if (hTrackerPxBToF)
0923         hTrackerPxBToF->Fill(itHit->tof());
0924       if (hTrackerPxBR)
0925         hTrackerPxBR->Fill(bSurface.toGlobal(itHit->localPosition()).perp());
0926       if (hTrackerPxPhi)
0927         hTrackerPxPhi->Fill(bSurface.toGlobal(itHit->localPosition()).phi());
0928       if (hTrackerPxEta)
0929         hTrackerPxEta->Fill(bSurface.toGlobal(itHit->localPosition()).eta());
0930 
0931     } else {
0932       edm::LogWarning(MsgLoggerCat) << "PxlBrl PSimHit " << i << " is expected to be (det,subdet) = (" << dTrk << ","
0933                                     << sdPxlBrl << "); value returned is: (" << detector << "," << subdetector << ")";
0934       continue;
0935     }  // end detector type check
0936   }  // end loop through PxlBrl Hits
0937 
0938   if (verbosity > 1) {
0939     eventout += "\n          Number of Pixel Barrel Hits collected:..... ";
0940     eventout += j;
0941   }
0942 
0943   nPxlHits += j;
0944 
0945   /////////////////////////////////
0946   // get Pixel Forward information
0947   ////////////////////////////////
0948   edm::PSimHitContainer thePxlFwdHits;
0949   // extract low container
0950   edm::Handle<edm::PSimHitContainer> PxlFwdLowContainer;
0951   iEvent.getByToken(PxlFwdLowSrc_Token_, PxlFwdLowContainer);
0952   if (!PxlFwdLowContainer.isValid()) {
0953     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsPixelEndcapLowTof in event!";
0954     return;
0955   }
0956   // extract high container
0957   edm::Handle<edm::PSimHitContainer> PxlFwdHighContainer;
0958   iEvent.getByToken(PxlFwdHighSrc_Token_, PxlFwdHighContainer);
0959   if (!PxlFwdHighContainer.isValid()) {
0960     edm::LogWarning("GlobalHitsProdHist_fillTrk") << "Unable to find TrackerHitsPixelEndcapHighTof in event!";
0961     return;
0962   }
0963   // place both containers into new container
0964   thePxlFwdHits.insert(thePxlFwdHits.end(), PxlFwdLowContainer->begin(), PxlFwdLowContainer->end());
0965   thePxlFwdHits.insert(thePxlFwdHits.end(), PxlFwdHighContainer->begin(), PxlFwdHighContainer->end());
0966 
0967   // cycle through new container
0968   i = 0;
0969   j = 0;
0970   for (itHit = thePxlFwdHits.begin(); itHit != thePxlFwdHits.end(); ++itHit) {
0971     ++i;
0972 
0973     // create a DetId from the detUnitId
0974     DetId theDetUnitId(itHit->detUnitId());
0975     int detector = theDetUnitId.det();
0976     int subdetector = theDetUnitId.subdetId();
0977 
0978     // check that expected detector is returned
0979     if ((detector == dTrk) && (subdetector == sdPxlFwd)) {
0980       // get the GeomDetUnit from the geometry using theDetUnitID
0981       const GeomDetUnit *theDet = theTracker.idToDetUnit(theDetUnitId);
0982 
0983       if (!theDet) {
0984         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from PxlFwdHits for Hit " << i;
0985         ;
0986         continue;
0987       }
0988 
0989       ++j;
0990 
0991       // get the Surface of the hit (knows how to go from local <-> global)
0992       const BoundPlane &bSurface = theDet->surface();
0993 
0994       if (hTrackerPxFToF)
0995         hTrackerPxFToF->Fill(itHit->tof());
0996       if (hTrackerPxFZ)
0997         hTrackerPxFZ->Fill(bSurface.toGlobal(itHit->localPosition()).z());
0998       if (hTrackerPxPhi)
0999         hTrackerPxPhi->Fill(bSurface.toGlobal(itHit->localPosition()).phi());
1000       if (hTrackerPxEta)
1001         hTrackerPxEta->Fill(bSurface.toGlobal(itHit->localPosition()).eta());
1002 
1003     } else {
1004       edm::LogWarning(MsgLoggerCat) << "PxlFwd PSimHit " << i << " is expected to be (det,subdet) = (" << dTrk << ","
1005                                     << sdPxlFwd << "); value returned is: (" << detector << "," << subdetector << ")";
1006       continue;
1007     }  // end detector type check
1008   }  // end loop through PxlFwd Hits
1009 
1010   if (verbosity > 1) {
1011     eventout += "\n          Number of Pixel Forward Hits collected:.... ";
1012     eventout += j;
1013   }
1014 
1015   nPxlHits += j;
1016 
1017   if (hTrackerPx[0])
1018     hTrackerPx[0]->Fill((float)nPxlHits);
1019   if (hTrackerPx[1])
1020     hTrackerPx[1]->Fill((float)nPxlHits);
1021 
1022   ///////////////////////////////////
1023   // get Silicon Barrel information
1024   //////////////////////////////////
1025   nSiHits = 0;
1026   edm::PSimHitContainer theSiBrlHits;
1027   // extract TIB low container
1028   edm::Handle<edm::PSimHitContainer> SiTIBLowContainer;
1029   iEvent.getByToken(SiTIBLowSrc_Token_, SiTIBLowContainer);
1030   if (!SiTIBLowContainer.isValid()) {
1031     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTIBLowTof in event!";
1032     return;
1033   }
1034   // extract TIB high container
1035   edm::Handle<edm::PSimHitContainer> SiTIBHighContainer;
1036   iEvent.getByToken(SiTIBHighSrc_Token_, SiTIBHighContainer);
1037   if (!SiTIBHighContainer.isValid()) {
1038     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTIBHighTof in event!";
1039     return;
1040   }
1041   // extract TOB low container
1042   edm::Handle<edm::PSimHitContainer> SiTOBLowContainer;
1043   iEvent.getByToken(SiTOBLowSrc_Token_, SiTOBLowContainer);
1044   if (!SiTOBLowContainer.isValid()) {
1045     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTOBLowTof in event!";
1046     return;
1047   }
1048   // extract TOB high container
1049   edm::Handle<edm::PSimHitContainer> SiTOBHighContainer;
1050   iEvent.getByToken(SiTOBHighSrc_Token_, SiTOBHighContainer);
1051   if (!SiTOBHighContainer.isValid()) {
1052     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTOBHighTof in event!";
1053     return;
1054   }
1055   // place all containers into new container
1056   theSiBrlHits.insert(theSiBrlHits.end(), SiTIBLowContainer->begin(), SiTIBLowContainer->end());
1057   theSiBrlHits.insert(theSiBrlHits.end(), SiTIBHighContainer->begin(), SiTIBHighContainer->end());
1058   theSiBrlHits.insert(theSiBrlHits.end(), SiTOBLowContainer->begin(), SiTOBLowContainer->end());
1059   theSiBrlHits.insert(theSiBrlHits.end(), SiTOBHighContainer->begin(), SiTOBHighContainer->end());
1060 
1061   // cycle through new container
1062   i = 0;
1063   j = 0;
1064   for (itHit = theSiBrlHits.begin(); itHit != theSiBrlHits.end(); ++itHit) {
1065     ++i;
1066 
1067     // create a DetId from the detUnitId
1068     DetId theDetUnitId(itHit->detUnitId());
1069     int detector = theDetUnitId.det();
1070     int subdetector = theDetUnitId.subdetId();
1071 
1072     // check that expected detector is returned
1073     if ((detector == dTrk) && ((subdetector == sdSiTIB) || (subdetector == sdSiTOB))) {
1074       // get the GeomDetUnit from the geometry using theDetUnitID
1075       const GeomDetUnit *theDet = theTracker.idToDetUnit(theDetUnitId);
1076 
1077       if (!theDet) {
1078         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from SiBrlHits for Hit " << i;
1079         continue;
1080       }
1081 
1082       ++j;
1083 
1084       // get the Surface of the hit (knows how to go from local <-> global)
1085       const BoundPlane &bSurface = theDet->surface();
1086 
1087       if (hTrackerSiBToF)
1088         hTrackerSiBToF->Fill(itHit->tof());
1089       if (hTrackerSiBR)
1090         hTrackerSiBR->Fill(bSurface.toGlobal(itHit->localPosition()).perp());
1091       if (hTrackerSiPhi)
1092         hTrackerSiPhi->Fill(bSurface.toGlobal(itHit->localPosition()).phi());
1093       if (hTrackerSiEta)
1094         hTrackerSiEta->Fill(bSurface.toGlobal(itHit->localPosition()).eta());
1095 
1096     } else {
1097       edm::LogWarning(MsgLoggerCat) << "SiBrl PSimHit " << i << " is expected to be (det,subdet) = (" << dTrk << ","
1098                                     << sdSiTIB << " || " << sdSiTOB << "); value returned is: (" << detector << ","
1099                                     << subdetector << ")";
1100       continue;
1101     }  // end detector type check
1102   }  // end loop through SiBrl Hits
1103 
1104   if (verbosity > 1) {
1105     eventout += "\n          Number of Silicon Barrel Hits collected:... ";
1106     eventout += j;
1107   }
1108 
1109   nSiHits += j;
1110 
1111   ///////////////////////////////////
1112   // get Silicon Forward information
1113   ///////////////////////////////////
1114   edm::PSimHitContainer theSiFwdHits;
1115   // extract TID low container
1116   edm::Handle<edm::PSimHitContainer> SiTIDLowContainer;
1117   iEvent.getByToken(SiTIDLowSrc_Token_, SiTIDLowContainer);
1118   if (!SiTIDLowContainer.isValid()) {
1119     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTIDLowTof in event!";
1120     return;
1121   }
1122   // extract TID high container
1123   edm::Handle<edm::PSimHitContainer> SiTIDHighContainer;
1124   iEvent.getByToken(SiTIDHighSrc_Token_, SiTIDHighContainer);
1125   if (!SiTIDHighContainer.isValid()) {
1126     edm::LogWarning("GlobalHitsProdHist_fillTrk") << "Unable to find TrackerHitsTIDHighTof in event!";
1127     return;
1128   }
1129   // extract TEC low container
1130   edm::Handle<edm::PSimHitContainer> SiTECLowContainer;
1131   iEvent.getByToken(SiTECLowSrc_Token_, SiTECLowContainer);
1132   if (!SiTECLowContainer.isValid()) {
1133     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTECLowTof in event!";
1134     return;
1135   }
1136   // extract TEC high container
1137   edm::Handle<edm::PSimHitContainer> SiTECHighContainer;
1138   iEvent.getByToken(SiTECHighSrc_Token_, SiTECHighContainer);
1139   if (!SiTECHighContainer.isValid()) {
1140     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTECHighTof in event!";
1141     return;
1142   }
1143   // place all containers into new container
1144   theSiFwdHits.insert(theSiFwdHits.end(), SiTIDLowContainer->begin(), SiTIDLowContainer->end());
1145   theSiFwdHits.insert(theSiFwdHits.end(), SiTIDHighContainer->begin(), SiTIDHighContainer->end());
1146   theSiFwdHits.insert(theSiFwdHits.end(), SiTECLowContainer->begin(), SiTECLowContainer->end());
1147   theSiFwdHits.insert(theSiFwdHits.end(), SiTECHighContainer->begin(), SiTECHighContainer->end());
1148 
1149   // cycle through container
1150   i = 0;
1151   j = 0;
1152   for (itHit = theSiFwdHits.begin(); itHit != theSiFwdHits.end(); ++itHit) {
1153     ++i;
1154 
1155     // create a DetId from the detUnitId
1156     DetId theDetUnitId(itHit->detUnitId());
1157     int detector = theDetUnitId.det();
1158     int subdetector = theDetUnitId.subdetId();
1159 
1160     // check that expected detector is returned
1161     if ((detector == dTrk) && ((subdetector == sdSiTID) || (subdetector == sdSiTEC))) {
1162       // get the GeomDetUnit from the geometry using theDetUnitID
1163       const GeomDetUnit *theDet = theTracker.idToDetUnit(theDetUnitId);
1164 
1165       if (!theDet) {
1166         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from SiFwdHits Hit " << i;
1167         return;
1168       }
1169 
1170       ++j;
1171 
1172       // get the Surface of the hit (knows how to go from local <-> global)
1173       const BoundPlane &bSurface = theDet->surface();
1174 
1175       if (hTrackerSiFToF)
1176         hTrackerSiFToF->Fill(itHit->tof());
1177       if (hTrackerSiFZ)
1178         hTrackerSiFZ->Fill(bSurface.toGlobal(itHit->localPosition()).z());
1179       if (hTrackerSiPhi)
1180         hTrackerSiPhi->Fill(bSurface.toGlobal(itHit->localPosition()).phi());
1181       if (hTrackerSiEta)
1182         hTrackerSiEta->Fill(bSurface.toGlobal(itHit->localPosition()).eta());
1183 
1184     } else {
1185       edm::LogWarning(MsgLoggerCat) << "SiFwd PSimHit " << i << " is expected to be (det,subdet) = (" << dTrk << ","
1186                                     << sdSiTOB << " || " << sdSiTEC << "); value returned is: (" << detector << ","
1187                                     << subdetector << ")";
1188       continue;
1189     }  // end check detector type
1190   }  // end loop through SiFwd Hits
1191 
1192   if (verbosity > 1) {
1193     eventout += "\n          Number of Silicon Forward Hits collected:.. ";
1194     eventout += j;
1195   }
1196 
1197   nSiHits += j;
1198 
1199   if (hTrackerSi[0])
1200     hTrackerSi[0]->Fill((float)nSiHits);
1201   if (hTrackerSi[1])
1202     hTrackerSi[1]->Fill((float)nSiHits);
1203 
1204   if (verbosity > 0)
1205     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
1206 
1207   return;
1208 }
1209 
1210 void GlobalHitsProdHist::fillMuon(edm::Event &iEvent, const edm::EventSetup &iSetup) {
1211   nMuonHits = 0;
1212   std::string MsgLoggerCat = "GlobalHitsProdHist_fillMuon";
1213 
1214   TString eventout;
1215   if (verbosity > 0)
1216     eventout = "\nGathering info:";
1217 
1218   // iterator to access containers
1219   edm::PSimHitContainer::const_iterator itHit;
1220 
1221   ///////////////////////
1222   // access the CSC Muon
1223   ///////////////////////
1224   // access the CSC Muon geometry
1225   const auto &theCSCGeometry = iSetup.getHandle(cscGeomToken_);
1226   if (!theCSCGeometry.isValid()) {
1227     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonGeometryRecord for the CSCGeometry in event!";
1228     return;
1229   }
1230   const CSCGeometry &theCSCMuon(*theCSCGeometry);
1231 
1232   // get Muon CSC information
1233   edm::Handle<edm::PSimHitContainer> MuonCSCContainer;
1234   iEvent.getByToken(MuonCscSrc_Token_, MuonCSCContainer);
1235   if (!MuonCSCContainer.isValid()) {
1236     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonCSCHits in event!";
1237     return;
1238   }
1239 
1240   // cycle through container
1241   int i = 0, j = 0;
1242   for (itHit = MuonCSCContainer->begin(); itHit != MuonCSCContainer->end(); ++itHit) {
1243     ++i;
1244 
1245     // create a DetId from the detUnitId
1246     DetId theDetUnitId(itHit->detUnitId());
1247     int detector = theDetUnitId.det();
1248     int subdetector = theDetUnitId.subdetId();
1249 
1250     // check that expected detector is returned
1251     if ((detector == dMuon) && (subdetector == sdMuonCSC)) {
1252       // get the GeomDetUnit from the geometry using theDetUnitID
1253       const GeomDetUnit *theDet = theCSCMuon.idToDetUnit(theDetUnitId);
1254 
1255       if (!theDet) {
1256         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from theCSCMuon for hit " << i;
1257         continue;
1258       }
1259 
1260       ++j;
1261 
1262       // get the Surface of the hit (knows how to go from local <-> global)
1263       const BoundPlane &bSurface = theDet->surface();
1264 
1265       if (hMuonCscToF[0])
1266         hMuonCscToF[0]->Fill(itHit->tof());
1267       if (hMuonCscToF[1])
1268         hMuonCscToF[1]->Fill(itHit->tof());
1269       if (hMuonCscZ)
1270         hMuonCscZ->Fill(bSurface.toGlobal(itHit->localPosition()).z());
1271       if (hMuonPhi)
1272         hMuonPhi->Fill(bSurface.toGlobal(itHit->localPosition()).phi());
1273       if (hMuonEta)
1274         hMuonEta->Fill(bSurface.toGlobal(itHit->localPosition()).eta());
1275 
1276     } else {
1277       edm::LogWarning(MsgLoggerCat) << "MuonCsc PSimHit " << i << " is expected to be (det,subdet) = (" << dMuon << ","
1278                                     << sdMuonCSC << "); value returned is: (" << detector << "," << subdetector << ")";
1279       continue;
1280     }  // end detector type check
1281   }  // end loop through CSC Hits
1282 
1283   if (verbosity > 1) {
1284     eventout += "\n          Number of CSC muon Hits collected:......... ";
1285     eventout += j;
1286   }
1287 
1288   nMuonHits += j;
1289 
1290   /////////////////////
1291   // access the DT Muon
1292   /////////////////////
1293   // access the DT Muon geometry
1294   const auto &theDTGeometry = iSetup.getHandle(dtGeomToken_);
1295   if (!theDTGeometry.isValid()) {
1296     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonGeometryRecord for the DTGeometry in event!";
1297     return;
1298   }
1299   const DTGeometry &theDTMuon(*theDTGeometry);
1300 
1301   // get Muon DT information
1302   edm::Handle<edm::PSimHitContainer> MuonDtContainer;
1303   iEvent.getByToken(MuonDtSrc_Token_, MuonDtContainer);
1304   if (!MuonDtContainer.isValid()) {
1305     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonDTHits in event!";
1306     return;
1307   }
1308 
1309   // cycle through container
1310   i = 0, j = 0;
1311   for (itHit = MuonDtContainer->begin(); itHit != MuonDtContainer->end(); ++itHit) {
1312     ++i;
1313 
1314     // create a DetId from the detUnitId
1315     DetId theDetUnitId(itHit->detUnitId());
1316     int detector = theDetUnitId.det();
1317     int subdetector = theDetUnitId.subdetId();
1318 
1319     // check that expected detector is returned
1320     if ((detector == dMuon) && (subdetector == sdMuonDT)) {
1321       // CSC uses wires and layers rather than the full detID
1322       // get the wireId
1323       DTWireId wireId(itHit->detUnitId());
1324 
1325       // get the DTLayer from the geometry using the wireID
1326       const DTLayer *theDet = theDTMuon.layer(wireId.layerId());
1327 
1328       if (!theDet) {
1329         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from theDtMuon for hit " << i;
1330         continue;
1331       }
1332 
1333       ++j;
1334 
1335       // get the Surface of the hit (knows how to go from local <-> global)
1336       const BoundPlane &bSurface = theDet->surface();
1337 
1338       if (hMuonDtToF[0])
1339         hMuonDtToF[0]->Fill(itHit->tof());
1340       if (hMuonDtToF[1])
1341         hMuonDtToF[1]->Fill(itHit->tof());
1342       if (hMuonDtR)
1343         hMuonDtR->Fill(bSurface.toGlobal(itHit->localPosition()).perp());
1344       if (hMuonPhi)
1345         hMuonPhi->Fill(bSurface.toGlobal(itHit->localPosition()).phi());
1346       if (hMuonEta)
1347         hMuonEta->Fill(bSurface.toGlobal(itHit->localPosition()).eta());
1348 
1349     } else {
1350       edm::LogWarning(MsgLoggerCat) << "MuonDt PSimHit " << i << " is expected to be (det,subdet) = (" << dMuon << ","
1351                                     << sdMuonDT << "); value returned is: (" << detector << "," << subdetector << ")";
1352       continue;
1353     }  // end detector type check
1354   }  // end loop through DT Hits
1355 
1356   if (verbosity > 1) {
1357     eventout += "\n          Number of DT muon Hits collected:.......... ";
1358     eventout += j;
1359   }
1360 
1361   nMuonHits += j;
1362 
1363   // int RPCBrl = 0, RPCFwd = 0;
1364   ///////////////////////
1365   // access the RPC Muon
1366   ///////////////////////
1367   // access the RPC Muon geometry
1368   const auto &theRPCGeometry = iSetup.getHandle(rpcGeomToken_);
1369   if (!theRPCGeometry.isValid()) {
1370     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonGeometryRecord for the RPCGeometry in event!";
1371     return;
1372   }
1373   const RPCGeometry &theRPCMuon(*theRPCGeometry);
1374 
1375   // get Muon RPC information
1376   edm::Handle<edm::PSimHitContainer> MuonRPCContainer;
1377   iEvent.getByToken(MuonRpcSrc_Token_, MuonRPCContainer);
1378   if (!MuonRPCContainer.isValid()) {
1379     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonRPCHits in event!";
1380     return;
1381   }
1382 
1383   // cycle through container
1384   i = 0, j = 0;
1385   int RPCBrl = 0, RPCFwd = 0;
1386   for (itHit = MuonRPCContainer->begin(); itHit != MuonRPCContainer->end(); ++itHit) {
1387     ++i;
1388 
1389     // create a DetID from the detUnitId
1390     DetId theDetUnitId(itHit->detUnitId());
1391     int detector = theDetUnitId.det();
1392     int subdetector = theDetUnitId.subdetId();
1393 
1394     // check that expected detector is returned
1395     if ((detector == dMuon) && (subdetector == sdMuonRPC)) {
1396       // get an RPCDetID from the detUnitID
1397       RPCDetId RPCId(itHit->detUnitId());
1398 
1399       // find the region of the RPC hit
1400       int region = RPCId.region();
1401 
1402       // get the GeomDetUnit from the geometry using the RPCDetId
1403       const GeomDetUnit *theDet = theRPCMuon.idToDetUnit(theDetUnitId);
1404 
1405       if (!theDet) {
1406         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from theRPCMuon for hit " << i;
1407         continue;
1408       }
1409 
1410       ++j;
1411 
1412       // get the Surface of the hit (knows how to go from local <-> global)
1413       const BoundPlane &bSurface = theDet->surface();
1414 
1415       // gather necessary information
1416       if ((region == sdMuonRPCRgnFwdp) || (region == sdMuonRPCRgnFwdn)) {
1417         ++RPCFwd;
1418 
1419         if (hMuonRpcFToF[0])
1420           hMuonRpcFToF[0]->Fill(itHit->tof());
1421         if (hMuonRpcFToF[1])
1422           hMuonRpcFToF[1]->Fill(itHit->tof());
1423         if (hMuonRpcFZ)
1424           hMuonRpcFZ->Fill(bSurface.toGlobal(itHit->localPosition()).z());
1425         if (hMuonPhi)
1426           hMuonPhi->Fill(bSurface.toGlobal(itHit->localPosition()).phi());
1427         if (hMuonEta)
1428           hMuonEta->Fill(bSurface.toGlobal(itHit->localPosition()).eta());
1429 
1430       } else if (region == sdMuonRPCRgnBrl) {
1431         ++RPCBrl;
1432 
1433         if (hMuonRpcBToF[0])
1434           hMuonRpcBToF[0]->Fill(itHit->tof());
1435         if (hMuonRpcBToF[1])
1436           hMuonRpcBToF[1]->Fill(itHit->tof());
1437         if (hMuonRpcBR)
1438           hMuonRpcBR->Fill(bSurface.toGlobal(itHit->localPosition()).perp());
1439         if (hMuonPhi)
1440           hMuonPhi->Fill(bSurface.toGlobal(itHit->localPosition()).phi());
1441         if (hMuonEta)
1442           hMuonEta->Fill(bSurface.toGlobal(itHit->localPosition()).eta());
1443 
1444       } else {
1445         edm::LogWarning(MsgLoggerCat) << "Invalid region for RPC Muon hit" << i;
1446         continue;
1447       }  // end check of region
1448     } else {
1449       edm::LogWarning(MsgLoggerCat) << "MuonRpc PSimHit " << i << " is expected to be (det,subdet) = (" << dMuon << ","
1450                                     << sdMuonRPC << "); value returned is: (" << detector << "," << subdetector << ")";
1451       continue;
1452     }  // end detector type check
1453   }  // end loop through RPC Hits
1454 
1455   if (verbosity > 1) {
1456     eventout += "\n          Number of RPC muon Hits collected:......... ";
1457     eventout += j;
1458     eventout += "\n                    RPC Barrel muon Hits:............ ";
1459     eventout += RPCBrl;
1460     eventout += "\n                    RPC Forward muon Hits:........... ";
1461     eventout += RPCFwd;
1462   }
1463 
1464   nMuonHits += j;
1465 
1466   if (hMuon[0])
1467     hMuon[0]->Fill((float)nMuonHits);
1468   if (hMuon[1])
1469     hMuon[1]->Fill((float)nMuonHits);
1470 
1471   if (verbosity > 0)
1472     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
1473 
1474   return;
1475 }
1476 
1477 void GlobalHitsProdHist::fillECal(edm::Event &iEvent, const edm::EventSetup &iSetup) {
1478   std::string MsgLoggerCat = "GlobalHitsProdHist_fillECal";
1479 
1480   TString eventout;
1481   if (verbosity > 0)
1482     eventout = "\nGathering info:";
1483 
1484   // access the calorimeter geometry
1485   const auto &theCaloGeometry = iSetup.getHandle(caloGeomToken_);
1486   if (!theCaloGeometry.isValid()) {
1487     edm::LogWarning(MsgLoggerCat) << "Unable to find CaloGeometryRecord in event!";
1488     return;
1489   }
1490   const CaloGeometry &theCalo(*theCaloGeometry);
1491 
1492   // iterator to access containers
1493   edm::PCaloHitContainer::const_iterator itHit;
1494 
1495   ///////////////////////////////
1496   // get  ECal information
1497   ///////////////////////////////
1498   edm::PCaloHitContainer theECalHits;
1499   // extract EB container
1500   edm::Handle<edm::PCaloHitContainer> EBContainer;
1501   iEvent.getByToken(ECalEBSrc_Token_, EBContainer);
1502   if (!EBContainer.isValid()) {
1503     edm::LogWarning(MsgLoggerCat) << "Unable to find EcalHitsEB in event!";
1504     return;
1505   }
1506   // extract EE container
1507   edm::Handle<edm::PCaloHitContainer> EEContainer;
1508   iEvent.getByToken(ECalEESrc_Token_, EEContainer);
1509   if (!EEContainer.isValid()) {
1510     edm::LogWarning(MsgLoggerCat) << "Unable to find EcalHitsEE in event!";
1511     return;
1512   }
1513   // place both containers into new container
1514   theECalHits.insert(theECalHits.end(), EBContainer->begin(), EBContainer->end());
1515   theECalHits.insert(theECalHits.end(), EEContainer->begin(), EEContainer->end());
1516 
1517   // cycle through new container
1518   int i = 0, j = 0;
1519   for (itHit = theECalHits.begin(); itHit != theECalHits.end(); ++itHit) {
1520     ++i;
1521 
1522     // create a DetId from the detUnitId
1523     DetId theDetUnitId(itHit->id());
1524     int detector = theDetUnitId.det();
1525     int subdetector = theDetUnitId.subdetId();
1526 
1527     // check that expected detector is returned
1528     if ((detector == dEcal) && ((subdetector == sdEcalBrl) || (subdetector == sdEcalFwd))) {
1529       // get the Cell geometry
1530       auto theDet = (theCalo.getSubdetectorGeometry(theDetUnitId))->getGeometry(theDetUnitId);
1531 
1532       if (!theDet) {
1533         edm::LogWarning(MsgLoggerCat) << "Unable to get CaloCellGeometry from ECalHits for Hit " << i;
1534         continue;
1535       }
1536 
1537       ++j;
1538 
1539       // get the global position of the cell
1540       const GlobalPoint &globalposition = theDet->getPosition();
1541 
1542       if (hCaloEcalE[0])
1543         hCaloEcalE[0]->Fill(itHit->energy());
1544       if (hCaloEcalE[1])
1545         hCaloEcalE[1]->Fill(itHit->energy());
1546       if (hCaloEcalToF[0])
1547         hCaloEcalToF[0]->Fill(itHit->time());
1548       if (hCaloEcalToF[1])
1549         hCaloEcalToF[1]->Fill(itHit->time());
1550       if (hCaloEcalPhi)
1551         hCaloEcalPhi->Fill(globalposition.phi());
1552       if (hCaloEcalEta)
1553         hCaloEcalEta->Fill(globalposition.eta());
1554 
1555     } else {
1556       edm::LogWarning(MsgLoggerCat) << "ECal PCaloHit " << i << " is expected to be (det,subdet) = (" << dEcal << ","
1557                                     << sdEcalBrl << " || " << sdEcalFwd << "); value returned is: (" << detector << ","
1558                                     << subdetector << ")";
1559       continue;
1560     }  // end detector type check
1561   }  // end loop through ECal Hits
1562 
1563   if (verbosity > 1) {
1564     eventout += "\n          Number of ECal Hits collected:............. ";
1565     eventout += j;
1566   }
1567 
1568   if (hCaloEcal[0])
1569     hCaloEcal[0]->Fill((float)j);
1570   if (hCaloEcal[1])
1571     hCaloEcal[1]->Fill((float)j);
1572 
1573   ////////////////////////////
1574   // Get Preshower information
1575   ////////////////////////////
1576   // extract PreShower container
1577   edm::Handle<edm::PCaloHitContainer> PreShContainer;
1578   iEvent.getByToken(ECalESSrc_Token_, PreShContainer);
1579   if (!PreShContainer.isValid()) {
1580     edm::LogWarning(MsgLoggerCat) << "Unable to find EcalHitsES in event!";
1581     return;
1582   }
1583 
1584   // cycle through container
1585   i = 0, j = 0;
1586   for (itHit = PreShContainer->begin(); itHit != PreShContainer->end(); ++itHit) {
1587     ++i;
1588 
1589     // create a DetId from the detUnitId
1590     DetId theDetUnitId(itHit->id());
1591     int detector = theDetUnitId.det();
1592     int subdetector = theDetUnitId.subdetId();
1593 
1594     // check that expected detector is returned
1595     if ((detector == dEcal) && (subdetector == sdEcalPS)) {
1596       // get the Cell geometry
1597       auto theDet = (theCalo.getSubdetectorGeometry(theDetUnitId))->getGeometry(theDetUnitId);
1598 
1599       if (!theDet) {
1600         edm::LogWarning(MsgLoggerCat) << "Unable to get CaloCellGeometry from PreShContainer for Hit " << i;
1601         continue;
1602       }
1603 
1604       ++j;
1605 
1606       // get the global position of the cell
1607       const GlobalPoint &globalposition = theDet->getPosition();
1608 
1609       if (hCaloPreShE[0])
1610         hCaloPreShE[0]->Fill(itHit->energy());
1611       if (hCaloPreShE[1])
1612         hCaloPreShE[1]->Fill(itHit->energy());
1613       if (hCaloPreShToF[0])
1614         hCaloPreShToF[0]->Fill(itHit->time());
1615       if (hCaloPreShToF[1])
1616         hCaloPreShToF[1]->Fill(itHit->time());
1617       if (hCaloPreShPhi)
1618         hCaloPreShPhi->Fill(globalposition.phi());
1619       if (hCaloPreShEta)
1620         hCaloPreShEta->Fill(globalposition.eta());
1621 
1622     } else {
1623       edm::LogWarning(MsgLoggerCat) << "PreSh PCaloHit " << i << " is expected to be (det,subdet) = (" << dEcal << ","
1624                                     << sdEcalPS << "); value returned is: (" << detector << "," << subdetector << ")";
1625       continue;
1626     }  // end detector type check
1627   }  // end loop through PreShower Hits
1628 
1629   if (verbosity > 1) {
1630     eventout += "\n          Number of PreSh Hits collected:............ ";
1631     eventout += j;
1632   }
1633 
1634   if (hCaloPreSh[0])
1635     hCaloPreSh[0]->Fill((float)j);
1636   if (hCaloPreSh[1])
1637     hCaloPreSh[1]->Fill((float)j);
1638 
1639   if (verbosity > 0)
1640     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
1641 
1642   return;
1643 }
1644 
1645 void GlobalHitsProdHist::fillHCal(edm::Event &iEvent, const edm::EventSetup &iSetup) {
1646   std::string MsgLoggerCat = "GlobalHitsProdHist_fillHCal";
1647 
1648   TString eventout;
1649   if (verbosity > 0)
1650     eventout = "\nGathering info:";
1651 
1652   // access the calorimeter geometry
1653   const auto &theCaloGeometry = iSetup.getHandle(caloGeomToken_);
1654   if (!theCaloGeometry.isValid()) {
1655     edm::LogWarning(MsgLoggerCat) << "Unable to find CaloGeometryRecord in event!";
1656     return;
1657   }
1658   const CaloGeometry &theCalo(*theCaloGeometry);
1659 
1660   // iterator to access containers
1661   edm::PCaloHitContainer::const_iterator itHit;
1662 
1663   ///////////////////////////////
1664   // get  HCal information
1665   ///////////////////////////////
1666   // extract HCal container
1667   edm::Handle<edm::PCaloHitContainer> HCalContainer;
1668   iEvent.getByToken(HCalSrc_Token_, HCalContainer);
1669   if (!HCalContainer.isValid()) {
1670     edm::LogWarning(MsgLoggerCat) << "Unable to find HCalHits in event!";
1671     return;
1672   }
1673 
1674   // cycle through container
1675   int i = 0, j = 0;
1676   for (itHit = HCalContainer->begin(); itHit != HCalContainer->end(); ++itHit) {
1677     ++i;
1678 
1679     // create a DetId from the detUnitId
1680     DetId theDetUnitId(itHit->id());
1681     int detector = theDetUnitId.det();
1682     int subdetector = theDetUnitId.subdetId();
1683 
1684     // check that expected detector is returned
1685     if ((detector == dHcal) && ((subdetector == sdHcalBrl) || (subdetector == sdHcalEC) || (subdetector == sdHcalOut) ||
1686                                 (subdetector == sdHcalFwd))) {
1687       // get the Cell geometry
1688       const HcalGeometry *theDet = dynamic_cast<const HcalGeometry *>(theCalo.getSubdetectorGeometry(theDetUnitId));
1689 
1690       if (!theDet) {
1691         edm::LogWarning(MsgLoggerCat) << "Unable to get HcalGeometry from HCalContainer for Hit " << i;
1692         continue;
1693       }
1694 
1695       ++j;
1696 
1697       // get the global position of the cell
1698       const GlobalPoint &globalposition = theDet->getPosition(theDetUnitId);
1699 
1700       if (hCaloHcalE[0])
1701         hCaloHcalE[0]->Fill(itHit->energy());
1702       if (hCaloHcalE[1])
1703         hCaloHcalE[1]->Fill(itHit->energy());
1704       if (hCaloHcalToF[0])
1705         hCaloHcalToF[0]->Fill(itHit->time());
1706       if (hCaloHcalToF[1])
1707         hCaloHcalToF[1]->Fill(itHit->time());
1708       if (hCaloHcalPhi)
1709         hCaloHcalPhi->Fill(globalposition.phi());
1710       if (hCaloHcalEta)
1711         hCaloHcalEta->Fill(globalposition.eta());
1712 
1713     } else {
1714       edm::LogWarning(MsgLoggerCat) << "HCal PCaloHit " << i << " is expected to be (det,subdet) = (" << dHcal << ","
1715                                     << sdHcalBrl << " || " << sdHcalEC << " || " << sdHcalOut << " || " << sdHcalFwd
1716                                     << "); value returned is: (" << detector << "," << subdetector << ")";
1717       continue;
1718     }  // end detector type check
1719   }  // end loop through HCal Hits
1720 
1721   if (verbosity > 1) {
1722     eventout += "\n          Number of HCal Hits collected:............. ";
1723     eventout += j;
1724   }
1725 
1726   if (hCaloHcal[0])
1727     hCaloHcal[0]->Fill((float)j);
1728   if (hCaloHcal[1])
1729     hCaloHcal[1]->Fill((float)j);
1730 
1731   if (verbosity > 0)
1732     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
1733 
1734   return;
1735 }
[ 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!