CMS 3D CMS Logo

List of all members | Public Member Functions | Protected Types | Protected Member Functions | Private Member Functions | Private Attributes
SensitiveDetector Class Referenceabstract

#include <SensitiveDetector.h>

Inheritance diagram for SensitiveDetector:
SensitiveCaloDetector SensitiveTkDetector CaloSD CaloTrkProcessing FiberSD HFChamberSD HFWedgeSD FP420SD MuonSensitiveDetector TimingSD TkAccumulatingSensitiveDetector TotemSD

Public Member Functions

virtual void clearHits ()=0
 
void EndOfEvent (G4HCofThisEvent *eventHC) override
 
const std::vector< std::string > & getNames () const
 
void Initialize (G4HCofThisEvent *eventHC) override
 
bool isCaloSD () const
 
G4bool ProcessHits (G4Step *step, G4TouchableHistory *tHistory) override=0
 
 SensitiveDetector (const std::string &iname, const DDCompactView &cpv, const SensitiveDetectorCatalog &, edm::ParameterSet const &p, bool calo)
 
virtual uint32_t setDetUnitId (const G4Step *step)=0
 
 ~SensitiveDetector () override
 

Protected Types

enum  coordinates { WorldCoordinates, LocalCoordinates }
 

Protected Member Functions

TrackInformationcmsTrackInformation (const G4Track *aTrack)
 
Local3DPoint ConvertToLocal3DPoint (const G4ThreeVector &point) const
 
Local3DPoint FinalStepPosition (const G4Step *step, coordinates) const
 
Local3DPoint InitialStepPosition (const G4Step *step, coordinates) const
 
Local3DPoint LocalPostStepPosition (const G4Step *step) const
 
Local3DPoint LocalPreStepPosition (const G4Step *step) const
 
void NaNTrap (const G4Step *step) const
 
void setNames (const std::vector< std::string > &)
 

Private Member Functions

void AssignSD (const std::string &vname)
 

Private Attributes

bool m_isCalo
 
std::vector< std::string > m_namesOfSD
 

Detailed Description

Definition at line 22 of file SensitiveDetector.h.

Member Enumeration Documentation

Enumerator
WorldCoordinates 
LocalCoordinates 

Definition at line 47 of file SensitiveDetector.h.

Constructor & Destructor Documentation

SensitiveDetector::SensitiveDetector ( const std::string &  iname,
const DDCompactView cpv,
const SensitiveDetectorCatalog clg,
edm::ParameterSet const &  p,
bool  calo 
)
explicit

Definition at line 18 of file SensitiveDetector.cc.

References AssignSD(), ecaldqm::collectionName, SensitiveDetectorCatalog::logicalNames(), and m_namesOfSD.

22  :
23  G4VSensitiveDetector(iname), m_isCalo(calo)
24 {
25  // for CMS hits
26  m_namesOfSD.push_back(iname);
27 
28  // Geant4 hit collection
29  collectionName.insert(iname);
30 
31  // register sensitive detector
32  G4SDManager * SDman = G4SDManager::GetSDMpointer();
33  SDman->AddNewDetector(this);
34 
35  const std::vector<std::string>& lvNames = clg.logicalNames(iname);
36  std::stringstream ss;
37  for (auto & lvname : lvNames) {
38  this->AssignSD(lvname);
39  ss << " " << lvname;
40  }
41  edm::LogVerbatim("SensitiveDetector")
42  << " <" << iname <<"> : Assigns SD to LVs " << ss.str();
43 }
const std::vector< std::string > & logicalNames(const std::string &readoutName) const
void AssignSD(const std::string &vname)
std::vector< std::string > m_namesOfSD
std::string const collectionName[nCollections]
Definition: Collections.h:47
SensitiveDetector::~SensitiveDetector ( )
override

Definition at line 45 of file SensitiveDetector.cc.

45 {}

Member Function Documentation

void SensitiveDetector::AssignSD ( const std::string &  vname)
private

Definition at line 51 of file SensitiveDetector.cc.

Referenced by ConvertToLocal3DPoint(), and SensitiveDetector().

52 {
53  G4LogicalVolumeStore * theStore = G4LogicalVolumeStore::GetInstance();
54  for (auto & lv : *theStore)
55  {
56  if (vname == lv->GetName()) {
57  lv->SetSensitiveDetector(this);
58  }
59  }
60 }
virtual void SensitiveDetector::clearHits ( )
pure virtual
TrackInformation * SensitiveDetector::cmsTrackInformation ( const G4Track *  aTrack)
protected

