Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-20 22:46:48

0001 #include "SimG4CMS/ShowerLibraryProducer/interface/HFWedgeSD.h"
0002 #include "DataFormats/Math/interface/Point3D.h"
0003 
0004 #include "G4VPhysicalVolume.hh"
0005 #include "G4PVPlacement.hh"
0006 #include "G4HCofThisEvent.hh"
0007 #include "G4TouchableHistory.hh"
0008 #include "G4Track.hh"
0009 #include "G4Step.hh"
0010 #include "G4VSolid.hh"
0011 #include "G4DynamicParticle.hh"
0012 #include "G4ParticleDefinition.hh"
0013 #include "G4SDManager.hh"
0014 #include "G4ios.hh"
0015 
0016 #include "G4PhysicalConstants.hh"
0017 #include "G4SystemOfUnits.hh"
0018 
0019 //#define EDM_ML_DEBUG
0020 
0021 HFWedgeSD::HFWedgeSD(const std::string& iname,
0022                            const SensitiveDetectorCatalog& clg,
0023                            const SimTrackManager* manager)
0024     : SensitiveCaloDetector(iname, clg), hcID(-1), theHC(nullptr), currentHit(nullptr) {
0025   edm::LogVerbatim("FiberSim") << "HFWedgeSD : Instantiated for " << iname;
0026 }
0027 
0028 HFWedgeSD::~HFWedgeSD() { delete theHC; }
0029 
0030 void HFWedgeSD::Initialize(G4HCofThisEvent* HCE) {
0031   edm::LogVerbatim("FiberSim") << "HFWedgeSD : Initialize called for " << GetName() << " in collection " << HCE;
0032   theHC = new HFShowerG4HitsCollection(GetName(), collectionName[0]);
0033   if (hcID < 0)
0034     hcID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
0035   HCE->AddHitsCollection(hcID, theHC);
0036   edm::LogVerbatim("FiberSim") << "HFWedgeSD : Add hit collectrion for " << collectionName[0] << ":" << hcID << ":" << theHC;
0037 
0038   clearHits();
0039 }
0040 
0041 G4bool HFWedgeSD::ProcessHits(G4Step* aStep, G4TouchableHistory*) {
0042   G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
0043   const G4VTouchable* touch = preStepPoint->GetTouchable();
0044   currentID = setDetUnitId(aStep);
0045   trackID = aStep->GetTrack()->GetTrackID();
0046   edep = aStep->GetTotalEnergyDeposit();
0047   time = (preStepPoint->GetGlobalTime()) / ns;
0048 
0049   globalPos = preStepPoint->GetPosition();
0050   localPos = touch->GetHistory()->GetTopTransform().TransformPoint(globalPos);
0051   const G4DynamicParticle* particle = aStep->GetTrack()->GetDynamicParticle();
0052   momDir = particle->GetMomentumDirection();
0053 
0054   if (hitExists() == false && edep > 0.)
0055     currentHit = createNewHit();
0056 
0057   return true;
0058 }
0059 
0060 void HFWedgeSD::EndOfEvent(G4HCofThisEvent* HCE) {
0061   edm::LogVerbatim("FiberSim") << "HFWedgeSD: Sees" << theHC->entries() << " hits";
0062   clear();
0063 }
0064 
0065 void HFWedgeSD::clear() {}
0066 
0067 void HFWedgeSD::DrawAll() {}
0068 
0069 void HFWedgeSD::PrintAll() {}
0070 
0071 G4bool HFWedgeSD::hitExists() {
0072   // Update if in the same detector, time-slice and for same track
0073   if (currentID == previousID) {
0074     updateHit(currentHit);
0075     return true;
0076   }
0077 
0078   std::map<int, HFShowerG4Hit*>::const_iterator it = hitMap.find(currentID);
0079   if (it != hitMap.end()) {
0080     updateHit(currentHit);
0081     return true;
0082   }
0083 
0084   return false;
0085 }
0086 
0087 HFShowerG4Hit* HFWedgeSD::createNewHit() {
0088 #ifdef EDM_ML_DEBUG
0089   edm::LogVerbatim("FiberSim") << "HFWedgeSD::CreateNewHit for ID " << currentID << " Track " << trackID
0090                                << " Edep: " << edep / CLHEP::MeV << " MeV; Time: " << time << " ns; Position (local) "
0091                                << localPos << " (global ) " << globalPos << " direction " << momDir;
0092 #endif
0093   HFShowerG4Hit* aHit = new HFShowerG4Hit;
0094   aHit->setHitId(currentID);
0095   aHit->setTrackId(trackID);
0096   aHit->setTime(time);
0097   aHit->setLocalPos(localPos);
0098   aHit->setGlobalPos(globalPos);
0099   aHit->setPrimMomDir(momDir);
0100   updateHit(aHit);
0101 
0102   theHC->insert(aHit);
0103   hitMap.insert(std::pair<int, HFShowerG4Hit*>(previousID, aHit));
0104 
0105   return aHit;
0106 }
0107 
0108 void HFWedgeSD::updateHit(HFShowerG4Hit* aHit) {
0109   if (edep != 0) {
0110     aHit->updateEnergy(edep);
0111 #ifdef EDM_ML_DEBUG
0112     edm::LogVerbatim("FiberSim") << "HFWedgeSD: Add energy deposit in " << currentID << " edep " << edep / CLHEP::MeV << " MeV";
0113 #endif
0114   }
0115   previousID = currentID;
0116 }
0117 
0118 void HFWedgeSD::clearHits() {
0119   hitMap.erase(hitMap.begin(), hitMap.end());
0120   previousID = -1;
0121 }
0122 
0123 uint32_t HFWedgeSD::setDetUnitId(const G4Step* aStep) {
0124   const G4VTouchable* touch = aStep->GetPreStepPoint()->GetTouchable();
0125   return (touch->GetReplicaNumber(0));
0126 }
0127 
0128 void HFWedgeSD::fillHits(edm::PCaloHitContainer&, const std::string&) {}