Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 /** \file GlobalHitsProducer.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/GlobalHitsProducer.h"
0012 
0013 GlobalHitsProducer::GlobalHitsProducer(const edm::ParameterSet &iPSet)
0014     : fName(""),
0015       verbosity(0),
0016       frequency(0),
0017       vtxunit(0),
0018       label(""),
0019       getAllProvenances(false),
0020       printProvenanceInfo(false),
0021       nRawGenPart(0),
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 = "GlobalHitsProducer_GlobalHitsProducer";
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   label = iPSet.getParameter<std::string>("Label");
0038   edm::ParameterSet m_Prov = iPSet.getParameter<edm::ParameterSet>("ProvenanceLookup");
0039   getAllProvenances = m_Prov.getUntrackedParameter<bool>("GetAllProvenances");
0040   printProvenanceInfo = m_Prov.getUntrackedParameter<bool>("PrintProvenanceInfo");
0041 
0042   // get Labels to use to extract information
0043   PxlBrlLowSrc_ = iPSet.getParameter<edm::InputTag>("PxlBrlLowSrc");
0044   PxlBrlHighSrc_ = iPSet.getParameter<edm::InputTag>("PxlBrlHighSrc");
0045   PxlFwdLowSrc_ = iPSet.getParameter<edm::InputTag>("PxlFwdLowSrc");
0046   PxlFwdHighSrc_ = iPSet.getParameter<edm::InputTag>("PxlFwdHighSrc");
0047 
0048   SiTIBLowSrc_ = iPSet.getParameter<edm::InputTag>("SiTIBLowSrc");
0049   SiTIBHighSrc_ = iPSet.getParameter<edm::InputTag>("SiTIBHighSrc");
0050   SiTOBLowSrc_ = iPSet.getParameter<edm::InputTag>("SiTOBLowSrc");
0051   SiTOBHighSrc_ = iPSet.getParameter<edm::InputTag>("SiTOBHighSrc");
0052   SiTIDLowSrc_ = iPSet.getParameter<edm::InputTag>("SiTIDLowSrc");
0053   SiTIDHighSrc_ = iPSet.getParameter<edm::InputTag>("SiTIDHighSrc");
0054   SiTECLowSrc_ = iPSet.getParameter<edm::InputTag>("SiTECLowSrc");
0055   SiTECHighSrc_ = iPSet.getParameter<edm::InputTag>("SiTECHighSrc");
0056 
0057   MuonCscSrc_ = iPSet.getParameter<edm::InputTag>("MuonCscSrc");
0058   MuonDtSrc_ = iPSet.getParameter<edm::InputTag>("MuonDtSrc");
0059   MuonRpcSrc_ = iPSet.getParameter<edm::InputTag>("MuonRpcSrc");
0060 
0061   ECalEBSrc_ = iPSet.getParameter<edm::InputTag>("ECalEBSrc");
0062   ECalEESrc_ = iPSet.getParameter<edm::InputTag>("ECalEESrc");
0063   ECalESSrc_ = iPSet.getParameter<edm::InputTag>("ECalESSrc");
0064 
0065   HCalSrc_ = iPSet.getParameter<edm::InputTag>("HCalSrc");
0066 
0067   // fix for consumes
0068   PxlBrlLowSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("PxlBrlLowSrc"));
0069   PxlBrlHighSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("PxlBrlHighSrc"));
0070   PxlFwdLowSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("PxlFwdLowSrc"));
0071   PxlFwdHighSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("PxlFwdHighSrc"));
0072 
0073   SiTIBLowSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("SiTIBLowSrc"));
0074   SiTIBHighSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("SiTIBHighSrc"));
0075   SiTOBLowSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("SiTOBLowSrc"));
0076   SiTOBHighSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("SiTOBHighSrc"));
0077   SiTIDLowSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("SiTIDLowSrc"));
0078   SiTIDHighSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("SiTIDHighSrc"));
0079   SiTECLowSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("SiTECLowSrc"));
0080   SiTECHighSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("SiTECHighSrc"));
0081 
0082   MuonCscSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("MuonCscSrc"));
0083   MuonDtSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("MuonDtSrc"));
0084   MuonRpcSrc_Token_ = consumes<edm::PSimHitContainer>(iPSet.getParameter<edm::InputTag>("MuonRpcSrc"));
0085 
0086   ECalEBSrc_Token_ = consumes<edm::PCaloHitContainer>(iPSet.getParameter<edm::InputTag>("ECalEBSrc"));
0087   ECalEESrc_Token_ = consumes<edm::PCaloHitContainer>(iPSet.getParameter<edm::InputTag>("ECalEESrc"));
0088   ECalESSrc_Token_ = consumes<edm::PCaloHitContainer>(iPSet.getParameter<edm::InputTag>("ECalESSrc"));
0089   HCalSrc_Token_ = consumes<edm::PCaloHitContainer>(iPSet.getParameter<edm::InputTag>("HCalSrc"));
0090 
0091   // use value of first digit to determine default output level (inclusive)
0092   // 0 is none, 1 is basic, 2 is fill output, 3 is gather output
0093   verbosity %= 10;
0094 
0095   // create persistent object
0096   produces<PGlobalSimHit>(label);
0097 
0098   // print out Parameter Set information being used
0099   if (verbosity >= 0) {
0100     edm::LogInfo(MsgLoggerCat)
0101         << "\n===============================\n"
0102         << "Initialized as EDProducer with parameter values:\n"
0103         << "    Name          = " << fName << "\n"
0104         << "    Verbosity     = " << verbosity << "\n"
0105         << "    Frequency     = " << frequency << "\n"
0106         << "    VtxUnit       = " << vtxunit << "\n"
0107         << "    Label         = " << label << "\n"
0108         << "    GetProv       = " << getAllProvenances << "\n"
0109         << "    PrintProv     = " << printProvenanceInfo << "\n"
0110         << "    PxlBrlLowSrc  = " << PxlBrlLowSrc_.label() << ":" << PxlBrlLowSrc_.instance() << "\n"
0111         << "    PxlBrlHighSrc = " << PxlBrlHighSrc_.label() << ":" << PxlBrlHighSrc_.instance() << "\n"
0112         << "    PxlFwdLowSrc  = " << PxlFwdLowSrc_.label() << ":" << PxlBrlLowSrc_.instance() << "\n"
0113         << "    PxlFwdHighSrc = " << PxlFwdHighSrc_.label() << ":" << PxlBrlHighSrc_.instance() << "\n"
0114         << "    SiTIBLowSrc   = " << SiTIBLowSrc_.label() << ":" << SiTIBLowSrc_.instance() << "\n"
0115         << "    SiTIBHighSrc  = " << SiTIBHighSrc_.label() << ":" << SiTIBHighSrc_.instance() << "\n"
0116         << "    SiTOBLowSrc   = " << SiTOBLowSrc_.label() << ":" << SiTOBLowSrc_.instance() << "\n"
0117         << "    SiTOBHighSrc  = " << SiTOBHighSrc_.label() << ":" << SiTOBHighSrc_.instance() << "\n"
0118         << "    SiTIDLowSrc   = " << SiTIDLowSrc_.label() << ":" << SiTIDLowSrc_.instance() << "\n"
0119         << "    SiTIDHighSrc  = " << SiTIDHighSrc_.label() << ":" << SiTIDHighSrc_.instance() << "\n"
0120         << "    SiTECLowSrc   = " << SiTECLowSrc_.label() << ":" << SiTECLowSrc_.instance() << "\n"
0121         << "    SiTECHighSrc  = " << SiTECHighSrc_.label() << ":" << SiTECHighSrc_.instance() << "\n"
0122         << "    MuonCscSrc    = " << MuonCscSrc_.label() << ":" << MuonCscSrc_.instance() << "\n"
0123         << "    MuonDtSrc     = " << MuonDtSrc_.label() << ":" << MuonDtSrc_.instance() << "\n"
0124         << "    MuonRpcSrc    = " << MuonRpcSrc_.label() << ":" << MuonRpcSrc_.instance() << "\n"
0125         << "    ECalEBSrc     = " << ECalEBSrc_.label() << ":" << ECalEBSrc_.instance() << "\n"
0126         << "    ECalEESrc     = " << ECalEESrc_.label() << ":" << ECalEESrc_.instance() << "\n"
0127         << "    ECalESSrc     = " << ECalESSrc_.label() << ":" << ECalESSrc_.instance() << "\n"
0128         << "    HCalSrc       = " << HCalSrc_.label() << ":" << HCalSrc_.instance() << "\n"
0129         << "===============================\n";
0130   }
0131 
0132   // migrated here from beginJob
0133   clear();
0134 }
0135 
0136 GlobalHitsProducer::~GlobalHitsProducer() {}
0137 
0138 void GlobalHitsProducer::beginJob(void) { return; }
0139 
0140 void GlobalHitsProducer::endJob() {
0141   std::string MsgLoggerCat = "GlobalHitsProducer_endJob";
0142   if (verbosity >= 0)
0143     edm::LogInfo(MsgLoggerCat) << "Terminating having processed " << count << " events.";
0144   return;
0145 }
0146 
0147 void GlobalHitsProducer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0148   std::string MsgLoggerCat = "GlobalHitsProducer_produce";
0149 
0150   // keep track of number of events processed
0151   ++count;
0152 
0153   // get event id information
0154   edm::RunNumber_t nrun = iEvent.id().run();
0155   edm::EventNumber_t nevt = iEvent.id().event();
0156 
0157   if (verbosity > 0) {
0158     edm::LogInfo(MsgLoggerCat) << "Processing run " << nrun << ", event " << nevt << " (" << count << " events total)";
0159   } else if (verbosity == 0) {
0160     if (nevt % frequency == 0 || nevt == 1) {
0161       edm::LogInfo(MsgLoggerCat) << "Processing run " << nrun << ", event " << nevt << " (" << count
0162                                  << " events total)";
0163     }
0164   }
0165 
0166   // clear event holders
0167   clear();
0168 
0169   // look at information available in the event
0170   if (getAllProvenances) {
0171     std::vector<const edm::StableProvenance *> AllProv;
0172     iEvent.getAllStableProvenance(AllProv);
0173 
0174     if (verbosity >= 0)
0175       edm::LogInfo(MsgLoggerCat) << "Number of Provenances = " << AllProv.size();
0176 
0177     if (printProvenanceInfo && (verbosity >= 0)) {
0178       TString eventout("\nProvenance info:\n");
0179 
0180       for (unsigned int i = 0; i < AllProv.size(); ++i) {
0181         eventout += "\n       ******************************";
0182         eventout += "\n       Module       : ";
0183         eventout += AllProv[i]->moduleLabel();
0184         eventout += "\n       ProductID    : ";
0185         eventout += AllProv[i]->productID().id();
0186         eventout += "\n       ClassName    : ";
0187         eventout += AllProv[i]->className();
0188         eventout += "\n       InstanceName : ";
0189         eventout += AllProv[i]->productInstanceName();
0190         eventout += "\n       BranchName   : ";
0191         eventout += AllProv[i]->branchName();
0192       }
0193       eventout += "\n       ******************************\n";
0194       edm::LogInfo(MsgLoggerCat) << eventout << "\n";
0195       printProvenanceInfo = false;
0196     }
0197     getAllProvenances = false;
0198   }
0199 
0200   // call fill functions
0201   // gather G4MC information from event
0202   fillG4MC(iEvent);
0203   // gather Tracker information from event
0204   fillTrk(iEvent, iSetup);
0205   // gather muon information from event
0206   fillMuon(iEvent, iSetup);
0207   // gather Ecal information from event
0208   fillECal(iEvent, iSetup);
0209   // gather Hcal information from event
0210   fillHCal(iEvent, iSetup);
0211 
0212   if (verbosity > 0)
0213     edm::LogInfo(MsgLoggerCat) << "Done gathering data from event.";
0214 
0215   // produce object to put into event
0216   std::unique_ptr<PGlobalSimHit> pOut(new PGlobalSimHit);
0217 
0218   if (verbosity > 2)
0219     edm::LogInfo(MsgLoggerCat) << "Saving event contents:";
0220 
0221   // call store functions
0222   // store G4MC information in product
0223   storeG4MC(*pOut);
0224   // store Tracker information in produce
0225   storeTrk(*pOut);
0226   // store Muon information in produce
0227   storeMuon(*pOut);
0228   // store ECal information in produce
0229   storeECal(*pOut);
0230   // store HCal information in produce
0231   storeHCal(*pOut);
0232 
0233   // store information in event
0234   iEvent.put(std::move(pOut), label);
0235 
0236   return;
0237 }
0238 
0239 //==================fill and store functions================================
0240 void GlobalHitsProducer::fillG4MC(edm::Event &iEvent) {
0241   std::string MsgLoggerCat = "GlobalHitsProducer_fillG4MC";
0242 
0243   TString eventout;
0244   if (verbosity > 0)
0245     eventout = "\nGathering info:";
0246 
0247   //////////////////////
0248   // get MC information
0249   /////////////////////
0250   edm::Handle<edm::HepMCProduct> HepMCEvt;
0251   std::vector<edm::Handle<edm::HepMCProduct>> AllHepMCEvt;
0252 
0253   //iEvent.getManyByType(AllHepMCEvt);
0254   throw cms::Exception("UnsupportedFunction") << "GlobalHitsProducer::fillG4MC: "
0255                                               << "getManyByType has not been supported by the Framework since 2015. "
0256                                               << "This module has been broken since then. Maybe it should be deleted. "
0257                                               << "Another possibility is to upgrade to use GetterOfProducts instead.";
0258 
0259   // loop through products and extract VtxSmearing if available. Any of them
0260   // should have the information needed
0261   for (unsigned int i = 0; i < AllHepMCEvt.size(); ++i) {
0262     HepMCEvt = AllHepMCEvt[i];
0263     if ((HepMCEvt.provenance()->branchDescription()).moduleLabel() == "generatorSmeared")
0264       break;
0265   }
0266 
0267   if (!HepMCEvt.isValid()) {
0268     edm::LogWarning(MsgLoggerCat) << "Unable to find HepMCProduct in event!";
0269     return;
0270   } else {
0271     eventout += "\n          Using HepMCProduct: ";
0272     eventout += (HepMCEvt.provenance()->branchDescription()).moduleLabel();
0273   }
0274   const HepMC::GenEvent *MCEvt = HepMCEvt->GetEvent();
0275   nRawGenPart = MCEvt->particles_size();
0276 
0277   if (verbosity > 1) {
0278     eventout += "\n          Number of Raw Particles collected:......... ";
0279     eventout += nRawGenPart;
0280   }
0281 
0282   ////////////////////////////
0283   // get G4Vertex information
0284   ////////////////////////////
0285   // convert unit stored in SimVertex to mm
0286   float unit = 0.;
0287   if (vtxunit == 0)
0288     unit = 1.;  // already in mm
0289   if (vtxunit == 1)
0290     unit = 10.;  // stored in cm, convert to mm
0291 
0292   edm::Handle<edm::SimVertexContainer> G4VtxContainer;
0293   iEvent.getByToken(G4VtxSrc_Token_, G4VtxContainer);
0294   if (!G4VtxContainer.isValid()) {
0295     edm::LogWarning(MsgLoggerCat) << "Unable to find SimVertex in event!";
0296     return;
0297   }
0298   int i = 0;
0299   edm::SimVertexContainer::const_iterator itVtx;
0300   for (itVtx = G4VtxContainer->begin(); itVtx != G4VtxContainer->end(); ++itVtx) {
0301     ++i;
0302 
0303     const math::XYZTLorentzVector G4Vtx1(
0304         itVtx->position().x(), itVtx->position().y(), itVtx->position().z(), itVtx->position().e());
0305     double G4Vtx[4];
0306     G4Vtx1.GetCoordinates(G4Vtx);
0307 
0308     G4VtxX.push_back((G4Vtx[0] * unit) / micrometer);
0309     G4VtxY.push_back((G4Vtx[1] * unit) / micrometer);
0310     G4VtxZ.push_back((G4Vtx[2] * unit) / millimeter);
0311   }
0312 
0313   if (verbosity > 1) {
0314     eventout += "\n          Number of G4Vertices collected:............ ";
0315     eventout += i;
0316   }
0317 
0318   ///////////////////////////
0319   // get G4Track information
0320   ///////////////////////////
0321   edm::Handle<edm::SimTrackContainer> G4TrkContainer;
0322   iEvent.getByToken(G4TrkSrc_Token_, G4TrkContainer);
0323 
0324   if (!G4TrkContainer.isValid()) {
0325     edm::LogWarning(MsgLoggerCat) << "Unable to find SimTrack in event!";
0326     return;
0327   }
0328   i = 0;
0329   edm::SimTrackContainer::const_iterator itTrk;
0330   for (itTrk = G4TrkContainer->begin(); itTrk != G4TrkContainer->end(); ++itTrk) {
0331     ++i;
0332 
0333     const math::XYZTLorentzVector G4Trk1(
0334         itTrk->momentum().x(), itTrk->momentum().y(), itTrk->momentum().z(), itTrk->momentum().e());
0335     double G4Trk[4];
0336     G4Trk1.GetCoordinates(G4Trk);
0337 
0338     G4TrkPt.push_back(sqrt(G4Trk[0] * G4Trk[0] + G4Trk[1] * G4Trk[1]));  // GeV
0339     G4TrkE.push_back(G4Trk[3]);                                          // GeV
0340   }
0341 
0342   if (verbosity > 1) {
0343     eventout += "\n          Number of G4Tracks collected:.............. ";
0344     eventout += i;
0345   }
0346 
0347   if (verbosity > 0)
0348     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
0349 
0350   return;
0351 }
0352 
0353 void GlobalHitsProducer::storeG4MC(PGlobalSimHit &product) {
0354   std::string MsgLoggerCat = "GlobalHitsProducer_storeG4MC";
0355 
0356   if (verbosity > 2) {
0357     TString eventout("\n       nRawGenPart        = ");
0358     eventout += nRawGenPart;
0359     eventout += "\n       nG4Vtx             = ";
0360     eventout += G4VtxX.size();
0361     for (unsigned int i = 0; i < G4VtxX.size(); ++i) {
0362       eventout += "\n          (x,y,z)         = (";
0363       eventout += G4VtxX[i];
0364       eventout += ", ";
0365       eventout += G4VtxY[i];
0366       eventout += ", ";
0367       eventout += G4VtxZ[i];
0368       eventout += ")";
0369     }
0370     eventout += "\n       nG4Trk             = ";
0371     eventout += G4TrkPt.size();
0372     for (unsigned int i = 0; i < G4TrkPt.size(); ++i) {
0373       eventout += "\n          (pt,e)          = (";
0374       eventout += G4TrkPt[i];
0375       eventout += ", ";
0376       eventout += G4TrkE[i];
0377       eventout += ")";
0378     }
0379     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
0380   }  // end verbose output
0381 
0382   product.putRawGenPart(nRawGenPart);
0383   product.putG4Vtx(G4VtxX, G4VtxY, G4VtxZ);
0384   product.putG4Trk(G4TrkPt, G4TrkE);
0385 
0386   return;
0387 }
0388 
0389 void GlobalHitsProducer::fillTrk(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0390   std::string MsgLoggerCat = "GlobalHitsProducer_fillTrk";
0391 
0392   TString eventout;
0393   if (verbosity > 0)
0394     eventout = "\nGathering info:";
0395 
0396   // access the tracker geometry
0397   const auto &theTrackerGeometry = iSetup.getHandle(tGeomToken_);
0398   if (!theTrackerGeometry.isValid()) {
0399     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerDigiGeometryRecord in event!";
0400     return;
0401   }
0402   const TrackerGeometry &theTracker(*theTrackerGeometry);
0403 
0404   // iterator to access containers
0405   edm::PSimHitContainer::const_iterator itHit;
0406 
0407   ///////////////////////////////
0408   // get Pixel Barrel information
0409   ///////////////////////////////
0410   edm::PSimHitContainer thePxlBrlHits;
0411   // extract low container
0412   edm::Handle<edm::PSimHitContainer> PxlBrlLowContainer;
0413   iEvent.getByToken(PxlBrlLowSrc_Token_, PxlBrlLowContainer);
0414   if (!PxlBrlLowContainer.isValid()) {
0415     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsPixelBarrelLowTof in event!";
0416     return;
0417   }
0418   // extract high container
0419   edm::Handle<edm::PSimHitContainer> PxlBrlHighContainer;
0420   iEvent.getByToken(PxlBrlHighSrc_Token_, PxlBrlHighContainer);
0421   if (!PxlBrlHighContainer.isValid()) {
0422     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsPixelBarrelHighTof in event!";
0423     return;
0424   }
0425   // place both containers into new container
0426   thePxlBrlHits.insert(thePxlBrlHits.end(), PxlBrlLowContainer->begin(), PxlBrlLowContainer->end());
0427   thePxlBrlHits.insert(thePxlBrlHits.end(), PxlBrlHighContainer->begin(), PxlBrlHighContainer->end());
0428 
0429   // cycle through new container
0430   int i = 0, j = 0;
0431   for (itHit = thePxlBrlHits.begin(); itHit != thePxlBrlHits.end(); ++itHit) {
0432     ++i;
0433 
0434     // create a DetId from the detUnitId
0435     DetId theDetUnitId(itHit->detUnitId());
0436     int detector = theDetUnitId.det();
0437     int subdetector = theDetUnitId.subdetId();
0438 
0439     // check that expected detector is returned
0440     if ((detector == dTrk) && (subdetector == sdPxlBrl)) {
0441       // get the GeomDetUnit from the geometry using theDetUnitID
0442       const GeomDetUnit *theDet = theTracker.idToDetUnit(theDetUnitId);
0443 
0444       if (!theDet) {
0445         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from PxlBrlHits for Hit " << i;
0446         continue;
0447       }
0448 
0449       ++j;
0450 
0451       // get the Surface of the hit (knows how to go from local <-> global)
0452       const BoundPlane &bSurface = theDet->surface();
0453 
0454       // gather necessary information
0455       PxlBrlToF.push_back(itHit->tof());
0456       PxlBrlR.push_back(bSurface.toGlobal(itHit->localPosition()).perp());
0457       PxlBrlPhi.push_back(bSurface.toGlobal(itHit->localPosition()).phi());
0458       PxlBrlEta.push_back(bSurface.toGlobal(itHit->localPosition()).eta());
0459 
0460     } else {
0461       edm::LogWarning(MsgLoggerCat) << "PxlBrl PSimHit " << i << " is expected to be (det,subdet) = (" << dTrk << ","
0462                                     << sdPxlBrl << "); value returned is: (" << detector << "," << subdetector << ")";
0463       continue;
0464     }  // end detector type check
0465   }    // end loop through PxlBrl Hits
0466 
0467   if (verbosity > 1) {
0468     eventout += "\n          Number of Pixel Barrel Hits collected:..... ";
0469     eventout += j;
0470   }
0471 
0472   /////////////////////////////////
0473   // get Pixel Forward information
0474   ////////////////////////////////
0475   edm::PSimHitContainer thePxlFwdHits;
0476   // extract low container
0477   edm::Handle<edm::PSimHitContainer> PxlFwdLowContainer;
0478   iEvent.getByToken(PxlFwdLowSrc_Token_, PxlFwdLowContainer);
0479   if (!PxlFwdLowContainer.isValid()) {
0480     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsPixelEndcapLowTof in event!";
0481     return;
0482   }
0483   // extract high container
0484   edm::Handle<edm::PSimHitContainer> PxlFwdHighContainer;
0485   iEvent.getByToken(PxlFwdHighSrc_Token_, PxlFwdHighContainer);
0486   if (!PxlFwdHighContainer.isValid()) {
0487     edm::LogWarning("GlobalHitsProducer_fillTrk") << "Unable to find TrackerHitsPixelEndcapHighTof in event!";
0488     return;
0489   }
0490   // place both containers into new container
0491   thePxlFwdHits.insert(thePxlFwdHits.end(), PxlFwdLowContainer->begin(), PxlFwdLowContainer->end());
0492   thePxlFwdHits.insert(thePxlFwdHits.end(), PxlFwdHighContainer->begin(), PxlFwdHighContainer->end());
0493 
0494   // cycle through new container
0495   i = 0;
0496   j = 0;
0497   for (itHit = thePxlFwdHits.begin(); itHit != thePxlFwdHits.end(); ++itHit) {
0498     ++i;
0499 
0500     // create a DetId from the detUnitId
0501     DetId theDetUnitId(itHit->detUnitId());
0502     int detector = theDetUnitId.det();
0503     int subdetector = theDetUnitId.subdetId();
0504 
0505     // check that expected detector is returned
0506     if ((detector == dTrk) && (subdetector == sdPxlFwd)) {
0507       // get the GeomDetUnit from the geometry using theDetUnitID
0508       const GeomDetUnit *theDet = theTracker.idToDetUnit(theDetUnitId);
0509 
0510       if (!theDet) {
0511         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from PxlFwdHits for Hit " << i;
0512         ;
0513         continue;
0514       }
0515 
0516       ++j;
0517 
0518       // get the Surface of the hit (knows how to go from local <-> global)
0519       const BoundPlane &bSurface = theDet->surface();
0520 
0521       // gather necessary information
0522       PxlFwdToF.push_back(itHit->tof());
0523       PxlFwdZ.push_back(bSurface.toGlobal(itHit->localPosition()).z());
0524       PxlFwdPhi.push_back(bSurface.toGlobal(itHit->localPosition()).phi());
0525       PxlFwdEta.push_back(bSurface.toGlobal(itHit->localPosition()).eta());
0526     } else {
0527       edm::LogWarning(MsgLoggerCat) << "PxlFwd PSimHit " << i << " is expected to be (det,subdet) = (" << dTrk << ","
0528                                     << sdPxlFwd << "); value returned is: (" << detector << "," << subdetector << ")";
0529       continue;
0530     }  // end detector type check
0531   }    // end loop through PxlFwd Hits
0532 
0533   if (verbosity > 1) {
0534     eventout += "\n          Number of Pixel Forward Hits collected:.... ";
0535     eventout += j;
0536   }
0537 
0538   ///////////////////////////////////
0539   // get Silicon Barrel information
0540   //////////////////////////////////
0541   edm::PSimHitContainer theSiBrlHits;
0542   // extract TIB low container
0543   edm::Handle<edm::PSimHitContainer> SiTIBLowContainer;
0544   iEvent.getByToken(SiTIBLowSrc_Token_, SiTIBLowContainer);
0545   if (!SiTIBLowContainer.isValid()) {
0546     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTIBLowTof in event!";
0547     return;
0548   }
0549   // extract TIB high container
0550   edm::Handle<edm::PSimHitContainer> SiTIBHighContainer;
0551   iEvent.getByToken(SiTIBHighSrc_Token_, SiTIBHighContainer);
0552   if (!SiTIBHighContainer.isValid()) {
0553     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTIBHighTof in event!";
0554     return;
0555   }
0556   // extract TOB low container
0557   edm::Handle<edm::PSimHitContainer> SiTOBLowContainer;
0558   iEvent.getByToken(SiTOBLowSrc_Token_, SiTOBLowContainer);
0559   if (!SiTOBLowContainer.isValid()) {
0560     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTOBLowTof in event!";
0561     return;
0562   }
0563   // extract TOB high container
0564   edm::Handle<edm::PSimHitContainer> SiTOBHighContainer;
0565   iEvent.getByToken(SiTOBHighSrc_Token_, SiTOBHighContainer);
0566   if (!SiTOBHighContainer.isValid()) {
0567     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTOBHighTof in event!";
0568     return;
0569   }
0570   // place all containers into new container
0571   theSiBrlHits.insert(theSiBrlHits.end(), SiTIBLowContainer->begin(), SiTIBLowContainer->end());
0572   theSiBrlHits.insert(theSiBrlHits.end(), SiTIBHighContainer->begin(), SiTIBHighContainer->end());
0573   theSiBrlHits.insert(theSiBrlHits.end(), SiTOBLowContainer->begin(), SiTOBLowContainer->end());
0574   theSiBrlHits.insert(theSiBrlHits.end(), SiTOBHighContainer->begin(), SiTOBHighContainer->end());
0575 
0576   // cycle through new container
0577   i = 0;
0578   j = 0;
0579   for (itHit = theSiBrlHits.begin(); itHit != theSiBrlHits.end(); ++itHit) {
0580     ++i;
0581 
0582     // create a DetId from the detUnitId
0583     DetId theDetUnitId(itHit->detUnitId());
0584     int detector = theDetUnitId.det();
0585     int subdetector = theDetUnitId.subdetId();
0586 
0587     // check that expected detector is returned
0588     if ((detector == dTrk) && ((subdetector == sdSiTIB) || (subdetector == sdSiTOB))) {
0589       // get the GeomDetUnit from the geometry using theDetUnitID
0590       const GeomDetUnit *theDet = theTracker.idToDetUnit(theDetUnitId);
0591 
0592       if (!theDet) {
0593         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from SiBrlHits for Hit " << i;
0594         continue;
0595       }
0596 
0597       ++j;
0598 
0599       // get the Surface of the hit (knows how to go from local <-> global)
0600       const BoundPlane &bSurface = theDet->surface();
0601 
0602       // gather necessary information
0603       SiBrlToF.push_back(itHit->tof());
0604       SiBrlR.push_back(bSurface.toGlobal(itHit->localPosition()).perp());
0605       SiBrlPhi.push_back(bSurface.toGlobal(itHit->localPosition()).phi());
0606       SiBrlEta.push_back(bSurface.toGlobal(itHit->localPosition()).eta());
0607     } else {
0608       edm::LogWarning(MsgLoggerCat) << "SiBrl PSimHit " << i << " is expected to be (det,subdet) = (" << dTrk << ","
0609                                     << sdSiTIB << " || " << sdSiTOB << "); value returned is: (" << detector << ","
0610                                     << subdetector << ")";
0611       continue;
0612     }  // end detector type check
0613   }    // end loop through SiBrl Hits
0614 
0615   if (verbosity > 1) {
0616     eventout += "\n          Number of Silicon Barrel Hits collected:... ";
0617     eventout += j;
0618   }
0619 
0620   ///////////////////////////////////
0621   // get Silicon Forward information
0622   ///////////////////////////////////
0623   edm::PSimHitContainer theSiFwdHits;
0624   // extract TID low container
0625   edm::Handle<edm::PSimHitContainer> SiTIDLowContainer;
0626   iEvent.getByToken(SiTIDLowSrc_Token_, SiTIDLowContainer);
0627   if (!SiTIDLowContainer.isValid()) {
0628     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTIDLowTof in event!";
0629     return;
0630   }
0631   // extract TID high container
0632   edm::Handle<edm::PSimHitContainer> SiTIDHighContainer;
0633   iEvent.getByToken(SiTIDHighSrc_Token_, SiTIDHighContainer);
0634   if (!SiTIDHighContainer.isValid()) {
0635     edm::LogWarning("GlobalHitsProducer_fillTrk") << "Unable to find TrackerHitsTIDHighTof in event!";
0636     return;
0637   }
0638   // extract TEC low container
0639   edm::Handle<edm::PSimHitContainer> SiTECLowContainer;
0640   iEvent.getByToken(SiTECLowSrc_Token_, SiTECLowContainer);
0641   if (!SiTECLowContainer.isValid()) {
0642     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTECLowTof in event!";
0643     return;
0644   }
0645   // extract TEC high container
0646   edm::Handle<edm::PSimHitContainer> SiTECHighContainer;
0647   iEvent.getByToken(SiTECHighSrc_Token_, SiTECHighContainer);
0648   if (!SiTECHighContainer.isValid()) {
0649     edm::LogWarning(MsgLoggerCat) << "Unable to find TrackerHitsTECHighTof in event!";
0650     return;
0651   }
0652   // place all containers into new container
0653   theSiFwdHits.insert(theSiFwdHits.end(), SiTIDLowContainer->begin(), SiTIDLowContainer->end());
0654   theSiFwdHits.insert(theSiFwdHits.end(), SiTIDHighContainer->begin(), SiTIDHighContainer->end());
0655   theSiFwdHits.insert(theSiFwdHits.end(), SiTECLowContainer->begin(), SiTECLowContainer->end());
0656   theSiFwdHits.insert(theSiFwdHits.end(), SiTECHighContainer->begin(), SiTECHighContainer->end());
0657 
0658   // cycle through container
0659   i = 0;
0660   j = 0;
0661   for (itHit = theSiFwdHits.begin(); itHit != theSiFwdHits.end(); ++itHit) {
0662     ++i;
0663 
0664     // create a DetId from the detUnitId
0665     DetId theDetUnitId(itHit->detUnitId());
0666     int detector = theDetUnitId.det();
0667     int subdetector = theDetUnitId.subdetId();
0668 
0669     // check that expected detector is returned
0670     if ((detector == dTrk) && ((subdetector == sdSiTID) || (subdetector == sdSiTEC))) {
0671       // get the GeomDetUnit from the geometry using theDetUnitID
0672       const GeomDetUnit *theDet = theTracker.idToDetUnit(theDetUnitId);
0673 
0674       if (!theDet) {
0675         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from SiFwdHits Hit " << i;
0676         return;
0677       }
0678 
0679       ++j;
0680 
0681       // get the Surface of the hit (knows how to go from local <-> global)
0682       const BoundPlane &bSurface = theDet->surface();
0683 
0684       // gather necessary information
0685       SiFwdToF.push_back(itHit->tof());
0686       SiFwdZ.push_back(bSurface.toGlobal(itHit->localPosition()).z());
0687       SiFwdPhi.push_back(bSurface.toGlobal(itHit->localPosition()).phi());
0688       SiFwdEta.push_back(bSurface.toGlobal(itHit->localPosition()).eta());
0689     } else {
0690       edm::LogWarning(MsgLoggerCat) << "SiFwd PSimHit " << i << " is expected to be (det,subdet) = (" << dTrk << ","
0691                                     << sdSiTOB << " || " << sdSiTEC << "); value returned is: (" << detector << ","
0692                                     << subdetector << ")";
0693       continue;
0694     }  // end check detector type
0695   }    // end loop through SiFwd Hits
0696 
0697   if (verbosity > 1) {
0698     eventout += "\n          Number of Silicon Forward Hits collected:.. ";
0699     eventout += j;
0700   }
0701 
0702   if (verbosity > 0)
0703     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
0704 
0705   return;
0706 }
0707 
0708 void GlobalHitsProducer::storeTrk(PGlobalSimHit &product) {
0709   std::string MsgLoggerCat = "GlobalHitsProducer_storeTrk";
0710 
0711   if (verbosity > 2) {
0712     TString eventout("\n       nPxlBrlHits        = ");
0713     eventout += PxlBrlToF.size();
0714     for (unsigned int i = 0; i < PxlBrlToF.size(); ++i) {
0715       eventout += "\n          (tof,r,phi,eta) = (";
0716       eventout += PxlBrlToF[i];
0717       eventout += ", ";
0718       eventout += PxlBrlR[i];
0719       eventout += ", ";
0720       eventout += PxlBrlPhi[i];
0721       eventout += ", ";
0722       eventout += PxlBrlEta[i];
0723       eventout += ")";
0724     }  // end PxlBrl output
0725     eventout += "\n       nPxlFwdHits        = ";
0726     eventout += PxlFwdToF.size();
0727     for (unsigned int i = 0; i < PxlFwdToF.size(); ++i) {
0728       eventout += "\n          (tof,z,phi,eta) = (";
0729       eventout += PxlFwdToF[i];
0730       eventout += ", ";
0731       eventout += PxlFwdZ[i];
0732       eventout += ", ";
0733       eventout += PxlFwdPhi[i];
0734       eventout += ", ";
0735       eventout += PxlFwdEta[i];
0736       eventout += ")";
0737     }  // end PxlFwd output
0738     eventout += "\n       nSiBrlHits         = ";
0739     eventout += SiBrlToF.size();
0740     for (unsigned int i = 0; i < SiBrlToF.size(); ++i) {
0741       eventout += "\n          (tof,r,phi,eta) = (";
0742       eventout += SiBrlToF[i];
0743       eventout += ", ";
0744       eventout += SiBrlR[i];
0745       eventout += ", ";
0746       eventout += SiBrlPhi[i];
0747       eventout += ", ";
0748       eventout += SiBrlEta[i];
0749       eventout += ")";
0750     }  // end SiBrl output
0751     eventout += "\n       nSiFwdHits         = ";
0752     eventout += SiFwdToF.size();
0753     for (unsigned int i = 0; i < SiFwdToF.size(); ++i) {
0754       eventout += "\n          (tof,z,phi,eta) = (";
0755       eventout += SiFwdToF[i];
0756       eventout += ", ";
0757       eventout += SiFwdZ[i];
0758       eventout += ", ";
0759       eventout += SiFwdPhi[i];
0760       eventout += ", ";
0761       eventout += SiFwdEta[i];
0762       eventout += ")";
0763     }  // end SiFwd output
0764     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
0765   }  // end verbose output
0766 
0767   product.putPxlBrlHits(PxlBrlToF, PxlBrlR, PxlBrlPhi, PxlBrlEta);
0768   product.putPxlFwdHits(PxlFwdToF, PxlFwdZ, PxlFwdPhi, PxlFwdEta);
0769   product.putSiBrlHits(SiBrlToF, SiBrlR, SiBrlPhi, SiBrlEta);
0770   product.putSiFwdHits(SiFwdToF, SiFwdZ, SiFwdPhi, SiFwdEta);
0771 
0772   return;
0773 }
0774 
0775 void GlobalHitsProducer::fillMuon(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0776   std::string MsgLoggerCat = "GlobalHitsProducer_fillMuon";
0777 
0778   TString eventout;
0779   if (verbosity > 0)
0780     eventout = "\nGathering info:";
0781 
0782   // iterator to access containers
0783   edm::PSimHitContainer::const_iterator itHit;
0784 
0785   // int i = 0, j = 0;
0786   ///////////////////////
0787   // access the CSC Muon
0788   ///////////////////////
0789   // access the CSC Muon geometry
0790   const auto &theCSCGeometry = iSetup.getHandle(cscGeomToken_);
0791   if (!theCSCGeometry.isValid()) {
0792     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonGeometryRecord for the CSCGeometry in event!";
0793     return;
0794   }
0795   const CSCGeometry &theCSCMuon(*theCSCGeometry);
0796 
0797   // get Muon CSC information
0798   edm::Handle<edm::PSimHitContainer> MuonCSCContainer;
0799   iEvent.getByToken(MuonCscSrc_Token_, MuonCSCContainer);
0800   if (!MuonCSCContainer.isValid()) {
0801     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonCSCHits in event!";
0802     return;
0803   }
0804 
0805   // cycle through container
0806   int i = 0, j = 0;
0807   for (itHit = MuonCSCContainer->begin(); itHit != MuonCSCContainer->end(); ++itHit) {
0808     ++i;
0809 
0810     // create a DetId from the detUnitId
0811     DetId theDetUnitId(itHit->detUnitId());
0812     int detector = theDetUnitId.det();
0813     int subdetector = theDetUnitId.subdetId();
0814 
0815     // check that expected detector is returned
0816     if ((detector == dMuon) && (subdetector == sdMuonCSC)) {
0817       // get the GeomDetUnit from the geometry using theDetUnitID
0818       const GeomDetUnit *theDet = theCSCMuon.idToDetUnit(theDetUnitId);
0819 
0820       if (!theDet) {
0821         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from theCSCMuon for hit " << i;
0822         continue;
0823       }
0824 
0825       ++j;
0826 
0827       // get the Surface of the hit (knows how to go from local <-> global)
0828       const BoundPlane &bSurface = theDet->surface();
0829 
0830       // gather necessary information
0831       MuonCscToF.push_back(itHit->tof());
0832       MuonCscZ.push_back(bSurface.toGlobal(itHit->localPosition()).z());
0833       MuonCscPhi.push_back(bSurface.toGlobal(itHit->localPosition()).phi());
0834       MuonCscEta.push_back(bSurface.toGlobal(itHit->localPosition()).eta());
0835     } else {
0836       edm::LogWarning(MsgLoggerCat) << "MuonCsc PSimHit " << i << " is expected to be (det,subdet) = (" << dMuon << ","
0837                                     << sdMuonCSC << "); value returned is: (" << detector << "," << subdetector << ")";
0838       continue;
0839     }  // end detector type check
0840   }    // end loop through CSC Hits
0841 
0842   if (verbosity > 1) {
0843     eventout += "\n          Number of CSC muon Hits collected:......... ";
0844     eventout += j;
0845   }
0846 
0847   // i = 0, j = 0;
0848   /////////////////////
0849   // access the DT Muon
0850   /////////////////////
0851   // access the DT Muon geometry
0852   const auto &theDTGeometry = iSetup.getHandle(dtGeomToken_);
0853   if (!theDTGeometry.isValid()) {
0854     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonGeometryRecord for the DTGeometry in event!";
0855     return;
0856   }
0857   const DTGeometry &theDTMuon(*theDTGeometry);
0858 
0859   // get Muon DT information
0860   edm::Handle<edm::PSimHitContainer> MuonDtContainer;
0861   iEvent.getByToken(MuonDtSrc_Token_, MuonDtContainer);
0862   if (!MuonDtContainer.isValid()) {
0863     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonDTHits in event!";
0864     return;
0865   }
0866 
0867   // cycle through container
0868   i = 0, j = 0;
0869   for (itHit = MuonDtContainer->begin(); itHit != MuonDtContainer->end(); ++itHit) {
0870     ++i;
0871 
0872     // create a DetId from the detUnitId
0873     DetId theDetUnitId(itHit->detUnitId());
0874     int detector = theDetUnitId.det();
0875     int subdetector = theDetUnitId.subdetId();
0876 
0877     // check that expected detector is returned
0878     if ((detector == dMuon) && (subdetector == sdMuonDT)) {
0879       // CSC uses wires and layers rather than the full detID
0880       // get the wireId
0881       DTWireId wireId(itHit->detUnitId());
0882 
0883       // get the DTLayer from the geometry using the wireID
0884       const DTLayer *theDet = theDTMuon.layer(wireId.layerId());
0885 
0886       if (!theDet) {
0887         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from theDtMuon for hit " << i;
0888         continue;
0889       }
0890 
0891       ++j;
0892 
0893       // get the Surface of the hit (knows how to go from local <-> global)
0894       const BoundPlane &bSurface = theDet->surface();
0895 
0896       // gather necessary information
0897       MuonDtToF.push_back(itHit->tof());
0898       MuonDtR.push_back(bSurface.toGlobal(itHit->localPosition()).perp());
0899       MuonDtPhi.push_back(bSurface.toGlobal(itHit->localPosition()).phi());
0900       MuonDtEta.push_back(bSurface.toGlobal(itHit->localPosition()).eta());
0901     } else {
0902       edm::LogWarning(MsgLoggerCat) << "MuonDt PSimHit " << i << " is expected to be (det,subdet) = (" << dMuon << ","
0903                                     << sdMuonDT << "); value returned is: (" << detector << "," << subdetector << ")";
0904       continue;
0905     }  // end detector type check
0906   }    // end loop through DT Hits
0907 
0908   if (verbosity > 1) {
0909     eventout += "\n          Number of DT muon Hits collected:.......... ";
0910     eventout += j;
0911   }
0912 
0913   // i = 0, j = 0;
0914   // int RPCBrl = 0, RPCFwd = 0;
0915   ///////////////////////
0916   // access the RPC Muon
0917   ///////////////////////
0918   // access the RPC Muon geometry
0919   const auto &theRPCGeometry = iSetup.getHandle(rpcGeomToken_);
0920   if (!theRPCGeometry.isValid()) {
0921     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonGeometryRecord for the RPCGeometry in event!";
0922     return;
0923   }
0924   const RPCGeometry &theRPCMuon(*theRPCGeometry);
0925 
0926   // get Muon RPC information
0927   edm::Handle<edm::PSimHitContainer> MuonRPCContainer;
0928   iEvent.getByToken(MuonRpcSrc_Token_, MuonRPCContainer);
0929   if (!MuonRPCContainer.isValid()) {
0930     edm::LogWarning(MsgLoggerCat) << "Unable to find MuonRPCHits in event!";
0931     return;
0932   }
0933 
0934   // cycle through container
0935   i = 0, j = 0;
0936   int RPCBrl = 0, RPCFwd = 0;
0937   for (itHit = MuonRPCContainer->begin(); itHit != MuonRPCContainer->end(); ++itHit) {
0938     ++i;
0939 
0940     // create a DetID from the detUnitId
0941     DetId theDetUnitId(itHit->detUnitId());
0942     int detector = theDetUnitId.det();
0943     int subdetector = theDetUnitId.subdetId();
0944 
0945     // check that expected detector is returned
0946     if ((detector == dMuon) && (subdetector == sdMuonRPC)) {
0947       // get an RPCDetID from the detUnitID
0948       RPCDetId RPCId(itHit->detUnitId());
0949 
0950       // find the region of the RPC hit
0951       int region = RPCId.region();
0952 
0953       // get the GeomDetUnit from the geometry using the RPCDetId
0954       const GeomDetUnit *theDet = theRPCMuon.idToDetUnit(theDetUnitId);
0955 
0956       if (!theDet) {
0957         edm::LogWarning(MsgLoggerCat) << "Unable to get GeomDetUnit from theRPCMuon for hit " << i;
0958         continue;
0959       }
0960 
0961       ++j;
0962 
0963       // get the Surface of the hit (knows how to go from local <-> global)
0964       const BoundPlane &bSurface = theDet->surface();
0965 
0966       // gather necessary information
0967       if ((region == sdMuonRPCRgnFwdp) || (region == sdMuonRPCRgnFwdn)) {
0968         ++RPCFwd;
0969 
0970         MuonRpcFwdToF.push_back(itHit->tof());
0971         MuonRpcFwdZ.push_back(bSurface.toGlobal(itHit->localPosition()).z());
0972         MuonRpcFwdPhi.push_back(bSurface.toGlobal(itHit->localPosition()).phi());
0973         MuonRpcFwdEta.push_back(bSurface.toGlobal(itHit->localPosition()).eta());
0974       } else if (region == sdMuonRPCRgnBrl) {
0975         ++RPCBrl;
0976 
0977         MuonRpcBrlToF.push_back(itHit->tof());
0978         MuonRpcBrlR.push_back(bSurface.toGlobal(itHit->localPosition()).perp());
0979         MuonRpcBrlPhi.push_back(bSurface.toGlobal(itHit->localPosition()).phi());
0980         MuonRpcBrlEta.push_back(bSurface.toGlobal(itHit->localPosition()).eta());
0981       } else {
0982         edm::LogWarning(MsgLoggerCat) << "Invalid region for RPC Muon hit" << i;
0983         continue;
0984       }  // end check of region
0985     } else {
0986       edm::LogWarning(MsgLoggerCat) << "MuonRpc PSimHit " << i << " is expected to be (det,subdet) = (" << dMuon << ","
0987                                     << sdMuonRPC << "); value returned is: (" << detector << "," << subdetector << ")";
0988       continue;
0989     }  // end detector type check
0990   }    // end loop through RPC Hits
0991 
0992   if (verbosity > 1) {
0993     eventout += "\n          Number of RPC muon Hits collected:......... ";
0994     eventout += j;
0995     eventout += "\n                    RPC Barrel muon Hits:............ ";
0996     eventout += RPCBrl;
0997     eventout += "\n                    RPC Forward muon Hits:........... ";
0998     eventout += RPCFwd;
0999   }
1000 
1001   if (verbosity > 0)
1002     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
1003 
1004   return;
1005 }
1006 
1007 void GlobalHitsProducer::storeMuon(PGlobalSimHit &product) {
1008   std::string MsgLoggerCat = "GlobalHitsProducer_storeMuon";
1009 
1010   if (verbosity > 2) {
1011     TString eventout("\n       nMuonCSCHits       = ");
1012     eventout += MuonCscToF.size();
1013     for (unsigned int i = 0; i < MuonCscToF.size(); ++i) {
1014       eventout += "\n          (tof,z,phi,eta) = (";
1015       eventout += MuonCscToF[i];
1016       eventout += ", ";
1017       eventout += MuonCscZ[i];
1018       eventout += ", ";
1019       eventout += MuonCscPhi[i];
1020       eventout += ", ";
1021       eventout += MuonCscEta[i];
1022       eventout += ")";
1023     }  // end MuonCsc output
1024     eventout += "\n       nMuonDtHits        = ";
1025     eventout += MuonDtToF.size();
1026     for (unsigned int i = 0; i < MuonDtToF.size(); ++i) {
1027       eventout += "\n          (tof,r,phi,eta) = (";
1028       eventout += MuonDtToF[i];
1029       eventout += ", ";
1030       eventout += MuonDtR[i];
1031       eventout += ", ";
1032       eventout += MuonDtPhi[i];
1033       eventout += ", ";
1034       eventout += MuonDtEta[i];
1035       eventout += ")";
1036     }  // end MuonDt output
1037     eventout += "\n       nMuonRpcBrlHits    = ";
1038     eventout += MuonRpcBrlToF.size();
1039     for (unsigned int i = 0; i < MuonRpcBrlToF.size(); ++i) {
1040       eventout += "\n          (tof,r,phi,eta) = (";
1041       eventout += MuonRpcBrlToF[i];
1042       eventout += ", ";
1043       eventout += MuonRpcBrlR[i];
1044       eventout += ", ";
1045       eventout += MuonRpcBrlPhi[i];
1046       eventout += ", ";
1047       eventout += MuonRpcBrlEta[i];
1048       eventout += ")";
1049     }  // end MuonRpcBrl output
1050     eventout += "\n       nMuonRpcFwdHits    = ";
1051     eventout += MuonRpcFwdToF.size();
1052     for (unsigned int i = 0; i < MuonRpcFwdToF.size(); ++i) {
1053       eventout += "\n          (tof,z,phi,eta) = (";
1054       eventout += MuonRpcFwdToF[i];
1055       eventout += ", ";
1056       eventout += MuonRpcFwdZ[i];
1057       eventout += ", ";
1058       eventout += MuonRpcFwdPhi[i];
1059       eventout += ", ";
1060       eventout += MuonRpcFwdEta[i];
1061       eventout += ")";
1062     }  // end MuonRpcFwd output
1063     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
1064   }  // end verbose output
1065 
1066   product.putMuonCscHits(MuonCscToF, MuonCscZ, MuonCscPhi, MuonCscEta);
1067   product.putMuonDtHits(MuonDtToF, MuonDtR, MuonDtPhi, MuonDtEta);
1068   product.putMuonRpcBrlHits(MuonRpcBrlToF, MuonRpcBrlR, MuonRpcBrlPhi, MuonRpcBrlEta);
1069   product.putMuonRpcFwdHits(MuonRpcFwdToF, MuonRpcFwdZ, MuonRpcFwdPhi, MuonRpcFwdEta);
1070 
1071   return;
1072 }
1073 
1074 void GlobalHitsProducer::fillECal(edm::Event &iEvent, const edm::EventSetup &iSetup) {
1075   std::string MsgLoggerCat = "GlobalHitsProducer_fillECal";
1076 
1077   TString eventout;
1078   if (verbosity > 0)
1079     eventout = "\nGathering info:";
1080 
1081   // access the calorimeter geometry
1082   const auto &theCaloGeometry = iSetup.getHandle(caloGeomToken_);
1083   if (!theCaloGeometry.isValid()) {
1084     edm::LogWarning(MsgLoggerCat) << "Unable to find CaloGeometryRecord in event!";
1085     return;
1086   }
1087   const CaloGeometry &theCalo(*theCaloGeometry);
1088 
1089   // iterator to access containers
1090   edm::PCaloHitContainer::const_iterator itHit;
1091 
1092   ///////////////////////////////
1093   // get  ECal information
1094   ///////////////////////////////
1095   edm::PCaloHitContainer theECalHits;
1096   // extract EB container
1097   edm::Handle<edm::PCaloHitContainer> EBContainer;
1098   iEvent.getByToken(ECalEBSrc_Token_, EBContainer);
1099   if (!EBContainer.isValid()) {
1100     edm::LogWarning(MsgLoggerCat) << "Unable to find EcalHitsEB in event!";
1101     return;
1102   }
1103   // extract EE container
1104   edm::Handle<edm::PCaloHitContainer> EEContainer;
1105   iEvent.getByToken(ECalEESrc_Token_, EEContainer);
1106   if (!EEContainer.isValid()) {
1107     edm::LogWarning(MsgLoggerCat) << "Unable to find EcalHitsEE in event!";
1108     return;
1109   }
1110   // place both containers into new container
1111   theECalHits.insert(theECalHits.end(), EBContainer->begin(), EBContainer->end());
1112   theECalHits.insert(theECalHits.end(), EEContainer->begin(), EEContainer->end());
1113 
1114   // cycle through new container
1115   int i = 0, j = 0;
1116   for (itHit = theECalHits.begin(); itHit != theECalHits.end(); ++itHit) {
1117     ++i;
1118 
1119     // create a DetId from the detUnitId
1120     DetId theDetUnitId(itHit->id());
1121     int detector = theDetUnitId.det();
1122     int subdetector = theDetUnitId.subdetId();
1123 
1124     // check that expected detector is returned
1125     if ((detector == dEcal) && ((subdetector == sdEcalBrl) || (subdetector == sdEcalFwd))) {
1126       // get the Cell geometry
1127       auto theDet = (theCalo.getSubdetectorGeometry(theDetUnitId))->getGeometry(theDetUnitId);
1128 
1129       if (!theDet) {
1130         edm::LogWarning(MsgLoggerCat) << "Unable to get CaloCellGeometry from ECalHits for Hit " << i;
1131         continue;
1132       }
1133 
1134       ++j;
1135 
1136       // get the global position of the cell
1137       const GlobalPoint &globalposition = theDet->getPosition();
1138 
1139       // gather necessary information
1140       ECalE.push_back(itHit->energy());
1141       ECalToF.push_back(itHit->time());
1142       ECalPhi.push_back(globalposition.phi());
1143       ECalEta.push_back(globalposition.eta());
1144 
1145     } else {
1146       edm::LogWarning(MsgLoggerCat) << "ECal PCaloHit " << i << " is expected to be (det,subdet) = (" << dEcal << ","
1147                                     << sdEcalBrl << " || " << sdEcalFwd << "); value returned is: (" << detector << ","
1148                                     << subdetector << ")";
1149       continue;
1150     }  // end detector type check
1151   }    // end loop through ECal Hits
1152 
1153   if (verbosity > 1) {
1154     eventout += "\n          Number of ECal Hits collected:............. ";
1155     eventout += j;
1156   }
1157 
1158   ////////////////////////////
1159   // Get Preshower information
1160   ////////////////////////////
1161   // extract PreShower container
1162   edm::Handle<edm::PCaloHitContainer> PreShContainer;
1163   iEvent.getByToken(ECalESSrc_Token_, PreShContainer);
1164   if (!PreShContainer.isValid()) {
1165     edm::LogWarning(MsgLoggerCat) << "Unable to find EcalHitsES in event!";
1166     return;
1167   }
1168 
1169   // cycle through container
1170   i = 0, j = 0;
1171   for (itHit = PreShContainer->begin(); itHit != PreShContainer->end(); ++itHit) {
1172     ++i;
1173 
1174     // create a DetId from the detUnitId
1175     DetId theDetUnitId(itHit->id());
1176     int detector = theDetUnitId.det();
1177     int subdetector = theDetUnitId.subdetId();
1178 
1179     // check that expected detector is returned
1180     if ((detector == dEcal) && (subdetector == sdEcalPS)) {
1181       // get the Cell geometry
1182       auto theDet = (theCalo.getSubdetectorGeometry(theDetUnitId))->getGeometry(theDetUnitId);
1183 
1184       if (!theDet) {
1185         edm::LogWarning(MsgLoggerCat) << "Unable to get CaloCellGeometry from PreShContainer for Hit " << i;
1186         continue;
1187       }
1188 
1189       ++j;
1190 
1191       // get the global position of the cell
1192       const GlobalPoint &globalposition = theDet->getPosition();
1193 
1194       // gather necessary information
1195       PreShE.push_back(itHit->energy());
1196       PreShToF.push_back(itHit->time());
1197       PreShPhi.push_back(globalposition.phi());
1198       PreShEta.push_back(globalposition.eta());
1199 
1200     } else {
1201       edm::LogWarning(MsgLoggerCat) << "PreSh PCaloHit " << i << " is expected to be (det,subdet) = (" << dEcal << ","
1202                                     << sdEcalPS << "); value returned is: (" << detector << "," << subdetector << ")";
1203       continue;
1204     }  // end detector type check
1205   }    // end loop through PreShower Hits
1206 
1207   if (verbosity > 1) {
1208     eventout += "\n          Number of PreSh Hits collected:............ ";
1209     eventout += j;
1210   }
1211 
1212   if (verbosity > 0)
1213     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
1214 
1215   return;
1216 }
1217 
1218 void GlobalHitsProducer::storeECal(PGlobalSimHit &product) {
1219   std::string MsgLoggerCat = "GlobalHitsProducer_storeECal";
1220 
1221   if (verbosity > 2) {
1222     TString eventout("\n       nECalHits          = ");
1223     eventout += ECalE.size();
1224     for (unsigned int i = 0; i < ECalE.size(); ++i) {
1225       eventout += "\n          (e,tof,phi,eta) = (";
1226       eventout += ECalE[i];
1227       eventout += ", ";
1228       eventout += ECalToF[i];
1229       eventout += ", ";
1230       eventout += ECalPhi[i];
1231       eventout += ", ";
1232       eventout += ECalEta[i];
1233       eventout += ")";
1234     }  // end ECal output
1235     eventout += "\n       nPreShHits         = ";
1236     eventout += PreShE.size();
1237     for (unsigned int i = 0; i < PreShE.size(); ++i) {
1238       eventout += "\n          (e,tof,phi,eta) = (";
1239       eventout += PreShE[i];
1240       eventout += ", ";
1241       eventout += PreShToF[i];
1242       eventout += ", ";
1243       eventout += PreShPhi[i];
1244       eventout += ", ";
1245       eventout += PreShEta[i];
1246       eventout += ")";
1247     }  // end PreShower output
1248     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
1249   }  // end verbose output
1250 
1251   product.putECalHits(ECalE, ECalToF, ECalPhi, ECalEta);
1252   product.putPreShHits(PreShE, PreShToF, PreShPhi, PreShEta);
1253 
1254   return;
1255 }
1256 
1257 void GlobalHitsProducer::fillHCal(edm::Event &iEvent, const edm::EventSetup &iSetup) {
1258   std::string MsgLoggerCat = "GlobalHitsProducer_fillHCal";
1259 
1260   TString eventout;
1261   if (verbosity > 0)
1262     eventout = "\nGathering info:";
1263 
1264   // access the calorimeter geometry
1265   const auto &theCaloGeometry = iSetup.getHandle(caloGeomToken_);
1266   if (!theCaloGeometry.isValid()) {
1267     edm::LogWarning(MsgLoggerCat) << "Unable to find CaloGeometryRecord in event!";
1268     return;
1269   }
1270   const CaloGeometry &theCalo(*theCaloGeometry);
1271 
1272   // iterator to access containers
1273   edm::PCaloHitContainer::const_iterator itHit;
1274 
1275   ///////////////////////////////
1276   // get  HCal information
1277   ///////////////////////////////
1278   // extract HCal container
1279   edm::Handle<edm::PCaloHitContainer> HCalContainer;
1280   iEvent.getByToken(HCalSrc_Token_, HCalContainer);
1281   if (!HCalContainer.isValid()) {
1282     edm::LogWarning(MsgLoggerCat) << "Unable to find HCalHits in event!";
1283     return;
1284   }
1285 
1286   // cycle through container
1287   int i = 0, j = 0;
1288   for (itHit = HCalContainer->begin(); itHit != HCalContainer->end(); ++itHit) {
1289     ++i;
1290 
1291     // create a DetId from the detUnitId
1292     DetId theDetUnitId(itHit->id());
1293     int detector = theDetUnitId.det();
1294     int subdetector = theDetUnitId.subdetId();
1295 
1296     // check that expected detector is returned
1297     if ((detector == dHcal) && ((subdetector == sdHcalBrl) || (subdetector == sdHcalEC) || (subdetector == sdHcalOut) ||
1298                                 (subdetector == sdHcalFwd))) {
1299       // get the Cell geometry
1300       const HcalGeometry *theDet = dynamic_cast<const HcalGeometry *>(theCalo.getSubdetectorGeometry(theDetUnitId));
1301 
1302       if (!theDet) {
1303         edm::LogWarning(MsgLoggerCat) << "Unable to get CaloCellGeometry from HCalContainer for Hit " << i;
1304         continue;
1305       }
1306 
1307       ++j;
1308 
1309       // get the global position of the cell
1310       const GlobalPoint &globalposition = theDet->getPosition(theDetUnitId);
1311 
1312       // gather necessary information
1313       HCalE.push_back(itHit->energy());
1314       HCalToF.push_back(itHit->time());
1315       HCalPhi.push_back(globalposition.phi());
1316       HCalEta.push_back(globalposition.eta());
1317 
1318     } else {
1319       edm::LogWarning(MsgLoggerCat) << "HCal PCaloHit " << i << " is expected to be (det,subdet) = (" << dHcal << ","
1320                                     << sdHcalBrl << " || " << sdHcalEC << " || " << sdHcalOut << " || " << sdHcalFwd
1321                                     << "); value returned is: (" << detector << "," << subdetector << ")";
1322       continue;
1323     }  // end detector type check
1324   }    // end loop through HCal Hits
1325 
1326   if (verbosity > 1) {
1327     eventout += "\n          Number of HCal Hits collected:............. ";
1328     eventout += j;
1329   }
1330 
1331   if (verbosity > 0)
1332     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
1333 
1334   return;
1335 }
1336 
1337 void GlobalHitsProducer::storeHCal(PGlobalSimHit &product) {
1338   std::string MsgLoggerCat = "GlobalHitsProducer_storeHCal";
1339 
1340   if (verbosity > 2) {
1341     TString eventout("\n       nHCalHits          = ");
1342     eventout += HCalE.size();
1343     for (unsigned int i = 0; i < HCalE.size(); ++i) {
1344       eventout += "\n          (e,tof,phi,eta) = (";
1345       eventout += HCalE[i];
1346       eventout += ", ";
1347       eventout += HCalToF[i];
1348       eventout += ", ";
1349       eventout += HCalPhi[i];
1350       eventout += ", ";
1351       eventout += HCalEta[i];
1352       eventout += ")";
1353     }  // end HCal output
1354     edm::LogInfo(MsgLoggerCat) << eventout << "\n";
1355   }  // end verbose output
1356 
1357   product.putHCalHits(HCalE, HCalToF, HCalPhi, HCalEta);
1358 
1359   return;
1360 }
1361 
1362 void GlobalHitsProducer::clear() {
1363   std::string MsgLoggerCat = "GlobalHitsProducer_clear";
1364 
1365   if (verbosity > 0)
1366     edm::LogInfo(MsgLoggerCat) << "Clearing event holders";
1367 
1368   // reset G4MC info
1369   nRawGenPart = 0;
1370   G4VtxX.clear();
1371   G4VtxY.clear();
1372   G4VtxZ.clear();
1373   G4TrkPt.clear();
1374   G4TrkE.clear();
1375 
1376   // reset electromagnetic info
1377   // reset ECal info
1378   ECalE.clear();
1379   ECalToF.clear();
1380   ECalPhi.clear();
1381   ECalEta.clear();
1382   // reset Preshower info
1383   PreShE.clear();
1384   PreShToF.clear();
1385   PreShPhi.clear();
1386   PreShEta.clear();
1387 
1388   // reset hadronic info
1389   // reset HCal info
1390   HCalE.clear();
1391   HCalToF.clear();
1392   HCalPhi.clear();
1393   HCalEta.clear();
1394 
1395   // reset tracker info
1396   // reset Pixel info
1397   PxlBrlToF.clear();
1398   PxlBrlR.clear();
1399   PxlBrlPhi.clear();
1400   PxlBrlEta.clear();
1401   PxlFwdToF.clear();
1402   PxlFwdZ.clear();
1403   PxlFwdPhi.clear();
1404   PxlFwdEta.clear();
1405   // reset strip info
1406   SiBrlToF.clear();
1407   SiBrlR.clear();
1408   SiBrlPhi.clear();
1409   SiBrlEta.clear();
1410   SiFwdToF.clear();
1411   SiFwdZ.clear();
1412   SiFwdPhi.clear();
1413   SiFwdEta.clear();
1414 
1415   // reset muon info
1416   // reset muon DT info
1417   MuonDtToF.clear();
1418   MuonDtR.clear();
1419   MuonDtPhi.clear();
1420   MuonDtEta.clear();
1421   // reset muon CSC info
1422   MuonCscToF.clear();
1423   MuonCscZ.clear();
1424   MuonCscPhi.clear();
1425   MuonCscEta.clear();
1426   // rest muon RPC info
1427   MuonRpcBrlToF.clear();
1428   MuonRpcBrlR.clear();
1429   MuonRpcBrlPhi.clear();
1430   MuonRpcBrlEta.clear();
1431   MuonRpcFwdToF.clear();
1432   MuonRpcFwdZ.clear();
1433   MuonRpcFwdPhi.clear();
1434   MuonRpcFwdEta.clear();
1435 
1436   return;
1437 }