Definition at line 100 of file SensitiveDetector.cc.

References info().

Referenced by ConvertToLocal3DPoint(), TkAccumulatingSensitiveDetector::createHit(), MuonSensitiveDetector::createHit(), CaloSD::createNewHit(), CaloSD::getResponseWt(), TimingSD::getStepInfo(), CaloSD::getTrackID(), CaloSD::setTrackID(), TkAccumulatingSensitiveDetector::update(), and CaloSD::update().

101 {
102  TrackInformation* info = (TrackInformation*)(aTrack->GetUserInformation());
103  if (!info) {
104  edm::LogWarning("SensitiveDetector")
105  << " no TrackInformation available for trackID= "
106  << aTrack->GetTrackID();
107  throw SimG4Exception("SimG4CoreSensitiveDetector: cannot handle hits for "
108  + GetName());
109  }
110  return info;
111 }
static const TGPicture * info(bool iBackgroundIsBlack)
Local3DPoint SensitiveDetector::ConvertToLocal3DPoint ( const G4ThreeVector &  point) const
inlineprotected

Definition at line 55 of file SensitiveDetector.h.

References AssignSD(), cmsTrackInformation(), NaNTrap(), setNames(), and AlCaHLTBitMon_QueryRunRegistry::string.

Referenced by TkAccumulatingSensitiveDetector::createHit(), MuonSensitiveDetector::createHit(), TimingSD::EndOfEvent(), FinalStepPosition(), MuonSensitiveDetector::FinalStepPositionVsParent(), InitialStepPosition(), MuonSensitiveDetector::InitialStepPositionVsParent(), LocalPostStepPosition(), LocalPreStepPosition(), and TkAccumulatingSensitiveDetector::updateHit().

56  { return Local3DPoint(point.x(),point.y(),point.z()); }
Point3DBase< float, LocalTag > Local3DPoint
Definition: LocalPoint.h:9
*vegas h *****************************************************used in the default bin number in original ***version of VEGAS is ***a higher bin number might help to derive a more precise ***grade subtle point
Definition: invegas.h:5
void SensitiveDetector::EndOfEvent ( G4HCofThisEvent *  eventHC)
override

Definition at line 49 of file SensitiveDetector.cc.

49 {}
Local3DPoint SensitiveDetector::FinalStepPosition ( const G4Step *  step,
coordinates  cd 
) const
protected

Definition at line 73 of file SensitiveDetector.cc.

References ConvertToLocal3DPoint(), and WorldCoordinates.

Referenced by MuonSensitiveDetector::createHit(), and MuonSensitiveDetector::updateHit().

74 {
75  const G4StepPoint * postStepPoint = step->GetPostStepPoint();
76  const G4ThreeVector& globalCoordinates = postStepPoint->GetPosition();
77  if (cd == WorldCoordinates) { return ConvertToLocal3DPoint(globalCoordinates); }
78  const G4StepPoint * preStepPoint = step->GetPreStepPoint();
79  const G4ThreeVector localCoordinates = preStepPoint->GetTouchable()->GetHistory()
80  ->GetTopTransform().TransformPoint(globalCoordinates);
81  return ConvertToLocal3DPoint(localCoordinates);
82 }
Local3DPoint ConvertToLocal3DPoint(const G4ThreeVector &point) const
step
Definition: StallMonitor.cc:94
const std::vector<std::string>& SensitiveDetector::getNames ( ) const
inline
void SensitiveDetector::Initialize ( G4HCofThisEvent *  eventHC)
override

Definition at line 47 of file SensitiveDetector.cc.

47 {}
Local3DPoint SensitiveDetector::InitialStepPosition ( const G4Step *  step,
coordinates  cd 
) const
protected

Definition at line 62 of file SensitiveDetector.cc.

References ConvertToLocal3DPoint(), and WorldCoordinates.

Referenced by MuonSensitiveDetector::createHit(), and MuonSensitiveDetector::ProcessHits().

