Belle II Software  release-05-01-25
SensitiveDetector.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010-2011 Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Hiroshi Nakano, Andreas Moll *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <ir/simulation/SensitiveDetector.h>
12 #include <ir/dataobjects/IRSimHit.h>
13 
14 #include <framework/logging/Logger.h>
15 #include <framework/datastore/StoreArray.h>
16 #include <framework/datastore/RelationArray.h>
17 #include <framework/gearbox/Unit.h>
18 
19 #include <TVector3.h>
20 
21 // Geant4
22 #include <G4Types.hh>
23 #include <G4ThreeVector.hh>
24 #include <G4Track.hh>
25 #include <G4VPhysicalVolume.hh>
26 #include <G4Step.hh>
27 #include <G4UserLimits.hh>
28 
29 #include <string>
30 
31 using namespace std;
32 
33 namespace Belle2 {
38  namespace ir {
39 
40  SensitiveDetector::SensitiveDetector() :
41  Simulation::SensitiveDetectorBase("IR ", Const::IR)
42  {
43  StoreArray<MCParticle> mcParticles;
44  StoreArray<IRSimHit> irSimHits;
45  RelationArray irSimHitRel(mcParticles, irSimHits);
46  registerMCParticleRelation(irSimHitRel);
47  }
48 
49  G4bool SensitiveDetector::step(G4Step* step, G4TouchableHistory*)
50  {
51  const G4Track& track = *step->GetTrack();
52  const G4int trackID = track.GetTrackID();
53  const G4int partPDGCode = track.GetDefinition()->GetPDGEncoding();
54 
55  const G4VPhysicalVolume& g4Volume = *track.GetVolume();
56  string Volname = g4Volume.GetName();
57 
58  const G4double depEnergy = step->GetTotalEnergyDeposit() * Unit::MeV;
59 
60  const G4StepPoint& preStep = *step->GetPreStepPoint();
61  const G4StepPoint& postStep = *step->GetPostStepPoint();
62  const G4ThreeVector& preStepPos = preStep.GetPosition();
63  const G4ThreeVector& postStepPos = postStep.GetPosition();
64  const G4ThreeVector momIn(preStep.GetMomentum());
65  const G4ThreeVector momOut(postStep.GetMomentum());
66  const G4ThreeVector preStepPosLocal = preStep.GetTouchableHandle()->GetHistory()->GetTopTransform().TransformPoint(preStepPos);
67  const G4ThreeVector postStepPosLocal = postStep.GetTouchableHandle()->GetHistory()->GetTopTransform().TransformPoint(postStepPos);
68  const G4ThreeVector momInLocal = preStep.GetTouchableHandle()->GetHistory()->GetTopTransform().TransformAxis(momIn);
69  const G4ThreeVector momOutLocal = preStep.GetTouchableHandle()->GetHistory()->GetTopTransform().TransformAxis(momOut);
70  TVector3 posInVec(preStepPosLocal.x() * Unit::mm, preStepPosLocal.y() * Unit::mm, preStepPosLocal.z() * Unit::mm);
71  TVector3 posOutVec(postStepPosLocal.x() * Unit::mm, postStepPosLocal.y() * Unit::mm, postStepPosLocal.z() * Unit::mm);
72  TVector3 momInVec(momInLocal.x() * Unit::MeV, momInLocal.y() * Unit::MeV, momInLocal.z() * Unit::MeV);
73  TVector3 momOutVec(momOutLocal.x() * Unit::MeV, momOutLocal.y() * Unit::MeV, momOutLocal.z() * Unit::MeV);
74 
75  // B2INFO("Step in volume: " << g4Volume.GetName())
76  // check that user limits are set properly
77  G4UserLimits* userLimits = g4Volume.GetLogicalVolume()->GetUserLimits();
78  if (userLimits) {
79  B2DEBUG(100, "Volume " << g4Volume.GetName() << ": max. allowed step set to " << userLimits->GetMaxAllowedStep(track));
80  }
81 
82  StoreArray<IRSimHit> irSimHits;
83  IRSimHit* simHit = irSimHits.appendNew(posInVec, momInVec, posOutVec, momOutVec, partPDGCode, depEnergy, Volname);
84 
85  // add relation to MCParticles
86  StoreArray<MCParticle> mcParticles;
87  RelationArray irSimHitRel(mcParticles, irSimHits);
88  irSimHitRel.add(trackID, simHit->getArrayIndex());
89 
90  return true;
91  }
92 
93  }
95 }
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::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::Unit::MeV
static const double MeV
[megaelectronvolt]
Definition: Unit.h:124
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::RelationsInterface::getArrayIndex
int getArrayIndex() const
Returns this object's array index (in StoreArray), or -1 if not found.
Definition: RelationsObject.h:387
Belle2::Unit::mm
static const double mm
[millimeters]
Definition: Unit.h:80
Belle2::StoreArray< MCParticle >
Belle2::IRSimHit
Class for saving raw hit data of the IR.
Definition: IRSimHit.h:33
Belle2::Const
This class provides a set of constants for the framework.
Definition: Const.h:36
Belle2::ir::SensitiveDetector::step
bool step(G4Step *aStep, G4TouchableHistory *) override
Process each step and calculate variables defined in IRSimHit.
Definition: SensitiveDetector.cc:49