25#include "FairDetector.h"
26#include <fairlogger/Logger.h>
27#include "FairRootManager.h"
29#include "FairRuntimeDb.h"
30#include "FairVolume.h"
31#include "FairRootManager.h"
33#include "TGeoManager.h"
36#include "TGeoVolume.h"
38#include "TVirtualMC.h"
39#include "TVirtualMCStack.h"
88 LOG(info) <<
"Building FT3 Detector: From file";
89 LOG(info) <<
" FT3 detector configuration: " << configFileName;
90 std::ifstream ifs(configFileName.c_str());
92 LOG(fatal) <<
" Invalid FT3Base.configFile!";
95 float z_layer, r_in, r_out, Layerx2X0;
98 while (std::getline(ifs, tempstr)) {
99 if (tempstr[0] ==
'#') {
100 LOG(info) <<
" Comment: " << tempstr;
103 std::istringstream iss(tempstr);
118 LOG(info) <<
"Adding Layer " << layerName <<
" at z = " << z_layer <<
" ; direction = " << direction <<
" ; r_in = " << r_in <<
" ; r_out = " << r_out <<
" x/X0 = " << Layerx2X0;
119 auto& thisLayer = mLayers[0].emplace_back(direction, layerNumber, layerName, z_layer, r_in, r_out, Layerx2X0);
133 std::string configFileName =
"FT3_layout.cfg";
135 LOG(info) <<
"Exporting FT3 Detector layout to " << configFileName;
137 std::ofstream fOut(configFileName.c_str(), std::ios::out);
139 printf(
"Cannot open file\n");
142 fOut <<
"# z_layer r_in r_out Layerx2X0" << std::endl;
143 for (
auto layers_dir : mLayers) {
144 for (
auto layer : layers_dir) {
145 fOut <<
layer.getZ() <<
" " <<
layer.getInnerRadius() <<
" " <<
layer.getOuterRadius() <<
" " <<
layer.getx2X0() << std::endl;
156 LOG(info) <<
"Building FT3 Detector: Conical Telescope";
158 auto z_first =
param.z0;
159 auto z_length =
param.zLength;
160 auto etaIn =
param.etaIn;
161 auto etaOut =
param.etaOut;
162 auto Layerx2X0 =
param.Layerx2X0;
170 for (
int direction : {0, 1}) {
171 for (
int layerNumber = 0; layerNumber <
mNumberOfLayers; layerNumber++) {
173 mLayerName[direction][layerNumber] = layerName;
176 float layerZ = z_first + (layerNumber * z_length / (
mNumberOfLayers - 1)) * std::copysign(1, z_first);
177 float rIn = std::abs(layerZ * std::tan(2.f * std::atan(std::exp(-etaIn))));
178 float rOut = std::abs(layerZ * std::tan(2.f * std::atan(std::exp(-etaOut))));
179 auto& thisLayer = mLayers[direction].emplace_back(direction, layerNumber, layerName, layerZ, rIn, rOut, Layerx2X0);
190 LOG(info) <<
"Building FT3 Detector: V1";
193 float sensorThickness = 30.e-4;
194 float layersx2X0 = 1.e-2;
195 std::vector<std::array<float, 5>> layersConfig{
196 {26., .5, 3., 0.1f * layersx2X0},
197 {30., .5, 3., 0.1f * layersx2X0},
198 {34., .5, 3., 0.1f * layersx2X0},
199 {77., 3.5, 35., layersx2X0},
200 {100., 3.5, 35., layersx2X0},
201 {122., 3.5, 35., layersx2X0},
202 {150., 3.5, 80.f, layersx2X0},
203 {180., 3.5, 80.f, layersx2X0},
204 {220., 3.5, 80.f, layersx2X0},
205 {279., 3.5, 80.f, layersx2X0}};
213 for (
auto direction : {0, 1}) {
214 for (
int layerNumber = 0; layerNumber <
mNumberOfLayers; layerNumber++) {
217 mLayerName[direction][layerNumber] = layerName;
218 auto&
z = layersConfig[layerNumber][0];
220 auto& rIn = layersConfig[layerNumber][1];
221 auto& rOut = layersConfig[layerNumber][2];
222 auto&
x0 = layersConfig[layerNumber][3];
224 LOG(info) <<
"Adding Layer " << layerName <<
" at z = " <<
z;
226 auto& thisLayer = mLayers[direction].emplace_back(direction, layerNumber, layerName,
z, rIn, rOut,
x0);
237 LOG(info) <<
"Building FT3 Detector: V3b";
240 float sensorThickness = 30.e-4;
241 float layersx2X0 = 1.e-2;
242 std::vector<std::array<float, 5>> layersConfig{
243 {26., .5, 3., 0.1f * layersx2X0},
244 {30., .5, 3., 0.1f * layersx2X0},
245 {34., .5, 3., 0.1f * layersx2X0},
246 {77., 5.0, 35., layersx2X0},
247 {100., 5.0, 35., layersx2X0},
248 {122., 5.0, 35., layersx2X0},
249 {150., 5.5, 80.f, layersx2X0},
250 {180., 6.6, 80.f, layersx2X0},
251 {220., 8.1, 80.f, layersx2X0},
252 {279., 10.2, 80.f, layersx2X0},
253 {340., 12.5, 80.f, layersx2X0},
254 {400., 14.7, 80.f, layersx2X0}};
262 for (
auto direction : {0, 1}) {
263 for (
int layerNumber = 0; layerNumber <
mNumberOfLayers; layerNumber++) {
266 mLayerName[direction][layerNumber] = layerName;
267 auto&
z = layersConfig[layerNumber][0];
269 auto& rIn = layersConfig[layerNumber][1];
270 auto& rOut = layersConfig[layerNumber][2];
271 auto&
x0 = layersConfig[layerNumber][3];
273 LOG(info) <<
"Adding Layer " << layerName <<
" at z = " <<
z;
275 auto& thisLayer = mLayers[direction].emplace_back(direction, layerNumber, layerName,
z, rIn, rOut,
x0);
289 LOG(info) <<
"Building FT3 Detector: After Upgrade Days March 2024 version";
292 float sensorThickness = 30.e-4;
293 float layersx2X0 = 1.e-2;
294 std::vector<std::array<float, 5>> layersConfigCSide{
295 {77., 7.0, 35., layersx2X0},
296 {100., 7.0, 35., layersx2X0},
297 {122., 7.0, 35., layersx2X0},
298 {150., 7.0, 68.f, layersx2X0},
299 {180., 7.0, 68.f, layersx2X0},
300 {220., 7.0, 68.f, layersx2X0},
301 {260., 7.0, 68.f, layersx2X0},
302 {300., 7.0, 68.f, layersx2X0},
303 {350., 7.0, 68.f, layersx2X0}};
305 std::vector<std::array<float, 5>> layersConfigASide{
306 {77., 5.0, 35., layersx2X0},
307 {100., 5.0, 35., layersx2X0},
308 {122., 5.0, 35., layersx2X0},
309 {150., 5.0, 68.f, layersx2X0},
310 {180., 5.0, 68.f, layersx2X0},
311 {220., 5.0, 68.f, layersx2X0},
312 {260., 5.0, 68.f, layersx2X0},
313 {300., 5.0, 68.f, layersx2X0},
314 {350., 5.0, 68.f, layersx2X0}};
322 for (
auto direction : {0, 1}) {
323 for (
int layerNumber = 0; layerNumber <
mNumberOfLayers; layerNumber++) {
326 mLayerName[direction][layerNumber] = layerName;
327 float z, rIn, rOut,
x0;
328 if (direction == 0) {
329 z = layersConfigCSide[layerNumber][0];
330 rIn = layersConfigCSide[layerNumber][1];
331 rOut = layersConfigCSide[layerNumber][2];
332 x0 = layersConfigCSide[layerNumber][3];
333 }
else if (direction == 1) {
334 z = layersConfigASide[layerNumber][0];
335 rIn = layersConfigASide[layerNumber][1];
336 rOut = layersConfigASide[layerNumber][2];
337 x0 = layersConfigASide[layerNumber][3];
340 LOG(info) <<
"Adding Layer " << layerName <<
" at z = " <<
z;
342 auto& thisLayer = mLayers[direction].emplace_back(direction, layerNumber, layerName,
z, rIn, rOut,
x0);
352 LOG(info) <<
"Building FT3 Detector: Scoping document version";
355 float sensorThickness = 30.e-4;
356 float layersx2X0 = 1.e-2;
357 std::vector<std::array<float, 5>> layersConfig{
358 {26., .5, 2.5, 0.1f * layersx2X0},
359 {30., .5, 2.5, 0.1f * layersx2X0},
360 {34., .5, 2.5, 0.1f * layersx2X0},
361 {77., 5.0, 35., layersx2X0},
362 {100., 5.0, 35., layersx2X0},
363 {122., 5.0, 35., layersx2X0},
364 {150., 5.0, 68.f, layersx2X0},
365 {180., 5.0, 68.f, layersx2X0},
366 {220., 5.0, 68.f, layersx2X0},
367 {260., 5.0, 68.f, layersx2X0},
368 {300., 5.0, 68.f, layersx2X0},
369 {350., 5.0, 68.f, layersx2X0}};
377 for (
auto direction : {0, 1}) {
378 for (
int layerNumber = 0; layerNumber <
mNumberOfLayers; layerNumber++) {
381 mLayerName[direction][layerNumber] = layerName;
382 auto&
z = layersConfig[layerNumber][0];
384 auto& rIn = layersConfig[layerNumber][1];
385 auto& rOut = layersConfig[layerNumber][2];
386 auto&
x0 = layersConfig[layerNumber][3];
388 LOG(info) <<
"Adding Layer " << layerName <<
" at z = " <<
z;
390 auto& thisLayer = mLayers[direction].emplace_back(direction, layerNumber, layerName,
z, rIn, rOut,
x0);
403 auto& ft3BaseParam = FT3BaseParam::Instance();
405 if (ft3BaseParam.configFile !=
"") {
406 LOG(info) <<
"FT3 Geometry configuration file provided. Overriding FT3Base.geoModel configuration.";
407 buildFT3FromFile(ft3BaseParam.configFile);
410 switch (ft3BaseParam.geoModel) {
412 buildFT3NewVacuumVessel();
415 buildBasicFT3(ft3BaseParam);
418 LOG(fatal) <<
"Invalid Geometry.\n";
433 mLayerID = rhs.mLayerID;
434 mLayerName = rhs.mLayerName;
435 mNumberOfLayers = rhs.mNumberOfLayers;
466 mLayerID = rhs.mLayerID;
467 mLayerName = rhs.mLayerName;
468 mNumberOfLayers = rhs.mNumberOfLayers;
469 mLayers = rhs.mLayers;
470 mTrackData = rhs.mTrackData;
482 LOG(info) <<
"Initialize FT3 O2Detector";
486 defineSensitiveVolumes();
493 if (!(fMC->TrackCharge())) {
497 int lay = 0, volID = vol->getMCid();
498 while ((lay <= mLayerID.size()) && (volID != mLayerID[lay])) {
504 bool startHit =
false, stopHit =
false;
505 unsigned char status = 0;
506 if (fMC->IsTrackEntering()) {
507 status |= Hit::kTrackEntering;
509 if (fMC->IsTrackInside()) {
510 status |= Hit::kTrackInside;
512 if (fMC->IsTrackExiting()) {
513 status |= Hit::kTrackExiting;
515 if (fMC->IsTrackOut()) {
516 status |= Hit::kTrackOut;
518 if (fMC->IsTrackStop()) {
519 status |= Hit::kTrackStopped;
521 if (fMC->IsTrackAlive()) {
522 status |= Hit::kTrackAlive;
526 if ((status & Hit::kTrackEntering) || (status & Hit::kTrackInside && !mTrackData.mHitStarted)) {
528 }
else if ((status & (Hit::kTrackExiting | Hit::kTrackOut | Hit::kTrackStopped))) {
534 mTrackData.mEnergyLoss += fMC->Edep();
536 if (!(startHit | stopHit)) {
540 mTrackData.mEnergyLoss = 0.;
541 fMC->TrackMomentum(mTrackData.mMomentumStart);
542 fMC->TrackPosition(mTrackData.mPositionStart);
543 mTrackData.mTrkStatusStart = status;
544 mTrackData.mHitStarted =
true;
547 TLorentzVector positionStop;
548 fMC->TrackPosition(positionStop);
552 Hit* p = addHit(
stack->GetCurrentTrackNumber(), chipindex, mTrackData.mPositionStart.Vect(), positionStop.Vect(),
553 mTrackData.mMomentumStart.Vect(), mTrackData.mMomentumStart.E(), positionStop.T(),
554 mTrackData.mEnergyLoss, mTrackData.mTrkStatusStart, status);
559 stack->addHit(GetDetId());
566void Detector::createMaterials()
572 float tmaxfdSi = 0.1;
573 float stemaxSi = 0.0075;
574 float deemaxSi = 0.1;
575 float epsilSi = 1.0E-4;
578 float tmaxfdAir = 0.1;
579 float stemaxAir = .10000E+01;
580 float deemaxAir = 0.1;
581 float epsilAir = 1.0E-4;
582 float stminAir = 0.0;
585 float aAir[4] = {12.0107, 14.0067, 15.9994, 39.948};
586 float zAir[4] = {6., 7., 8., 18.};
587 float wAir[4] = {0.000124, 0.755267, 0.231781, 0.012827};
588 float dAir = 1.20479E-3;
591 o2::base::Detector::Medium(1,
"AIR$", 1, 0, ifield, fieldm, tmaxfdAir, stemaxAir, deemaxAir, epsilAir, stminAir);
594 o2::base::Detector::Medium(3,
"SILICON$", 3, 0, ifield, fieldm, tmaxfdSi, stemaxSi, deemaxSi, epsilSi, stminSi);
607 if (FairRootManager::Instance()) {
608 FairRootManager::Instance()->RegisterAny(
addNameTo(
"Hit").
data(), mHits, kTRUE);
631void Detector::createGeometry()
641 TGeoVolume* vALIC = gGeoManager->GetVolume(
"barrel");
643 LOG(fatal) <<
"Could not find the top volume";
646 TGeoVolume* A3IPvac = gGeoManager->GetVolume(
"OUT_PIPEVACUUM");
648 LOG(info) <<
"Running simulation with no beam pipe.";
651 LOG(
debug) <<
"FT3 createGeometry: "
652 << Form(
"gGeoManager name is %s title is %s", gGeoManager->GetName(), gGeoManager->GetTitle());
654 if (mLayers.size() == 2) {
656 for (
int direction : {0, 1}) {
657 std::string directionString = direction ?
"Forward" :
"Backward";
658 LOG(info) <<
"Creating FT3 " << directionString <<
" layers:";
659 for (
int iLayer = 0; iLayer < mLayers[direction].size(); iLayer++) {
660 mLayers[direction][iLayer].createLayer(volFT3);
663 vALIC->AddNode(volFT3, 2,
new TGeoTranslation(0., 30., 0.));
665 for (
int direction : {0, 1}) {
666 std::string directionString = direction ?
"Forward" :
"Backward";
667 LOG(info) <<
"Creating FT3 " << directionString <<
" layers:";
668 for (
int iLayer = 0; iLayer < mLayers[direction].size(); iLayer++) {
670 mLayers[direction][iLayer].createLayer(volIFT3);
672 mLayers[direction][iLayer].createLayer(volFT3);
676 A3IPvac->AddNode(volIFT3, 2,
new TGeoTranslation(0., 0., 0.));
677 vALIC->AddNode(volFT3, 2,
new TGeoTranslation(0., 30., 0.));
680 for (
auto direction : {0, 1}) {
681 std::string directionString = direction ?
"Forward" :
"Backward";
682 LOG(info) <<
"Registering FT3 " << directionString <<
" LayerIDs:";
683 for (
int iLayer = 0; iLayer < mLayers[direction].size(); iLayer++) {
685 mLayerID.push_back(layerID);
686 LOG(info) <<
" " << directionString <<
" layer " << iLayer <<
" LayerID " << layerID;
691 if (mLayers.size() == 1) {
692 LOG(info) <<
"Creating FT3 layers:";
694 for (
int iLayer = 0; iLayer < mLayers[0].size(); iLayer++) {
695 if (std::abs(mLayers[0][iLayer].getZ()) < 25) {
696 mLayers[0][iLayer].createLayer(volIFT3);
698 mLayers[0][iLayer].createLayer(volFT3);
701 A3IPvac->AddNode(volIFT3, 2,
new TGeoTranslation(0., 0., 0.));
702 vALIC->AddNode(volFT3, 2,
new TGeoTranslation(0., 30., 0.));
704 for (
int iLayer = 0; iLayer < mLayers[0].size(); iLayer++) {
705 mLayers[0][iLayer].createLayer(volFT3);
707 vALIC->AddNode(volFT3, 2,
new TGeoTranslation(0., 30., 0.));
709 LOG(info) <<
"Registering FT3 LayerIDs:";
710 for (
int iLayer = 0; iLayer < mLayers[0].size(); iLayer++) {
712 mLayerID.push_back(layerID);
713 LOG(info) <<
" mLayerID[" << iLayer <<
"] = " << layerID;
719void Detector::defineSensitiveVolumes()
721 TGeoManager* geoManager = gGeoManager;
725 LOG(info) <<
"Adding FT3 Sensitive Volumes";
728 if (mLayers.size() == 2) {
729 for (
int direction : {0, 1}) {
730 for (
int iLayer = 0; iLayer < mNumberOfLayers; iLayer++) {
733 LOG(info) <<
"Adding FT3 Sensitive Volume => " <<
v->GetName();
734 AddSensitiveVolume(
v);
739 if (mLayers.size() == 1) {
740 for (
int iLayer = 0; iLayer < mLayers[0].size(); iLayer++) {
743 LOG(info) <<
"Adding FT3 Sensitive Volume => " <<
v->GetName();
744 AddSensitiveVolume(
v);
751 const TVector3& startMom,
double startE,
double endTime,
double eLoss,
unsigned char startStatus,
752 unsigned char endStatus)
754 mHits->emplace_back(trackID, detID, startPos, endPos, startMom, startE, endTime, eLoss, startStatus, endStatus);
755 return &(mHits->back());
Definition of the Stack class.
Definition of the ITSMFT Hit class.
Definition of the FT3Layer class.
Definition of the GeometryTGeo class.
Definition of the Detector class.
Detector & operator=(const Detector &)
void Mixture(Int_t imat, const char *name, Float_t *a, Float_t *z, Float_t dens, Int_t nlmat, Float_t *wmat)
void Medium(Int_t numed, const char *name, Int_t nmat, Int_t isvol, Int_t ifield, Float_t fieldm, Float_t tmaxfd, Float_t stemax, Float_t deemax, Float_t epsil, Float_t stmin, Float_t *ubuf=nullptr, Int_t nbuf=0)
static void initFieldTrackingParams(int &mode, float &maxfield)
virtual void InitializeO2Detector()=0
void Material(Int_t imat, const char *name, Float_t a, Float_t z, Float_t dens, Float_t radl, Float_t absl, Float_t *buf=nullptr, Int_t nwbuf=0)
std::string addNameTo(const char *ext) const
std::vector< std::vector< TString > > mLayerName
void EndOfEvent() override
o2::itsmft::Hit * addHit(int trackID, int detID, const TVector3 &startPos, const TVector3 &endPos, const TVector3 &startMom, double startE, double endTime, double eLoss, unsigned char startStatus, unsigned char endStatus)
This method is an example of how to add your own point of type Hit to the clones array.
void Register() override
Registers the produced collections in FAIRRootManager.
void buildFT3FromFile(std::string)
void ConstructGeometry() override
Base class to create the detector geometry.
Detector()
Default constructor.
Bool_t ProcessHits(FairVolume *v=nullptr) override
This method is called for each step during simulation (see FairMCApplication::Stepping())
std::vector< Int_t > mLayerID
void buildBasicFT3(const FT3BaseParam ¶m)
void buildFT3NewVacuumVessel()
void Reset() override
Has to be called after each event to reset the containers.
void exportLayout()
access to geometry details
static const char * getFT3SensorPattern()
static const char * getFT3LayerPattern()
static GeometryTGeo * Instance()
static const char * getFT3VolPattern()
static const char * getFT3InnerVolPattern()
static ShmManager & Instance()
GLenum GLuint GLint GLint layer
GLdouble GLdouble GLdouble z
void freeSimVector(std::vector< T > *ptr)
a couple of static helper functions to create timestamp values for CCDB queries or override obsolete ...
std::string to_string(gsl::span< T, Size > span)
Common utility functions.
LOG(info)<< "Compressed in "<< sw.CpuTime()<< " s"