63 {
64  const G4StepPoint * preStepPoint = step->GetPreStepPoint();
65  const G4ThreeVector& globalCoordinates = preStepPoint->GetPosition();
66  if (cd == WorldCoordinates) { return ConvertToLocal3DPoint(globalCoordinates); }
67  const G4TouchableHistory * theTouchable=static_cast<const G4TouchableHistory *>(preStepPoint->GetTouchable());
68  const G4ThreeVector localCoordinates = theTouchable->GetHistory()
69  ->GetTopTransform().TransformPoint(globalCoordinates);
70  return ConvertToLocal3DPoint(localCoordinates);
71 }
Local3DPoint ConvertToLocal3DPoint(const G4ThreeVector &point) const
step
Definition: StallMonitor.cc:94
bool SensitiveDetector::isCaloSD ( ) const
inline

Definition at line 42 of file SensitiveDetector.h.

References m_isCalo.

42 { return m_isCalo; }
Local3DPoint SensitiveDetector::LocalPostStepPosition ( const G4Step *  step) const
protected

Definition at line 92 of file SensitiveDetector.cc.

References ConvertToLocal3DPoint().

Referenced by TkAccumulatingSensitiveDetector::createHit(), and TkAccumulatingSensitiveDetector::updateHit().

93 {
94  const G4ThreeVector& globalCoordinates = step->GetPostStepPoint()->GetPosition();
95  G4ThreeVector localCoordinates = step->GetPreStepPoint()->GetTouchable()->GetHistory()
96  ->GetTopTransform().TransformPoint(globalCoordinates);
97  return ConvertToLocal3DPoint(localCoordinates);
98 }
Local3DPoint ConvertToLocal3DPoint(const G4ThreeVector &point) const
step
Definition: StallMonitor.cc:94
Local3DPoint SensitiveDetector::LocalPreStepPosition ( const G4Step *  step) const
protected

Definition at line 84 of file SensitiveDetector.cc.

References ConvertToLocal3DPoint().

Referenced by TkAccumulatingSensitiveDetector::closeHit(), and TkAccumulatingSensitiveDetector::createHit().

85 {
86  const G4StepPoint * preStepPoint = step->GetPreStepPoint();
87  G4ThreeVector localCoordinates = preStepPoint->GetTouchable()->GetHistory()
88  ->GetTopTransform().TransformPoint(preStepPoint->GetPosition());
89  return ConvertToLocal3DPoint(localCoordinates);
90 }
Local3DPoint ConvertToLocal3DPoint(const G4ThreeVector &point) const
step
Definition: StallMonitor.cc:94
void SensitiveDetector::NaNTrap ( const G4Step *  step) const
protected

Definition at line 119 of file SensitiveDetector.cc.

References edm::isNotFinite().

Referenced by ConvertToLocal3DPoint(), and CaloSD::ProcessHits().

120 {
121  const G4Track* CurrentTrk = aStep->GetTrack();
122  const G4ThreeVector& CurrentPos = CurrentTrk->GetPosition();
123  double xyz = CurrentPos.x() + CurrentPos.y() + CurrentPos.z();
124  const G4ThreeVector& CurrentMom = CurrentTrk->GetMomentum();
125  xyz += CurrentMom.x() + CurrentMom.y() + CurrentMom.z();
126 
127  if( edm::isNotFinite(xyz) ) {
128 
129  const G4VPhysicalVolume* pCurrentVol = aStep->GetPreStepPoint()->GetPhysicalVolume();
130  const G4String& NameOfVol = pCurrentVol->GetName();
131  throw SimG4Exception("SimG4CoreSensitiveDetector: Corrupted Event - NaN detected in volume "
132  + NameOfVol);
133  }
134 }
constexpr bool isNotFinite(T x)
Definition: isFinite.h:9
G4bool SensitiveDetector::ProcessHits ( G4Step *  step,
G4TouchableHistory *  tHistory 
)
overridepure virtual
virtual uint32_t SensitiveDetector::setDetUnitId ( const G4Step *  step)
pure virtual
void SensitiveDetector::setNames ( const std::vector< std::string > &  hnames)
protected

Definition at line 113 of file SensitiveDetector.cc.

References m_namesOfSD.

Referenced by ConvertToLocal3DPoint(), and TkAccumulatingSensitiveDetector::TkAccumulatingSensitiveDetector().

114 {
115  m_namesOfSD.clear();
116  m_namesOfSD = hnames;
117 }
std::vector< std::string > m_namesOfSD

Member Data Documentation

bool SensitiveDetector::m_isCalo
private

Definition at line 68 of file SensitiveDetector.h.

Referenced by isCaloSD().

std::vector<std::string> SensitiveDetector::m_namesOfSD
private

Definition at line 67 of file SensitiveDetector.h.

Referenced by getNames(), SensitiveDetector(), and setNames().