Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-05-10 02:21:34 

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

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!