Belle II Software  release-05-02-19
SensitiveAero.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Luka Santelj *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <arich/simulation/SensitiveAero.h>
12 #include <arich/dataobjects/ARICHAeroHit.h>
13 
14 #include <G4Step.hh>
15 #include <G4Track.hh>
16 
17 #include <framework/datastore/StoreArray.h>
18 #include <framework/datastore/RelationArray.h>
19 #include <framework/gearbox/Unit.h>
20 
21 #include <TVector3.h>
22 
23 using namespace std;
24 
25 namespace Belle2 {
30  namespace arich {
31 
32  SensitiveAero::SensitiveAero():
33  Simulation::SensitiveDetectorBase("ARICH", Const::ARICH)
34  {
35 
36  // registration of store arrays and relations
37 
38  StoreArray<MCParticle> mcParticles;
39  StoreArray<ARICHAeroHit> aeroHits;
40  aeroHits.registerInDataStore();
41  mcParticles.registerRelationTo(aeroHits);
42 
43  // additional registration of MCParticle relation (required for correct relations)
44 
45  RelationArray relation(mcParticles, aeroHits);
47 
48  }
49 
50  bool SensitiveAero::step(G4Step* aStep, G4TouchableHistory*)
51  {
52  // Get track parameters
53 
54  G4Track* aTrack = aStep->GetTrack();
55 
56  G4StepPoint* PostPosition = aStep->GetPostStepPoint();
57  G4ThreeVector worldPosition = PostPosition->GetPosition();
58  G4ParticleDefinition* particle = aTrack->GetDefinition();
59  G4double PDGCharge = particle->GetPDGCharge();
60  G4ThreeVector momentum = PostPosition->GetMomentum();
61 
62  // Save only tracks of charged particles
63  if (PDGCharge == 0) return (true);
64 
65  // Track parameters are saved at the entrance in aerogel
66 
67  if ((PostPosition->GetStepStatus() == fGeomBoundary) && (momentum.z() > 0.0)) {
68 
69  /* B2INFO ("SensAero: " << aTrack->GetDefinition()->GetParticleName()
70  << " " << aTrack->GetTrackID()
71  << " " << aTrack->GetParentID()
72  << " " << G4BestUnit(worldPosition,"Length")
73  << " " << G4BestUnit(aTrack->GetMomentum(), "Energy")
74  << " " << G4BestUnit(aTrack->GetGlobalTime(), "Time")
75  << " Edep is " << G4BestUnit(aStep->GetTotalEnergyDeposit(),"Energy"));
76  */
77 
78  int trackID = aTrack->GetTrackID();
79  int PDGEncoding = particle->GetPDGEncoding();
80 
81  TVector3 TPosition(worldPosition.x() * Unit::mm, worldPosition.y() * Unit::mm, worldPosition.z() * Unit::mm);
82  TVector3 TMomentum(momentum.x() * Unit::MeV, momentum.y() * Unit::MeV , momentum.z() * Unit::MeV);
83 
84  // write the hit in datastore"
85  StoreArray<ARICHAeroHit> aeroHits;
86  ARICHAeroHit* aeroHit = aeroHits.appendNew(PDGEncoding, TPosition, TMomentum);
87 
88  // Create relation to MCParticle
89  StoreArray<MCParticle> mcParticles;
90  RelationArray rel(mcParticles, aeroHits);
91  rel.add(trackID, aeroHit->getArrayIndex());
92 
93  }
94 
95  return true;
96 
97  }
98 
99  } // end of namespace arich
101 } // end of namespace Belle2
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::StoreArray::registerRelationTo
bool registerRelationTo(const StoreArray< TO > &toArray, DataStore::EDurability durability=DataStore::c_Event, DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut, const std::string &namedRelation="") const
Register a relation to the given StoreArray.
Definition: StoreArray.h:150
Belle2::arich::SensitiveAero::step
bool step(G4Step *aStep, G4TouchableHistory *) override
Process each step and calculate variables defined in PXDSimHit.
Definition: SensitiveAero.cc:50
Belle2::ARICHAeroHit
Datastore class that holds information on track parameters at the entrance in aerogel.
Definition: ARICHAeroHit.h: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::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::c_deleteElement
@ c_deleteElement
Delete the whole relation element if the original element got re-attributed.
Definition: RelationArray.h:91
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::Const
This class provides a set of constants for the framework.
Definition: Const.h:36