Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 10:04:18

0001 #include "SimG4CMS/ShowerLibraryProducer/interface/HFChamberSD.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 HFChamberSD::HFChamberSD(const std::string& name, const SensitiveDetectorCatalog& clg, const SimTrackManager* manager)
0022     : SensitiveCaloDetector(name, clg), theHCID(-1), theHC(nullptr), theNSteps(0) {
0023   edm::LogVerbatim("FiberSim") << "HFChamberSD : Instantiated for " << name;
0024 }
0025 
0026 HFChamberSD::~HFChamberSD() { delete theHC; }
0027 
0028 void HFChamberSD::Initialize(G4HCofThisEvent* HCE) {
0029   edm::LogVerbatim("FiberSim") << "HFChamberSD : Initialize called for " << GetName() << " in collection " << HCE;
0030   theHC = new HFShowerG4HitsCollection(GetName(), collectionName[0]);
0031   if (theHCID < 0)
0032     theHCID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
0033   HCE->AddHitsCollection(theHCID, theHC);
0034   edm::LogVerbatim("FiberSim") << "HFChamberSD : Add hit collectrion for " << collectionName[0] << ":" << theHCID << ":"
0035                                << theHC;
0036 }
0037 
0038 G4bool HFChamberSD::ProcessHits(G4Step* aStep, G4TouchableHistory*) {
0039   //do not process hits other than primary particle hits:
0040   double charge = aStep->GetTrack()->GetDefinition()->GetPDGCharge();
0041   int trackID = aStep->GetTrack()->GetTrackID();
0042   if (charge == 0. || trackID != 1 || aStep->GetTrack()->GetParentID() != 0 ||
0043       aStep->GetTrack()->GetCreatorProcess() != nullptr)
0044     return false;
0045   ++theNSteps;
0046   //if(theNSteps>1)return false;
0047 
0048   G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
0049   const G4VTouchable* touch = preStepPoint->GetTouchable();
0050   int detID = setDetUnitId(aStep);
0051 
0052   double edep = aStep->GetTotalEnergyDeposit();
0053   double time = (preStepPoint->GetGlobalTime()) / ns;
0054 
0055   const G4ThreeVector& globalPos = preStepPoint->GetPosition();
0056   G4ThreeVector localPos = touch->GetHistory()->GetTopTransform().TransformPoint(globalPos);
0057   const G4DynamicParticle* particle = aStep->GetTrack()->GetDynamicParticle();
0058   const G4ThreeVector& momDir = particle->GetMomentumDirection();
0059 
0060   HFShowerG4Hit* aHit = new HFShowerG4Hit(detID, trackID, edep, time);
0061   aHit->setLocalPos(localPos);
0062   aHit->setGlobalPos(globalPos);
0063   aHit->setPrimMomDir(momDir);
0064 
0065 #ifdef EDM_ML_DEBUG
0066   edm::LogVerbatim("FiberSim") << "HFChamberSD: Hit created in (" << touch->GetVolume(0)->GetLogicalVolume()->GetName()
0067                                << ")  ID " << detID << " Track " << trackID << " Edep: " << edep / CLHEP::MeV
0068                                << " MeV; Time: " << time << " ns; Position (local) " << localPos << " (global ) "
0069                                << globalPos << " direction " << momDir;
0070 #endif
0071   theHC->insert(aHit);
0072   return true;
0073 }
0074 
0075 void HFChamberSD::EndOfEvent(G4HCofThisEvent* HCE) {
0076   edm::LogVerbatim("FiberSim") << "HFChamberSD: Finds " << theHC->entries() << " hits";
0077   clear();
0078 }
0079 
0080 void HFChamberSD::clear() { theNSteps = 0; }
0081 
0082 void HFChamberSD::DrawAll() {}
0083 
0084 void HFChamberSD::PrintAll() {}
0085 
0086 void HFChamberSD::clearHits() {}
0087 
0088 uint32_t HFChamberSD::setDetUnitId(const G4Step* aStep) {
0089   const G4VTouchable* touch = aStep->GetPreStepPoint()->GetTouchable();
0090   return (touch->GetReplicaNumber(0));
0091 }
0092 
0093 void HFChamberSD::fillHits(edm::PCaloHitContainer&, const std::string&) {}