21#include <TGeoVolume.h>
22#include <TVirtualMC.h>
23#include <TVirtualMCStack.h>
25#include <FairVolume.h>
39 return 2. * (10. + radius * std::cos(2 * std::atan(std::exp(-eta))));
43 :
o2::base::DetImpl<Detector>(
"TRK", true),
49Detector::Detector(
bool active)
50 :
o2::base::DetImpl<Detector>(
"TRK", true),
56 if (trkPars.configFile !=
"") {
57 configFromFile(trkPars.configFile);
64 LOGP(info,
"Summary of TRK configuration:");
65 for (
auto&
layer : mLayers) {
66 LOGP(info,
"Layer: {} name: {} r: {} cm | z: {} cm | thickness: {} cm",
layer->getNumber(),
layer->getName(),
layer->getInnerRadius(),
layer->getZ(),
layer->getChipThickness());
70Detector::Detector(
const Detector&
other)
71 :
o2::base::DetImpl<Detector>(
other),
84void Detector::ConstructGeometry()
90void Detector::configMLOT()
92 auto& trkPars = TRKBaseParam::Instance();
96 const std::vector<float> rInn{7.f, 9.f, 12.f, 20.f, 30.f, 45.f, 60.f, 80.f};
97 const float thick = 100.e-3;
99 switch (trkPars.layoutMLOT) {
101 const std::vector<float>
length{128.35f, 128.35f, 128.35f, 128.35f, 128.35f, 256.7f, 256.7f, 256.7f};
102 LOGP(warning,
"Loading cylindrical configuration for ALICE3 TRK");
103 for (
int i{0};
i < constants::ML::nLayers + constants::OT::nLayers; ++
i) {
105 mLayers.push_back(std::make_unique<TRKCylindricalLayer>(
i,
name, rInn[
i],
length[
i], thick, MatBudgetParamMode::Thickness));
110 const std::vector<float> tiltAngles{11.2f, 11.9f, 11.4f, 0.f, 0.f, 0.f, 0.f, 0.f};
112 const std::vector<int> nStaves{10, 14, 18, 26, 38, 32, 42, 56};
114 const std::vector<int> nMods{10, 10, 10, 10, 10, 20, 20, 20};
116 const std::vector<float> stagOffsets{0.f, 0.f, 0.f, 1.17f, 0.89f};
118 LOGP(warning,
"Loading segmented configuration for ALICE3 TRK");
119 for (
int i{0};
i < constants::ML::nLayers + constants::OT::nLayers; ++
i) {
121 if (
i < constants::ML::nLayers) {
122 mLayers.push_back(std::make_unique<TRKMLLayer>(
i,
name, rInn[
i], stagOffsets[
i], tiltAngles[
i], nStaves[
i], nMods[
i], thick, MatBudgetParamMode::Thickness));
124 mLayers.push_back(std::make_unique<TRKOTLayer>(
i,
name, rInn[
i], tiltAngles[
i], nStaves[
i], nMods[
i], thick, MatBudgetParamMode::Thickness));
130 LOGP(fatal,
"Unknown option {} for configMLOT",
static_cast<int>(trkPars.layoutMLOT));
135void Detector::configFromFile(std::string fileName)
138 std::ifstream confFile(fileName);
139 if (!confFile.good()) {
140 LOGP(fatal,
"File {} not found, aborting.", fileName);
143 auto& trkPars = TRKBaseParam::Instance();
147 LOGP(info,
"Overriding geometry of ALICE3 TRK using {} file.", fileName);
150 std::vector<float> tmpBuff;
152 while (std::getline(confFile, line)) {
153 if (line[0] ==
'/') {
157 std::stringstream ss(line);
160 while (getline(ss, substr,
'\t')) {
161 tmpBuff.push_back(std::stof(substr));
166 switch (trkPars.layoutMLOT) {
169 if (tmpBuff.size() < 3) {
170 LOGP(fatal,
"Invalid configuration for cylindrical layer {}: insufficient parameters.", layerCount);
175 if (tmpBuff.size() >= 4) {
179 mLayers.push_back(std::make_unique<TRKCylindricalLayer>(layerCount,
name, tmpBuff[0], tmpBuff[1], tmpBuff[2],
mode));
193 if (tmpBuff.size() < 5) {
194 LOGP(fatal,
"Invalid configuration for segmented layer {}: missing base parameters.", layerCount);
197 float rInn = tmpBuff[0];
198 float thick = tmpBuff[1];
199 float tiltAngle = tmpBuff[2];
200 int nStaves =
static_cast<int>(tmpBuff[3]);
201 int nMods =
static_cast<int>(tmpBuff[4]);
206 if (layerCount < constants::ML::nLayers) {
208 if (tmpBuff.size() < 6) {
209 LOGP(fatal,
"Invalid configuration for ML layer {}: stagOffset is missing.", layerCount);
211 float stagOffset = tmpBuff[5];
213 if (tmpBuff.size() >= 7) {
217 mLayers.push_back(std::make_unique<TRKMLLayer>(layerCount,
name, rInn, stagOffset, tiltAngle, nStaves, nMods, thick,
mode));
220 if (tmpBuff.size() >= 6) {
224 mLayers.push_back(std::make_unique<TRKOTLayer>(layerCount,
name, rInn, tiltAngle, nStaves, nMods, thick,
mode));
229 LOGP(fatal,
"Unknown option {} for configMLOT",
static_cast<int>(trkPars.layoutMLOT));
237void Detector::configToFile(std::string fileName)
239 LOGP(info,
"Exporting TRK Detector layout to {}", fileName);
240 std::ofstream conFile(fileName.c_str(), std::ios::out);
241 conFile <<
"/// TRK configuration file: inn_radius z_length lay_thickness" << std::endl;
242 for (
const auto&
layer : mLayers) {
243 conFile <<
layer->getInnerRadius() <<
"\t" <<
layer->getZ() <<
"\t" <<
layer->getChipThickness() << std::endl;
247void Detector::configServices()
252void Detector::createMaterials()
258 float tmaxfdSi = 0.1;
259 float stemaxSi = 0.0075;
260 float deemaxSi = 0.1;
261 float epsilSi = 1.0E-4;
264 float tmaxfdAir = 0.1;
265 float stemaxAir = .10000E+01;
266 float deemaxAir = 0.1;
267 float epsilAir = 1.0E-4;
268 float stminAir = 0.0;
270 float tmaxfdCer = 0.1;
271 float stemaxCer = .10000E+01;
272 float deemaxCer = 0.1;
273 float epsilCer = 1.0E-4;
274 float stminCer = 0.0;
277 float aAir[4] = {12.0107, 14.0067, 15.9994, 39.948};
278 float zAir[4] = {6., 7., 8., 18.};
279 float wAir[4] = {0.000124, 0.755267, 0.231781, 0.012827};
280 float dAir = 1.20479E-3;
283 float aCf[2] = {12.0107, 1.00794};
284 float zCf[2] = {6., 1.};
287 o2::base::Detector::Medium(1,
"AIR$", 1, 0, ifield, fieldm, tmaxfdAir, stemaxAir, deemaxAir, epsilAir, stminAir);
290 o2::base::Detector::Medium(3,
"SILICON$", 3, 0, ifield, fieldm, tmaxfdSi, stemaxSi, deemaxSi, epsilSi, stminSi);
293void Detector::createGeometry()
295 TGeoManager* geoManager = gGeoManager;
296 TGeoVolume* vALIC = geoManager->GetVolume(
"barrel");
298 LOGP(fatal,
"Could not find barrel volume while constructing TRK geometry");
300 new TGeoVolumeAssembly(GeometryTGeo::getTRKVolPattern());
301 TGeoVolume* vTRK = geoManager->GetVolume(GeometryTGeo::getTRKVolPattern());
302 vALIC->AddNode(vTRK, 2,
new TGeoTranslation(0, 30., 0));
304 char vstrng[100] =
"TRKVol";
305 vTRK->SetTitle(vstrng);
307 for (
auto&
layer : mLayers) {
308 layer->createLayer(vTRK);
312 mServices.createServices(vTRK);
316 auto& trkPars = TRKBaseParam::Instance();
320 switch (trkPars.layoutVD) {
322 LOG(info) <<
"Building VD with IRIS4 layout";
326 LOG(info) <<
"Building VD with IRIS fully cylindrical layout";
330 LOG(info) <<
"Building VD with IRIS fully cylindrical layout with 3 inclined walls";
334 LOG(info) <<
"Building VD with IRIS5 layout";
338 LOG(info) <<
"Building VD with IRIS4a layout";
342 LOG(fatal) <<
"Unknown VD layout option: " <<
static_cast<int>(trkPars.layoutVD);
348 mNumberOfVolumesVD =
static_cast<int>(regs.size());
349 mNumberOfVolumes = mNumberOfVolumesVD + mLayers.size();
350 mSensorName.resize(mNumberOfVolumes);
354 for (
const auto& sensor : regs) {
355 mSensorName[VDvolume] = sensor.name;
360 for (
int i = 0;
i < mLayers.size();
i++) {
361 mSensorName[VDvolume++].Form(
"%s%d", GeometryTGeo::getTRKSensorPattern(),
i);
364 for (
auto vd : mSensorName) {
365 std::cout <<
"Volume name: " << vd << std::endl;
368 mServices.excavateFromVacuum(
"IRIS_CUTOUTsh");
369 mServices.registerVacuum(vTRK);
372void Detector::InitializeO2Detector()
374 LOG(info) <<
"Initialize TRK O2Detector";
375 mGeometryTGeo = GeometryTGeo::Instance();
376 defineSensitiveVolumes();
378 mSensorID.resize(mNumberOfVolumes);
379 for (
int i = 0;
i < mNumberOfVolumes;
i++) {
380 mSensorID[
i] = gMC ? TVirtualMC::GetMC()->VolId(mSensorName[
i]) : 0;
381 LOGP(info,
"{}: mSensorID={}, mSensorName={}",
i, mSensorID[
i], mSensorName[
i].Data());
385void Detector::defineSensitiveVolumes()
387 TGeoManager* geoManager = gGeoManager;
391 LOGP(info,
"Adding TRK Sensitive Volumes");
394 for (
const auto& s :
o2::trk::vdSensorRegistry()) {
395 TGeoVolume*
v = gGeoManager->GetVolume(s.name.c_str());
397 LOGP(warning,
"VD sensor volume '{}' not found", s.name);
400 LOGP(info,
"Adding VD Sensitive Volume {}",
v->GetName());
401 AddSensitiveVolume(
v);
404 mFirstOrLastLayers.push_back(
s.name);
409 for (
int j{0};
j < mLayers.size();
j++) {
410 volumeName = GeometryTGeo::getTRKSensorPattern() + TString::Itoa(
j, 10);
411 if (
j == mLayers.size() - 1) {
412 mFirstOrLastLayers.push_back(volumeName.Data());
414 LOGP(info,
"Trying {}", volumeName.Data());
415 v = geoManager->GetVolume(volumeName.Data());
416 LOGP(info,
"Adding TRK Sensitive Volume {}",
v->GetName());
417 AddSensitiveVolume(
v);
421void Detector::EndOfEvent() { Reset(); }
423void Detector::Register()
429 if (FairRootManager::Instance()) {
430 FairRootManager::Instance()->RegisterAny(addNameTo(
"Hit").
data(), mHits,
true);
434void Detector::Reset()
441bool Detector::InsideFirstOrLastLayer(std::string layerName)
444 for (
auto& firstOrLastLayer : mFirstOrLastLayers) {
445 if (firstOrLastLayer == layerName) {
453bool Detector::ProcessHits(FairVolume* vol)
456 if (!(fMC->TrackCharge())) {
463 int volID = vol->getMCid();
465 bool notSens =
false;
466 while ((volume < mNumberOfVolumes) && (notSens = (volID != mSensorID[volume]))) {
474 if (volume < mNumberOfVolumesVD) {
478 layer = volume - mNumberOfVolumesVD;
484 if (fMC->IsTrackExiting() && InsideFirstOrLastLayer(vol->GetName())) {
487 tr.setTrackID(
stack->GetCurrentTrackNumber());
488 tr.setUserId(volume);
489 stack->addTrackReference(tr);
491 bool startHit =
false, stopHit =
false;
492 unsigned char status = 0;
493 if (fMC->IsTrackEntering()) {
494 status |= Hit::kTrackEntering;
496 if (fMC->IsTrackInside()) {
497 status |= Hit::kTrackInside;
499 if (fMC->IsTrackExiting()) {
500 status |= Hit::kTrackExiting;
502 if (fMC->IsTrackOut()) {
503 status |= Hit::kTrackOut;
505 if (fMC->IsTrackStop()) {
506 status |= Hit::kTrackStopped;
508 if (fMC->IsTrackAlive()) {
509 status |= Hit::kTrackAlive;
513 if ((status & Hit::kTrackEntering) || (status & Hit::kTrackInside && !mTrackData.mHitStarted)) {
515 }
else if ((status & (Hit::kTrackExiting | Hit::kTrackOut | Hit::kTrackStopped))) {
521 mTrackData.mEnergyLoss += fMC->Edep();
523 if (!(startHit | stopHit)) {
528 mTrackData.mEnergyLoss = 0.;
529 fMC->TrackMomentum(mTrackData.mMomentumStart);
530 fMC->TrackPosition(mTrackData.mPositionStart);
531 mTrackData.mTrkStatusStart = status;
532 mTrackData.mHitStarted =
true;
535 TLorentzVector positionStop;
536 fMC->TrackPosition(positionStop);
539 int stave(0), halfstave(0), mod(0), chip(0);
541 auto& trkPars = TRKBaseParam::Instance();
545 fMC->CurrentVolOffID(1, chip);
546 fMC->CurrentVolOffID(2, mod);
547 if (mGeometryTGeo->getNumberOfHalfStaves(
layer) == 2) {
548 fMC->CurrentVolOffID(3, halfstave);
549 fMC->CurrentVolOffID(4, stave);
550 }
else if (mGeometryTGeo->getNumberOfHalfStaves(
layer) == 1) {
551 fMC->CurrentVolOffID(3, stave);
553 LOGP(fatal,
"Wrong number of halfstaves for layer {}",
layer);
558 unsigned short chipID = mGeometryTGeo->getChipIndex(subDetID, volume,
layer, stave, halfstave, mod, chip);
564 Hit*
p = addHit(
stack->GetCurrentTrackNumber(), chipID, mTrackData.mPositionStart.Vect(), positionStop.Vect(),
565 mTrackData.mMomentumStart.Vect(), mTrackData.mMomentumStart.E(), positionStop.T(),
566 mTrackData.mEnergyLoss, mTrackData.mTrkStatusStart, status);
571 stack->addHit(GetDetId());
577o2::trk::Hit* Detector::addHit(
int trackID,
unsigned short detID,
const TVector3& startPos,
const TVector3& endPos,
578 const TVector3& startMom,
double startE,
double endTime,
double eLoss,
unsigned char startStatus,
579 unsigned char endStatus)
581 mHits->emplace_back(trackID, detID, startPos, endPos, startMom, startE, endTime, eLoss, startStatus, endStatus);
582 return &(mHits->back());
585void Detector::Print(FairVolume* vol,
int volume,
int subDetID,
int layer,
int stave,
int halfstave,
int mod,
int chip,
int chipID)
const
588 LOG(info) <<
"Current volume name: " << fMC->CurrentVolName() <<
" and ID " << fMC->CurrentVolID(currentVol);
589 LOG(info) <<
"volume: " << volume <<
"/" << mNumberOfVolumes - 1;
591 auto& trkPars = TRKBaseParam::Instance();
595 LOG(info) <<
"off volume name 1 " << fMC->CurrentVolOffName(1) <<
" chip: " << chip;
596 LOG(info) <<
"SubDetector ID: " << subDetID <<
" Layer: " <<
layer <<
" Chip ID: " << chipID;
598 LOG(info) <<
"off volume name 1 " << fMC->CurrentVolOffName(1) <<
" chip: " << chip;
599 LOG(info) <<
"off volume name 2 " << fMC->CurrentVolOffName(2) <<
" module: " << mod;
600 if (mGeometryTGeo->getNumberOfHalfStaves(
layer) == 2) {
601 LOG(info) <<
"off volume name 3 " << fMC->CurrentVolOffName(3) <<
" halfstave: " << halfstave;
602 LOG(info) <<
"off volume name 4 " << fMC->CurrentVolOffName(4) <<
" stave: " << stave;
603 LOG(info) <<
"SubDetector ID: " << subDetID <<
" Layer: " <<
layer <<
" staveinLayer: " << stave <<
" Chip ID: " << chipID;
604 }
else if (mGeometryTGeo->getNumberOfHalfStaves(
layer) == 1) {
605 LOG(info) <<
"off volume name 3 " << fMC->CurrentVolOffName(3) <<
" stave: " << stave;
606 LOG(info) <<
"SubDetector ID: " << subDetID <<
" Layer: " <<
layer <<
" staveinLayer: " << stave <<
" Chip ID: " << chipID;
611 LOG(info) <<
"SubDetector ID: " << subDetID <<
" Chip ID: " << chipID;
Definition of the Stack class.
Definition of the TRK Hit class.
ClassImp(o2::trk::Detector)
o2::base::Detector * create_detector_trk(bool active)
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)
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)
static const TRKBaseParam & Instance()
static o2::base::Detector * create(bool active)
static ShmManager & Instance()
GLuint const GLchar * name
GLuint GLsizei GLsizei * length
GLenum GLuint GLint GLint layer
void createGeometry(TGeoManager &geom, TGeoVolume &topVolume)
std::vector< VDSensorDesc > & vdSensorRegistry()
float getDetLengthFromEta(const float eta, const float radius)
void createIRISGeometry3InclinedWalls(TGeoVolume *motherVolume)
void createIRISGeometryFullCyl(TGeoVolume *motherVolume)
void createIRIS4aGeometry(TGeoVolume *motherVolume)
void createIRIS4Geometry(TGeoVolume *motherVolume)
void createIRIS5Geometry(TGeoVolume *motherVolume)
@ kIRISFullCyl3InclinedWalls
void clearVDSensorRegistry()
void freeSimVector(std::vector< T > *ptr)
std::vector< T > * createSimVector()
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.
VectorOfTObjectPtrs other
LOG(info)<< "Compressed in "<< sw.CpuTime()<< " s"