diff --git a/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/GeometryTGeo.h b/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/GeometryTGeo.h index 0e9ff8727a977..f900065ad738a 100644 --- a/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/GeometryTGeo.h +++ b/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/GeometryTGeo.h @@ -48,8 +48,13 @@ class GeometryTGeo : public o2::detectors::DetMatrixCache static const char* getTRKPetalDiskPattern() { return sPetalDiskName.c_str(); } static const char* getTRKPetalLayerPattern() { return sPetalLayerName.c_str(); } static const char* getTRKStavePattern() { return sStaveName.c_str(); } + static const char* getTRKHalfStavePattern() { return sHalfStaveName.c_str(); } + static const char* getTRKModulePattern() { return sModuleName.c_str(); } static const char* getTRKChipPattern() { return sChipName.c_str(); } static const char* getTRKSensorPattern() { return sSensorName.c_str(); } + static const char* getTRKDeadzonePattern() { return sDeadzoneName.c_str(); } + static const char* getTRKMetalStackPattern() { return sMetalStackName.c_str(); } + static const char* getTRKWrapVolPattern() { return sWrapperVolumeName.c_str(); } int getNumberOfChips() const { return mSize; } @@ -63,6 +68,8 @@ class GeometryTGeo : public o2::detectors::DetMatrixCache int extractNumberOfChipsPerPetalVD() const; int extractNumberOfStavesMLOT(int lay) const; int extractNumberOfHalfStavesMLOT(int lay) const; + int extractNumberOfModulesMLOT(int lay) const; + int extractNumberOfChipsMLOT(int lay) const; /// Extract number following the prefix in the name string int extractVolumeCopy(const char* name, const char* prefix) const; @@ -75,33 +82,39 @@ class GeometryTGeo : public o2::detectors::DetMatrixCache void setOwner(bool v) { mOwner = v; } void Print(Option_t* opt = "") const; - void PrintChipID(int index, int subDetID, int petalcase, int disk, int lay, int stave, int halfstave) const; + void PrintChipID(int index, int subDetID, int petalcase, int disk, int lay, int stave, int halfstave, int mod, int chip) const; - int getLayer(int index) const; - int getStave(int index) const; - int getHalfStave(int index) const; int getSubDetID(int index) const; int getPetalCase(int index) const; int getDisk(int index) const; + int getLayer(int index) const; + int getStave(int index) const; + int getHalfStave(int index) const; + int getModule(int index) const; + int getChip(int index) const; - /// This routine computes the chip index number from the subDetID, petal, disk, layer, stave /// TODO: retrieve also from chip when chips will be available + /// This routine computes the chip index number from the subDetID, petal, disk, layer, stave, half stave, module, chip /// \param int subDetID The subdetector ID, 0 for VD, 1 for MLOT /// \param int petalcase The petal case number for VD, from 0 to 3 /// \param int disk The disk number for VD, from 0 to 5 /// \param int lay The layer number. Starting from 0 both for VD and MLOT /// \param int stave The stave number for MLOT. Starting from 0 /// \param int halfstave The half stave number for MLOT. Can be 0 or 1 - int getChipIndex(int subDetID, int petalcase, int disk, int lay, int stave, int halfstave) const; + /// \param int module The module number for MLOT, from 0 to 10 (or 20) + /// \param int chip The chip number for MLOT, from 0 to 8 + unsigned short getChipIndex(int subDetID, int petalcase, int disk, int lay, int stave, int halfstave, int mod, int chip) const; - /// This routine computes the chip index number from the subDetID, volume, layer, stave /// TODO: retrieve also from chip when chips will be available + /// This routine computes the chip index number from the subDetID, volume, layer, stave, half stave, module, chip /// \param int subDetID The subdetector ID, 0 for VD, 1 for MLOT /// \param int volume is needed only with the current configuration for VD where each single element is a volume. // TODO: when the geometry naming scheme will be changed, change this method /// \param int lay The layer number for the MLOT. In the current configuration for VD this is not needed. // TODO: when the geometry naming scheme will be changed, change this method /// \param int stave The stave number in each layer for MLOT. Starting from 0. /// \param int halfstave The half stave number for MLOT. Can be 0 or 1 - int getChipIndex(int subDetID, int volume, int lay, int stave, int halfstave) const; + /// \param int module The module number for MLOT, from 0 to 10 (or 20) + /// \param int chip The chip number for MLOT, from 0 to 8 + unsigned short getChipIndex(int subDetID, int volume, int lay, int stave, int halfstave, int mod, int chip) const; - /// This routine computes subDetID, petal, disk, layer, stave given the chip index number /// TODO: copute also from chip when chips will be available + /// This routine computes subDetID, petal, disk, layer, stave, half stave, module, chip, given the chip index number /// \param int index The chip index number, starting from 0 /// \param int subDetID The subdetector ID, 0 for VD, 1 for MLOT /// \param int petalcase The petal case number for VD, from 0 to 3 @@ -109,10 +122,12 @@ class GeometryTGeo : public o2::detectors::DetMatrixCache /// \param int lay The layer number. Starting from 0 both for VD and MLOT /// \param int stave The stave number for MLOT. Starting from 0 /// \param int halfstave The half stave number for MLOT. Can be 0 or 1 - bool getChipID(int index, int& subDetID, int& petalcase, int& disk, int& lay, int& stave, int& halfstave) const; + /// \param int module The module number for MLOT, from 0 to 10 (or 20) + /// \param int chip The chip number for MLOT, from 0 to 8 + bool getChipID(int index, int& subDetID, int& petalcase, int& disk, int& lay, int& stave, int& halfstave, int& mod, int& chip) const; - int getLastChipIndex(int lay) const { return mLastChipIndex[lay]; } - int getFirstChipIndex(int lay, int petalcase, int subDetID) const + unsigned short getLastChipIndex(int lay) const { return mLastChipIndex[lay]; } + unsigned short getFirstChipIndex(int lay, int petalcase, int subDetID) const { /// Get the first chip index of the active petal (VD) or layer (MLOT) if (subDetID == 0) { // VD @@ -138,7 +153,8 @@ class GeometryTGeo : public o2::detectors::DetMatrixCache static const char* composeSymNameLayer(int d, int layer); static const char* composeSymNameStave(int d, int layer); - static const char* composeSymNameChip(int d, int lr); + static const char* composeSymNameModule(int d, int layer); + static const char* composeSymNameChip(int d, int layer); static const char* composeSymNameSensor(int d, int layer); protected: @@ -151,25 +167,36 @@ class GeometryTGeo : public o2::detectors::DetMatrixCache static std::string sPetalDiskName; static std::string sPetalLayerName; static std::string sStaveName; + static std::string sHalfStaveName; + static std::string sModuleName; static std::string sChipName; static std::string sSensorName; - static std::string sWrapperVolumeName; ///< Wrapper volume name + static std::string sDeadzoneName; + static std::string sMetalStackName; + + static std::string sWrapperVolumeName; ///< Wrapper volume name, not implemented at the moment Int_t mNumberOfLayersMLOT; ///< number of layers Int_t mNumberOfActivePartsVD; ///< number of layers Int_t mNumberOfLayersVD; ///< number of layers Int_t mNumberOfPetalsVD; ///< number of Petals = chip in each VD layer Int_t mNumberOfDisksVD; ///< number of Disks = 6 - std::vector mLastChipIndex; ///< max ID of the detctor in the petal(VD) or layer(MLOT) - std::vector mLastChipIndexVD; ///< max ID of the detctor in the layer for the VD - std::vector mLastChipIndexMLOT; ///< max ID of the detctor in the layer for the MLOT + std::vector mNumberOfStaves; ///< Number Of Staves per layer in ML/OT + std::vector mNumberOfHalfStaves; ///< Number Of Half staves in each stave of the layer in ML/OT + std::vector mNumberOfModules; ///< Number Of Modules per stave (half stave) in ML/OT + std::vector mNumberOfChips; ///< number of chips per module in ML/OT std::vector mNumberOfChipsPerLayerVD; ///< number of chips per layer VD ( = number of petals) - std::vector mNumberOfChipsPerLayerMLOT; ///< number of chips per layer MLOT ( = 1 for the moment) + std::vector mNumberOfChipsPerLayerMLOT; ///< number of chips per layer MLOT std::vector mNumbersOfChipPerDiskVD; ///< numbersOfChipPerDiskVD std::vector mNumberOfChipsPerPetalVD; ///< numbersOfChipPerPetalVD - std::vector mNumberOfStaves; ///< Number Of Staves per layer in ML/OT - std::vector mNumberOfHalfStaves; ///< Number Of Staves in each stave of the layer in ML/OT - std::array mLayerToWrapper; ///< Layer to wrapper correspondence + // std::vector mNumberOfChipsPerStave; ///< number of chips per stave in ML/OT + // std::vector mNumberOfChipsPerHalfStave; ///< number of chips per half stave in ML/OT + // std::vector mNumberOfChipsPerModule; ///< number of chips per module in ML/OT + std::vector mLastChipIndex; ///< max ID of the detctor in the petal(VD) or layer(MLOT) + std::vector mLastChipIndexVD; ///< max ID of the detctor in the layer for the VD + std::vector mLastChipIndexMLOT; ///< max ID of the detctor in the layer for the MLOT + + std::array mLayerToWrapper; ///< Layer to wrapper correspondence, not implemented yet bool mOwner = true; //! is it owned by the singleton? diff --git a/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/Specs.h b/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/Specs.h index bd95e5207b7ee..559d8f6154c59 100644 --- a/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/Specs.h +++ b/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/Specs.h @@ -22,7 +22,7 @@ // Each TGeoShape has the following properties // length: dimension in z-axis // width: dimension in xy-axes -// color: for visulisation +// color: for visualisation namespace o2::trk::constants { // Default unit of TGeo = cm @@ -84,13 +84,11 @@ constexpr int nRows{static_cast(width / pitchX)}; // constexpr int nCols{static_cast(length / pitchZ)}; // number of rows in the chip constexpr double totalThickness{silicon::thickness + metalstack::thickness}; // total thickness of the chip /// Set to 0 for the moment, to be adjusted with the actual design of the chip if needed -static constexpr float PassiveEdgeReadOut = 0.f; // width of the readout edge (Passive bottom) -static constexpr float PassiveEdgeTop = 0.f; // Passive area on top -static constexpr float PassiveEdgeSide = 0.f; // width of Passive area on left/right of the sensor +static constexpr float passiveEdgeReadOut{1.5 * mm}; // width of the readout edge } // namespace chip namespace gaps { -constexpr double interChips{0.2 * mm}; // gap between the chips +constexpr double interChips{50 * mu}; // gap between the chips constexpr double outerEdgeLongSide{1 * mm}; // gap between the chips and the outer edges (long side) constexpr double outerEdgeShortSide{0.1 * mm}; // gap between the chips and the outer edges (short side) } // namespace gaps diff --git a/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/TRKBaseParam.h b/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/TRKBaseParam.h index 63c95b1e6b2f6..3f3f656c4b417 100644 --- a/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/TRKBaseParam.h +++ b/Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/TRKBaseParam.h @@ -20,6 +20,11 @@ namespace o2 namespace trk { +enum eOverallGeom { + kDefaultRadii = 0, // After Upgrade Days March 2024 + kModRadii, +}; + enum eLayout { kCylinder = 0, kTurboStaves, @@ -31,6 +36,8 @@ struct TRKBaseParam : public o2::conf::ConfigurableParamHelper { float serviceTubeX0 = 0.02f; // X0 Al2O3 Bool_t irisOpen = false; + eOverallGeom overallGeom = kDefaultRadii; // Overall geometry option, to be used in Detector::buildTRKMiddleOuterLayers + eLayout layoutML = kCylinder; // Type of segmentation for the middle layers eLayout layoutOL = kCylinder; // Type of segmentation for the outer layers diff --git a/Detectors/Upgrades/ALICE3/TRK/base/src/GeometryTGeo.cxx b/Detectors/Upgrades/ALICE3/TRK/base/src/GeometryTGeo.cxx index 9325f5079375d..7fe9717313343 100644 --- a/Detectors/Upgrades/ALICE3/TRK/base/src/GeometryTGeo.cxx +++ b/Detectors/Upgrades/ALICE3/TRK/base/src/GeometryTGeo.cxx @@ -13,6 +13,8 @@ #include #include "TRKBase/SegmentationChip.h" +#include + using Segmentation = o2::trk::SegmentationChip; namespace o2 @@ -29,8 +31,12 @@ std::string GeometryTGeo::sPetalName = "PETALCASE"; std::string GeometryTGeo::sPetalDiskName = "DISK"; std::string GeometryTGeo::sPetalLayerName = "LAYER"; std::string GeometryTGeo::sStaveName = "TRKStave"; +std::string GeometryTGeo::sHalfStaveName = "TRKHalfStave"; +std::string GeometryTGeo::sModuleName = "TRKModule"; std::string GeometryTGeo::sChipName = "TRKChip"; std::string GeometryTGeo::sSensorName = "TRKSensor"; +std::string GeometryTGeo::sDeadzoneName = "TRKDeadzone"; +std::string GeometryTGeo::sMetalStackName = "TRKMetalStack"; std::string GeometryTGeo::sWrapperVolumeName = "TRKUWrapVol"; ///< Wrapper volume name, not implemented at the moment @@ -77,18 +83,24 @@ void GeometryTGeo::Build(int loadTrans) mNumberOfStaves.resize(mNumberOfLayersMLOT); mNumberOfHalfStaves.resize(mNumberOfLayersMLOT); - mLastChipIndex.resize(mNumberOfPetalsVD + mNumberOfLayersMLOT); - mLastChipIndexVD.resize(mNumberOfPetalsVD); - mLastChipIndexMLOT.resize(mNumberOfLayersMLOT); /// ML and OT are part of TRK as the same detector, without disks + mNumberOfModules.resize(mNumberOfLayersMLOT); + mNumberOfChips.resize(mNumberOfLayersMLOT); + mNumberOfChipsPerLayerVD.resize(mNumberOfLayersVD); mNumberOfChipsPerLayerMLOT.resize(mNumberOfLayersMLOT); mNumbersOfChipPerDiskVD.resize(mNumberOfDisksVD); mNumberOfChipsPerPetalVD.resize(mNumberOfPetalsVD); + mLastChipIndex.resize(mNumberOfPetalsVD + mNumberOfLayersMLOT); + mLastChipIndexVD.resize(mNumberOfPetalsVD); + mLastChipIndexMLOT.resize(mNumberOfLayersMLOT); /// ML and OT are part of TRK as the same detector, without disks + for (int i = 0; i < mNumberOfLayersMLOT; i++) { std::cout << "Layer MLOT: " << i << std::endl; mNumberOfStaves[i] = extractNumberOfStavesMLOT(i); mNumberOfHalfStaves[i] = extractNumberOfHalfStavesMLOT(i); + mNumberOfModules[i] = extractNumberOfModulesMLOT(i); + mNumberOfChips[i] = extractNumberOfChipsMLOT(i); } int numberOfChipsTotal = 0; @@ -103,13 +115,13 @@ void GeometryTGeo::Build(int loadTrans) /// filling the information for the MLOT for (int i = 0; i < mNumberOfLayersMLOT; i++) { - mNumberOfChipsPerLayerMLOT[i] = extractNumberOfStavesMLOT(i) * extractNumberOfHalfStavesMLOT(i); // for the moment, considering 1 half stave = 1 chip. TODO: add the final segmentation in chips + mNumberOfChipsPerLayerMLOT[i] = mNumberOfStaves[i] * mNumberOfHalfStaves[i] * mNumberOfModules[i] * mNumberOfChips[i]; numberOfChipsTotal += mNumberOfChipsPerLayerMLOT[i]; mLastChipIndex[i + mNumberOfPetalsVD] = numberOfChipsTotal - 1; mLastChipIndexMLOT[i] = numberOfChipsTotal - 1; } - setSize(numberOfChipsTotal); /// temporary, number of chips = number of staves and active parts + setSize(numberOfChipsTotal); fillMatrixCache(loadTrans); } @@ -132,9 +144,7 @@ int GeometryTGeo::getPetalCase(int index) const int subDetID = getSubDetID(index); if (subDetID == 1) { return -1; - } - - else if (index <= mLastChipIndexVD[mNumberOfPetalsVD - 1]) { + } else if (index <= mLastChipIndexVD[mNumberOfPetalsVD - 1]) { while (index > mLastChipIndexVD[petalcase]) { petalcase++; } @@ -142,6 +152,22 @@ int GeometryTGeo::getPetalCase(int index) const return petalcase; } +//__________________________________________________________________________ +int GeometryTGeo::getDisk(int index) const +{ + int subDetID = getSubDetID(index); + int petalcase = getPetalCase(index); + + if (subDetID == 0) { /// VD + if (index % mNumberOfChipsPerPetalVD[petalcase] < mNumberOfLayersVD) { + return -1; /// layers + } + return (index % mNumberOfChipsPerPetalVD[petalcase]) - mNumberOfLayersVD; + } + + return -1; /// not found or ML/OT +} + //__________________________________________________________________________ int GeometryTGeo::getLayer(int index) const { @@ -158,7 +184,7 @@ int GeometryTGeo::getLayer(int index) const while (index > mLastChipIndex[lay]) { lay++; } - return lay - mNumberOfPetalsVD; /// numeration of MLOT layesrs starting from 0 + return lay - mNumberOfPetalsVD; /// numeration of MLOT layers starting from 0 } return -1; /// -1 if not found } @@ -174,9 +200,23 @@ int GeometryTGeo::getStave(int index) const } else if (subDetID == 1) { /// MLOT int lay = getLayer(index); index -= getFirstChipIndex(lay, petalcase, subDetID); // get the index of the sensing element in the layer - return index / mNumberOfHalfStaves[lay]; + + const int Nhs = mNumberOfHalfStaves[lay]; + const int Nmod = mNumberOfModules[lay]; + const int Nchip = mNumberOfChips[lay]; + + if (Nhs == 2) { + int chipsPerModule = Nchip; + int chipsPerHalfStave = Nmod * chipsPerModule; + int chipsPerStave = Nhs * chipsPerHalfStave; + return index / chipsPerStave; + } else if (Nhs == 1) { + int chipsPerModule = Nchip; + int chipsPerStave = Nmod * chipsPerModule; + return index / chipsPerStave; + } } - return -1; /// not found + return -1; } //__________________________________________________________________________ @@ -185,36 +225,89 @@ int GeometryTGeo::getHalfStave(int index) const int subDetID = getSubDetID(index); int lay = getLayer(index); int petalcase = getPetalCase(index); - int stave = getStave(index); if (subDetID == 0) { /// VD return -1; } else if (subDetID == 1) { /// MLOT int lay = getLayer(index); index -= getFirstChipIndex(lay, petalcase, subDetID); // get the index of the sensing element in the layer - return index % 2; /// 0 = half stave left, 1 = half stave right, as geometry is filled /// TODO: generalize once chips will be in place. Can it be working also with chips? + + const int Nhs = mNumberOfHalfStaves[lay]; + const int Nmod = mNumberOfModules[lay]; + const int Nchip = mNumberOfChips[lay]; + + int chipsPerModule = Nchip; + int chipsPerHalfStave = Nmod * chipsPerModule; + int chipsPerStave = Nhs * chipsPerHalfStave; + + int rem = index % chipsPerStave; + return rem / chipsPerHalfStave; // 0 = left, 1 = right } - return -1; /// not found + return -1; } //__________________________________________________________________________ -int GeometryTGeo::getDisk(int index) const +int GeometryTGeo::getModule(int index) const { int subDetID = getSubDetID(index); + int lay = getLayer(index); int petalcase = getPetalCase(index); if (subDetID == 0) { /// VD - if (index % mNumberOfChipsPerPetalVD[petalcase] < mNumberOfLayersVD) { - return -1; /// layers + return -1; + } else if (subDetID == 1) { /// MLOT + int lay = getLayer(index); + index -= getFirstChipIndex(lay, petalcase, subDetID); // get the index of the sensing element in the layer + + const int Nhs = mNumberOfHalfStaves[lay]; + const int Nmod = mNumberOfModules[lay]; + const int Nchip = mNumberOfChips[lay]; + + if (Nhs == 2) { + int chipsPerModule = Nchip; + int chipsPerHalfStave = Nmod * chipsPerModule; + int rem = index % (Nhs * chipsPerHalfStave); + rem = rem % chipsPerHalfStave; + return rem / chipsPerModule; + } else if (Nhs == 1) { + int chipsPerModule = Nchip; + int rem = index % (Nmod * chipsPerModule); + return rem / chipsPerModule; } - return (index % mNumberOfChipsPerPetalVD[petalcase]) - mNumberOfLayersVD; } + return -1; +} - return -1; /// not found or ML/OT +//__________________________________________________________________________ +int GeometryTGeo::getChip(int index) const +{ + int subDetID = getSubDetID(index); + int lay = getLayer(index); + int petalcase = getPetalCase(index); + + if (subDetID == 0) { /// VD + return -1; + } else if (subDetID == 1) { /// MLOT + int lay = getLayer(index); + index -= getFirstChipIndex(lay, petalcase, subDetID); // get the index of the sensing element in the layer + + const int Nhs = mNumberOfHalfStaves[lay]; + const int Nmod = mNumberOfModules[lay]; + const int Nchip = mNumberOfChips[lay]; + + if (Nhs == 2) { + int chipsPerModule = Nchip; + return index % chipsPerModule; + } else if (Nhs == 1) { + int chipsPerModule = Nchip; + return index % chipsPerModule; + } + } + return -1; } //__________________________________________________________________________ -int GeometryTGeo::getChipIndex(int subDetID, int petalcase, int disk, int lay, int stave, int halfstave) const +unsigned short GeometryTGeo::getChipIndex(int subDetID, int petalcase, int disk, int lay, int stave, int halfstave, int mod, int chip) const { if (subDetID == 0) { // VD if (lay == -1) { // disk @@ -222,41 +315,70 @@ int GeometryTGeo::getChipIndex(int subDetID, int petalcase, int disk, int lay, i } else { // layer return getFirstChipIndex(lay, petalcase, subDetID) + lay; } - } else if (subDetID == 1) { // MLOT - if (mNumberOfHalfStaves[lay] == 2) { // staggered geometry - return getFirstChipIndex(lay, petalcase, subDetID) + stave * mNumberOfHalfStaves[lay] + halfstave; - } else if (mNumberOfHalfStaves[lay] == 1) { // turbo geometry - return getFirstChipIndex(lay, petalcase, subDetID) + stave; + } else if (subDetID == 1) { // MLOT + const int Nhs = mNumberOfHalfStaves[lay]; // 1 or 2 + const int Nmod = mNumberOfModules[lay]; // module per half-stave (per stave if Nhs==1) + const int Nchip = mNumberOfChips[lay]; // chips per module + + if (Nhs == 2) { // staggered geometry: layer -> stave -> halfstave -> mod -> chip + int chipsPerModule = Nchip; + int chipsPerHalfStave = Nmod * chipsPerModule; + int chipsPerStave = Nhs * chipsPerHalfStave; + return getFirstChipIndex(lay, petalcase, subDetID) + stave * chipsPerStave + halfstave * chipsPerHalfStave + mod * chipsPerModule + chip; + } else if (Nhs == 1) { // turbo geometry: layer -> stave -> mod -> chip (no halfstave) + int chipsPerModule = Nchip; + int chipsPerStave = Nmod * chipsPerModule; + return getFirstChipIndex(lay, petalcase, subDetID) + stave * chipsPerStave + mod * chipsPerModule + chip; } } - return -1; // not found + + LOGP(warning, "Chip index not found for subDetID %d, petalcase %d, disk %d, layer %d, stave %d, halfstave %d, module %d, chip %d, returning numeric limit", subDetID, petalcase, disk, lay, stave, halfstave, mod, chip); + return std::numeric_limits::max(); // not found } //__________________________________________________________________________ -int GeometryTGeo::getChipIndex(int subDetID, int volume, int lay, int stave, int halfstave) const +unsigned short GeometryTGeo::getChipIndex(int subDetID, int volume, int lay, int stave, int halfstave, int mod, int chip) const { if (subDetID == 0) { // VD return volume; /// In the current configuration for VD, each volume is the sensor element = chip. // TODO: when the geometry naming scheme will be changed, change this method - } else if (subDetID == 1) { // MLOT - if (mNumberOfHalfStaves[lay] == 2) { // staggered geometry - return getFirstChipIndex(lay, -1, subDetID) + stave * mNumberOfHalfStaves[lay] + halfstave; - } else if (mNumberOfHalfStaves[lay] == 1) { // turbo geometry - return getFirstChipIndex(lay, -1, subDetID) + stave; + } else if (subDetID == 1) { // MLOT + const int Nhs = mNumberOfHalfStaves[lay]; // 1 or 2 + const int Nmod = mNumberOfModules[lay]; // module per half-stave (per stave if Nhs==1) + const int Nchip = mNumberOfChips[lay]; // chips per module + + if (Nhs == 2) { // staggered geometry: layer -> stave -> halfstave -> mod -> chip + int chipsPerModule = Nchip; + int chipsPerHalfStave = Nmod * chipsPerModule; + int chipsPerStave = Nhs * chipsPerHalfStave; + return getFirstChipIndex(lay, -1, subDetID) + stave * chipsPerStave + halfstave * chipsPerHalfStave + mod * chipsPerModule + chip; + } else if (Nhs == 1) { // turbo geometry: layer -> stave -> mod -> chip (no halfstave) + int chipsPerModule = Nchip; + int chipsPerStave = Nmod * chipsPerModule; + return getFirstChipIndex(lay, -1, subDetID) + stave * chipsPerStave + mod * chipsPerModule + chip; } } - return -1; // not found + + LOGP(warning, "Chip index not found for subDetID %d, volume %d, layer %d, stave %d, halfstave %d, module %d, chip %d, returning numeric limit", subDetID, volume, lay, stave, halfstave, mod, chip); + return std::numeric_limits::max(); // not found } //__________________________________________________________________________ -bool GeometryTGeo::getChipID(int index, int& subDetID, int& petalcase, int& disk, int& lay, int& stave, int& halfstave) const +bool GeometryTGeo::getChipID(int index, int& subDetID, int& petalcase, int& disk, int& lay, int& stave, int& halfstave, int& mod, int& chip) const { subDetID = getSubDetID(index); petalcase = getPetalCase(index); disk = getDisk(index); lay = getLayer(index); stave = getStave(index); + if (mNumberOfHalfStaves[lay] == 2) { + halfstave = getHalfStave(index); + } else { + halfstave = 0; // if not staggered geometry, return 0 + } halfstave = getHalfStave(index); + mod = getModule(index); + chip = getChip(index); return kTRUE; } @@ -265,10 +387,10 @@ bool GeometryTGeo::getChipID(int index, int& subDetID, int& petalcase, int& disk TString GeometryTGeo::getMatrixPath(int index) const { - int subDetID, petalcase, disk, layer, stave, halfstave; //// TODO: add chips in a second step - getChipID(index, subDetID, petalcase, disk, layer, stave, halfstave); + int subDetID, petalcase, disk, layer, stave, halfstave, mod, chip; + getChipID(index, subDetID, petalcase, disk, layer, stave, halfstave, mod, chip); - // PrintChipID(index, subDetID, petalcase, disk, layer, stave, halfstave); + // PrintChipID(index, subDetID, petalcase, disk, layer, stave, halfstave, mod, chip); // TString path = "/cave_1/barrel_1/TRKV_2/TRKLayer0_1/TRKStave0_1/TRKChip0_1/TRKSensor0_1/"; /// dummy path, to be used for tests TString path = Form("/cave_1/barrel_1/%s_2/", GeometryTGeo::getTRKVolPattern()); @@ -284,15 +406,15 @@ TString GeometryTGeo::getMatrixPath(int index) const path += Form("%s%d_%s%d_%s%d_1/", getTRKPetalPattern(), petalcase, getTRKPetalLayerPattern(), layer, getTRKChipPattern(), layer); // PETALCASEx_LAYERy_TRKChipy_1 path += Form("%s%d_%s%d_%s%d_1/", getTRKPetalPattern(), petalcase, getTRKPetalLayerPattern(), layer, getTRKSensorPattern(), layer); // PETALCASEx_LAYERy_TRKSensory_1 } - } else if (subDetID == 1) { // MLOT - path += Form("%s%d_1/", getTRKLayerPattern(), layer); // TRKLayerx_1 - path += Form("%s%d_%d/", getTRKStavePattern(), layer, stave); // TRKStavex_y - if (mNumberOfHalfStaves[layer] == 2) { // staggered geometry - path += Form("%s%d_%d/", getTRKChipPattern(), layer, halfstave); // TRKChipx_0/1 - } else if (mNumberOfHalfStaves[layer] == 1) { // turbo geometry - path += Form("%s%d_1/", getTRKChipPattern(), layer); // TRKChipx_1 + } else if (subDetID == 1) { // MLOT + path += Form("%s%d_1/", getTRKLayerPattern(), layer); // TRKLayerx_1 + path += Form("%s%d_%d/", getTRKStavePattern(), layer, stave); // TRKStavex_y + if (mNumberOfHalfStaves[layer] == 2) { // staggered geometry + path += Form("%s%d_%d/", getTRKHalfStavePattern(), layer, halfstave); // TRKHalfStavex_y } - path += Form("%s%d_1/", getTRKSensorPattern(), layer); // TRKSensorx_1 + path += Form("%s%d_%d/", getTRKModulePattern(), layer, mod); // TRKModulx_y + path += Form("%s%d_%d/", getTRKChipPattern(), layer, chip); // TRKChipx_y + path += Form("%s%d_1/", getTRKSensorPattern(), layer); // TRKSensorx_1 } return path; } @@ -370,25 +492,30 @@ void GeometryTGeo::fillMatrixCache(int mask) //__________________________________________________________________________ #ifdef ENABLE_UPGRADES -const char* GeometryTGeo::composeSymNameLayer(int d, int lr) +const char* GeometryTGeo::composeSymNameLayer(int d, int layer) { - return Form("%s/%s%d", composeSymNameTRK(d), getTRKLayerPattern(), lr); + return Form("%s/%s%d", composeSymNameTRK(d), getTRKLayerPattern(), layer); } #endif -const char* GeometryTGeo::composeSymNameStave(int d, int lr) +const char* GeometryTGeo::composeSymNameStave(int d, int layer) { - return Form("%s/%s%d", composeSymNameLayer(d, lr), getTRKStavePattern(), lr); + return Form("%s/%s%d", composeSymNameLayer(d, layer), getTRKStavePattern(), layer); } -const char* GeometryTGeo::composeSymNameChip(int d, int lr) +const char* GeometryTGeo::composeSymNameModule(int d, int layer) { - return Form("%s/%s%d", composeSymNameStave(d, lr), getTRKChipPattern(), lr); + return Form("%s/%s%d", composeSymNameStave(d, layer), getTRKModulePattern(), layer); } -const char* GeometryTGeo::composeSymNameSensor(int d, int lr) +const char* GeometryTGeo::composeSymNameChip(int d, int layer) { - return Form("%s/%s%d", composeSymNameChip(d, lr), getTRKSensorPattern(), lr); + return Form("%s/%s%d", composeSymNameStave(d, layer), getTRKChipPattern(), layer); +} + +const char* GeometryTGeo::composeSymNameSensor(int d, int layer) +{ + return Form("%s/%s%d", composeSymNameChip(d, layer), getTRKSensorPattern(), layer); } //__________________________________________________________________________ @@ -845,15 +972,71 @@ int GeometryTGeo::extractNumberOfHalfStavesMLOT(int lay) const for (int j = 0; j < nNodes; j++) { auto nd = dynamic_cast(nodes->At(j)); /// layer node const char* name = nd->GetName(); - if (strstr(name, getTRKChipPattern()) != nullptr) { + if (strstr(name, getTRKHalfStavePattern()) != nullptr) { numberOfHalfStaves++; } } + + if (numberOfHalfStaves == 0) { + numberOfHalfStaves = 1; /// in case of turbo geometry, there is no half stave volume, but only stave volume + } return numberOfHalfStaves; } //__________________________________________________________________________ -void GeometryTGeo::PrintChipID(int index, int subDetID, int petalcase, int disk, int lay, int stave, int halfstave) const +int GeometryTGeo::extractNumberOfModulesMLOT(int lay) const +{ + int numberOfModules = 0; + + std::string staveName = Form("%s%d", (mNumberOfHalfStaves[lay] == 2 ? getTRKHalfStavePattern() : getTRKStavePattern()), lay); + TGeoVolume* staveV = gGeoManager->GetVolume(staveName.c_str()); + + if (staveV == nullptr) { + LOG(fatal) << getName() << " volume " << (mNumberOfHalfStaves[lay] == 2 ? getTRKHalfStavePattern() : getTRKStavePattern()) << " is not in the geometry"; + } + + // Loop on all staveV nodes, count Module volumes by checking names + TObjArray* nodes = staveV->GetNodes(); + int nNodes = nodes->GetEntriesFast(); + + for (int j = 0; j < nNodes; j++) { + auto nd = dynamic_cast(nodes->At(j)); /// stave node + const char* name = nd->GetName(); + if (strstr(name, getTRKModulePattern()) != nullptr) { + numberOfModules++; + } + } + return numberOfModules; +} + +//__________________________________________________________________________ +int GeometryTGeo::extractNumberOfChipsMLOT(int lay) const +{ + int numberOfChips = 0; + + std::string moduleName = Form("%s%d", getTRKModulePattern(), lay); + TGeoVolume* moduleV = gGeoManager->GetVolume(moduleName.c_str()); + + if (moduleV == nullptr) { + LOG(fatal) << getName() << " volume " << getTRKModulePattern() << " is not in the geometry"; + } + + // Loop on all moduleV nodes, count Chip volumes by checking names + TObjArray* nodes = moduleV->GetNodes(); + int nNodes = nodes->GetEntriesFast(); + + for (int j = 0; j < nNodes; j++) { + auto nd = dynamic_cast(nodes->At(j)); /// module node + const char* name = nd->GetName(); + if (strstr(name, getTRKChipPattern()) != nullptr) { + numberOfChips++; + } + } + return numberOfChips; +} + +//__________________________________________________________________________ +void GeometryTGeo::PrintChipID(int index, int subDetID, int petalcase, int disk, int lay, int stave, int halfstave, int mod, int chip) const { std::cout << "\nindex = " << index << std::endl; std::cout << "subDetID = " << subDetID << std::endl; @@ -863,6 +1046,8 @@ void GeometryTGeo::PrintChipID(int index, int subDetID, int petalcase, int disk, std::cout << "first chip index = " << getFirstChipIndex(lay, petalcase, subDetID) << std::endl; std::cout << "stave = " << stave << std::endl; std::cout << "halfstave = " << halfstave << std::endl; + std::cout << "module = " << mod << std::endl; + std::cout << "chip = " << chip << std::endl; } //__________________________________________________________________________ @@ -890,6 +1075,18 @@ void GeometryTGeo::Print(Option_t*) const mlot = (i < 4) ? "ML" : "OT"; LOGF(info, "Layer: %d, %s, %d staves, %d half staves per stave", i, mlot.c_str(), mNumberOfStaves[i], mNumberOfHalfStaves[i]); } + LOGF(info, "Number of modules per layer MLOT: "); + for (int i = 0; i < mNumberOfLayersMLOT; i++) { + LOGF(info, "%d", mNumberOfModules[i]); + } + LOGF(info, "Number of chips per module MLOT: "); + for (int i = 0; i < mNumberOfLayersMLOT; i++) { + LOGF(info, "%d", mNumberOfChips[i]); + } + LOGF(info, "Number of chips per layer MLOT: "); + for (int i = 0; i < mNumberOfLayersMLOT; i++) { + LOGF(info, "%d", mNumberOfChipsPerLayerMLOT[i]); + } LOGF(info, "Total number of chips: %d", getNumberOfChips()); std::cout << "mLastChipIndex = ["; diff --git a/Detectors/Upgrades/ALICE3/TRK/simulation/CMakeLists.txt b/Detectors/Upgrades/ALICE3/TRK/simulation/CMakeLists.txt index 0c3c35d49f722..10f117750d793 100644 --- a/Detectors/Upgrades/ALICE3/TRK/simulation/CMakeLists.txt +++ b/Detectors/Upgrades/ALICE3/TRK/simulation/CMakeLists.txt @@ -10,7 +10,8 @@ # or submit itself to any jurisdiction. o2_add_library(TRKSimulation - SOURCES src/TRKLayer.cxx + SOURCES src/Hit.cxx + src/TRKLayer.cxx src/ChipDigitsContainer.cxx src/ChipSimResponse.cxx src/Detector.cxx @@ -27,7 +28,8 @@ o2_add_library(TRKSimulation O2::SimulationDataFormat) o2_target_root_dictionary(TRKSimulation - HEADERS include/TRKSimulation/ChipDigitsContainer.h + HEADERS include/TRKSimulation/Hit.h + include/TRKSimulation/ChipDigitsContainer.h include/TRKSimulation/ChipSimResponse.h include/TRKSimulation/DigiParams.h include/TRKSimulation/Digitizer.h diff --git a/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/Detector.h b/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/Detector.h index 92cebd681176d..32bdc89109269 100644 --- a/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/Detector.h +++ b/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/Detector.h @@ -13,7 +13,7 @@ #define ALICEO2_TRK_DETECTOR_H #include "DetectorsBase/Detector.h" -#include "ITSMFTSimulation/Hit.h" +#include "TRKSimulation/Hit.h" #include "TRKSimulation/TRKLayer.h" #include "TRKSimulation/TRKServices.h" @@ -42,9 +42,9 @@ class Detector : public o2::base::DetImpl void ConstructGeometry() 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); + o2::trk::Hit* addHit(int trackID, unsigned short detID, const TVector3& startPos, const TVector3& endPos, + const TVector3& startMom, double startE, double endTime, double eLoss, + unsigned char startStatus, unsigned char endStatus); // Mandatory overrides void BeginPrimary() override { ; } @@ -57,8 +57,8 @@ class Detector : public o2::base::DetImpl void Register() override; void Reset() override; - // Custom memer functions - std::vector* getHits(int iColl) const + // Custom member functions + std::vector* getHits(int iColl) const { if (!iColl) { return mHits; @@ -81,14 +81,14 @@ class Detector : public o2::base::DetImpl // Transient data about track passing the sensor struct TrackData { - bool mHitStarted; // hit creation started - unsigned char mTrkStatusStart; // track status flag - TLorentzVector mPositionStart; // position at entrance - TLorentzVector mMomentumStart; // momentum - double mEnergyLoss; // energy loss - } mTrackData; //! transient data - GeometryTGeo* mGeometryTGeo; //! - std::vector* mHits; // ITSMFT ones for the moment + bool mHitStarted; // hit creation started + unsigned char mTrkStatusStart; // track status flag + TLorentzVector mPositionStart; // position at entrance + TLorentzVector mMomentumStart; // momentum + double mEnergyLoss; // energy loss + } mTrackData; //! transient data + GeometryTGeo* mGeometryTGeo; //! + std::vector* mHits; // ITSMFT ones for the moment std::vector mLayers; TRKServices mServices; // Houses the services of the TRK, but not the Iris tracker @@ -105,11 +105,11 @@ class Detector : public o2::base::DetImpl static constexpr Int_t sNumberVDPetalCases = 4; //! Number of VD petals int getNumberOfLayers() const { return mLayers.size(); } //! Number of TRK layers - void Print(FairVolume* vol, int volume, int subDetID, int layer, int stave, int halfstave, int chipID) const; + void Print(FairVolume* vol, int volume, int subDetID, int layer, int stave, int halfstave, int mod, int chip, int chipID) const; template friend class o2::base::DetImpl; - ClassDefOverride(Detector, 1); + ClassDefOverride(Detector, 2); }; } // namespace trk } // namespace o2 diff --git a/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/Hit.h b/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/Hit.h new file mode 100644 index 0000000000000..a178c30069f14 --- /dev/null +++ b/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/Hit.h @@ -0,0 +1,149 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file Hit.h +/// \brief Definition of the TRK Hit class + +#ifndef ALICEO2_TRK_POINT_H_ +#define ALICEO2_TRK_POINT_H_ + +#include "SimulationDataFormat/BaseHits.h" // for BasicXYZEHit +#include "Rtypes.h" // for Bool_t, Double_t, Int_t, Double32_t, etc +#include "TVector3.h" // for TVector3 +#include +#include "CommonUtils/ShmAllocator.h" + +namespace o2 +{ +namespace trk +{ + +class Hit : public o2::BasicXYZEHit +{ + + public: + enum HitStatus_t { + kTrackEntering = 0x1, + kTrackInside = 0x1 << 1, + kTrackExiting = 0x1 << 2, + kTrackOut = 0x1 << 3, + kTrackStopped = 0x1 << 4, + kTrackAlive = 0x1 << 5 + }; + + /// Default constructor + Hit() = default; + + /// Class Constructor + /// \param trackID Index of MCTrack + /// \param detID Detector ID + /// \param startPos Coordinates at entrance to active volume [cm] + /// \param pos Coordinates to active volume [cm] + /// \param mom Momentum of track at entrance [GeV] + /// \param endTime Time at entrance [ns] + /// \param time Time since event start [ns] + /// \param eLoss Energy deposit [GeV] + /// \param startStatus: status at entrance + /// \param endStatus: status at exit + inline Hit(int trackID, unsigned short detID, const TVector3& startPos, const TVector3& pos, const TVector3& mom, double startE, + double endTime, double eLoss, unsigned char statusStart, unsigned char status); + + // Entrance position getters + math_utils::Point3D GetPosStart() const { return mPosStart; } + Float_t GetStartX() const { return mPosStart.X(); } + Float_t GetStartY() const { return mPosStart.Y(); } + Float_t GetStartZ() const { return mPosStart.Z(); } + template + void GetStartPosition(F& x, F& y, F& z) const + { + x = GetStartX(); + y = GetStartY(); + z = GetStartZ(); + } + // momentum getters + math_utils::Vector3D GetMomentum() const { return mMomentum; } + math_utils::Vector3D& GetMomentum() { return mMomentum; } + Float_t GetPx() const { return mMomentum.X(); } + Float_t GetPy() const { return mMomentum.Y(); } + Float_t GetPz() const { return mMomentum.Z(); } + Float_t GetE() const { return mE; } + Float_t GetTotalEnergy() const { return GetE(); } + + UChar_t GetStatusEnd() const { return mTrackStatusEnd; } + UChar_t GetStatusStart() const { return mTrackStatusStart; } + + Bool_t IsEntering() const { return mTrackStatusEnd & kTrackEntering; } + Bool_t IsInside() const { return mTrackStatusEnd & kTrackInside; } + Bool_t IsExiting() const { return mTrackStatusEnd & kTrackExiting; } + Bool_t IsOut() const { return mTrackStatusEnd & kTrackOut; } + Bool_t IsStopped() const { return mTrackStatusEnd & kTrackStopped; } + Bool_t IsAlive() const { return mTrackStatusEnd & kTrackAlive; } + + Bool_t IsEnteringStart() const { return mTrackStatusStart & kTrackEntering; } + Bool_t IsInsideStart() const { return mTrackStatusStart & kTrackInside; } + Bool_t IsExitingStart() const { return mTrackStatusStart & kTrackExiting; } + Bool_t IsOutStart() const { return mTrackStatusStart & kTrackOut; } + Bool_t IsStoppedStart() const { return mTrackStatusStart & kTrackStopped; } + Bool_t IsAliveStart() const { return mTrackStatusStart & kTrackAlive; } + + // Entrance position setter + void SetPosStart(const math_utils::Point3D& p) { mPosStart = p; } + + /// Output to screen + void Print(const Option_t* opt) const; + friend std::ostream& operator<<(std::ostream& of, const Hit& point) + { + of << "-I- Hit: O2its point for track " << point.GetTrackID() << " in detector " << point.GetDetectorID() << std::endl; + /* + of << " Position (" << point.fX << ", " << point.fY << ", " << point.fZ << ") cm" << std::endl; + of << " Momentum (" << point.fPx << ", " << point.fPy << ", " << point.fPz << ") GeV" << std::endl; + of << " Time " << point.fTime << " ns, Length " << point.fLength << " cm, Energy loss " + << point.fELoss * 1.0e06 << " keV" << std::endl; + */ + return of; + } + + private: + math_utils::Vector3D mMomentum; ///< momentum at entrance + math_utils::Point3D mPosStart; ///< position at entrance (base mPos give position on exit) + Float_t mE; ///< total energy at entrance + UChar_t mTrackStatusEnd; ///< MC status flag at exit + UChar_t mTrackStatusStart; ///< MC status at starting point + + ClassDefNV(Hit, 1); +}; + +Hit::Hit(int trackID, unsigned short detID, const TVector3& startPos, const TVector3& endPos, const TVector3& startMom, + double startE, double endTime, double eLoss, unsigned char startStatus, unsigned char endStatus) + : BasicXYZEHit(endPos.X(), endPos.Y(), endPos.Z(), endTime, eLoss, trackID, detID), + mMomentum(startMom.Px(), startMom.Py(), startMom.Pz()), + mPosStart(startPos.X(), startPos.Y(), startPos.Z()), + mE(startE), + mTrackStatusEnd(endStatus), + mTrackStatusStart(startStatus) +{ +} + +} // namespace trk +} // namespace o2 + +#ifdef USESHM +namespace std +{ +template <> +class allocator : public o2::utils::ShmAllocator +{ +}; +} // namespace std + +#endif + +#endif diff --git a/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/TRKLayer.h b/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/TRKLayer.h index ba894f6d7a92b..0a7a45e87bfd8 100644 --- a/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/TRKLayer.h +++ b/Detectors/Upgrades/ALICE3/TRK/simulation/include/TRKSimulation/TRKLayer.h @@ -16,6 +16,7 @@ #include #include "TRKBase/TRKBaseParam.h" +#include "TRKBase/Specs.h" namespace o2 { @@ -25,23 +26,27 @@ class TRKLayer { public: TRKLayer() = default; - TRKLayer(int layerNumber, std::string layerName, float rInn, float rOut, float zLength, float layerX2X0); - TRKLayer(int layerNumber, std::string layerName, float rInn, float zLength, float thick); + TRKLayer(int layerNumber, std::string layerName, float rInn, float rOut, int numberOfModules, float layerX2X0); + TRKLayer(int layerNumber, std::string layerName, float rInn, int numberOfModules, float thick); ~TRKLayer() = default; void setLayout(eLayout layout) { mLayout = layout; }; auto getInnerRadius() const { return mInnerRadius; } auto getOuterRadius() const { return mOuterRadius; } - auto getZ() const { return mZ; } + auto getZ() const { return constants::moduleMLOT::length * mNumberOfModules; } auto getx2X0() const { return mX2X0; } auto getChipThickness() const { return mChipThickness; } auto getNumber() const { return mLayerNumber; } auto getName() const { return mLayerName; } - TGeoVolume* createSensor(std::string type, double width = -1); - TGeoVolume* createChip(std::string type, double width = -1); - TGeoVolume* createStave(std::string type, double width = -1); + TGeoVolume* createSensor(std::string type); + TGeoVolume* createDeadzone(std::string type); + TGeoVolume* createMetalStack(std::string type); + TGeoVolume* createChip(std::string type); + TGeoVolume* createModule(std::string type); + TGeoVolume* createStave(std::string type); + TGeoVolume* createHalfStave(std::string type); void createLayer(TGeoVolume* motherVolume); private: @@ -49,16 +54,20 @@ class TRKLayer static constexpr float mLogicalVolumeThickness = 1; int mLayerNumber; + eLayout mLayout; std::string mLayerName; float mInnerRadius; float mOuterRadius; - float mZ; + int mNumberOfModules; float mX2X0; + float mChipWidth; + float mChipLength; float mChipThickness; - float mModuleWidth; // u.m. = cm - eLayout mLayout; + float mDeadzoneWidth; + float mSensorThickness; + int mHalfNumberOfChips; - ClassDef(TRKLayer, 1); + ClassDef(TRKLayer, 2); }; } // namespace trk diff --git a/Detectors/Upgrades/ALICE3/TRK/simulation/src/Detector.cxx b/Detectors/Upgrades/ALICE3/TRK/simulation/src/Detector.cxx index a4d99ccf9f79f..4d7e560d50dc2 100644 --- a/Detectors/Upgrades/ALICE3/TRK/simulation/src/Detector.cxx +++ b/Detectors/Upgrades/ALICE3/TRK/simulation/src/Detector.cxx @@ -16,13 +16,16 @@ #include #include "DetectorsBase/Stack.h" -#include "ITSMFTSimulation/Hit.h" +#include "TRKSimulation/Hit.h" #include "TRKSimulation/Detector.h" #include "TRKBase/TRKBaseParam.h" #include "TRKSimulation/VDGeometryBuilder.h" #include "TRKSimulation/VDSensorRegistry.h" -using o2::itsmft::Hit; +#include +#include + +using o2::trk::Hit; namespace o2 { @@ -37,14 +40,14 @@ float getDetLengthFromEta(const float eta, const float radius) Detector::Detector() : o2::base::DetImpl("TRK", true), mTrackData(), - mHits(o2::utils::createSimVector()) + mHits(o2::utils::createSimVector()) { } Detector::Detector(bool active) : o2::base::DetImpl("TRK", true), mTrackData(), - mHits(o2::utils::createSimVector()) + mHits(o2::utils::createSimVector()) { auto& trkPars = TRKBaseParam::Instance(); @@ -83,40 +86,58 @@ void Detector::configDefault() mLayers.clear(); LOGP(warning, "Loading Scoping Document configuration for ALICE3 TRK"); - // mLayers.emplace_back(0, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(0)}, 0.5f, 50.f, 100.e-4); - // mLayers.emplace_back(1, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(1)}, 1.2f, 50.f, 100.e-4); - // mLayers.emplace_back(2, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(2)}, 2.5f, 50.f, 100.e-4); - mLayers.emplace_back(0, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(0)}, 3.78f, 124.f, 100.e-3); - mLayers.emplace_back(1, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(1)}, 7.f, 124.f, 100.e-3); - mLayers.emplace_back(2, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(2)}, 12.f, 124.f, 100.e-3); - mLayers.emplace_back(3, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(3)}, 20.f, 124.f, 100.e-3); - mLayers.emplace_back(4, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(4)}, 30.f, 124.f, 100.e-3); - mLayers.emplace_back(5, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(5)}, 45.f, 258.f, 100.e-3); - mLayers.emplace_back(6, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(6)}, 60.f, 258.f, 100.e-3); - mLayers.emplace_back(7, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(7)}, 80.f, 258.f, 100.e-3); + mLayers.emplace_back(0, GeometryTGeo::getTRKLayerPattern() + std::to_string(0), 3.78f, 10, 100.e-3); + mLayers.emplace_back(1, GeometryTGeo::getTRKLayerPattern() + std::to_string(1), 7.f, 10, 100.e-3); + mLayers.emplace_back(2, GeometryTGeo::getTRKLayerPattern() + std::to_string(2), 12.f, 10, 100.e-3); + mLayers.emplace_back(3, GeometryTGeo::getTRKLayerPattern() + std::to_string(3), 20.f, 10, 100.e-3); + mLayers.emplace_back(4, GeometryTGeo::getTRKLayerPattern() + std::to_string(4), 30.f, 10, 100.e-3); + mLayers.emplace_back(5, GeometryTGeo::getTRKLayerPattern() + std::to_string(5), 45.f, 20, 100.e-3); + mLayers.emplace_back(6, GeometryTGeo::getTRKLayerPattern() + std::to_string(6), 60.f, 20, 100.e-3); + mLayers.emplace_back(7, GeometryTGeo::getTRKLayerPattern() + std::to_string(7), 80.f, 20, 100.e-3); } void Detector::buildTRKMiddleOuterLayers() { - // Build the TRK detector according to changes proposed during - // https://indico.cern.ch/event/1407704/ - // to adhere to the changes that were presented at the ALICE 3 Upgrade days in March 2024 - // L3 -> 7 cm, L4 -> 9 cm + auto& trkPars = TRKBaseParam::Instance(); mLayers.clear(); - LOGP(warning, "Loading \"After Upgrade Days March 2024\" configuration for ALICE3 TRK"); - mLayers.emplace_back(0, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(0)}, 7.f, 124.f, 100.e-3); - LOGP(info, "TRKLayer created. Name: {}", std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(0)}); - mLayers.emplace_back(1, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(1)}, 9.f, 124.f, 100.e-3); - mLayers.emplace_back(2, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(2)}, 12.f, 124.f, 100.e-3); - mLayers.emplace_back(3, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(3)}, 20.f, 124.f, 100.e-3); - mLayers.emplace_back(4, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(4)}, 30.f, 124.f, 100.e-3); - mLayers.emplace_back(5, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(5)}, 45.f, 258.f, 100.e-3); - mLayers.emplace_back(6, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(6)}, 60.f, 258.f, 100.e-3); - mLayers.emplace_back(7, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(7)}, 80.f, 258.f, 100.e-3); - - auto& trkPars = TRKBaseParam::Instance(); + switch (trkPars.overallGeom) { + case kDefaultRadii: + // Build the TRK detector according to changes proposed during + // https://indico.cern.ch/event/1407704/ + // to adhere to the changes that were presented at the ALICE 3 Upgrade days in March 2024 + // L3 -> 7 cm, L4 -> 9 cm, L5 -> 12 cm, L6 -> 20 cm + + LOGP(warning, "Loading \"After Upgrade Days March 2024\" configuration for ALICE3 TRK"); + LOGP(warning, "Building TRK with new vacuum vessel and L3 at 7 cm, L4 at 9 cm, L5 at 12 cm, L6 at 20 cm"); + mLayers.emplace_back(0, GeometryTGeo::getTRKLayerPattern() + std::to_string(0), 7.f, 10, 100.e-3); + LOGP(info, "TRKLayer created. Name: {}", GeometryTGeo::getTRKLayerPattern() + std::to_string(0)); + mLayers.emplace_back(1, GeometryTGeo::getTRKLayerPattern() + std::to_string(1), 9.f, 10, 100.e-3); + mLayers.emplace_back(2, GeometryTGeo::getTRKLayerPattern() + std::to_string(2), 12.f, 10, 100.e-3); + mLayers.emplace_back(3, GeometryTGeo::getTRKLayerPattern() + std::to_string(3), 20.f, 10, 100.e-3); + mLayers.emplace_back(4, GeometryTGeo::getTRKLayerPattern() + std::to_string(4), 30.f, 10, 100.e-3); + mLayers.emplace_back(5, GeometryTGeo::getTRKLayerPattern() + std::to_string(5), 45.f, 20, 100.e-3); + mLayers.emplace_back(6, GeometryTGeo::getTRKLayerPattern() + std::to_string(6), 60.f, 20, 100.e-3); + mLayers.emplace_back(7, GeometryTGeo::getTRKLayerPattern() + std::to_string(7), 80.f, 20, 100.e-3); + break; + case kModRadii: + LOGP(warning, "Loading \"Alternative\" configuration for ALICE3 TRK"); + LOGP(warning, "Building TRK with new vacuum vessel and L3 at 7 cm, L4 at 11 cm, L5 at 15 cm, L6 at 19 cm"); + mLayers.emplace_back(0, GeometryTGeo::getTRKLayerPattern() + std::to_string(0), 7.f, 10, 100.e-3); + LOGP(info, "TRKLayer created. Name: {}", GeometryTGeo::getTRKLayerPattern() + std::to_string(0)); + mLayers.emplace_back(1, GeometryTGeo::getTRKLayerPattern() + std::to_string(1), 11.f, 10, 100.e-3); + mLayers.emplace_back(2, GeometryTGeo::getTRKLayerPattern() + std::to_string(2), 15.f, 10, 100.e-3); + mLayers.emplace_back(3, GeometryTGeo::getTRKLayerPattern() + std::to_string(3), 19.f, 10, 100.e-3); + mLayers.emplace_back(4, GeometryTGeo::getTRKLayerPattern() + std::to_string(4), 30.f, 10, 100.e-3); + mLayers.emplace_back(5, GeometryTGeo::getTRKLayerPattern() + std::to_string(5), 45.f, 20, 100.e-3); + mLayers.emplace_back(6, GeometryTGeo::getTRKLayerPattern() + std::to_string(6), 60.f, 20, 100.e-3); + mLayers.emplace_back(7, GeometryTGeo::getTRKLayerPattern() + std::to_string(7), 80.f, 20, 100.e-3); + break; + default: + LOGP(fatal, "Unknown option {} for buildTRKMiddleOuterLayers", static_cast(trkPars.overallGeom)); + break; + } // Middle layers mLayers[0].setLayout(trkPars.layoutML); @@ -157,7 +178,7 @@ void Detector::configFromFile(std::string fileName) while (getline(ss, substr, '\t')) { tmpBuff.push_back(std::stof(substr)); } - mLayers.emplace_back(layerCount, std::string{GeometryTGeo::getTRKLayerPattern() + std::to_string(layerCount)}, tmpBuff[0], tmpBuff[1], tmpBuff[2]); + mLayers.emplace_back(layerCount, GeometryTGeo::getTRKLayerPattern() + std::to_string(layerCount), tmpBuff[0], tmpBuff[1], tmpBuff[2]); ++layerCount; } } @@ -364,7 +385,6 @@ bool Detector::ProcessHits(FairVolume* vol) int subDetID = -1; int layer = -1; int volume = 0; - int stave = -1; int volID = vol->getMCid(); bool notSens = false; @@ -440,15 +460,23 @@ bool Detector::ProcessHits(FairVolume* vol) TLorentzVector positionStop; fMC->TrackPosition(positionStop); // Retrieve the indices with the volume path - int stave(0), halfstave(0); + int stave(0), halfstave(0), mod(0), chip(0); if (subDetID == 1) { - fMC->CurrentVolOffID(1, halfstave); - fMC->CurrentVolOffID(2, stave); + fMC->CurrentVolOffID(1, chip); + fMC->CurrentVolOffID(2, mod); + if (mGeometryTGeo->getNumberOfHalfStaves(layer) == 2) { + fMC->CurrentVolOffID(3, halfstave); + fMC->CurrentVolOffID(4, stave); + } else if (mGeometryTGeo->getNumberOfHalfStaves(layer) == 1) { + fMC->CurrentVolOffID(3, stave); + } else { + LOGP(fatal, "Wrong number of halfstaves for layer {}", layer); + } } /// if VD, for the moment the volume is the "chipID" so no need to retrieve other elments - int chipID = mGeometryTGeo->getChipIndex(subDetID, volume, layer, stave, halfstave); + unsigned short chipID = mGeometryTGeo->getChipIndex(subDetID, volume, layer, stave, halfstave, mod, chip); - Print(vol, volume, subDetID, layer, stave, halfstave, chipID); + Print(vol, volume, subDetID, layer, stave, halfstave, mod, chip, chipID); mGeometryTGeo->Print(); @@ -465,25 +493,27 @@ bool Detector::ProcessHits(FairVolume* vol) return true; } -o2::itsmft::Hit* Detector::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) +o2::trk::Hit* Detector::addHit(int trackID, unsigned short detID, const TVector3& startPos, const TVector3& endPos, + const TVector3& startMom, double startE, double endTime, double eLoss, unsigned char startStatus, + unsigned char endStatus) { mHits->emplace_back(trackID, detID, startPos, endPos, startMom, startE, endTime, eLoss, startStatus, endStatus); return &(mHits->back()); } -void Detector::Print(FairVolume* vol, int volume, int subDetID, int layer, int stave, int halfstave, int chipID) const +void Detector::Print(FairVolume* vol, int volume, int subDetID, int layer, int stave, int halfstave, int mod, int chip, int chipID) const { int currentVol(0); LOG(info) << "Current volume name: " << fMC->CurrentVolName() << " and ID " << fMC->CurrentVolID(currentVol); LOG(info) << "volume: " << volume << "/" << mNumberOfVolumes - 1; + LOG(info) << "off volume name 1 " << fMC->CurrentVolOffName(1) << " chip: " << chip; + LOG(info) << "off volume name 2 " << fMC->CurrentVolOffName(2) << " module: " << mod; if (subDetID == 1 && mGeometryTGeo->getNumberOfHalfStaves(layer) == 2) { // staggered geometry - LOG(info) << "off volume name 1 " << fMC->CurrentVolOffName(1) << " halfstave: " << halfstave; - LOG(info) << "off volume name 2 " << fMC->CurrentVolOffName(2) << " stave: " << stave; + LOG(info) << "off volume name 3 " << fMC->CurrentVolOffName(3) << " halfstave: " << halfstave; + LOG(info) << "off volume name 4 " << fMC->CurrentVolOffName(4) << " stave: " << stave; LOG(info) << "SubDetector ID: " << subDetID << " Layer: " << layer << " staveinLayer: " << stave << " Chip ID: " << chipID; } else if (subDetID == 1 && mGeometryTGeo->getNumberOfHalfStaves(layer) == 1) { // turbo geometry - LOG(info) << "off volume name 2 " << fMC->CurrentVolOffName(2) << " stave: " << stave; + LOG(info) << "off volume name 3 " << fMC->CurrentVolOffName(3) << " stave: " << stave; LOG(info) << "SubDetector ID: " << subDetID << " Layer: " << layer << " staveinLayer: " << stave << " Chip ID: " << chipID; } else { LOG(info) << "SubDetector ID: " << subDetID << " Chip ID: " << chipID; diff --git a/Detectors/Upgrades/ALICE3/TRK/simulation/src/Hit.cxx b/Detectors/Upgrades/ALICE3/TRK/simulation/src/Hit.cxx new file mode 100644 index 0000000000000..fe496bc59692f --- /dev/null +++ b/Detectors/Upgrades/ALICE3/TRK/simulation/src/Hit.cxx @@ -0,0 +1,34 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file Hit.cxx +/// \brief Implementation of the Hit class + +#include "TRKSimulation/Hit.h" + +#include +#include + +ClassImp(o2::trk::Hit); + +using std::cout; +using std::endl; +using namespace o2::trk; +using namespace o2; //::base; + +void Hit::Print(const Option_t* opt) const +{ + printf( + "Det: %5d Track: %6d E.loss: %.3e P: %+.3e %+.3e %+.3e\n" + "PosIn: %+.3e %+.3e %+.3e PosOut: %+.3e %+.3e %+.3e\n", + GetDetectorID(), GetTrackID(), GetEnergyLoss(), GetPx(), GetPy(), GetPz(), + GetStartX(), GetStartY(), GetStartZ(), GetX(), GetY(), GetZ()); +} diff --git a/Detectors/Upgrades/ALICE3/TRK/simulation/src/TRKLayer.cxx b/Detectors/Upgrades/ALICE3/TRK/simulation/src/TRKLayer.cxx index a95418afbba25..223c8b5c477a1 100644 --- a/Detectors/Upgrades/ALICE3/TRK/simulation/src/TRKLayer.cxx +++ b/Detectors/Upgrades/ALICE3/TRK/simulation/src/TRKLayer.cxx @@ -11,6 +11,7 @@ #include "TRKSimulation/TRKLayer.h" #include "TRKBase/GeometryTGeo.h" +#include "TRKBase/Specs.h" #include "Framework/Logger.h" @@ -24,37 +25,34 @@ namespace o2 { namespace trk { -TRKLayer::TRKLayer(int layerNumber, std::string layerName, float rInn, float rOut, float zLength, float layerX2X0) - : mLayerNumber(layerNumber), mLayerName(layerName), mInnerRadius(rInn), mOuterRadius(rOut), mZ(zLength), mX2X0(layerX2X0), mModuleWidth(4.54), mLayout(kCylinder) +TRKLayer::TRKLayer(int layerNumber, std::string layerName, float rInn, float rOut, int numberOfModules, float layerX2X0) + : mLayerNumber(layerNumber), mLayout(kCylinder), mLayerName(layerName), mInnerRadius(rInn), mOuterRadius(rOut), mNumberOfModules(numberOfModules), mX2X0(layerX2X0), mChipWidth(constants::moduleMLOT::chip::width), mChipLength(constants::moduleMLOT::chip::length), mDeadzoneWidth(constants::moduleMLOT::chip::passiveEdgeReadOut), mSensorThickness(constants::moduleMLOT::silicon::thickness), mHalfNumberOfChips(4) { float Si_X0 = 9.5f; mChipThickness = mX2X0 * Si_X0; - LOGP(info, "Creating layer: id: {} rInner: {} rOuter: {} zLength: {} x2X0: {}", mLayerNumber, mInnerRadius, mOuterRadius, mZ, mX2X0); + LOGP(info, "Creating layer: id: {} rInner: {} rOuter: {} zLength: {} x2X0: {}", mLayerNumber, mInnerRadius, mOuterRadius, getZ(), mX2X0); } -TRKLayer::TRKLayer(int layerNumber, std::string layerName, float rInn, float zLength, float thick) - : mLayerNumber(layerNumber), mLayerName(layerName), mInnerRadius(rInn), mZ(zLength), mChipThickness(thick), mModuleWidth(4.54), mLayout(kCylinder) +TRKLayer::TRKLayer(int layerNumber, std::string layerName, float rInn, int numberOfModules, float thick) + : mLayerNumber(layerNumber), mLayout(kCylinder), mLayerName(layerName), mInnerRadius(rInn), mNumberOfModules(numberOfModules), mChipThickness(thick), mChipWidth(constants::moduleMLOT::chip::width), mChipLength(constants::moduleMLOT::chip::length), mDeadzoneWidth(constants::moduleMLOT::chip::passiveEdgeReadOut), mSensorThickness(constants::moduleMLOT::silicon::thickness), mHalfNumberOfChips(4) { float Si_X0 = 9.5f; mOuterRadius = rInn + thick; mX2X0 = mChipThickness / Si_X0; - LOGP(info, "Creating layer: id: {} rInner: {} rOuter: {} zLength: {} x2X0: {}", mLayerNumber, mInnerRadius, mOuterRadius, mZ, mX2X0); + LOGP(info, "Creating layer: id: {} rInner: {} rOuter: {} zLength: {} x2X0: {}", mLayerNumber, mInnerRadius, mOuterRadius, getZ(), mX2X0); } -TGeoVolume* TRKLayer::createSensor(std::string type, double width) +TGeoVolume* TRKLayer::createSensor(std::string type) { TGeoMedium* medSi = gGeoManager->GetMedium("TRK_SILICON$"); - std::string sensName = Form("%s%d", GeometryTGeo::getTRKSensorPattern(), this->mLayerNumber); + std::string sensName = GeometryTGeo::getTRKSensorPattern() + std::to_string(mLayerNumber); TGeoShape* sensor; if (type == "cylinder") { - sensor = new TGeoTube(mInnerRadius, mInnerRadius + mChipThickness, mZ / 2); + sensor = new TGeoTube(mInnerRadius, mInnerRadius + mSensorThickness, mChipLength / 2); // TO BE CHECKED !!! } else if (type == "flat") { - if (width < 0) { - LOGP(fatal, "Attempting to create sensor with invalid width"); - } - sensor = new TGeoBBox(width / 2, mChipThickness / 2, mZ / 2); + sensor = new TGeoBBox((mChipWidth - mDeadzoneWidth) / 2, mSensorThickness / 2, mChipLength / 2); // TO BE CHECKED !!! } else { LOGP(fatal, "Sensor of type '{}' is not implemented", type); } @@ -65,75 +63,262 @@ TGeoVolume* TRKLayer::createSensor(std::string type, double width) return sensVol; }; -TGeoVolume* TRKLayer::createChip(std::string type, double width) +TGeoVolume* TRKLayer::createDeadzone(std::string type) { TGeoMedium* medSi = gGeoManager->GetMedium("TRK_SILICON$"); - std::string chipName = o2::trk::GeometryTGeo::getTRKChipPattern() + std::to_string(mLayerNumber); + std::string deadName = GeometryTGeo::getTRKDeadzonePattern() + std::to_string(mLayerNumber); + + TGeoShape* deadzone; + + if (type == "cylinder") { + deadzone = new TGeoTube(mInnerRadius, mInnerRadius + mSensorThickness, mChipLength / 2); // TO BE CHECKED !!! + } else if (type == "flat") { + deadzone = new TGeoBBox(mDeadzoneWidth / 2, mSensorThickness / 2, mChipLength / 2); // TO BE CHECKED !!! + } else { + LOGP(fatal, "Deadzone of type '{}' is not implemented", type); + } + + TGeoVolume* deadVol = new TGeoVolume(deadName.c_str(), deadzone, medSi); + deadVol->SetLineColor(kGray); + + return deadVol; +}; + +TGeoVolume* TRKLayer::createMetalStack(std::string type) +{ + TGeoMedium* medSi = gGeoManager->GetMedium("TRK_SILICON$"); + std::string metalName = GeometryTGeo::getTRKMetalStackPattern() + std::to_string(mLayerNumber); + + TGeoShape* metalStack; + + if (type == "cylinder") { + metalStack = new TGeoTube(mInnerRadius + mSensorThickness, mInnerRadius + mChipThickness, mChipLength / 2); // TO BE CHECKED !!! + } else if (type == "flat") { + metalStack = new TGeoBBox(mChipWidth / 2, mChipThickness - mSensorThickness / 2, mChipLength / 2); // TO BE CHECKED !!! + } else { + LOGP(fatal, "Metal stack of type '{}' is not implemented", type); + } + + TGeoVolume* metalVol = new TGeoVolume(metalName.c_str(), metalStack, medSi); + metalVol->SetLineColor(kGray); + + return metalVol; +}; + +TGeoVolume* TRKLayer::createChip(std::string type) +{ + TGeoMedium* medSi = gGeoManager->GetMedium("TRK_SILICON$"); + std::string chipName = GeometryTGeo::getTRKChipPattern() + std::to_string(mLayerNumber); TGeoShape* chip; + TGeoVolume* chipVol; + TGeoVolume* sensVol; + TGeoVolume* deadVol; + TGeoVolume* metalVol; if (type == "cylinder") { - chip = new TGeoTube(mInnerRadius, mInnerRadius + mChipThickness, mZ / 2); + chip = new TGeoTube(mInnerRadius, mInnerRadius + mChipThickness, mChipLength / 2); + chipVol = new TGeoVolume(chipName.c_str(), chip, medSi); + sensVol = createSensor("cylinder"); + metalVol = createMetalStack("cylinder"); + + TGeoCombiTrans* transSens = new TGeoCombiTrans(); + transSens->SetTranslation(0, -(mChipThickness - mSensorThickness) / 2, 0); // TO BE CHECKED !!! + LOGP(info, "Inserting {} in {} ", sensVol->GetName(), chipVol->GetName()); + chipVol->AddNode(sensVol, 1, transSens); + + TGeoCombiTrans* transMetal = new TGeoCombiTrans(); + transMetal->SetTranslation(0, mSensorThickness / 2, 0); // TO BE CHECKED !!! + LOGP(info, "Inserting {} in {} ", metalVol->GetName(), chipVol->GetName()); + chipVol->AddNode(metalVol, 1, transMetal); + + // deadVol = createDeadzone("cylinder"); } else if (type == "flat") { - if (width < 0) { - LOGP(fatal, "Attempting to create chip with invalid width"); - } - chip = new TGeoBBox(width / 2, mChipThickness / 2, mZ / 2); - sensVol = createSensor("flat", width); + chip = new TGeoBBox(mChipWidth / 2, mChipThickness / 2, mChipLength / 2); // TO BE CHECKED !!! + chipVol = new TGeoVolume(chipName.c_str(), chip, medSi); + + sensVol = createSensor("flat"); + deadVol = createDeadzone("flat"); + metalVol = createMetalStack("flat"); + + TGeoCombiTrans* transSens = new TGeoCombiTrans(); + transSens->SetTranslation(-mDeadzoneWidth / 2, -(mChipThickness - mSensorThickness) / 2, 0); // TO BE CHECKED !!! + LOGP(info, "Inserting {} in {} ", sensVol->GetName(), chipVol->GetName()); + chipVol->AddNode(sensVol, 1, transSens); + + TGeoCombiTrans* transDead = new TGeoCombiTrans(); + transDead->SetTranslation((mChipWidth - mDeadzoneWidth) / 2, -(mChipThickness - mSensorThickness) / 2, 0); // TO BE CHECKED !!! + LOGP(info, "Inserting {} in {} ", deadVol->GetName(), chipVol->GetName()); + chipVol->AddNode(deadVol, 1, transDead); + + TGeoCombiTrans* transMetal = new TGeoCombiTrans(); + transMetal->SetTranslation(0, mSensorThickness / 2, 0); // TO BE CHECKED !!! + LOGP(info, "Inserting {} in {} ", metalVol->GetName(), chipVol->GetName()); + chipVol->AddNode(metalVol, 1, transMetal); } else { LOGP(fatal, "Sensor of type '{}' is not implemented", type); } - TGeoVolume* chipVol = new TGeoVolume(chipName.c_str(), chip, medSi); - LOGP(info, "Inserting {} in {} ", sensVol->GetName(), chipVol->GetName()); - chipVol->AddNode(sensVol, 1, nullptr); chipVol->SetLineColor(kYellow); return chipVol; } -TGeoVolume* TRKLayer::createStave(std::string type, double width) +TGeoVolume* TRKLayer::createModule(std::string type) { TGeoMedium* medAir = gGeoManager->GetMedium("TRK_AIR$"); - std::string staveName = o2::trk::GeometryTGeo::getTRKStavePattern() + std::to_string(mLayerNumber); + std::string moduleName = GeometryTGeo::getTRKModulePattern() + std::to_string(mLayerNumber); + + TGeoShape* module; + TGeoVolume* moduleVol; + + if (type == "cylinder") { + module = new TGeoTube(mInnerRadius, mInnerRadius + mChipThickness, mChipLength / 2); + moduleVol = new TGeoVolume(moduleName.c_str(), module, medAir); + + TGeoVolume* chipVol = createChip("cylinder"); + LOGP(info, "Inserting {} in {} ", chipVol->GetName(), moduleVol->GetName()); + moduleVol->AddNode(chipVol, 1, nullptr); + } else if (type == "flat") { + double moduleWidth = constants::moduleMLOT::width; + double moduleLength = constants::moduleMLOT::length; + + module = new TGeoBBox(moduleWidth / 2, mChipThickness / 2, moduleLength / 2); // TO BE CHECKED !!! + moduleVol = new TGeoVolume(moduleName.c_str(), module, medAir); + + for (int iChip = 0; iChip < mHalfNumberOfChips; iChip++) { + TGeoVolume* chipVolLeft = createChip("flat"); + TGeoVolume* chipVolRight = createChip("flat"); + + // Put the chips in the correct position + double xLeft = -moduleWidth / 2 + constants::moduleMLOT::gaps::outerEdgeLongSide + constants::moduleMLOT::chip::width / 2; + double zLeft = -moduleLength / 2 + constants::moduleMLOT::gaps::outerEdgeShortSide + iChip * (constants::moduleMLOT::chip::length + constants::moduleMLOT::gaps::interChips) + constants::moduleMLOT::chip::length / 2; + + TGeoCombiTrans* transLeft = new TGeoCombiTrans(); + transLeft->SetTranslation(xLeft, 0, zLeft); // TO BE CHECKED !!! + LOGP(info, "Inserting {} in {} ", chipVolLeft->GetName(), moduleVol->GetName()); + moduleVol->AddNode(chipVolLeft, iChip * 2, transLeft); + + double xRight = +moduleWidth / 2 - constants::moduleMLOT::gaps::outerEdgeLongSide - constants::moduleMLOT::chip::width / 2; + double zRight = -moduleLength / 2 + constants::moduleMLOT::gaps::outerEdgeShortSide + iChip * (constants::moduleMLOT::chip::length + constants::moduleMLOT::gaps::interChips) + constants::moduleMLOT::chip::length / 2; + + TGeoCombiTrans* transRight = new TGeoCombiTrans(); + transRight->SetTranslation(xRight, 0, zRight); // TO BE CHECKED !!! + TGeoRotation* rot = new TGeoRotation(); + rot->RotateY(180); + transRight->SetRotation(rot); + LOGP(info, "Inserting {} in {} ", chipVolRight->GetName(), moduleVol->GetName()); + moduleVol->AddNode(chipVolRight, iChip * 2 + 1, transRight); + } + } else { + LOGP(fatal, "Chip of type '{}' is not implemented", type); + } + + moduleVol->SetLineColor(kYellow); + + return moduleVol; +} + +TGeoVolume* TRKLayer::createHalfStave(std::string type) +{ + TGeoMedium* medAir = gGeoManager->GetMedium("TRK_AIR$"); + std::string halfStaveName = GeometryTGeo::getTRKHalfStavePattern() + std::to_string(mLayerNumber); + + TGeoShape* halfStave; + TGeoVolume* halfStaveVol; + + if (type == "cylinder") { + halfStave = new TGeoTube(mInnerRadius, mInnerRadius + mChipThickness, mChipLength / 2); + halfStaveVol = new TGeoVolume(halfStaveName.c_str(), halfStave, medAir); + + TGeoVolume* moduleVol = createModule("cylinder"); + LOGP(info, "Inserting {} in {} ", moduleVol->GetName(), halfStaveVol->GetName()); + halfStaveVol->AddNode(moduleVol, 1, nullptr); + } else if (type == "flat") { + double moduleLength = constants::moduleMLOT::length; + double halfStaveWidth = constants::OT::halfstave::width; + double halfStaveLength = constants::moduleMLOT::length * mNumberOfModules; + + halfStave = new TGeoBBox(halfStaveWidth / 2, mChipThickness / 2, halfStaveLength / 2); + halfStaveVol = new TGeoVolume(halfStaveName.c_str(), halfStave, medAir); + + for (int iModule = 0; iModule < mNumberOfModules; iModule++) { + TGeoVolume* moduleVol = createModule("flat"); + + // Put the modules in the correct position + double zPos = -0.5 * mNumberOfModules * moduleLength + (iModule + 0.5) * moduleLength; + + TGeoCombiTrans* trans = new TGeoCombiTrans(); + trans->SetTranslation(0, 0, zPos); // TO BE CHECKED !!! + + LOGP(info, "Inserting {} in {} ", moduleVol->GetName(), halfStaveVol->GetName()); + halfStaveVol->AddNode(moduleVol, iModule, trans); + } + } + return halfStaveVol; +} + +TGeoVolume* TRKLayer::createStave(std::string type) +{ + TGeoMedium* medAir = gGeoManager->GetMedium("TRK_AIR$"); + std::string staveName = GeometryTGeo::getTRKStavePattern() + std::to_string(mLayerNumber); TGeoShape* stave; TGeoVolume* staveVol; - TGeoVolume* chipVol; if (type == "cylinder") { - stave = new TGeoTube(mInnerRadius, mInnerRadius + mChipThickness, mZ / 2); - chipVol = createChip("cylinder"); + stave = new TGeoTube(mInnerRadius, mInnerRadius + mChipThickness, mChipLength / 2); staveVol = new TGeoVolume(staveName.c_str(), stave, medAir); - LOGP(info, "Inserting {} in {} ", chipVol->GetName(), staveVol->GetName()); - staveVol->AddNode(chipVol, 1, nullptr); + + TGeoVolume* moduleVol = createModule("cylinder"); + LOGP(info, "Inserting {} in {} ", moduleVol->GetName(), staveVol->GetName()); + staveVol->AddNode(moduleVol, 1, nullptr); } else if (type == "flat") { - if (width < 0) { - LOGP(fatal, "Attempting to create stave with invalid width"); - } - stave = new TGeoBBox(width / 2, mChipThickness / 2, mZ / 2); - chipVol = createChip("flat", width); + double moduleLength = constants::moduleMLOT::length; + double staveWidth = constants::ML::width; + double staveLength = constants::moduleMLOT::length * mNumberOfModules; + + stave = new TGeoBBox(staveWidth / 2, mChipThickness / 2, staveLength / 2); staveVol = new TGeoVolume(staveName.c_str(), stave, medAir); - LOGP(info, "Inserting {} in {} ", chipVol->GetName(), staveVol->GetName()); - staveVol->AddNode(chipVol, 1, nullptr); + + for (int iModule = 0; iModule < mNumberOfModules; iModule++) { + TGeoVolume* moduleVol = createModule("flat"); + + // Put the modules in the correct position + double zPos = -0.5 * mNumberOfModules * moduleLength + (iModule + 0.5) * moduleLength; + + TGeoCombiTrans* trans = new TGeoCombiTrans(); + trans->SetTranslation(0, 0, zPos); // TO BE CHECKED !!! + + LOGP(info, "Inserting {} in {} ", moduleVol->GetName(), staveVol->GetName()); + staveVol->AddNode(moduleVol, iModule, trans); + } } else if (type == "staggered") { - double width = mModuleWidth * 2; // Each stave has two modules (based on the LOI design) - stave = new TGeoBBox(width / 2, mLogicalVolumeThickness / 2, mZ / 2); - TGeoVolume* chipVolLeft = createChip("flat", mModuleWidth); - TGeoVolume* chipVolRight = createChip("flat", mModuleWidth); + /*double moduleWidth = constants::moduleMLOT::width; + double moduleLength = constants::moduleMLOT::length;*/ + + double halfstaveWidth = constants::ML::width; + double staveWidth = constants::OT::width; // Each stave has two modules (based on the LOI design) + double staveLength = constants::moduleMLOT::length * mNumberOfModules; + + stave = new TGeoBBox(staveWidth / 2, mLogicalVolumeThickness / 2, staveLength / 2); staveVol = new TGeoVolume(staveName.c_str(), stave, medAir); + // Put the half staves in the correct position + TGeoVolume* halfStaveVolLeft = createHalfStave("flat"); + TGeoVolume* halfStaveVolRight = createHalfStave("flat"); + TGeoCombiTrans* transLeft = new TGeoCombiTrans(); - transLeft->SetTranslation(-mModuleWidth / 2 + 0.05, 0, 0); // 1mm overlap between the modules - LOGP(info, "Inserting {} in {} ", chipVolLeft->GetName(), staveVol->GetName()); - staveVol->AddNode(chipVolLeft, 0, transLeft); + transLeft->SetTranslation(-halfstaveWidth / 2 + 0.05, 0, 0); // TO BE CHECKED !!! 1mm overlap between the modules + LOGP(info, "Inserting {} in {} ", halfStaveVolLeft->GetName(), staveVol->GetName()); + staveVol->AddNode(halfStaveVolLeft, 0, transLeft); TGeoCombiTrans* transRight = new TGeoCombiTrans(); - transRight->SetTranslation(mModuleWidth / 2 - 0.05, 0.2, 0); - LOGP(info, "Inserting {} in {} ", chipVolRight->GetName(), staveVol->GetName()); - staveVol->AddNode(chipVolRight, 1, transRight); + transRight->SetTranslation(halfstaveWidth / 2 - 0.05, 0.2, 0); // TO BE CHECKED !!! 1mm overlap between the modules + LOGP(info, "Inserting {} in {} ", halfStaveVolRight->GetName(), staveVol->GetName()); + staveVol->AddNode(halfStaveVolRight, 1, transRight); } else { LOGP(fatal, "Chip of type '{}' is not implemented", type); } @@ -145,47 +330,49 @@ TGeoVolume* TRKLayer::createStave(std::string type, double width) void TRKLayer::createLayer(TGeoVolume* motherVolume) { - TGeoMedium* medSi = gGeoManager->GetMedium("TRK_SILICON$"); TGeoMedium* medAir = gGeoManager->GetMedium("TRK_AIR$"); - std::string staveName = o2::trk::GeometryTGeo::getTRKStavePattern() + std::to_string(mLayerNumber), - chipName = o2::trk::GeometryTGeo::getTRKChipPattern() + std::to_string(mLayerNumber), - sensName = Form("%s%d", GeometryTGeo::getTRKSensorPattern(), mLayerNumber); - double layerThickness = mChipThickness; if (mLayout != eLayout::kCylinder) { layerThickness = mLogicalVolumeThickness; } - TGeoTube* layer = new TGeoTube(mInnerRadius - 0.333 * layerThickness, mInnerRadius + 0.667 * layerThickness, mZ / 2); - TGeoVolume* layerVol = new TGeoVolume(mLayerName.c_str(), layer, medAir); - layerVol->SetLineColor(kYellow); + TGeoTube* layer; + TGeoVolume* layerVol; if (mLayout == eLayout::kCylinder) { - auto staveVol = createStave("cylinder"); + layer = new TGeoTube(mInnerRadius - 0.333 * layerThickness, mInnerRadius + 0.667 * layerThickness, mChipLength / 2); + layerVol = new TGeoVolume(mLayerName.c_str(), layer, medAir); + + TGeoVolume* staveVol = createStave("cylinder"); LOGP(info, "Inserting {} in {} ", staveVol->GetName(), layerVol->GetName()); layerVol->AddNode(staveVol, 1, nullptr); } else if (mLayout == eLayout::kTurboStaves) { - // Compute the number of staves - double width = mModuleWidth; // Each stave has two modules (based on the LOI design) + double layerLength = constants::moduleMLOT::length * mNumberOfModules; + double staveWidth = constants::ML::width; // Each stave has two modules (based on the LOI design) + if (mInnerRadius > 25) { - width *= 2; // Outer layers have two modules per stave + staveWidth = constants::OT::width; // Outer layers have two modules per stave } - int nStaves = (int)std::ceil(mInnerRadius * 2 * TMath::Pi() / width); + layer = new TGeoTube(mInnerRadius - 0.333 * layerThickness, mInnerRadius + 0.667 * layerThickness, layerLength / 2); + layerVol = new TGeoVolume(mLayerName.c_str(), layer, medAir); + + // Compute the number of staves + int nStaves = (int)std::ceil(mInnerRadius * 2 * TMath::Pi() / staveWidth); nStaves += nStaves % 2; // Require an even number of staves // Compute the size of the overlap region double theta = 2 * TMath::Pi() / nStaves; - double theta1 = std::atan(width / 2 / mInnerRadius); + double theta1 = std::atan(staveWidth / 2 / mInnerRadius); double st = std::sin(theta); double ct = std::cos(theta); - double theta2 = std::atan((mInnerRadius * st - width / 2 * ct) / (mInnerRadius * ct + width / 2 * st)); + double theta2 = std::atan((mInnerRadius * st - staveWidth / 2 * ct) / (mInnerRadius * ct + staveWidth / 2 * st)); double overlap = (theta1 - theta2) * mInnerRadius; LOGP(info, "Creating a layer with {} staves and {} mm overlap", nStaves, overlap * 10); for (int iStave = 0; iStave < nStaves; iStave++) { - TGeoVolume* staveVol = createStave("flat", width); + TGeoVolume* staveVol = createStave("flat"); // Put the staves in the correct position and orientation TGeoCombiTrans* trans = new TGeoCombiTrans(); @@ -198,17 +385,22 @@ void TRKLayer::createLayer(TGeoVolume* motherVolume) layerVol->AddNode(staveVol, iStave, trans); } } else if (mLayout == kStaggered) { + double layerLength = constants::moduleMLOT::length * mNumberOfModules; + + layer = new TGeoTube(mInnerRadius - 0.333 * layerThickness, mInnerRadius + 0.667 * layerThickness, layerLength / 2); + layerVol = new TGeoVolume(mLayerName.c_str(), layer, medAir); + // Compute the number of staves - double width = mModuleWidth * 2; // Each stave has two modules (based on the LOI design) - int nStaves = (int)std::ceil(mInnerRadius * 2 * TMath::Pi() / width); + double staveWidth = constants::OT::width; // Each stave has two modules (based on the LOI design) + int nStaves = (int)std::ceil(mInnerRadius * 2 * TMath::Pi() / staveWidth); nStaves += nStaves % 2; // Require an even number of staves // Compute the size of the overlap region double theta = 2 * TMath::Pi() / nStaves; - double theta1 = std::atan(width / 2 / mInnerRadius); + double theta1 = std::atan(staveWidth / 2 / mInnerRadius); double st = std::sin(theta); double ct = std::cos(theta); - double theta2 = std::atan((mInnerRadius * st - width / 2 * ct) / (mInnerRadius * ct + width / 2 * st)); + double theta2 = std::atan((mInnerRadius * st - staveWidth / 2 * ct) / (mInnerRadius * ct + staveWidth / 2 * st)); double overlap = (theta1 - theta2) * mInnerRadius; LOGP(info, "Creating a layer with {} staves and {} mm overlap", nStaves, overlap * 10); @@ -228,6 +420,8 @@ void TRKLayer::createLayer(TGeoVolume* motherVolume) } else { LOGP(fatal, "Layout not implemented"); } + layerVol->SetLineColor(kYellow); + LOGP(info, "Inserting {} in {} ", layerVol->GetName(), motherVolume->GetName()); motherVolume->AddNode(layerVol, 1, nullptr); } diff --git a/Detectors/Upgrades/ALICE3/TRK/simulation/src/TRKSimulationLinkDef.h b/Detectors/Upgrades/ALICE3/TRK/simulation/src/TRKSimulationLinkDef.h index 1a2e93636491c..9af868a2de44c 100644 --- a/Detectors/Upgrades/ALICE3/TRK/simulation/src/TRKSimulationLinkDef.h +++ b/Detectors/Upgrades/ALICE3/TRK/simulation/src/TRKSimulationLinkDef.h @@ -15,6 +15,9 @@ #pragma link off all classes; #pragma link off all functions; +#pragma link C++ class o2::trk::Hit + ; +#pragma link C++ class std::vector < o2::trk::Hit> + ; + #pragma link C++ class o2::trk::TRKLayer + ; #pragma link C++ class o2::trk::VDLayer + ; #pragma link C++ class o2::trk::TRKServices + ;