File indexing completed on 2023-03-17 11:25:05
0001 #include "SimG4Core/SensitiveDetector/interface/SensitiveDetector.h"
0002
0003 #include "FWCore/Utilities/interface/isFinite.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005
0006 #include "G4SDManager.hh"
0007 #include "G4Step.hh"
0008 #include "G4Track.hh"
0009 #include "G4StepPoint.hh"
0010 #include "G4Transform3D.hh"
0011 #include "G4LogicalVolumeStore.hh"
0012 #include "G4TouchableHistory.hh"
0013 #include "G4VUserTrackInformation.hh"
0014
0015 #include <sstream>
0016
0017 SensitiveDetector::SensitiveDetector(const std::string& iname, const SensitiveDetectorCatalog& clg, bool calo)
0018 : G4VSensitiveDetector(iname), m_isCalo(calo) {
0019
0020 m_namesOfSD.push_back(iname);
0021
0022
0023 collectionName.insert(iname);
0024
0025
0026 G4SDManager* SDman = G4SDManager::GetSDMpointer();
0027 SDman->AddNewDetector(this);
0028
0029 const std::vector<std::string_view>& lvNames = clg.logicalNames(iname);
0030 std::stringstream ss;
0031 for (auto& lvname : lvNames) {
0032 this->AssignSD({lvname.data(), lvname.size()});
0033 ss << " " << lvname;
0034 }
0035 edm::LogVerbatim("SensitiveDetector") << " <" << iname << "> : Assigns SD to LVs " << ss.str();
0036 }
0037
0038 SensitiveDetector::~SensitiveDetector() {}
0039
0040 void SensitiveDetector::Initialize(G4HCofThisEvent* eventHC) {}
0041
0042 void SensitiveDetector::EndOfEvent(G4HCofThisEvent* eventHC) {}
0043
0044 void SensitiveDetector::AssignSD(const std::string& vname) {
0045 G4LogicalVolumeStore* theStore = G4LogicalVolumeStore::GetInstance();
0046 for (auto& lv : *theStore) {
0047 if (vname == lv->GetName()) {
0048 lv->SetSensitiveDetector(this);
0049 }
0050 }
0051 }
0052
0053 Local3DPoint SensitiveDetector::InitialStepPosition(const G4Step* step, coordinates cd) const {
0054 const G4StepPoint* preStepPoint = step->GetPreStepPoint();
0055 const G4ThreeVector& globalCoordinates = preStepPoint->GetPosition();
0056 if (cd == WorldCoordinates) {
0057 return ConvertToLocal3DPoint(globalCoordinates);
0058 }
0059 const G4TouchableHistory* theTouchable = static_cast<const G4TouchableHistory*>(preStepPoint->GetTouchable());
0060 const G4ThreeVector localCoordinates =
0061 theTouchable->GetHistory()->GetTopTransform().TransformPoint(globalCoordinates);
0062 return ConvertToLocal3DPoint(localCoordinates);
0063 }
0064
0065 Local3DPoint SensitiveDetector::FinalStepPosition(const G4Step* step, coordinates cd) const {
0066 const G4StepPoint* postStepPoint = step->GetPostStepPoint();
0067 const G4ThreeVector& globalCoordinates = postStepPoint->GetPosition();
0068 if (cd == WorldCoordinates) {
0069 return ConvertToLocal3DPoint(globalCoordinates);
0070 }
0071 const G4StepPoint* preStepPoint = step->GetPreStepPoint();
0072 const G4ThreeVector localCoordinates =
0073 preStepPoint->GetTouchable()->GetHistory()->GetTopTransform().TransformPoint(globalCoordinates);
0074 return ConvertToLocal3DPoint(localCoordinates);
0075 }
0076
0077 Local3DPoint SensitiveDetector::LocalPreStepPosition(const G4Step* step) const {
0078 const G4StepPoint* preStepPoint = step->GetPreStepPoint();
0079 G4ThreeVector localCoordinates =
0080 preStepPoint->GetTouchable()->GetHistory()->GetTopTransform().TransformPoint(preStepPoint->GetPosition());
0081 return ConvertToLocal3DPoint(localCoordinates);
0082 }
0083
0084 Local3DPoint SensitiveDetector::LocalPostStepPosition(const G4Step* step) const {
0085 const G4ThreeVector& globalCoordinates = step->GetPostStepPoint()->GetPosition();
0086 G4ThreeVector localCoordinates =
0087 step->GetPreStepPoint()->GetTouchable()->GetHistory()->GetTopTransform().TransformPoint(globalCoordinates);
0088 return ConvertToLocal3DPoint(localCoordinates);
0089 }
0090
0091 TrackInformation* SensitiveDetector::cmsTrackInformation(const G4Track* aTrack) {
0092 TrackInformation* info = (TrackInformation*)(aTrack->GetUserInformation());
0093 if (nullptr == info) {
0094 edm::LogWarning("SensitiveDetector") << " no TrackInformation available for trackID= " << aTrack->GetTrackID()
0095 << " inside SD " << GetName();
0096 G4Exception(
0097 "SensitiveDetector::cmsTrackInformation()", "sd01", FatalException, "cannot handle hits without trackinfo");
0098 }
0099 return info;
0100 }
0101
0102 void SensitiveDetector::setNames(const std::vector<std::string>& hnames) {
0103 m_namesOfSD.clear();
0104 m_namesOfSD = hnames;
0105 }
0106
0107 void SensitiveDetector::NaNTrap(const G4Step* aStep) const {
0108 G4Track* currentTrk = aStep->GetTrack();
0109 double ekin = currentTrk->GetKineticEnergy();
0110 if (ekin < 0.0) {
0111 const G4VPhysicalVolume* pCurrentVol = aStep->GetPreStepPoint()->GetPhysicalVolume();
0112 edm::LogWarning("SensitiveDetector") << "Negative kinetic energy Ekin(MeV)=" << ekin / CLHEP::MeV << " of "
0113 << currentTrk->GetDefinition()->GetParticleName()
0114 << " trackID= " << currentTrk->GetTrackID() << " inside "
0115 << pCurrentVol->GetName();
0116 currentTrk->SetKineticEnergy(0.0);
0117 }
0118 const G4ThreeVector& currentPos = currentTrk->GetPosition();
0119 double xyz = currentPos.x() + currentPos.y() + currentPos.z();
0120 const G4ThreeVector& currentMom = currentTrk->GetMomentum();
0121 xyz += currentMom.x() + currentMom.y() + currentMom.z();
0122
0123 if (edm::isNotFinite(xyz)) {
0124 const G4VPhysicalVolume* pCurrentVol = aStep->GetPreStepPoint()->GetPhysicalVolume();
0125 edm::LogWarning("SensitiveDetector") << "NaN detected for trackID= " << currentTrk->GetTrackID() << " inside "
0126 << pCurrentVol->GetName();
0127 G4Exception("SensitiveDetector::NaNTrap()", "sd01", FatalException, "corrupted event or step");
0128 }
0129 }