Belle II Software  release-05-01-25
SensitiveDetector.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Martin Ritter, Igal Jaegle *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <beast/claw/simulation/SensitiveDetector.h>
12 #include <beast/claw/dataobjects/ClawSimHit.h>
13 
14 #include <framework/datastore/StoreArray.h>
15 #include <framework/datastore/RelationArray.h>
16 
17 #include <G4Track.hh>
18 #include <G4Step.hh>
19 #include <G4Version.hh>
20 
21 namespace Belle2 {
27  namespace claw {
28 
30  Simulation::SensitiveDetectorBase("ClawSensitiveDetector", Const::invalidDetector)
31  {
32  m_hitNum = 0;
33  m_EvnetNumber = 0;
34  m_oldEvnetNumber = 0;
35  m_trackID = 0;
36  m_startTime = 0;
37  m_endTime = 0;
38  m_WightedTime = 0;
39  m_startEnergy = 0;
40  m_energyDeposit = 0;
41  m_trackLength = 0;
42  iECLCell = 0;
43  TimeIndex = 0;
44  local_pos = 0;
45  T_ave = 0;
46  firstcall = 0;
47  m_phiID = 0;
48  m_thetaID = 0;
49  m_cellID = 0;
50 
51  //Make sure all collections are registered
52  StoreArray<MCParticle> mcParticles;
53  StoreArray<ClawSimHit> simHits;
54  RelationArray relMCSimHit(mcParticles, simHits);
55 
56  //Register all collections we want to modify and require those we want to use
57  mcParticles.registerInDataStore();
58  simHits.registerInDataStore();
59  relMCSimHit.registerInDataStore();
60 
61  //Register the Relation so that the TrackIDs get replaced by the actual
62  //MCParticle indices after simulating the events. This is needed as
63  //secondary particles might not be stored so everything relating to those
64  //particles will be attributed to the last saved mother particle
65  registerMCParticleRelation(relMCSimHit);
66 #if G4VERSION_NUMBER < 1010
67  saturationEngine = new G4EmSaturation();
68 #else
69  saturationEngine = new G4EmSaturation(1); // verbose level
70 #endif
71 
72  }
73 
74  SensitiveDetector::~SensitiveDetector()
75  {
76 
77  }
78 
79  bool SensitiveDetector::step(G4Step* step, G4TouchableHistory*)
80  {
81  const G4StepPoint& preStep = *step->GetPreStepPoint();
82  const G4StepPoint& postStep = *step->GetPostStepPoint();
83 
84  G4Track& track = *step->GetTrack();
85  if (m_trackID != track.GetTrackID()) {
86  //TrackID changed, store track informations
87  m_trackID = track.GetTrackID();
88  //Get momentum
89  m_momentum = preStep.GetMomentum() ;
90  //Get energy
91  m_startEnergy = preStep.GetKineticEnergy() ;
92  //Reset energy deposit;
93  m_energyDeposit = 0;
94  //Reset Wighted Time;
95  m_WightedTime = 0;
96  //Reset m_WightedPos;
97  m_WightedPos.SetXYZ(0, 0, 0);
98 
99  }
100  //Update energy deposit
101  const double depEnergy = step->GetTotalEnergyDeposit() ;
102  m_energyDeposit += saturationEngine->VisibleEnergyDeposition(track.GetDefinition(), track.GetMaterialCutsCouple(),
103  step->GetStepLength(), depEnergy, 0.);
104  //step->GetTotalEnergyDeposit() ;
105 
106  m_startTime = preStep.GetGlobalTime();
107  m_endTime = postStep.GetGlobalTime();
108  m_WightedTime += (m_startTime + m_endTime) / 2 * (step->GetTotalEnergyDeposit());
109 
110  m_startPos = preStep.GetPosition();
111  m_endPos = postStep.GetPosition();
112  TVector3 position((m_startPos.getX() + m_endPos.getX()) / 2 / CLHEP::cm, (m_startPos.getY() + m_endPos.getY()) / 2 / CLHEP::cm,
113  (m_startPos.getZ() + m_endPos.getZ()) / 2 / CLHEP::cm);
114  m_WightedPos += position * (step->GetTotalEnergyDeposit());
115 
116  //Save Hit if track leaves volume or is killed
117  if (track.GetNextVolume() != track.GetVolume() || track.GetTrackStatus() >= fStopAndKill) {
118  int pdgCode = track.GetDefinition()->GetPDGEncoding();
119 
120  const G4VPhysicalVolume& v = * track.GetVolume();
121  G4ThreeVector posCell = v.GetTranslation();
122  // Get layer ID
123 
124  //if (v.GetName().find("Crystal") != std::string::npos) {
125  //CsiGeometryPar* eclp = CsiGeometryPar::Instance();
126  m_cellID = track.GetVolume()->GetCopyNo();
127 
128  double dTotalEnergy = 1 / m_energyDeposit; //avoid the error no match for 'operator/'
129  if (m_energyDeposit > 0.) {
130  saveSimHit(m_cellID, m_trackID, pdgCode, m_WightedTime / m_energyDeposit,
131  m_energyDeposit, m_momentum, m_WightedPos * dTotalEnergy);
132  //}
133  }
134 
135  //Reset TrackID
136  m_trackID = 0;
137  }
138 
139  return true;
140  }
141 
142  int SensitiveDetector::saveSimHit(
143  const G4int cellId,
144  const G4int trackID,
145  const G4int pid,
146  const G4double tof,
147  const G4double edep,
148  G4ThreeVector mom,
149  TVector3 posAve)
150  {
151 
152  //Get the datastore arraus
153  StoreArray<MCParticle> mcParticles;
154  StoreArray<ClawSimHit> simHits;
155  RelationArray relMClawmHit(mcParticles, simHits);
156 
157  StoreArray<ClawSimHit> ClawHits;
158  RelationArray clawSimHitRel(mcParticles, ClawHits);
159  TVector3 momentum(mom.getX() / CLHEP::GeV, mom.getY() / CLHEP::GeV, mom.getZ() / CLHEP::GeV);
160  ClawHits.appendNew(cellId, trackID, pid, tof / CLHEP::ns, edep / CLHEP::GeV, momentum, posAve);
161  int simhitNumber = ClawHits.getEntries() - 1;
162  B2DEBUG(150, "HitNumber: " << simhitNumber);
163  clawSimHitRel.add(trackID, simhitNumber);
164  return (simhitNumber);
165  }//saveSimHit
166 
167 
168  } //claw namespace
170 } //Belle2 namespace
Belle2::StoreArray::appendNew
T * appendNew()
Construct a new T object at the end of the array.
Definition: StoreArray.h:256
Belle2::RelationArray
Low-level class to create/modify relations between StoreArrays.
Definition: RelationArray.h:72
Belle2::claw::SensitiveDetector::SensitiveDetector
SensitiveDetector()
Constructor.
Definition: SensitiveDetector.cc:37
Belle2::Simulation::SensitiveDetectorBase::registerMCParticleRelation
static void registerMCParticleRelation(const std::string &name, RelationArray::EConsolidationAction ignoreAction=RelationArray::c_negativeWeight)
Register an relation involving MCParticles.
Definition: SensitiveDetectorBase.cc:24
Belle2::cave::SensitiveDetector::step
bool step(G4Step *step, G4TouchableHistory *) override
Step processing method.
Definition: SensitiveDetector.cc:56
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::RelationArray::add
void add(index_type from, index_type to, weight_type weight=1.0)
Add a new element to the relation.
Definition: RelationArray.h:291
Belle2::StoreArray< MCParticle >
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226