Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-04-28 04:48:45

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, const SensitiveDetectorCatalog& clg, const SimTrackManager* manager)
0022     : SensitiveCaloDetector(iname, clg), hcID(-1), theHC(nullptr), currentHit(nullptr) {
0023   edm::LogVerbatim("FiberSim") << "HFWedgeSD : Instantiated for " << iname;
0024 }
0025 
0026 HFWedgeSD::~HFWedgeSD() { delete theHC; }
0027 
0028 void HFWedgeSD::Initialize(G4HCofThisEvent* HCE) {
0029   edm::LogVerbatim("FiberSim") << "HFWedgeSD : Initialize called for " << GetName() << " in collection " << HCE;
0030   theHC = new HFShowerG4HitsCollection(GetName(), collectionName[0]);
0031   if (hcID < 0)
0032     hcID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
0033   HCE->AddHitsCollection(hcID, theHC);
0034   edm::LogVerbatim("FiberSim") << "HFWedgeSD : Add hit collectrion for " << collectionName[0] << ":" << hcID << ":"
0035                                << theHC;
0036 
0037   clearHits();
0038 }
0039 
0040 G4bool HFWedgeSD::ProcessHits(G4Step* aStep, G4TouchableHistory*) {
0041   G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
0042   const G4VTouchable* touch = preStepPoint->GetTouchable();
0043   currentID = setDetUnitId(aStep);
0044   trackID = aStep->GetTrack()->GetTrackID();
0045   edep = aStep->GetTotalEnergyDeposit();
0046   time = (preStepPoint->GetGlobalTime()) / ns;
0047 
0048   globalPos = preStepPoint->GetPosition();
0049   localPos = touch->GetHistory()->GetTopTransform().TransformPoint(globalPos);
0050   const G4DynamicParticle* particle = aStep->GetTrack()->GetDynamicParticle();
0051   momDir = particle->GetMomentumDirection();
0052 
0053   if (hitExists() == false && edep > 0.)
0054     currentHit = createNewHit();
0055 
0056   return true;
0057 }
0058 
0059 void HFWedgeSD::EndOfEvent(G4HCofThisEvent* HCE) {
0060   edm::LogVerbatim("FiberSim") << "HFWedgeSD: Sees" << theHC->entries() << " hits";
0061   clear();
0062 }
0063 
0064 void HFWedgeSD::clear() {}
0065 
0066 void HFWedgeSD::DrawAll() {}
0067 
0068 void HFWedgeSD::PrintAll() {}
0069 
0070 G4bool HFWedgeSD::hitExists() {
0071   // Update if in the same detector, time-slice and for same track
0072   if (currentID == previousID) {
0073     updateHit(currentHit);
0074     return true;
0075   }
0076 
0077   std::map<int, HFShowerG4Hit*>::const_iterator it = hitMap.find(currentID);
0078   if (it != hitMap.end()) {
0079     updateHit(currentHit);
0080     return true;
0081   }
0082 
0083   return false;
0084 }
0085 
0086 HFShowerG4Hit* HFWedgeSD::createNewHit() {
0087 #ifdef EDM_ML_DEBUG
0088   edm::LogVerbatim("FiberSim") << "HFWedgeSD::CreateNewHit for ID " << currentID << " Track " << trackID
0089                                << " Edep: " << edep / CLHEP::MeV << " MeV; Time: " << time << " ns; Position (local) "
0090                                << localPos << " (global ) " << globalPos << " direction " << momDir;
0091 #endif
0092   HFShowerG4Hit* aHit = new HFShowerG4Hit;
0093   aHit->setHitId(currentID);
0094   aHit->setTrackId(trackID);
0095   aHit->setTime(time);
0096   aHit->setLocalPos(localPos);
0097   aHit->setGlobalPos(globalPos);
0098   aHit->setPrimMomDir(momDir);
0099   updateHit(aHit);
0100 
0101   theHC->insert(aHit);
0102   hitMap.insert(std::pair<int, HFShowerG4Hit*>(previousID, aHit));
0103 
0104   return aHit;
0105 }
0106 
0107 void HFWedgeSD::updateHit(HFShowerG4Hit* aHit) {
0108   if (edep != 0) {
0109     aHit->updateEnergy(edep);
0110 #ifdef EDM_ML_DEBUG
0111     edm::LogVerbatim("FiberSim") << "HFWedgeSD: Add energy deposit in " << currentID << " edep " << edep / CLHEP::MeV
0112                                  << " 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&) {}