Skip to content

Commit a84a262

Browse files
Andreas Mathiswiechula
Andreas Mathis
authored and
wiechula
committed
Use new AliceO2 hits class + use reference for TVirtualMC
1 parent d3665c8 commit a84a262

7 files changed

Lines changed: 63 additions & 90 deletions

File tree

‎Detectors/TPC/simulation/include/TPCSimulation/Detector.h‎

Lines changed: 12 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -5,18 +5,17 @@
55
#include "DetectorsBase/Detector.h" // for Detector
66
#include "Rtypes.h" // for Int_t, Double32_t, Double_t, Bool_t, etc
77
#include "TLorentzVector.h" // for TLorentzVector
8-
#include "TVector3.h" // for TVector3
8+
#include "TClonesArray.h"
99
#include "TString.h"
1010

11+
#include "TPCSimulation/Point.h"
12+
1113
class FairVolume; // lines 10-10
12-
class TClonesArray; // lines 11-11
13-
namespace AliceO2 { namespace TPC { class Point; } } // lines 15-15
1414

1515
class AliTPCParam;
1616

1717
namespace AliceO2 {
1818
namespace TPC {
19-
class Point;
2019

2120
class Detector: public AliceO2::Base::Detector {
2221

@@ -63,10 +62,7 @@ class Detector: public AliceO2::Base::Detector {
6362
/** This method is an example of how to add your own point
6463
* of type DetectorPoint to the clones array
6564
*/
66-
Point* AddHit(Int_t trackID, Int_t detID,
67-
TVector3 pos, TVector3 mom,
68-
Double_t time, Double_t length,
69-
Double_t eLoss);
65+
Point* addHit(float x, float y, float z, float time, float nElectrons, float trackID, float detID);
7066

7167

7268
/// Copied from AliRoot - should go to someplace else
@@ -131,6 +127,14 @@ class Detector: public AliceO2::Base::Detector {
131127

132128
ClassDef(Detector,1)
133129
};
130+
131+
inline
132+
Point* Detector::addHit(float x, float y, float z, float time, float nElectrons, float trackID, float detID)
133+
{
134+
TClonesArray& clref = *mPointCollection;
135+
Int_t size = clref.GetEntriesFast();
136+
return new(clref[size]) Point(x, y, z, time, nElectrons, trackID, detID);
137+
}
134138
}
135139
}
136140

‎Detectors/TPC/simulation/include/TPCSimulation/Point.h‎

Lines changed: 13 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -3,20 +3,18 @@
33
#ifndef ALICEO2_TPC_POINT_H
44
#define ALICEO2_TPC_POINT_H
55

6+
#include "SimulationDataFormat/BaseHits.h"
67

7-
#include "FairMCPoint.h" // for FairMCPoint
8-
#include "Rtypes.h" // for Double_t, Int_t, Point::Class, ClassDef, etc
9-
#include "TVector3.h" // for TVector3
108
namespace AliceO2 {
119
namespace TPC {
1210

13-
class Point : public FairMCPoint
11+
class Point : public AliceO2::BasicXYZEHit<float>
1412
{
1513

1614
public:
1715

1816
/// Default constructor
19-
Point();
17+
Point() = default;
2018

2119
/// Constructor with arguments
2220
/// @param trackID Index of MCTrack
@@ -26,21 +24,27 @@ class Point : public FairMCPoint
2624
/// @param tof Time since event start [ns]
2725
/// @param length Track length since creation [cm]
2826
/// @param eLoss Energy deposit [GeV]
29-
Point(Int_t trackID, Int_t detID, TVector3 pos, TVector3 mom, Double_t tof, Double_t length, Double_t eLoss);
27+
Point(float x, float y, float z, float time, float nElectrons, float trackID, float detID);
3028

3129
/// Destructor
32-
virtual ~Point();
30+
virtual ~Point() = default;
3331

3432
/// Output to screen
35-
virtual void Print(const Option_t* opt) const;
33+
virtual void Print(const Option_t* opt) const override;
3634

3735
private:
3836
/// Copy constructor
3937
Point(const Point& point);
4038
Point operator=(const Point& point);
4139

42-
ClassDef(AliceO2::TPC::Point,1)
40+
ClassDefOverride(AliceO2::TPC::Point,1)
4341
};
42+
43+
inline
44+
Point::Point(float x, float y, float z, float time, float nElectrons, float trackID, float detID)
45+
: BasicXYZEHit<float>(x, y, z, time, nElectrons, trackID, detID)
46+
{}
47+
4448
}
4549
}
4650

‎Detectors/TPC/simulation/src/Detector.cxx‎

Lines changed: 20 additions & 34 deletions
Original file line numberDiff line numberDiff line change
@@ -24,7 +24,6 @@
2424
#include "TSystem.h"
2525
#include "TClonesArray.h"
2626
#include "TVirtualMC.h"
27-
#include "TVectorD.h"
2827

2928
#include "TFile.h"
3029

@@ -176,12 +175,14 @@ void Detector::SetSpecialPhysicsCuts()
176175

177176
Bool_t Detector::ProcessHits(FairVolume* vol)
178177
{
178+
static auto *refMC = TVirtualMC::GetMC();
179+
179180
/* This method is called from the MC stepping for the sensitive volume only */
180181
// LOG(INFO) << "TPC::ProcessHits" << FairLogger::endl;
181-
if(static_cast<int>(TVirtualMC::GetMC()->TrackCharge()) == 0) {
182+
if(static_cast<int>(refMC->TrackCharge()) == 0) {
182183

183184
// set a very large step size for neutral particles
184-
TVirtualMC::GetMC()->SetMaxStep(1.e10);
185+
refMC->SetMaxStep(1.e10);
185186
return kFALSE; // take only charged particles
186187
}
187188

@@ -197,21 +198,21 @@ Bool_t Detector::ProcessHits(FairVolume* vol)
197198
// https://indico.cern.ch/event/316891/contributions/732168/
198199

199200
TLorentzVector momentum;
200-
TVirtualMC::GetMC()->TrackMomentum(momentum);
201-
const Double_t rnd = TVirtualMC::GetMC()->GetRandom()->Rndm();
201+
refMC->TrackMomentum(momentum);
202+
const Double_t rnd = refMC->GetRandom()->Rndm();
202203
if(mSimulationType == SimulationType::GEANT3) {
203204

204205
// betagamma = p/m
205-
Float_t betaGamma = momentum.P()/TVirtualMC::GetMC()->TrackMass();
206+
Float_t betaGamma = momentum.P()/refMC->TrackMass();
206207
betaGamma = TMath::Max(betaGamma, static_cast<Float_t>(7.e-3)); // protection against too small bg
207208

208209
// NPRIM etc. are defined in "TPCSimulation/Constants.h"
209210
const Float_t pp = NPRIM * BetheBlochAleph(betaGamma, BBPARAM[0], BBPARAM[1], BBPARAM[2], BBPARAM[3], BBPARAM[4]);
210211

211-
TVirtualMC::GetMC()->SetMaxStep(-TMath::Log(rnd)/pp);
212+
refMC->SetMaxStep(-TMath::Log(rnd)/pp);
212213
} else {
213214

214-
TVirtualMC::GetMC()->SetMaxStep(0.2+(2.*rnd-1.)*0.05); // 2 mm +- rndm*0.5mm step
215+
refMC->SetMaxStep(0.2+(2.*rnd-1.)*0.05); // 2 mm +- rndm*0.5mm step
215216
}
216217

217218
// CONVERT THE ENERGY LOSS TO IONIZATION
@@ -227,15 +228,15 @@ Bool_t Detector::ProcessHits(FairVolume* vol)
227228
Int_t nel=0;
228229
if(mSimulationType == SimulationType::GEANT3) {
229230

230-
nel = 1 + static_cast<int>((TVirtualMC::GetMC()->Edep()-IPOT) / WION);
231+
nel = 1 + static_cast<int>((refMC->Edep()-IPOT) / WION);
231232
// LOG(INFO) << "TPC::AddHit" << FairLogger::endl << "GEANT3: Nelectrons: " << nel << FairLogger::endl;
232233
} else {
233234

234-
const Double_t meanIon = TVirtualMC::GetMC()->Edep() / (WION*SCALEWIONG4);
235+
const Double_t meanIon = refMC->Edep() / (WION*SCALEWIONG4);
235236
if(meanIon > 0)
236237
nel = static_cast<int>(FANOFACTORG4 * Gamma(meanIon/FANOFACTORG4));
237238
// LOG(INFO) << "TPC::AddHit" << FairLogger::endl << "GEANT4: Eloss: "
238-
// << TVirtualMC::GetMC()->Edep() << ", Nelectrons: "
239+
// << refMC->Edep() << ", Nelectrons: "
239240
// << nel << FairLogger::endl;
240241
}
241242

@@ -246,14 +247,12 @@ Bool_t Detector::ProcessHits(FairVolume* vol)
246247

247248
// ADD HIT
248249
TLorentzVector position;
249-
TVirtualMC::GetMC()->TrackPosition(position);
250-
Double32_t time = TVirtualMC::GetMC()->TrackTime() * 1.0e09;
251-
Double32_t length = TVirtualMC::GetMC()->TrackLength();
252-
int trackNumberID = TVirtualMC::GetMC()->GetStack()->GetCurrentTrackNumber();
253-
int volumeID = vol->getMCid();
254-
AddHit(trackNumberID, volumeID, TVector3(position.X(), position.Y(), position.Z()),
255-
TVector3(momentum.Px(), momentum.Py(), momentum.Pz()), time, length,
256-
nel);
250+
refMC->TrackPosition(position);
251+
float time = refMC->TrackTime() * 1.0e09;
252+
int trackID = refMC->GetStack()->GetCurrentTrackNumber();
253+
int detID = vol->getMCid();
254+
addHit(position.X(), position.Y(), position.Z(), time, nel, trackID, detID);
255+
257256
//LOG(INFO) << "TPC::AddHit" << FairLogger::endl
258257
//<< " -- " << trackNumberID <<"," << volumeID << " " << vol->GetName()
259258
//<< ", Pos: (" << position.X() << ", " << position.Y() <<", "<< position.Z()<< ", " << r << ") "
@@ -262,7 +261,7 @@ Bool_t Detector::ProcessHits(FairVolume* vol)
262261
//nel << FairLogger::endl;
263262

264263
// Increment number of Detector det points in TParticle
265-
AliceO2::Data::Stack* stack = (AliceO2::Data::Stack*)TVirtualMC::GetMC()->GetStack();
264+
AliceO2::Data::Stack* stack = (AliceO2::Data::Stack*)refMC->GetStack();
266265
stack->AddPoint(kAliTpc);
267266

268267
return kTRUE;
@@ -288,7 +287,7 @@ void Detector::Register()
288287
only during the simulation.
289288
*/
290289

291-
FairGenericRootManager::Instance()->Register("TPCPoint", "TPC",mPointCollection, kTRUE);
290+
FairRootManager::Instance()->Register("TPCPoint", "TPC",mPointCollection, kTRUE);
292291

293292
}
294293

@@ -3041,19 +3040,6 @@ void Detector::DefineSensitiveVolumes()
30413040
}
30423041
}
30433042

3044-
3045-
Point* Detector::AddHit(Int_t trackID, Int_t detID,
3046-
TVector3 pos, TVector3 mom,
3047-
Double_t time, Double_t length,
3048-
Double_t eLoss)
3049-
{
3050-
TClonesArray& clref = *mPointCollection;
3051-
Int_t size = clref.GetEntriesFast();
3052-
return new(clref[size]) Point(trackID, detID, pos, mom,
3053-
time, length, eLoss);
3054-
}
3055-
3056-
30573043
Double_t Detector::BetheBlochAleph(Double_t bg, Double_t kp1, Double_t kp2, Double_t kp3, Double_t kp4, Double_t kp5){
30583044
Double_t beta = bg/TMath::Sqrt(1.+ bg*bg);
30593045

‎Detectors/TPC/simulation/src/Digitizer.cxx‎

Lines changed: 12 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -44,12 +44,17 @@ DigitContainer *Digitizer::Process(TClonesArray *points)
4444
{
4545
mDigitContainer->reset();
4646
const Mapper& mapper = Mapper::instance();
47+
FairRootManager *mgr = FairRootManager::Instance();
48+
49+
const float eventTime = mgr->GetEventTime() * 0.001; /// transform in us
4750

4851
/// \todo static_thread for thread savety?
4952
static GEMAmplification gemAmplification;
5053
static ElectronTransport electronTransport;
5154
static PadResponse padResponse;
5255

56+
const int MCEventID = mgr->GetEntryNr();
57+
5358
static std::array<float, mNShapedPoints> signalArray;
5459

5560
for(auto pointObject : *points) {
@@ -58,10 +63,9 @@ DigitContainer *Digitizer::Process(TClonesArray *points)
5863
const GlobalPosition3D posEle(inputpoint->GetX(), inputpoint->GetY(), inputpoint->GetZ());
5964

6065
// The energy loss stored is really nElectrons
61-
int nPrimaryElectrons = static_cast<int>(inputpoint->GetEnergyLoss());
66+
const int nPrimaryElectrons = static_cast<int>(inputpoint->GetEnergyLoss());
6267

63-
int MCEventID = inputpoint->GetEventID();
64-
int MCTrackID = inputpoint->GetTrackID();
68+
const int MCTrackID = inputpoint->GetTrackID();
6569

6670
/// Loop over electrons
6771
/// \todo can be vectorized?
@@ -72,7 +76,8 @@ DigitContainer *Digitizer::Process(TClonesArray *points)
7276
const GlobalPosition3D posEleDiff = electronTransport.getElectronDrift(posEle);
7377

7478
/// \todo Time management in continuous mode (adding the time of the event?)
75-
const float driftTime = getTime(posEleDiff.getZ());
79+
const float driftTime = getTime(posEleDiff.getZ()) + inputpoint->GetTime() * 0.001; /// in us
80+
const float absoluteTime = driftTime + eventTime;
7681

7782
/// Attachment
7883
if(electronTransport.isElectronAttachment(driftTime)) continue;
@@ -108,18 +113,18 @@ DigitContainer *Digitizer::Process(TClonesArray *points)
108113
if(mDebugFlagPRF) {
109114
/// \todo Write out the debug output
110115
GEMresponse.CRU = digiPos.getCRU().number();
111-
GEMresponse.time = driftTime;
116+
GEMresponse.time = absoluteTime;
112117
GEMresponse.row = row;
113118
GEMresponse.pad = pad;
114119
GEMresponse.nElectrons = nElectronsGEM * normalizedPadResponse;
115120
//mDebugTreePRF->Fill();
116121
}
117122

118123
const float ADCsignal = SAMPAProcessing::getADCvalue(nElectronsGEM * normalizedPadResponse);
119-
SAMPAProcessing::getShapedSignal(ADCsignal, driftTime, signalArray);
124+
SAMPAProcessing::getShapedSignal(ADCsignal, absoluteTime, signalArray);
120125

121126
for(float i=0; i<mNShapedPoints; ++i) {
122-
float time = driftTime + i * ZBINWIDTH;
127+
const float time = absoluteTime + i * ZBINWIDTH;
123128
mDigitContainer->addDigit(MCEventID, MCTrackID, digiPos.getCRU().number(), getTimeBinFromTime(time), row, pad, signalArray[i]);
124129
}
125130

‎Detectors/TPC/simulation/src/DigitizerTask.cxx‎

Lines changed: 1 addition & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -132,8 +132,7 @@ void DigitizerTask::fillHitArrayFromFile()
132132
// printf("Filling hit %d (event %d)\n", ihit, fEvent);
133133

134134
const int size = dummy.GetEntriesFast();
135-
new(dummy[size]) Point(fTrack, 98, TVector3(fX, fY, fZ), TVector3(0,0,0),
136-
fTime, 0., fQ*WION);
135+
new(dummy[size]) Point(fX, fY, fZ, fTime, fQ, fTrack, 98);
137136
}
138137

139138
printf("Converted hits: %d\n", dummy.GetEntriesFast());

‎Detectors/TPC/simulation/src/Point.cxx‎

Lines changed: 4 additions & 30 deletions
Original file line numberDiff line numberDiff line change
@@ -5,39 +5,13 @@ using std::endl;
55

66
using namespace AliceO2::TPC;
77

8-
// ----- Default constructor -------------------------------------------
9-
Point::Point()
10-
: FairMCPoint()
11-
{
12-
}
13-
// -------------------------------------------------------------------------
14-
15-
// ----- Standard constructor ------------------------------------------
16-
Point::Point(Int_t trackID, Int_t detID,
17-
TVector3 pos, TVector3 mom,
18-
Double_t tof, Double_t length,
19-
Double_t eLoss)
20-
: FairMCPoint(trackID, detID, pos, mom, tof, length, eLoss)
21-
{
22-
}
23-
// -------------------------------------------------------------------------
24-
25-
// ----- Destructor ----------------------------------------------------
26-
Point::~Point() = default;
27-
// -------------------------------------------------------------------------
28-
29-
// ----- Public method Print -------------------------------------------
308
void Point::Print(const Option_t* opt) const
319
{
32-
cout << "-I- Point: O2tpc point for track " << fTrackID
33-
<< " in detector " << fDetectorID << endl;
34-
cout << " Position (" << fX << ", " << fY << ", " << fZ
10+
cout << "-I- Point: O2tpc point for track " << GetTrackID()
11+
<< " in detector " << GetDetectorID() << endl;
12+
cout << " Position (" << GetX() << ", " << GetY() << ", " << GetZ()
3513
<< ") cm" << endl;
36-
cout << " Momentum (" << fPx << ", " << fPy << ", " << fPz
37-
<< ") GeV" << endl;
38-
cout << " Time " << fTime << " ns, Length " << fLength
39-
<< " cm, Energy loss " << fELoss*1.0e06 << " keV" << endl;
14+
cout << " Time " << GetTime() << " ns, n electrons " << GetEnergyLoss() << endl;
4015
}
41-
// -------------------------------------------------------------------------
4216

4317
ClassImp(Point)

‎cmake/O2Dependencies.cmake‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -577,6 +577,7 @@ o2_define_bucket(
577577
${CMAKE_SOURCE_DIR}/Detectors/Base/include
578578
${CMAKE_SOURCE_DIR}/Detectors/Passive/include
579579
${CMAKE_SOURCE_DIR}/Detectors/TPC/base/include
580+
${CMAKE_SOURCE_DIR}/DataFormats/simulation/include
580581
)
581582

582583

0 commit comments

Comments
 (0)