CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
HGCSD.cc
Go to the documentation of this file.
1 // File: HGCSD.cc
3 // Description: Sensitive Detector class for Combined Forward Calorimeter
5 
7 
16 
17 #include "G4LogicalVolumeStore.hh"
18 #include "G4LogicalVolume.hh"
19 #include "G4Step.hh"
20 #include "G4Track.hh"
21 #include "G4ParticleTable.hh"
22 #include "G4VProcess.hh"
23 #include "G4Trap.hh"
24 
25 #include <iostream>
26 #include <fstream>
27 #include <iomanip>
28 
29 //#define DebugLog
30 
31 HGCSD::HGCSD(G4String name, const DDCompactView & cpv,
32  const SensitiveDetectorCatalog & clg,
33  edm::ParameterSet const & p, const SimTrackManager* manager) :
34  CaloSD(name, cpv, clg, p, manager,
35  (float)(p.getParameter<edm::ParameterSet>("HGCSD").getParameter<double>("TimeSliceUnit")),
36  p.getParameter<edm::ParameterSet>("HGCSD").getParameter<bool>("IgnoreTrackID")),
37  numberingScheme(0) {
38 
40  eminHit = m_HGC.getParameter<double>("EminHit")*CLHEP::MeV;
41  bool checkID = m_HGC.getUntrackedParameter<bool>("CheckID", false);
42  verbosity = m_HGC.getUntrackedParameter<int>("Verbosity",0);
43 
44  //this is defined in the hgcsens.xml
45  G4String myName(this->nameOfSD());
47  std::string nameX("HGCal");
48  if (myName.find("HitsEE")!=std::string::npos) {
50  nameX = "HGCalEESensitive";
51  } else if (myName.find("HitsHEfront")!=std::string::npos) {
53  nameX = "HGCalHESiliconSensitive";
54  } else if (myName.find("HitsHEback")!=std::string::npos) {
56  nameX = "HGCalHEScintillatorSensitive";
57  }
58 
59 #ifdef DebugLog
60  LogDebug("HGCSim") << "**************************************************"
61  << "\n"
62  << "* *"
63  << "\n"
64  << "* Constructing a HGCSD with name " << name << "\n"
65  << "* *"
66  << "\n"
67  << "**************************************************";
68 #endif
69  edm::LogInfo("HGCSim") << "HGCSD:: Threshold for storing hits: " << eminHit;
70 
71  numberingScheme = new HGCNumberingScheme(cpv,nameX,checkID,verbosity);
72 }
73 
75  if (numberingScheme) delete numberingScheme;
76 }
77 
78 bool HGCSD::ProcessHits(G4Step * aStep, G4TouchableHistory * ) {
79 
80  NaNTrap( aStep ) ;
81 
82  if (aStep == NULL) {
83  return true;
84  } else {
85 #ifdef DebugLog
86  G4int parCode = aStep->GetTrack()->GetDefinition()->GetPDGEncoding();
87  bool notaMuon = (parCode == mupPDG || parCode == mumPDG ) ? false : true;
88  G4LogicalVolume* lv =
89  aStep->GetPreStepPoint()->GetPhysicalVolume()->GetLogicalVolume();
90  edm::LogInfo("HGCSim") << "HGCSD: Hit from standard path from "
91  << lv->GetName() << " for Track "
92  << aStep->GetTrack()->GetTrackID() << " ("
93  << aStep->GetTrack()->GetDefinition()->GetParticleName()
94  << ":" << notaMuon << ")";
95 #endif
96  if (getStepInfo(aStep)) {
97  if (hitExists() == false && edepositEM+edepositHAD>0.) currentHit = createNewHit();
98  }
99  return true;
100  }
101 }
102 
103 double HGCSD::getEnergyDeposit(G4Step* aStep) {
104  double destep = aStep->GetTotalEnergyDeposit();
105  return destep;
106 }
107 
108 uint32_t HGCSD::setDetUnitId(G4Step * aStep) {
109 
110  G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
111  const G4VTouchable* touch = preStepPoint->GetTouchable();
112 
113  //determine the exact position in global coordinates in the mass geometry
114  G4ThreeVector hitPoint = preStepPoint->GetPosition();
115  float globalZ=touch->GetTranslation(0).z();
116  int iz( globalZ>0 ? 1 : -1);
117 
118  //convert to local coordinates (=local to the current volume):
119  G4ThreeVector localpos = touch->GetHistory()->GetTopTransform().TransformPoint(hitPoint);
120 
121  //get the det unit id with
123 
124  int layer = touch->GetReplicaNumber(0);
125  int module = touch->GetReplicaNumber(1);
126  if (verbosity > 0)
127  std::cout << "HGCSD::Global " << hitPoint << " local " << localpos
128  << std::endl;
129  return setDetUnitId (subdet, layer, module, iz, localpos);
130 }
131 
133  G4ParticleTable * theParticleTable = G4ParticleTable::GetParticleTable();
134  G4String particleName;
135  mumPDG = theParticleTable->FindParticle(particleName="mu-")->GetPDGEncoding();
136  mupPDG = theParticleTable->FindParticle(particleName="mu+")->GetPDGEncoding();
137 #ifdef DebugLog
138  LogDebug("HGCSim") << "HGCSD: Particle code for mu- = " << mumPDG
139  << " for mu+ = " << mupPDG;
140 #endif
141 }
142 
143 bool HGCSD::filterHit(CaloG4Hit* aHit, double time) {
144  return ((time <= tmaxHit) && (aHit->getEnergyDeposit() > eminHit));
145 }
146 
147 
148 //
149 uint32_t HGCSD::setDetUnitId (ForwardSubdetector &subdet, int &layer, int &module, int &iz, G4ThreeVector &pos) {
150  return (numberingScheme ? numberingScheme->getUnitID(subdet, layer, module, iz, pos) : 0);
151 }
152 
153 //
154 int HGCSD::setTrackID (G4Step* aStep) {
155  theTrack = aStep->GetTrack();
156 
157  double etrack = preStepPoint->GetKineticEnergy();
158  TrackInformation * trkInfo = (TrackInformation *)(theTrack->GetUserInformation());
159  int primaryID = trkInfo->getIDonCaloSurface();
160  if (primaryID == 0) {
161 #ifdef DebugLog
162  edm::LogInfo("HGCSim") << "HGCSD: Problem with primaryID **** set by "
163  << "force to TkID **** " <<theTrack->GetTrackID();
164 #endif
165  primaryID = theTrack->GetTrackID();
166  }
167 
168  if (primaryID != previousID.trackID())
169  resetForNewPrimary(preStepPoint->GetPosition(), etrack);
170 
171  return primaryID;
172 }
#define LogDebug(id)
float edepositEM
Definition: CaloSD.h:121
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
virtual ~HGCSD()
Definition: HGCSD.cc:74
HGCSD(G4String, const DDCompactView &, const SensitiveDetectorCatalog &, edm::ParameterSet const &, const SimTrackManager *)
Definition: HGCSD.cc:31
std::string nameOfSD()
HGCNumberingScheme * numberingScheme
Definition: HGCSD.h:46
int getIDonCaloSurface() const
Definition: CaloSD.h:42
virtual bool ProcessHits(G4Step *, G4TouchableHistory *)
Definition: HGCSD.cc:78
#define NULL
Definition: scimark2.h:8
virtual uint32_t getUnitID(ForwardSubdetector subdet, int layer, int module, int iz, const G4ThreeVector &pos)
assigns the det id to a hit
ForwardSubdetector
type of data representation of DDCompactView
Definition: DDCompactView.h:77
static TrackerG4SimHitNumberingScheme & numberingScheme(const DDCompactView &cpv, const GeometricDet &det)
G4int mupPDG
Definition: HGCSD.h:47
double eminHit
Definition: HGCSD.h:48
const double MeV
void resetForNewPrimary(const G4ThreeVector &, double)
Definition: CaloSD.cc:454
virtual void initRun()
Definition: HGCSD.cc:132
float edepositHAD
Definition: CaloSD.h:121
int trackID() const
Definition: CaloHitID.h:25
CaloHitID previousID
Definition: CaloSD.h:117
CaloG4Hit * currentHit
Definition: CaloSD.h:128
void NaNTrap(G4Step *step)
G4Track * theTrack
Definition: CaloSD.h:118
virtual G4bool getStepInfo(G4Step *aStep)
Definition: CaloSD.cc:250
double tmaxHit
Definition: CaloSD.h:123
G4int verbosity
Definition: HGCSD.h:47
G4StepPoint * preStepPoint
Definition: CaloSD.h:120
virtual double getEnergyDeposit(G4Step *)
Definition: HGCSD.cc:103
virtual bool filterHit(CaloG4Hit *, double)
Definition: HGCSD.cc:143
int setTrackID(G4Step *step)
Definition: HGCSD.cc:154
G4bool hitExists()
Definition: CaloSD.cc:310
tuple cout
Definition: gather_cfg.py:121
ForwardSubdetector myFwdSubdet_
Definition: HGCSD.h:49
Definition: vlib.h:208
virtual uint32_t setDetUnitId(G4Step *step)
Definition: HGCSD.cc:108
CaloG4Hit * createNewHit()
Definition: CaloSD.cc:363
G4int mumPDG
Definition: HGCSD.h:47
double getEnergyDeposit() const
Definition: CaloG4Hit.h:81