Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-13 03:16:20

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