diff --git a/CODEOWNERS b/CODEOWNERS index 4ddd4cdf14692..49d64d1913892 100644 --- a/CODEOWNERS +++ b/CODEOWNERS @@ -41,7 +41,7 @@ /DataFormats/Detectors/TPC @davidrohr @wiechula @shahor02 /DataFormats/Detectors/TRD @f3sch @bazinski @tdietel @martenole /DataFormats/Detectors/Upgrades @qgp @marcovanleeuwen @mconcas -/DataFormats/Detectors/Upgrades/IT3 @fgrosa @arossi81 +/DataFormats/Detectors/Upgrades/ITS3 @fgrosa @arossi81 /DataFormats/Detectors/ZDC @coppedis #/DataFormats/Headers diff --git a/Detectors/ITSMFT/ITS/base/src/GeometryTGeo.cxx b/Detectors/ITSMFT/ITS/base/src/GeometryTGeo.cxx index 3431b1e3796b6..949e875011210 100644 --- a/Detectors/ITSMFT/ITS/base/src/GeometryTGeo.cxx +++ b/Detectors/ITSMFT/ITS/base/src/GeometryTGeo.cxx @@ -374,8 +374,11 @@ TGeoHMatrix* GeometryTGeo::extractMatrixSensor(int index) const if (mNumberOfModules[lay] > 0) { path += Form("%s%d_%d/", mIsLayerITS3[lay] ? getITS3ModulePattern() : getITSModulePattern(), lay, mod); } - path += - Form("%s%d_%d/%s%d_1", mIsLayerITS3[lay] ? getITS3ChipPattern() : getITSChipPattern(), lay, chipInMod, mIsLayerITS3[lay] ? getITS3SensorPattern() : getITSSensorPattern(), lay); + if (!mIsLayerITS3[lay]) { + path += Form("%s%d_%d/%s%d_1", getITSChipPattern(), lay, chipInMod, getITSSensorPattern(), lay); + } else { + path += Form("%s%d_%d", getITS3ChipPattern(), lay, chipInMod); // for ITS3 currently we might have more sensors than chips, so we have to take the chip to avoid mismatches with chipId + } static TGeoHMatrix matTmp; gGeoManager->PushPath(); @@ -388,8 +391,8 @@ TGeoHMatrix* GeometryTGeo::extractMatrixSensor(int index) const matTmp = *gGeoManager->GetCurrentMatrix(); // matrix may change after cd // RSS - // printf("%d/%d/%d %s\n",lay,stav,detInSta,path.Data()); - // mat->Print(); + // printf("%d/%d/%d %s\n", lay, stav, detInSta, path.Data()); + // matTmp.Print(); // Restore the modeler state. gGeoManager->PopPath(); @@ -397,7 +400,7 @@ TGeoHMatrix* GeometryTGeo::extractMatrixSensor(int index) const double delta = 0.; if (mIsLayerITS3[lay]) { #ifdef ENABLE_UPGRADES - delta = SegmentationITS3::SensorLayerThickness - SegmentationITS3::SensorLayerThicknessEff; + delta = SegmentationITS3::mSensorLayerThickness - SegmentationITS3::mSensorLayerThicknessEff; #endif } else { delta = Segmentation::SensorLayerThickness - Segmentation::SensorLayerThicknessEff; diff --git a/Detectors/Upgrades/CMakeLists.txt b/Detectors/Upgrades/CMakeLists.txt index 291c1bc465329..0123c1b3de901 100644 --- a/Detectors/Upgrades/CMakeLists.txt +++ b/Detectors/Upgrades/CMakeLists.txt @@ -10,5 +10,5 @@ # or submit itself to any jurisdiction. message(STATUS "Building detectors for upgrades") -add_subdirectory(IT3) +add_subdirectory(ITS3) add_subdirectory(ALICE3) \ No newline at end of file diff --git a/Detectors/Upgrades/IT3/base/include/ITS3Base/SegmentationSuperAlpide.h b/Detectors/Upgrades/IT3/base/include/ITS3Base/SegmentationSuperAlpide.h deleted file mode 100644 index 5605c4f8ec7e9..0000000000000 --- a/Detectors/Upgrades/IT3/base/include/ITS3Base/SegmentationSuperAlpide.h +++ /dev/null @@ -1,181 +0,0 @@ -// 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 SegmentationSuperAlpide.h -/// \brief Definition of the SegmentationSuperAlpide class - -#ifndef ALICEO2_ITS3_SEGMENTATIONSUPERALPIDE_H_ -#define ALICEO2_ITS3_SEGMENTATIONSUPERALPIDE_H_ - -#include -#include "MathUtils/Cartesian.h" -#include "CommonConstants/MathConstants.h" - -namespace o2 -{ -namespace its3 -{ - -/// Segmentation and response for pixels in ITS3 upgrade -class SegmentationSuperAlpide -{ - public: - SegmentationSuperAlpide(int layer = 0) : mLayer{layer}, - NPixels{NRows * NCols}, - NRows{static_cast(double(Radii[layer]) * double(constants::math::TwoPI) / double(PitchRow) + 1)}, - ActiveMatrixSizeRows{PitchRow * NRows}, - SensorSizeRows{ActiveMatrixSizeRows + PassiveEdgeTop + PassiveEdgeReadOut} - { - } - int mLayer; - static constexpr int NLayers = 4; - static constexpr float Length = 27.15f; - static constexpr float Radii[NLayers] = {1.8f, 2.4f, 3.0f, 7.0f}; - static constexpr float PitchCol = 29.24e-4; - static constexpr float PitchRow = 26.88e-4; - static constexpr int NCols = Length / PitchCol; - int NRows; - int NPixels; - static constexpr float PassiveEdgeReadOut = 0.; // width of the readout edge (Passive bottom) - static constexpr float PassiveEdgeTop = 0.; // Passive area on top - static constexpr float PassiveEdgeSide = 0.; // width of Passive area on left/right of the sensor - static constexpr float ActiveMatrixSizeCols = PitchCol * NCols; // Active size along columns - float ActiveMatrixSizeRows; // Active size along rows - - // effective thickness of sensitive layer, accounting for charge collection non-unifoemity, https://alice.its.cern.ch/jira/browse/AOC-46 - static constexpr float SensorLayerThicknessEff = 28.e-4; - static constexpr float SensorLayerThickness = 30.e-4; // physical thickness of sensitive part - static constexpr float SensorSizeCols = ActiveMatrixSizeCols + PassiveEdgeSide + PassiveEdgeSide; // SensorSize along columns - float SensorSizeRows; // SensorSize along rows - - ~SegmentationSuperAlpide() = default; - - /// Transformation from Geant detector centered local coordinates (cm) to - /// Pixel cell numbers iRow and iCol. - /// Returns kTRUE if point x,z is inside sensitive volume, kFALSE otherwise. - /// A value of -1 for iRow or iCol indicates that this point is outside of the - /// detector segmentation as defined. - /// \param float x Detector local coordinate x in cm with respect to - /// the center of the sensitive volume. - /// \param float z Detector local coordinate z in cm with respect to - /// the center of the sensitive volulme. - /// \param int iRow Detector x cell coordinate. Has the range 0 <= iRow < mNumberOfRows - /// \param int iCol Detector z cell coordinate. Has the range 0 <= iCol < mNumberOfColumns - bool localToDetector(float x, float z, int& iRow, int& iCol); - /// same but w/o check for row/column range - void localToDetectorUnchecked(float xRow, float zCol, int& iRow, int& iCol); - - /// Transformation from Detector cell coordiantes to Geant detector centered - /// local coordinates (cm) - /// \param int iRow Detector x cell coordinate. Has the range 0 <= iRow < mNumberOfRows - /// \param int iCol Detector z cell coordinate. Has the range 0 <= iCol < mNumberOfColumns - /// \param float x Detector local coordinate x in cm with respect to the - /// center of the sensitive volume. - /// \param float z Detector local coordinate z in cm with respect to the - /// center of the sensitive volulme. - /// If iRow and or iCol is outside of the segmentation range a value of -0.5*Dx() - /// or -0.5*Dz() is returned. - bool detectorToLocal(int iRow, int iCol, float& xRow, float& zCol); - bool detectorToLocal(float row, float col, float& xRow, float& zCol); - bool detectorToLocal(float row, float col, math_utils::Point3D& loc); - - // same but w/o check for row/col range - void detectorToLocalUnchecked(int iRow, int iCol, float& xRow, float& zCol); - void detectorToLocalUnchecked(float row, float col, float& xRow, float& zCol); - void detectorToLocalUnchecked(float row, float col, math_utils::Point3D& loc); - - float getFirstRowCoordinate() - { - return 0.5 * ((ActiveMatrixSizeRows - PassiveEdgeTop + PassiveEdgeReadOut) - PitchRow); - } - static constexpr float getFirstColCoordinate() { return 0.5 * (PitchCol - ActiveMatrixSizeCols); } - - void print(); - - ClassDefNV(SegmentationSuperAlpide, 1); // Segmentation class upgrade pixels -}; - -inline void SegmentationSuperAlpide::localToDetectorUnchecked(float xRow, float zCol, int& iRow, int& iCol) -{ - // convert to row/col w/o over/underflow check - xRow = 0.5 * (ActiveMatrixSizeRows - PassiveEdgeTop + PassiveEdgeReadOut) - xRow; // coordinate wrt top edge of Active matrix - zCol += 0.5 * ActiveMatrixSizeCols; // coordinate wrt left edge of Active matrix - iRow = int(xRow / PitchRow); - iCol = int(zCol / PitchCol); - if (xRow < 0) { - iRow -= 1; - } - if (zCol < 0) { - iCol -= 1; - } -} - -inline bool SegmentationSuperAlpide::localToDetector(float xRow, float zCol, int& iRow, int& iCol) -{ - // convert to row/col - xRow = 0.5 * (ActiveMatrixSizeRows - PassiveEdgeTop + PassiveEdgeReadOut) - xRow; // coordinate wrt left edge of Active matrix - zCol += 0.5 * ActiveMatrixSizeCols; // coordinate wrt bottom edge of Active matrix - if (xRow < 0 || xRow >= ActiveMatrixSizeRows || zCol < 0 || zCol >= ActiveMatrixSizeCols) { - iRow = iCol = -1; - return false; - } - iRow = int(xRow / PitchRow); - iCol = int(zCol / PitchCol); - return true; -} - -inline void SegmentationSuperAlpide::detectorToLocalUnchecked(int iRow, int iCol, float& xRow, float& zCol) -{ - xRow = getFirstRowCoordinate() - iRow * PitchRow; - zCol = iCol * PitchCol + getFirstColCoordinate(); -} - -inline void SegmentationSuperAlpide::detectorToLocalUnchecked(float row, float col, float& xRow, float& zCol) -{ - xRow = getFirstRowCoordinate() - row * PitchRow; - zCol = col * PitchCol + getFirstColCoordinate(); -} - -inline void SegmentationSuperAlpide::detectorToLocalUnchecked(float row, float col, math_utils::Point3D& loc) -{ - loc.SetCoordinates(getFirstRowCoordinate() - row * PitchRow, 0.f, col * PitchCol + getFirstColCoordinate()); -} - -inline bool SegmentationSuperAlpide::detectorToLocal(int iRow, int iCol, float& xRow, float& zCol) -{ - if (iRow < 0 || iRow >= NRows || iCol < 0 || iCol >= NCols) { - return false; - } - detectorToLocalUnchecked(iRow, iCol, xRow, zCol); - return true; -} - -inline bool SegmentationSuperAlpide::detectorToLocal(float row, float col, float& xRow, float& zCol) -{ - if (row < 0 || row >= NRows || col < 0 || col >= NCols) { - return false; - } - detectorToLocalUnchecked(row, col, xRow, zCol); - return true; -} - -inline bool SegmentationSuperAlpide::detectorToLocal(float row, float col, math_utils::Point3D& loc) -{ - if (row < 0 || row >= NRows || col < 0 || col >= NCols) { - return false; - } - detectorToLocalUnchecked(row, col, loc); - return true; -} -} // namespace its3 -} // namespace o2 - -#endif diff --git a/Detectors/Upgrades/IT3/macros/test/CMakeLists.txt b/Detectors/Upgrades/IT3/macros/test/CMakeLists.txt deleted file mode 100644 index db0b8a378d221..0000000000000 --- a/Detectors/Upgrades/IT3/macros/test/CMakeLists.txt +++ /dev/null @@ -1,33 +0,0 @@ -# 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. - -# o2_add_test_root_macro(CheckDigitsITS3.C -# PUBLIC_LINK_LIBRARIES O2::ITSBase -# O2::ITS3Base -# O2::ITSMFTBase -# O2::ITSMFTSimulation -# O2::ITS3Simulation -# O2::MathUtils -# O2::SimulationDataFormat -# O2::DetectorsBase -# LABELS its3) - -# o2_add_test_root_macro(CheckClustersITS3.C -# PUBLIC_LINK_LIBRARIES O2::ITSBase -# O2::ITS3Base -# O2::ITSMFTBase -# O2::ITSMFTSimulation -# O2::ITS3Simulation -# O2::ITS3Reconstruction -# O2::MathUtils -# O2::SimulationDataFormat -# O2::DetectorsBase -# LABELS its3) \ No newline at end of file diff --git a/Detectors/Upgrades/IT3/CMakeLists.txt b/Detectors/Upgrades/ITS3/CMakeLists.txt similarity index 87% rename from Detectors/Upgrades/IT3/CMakeLists.txt rename to Detectors/Upgrades/ITS3/CMakeLists.txt index 2d9ff8a9ac78e..0c91f53afb082 100644 --- a/Detectors/Upgrades/IT3/CMakeLists.txt +++ b/Detectors/Upgrades/ITS3/CMakeLists.txt @@ -9,5 +9,8 @@ # granted to it by virtue of its status as an Intergovernmental Organization # or submit itself to any jurisdiction. +add_subdirectory(macros) add_subdirectory(simulation) add_subdirectory(base) +add_subdirectory(workflow) +add_subdirectory(reconstruction) \ No newline at end of file diff --git a/Detectors/Upgrades/IT3/README.md b/Detectors/Upgrades/ITS3/README.md similarity index 99% rename from Detectors/Upgrades/IT3/README.md rename to Detectors/Upgrades/ITS3/README.md index 1a844ff7e8871..1811bdba48086 100644 --- a/Detectors/Upgrades/IT3/README.md +++ b/Detectors/Upgrades/ITS3/README.md @@ -2,7 +2,7 @@ \page refDetectorsUpgradesIT3 UpgradesIT3 /doxy --> -# IT3 +# ITS3 Upgraded version of the ITS that includes upgraded truly-cylindrical inner barrel. # Run the full simulation diff --git a/Detectors/Upgrades/IT3/base/CMakeLists.txt b/Detectors/Upgrades/ITS3/base/CMakeLists.txt similarity index 100% rename from Detectors/Upgrades/IT3/base/CMakeLists.txt rename to Detectors/Upgrades/ITS3/base/CMakeLists.txt diff --git a/Detectors/Upgrades/IT3/base/include/ITS3Base/MisalignmentParameter.h b/Detectors/Upgrades/ITS3/base/include/ITS3Base/MisalignmentParameter.h similarity index 100% rename from Detectors/Upgrades/IT3/base/include/ITS3Base/MisalignmentParameter.h rename to Detectors/Upgrades/ITS3/base/include/ITS3Base/MisalignmentParameter.h diff --git a/Detectors/Upgrades/ITS3/base/include/ITS3Base/SegmentationSuperAlpide.h b/Detectors/Upgrades/ITS3/base/include/ITS3Base/SegmentationSuperAlpide.h new file mode 100644 index 0000000000000..340e95d557bdd --- /dev/null +++ b/Detectors/Upgrades/ITS3/base/include/ITS3Base/SegmentationSuperAlpide.h @@ -0,0 +1,224 @@ +// 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 SegmentationSuperAlpide.h +/// \brief Definition of the SegmentationSuperAlpide class + +#ifndef ALICEO2_ITS3_SEGMENTATIONSUPERALPIDE_H_ +#define ALICEO2_ITS3_SEGMENTATIONSUPERALPIDE_H_ + +#include +#include "MathUtils/Cartesian.h" +#include "CommonConstants/MathConstants.h" +#include + +namespace o2 +{ +namespace its3 +{ + +/// Segmentation and response for pixels in ITS3 upgrade +class SegmentationSuperAlpide +{ + public: + SegmentationSuperAlpide(int layer = 0) : mLayer{layer}, + mNPixels{mNRows * mNCols}, + mNRows{static_cast(double(mRadii[layer] + mSensorLayerThickness / 2) * double(constants::math::PI) / double(mPitchRow) + 1)}, + mActiveMatrixSizeRows{mPitchRow * mNRows}, + mSensorSizeRows{mActiveMatrixSizeRows + mPassiveEdgeTop + mPassiveEdgeReadOut} + { + LOGP(info, "rows: {} cols: {} npixels: {}", mNRows, mNCols, mNPixels); + LOGP(info, "SegmentationSuperAlpide: layer {} ActiveMatrixSizeRows: {} ActiveMatrixSizeCols: {}", mLayer, mActiveMatrixSizeCols, mActiveMatrixSizeRows); + } + static constexpr std::array mRadii = {1.8f, 2.4f, 3.0f, 7.0f, 10.f}; + static constexpr float mLength = 27.15f; + static constexpr float mPitchCol = 20.e-4; + static constexpr float mPitchRow = 20.e-4; + static constexpr int mNCols = mLength / mPitchCol; + int mNRows; + int mNPixels; + int mLayer; + static constexpr float mPassiveEdgeReadOut = 0.; // width of the readout edge (Passive bottom) + static constexpr float mPassiveEdgeTop = 0.; // Passive area on top + static constexpr float mPassiveEdgeSide = 0.; // width of Passive area on left/right of the sensor + static constexpr float mActiveMatrixSizeCols = mPitchCol * mNCols; // Active size along columns + float mActiveMatrixSizeRows; // Active size along rows + + // effective thickness of sensitive layer, accounting for charge collection non-unifoemity, https://alice.its.cern.ch/jira/browse/AOC-46 + static constexpr float mSensorLayerThicknessEff = 28.e-4; + static constexpr float mSensorLayerThickness = 30.e-4; // physical thickness of sensitive part + static constexpr float mSensorSizeCols = mActiveMatrixSizeCols + mPassiveEdgeSide + mPassiveEdgeSide; // SensorSize along columns + float mSensorSizeRows; // SensorSize along rows + + ~SegmentationSuperAlpide() = default; + + /// Transformation from the curved surface to a flat surface + /// It works only if the detector is not rototraslated + /// \param xCurved Detector local curved coordinate x in cm with respect to + /// the center of the sensitive volume. + /// \param yCurved Detector local curved coordinate y in cm with respect to + /// the center of the sensitive volulme. + /// \param xFlat Detector local flat coordinate x in cm with respect to + /// the center of the sensitive volume. + /// \param yFlat Detector local flat coordinate y in cm with respect to + /// the center of the sensitive volulme. + void curvedToFlat(float xCurved, float yCurved, float& xFlat, float& yFlat); + + /// Transformation from the flat surface to a curved surface + /// It works only if the detector is not rototraslated + /// \param xFlat Detector local flat coordinate x in cm with respect to + /// the center of the sensitive volume. + /// \param yFlat Detector local flat coordinate y in cm with respect to + /// the center of the sensitive volulme. + /// \param xCurved Detector local curved coordinate x in cm with respect to + /// the center of the sensitive volume. + /// \param yCurved Detector local curved coordinate y in cm with respect to + /// the center of the sensitive volulme. + void flatToCurved(float xFlat, float yFlat, float& xCurved, float& yCurved); + + /// Transformation from Geant detector centered local coordinates (cm) to + /// Pixel cell numbers iRow and iCol. + /// Returns kTRUE if point x,z is inside sensitive volume, kFALSE otherwise. + /// A value of -1 for iRow or iCol indicates that this point is outside of the + /// detector segmentation as defined. + /// \param float x Detector local coordinate x in cm with respect to + /// the center of the sensitive volume. + /// \param float z Detector local coordinate z in cm with respect to + /// the center of the sensitive volulme. + /// \param int iRow Detector x cell coordinate. Has the range 0 <= iRow < mNumberOfRows + /// \param int iCol Detector z cell coordinate. Has the range 0 <= iCol < mNumberOfColumns + bool localToDetector(float x, float z, int& iRow, int& iCol); + /// same but w/o check for row/column range + void localToDetectorUnchecked(float xRow, float zCol, int& iRow, int& iCol); + + /// Transformation from Detector cell coordiantes to Geant detector centered + /// local coordinates (cm) + /// \param int iRow Detector x cell coordinate. Has the range 0 <= iRow < mNumberOfRows + /// \param int iCol Detector z cell coordinate. Has the range 0 <= iCol < mNumberOfColumns + /// \param float x Detector local coordinate x in cm with respect to the + /// center of the sensitive volume. + /// \param float z Detector local coordinate z in cm with respect to the + /// center of the sensitive volulme. + /// If iRow and or iCol is outside of the segmentation range a value of -0.5*Dx() + /// or -0.5*Dz() is returned. + bool detectorToLocal(int iRow, int iCol, float& xRow, float& zCol); + bool detectorToLocal(float row, float col, float& xRow, float& zCol); + bool detectorToLocal(float row, float col, math_utils::Point3D& loc); + + // same but w/o check for row/col range + void detectorToLocalUnchecked(int iRow, int iCol, float& xRow, float& zCol); + void detectorToLocalUnchecked(float row, float col, float& xRow, float& zCol); + void detectorToLocalUnchecked(float row, float col, math_utils::Point3D& loc); + + float getFirstRowCoordinate() + { + return 0.5 * ((mActiveMatrixSizeRows - mPassiveEdgeTop + mPassiveEdgeReadOut) - mPitchRow); + } + static constexpr float getFirstColCoordinate() { return 0.5 * (mPitchCol - mActiveMatrixSizeCols); } + + void print(); + + ClassDefNV(SegmentationSuperAlpide, 1); // Segmentation class upgrade pixels +}; + +inline void SegmentationSuperAlpide::curvedToFlat(float xCurved, float yCurved, float& xFlat, float& yFlat) +{ + float dist = std::sqrt(xCurved * xCurved + yCurved * yCurved); + float phi = (double)constants::math::PI / 2 - std::atan2((double)yCurved, (double)xCurved); + xFlat = (mRadii[mLayer] + SegmentationSuperAlpide::mSensorLayerThickness / 2) * phi; + yFlat = dist - (mRadii[mLayer] + SegmentationSuperAlpide::mSensorLayerThickness / 2); +} + +inline void SegmentationSuperAlpide::flatToCurved(float xFlat, float yFlat, float& xCurved, float& yCurved) +{ + float phi = xFlat / (mRadii[mLayer] + SegmentationSuperAlpide::mSensorLayerThickness / 2); + float dist = yFlat + (mRadii[mLayer] + SegmentationSuperAlpide::mSensorLayerThickness / 2); + float tang = std::tan((double)constants::math::PI / 2 - (double)phi); + xCurved = (xFlat > 0 ? 1.f : -1.f) * dist / std::sqrt(1 + tang * tang); + yCurved = xCurved * tang; +} + +inline void SegmentationSuperAlpide::localToDetectorUnchecked(float xRow, float zCol, int& iRow, int& iCol) +{ + // convert to row/col w/o over/underflow check + xRow = 0.5 * (mActiveMatrixSizeRows - mPassiveEdgeTop + mPassiveEdgeReadOut) - xRow; // coordinate wrt top edge of Active matrix + zCol += 0.5 * mActiveMatrixSizeCols; // coordinate wrt left edge of Active matrix + iRow = int(xRow / mPitchRow); + iCol = int(zCol / mPitchCol); + if (xRow < 0) { + iRow -= 1; + } + if (zCol < 0) { + iCol -= 1; + } +} + +inline bool SegmentationSuperAlpide::localToDetector(float xRow, float zCol, int& iRow, int& iCol) +{ + // convert to row/col + xRow = 0.5 * (mActiveMatrixSizeRows - mPassiveEdgeTop + mPassiveEdgeReadOut) - xRow; // coordinate wrt left edge of Active matrix + zCol += 0.5 * mActiveMatrixSizeCols; // coordinate wrt bottom edge of Active matrix + if (xRow < 0 || xRow >= mActiveMatrixSizeRows || zCol < 0 || zCol >= mActiveMatrixSizeCols) { + iRow = iCol = -1; + return false; + } + iRow = int(xRow / mPitchRow); + iCol = int(zCol / mPitchCol); + return true; +} + +inline void SegmentationSuperAlpide::detectorToLocalUnchecked(int iRow, int iCol, float& xRow, float& zCol) +{ + xRow = getFirstRowCoordinate() - iRow * mPitchRow; + zCol = iCol * mPitchCol + getFirstColCoordinate(); +} + +inline void SegmentationSuperAlpide::detectorToLocalUnchecked(float row, float col, float& xRow, float& zCol) +{ + xRow = getFirstRowCoordinate() - row * mPitchRow; + zCol = col * mPitchCol + getFirstColCoordinate(); +} + +inline void SegmentationSuperAlpide::detectorToLocalUnchecked(float row, float col, math_utils::Point3D& loc) +{ + loc.SetCoordinates(getFirstRowCoordinate() - row * mPitchRow, 0.f, col * mPitchCol + getFirstColCoordinate()); +} + +inline bool SegmentationSuperAlpide::detectorToLocal(int iRow, int iCol, float& xRow, float& zCol) +{ + if (iRow < 0 || iRow >= mNRows || iCol < 0 || iCol >= mNCols) { + return false; + } + detectorToLocalUnchecked(iRow, iCol, xRow, zCol); + return true; +} + +inline bool SegmentationSuperAlpide::detectorToLocal(float row, float col, float& xRow, float& zCol) +{ + if (row < 0 || row >= mNRows || col < 0 || col >= mNCols) { + return false; + } + detectorToLocalUnchecked(row, col, xRow, zCol); + return true; +} + +inline bool SegmentationSuperAlpide::detectorToLocal(float row, float col, math_utils::Point3D& loc) +{ + if (row < 0 || row >= mNRows || col < 0 || col >= mNCols) { + return false; + } + detectorToLocalUnchecked(row, col, loc); + return true; +} +} // namespace its3 +} // namespace o2 + +#endif diff --git a/Detectors/Upgrades/IT3/base/src/ITS3BaseLinkDef.h b/Detectors/Upgrades/ITS3/base/src/ITS3BaseLinkDef.h similarity index 99% rename from Detectors/Upgrades/IT3/base/src/ITS3BaseLinkDef.h rename to Detectors/Upgrades/ITS3/base/src/ITS3BaseLinkDef.h index 5af984aa0c3c8..9aeaec7d3bb76 100644 --- a/Detectors/Upgrades/IT3/base/src/ITS3BaseLinkDef.h +++ b/Detectors/Upgrades/ITS3/base/src/ITS3BaseLinkDef.h @@ -18,6 +18,4 @@ #pragma link C++ class o2::its3::SegmentationSuperAlpide + ; #pragma link C++ class o2::its3::MisalignmentParameter + ; - - #endif diff --git a/Detectors/Upgrades/IT3/base/src/MisalignmentParameter.cxx b/Detectors/Upgrades/ITS3/base/src/MisalignmentParameter.cxx similarity index 100% rename from Detectors/Upgrades/IT3/base/src/MisalignmentParameter.cxx rename to Detectors/Upgrades/ITS3/base/src/MisalignmentParameter.cxx diff --git a/Detectors/Upgrades/IT3/base/src/SegmentationSuperAlpide.cxx b/Detectors/Upgrades/ITS3/base/src/SegmentationSuperAlpide.cxx similarity index 78% rename from Detectors/Upgrades/IT3/base/src/SegmentationSuperAlpide.cxx rename to Detectors/Upgrades/ITS3/base/src/SegmentationSuperAlpide.cxx index 10ce7f9a98fb4..29fc0c11b739d 100644 --- a/Detectors/Upgrades/IT3/base/src/SegmentationSuperAlpide.cxx +++ b/Detectors/Upgrades/ITS3/base/src/SegmentationSuperAlpide.cxx @@ -21,10 +21,9 @@ using namespace o2::its3; void SegmentationSuperAlpide::print() { - printf("ITS3 segmentation for layer: %d \n", mLayer); - printf("Pixel size: %.2f (along %d rows) %.2f (along %d columns) microns\n", PitchRow * 1e4, NRows, PitchCol * 1e4, NCols); + printf("Pixel size: %.2f (along %d rows) %.2f (along %d columns) microns\n", mPitchRow * 1e4, mNRows, mPitchCol * 1e4, mNCols); printf("Passive edges: bottom: %.2f, top: %.2f, left/right: %.2f microns\n", - PassiveEdgeReadOut * 1e4, PassiveEdgeTop * 1e4, PassiveEdgeSide * 1e4); - printf("Active/Total size: %.6f/%.6f (rows) %.6f/%.6f (cols) cm\n", ActiveMatrixSizeRows, SensorSizeRows, - ActiveMatrixSizeCols, SensorSizeCols); + mPassiveEdgeReadOut * 1e4, mPassiveEdgeTop * 1e4, mPassiveEdgeSide * 1e4); + printf("Active/Total size: %.6f/%.6f (rows) %.6f/%.6f (cols) cm\n", mActiveMatrixSizeRows, mSensorSizeRows, + mActiveMatrixSizeCols, mSensorSizeCols); } diff --git a/Detectors/Upgrades/IT3/macros/CMakeLists.txt b/Detectors/Upgrades/ITS3/macros/CMakeLists.txt similarity index 100% rename from Detectors/Upgrades/IT3/macros/CMakeLists.txt rename to Detectors/Upgrades/ITS3/macros/CMakeLists.txt diff --git a/Detectors/Upgrades/ITS3/macros/test/CMakeLists.txt b/Detectors/Upgrades/ITS3/macros/test/CMakeLists.txt new file mode 100644 index 0000000000000..87ee24626121d --- /dev/null +++ b/Detectors/Upgrades/ITS3/macros/test/CMakeLists.txt @@ -0,0 +1,45 @@ +# 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. + +o2_add_test_root_macro(CheckDigitsITS3.C + PUBLIC_LINK_LIBRARIES O2::ITSBase + O2::ITS3Base + O2::ITSMFTBase + O2::ITSMFTSimulation + O2::ITS3Simulation + O2::MathUtils + O2::SimulationDataFormat + O2::DetectorsBase + LABELS its3) + +o2_add_test_root_macro(CheckClustersITS3.C + PUBLIC_LINK_LIBRARIES O2::ITSBase + O2::ITS3Base + O2::ITSMFTBase + O2::ITSMFTSimulation + O2::ITS3Simulation + O2::ITS3Reconstruction + O2::MathUtils + O2::SimulationDataFormat + O2::DetectorsBase + LABELS its3) + + # o2_add_test_root_macro(CheckSquasherITS3.C +# PUBLIC_LINK_LIBRARIES O2::ITSBase +# O2::ITS3Base +# O2::ITSMFTBase +# O2::ITSMFTSimulation +# O2::ITS3Simulation +# O2::ITS3Reconstruction +# O2::MathUtils +# O2::SimulationDataFormat +# O2::DetectorsBase +# LABELS its3) \ No newline at end of file diff --git a/Detectors/Upgrades/IT3/macros/test/CheckClustersITS3.notest b/Detectors/Upgrades/ITS3/macros/test/CheckClustersITS3.C similarity index 71% rename from Detectors/Upgrades/IT3/macros/test/CheckClustersITS3.notest rename to Detectors/Upgrades/ITS3/macros/test/CheckClustersITS3.C index 7c8d55eaa7b51..009b42ec7c413 100644 --- a/Detectors/Upgrades/IT3/macros/test/CheckClustersITS3.notest +++ b/Detectors/Upgrades/ITS3/macros/test/CheckClustersITS3.C @@ -1,3 +1,14 @@ +// 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 CheckClusters.C /// \brief Simple macro to check ITSU clusters @@ -9,6 +20,8 @@ #include #include #include +#include +#include #include "DetectorsCommonDataFormats/DetID.h" #include "ITSMFTBase/SegmentationAlpide.h" @@ -26,10 +39,10 @@ #include "DetectorsCommonDataFormats/DetectorNameConf.h" #endif -void CheckClustersITS3(std::string clusfile = "o2clus_it3.root", std::string hitfile = "o2sim_HitsIT3.root", - std::string inputGeom = "", std::string dictfile = "") +void CheckClustersITS3(int nITS3layers = 3, std::string clusfile = "o2clus_it3.root", std::string hitfile = "o2sim_HitsIT3.root", + std::string inputGeom = "", std::string dictfile = "", bool batch = true) { - const int QEDSourceID = 99; // Clusters from this MC source correspond to QED electrons + gROOT->SetBatch(batch); using namespace o2::base; using namespace o2::its; @@ -46,6 +59,17 @@ void CheckClustersITS3(std::string clusfile = "o2clus_it3.root", std::string hit std::vector hitVecPool; std::vector mc2hitVec; + // we assume that we have 2 chips per layer + const int nChipsPerLayer = 2; + std::vector segs{}; + for (int iLayer{0}; iLayer < nITS3layers; ++iLayer) { + for (int iChip{0}; iChip < nChipsPerLayer; ++iChip) { + segs.push_back(SegmentationSuperAlpide(iLayer)); + } + } + + const int QEDSourceID = 99; // Clusters from this MC source correspond to QED electrons + TFile fout("CheckClusters.root", "recreate"); TNtuple nt("ntc", "cluster ntuple", "ev:lab:hlx:hlz:tx:tz:cgx:cgy:cgz:dx:dy:dz:ex:ez:patid:rof:npx:id"); @@ -158,7 +182,7 @@ void CheckClustersITS3(std::string clusfile = "o2clus_it3.root", std::string hit auto pattID = cluster.getPatternID(); o2::math_utils::Point3D locC; auto chipID = cluster.getSensorID(); - if (pattID == o2::itsmft::CompCluster::InvalidPatternID || dict.isGroup(pattID)) { + if (pattID == o2::its3::CompCluster::InvalidPatternID || dict.isGroup(pattID)) { o2::itsmft::ClusterPattern patt(pattIt); locC = dict.getClusterCoordinates(chipID, cluster, patt, false); } else { @@ -166,15 +190,16 @@ void CheckClustersITS3(std::string clusfile = "o2clus_it3.root", std::string hit errX = dict.getErrX(pattID); errZ = dict.getErrZ(pattID); npix = dict.getNpixels(pattID); + LOGP(info, "I am invalid and I am on chip {}", chipID); } - // Transformation to the local --> global - auto gloC = gman->getMatrixL2G(chipID) * locC; - if (chipID < o2::its3::SegmentationSuperAlpide::NLayers) { - double radius = SegmentationSuperAlpide::Radii[chipID]; - double phi = locC.X() / radius; - gloC.SetXYZ(radius * std::cos(phi), radius * std::sin(phi), locC.Z()); + if (chipID / nChipsPerLayer < nITS3layers) { + float xCurved{0.f}, yCurved{0.f}; + segs[chipID].flatToCurved(locC.X(), locC.Y(), xCurved, yCurved); + locC.SetXYZ(xCurved, yCurved, locC.Z()); } + // Transformation to the local --> global + auto gloC = gman->getMatrixL2G(chipID)(locC); const auto& lab = (clusLabArr->getLabels(clEntry))[0]; if (!lab.isValid() || lab.getSourceID() == QEDSourceID) @@ -195,24 +220,24 @@ void CheckClustersITS3(std::string clusfile = "o2clus_it3.root", std::string hit float dx = 0, dz = 0; int ievH = lab.getEventID(); o2::math_utils::Point3D locH, locHsta; + o2::math_utils::Point3D gloH, gloHsta; // mean local position of the hit locH = gman->getMatrixL2G(chipID) ^ (hit.GetPos()); // inverse conversion from global to local locHsta = gman->getMatrixL2G(chipID) ^ (hit.GetPosStart()); - if (chipID < 4) { - float startPhi{std::atan2(-hit.GetPosStart().Y(), -hit.GetPosStart().X())}; - float endPhi{std::atan2(-hit.GetPos().Y(), -hit.GetPos().X())}; - locH.SetXYZ(SegmentationSuperAlpide::Radii[chipID] * endPhi, 0.f, hit.GetPos().Z()); - locHsta.SetXYZ(SegmentationSuperAlpide::Radii[chipID] * startPhi, 0.f, hit.GetPosStart().Z()); - } - auto x0 = locHsta.X(), dltx = locH.X() - x0; auto y0 = locHsta.Y(), dlty = locH.Y() - y0; auto z0 = locHsta.Z(), dltz = locH.Z() - z0; - auto r = (0.5 * (Segmentation::SensorLayerThickness - Segmentation::SensorLayerThicknessEff) - y0) / dlty; - locH.SetXYZ(x0 + r * dltx, y0 + r * dlty, z0 + r * dltz); - //locH.SetXYZ(0.5 * (locH.X() + locHsta.X()), 0.5 * (locH.Y() + locHsta.Y()), 0.5 * (locH.Z() + locHsta.Z())); + + if (chipID / nChipsPerLayer >= nITS3layers) { + auto r = (0.5 * (Segmentation::SensorLayerThickness - Segmentation::SensorLayerThicknessEff) - y0) / dlty; + locH.SetXYZ(x0 + r * dltx, y0 + r * dlty, z0 + r * dltz); + } else { + // not really precise, but okish + locH.SetXYZ(0.5 * (locH.X() + locHsta.X()), 0.5 * (locH.Y() + locHsta.Y()), 0.5 * (locH.Z() + locHsta.Z())); + } + std::array data = {(float)lab.getEventID(), (float)trID, locH.X(), locH.Z(), dltx / dlty, dltz / dlty, gloC.X(), gloC.Y(), gloC.Z(), @@ -223,30 +248,34 @@ void CheckClustersITS3(std::string clusfile = "o2clus_it3.root", std::string hit } } - new TCanvas; - nt.Draw("cgy:cgx"); - new TCanvas; - nt.Draw("dz:dx", "abs(dz)<0.01 && abs(dx)<0.01"); - new TCanvas; - nt.Draw("dz:tz", "abs(dz)<0.005 && abs(tz)<2"); - - auto c1 = new TCanvas("p1", "pullX"); - c1->cd(); - c1->SetLogy(); - nt.Draw("dx/ex", "abs(dx/ex)<10&&patid<10"); - auto c2 = new TCanvas("p2", "pullZ"); - c2->cd(); - c2->SetLogy(); - nt.Draw("dz/ez", "abs(dz/ez)<10&&patid<10"); - - auto d1 = new TCanvas("d1", "deltaX"); - d1->cd(); - d1->SetLogy(); - nt.Draw("dx", "abs(dx)<5"); - auto d2 = new TCanvas("d2", "deltaZ"); - d2->cd(); - d2->SetLogy(); - nt.Draw("dz", "abs(dz)<5"); + auto canvCgXCgY = new TCanvas("canvCgXCgY", "", 1600, 1600); + canvCgXCgY->Divide(2, 2); + canvCgXCgY->cd(1); + nt.Draw("cgy:cgx>>h_cgy_vs_cgx_IB(1000, -10, 10, 1000, -10, 10)", "id < 6", "colz"); + canvCgXCgY->cd(2); + nt.Draw("cgy:cgz>>h_cgy_vs_cgz_IB(1000, -15, 15, 1000, -10, 10)", "id < 6", "colz"); + canvCgXCgY->cd(3); + nt.Draw("cgy:cgx>>h_cgy_vs_cgx_OB(1000, -50, 50, 1000, -50, 50)", "id >= 6", "colz"); + canvCgXCgY->cd(4); + nt.Draw("cgy:cgz>>h_cgy_vs_cgz_OB(1000, -100, 100, 1000, -50, 50)", "id >= 6", "colz"); + canvCgXCgY->SaveAs("it3clusters_y_vs_x_vs_z.pdf"); + + auto canvdXdZ = new TCanvas("canvdXdZ", "", 1600, 800); + canvdXdZ->Divide(2, 1); + canvdXdZ->cd(1); + nt.Draw("dx:dz>>h_dx_vs_dz_IB(1000, -0.0025, 0.0025, 1000, -0.0025, 0.0025)", "id < 6", "colz"); + canvdXdZ->cd(2); + nt.Draw("dx:dz>>h_dx_vs_dz_OB(1000, -0.0025, 0.0025, 1000, -0.0025, 0.0025)", "id >= 6", "colz"); + canvdXdZ->SaveAs("it3clusters_dx_vs_dz.pdf"); + + // auto c1 = new TCanvas("p1", "pullX"); + // c1->cd(); + // c1->SetLogy(); + // nt.Draw("dx/ex", "abs(dx/ex)<10&&patid<10"); + // auto c2 = new TCanvas("p2", "pullZ"); + // c2->cd(); + // c2->SetLogy(); + // nt.Draw("dz/ez", "abs(dz/ez)<10&&patid<10"); fout.cd(); nt.Write(); diff --git a/Detectors/Upgrades/IT3/macros/test/CheckDigitsITS3.notest b/Detectors/Upgrades/ITS3/macros/test/CheckDigitsITS3.C similarity index 67% rename from Detectors/Upgrades/IT3/macros/test/CheckDigitsITS3.notest rename to Detectors/Upgrades/ITS3/macros/test/CheckDigitsITS3.C index bac555eec6e2e..abb44876d127d 100644 --- a/Detectors/Upgrades/IT3/macros/test/CheckDigitsITS3.notest +++ b/Detectors/Upgrades/ITS3/macros/test/CheckDigitsITS3.C @@ -1,3 +1,14 @@ +// 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 CheckDigitsITS3.C /// \brief Simple macro to check ITS3 digits @@ -9,9 +20,10 @@ #include #include #include +#include #include -#include "ITS3Base/GeometryTGeo.h" +#include "ITSBase/GeometryTGeo.h" #include "DataFormatsITSMFT/Digit.h" #include "ITS3Base/SegmentationSuperAlpide.h" #include "ITSMFTBase/SegmentationAlpide.h" @@ -26,8 +38,9 @@ #endif -void CheckDigitsITS3(std::string digifile = "it3digits.root", std::string hitfile = "o2sim_HitsIT3.root", std::string inputGeom = "", std::string paramfile = "o2sim_par.root") +void CheckDigitsITS3(int nITS3layers = 3, std::string digifile = "it3digits.root", std::string hitfile = "o2sim_HitsIT3.root", std::string inputGeom = "", std::string paramfile = "o2sim_par.root", bool batch = true) { + gROOT->SetBatch(batch); using namespace o2::base; @@ -38,16 +51,22 @@ void CheckDigitsITS3(std::string digifile = "it3digits.root", std::string hitfil using o2::itsmft::SegmentationAlpide; TFile* f = TFile::Open("CheckDigits.root", "recreate"); - TNtuple* nt = new TNtuple("ntd", "digit ntuple", "id:x:y:z:rowD:colD:rowH:colH:xlH:zlH:xlcH:zlcH:dx:dz"); // Geometry o2::base::GeometryManager::loadGeometry(inputGeom); - auto* gman = o2::its3::GeometryTGeo::Instance(); + auto* gman = o2::its::GeometryTGeo::Instance(); gman->fillMatrixCache(o2::math_utils::bit2Mask(o2::math_utils::TransformType::L2G)); - SegmentationSuperAlpide segs[4]{SegmentationSuperAlpide(0), SegmentationSuperAlpide(1), SegmentationSuperAlpide(2), SegmentationSuperAlpide(3)}; - SegmentationSuperAlpide& seg = segs[0]; + // we assume that we have 2 chips per layer + const int nChipsPerLayer = 2; + + std::vector segs{}; + for (int iLayer{0}; iLayer < nITS3layers; ++iLayer) { + for (int iChip{0}; iChip < nChipsPerLayer; ++iChip) { + segs.push_back(SegmentationSuperAlpide(iLayer)); + } + } // Hits TFile* hitFile = TFile::Open(hitfile.data()); @@ -148,21 +167,21 @@ void CheckDigitsITS3(std::string digifile = "it3digits.root", std::string hitfil // LOOP on : digits array for (unsigned int iDigit = rofIndex; iDigit < rofIndex + rofNEntries; iDigit++) { - Int_t ix = (*digArr)[iDigit].getRow(), iz = (*digArr)[iDigit].getColumn(); - Float_t x = 0.f, z = 0.f; + int ix = (*digArr)[iDigit].getRow(), iz = (*digArr)[iDigit].getColumn(); + float x{0.f}, y{0.f}, z{0.f}; - Int_t chipID = (*digArr)[iDigit].getChipIndex(); + int chipID = (*digArr)[iDigit].getChipIndex(); - if (chipID < 4) { - segs[chipID].detectorToLocal(ix, iz, x, z); + if (chipID / nChipsPerLayer < nITS3layers) { + float xFlat{0.f}; + segs[chipID].detectorToLocal(ix, iz, xFlat, z); + segs[chipID].flatToCurved(xFlat, 0., x, y); } else { SegmentationAlpide::detectorToLocal(ix, iz, x, z); } - const o2::math_utils::Point3D locD(x, 0., z); - + const o2::math_utils::Point3D locD(x, y, z); auto lab = (labels.getLabels(iDigit))[0]; - int trID = lab.getTrackID(); if (lab.isValid()) { // not a noise @@ -170,16 +189,6 @@ void CheckDigitsITS3(std::string digifile = "it3digits.root", std::string hitfil nDigitRead++; auto gloD = gman->getMatrixL2G(chipID)(locD); // convert to global - if (chipID < 4) { - // - // invert - //xyzLocS = {SegmentationSuperAlpide::Radii[detID] * startPhi, 0.f, startPos.Z()}; - //xyzLocE = {SegmentationSuperAlpide::Radii[detID] * endPhi, 0.f, endPos.Z()}; - // - double radius = SegmentationSuperAlpide::Radii[chipID]; - double phi = locD.X() / radius; - gloD.SetXYZ(radius * std::cos(phi), radius * std::sin(phi), locD.Z()); - } float dx = 0., dz = 0.; std::unordered_map* mc2hit = &mc2hitVec[lab.getEventID()]; @@ -194,26 +203,24 @@ void CheckDigitsITS3(std::string digifile = "it3digits.root", std::string hitfil continue; } + ////// HITS Hit& hit = (*hitArray[lab.getEventID()])[hitEntry->second]; - // FIXME FOR INNER BARREL auto locH = gman->getMatrixL2G(chipID) ^ (hit.GetPos()); // inverse conversion from global to local auto locHsta = gman->getMatrixL2G(chipID) ^ (hit.GetPosStart()); - if (chipID < 4) { - float startPhi{std::atan2(-hit.GetPosStart().Y(), -hit.GetPosStart().X())}; - float endPhi{std::atan2(-hit.GetPos().Y(), -hit.GetPos().X())}; - locH.SetXYZ(SegmentationSuperAlpide::Radii[chipID] * endPhi, 0.f, hit.GetPos().Z()); - locHsta.SetXYZ(SegmentationSuperAlpide::Radii[chipID] * startPhi, 0.f, hit.GetPosStart().Z()); - } - locH.SetXYZ(0.5 * (locH.X() + locHsta.X()), 0.5 * (locH.Y() + locHsta.Y()), 0.5 * (locH.Z() + locHsta.Z())); int row, col; float xlc = 0., zlc = 0.; - segs[chipID < 4 ? chipID : 0].localToDetector(locH.X(), locH.Z(), row, col); - segs[chipID < 4 ? chipID : 0].detectorToLocal(row, col, xlc, zlc); + if (chipID / nChipsPerLayer < nITS3layers) { + segs[chipID].localToDetector(locH.X(), locH.Z(), row, col); + segs[chipID].detectorToLocal(row, col, xlc, zlc); + } else { + SegmentationAlpide::localToDetector(locH.X(), locH.Z(), row, col); + SegmentationAlpide::detectorToLocal(row, col, xlc, zlc); + } nt->Fill(chipID, gloD.X(), gloD.Y(), gloD.Z(), ix, iz, row, col, locH.X(), locH.Z(), xlc, zlc, locH.X() - locD.X(), locH.Z() - locD.Z()); @@ -224,10 +231,25 @@ void CheckDigitsITS3(std::string digifile = "it3digits.root", std::string hitfil } // end loop on ROFRecords array - new TCanvas; - nt->Draw("y:x"); - new TCanvas; - nt->Draw("dx:dz", "abs(dx)<0.02 && abs(dz)<0.02"); + auto canvXY = new TCanvas("canvXY", "", 1600, 1600); + canvXY->Divide(2, 2); + canvXY->cd(1); + nt->Draw("y:x>>h_y_vs_x_IB(1000, -10, 10, 1000, -10, 10)", "id < 6", "colz"); + canvXY->cd(2); + nt->Draw("y:z>>h_y_vs_z_IB(1000, -15, 15, 1000, -10, 10)", "id < 6", "colz"); + canvXY->cd(3); + nt->Draw("y:x>>h_y_vs_x_OB(1000, -50, 50, 1000, -50, 50)", "id >= 6", "colz"); + canvXY->cd(4); + nt->Draw("y:z>>h_y_vs_z_OB(1000, -100, 100, 1000, -50, 50)", "id >= 6", "colz"); + canvXY->SaveAs("it3digits_y_vs_x_vs_z.pdf"); + + auto canvdXdZ = new TCanvas("canvdXdZ", "", 1600, 800); + canvdXdZ->Divide(2, 1); + canvdXdZ->cd(1); + nt->Draw("dx:dz>>h_dx_vs_dz_IB(1000, -0.025, 0.025, 1000, -0.025, 0.025)", "id < 6", "colz"); + canvdXdZ->cd(2); + nt->Draw("dx:dz>>h_dx_vs_dz_OB(1000, -0.025, 0.025, 1000, -0.025, 0.025)", "id >= 6", "colz"); + canvdXdZ->SaveAs("it3digits_dx_vs_dz.pdf"); f->Write(); f->Close(); diff --git a/Detectors/Upgrades/ITS3/macros/test/CheckSquasherITS3.notest b/Detectors/Upgrades/ITS3/macros/test/CheckSquasherITS3.notest new file mode 100644 index 0000000000000..b36596abfab21 --- /dev/null +++ b/Detectors/Upgrades/ITS3/macros/test/CheckSquasherITS3.notest @@ -0,0 +1,183 @@ +// 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. + +#if !defined(__CLING__) || defined(__ROOTCLING__) +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#endif + +void getClusterPatterns(std::vector& pattVec, std::vector* ITSclus, std::vector* ITSpatt, o2::itsmft::TopologyDictionary& mdict); + +void CheckSquasherITS3(const uint chipId = 0, const uint startingROF = 0, const unsigned int nRofs = 3, const string fname = "o2clus_it3.root") +{ + // TColor::InvertPalette(); + gStyle->SetOptStat(0); + gStyle->SetPalette(kInvertedDarkBodyRadiator); + // Geometry + o2::base::GeometryManager::loadGeometry(""); + auto gman = o2::its::GeometryTGeo::Instance(); + gman->fillMatrixCache(o2::math_utils::bit2Mask(o2::math_utils::TransformType::T2L, o2::math_utils::TransformType::L2G)); + // Topology dictionary + auto& cc = o2::ccdb::BasicCCDBManager::instance(); + auto mdict = cc.get("ITS/Calib/ClusterDictionary"); + auto fITSclus = TFile::Open(fname.data(), "r"); + auto treeITSclus = (TTree*)fITSclus->Get("o2sim"); + + std::vector* ITSclus = nullptr; + std::vector* ITSrof = nullptr; + std::vector* ITSpatt = nullptr; + o2::dataformats::MCTruthContainer* clusLabArr = nullptr; + + treeITSclus->SetBranchAddress("IT3ClusterComp", &ITSclus); + treeITSclus->SetBranchAddress("IT3ClustersROF", &ITSrof); + treeITSclus->SetBranchAddress("IT3ClusterPatt", &ITSpatt); + // treeITSclus->SetBranchAddress("ITSClusterMCTruth", &clusLabArr); + + auto clSpan = gsl::span(ITSclus->data(), ITSclus->size()); + std::vector hHitMapsVsFrame(nRofs); + + treeITSclus->GetEvent(0); + LOGP(info, "there are {} rofs in this TF", ITSrof->size()); + + // Get patterns + std::vector pattVec; + getClusterPatterns(pattVec, ITSclus, ITSpatt, *mdict); + + for (unsigned int iR{0}; iR < nRofs; iR++) { + LOGP(info, "Processing rof {}", iR + startingROF); + switch(chipId/2) { + case 0: + { + hHitMapsVsFrame[iR] = new TH2D(Form("chip%i_rof%i", chipId, startingROF + iR), Form("chip %i rof %i; ; ; Counts", chipId, startingROF + iR), 13575, -0.5, 13574.5, 2828, -0.5, 2827.5); + break; + } + case 1: + { + hHitMapsVsFrame[iR] = new TH2D(Form("chip%i_rof%i", chipId, startingROF + iR), Form("chip %i rof %i; ; ; Counts", chipId, startingROF + iR), 13575, -0.5, 13574.5, 3770, -0.5, 3769.5); + break; + } + case 2: + { + hHitMapsVsFrame[iR] = new TH2D(Form("chip%i_rof%i", chipId, startingROF + iR), Form("chip %i rof %i; ; ; Counts", chipId, startingROF + iR), 13575, -0.5, 13574.5, 4713, -0.5, 4712.5); + break; + } + default: + { + hHitMapsVsFrame[iR] = new TH2D(Form("chip%i_rof%i", chipId, startingROF + iR), Form("chip %i rof %i; ; ; Counts", chipId, startingROF + iR), 1024, -0.5, 1023.5, 512, -0.5, 511.5); + break; + } + } + + // work on data + const auto& rof = (*ITSrof)[startingROF + iR]; + auto clustersInFrame = rof.getROFData(*ITSclus); + auto patternsInFrame = rof.getROFData(pattVec); + // auto pattIt = ITSpatt->cbegin(); + + for (unsigned int clusInd{0}; clusInd < clustersInFrame.size(); clusInd++) { + const auto& clus = clustersInFrame[clusInd]; + auto sID = clus.getSensorID(); + + if (sID == chipId) { + clus.print(); + + // auto labels = clusLabArr->getLabels(clusInd); + // extract pattern info + auto col = clus.getCol(); + auto row = clus.getRow(); + + std::cout << patternsInFrame[clusInd]; + + std::cout << std::endl; + int ic = 0, ir = 0; + + auto colSpan = patternsInFrame[clusInd].getColumnSpan(); + auto rowSpan = patternsInFrame[clusInd].getRowSpan(); + auto nBits = rowSpan * colSpan; + + for (int i = 2; i < patternsInFrame[clusInd].getUsedBytes() + 2; i++) { + unsigned char tempChar = patternsInFrame[clusInd].getByte(i); + int s = 128; // 0b10000000 + while (s > 0) { + if ((tempChar & s) != 0) // checking active pixels + { + hHitMapsVsFrame[iR]->Fill(col + ic, row + ir); + } + ic++; + s >>= 1; + if ((ir + 1) * ic == nBits) { + break; + } + if (ic == colSpan) { + ic = 0; + ir++; + } + if ((ir + 1) * ic == nBits) { + break; + } + } + } + } + } + } + auto canvas = new TCanvas(Form("chip%d", chipId), Form("chip%d", chipId), nRofs * 1000, 600); + + canvas->Divide(nRofs, 1); + for (unsigned int i{0}; i < nRofs; ++i) { + canvas->cd(i + 1); + gPad->SetGridx(); + gPad->SetGridy(); + hHitMapsVsFrame[i]->Draw("colz"); + } +} + +void getClusterPatterns(std::vector& pattVec, std::vector* ITSclus, std::vector* ITSpatt, o2::itsmft::TopologyDictionary& mdict) +{ + pattVec.reserve(ITSclus->size()); + auto pattIt = ITSpatt->cbegin(); + + for (unsigned int iClus{0}; iClus < ITSclus->size(); ++iClus) { + auto& clus = (*ITSclus)[iClus]; + + auto pattID = clus.getPatternID(); + int npix; + o2::itsmft::ClusterPattern patt; + + if (pattID == o2::its3::CompCluster::InvalidPatternID || mdict.isGroup(pattID)) { + patt.acquirePattern(pattIt); + npix = patt.getNPixels(); + } else { + + npix = mdict.getNpixels(pattID); + patt = mdict.getPattern(pattID); + } + + pattVec.push_back(patt); + } +} diff --git a/Detectors/Upgrades/IT3/reconstruction/CMakeLists.txt b/Detectors/Upgrades/ITS3/reconstruction/CMakeLists.txt similarity index 93% rename from Detectors/Upgrades/IT3/reconstruction/CMakeLists.txt rename to Detectors/Upgrades/ITS3/reconstruction/CMakeLists.txt index e32106dfe261f..38c25425c7ba3 100644 --- a/Detectors/Upgrades/IT3/reconstruction/CMakeLists.txt +++ b/Detectors/Upgrades/ITS3/reconstruction/CMakeLists.txt @@ -13,13 +13,13 @@ o2_add_library(ITS3Reconstruction TARGETVARNAME targetName SOURCES src/Clusterer.cxx src/TopologyDictionary.cxx - + # src/FastMultEst.cxx PUBLIC_LINK_LIBRARIES O2::ITSMFTBase O2::ITSMFTReconstruction O2::ITS3Base O2::CommonDataFormat O2::DetectorsRaw - O2::SimulationDataFormat + O2::SimulationDataFormat O2::DataFormatsITSMFT O2::DataFormatsITS3 O2::DPLUtils @@ -36,4 +36,3 @@ if (OpenMP_CXX_FOUND) target_compile_definitions(${targetName} PRIVATE WITH_OPENMP) target_link_libraries(${targetName} PRIVATE OpenMP::OpenMP_CXX) endif() - diff --git a/Detectors/Upgrades/IT3/reconstruction/include/ITS3Reconstruction/Clusterer.h b/Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/Clusterer.h similarity index 63% rename from Detectors/Upgrades/IT3/reconstruction/include/ITS3Reconstruction/Clusterer.h rename to Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/Clusterer.h index 9be735ec39e4b..0ce3a7e9ac6a1 100644 --- a/Detectors/Upgrades/IT3/reconstruction/include/ITS3Reconstruction/Clusterer.h +++ b/Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/Clusterer.h @@ -16,6 +16,10 @@ #define _PERFORM_TIMING_ +// uncomment this to not allow diagonal clusters, e.g. like |* | +// | *| +#define _ALLOW_DIAGONAL_ALPIDE_CLUSTERS_ + #include #include #include @@ -51,23 +55,60 @@ class MCTruthContainer; namespace its3 { -using CompClusCont = std::vector; +using CompClusCont = std::vector; using PatternCont = std::vector; using ROFRecCont = std::vector; -//template // container types (PMR or std::vectors) +// template // container types (PMR or std::vectors) class Clusterer { using PixelReader = o2::itsmft::PixelReader; using PixelData = o2::itsmft::PixelData; using ChipPixelData = o2::itsmft::ChipPixelData; + using CompCluster = o2::its3::CompCluster; + using CompClusterExt = o2::its3::CompClusterExt; using Label = o2::MCCompLabel; using MCTruth = o2::dataformats::MCTruthContainer; using ConstMCTruth = o2::dataformats::ConstMCTruthContainerView; public: static constexpr int MaxLabels = 10; + static constexpr int MaxHugeClusWarn = 5; // max number of warnings for HugeCluster + + struct BBox { + uint16_t chipID = 0xffff; + uint16_t rowMin = 0xffff; + uint16_t colMin = 0xffff; + uint16_t rowMax = 0; + uint16_t colMax = 0; + BBox(uint16_t c) : chipID(c) {} + bool isInside(uint16_t row, uint16_t col) const { return row >= rowMin && row <= rowMax && col >= colMin && col <= colMax; } + auto rowSpan() const { return rowMax - rowMin + 1; } + auto colSpan() const { return colMax - colMin + 1; } + bool isAcceptableSize() const { return colMax - colMin < o2::itsmft::ClusterPattern::MaxColSpan && rowMax - rowMin < o2::itsmft::ClusterPattern::MaxRowSpan; } + void clear() + { + rowMin = colMin = 0xffff; + rowMax = colMax = 0; + } + void adjust(uint16_t row, uint16_t col) + { + if (row < rowMin) { + rowMin = row; + } + if (row > rowMax) { + rowMax = row; + } + if (col < colMin) { + colMin = col; + } + if (col > colMax) { + colMax = col; + } + } + }; + //========================================================= /// methods and transient data used within a thread struct ThreadStat { @@ -81,7 +122,7 @@ class Clusterer }; struct ClustererThread { - + int id = -1; Clusterer* parent = nullptr; // parent clusterer // buffers for entries in preClusterIndices in 2 columns, to avoid boundary checks, we reserve // extra elements in the beginning and the end @@ -132,13 +173,8 @@ class Clusterer curr[row] = lastIndex; // store index of the new precluster in the current column buffer } - void streamCluster(const std::vector& pixbuf, uint16_t rowMin, uint16_t rowSpanW, uint16_t colMin, uint16_t colSpanW, - uint16_t chipID, - CompClusCont* compClusPtr, PatternCont* patternsPtr, - MCTruth* labelsClusPtr, int nlab, bool isHuge = false); - void fetchMCLabels(int digID, const ConstMCTruth* labelsDig, int& nfilled); - void initChip(const ChipPixelData* curChipData, uint32_t first, int chipID); + void initChip(const ChipPixelData* curChipData, uint32_t first); void updateChip(const ChipPixelData* curChipData, uint32_t ip); void finishChip(ChipPixelData* curChipData, CompClusCont* compClus, PatternCont* patterns, const ConstMCTruth* labelsDig, MCTruth* labelsClus); @@ -147,7 +183,6 @@ class Clusterer void process(uint16_t chip, uint16_t nChips, CompClusCont* compClusPtr, PatternCont* patternsPtr, const ConstMCTruth* labelsDigPtr, MCTruth* labelsClPtr, const itsmft::ROFRecord& rofPtr); - ClustererThread(Clusterer* par = nullptr) : parent(par) {} ~ClustererThread() { if (column1) { @@ -157,6 +192,11 @@ class Clusterer delete[] column2; } } + ClustererThread(Clusterer* par = nullptr, int _id = -1) : parent(par), id(_id), curr(column2 + 1), prev(column1 + 1) + { + // std::fill(std::begin(column1), std::end(column1), -1); + // std::fill(std::begin(column2), std::end(column2), -1); + } }; //========================================================= @@ -168,6 +208,10 @@ class Clusterer void process(int nThreads, PixelReader& r, CompClusCont* compClus, PatternCont* patterns, ROFRecCont* vecROFRec, MCTruth* labelsCl = nullptr); + template + static void streamCluster(const std::vector& pixbuf, const std::array* lblBuff, const BBox& bbox, const itsmft::LookUp& pattIdConverter, + VCLUS* compClusPtr, VPAT* patternsPtr, MCTruth* labelsClusPtr, int nlab, bool isHuge = false); + bool isContinuousReadOut() const { return mContinuousReadout; } void setContinuousReadOut(bool v) { mContinuousReadout = v; } @@ -177,6 +221,12 @@ class Clusterer int getMaxRowColDiffToMask() const { return mMaxRowColDiffToMask; } void setMaxRowColDiffToMask(int v) { mMaxRowColDiffToMask = v; } + int getMaxROFDepthToSquash() const { return mSquashingDepth; } + void setMaxROFDepthToSquash(int v) { mSquashingDepth = v; } + + int getMaxBCSeparationToSquash() const { return mMaxBCSeparationToSquash; } + void setMaxBCSeparationToSquash(int n) { mMaxBCSeparationToSquash = n; } + void print() const; void clear(); @@ -188,28 +238,12 @@ class Clusterer ///< load the dictionary of cluster topologies void loadDictionary(const std::string& fileName) { mPattIdConverter.loadDictionary(fileName); } + void setDictionary(const itsmft::TopologyDictionary* dict) { mPattIdConverter.setDictionary(dict); } TStopwatch& getTimer() { return mTimer; } // cannot be const TStopwatch& getTimerMerge() { return mTimerMerge; } // cannot be const private: - ///< recalculate min max row and column of the cluster accounting for the position of pix - static void adjustBoundingBox(uint16_t row, uint16_t col, uint16_t& rMin, uint16_t& rMax, uint16_t& cMin, uint16_t& cMax) - { - if (row < rMin) { - rMin = row; - } - if (row > rMax) { - rMax = row; - } - if (col < cMin) { - cMin = col; - } - if (col > cMax) { - cMax = col; - } - } - void flushClusters(CompClusCont* compClus, MCTruth* labels); // clusterization options @@ -218,6 +252,11 @@ class Clusterer ///< mask continuosly fired pixels in frames separated by less than this amount of BCs (fired from hit in prev. ROF) int mMaxBCSeparationToMask = 6000. / o2::constants::lhc::LHCBunchSpacingNS + 10; int mMaxRowColDiffToMask = 0; ///< provide their difference in col/row is <= than this + int mNHugeClus = 0; ///< number of encountered huge clusters + + ///< Squashing options + int mSquashingDepth = 0; ///< squashing is applied to next N rofs + int mMaxBCSeparationToSquash = 6000. / o2::constants::lhc::LHCBunchSpacingNS + 10; std::vector> mThreads; // buffers for threads std::vector mChips; // currently processed ROF's chips data @@ -230,6 +269,48 @@ class Clusterer TStopwatch mTimerMerge; }; +template +void Clusterer::streamCluster(const std::vector& pixbuf, const std::array* lblBuff, const Clusterer::BBox& bbox, const itsmft::LookUp& pattIdConverter, + VCLUS* compClusPtr, VPAT* patternsPtr, MCTruth* labelsClusPtr, int nlab, bool isHuge) +{ + if (labelsClusPtr && lblBuff) { // MC labels were requested + auto cnt = compClusPtr->size(); + for (int i = nlab; i--;) { + labelsClusPtr->addElement(cnt, (*lblBuff)[i]); + } + } + auto colSpanW = bbox.colSpan(); + auto rowSpanW = bbox.rowSpan(); + // add to compact clusters, which must be always filled + std::array patt{}; + for (const auto& pix : pixbuf) { + uint32_t ir = pix.getRowDirect() - bbox.rowMin, ic = pix.getCol() - bbox.colMin; + int nbits = ir * colSpanW + ic; + patt[nbits >> 3] |= (0x1 << (7 - (nbits % 8))); + } + uint16_t pattID = (isHuge || pattIdConverter.size() == 0) ? CompCluster::InvalidPatternID : pattIdConverter.findGroupID(rowSpanW, colSpanW, patt.data()); + uint16_t row = bbox.rowMin, col = bbox.colMin; + if (pattID == CompCluster::InvalidPatternID || pattIdConverter.isGroup(pattID)) { + if (pattID != CompCluster::InvalidPatternID) { + // For groupped topologies, the reference pixel is the COG pixel + float xCOG = 0., zCOG = 0.; + itsmft::ClusterPattern::getCOG(rowSpanW, colSpanW, patt.data(), xCOG, zCOG); + row += round(xCOG); + col += round(zCOG); + } + if (patternsPtr) { + patternsPtr->emplace_back((unsigned char)rowSpanW); + patternsPtr->emplace_back((unsigned char)colSpanW); + int nBytes = rowSpanW * colSpanW / 8; + if (((rowSpanW * colSpanW) % 8) != 0) { + nBytes++; + } + patternsPtr->insert(patternsPtr->end(), std::begin(patt), std::begin(patt) + nBytes); + } + } + compClusPtr->emplace_back(row, col, pattID, bbox.chipID); +} + } // namespace its3 } // namespace o2 #endif /* ALICEO2_ITS_CLUSTERER_H */ \ No newline at end of file diff --git a/Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/FastMultEst.h b/Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/FastMultEst.h new file mode 100644 index 0000000000000..e9da619c0efbf --- /dev/null +++ b/Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/FastMultEst.h @@ -0,0 +1,66 @@ +// 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. + +#ifndef ALICEO2_ITS3_FASTMULTEST_ +#define ALICEO2_ITS3_FASTMULTEST_ + +#include "ITSMFTReconstruction/ChipMappingITS.h" +#include "DataFormatsITSMFT/ROFRecord.h" +#include "DataFormatsITS3/CompCluster.h" +#include +#include "ITS3Reconstruction/FastMultEstConfig.h" +#include +#include + +namespace o2 +{ +namespace its3 +{ + +struct FastMultEst { + + static constexpr int NLayers = o2::itsmft::ChipMappingITS::NLayers; + + float mult = 0.; /// estimated signal clusters multipliciy at reference (1st?) layer + float noisePerChip = 0.; /// estimated or imposed noise per chip + float cov[3] = {0.}; /// covariance matrix of estimation + float chi2 = 0.; /// chi2 + int nLayersUsed = 0; /// number of layers actually used + uint32_t lastRandomSeed = 0; /// state of the gRandom before + + std::array nClPerLayer{0}; // measured N Cl per layer selectROFs + FastMultEst(); + + static uint32_t getCurrentRandomSeed(); + int selectROFs(const gsl::span rofs, const gsl::span clus, + const gsl::span trig, std::vector& sel); + + void fillNClPerLayer(const gsl::span& clusters); + float process(const std::array ncl) + { + return FastMultEstConfig::Instance().imposeNoisePerChip > 0 ? processNoiseImposed(ncl) : processNoiseFree(ncl); + } + float processNoiseFree(const std::array ncl); + float processNoiseImposed(const std::array ncl); + float process(const gsl::span& clusters) + { + fillNClPerLayer(clusters); + return process(nClPerLayer); + } + static bool sSeedSet; + + ClassDefNV(FastMultEst, 1); +}; + +} // namespace its3 +} // namespace o2 + +#endif diff --git a/Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/FastMultEstConfig.h b/Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/FastMultEstConfig.h new file mode 100644 index 0000000000000..1857176b19f1f --- /dev/null +++ b/Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/FastMultEstConfig.h @@ -0,0 +1,57 @@ +// 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 FastMultEstConfig.h +/// \brief Configuration parameters for ITS fast multiplicity estimator +/// \author ruben.shahoyan@cern.ch + +#ifndef ALICEO2_ITS_FASTMULTESTCONF_H_ +#define ALICEO2_ITS_FASTMULTESTCONF_H_ + +#include "CommonUtils/ConfigurableParam.h" +#include "CommonUtils/ConfigurableParamHelper.h" +#include "ITSMFTReconstruction/ChipMappingITS.h" + +namespace o2 +{ +namespace its +{ +struct FastMultEstConfig : public o2::conf::ConfigurableParamHelper { + static constexpr int NLayers = 7; // FIXME + + /// acceptance correction per layer (relative to 1st one) + float accCorr[NLayers] = {1.f, 0.895, 0.825, 0.803, 0.720, 0.962, 0.911}; + int firstLayer = 3; /// 1st layer to account + int lastLayer = 6; /// last layer to account + float imposeNoisePerChip = 1.e-7 * 1024 * 512; // assumed noise, free parameter if<0 + + // cuts to reject to low or too high mult events + float cutMultClusLow = 0; /// reject ROF with estimated cluster mult. below this value (no cut if <0) + float cutMultClusHigh = -1; /// reject ROF with estimated cluster mult. above this value (no cut if <0) + float cutMultVtxLow = -1; /// reject seed vertex if its multiplicity below this value (no cut if <0) + float cutMultVtxHigh = -1; /// reject seed vertex if its multiplicity above this value (no cut if <0) + float cutRandomFraction = -1.; /// apply random cut rejecting requested fraction + int randomSeed = 0; /// 0 - do not seet seed, >0 : set as is, <0 : use current time + bool preferTriggered = true; /// prefer ROFs with highest number of physics triggers + + bool isMultCutRequested() const { return cutMultClusLow >= 0.f && cutMultClusHigh > 0.f; }; + bool isVtxMultCutRequested() const { return cutMultVtxLow >= 0.f && cutMultVtxHigh > 0.f; }; + bool isPassingRandomRejection() const; + bool isPassingMultCut(float mult) const { return mult >= cutMultClusLow && (mult <= cutMultClusHigh || cutMultClusHigh <= 0.f); } + bool isPassingVtxMultCut(int mult) const { return mult >= cutMultVtxLow && (mult <= cutMultVtxHigh || cutMultVtxHigh <= 0.f); } + + O2ParamDef(FastMultEstConfig, "fastMultConfig"); +}; + +} // namespace its +} // namespace o2 + +#endif diff --git a/Detectors/Upgrades/IT3/reconstruction/include/ITS3Reconstruction/TopologyDictionary.h b/Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/TopologyDictionary.h similarity index 92% rename from Detectors/Upgrades/IT3/reconstruction/include/ITS3Reconstruction/TopologyDictionary.h rename to Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/TopologyDictionary.h index 03b84e5336c29..05aba4700985a 100644 --- a/Detectors/Upgrades/IT3/reconstruction/include/ITS3Reconstruction/TopologyDictionary.h +++ b/Detectors/Upgrades/ITS3/reconstruction/include/ITS3Reconstruction/TopologyDictionary.h @@ -29,9 +29,9 @@ class TopologyDictionary : public itsmft::TopologyDictionary public: TopologyDictionary(itsmft::TopologyDictionary top) : itsmft::TopologyDictionary{top} {} - ///Returns the local position of a compact cluster + /// Returns the local position of a compact cluster math_utils::Point3D getClusterCoordinates(int detID, const its3::CompCluster& cl) const; - ///Returns the local position of a compact cluster + /// Returns the local position of a compact cluster static math_utils::Point3D getClusterCoordinates(int detID, const its3::CompCluster& cl, const itsmft::ClusterPattern& patt, bool isGroup = true); }; } // namespace its3 diff --git a/Detectors/Upgrades/IT3/reconstruction/src/Clusterer.cxx b/Detectors/Upgrades/ITS3/reconstruction/src/Clusterer.cxx similarity index 76% rename from Detectors/Upgrades/IT3/reconstruction/src/Clusterer.cxx rename to Detectors/Upgrades/ITS3/reconstruction/src/Clusterer.cxx index 96822b01c178b..e643dc5182ec7 100644 --- a/Detectors/Upgrades/IT3/reconstruction/src/Clusterer.cxx +++ b/Detectors/Upgrades/ITS3/reconstruction/src/Clusterer.cxx @@ -23,7 +23,6 @@ #ifdef WITH_OPENMP #include #endif - using namespace o2::itsmft; namespace o2::its3 @@ -39,6 +38,7 @@ void Clusterer::process(int nThreads, PixelReader& reader, CompClusCont* compClu nThreads = 1; } auto autoDecode = reader.getDecodeNextAuto(); + int rofcount{0}; do { if (autoDecode) { reader.setDecodeNextAuto(false); // internally do not autodecode @@ -46,6 +46,9 @@ void Clusterer::process(int nThreads, PixelReader& reader, CompClusCont* compClu break; // on the fly decoding was requested, but there were no data left } } + if (reader.getInteractionRecord().isDummy()) { + continue; // No IR info was found + } // pre-fetch all non-empty chips of current ROF ChipPixelData* curChipData = nullptr; mFiredChipsPtr.clear(); @@ -55,7 +58,7 @@ void Clusterer::process(int nThreads, PixelReader& reader, CompClusCont* compClu nPix += curChipData->getData().size(); } - auto& rof = vecROFRec->emplace_back(reader.getInteractionRecord(), 0, compClus->size(), 0); // create new ROF + auto& rof = vecROFRec->emplace_back(reader.getInteractionRecord(), vecROFRec->size(), compClus->size(), 0); // create new ROF uint16_t nFired = mFiredChipsPtr.size(); if (!nFired) { @@ -67,22 +70,20 @@ void Clusterer::process(int nThreads, PixelReader& reader, CompClusCont* compClu if (nFired < nThreads) { nThreads = nFired; } -#ifdef WITH_OPENMP - omp_set_num_threads(nThreads); -#else +#ifndef WITH_OPENMP nThreads = 1; #endif - uint16_t chipStep = nThreads > 1 ? (nThreads == 2 ? 20 : (nThreads < 5 ? 5 : 1)) : nFired; + uint16_t chipStep = nThreads > 1 ? (nThreads == 2 ? 20 : 10) : nFired; int dynGrp = std::min(4, std::max(1, nThreads / 2)); if (nThreads > mThreads.size()) { int oldSz = mThreads.size(); mThreads.resize(nThreads); for (int i = oldSz; i < nThreads; i++) { - mThreads[i] = std::make_unique(this); + mThreads[i] = std::make_unique(this, i); } } #ifdef WITH_OPENMP -#pragma omp parallel for schedule(dynamic, dynGrp) +#pragma omp parallel for schedule(dynamic, dynGrp) num_threads(nThreads) //>> start of MT region for (uint16_t ic = 0; ic < nFired; ic += chipStep) { auto ith = omp_get_thread_num(); @@ -108,6 +109,7 @@ void Clusterer::process(int nThreads, PixelReader& reader, CompClusCont* compClu size_t nClTot = 0, nPattTot = 0; int chid = 0, thrStatIdx[nThreads]; for (int ith = 0; ith < nThreads; ith++) { + std::sort(mThreads[ith]->stats.begin(), mThreads[ith]->stats.end(), [](const ThreadStat& a, const ThreadStat& b) { return a.firstChip < b.firstChip; }); thrStatIdx[ith] = 0; nClTot += mThreads[ith]->compClusters.size(); nPattTot += mThreads[ith]->patterns.size(); @@ -162,10 +164,9 @@ void Clusterer::process(int nThreads, PixelReader& reader, CompClusCont* compClu void Clusterer::ClustererThread::process(uint16_t chip, uint16_t nChips, CompClusCont* compClusPtr, PatternCont* patternsPtr, const ConstMCTruth* labelsDigPtr, MCTruth* labelsClPtr, const ROFRecord& rofPtr) { - if (stats.empty() || stats.back().firstChip + stats.back().nChips < chip) { // there is a jump, register new block + if (stats.empty() || stats.back().firstChip + stats.back().nChips != chip) { // there is a jump, register new block stats.emplace_back(ThreadStat{chip, 0, uint32_t(compClusPtr->size()), patternsPtr ? uint32_t(patternsPtr->size()) : 0, 0, 0}); } - for (int ic = 0; ic < nChips; ic++) { auto* curChipData = parent->mFiredChipsPtr[chip + ic]; auto chipID = curChipData->getChipID(); @@ -175,6 +176,7 @@ void Clusterer::ClustererThread::process(uint16_t chip, uint16_t nChips, CompClu parent->mMaxRowColDiffToMask ? curChipData->maskFiredInSample(parent->mChipsOld[chipID], parent->mMaxRowColDiffToMask) : curChipData->maskFiredInSample(parent->mChipsOld[chipID]); } } + auto nclus0 = compClusPtr->size(); auto validPixID = curChipData->getFirstUnmasked(); auto npix = curChipData->getData().size(); if (validPixID < npix) { // chip data may have all of its pixels masked! @@ -182,7 +184,7 @@ void Clusterer::ClustererThread::process(uint16_t chip, uint16_t nChips, CompClu if (validPixID == npix) { // special case of a single pixel fired on the chip finishChipSingleHitFast(valp, curChipData, compClusPtr, patternsPtr, labelsDigPtr, labelsClPtr); } else { - initChip(curChipData, valp, chipID); + initChip(curChipData, valp); for (; validPixID < npix; validPixID++) { if (!curChipData->getData()[validPixID].isMasked()) { updateChip(curChipData, validPixID); @@ -205,15 +207,13 @@ void Clusterer::ClustererThread::process(uint16_t chip, uint16_t nChips, CompClu void Clusterer::ClustererThread::finishChip(ChipPixelData* curChipData, CompClusCont* compClusPtr, PatternCont* patternsPtr, const ConstMCTruth* labelsDigPtr, MCTruth* labelsClusPtr) { - auto clustersCount = compClusPtr->size(); const auto& pixData = curChipData->getData(); for (int i1 = 0; i1 < preClusterHeads.size(); ++i1) { auto ci = preClusterIndices[i1]; if (ci < 0) { continue; } - uint16_t rowMax = 0, rowMin = 65535; - uint16_t colMax = 0, colMin = 65535; + BBox bbox(curChipData->getChipID()); int nlab = 0; int next = preClusterHeads[i1]; pixArrBuff.clear(); @@ -221,9 +221,13 @@ void Clusterer::ClustererThread::finishChip(ChipPixelData* curChipData, CompClus const auto& pixEntry = pixels[next]; const auto pix = pixData[pixEntry.second]; pixArrBuff.push_back(pix); // needed for cluster topology - adjustBoundingBox(pix.getRowDirect(), pix.getCol(), rowMin, rowMax, colMin, colMax); - if (labelsClusPtr) { // the MCtruth for this pixel is at curChipData->startID+pixEntry.second - fetchMCLabels(pixEntry.second + curChipData->getStartID(), labelsDigPtr, nlab); + bbox.adjust(pix.getRowDirect(), pix.getCol()); + if (labelsClusPtr) { + if (parent->mSquashingDepth) { // the MCtruth for this pixel is stored in chip data: due to squashing we lose contiguity + fetchMCLabels(curChipData->getOrderedPixId(pixEntry.second), labelsDigPtr, nlab); + } else { // the MCtruth for this pixel is at curChipData->startID+pixEntry.second + fetchMCLabels(pixEntry.second + curChipData->getStartID(), labelsDigPtr, nlab); + } } next = pixEntry.first; } @@ -237,107 +241,60 @@ void Clusterer::ClustererThread::finishChip(ChipPixelData* curChipData, CompClus const auto& pixEntry = pixels[next]; const auto pix = pixData[pixEntry.second]; // PixelData pixArrBuff.push_back(pix); // needed for cluster topology - adjustBoundingBox(pix.getRowDirect(), pix.getCol(), rowMin, rowMax, colMin, colMax); - if (labelsClusPtr) { // the MCtruth for this pixel is at curChipData->startID+pixEntry.second - fetchMCLabels(pixEntry.second + curChipData->getStartID(), labelsDigPtr, nlab); + bbox.adjust(pix.getRowDirect(), pix.getCol()); + if (labelsClusPtr) { + if (parent->mSquashingDepth) { // the MCtruth for this pixel is stored in chip data: due to squashing we lose contiguity + fetchMCLabels(curChipData->getOrderedPixId(pixEntry.second), labelsDigPtr, nlab); + } else { // the MCtruth for this pixel is at curChipData->startID+pixEntry.second + fetchMCLabels(pixEntry.second + curChipData->getStartID(), labelsDigPtr, nlab); + } } next = pixEntry.first; } preClusterIndices[i2] = -1; } - - auto chipID = curChipData->getChipID(); - uint16_t colSpan = (colMax - colMin + 1); - uint16_t rowSpan = (rowMax - rowMin + 1); - if (colSpan <= o2::itsmft::ClusterPattern::MaxColSpan && - rowSpan <= o2::itsmft::ClusterPattern::MaxRowSpan) { - streamCluster(pixArrBuff, rowMin, rowSpan, colMin, colSpan, chipID, - compClusPtr, patternsPtr, labelsClusPtr, nlab); + if (bbox.isAcceptableSize()) { + parent->streamCluster(pixArrBuff, &labelsBuff, bbox, parent->mPattIdConverter, compClusPtr, patternsPtr, labelsClusPtr, nlab); } else { - LOG(alarm) << "Splitting a huge cluster ! ChipID: " << chipID; - - colSpan %= o2::itsmft::ClusterPattern::MaxColSpan; - if (colSpan == 0) { - colSpan = o2::itsmft::ClusterPattern::MaxColSpan; - } - - rowSpan %= o2::itsmft::ClusterPattern::MaxRowSpan; - if (rowSpan == 0) { - rowSpan = o2::itsmft::ClusterPattern::MaxRowSpan; + auto warnLeft = MaxHugeClusWarn - parent->mNHugeClus; + if (warnLeft > 0) { + LOGP(warn, "Splitting a huge cluster: chipID {}, rows {}:{} cols {}:{}{}", bbox.chipID, bbox.rowMin, bbox.rowMax, bbox.colMin, bbox.colMax, + warnLeft == 1 ? " (Further warnings will be muted)" : ""); +#ifdef WITH_OPENMP +#pragma omp critical +#endif + { + parent->mNHugeClus++; + } } - + BBox bboxT(bbox); // truncated box + std::vector pixbuf; do { - uint16_t r = rowMin, rsp = rowSpan; - - do { - // Select a subset of pixels fitting the reduced bounding box - std::vector pixbuf; - auto colMax = colMin + colSpan, rowMax = r + rsp; + bboxT.rowMin = bbox.rowMin; + bboxT.colMax = std::min(bbox.colMax, uint16_t(bboxT.colMin + o2::itsmft::ClusterPattern::MaxColSpan - 1)); + do { // Select a subset of pixels fitting the reduced bounding box + bboxT.rowMax = std::min(bbox.rowMax, uint16_t(bboxT.rowMin + o2::itsmft::ClusterPattern::MaxRowSpan - 1)); for (const auto& pix : pixArrBuff) { - if (pix.getRowDirect() >= r && pix.getRowDirect() < rowMax && - pix.getCol() >= colMin && pix.getCol() < colMax) { + if (bboxT.isInside(pix.getRowDirect(), pix.getCol())) { pixbuf.push_back(pix); } } - // Stream a piece of cluster only if the reduced bounding box is not empty - if (!pixbuf.empty()) { - streamCluster(pixbuf, r, rsp, colMin, colSpan, chipID, - compClusPtr, patternsPtr, labelsClusPtr, nlab, true); + if (!pixbuf.empty()) { // Stream a piece of cluster only if the reduced bounding box is not empty + parent->streamCluster(pixbuf, &labelsBuff, bboxT, parent->mPattIdConverter, compClusPtr, patternsPtr, labelsClusPtr, nlab, true); + pixbuf.clear(); } - r += rsp; - rsp = o2::itsmft::ClusterPattern::MaxRowSpan; - } while (r < rowMax); - - colMin += colSpan; - colSpan = o2::itsmft::ClusterPattern::MaxColSpan; - } while (colMin < colMax); + bboxT.rowMin = bboxT.rowMax + 1; + } while (bboxT.rowMin < bbox.rowMax); + bboxT.colMin = bboxT.colMax + 1; + } while (bboxT.colMin < bbox.colMax); } } } -void Clusterer::ClustererThread::streamCluster(const std::vector& pixbuf, uint16_t rowMin, uint16_t rowSpanW, uint16_t colMin, uint16_t colSpanW, uint16_t chipID, CompClusCont* compClusPtr, PatternCont* patternsPtr, MCTruth* labelsClusPtr, int nlab, bool isHuge) -{ - if (labelsClusPtr) { // MC labels were requested - auto cnt = compClusPtr->size(); - for (int i = nlab; i--;) { - labelsClusPtr->addElement(cnt, labelsBuff[i]); - } - } - - // add to compact clusters, which must be always filled - unsigned char patt[ClusterPattern::MaxPatternBytes] = {0}; // RSTODO FIX pattern filling - for (const auto& pix : pixbuf) { - unsigned short ir = pix.getRowDirect() - rowMin, ic = pix.getCol() - colMin; - int nbits = ir * colSpanW + ic; - patt[nbits >> 3] |= (0x1 << (7 - (nbits % 8))); - } - uint16_t pattID = (isHuge || parent->mPattIdConverter.size() == 0) ? CompCluster::InvalidPatternID : parent->mPattIdConverter.findGroupID(rowSpanW, colSpanW, patt); - if (pattID == CompCluster::InvalidPatternID || parent->mPattIdConverter.isGroup(pattID)) { - if (pattID != CompCluster::InvalidPatternID) { - //For groupped topologies, the reference pixel is the COG pixel - float xCOG = 0., zCOG = 0.; - ClusterPattern::getCOG(rowSpanW, colSpanW, patt, xCOG, zCOG); - rowMin += round(xCOG); - colMin += round(zCOG); - } - if (patternsPtr) { - patternsPtr->emplace_back((unsigned char)rowSpanW); - patternsPtr->emplace_back((unsigned char)colSpanW); - int nBytes = rowSpanW * colSpanW / 8; - if (((rowSpanW * colSpanW) % 8) != 0) { - nBytes++; - } - patternsPtr->insert(patternsPtr->end(), std::begin(patt), std::begin(patt) + nBytes); - } - } - compClusPtr->emplace_back(rowMin, colMin, pattID, chipID); -} - //__________________________________________________ void Clusterer::ClustererThread::finishChipSingleHitFast(uint32_t hit, ChipPixelData* curChipData, CompClusCont* compClusPtr, PatternCont* patternsPtr, const ConstMCTruth* labelsDigPtr, MCTruth* labelsClusPtr) { - auto clustersCount = compClusPtr->size(); auto pix = curChipData->getData()[hit]; uint16_t row = pix.getRowDirect(), col = pix.getCol(); @@ -373,13 +330,13 @@ Clusterer::Clusterer() : mPattIdConverter() } //__________________________________________________ -void Clusterer::ClustererThread::initChip(const ChipPixelData* curChipData, uint32_t first, int chip) +void Clusterer::ClustererThread::initChip(const ChipPixelData* curChipData, uint32_t first) { // init chip with the 1st unmasked pixel (entry "from" in the mChipData) size = itsmft::SegmentationAlpide::NRows + 2; - if (chip < SegmentationSuperAlpide::NLayers) { - SegmentationSuperAlpide seg(chip); - size = seg.NRows + 2; + if (curChipData->getChipID() < 6) { // TODO Fix for mutable layouts + SegmentationSuperAlpide seg(curChipData->getChipID() / 2); + size = seg.mNRows + 2; } if (column1) { delete[] column1; @@ -437,7 +394,11 @@ void Clusterer::ClustererThread::updateChip(const ChipPixelData* curChipData, ui return; } } else { +#ifdef _ALLOW_DIAGONAL_ALPIDE_CLUSTERS_ int neighbours[]{curr[row - 1], prev[row], prev[row + 1], prev[row - 1]}; +#else + int neighbours[]{curr[row - 1], prev[row]}; +#endif for (auto pci : neighbours) { if (pci < 0) { continue; @@ -499,8 +460,10 @@ void Clusterer::clear() void Clusterer::print() const { // print settings + LOGP(info, "Clusterizer squashes overflow pixels separated by {} BC and <= {} in row/col seeking down to {} neighbour ROFs", mMaxBCSeparationToSquash, mMaxRowColDiffToMask, mSquashingDepth); LOG(info) << "Clusterizer masks overflow pixels separated by < " << mMaxBCSeparationToMask << " BC and <= " << mMaxRowColDiffToMask << " in row/col"; + #ifdef _PERFORM_TIMING_ auto& tmr = const_cast(mTimer); // ugly but this is what root does internally auto& tmrm = const_cast(mTimerMerge); @@ -511,5 +474,4 @@ void Clusterer::print() const #endif } - } // namespace o2::its3 \ No newline at end of file diff --git a/Detectors/Upgrades/ITS3/reconstruction/src/FastMultEst.cxx b/Detectors/Upgrades/ITS3/reconstruction/src/FastMultEst.cxx new file mode 100644 index 0000000000000..fa2ce319328b5 --- /dev/null +++ b/Detectors/Upgrades/ITS3/reconstruction/src/FastMultEst.cxx @@ -0,0 +1,207 @@ +// 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. + +#include "ITS3Reconstruction/FastMultEst.h" +#include "ITSMFTBase/DPLAlpideParam.h" +#include "Framework/Logger.h" +#include +#include +#include + +using namespace o2::its; + +bool FastMultEst::sSeedSet = false; + +///______________________________________________________ +FastMultEst::FastMultEst() +{ + if (!sSeedSet && FastMultEstConfig::Instance().cutRandomFraction > 0.f) { + sSeedSet = true; + if (FastMultEstConfig::Instance().randomSeed > 0) { + gRandom->SetSeed(FastMultEstConfig::Instance().randomSeed); + } else if (FastMultEstConfig::Instance().randomSeed < 0) { + gRandom->SetSeed(std::time(nullptr) % 0xffff); + } + } +} + +///______________________________________________________ +/// find multiplicity for given set of clusters +void FastMultEst::fillNClPerLayer(const gsl::span& clusters) +{ + int lr = FastMultEst::NLayers - 1, nchAcc = o2::itsmft::ChipMappingITS::getNChips() - o2::itsmft::ChipMappingITS::getNChipsPerLr(lr); + std::memset(&nClPerLayer[0], 0, sizeof(int) * FastMultEst::NLayers); + for (int i = clusters.size(); i--;) { // profit from clusters being ordered in chip increasing order + while (clusters[i].getSensorID() < nchAcc) { + assert(lr >= 0); + nchAcc -= o2::itsmft::ChipMappingITS::getNChipsPerLr(--lr); + } + nClPerLayer[lr]++; + } +} + +///______________________________________________________ +/// find multiplicity for given number of clusters per layer +float FastMultEst::processNoiseFree(const std::array ncl) +{ + // we assume that on the used layers the observed number of clusters is defined by the + // the noise ~ nu * Nchips and contribution from the signal tracks Ntr*mAccCorr + const auto& conf = FastMultEstConfig::Instance(); + + float mat[3] = {0}, b[2] = {0}; + nLayersUsed = 0; + for (int il = conf.firstLayer; il <= conf.lastLayer; il++) { + if (ncl[il] > 0) { + int nch = o2::itsmft::ChipMappingITS::getNChipsPerLr(il); + float err2i = 1. / ncl[il]; + float m2n = nch * err2i; + mat[0] += err2i * conf.accCorr[il] * conf.accCorr[il]; + mat[2] += nch * m2n; + mat[1] += conf.accCorr[il] * m2n; // non-diagonal element + b[0] += conf.accCorr[il]; + b[1] += nch; + nLayersUsed++; + } + } + mult = noisePerChip = chi2 = -1; + float det = mat[0] * mat[2] - mat[1] * mat[1]; + if (nLayersUsed < 2 || std::abs(det) < 1e-15) { + return -1; + } + float detI = 1. / det; + mult = detI * (b[0] * mat[2] - b[1] * mat[1]); + noisePerChip = detI * (b[1] * mat[0] - b[0] * mat[1]); + cov[0] = mat[2] * detI; + cov[2] = mat[0] * detI; + cov[1] = -mat[1] * detI; + chi2 = 0.; + for (int il = conf.firstLayer; il <= conf.lastLayer; il++) { + if (ncl[il] > 0) { + int nch = o2::itsmft::ChipMappingITS::getNChipsPerLr(il); + float diff = mult * conf.accCorr[il] + nch * noisePerChip - ncl[il]; + chi2 += diff * diff / ncl[il]; + } + } + chi2 = nLayersUsed > 2 ? chi2 / (nLayersUsed - 2) : 0.; + return mult > 0 ? mult : 0; +} + +///______________________________________________________ +/// find multiplicity for given number of clusters per layer with mean noise imposed +float FastMultEst::processNoiseImposed(const std::array ncl) +{ + // we assume that on the used layers the observed number of clusters is defined by the + // the noise ~ nu * Nchips and contribution from the signal tracks Ntr*conf.accCorr + // + // minimize the form sum_lr (noise_i - mu nchips_i)^2 / (mu nchips_i) + lambda_i * (noise_i + mult*acc_i - ncl_i) + // whith noise_i being estimate of the noise clusters in nchips_i of layer i, mu is the mean noise per chip, + // mult is the number of signal clusters on the ref. (1st) layer and the acc_i is the acceptance of layer i wrt 1st. + // The lambda_i is hust a Lagrange multiplier. + + const auto& conf = FastMultEstConfig::Instance(); + float w2sum = 0., wnsum = 0., wsum = 0.; + nLayersUsed = 0; + for (int il = conf.firstLayer; il <= conf.lastLayer; il++) { + if (ncl[il] > 0) { + float nchInv = 1. / o2::itsmft::ChipMappingITS::getNChipsPerLr(il); + w2sum += conf.accCorr[il] * conf.accCorr[il] * nchInv; + wnsum += ncl[il] * nchInv * conf.accCorr[il]; + wsum += conf.accCorr[il]; + nLayersUsed++; + } + } + mult = 0; + chi2 = -1; + noisePerChip = conf.imposeNoisePerChip; + if (nLayersUsed < 1) { + return -1; + } + auto w2sumI = 1. / w2sum; + mult = (wnsum - noisePerChip * wsum) * w2sumI; + cov[0] = wnsum * w2sumI; + cov[2] = 0.; + cov[1] = 0.; + + chi2 = 0.; + for (int il = conf.firstLayer; il <= conf.lastLayer; il++) { + if (ncl[il] > 0) { + float noise = ncl[il] - mult * conf.accCorr[il], estNoise = o2::itsmft::ChipMappingITS::getNChipsPerLr(il) * noisePerChip; + float diff = noise - estNoise; + chi2 += diff * diff / estNoise; + } + } + chi2 = nLayersUsed > 2 ? chi2 / (nLayersUsed - 2) : 0.; + return mult > 0 ? mult : 0; +} + +int FastMultEst::selectROFs(const gsl::span rofs, const gsl::span clus, + const gsl::span trig, std::vector& sel) +{ + int nrof = rofs.size(), nsel = 0; + const auto& multEstConf = FastMultEstConfig::Instance(); // parameters for mult estimation and cuts + sel.clear(); + sel.resize(nrof, true); // by default select all + lastRandomSeed = gRandom->GetSeed(); + if (multEstConf.isMultCutRequested()) { + for (uint32_t irof = 0; irof < nrof; irof++) { + nsel += sel[irof] = multEstConf.isPassingMultCut(process(rofs[irof].getROFData(clus))); + } + } else { + nsel = nrof; + } + using IdNT = std::pair; + if (multEstConf.cutRandomFraction > 0.) { + int ntrig = trig.size(), currTrig = 0; + if (multEstConf.preferTriggered) { + const auto& alpParams = o2::itsmft::DPLAlpideParam::Instance(); + std::vector nTrigROF; + nTrigROF.reserve(nrof); + for (uint32_t irof = 0; irof < nrof; irof++) { + if (sel[irof]) { + if (nsel && gRandom->Rndm() < multEstConf.cutRandomFraction) { + nsel--; + } + auto irROF = rofs[irof].getBCData(); + while (currTrig < ntrig && trig[currTrig].ir < irROF) { // triggers are sorted, jump to 1st one not less than current ROF + currTrig++; + } + auto& trof = nTrigROF.emplace_back(irof, 0); + irROF += alpParams.roFrameLengthInBC; + while (currTrig < ntrig && trig[currTrig].ir < irROF) { + trof.second++; + currTrig++; + } + } + } + if (nsel > 0) { + sort(nTrigROF.begin(), nTrigROF.end(), [](const IdNT& a, const IdNT& b) { return a.second > b.second; }); // order in number of triggers + auto last = nTrigROF.begin() + nsel; + sort(nTrigROF.begin(), last, [](const IdNT& a, const IdNT& b) { return a.first < b.first; }); // order in ROF ID first nsel ROFs + } + for (int i = nsel; i < int(nTrigROF.size()); i++) { // reject ROFs in the tail + sel[nTrigROF[i].first] = false; + } + } else { // dummy random rejection + for (int irof = 0; irof < nrof; irof++) { + if (sel[irof]) { + float sr = gRandom->Rndm(); + if (gRandom->Rndm() < multEstConf.cutRandomFraction) { + sel[irof] = false; + nsel--; + } + } + } + } + } + LOGP(debug, "NSel = {} of {} rofs Seeds: before {} after {}", nsel, nrof, lastRandomSeed, gRandom->GetSeed()); + + return nsel; +} diff --git a/Detectors/Upgrades/ITS3/reconstruction/src/FastMultEstConfig.cxx b/Detectors/Upgrades/ITS3/reconstruction/src/FastMultEstConfig.cxx new file mode 100644 index 0000000000000..4f3a8a44b0391 --- /dev/null +++ b/Detectors/Upgrades/ITS3/reconstruction/src/FastMultEstConfig.cxx @@ -0,0 +1,22 @@ +// 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. + +#include "ITS3Reconstruction/FastMultEstConfig.h" +#include "TRandom.h" + +O2ParamImpl(o2::its::FastMultEstConfig); + +using namespace o2::its; + +bool FastMultEstConfig::isPassingRandomRejection() const +{ + return (cutRandomFraction <= 0. || gRandom->Rndm() > cutRandomFraction) ? true : false; +} diff --git a/Detectors/Upgrades/IT3/reconstruction/src/ITS3ReconstructionLinkDef.h b/Detectors/Upgrades/ITS3/reconstruction/src/ITS3ReconstructionLinkDef.h similarity index 100% rename from Detectors/Upgrades/IT3/reconstruction/src/ITS3ReconstructionLinkDef.h rename to Detectors/Upgrades/ITS3/reconstruction/src/ITS3ReconstructionLinkDef.h diff --git a/Detectors/Upgrades/IT3/reconstruction/src/TopologyDictionary.cxx b/Detectors/Upgrades/ITS3/reconstruction/src/TopologyDictionary.cxx similarity index 62% rename from Detectors/Upgrades/IT3/reconstruction/src/TopologyDictionary.cxx rename to Detectors/Upgrades/ITS3/reconstruction/src/TopologyDictionary.cxx index ea632ea39401e..51ae9d2699782 100644 --- a/Detectors/Upgrades/IT3/reconstruction/src/TopologyDictionary.cxx +++ b/Detectors/Upgrades/ITS3/reconstruction/src/TopologyDictionary.cxx @@ -13,6 +13,7 @@ #include "ITS3Reconstruction/TopologyDictionary.h" #include "ITS3Base/SegmentationSuperAlpide.h" +#include "ITSMFTBase/SegmentationAlpide.h" namespace o2::its3 { @@ -20,10 +21,14 @@ namespace o2::its3 math_utils::Point3D TopologyDictionary::getClusterCoordinates(int detID, const its3::CompCluster& cl) const { - static SegmentationSuperAlpide segmentations[SegmentationSuperAlpide::NLayers]{SegmentationSuperAlpide(0), SegmentationSuperAlpide(1), SegmentationSuperAlpide(2), SegmentationSuperAlpide(3)}; + static SegmentationSuperAlpide segmentations[6]{SegmentationSuperAlpide(0), + SegmentationSuperAlpide(0), + SegmentationSuperAlpide(1), + SegmentationSuperAlpide(1), + SegmentationSuperAlpide(2), + SegmentationSuperAlpide(2)}; // TODO: fix NLayers math_utils::Point3D locCl; - if (detID >= SegmentationSuperAlpide::NLayers) { - + if (detID >= 6) { // TODO: fix NLayers o2::itsmft::SegmentationAlpide::detectorToLocalUnchecked(cl.getRow(), cl.getCol(), locCl); } else { segmentations[detID].detectorToLocalUnchecked(cl.getRow(), cl.getCol(), locCl); @@ -35,7 +40,12 @@ math_utils::Point3D TopologyDictionary::getClusterCoordinates(int detID, math_utils::Point3D TopologyDictionary::getClusterCoordinates(int detID, const its3::CompCluster& cl, const itsmft::ClusterPattern& patt, bool isGroup) { - static SegmentationSuperAlpide segmentations[SegmentationSuperAlpide::NLayers]{SegmentationSuperAlpide(0), SegmentationSuperAlpide(1), SegmentationSuperAlpide(2), SegmentationSuperAlpide(3)}; + static SegmentationSuperAlpide segmentations[6]{SegmentationSuperAlpide(0), + SegmentationSuperAlpide(0), + SegmentationSuperAlpide(1), + SegmentationSuperAlpide(1), + SegmentationSuperAlpide(2), + SegmentationSuperAlpide(2)}; // TODO: fix NLayers auto refRow = cl.getRow(); auto refCol = cl.getCol(); @@ -46,7 +56,7 @@ math_utils::Point3D TopologyDictionary::getClusterCoordinates(int detID, refCol -= round(zCOG); } math_utils::Point3D locCl; - if (detID >= SegmentationSuperAlpide::NLayers) { + if (detID >= 6) { // TODO: fix NLayers o2::itsmft::SegmentationAlpide::detectorToLocalUnchecked(refRow + xCOG, refCol + zCOG, locCl); } else { segmentations[detID].detectorToLocalUnchecked(refRow + xCOG, refCol + zCOG, locCl); diff --git a/Detectors/Upgrades/IT3/simulation/CMakeLists.txt b/Detectors/Upgrades/ITS3/simulation/CMakeLists.txt similarity index 71% rename from Detectors/Upgrades/IT3/simulation/CMakeLists.txt rename to Detectors/Upgrades/ITS3/simulation/CMakeLists.txt index 0c8a2f70b8e93..7d652db91da59 100644 --- a/Detectors/Upgrades/IT3/simulation/CMakeLists.txt +++ b/Detectors/Upgrades/ITS3/simulation/CMakeLists.txt @@ -10,15 +10,19 @@ # or submit itself to any jurisdiction. o2_add_library(ITS3Simulation - SOURCES src/ITS3Layer.cxx src/DescriptorInnerBarrelITS3.cxx src/DescriptorInnerBarrelITS3Param.cxx - PUBLIC_LINK_LIBRARIES O2::SimulationDataFormat O2::ITSBase O2::ITSMFTSimulation + SOURCES src/ITS3Layer.cxx + src/DescriptorInnerBarrelITS3.cxx + src/DescriptorInnerBarrelITS3Param.cxx + src/Digitizer.cxx + PUBLIC_LINK_LIBRARIES O2::SimulationDataFormat + O2::ITSBase O2::ITSMFTSimulation ROOT::Physics) o2_target_root_dictionary(ITS3Simulation HEADERS include/ITS3Simulation/ITS3Layer.h include/ITS3Simulation/DescriptorInnerBarrelITS3.h include/ITS3Simulation/DescriptorInnerBarrelITS3Param.h - # include/ITS3Simulation/Digitizer.h + include/ITS3Simulation/Digitizer.h ) o2_data_file(COPY data DESTINATION Detectors/ITS3/simulation) \ No newline at end of file diff --git a/Detectors/Upgrades/IT3/simulation/data/simcuts.dat b/Detectors/Upgrades/ITS3/simulation/data/simcuts.dat similarity index 100% rename from Detectors/Upgrades/IT3/simulation/data/simcuts.dat rename to Detectors/Upgrades/ITS3/simulation/data/simcuts.dat diff --git a/Detectors/Upgrades/IT3/simulation/include/ITS3Simulation/DescriptorInnerBarrelITS3.h b/Detectors/Upgrades/ITS3/simulation/include/ITS3Simulation/DescriptorInnerBarrelITS3.h similarity index 100% rename from Detectors/Upgrades/IT3/simulation/include/ITS3Simulation/DescriptorInnerBarrelITS3.h rename to Detectors/Upgrades/ITS3/simulation/include/ITS3Simulation/DescriptorInnerBarrelITS3.h diff --git a/Detectors/Upgrades/IT3/simulation/include/ITS3Simulation/DescriptorInnerBarrelITS3Param.h b/Detectors/Upgrades/ITS3/simulation/include/ITS3Simulation/DescriptorInnerBarrelITS3Param.h similarity index 100% rename from Detectors/Upgrades/IT3/simulation/include/ITS3Simulation/DescriptorInnerBarrelITS3Param.h rename to Detectors/Upgrades/ITS3/simulation/include/ITS3Simulation/DescriptorInnerBarrelITS3Param.h diff --git a/Detectors/Upgrades/IT3/simulation/include/ITS3Simulation/Digitizer.h b/Detectors/Upgrades/ITS3/simulation/include/ITS3Simulation/Digitizer.h similarity index 95% rename from Detectors/Upgrades/IT3/simulation/include/ITS3Simulation/Digitizer.h rename to Detectors/Upgrades/ITS3/simulation/include/ITS3Simulation/Digitizer.h index 6b94a043dc116..b963fca38bad3 100644 --- a/Detectors/Upgrades/IT3/simulation/include/ITS3Simulation/Digitizer.h +++ b/Detectors/Upgrades/ITS3/simulation/include/ITS3Simulation/Digitizer.h @@ -25,7 +25,7 @@ #include "ITSMFTSimulation/AlpideSimResponse.h" #include "ITSMFTSimulation/DigiParams.h" #include "ITSMFTSimulation/Hit.h" -#include "ITS3Base/GeometryTGeo.h" +#include "ITSBase/GeometryTGeo.h" #include "ITS3Base/SegmentationSuperAlpide.h" #include "DataFormatsITSMFT/Digit.h" #include "DataFormatsITSMFT/ROFRecord.h" @@ -79,7 +79,7 @@ class Digitizer : public TObject const o2::itsmft::DigiParams& getDigitParams() const { return mParams; } // provide the common itsmft::GeometryTGeo to access matrices and segmentation - void setGeometry(const o2::its3::GeometryTGeo* gm) { mGeometry = gm; } + void setGeometry(const o2::its::GeometryTGeo* gm) { mGeometry = gm; } uint32_t getEventROFrameMin() const { return mEventROFrameMin; } uint32_t getEventROFrameMax() const { return mEventROFrameMax; } @@ -107,6 +107,7 @@ class Digitizer : public TObject } std::vector mSuperSegmentations; + std::vector mLayerID; static constexpr float sec2ns = 1e9; o2::itsmft::DigiParams mParams; ///< digitization parameters @@ -122,7 +123,7 @@ class Digitizer : public TObject o2::itsmft::AlpideSimResponse* mAlpSimResp = nullptr; // simulated response - const o2::its3::GeometryTGeo* mGeometry = nullptr; ///< ITS OR MFT upgrade geometry + const o2::its::GeometryTGeo* mGeometry = nullptr; ///< ITS OR MFT upgrade geometry std::vector mChips; ///< Array of chips digits containers std::deque> mExtraBuff; ///< burrer (per roFrame) for extra digits diff --git a/Detectors/Upgrades/IT3/simulation/include/ITS3Simulation/ITS3Layer.h b/Detectors/Upgrades/ITS3/simulation/include/ITS3Simulation/ITS3Layer.h similarity index 100% rename from Detectors/Upgrades/IT3/simulation/include/ITS3Simulation/ITS3Layer.h rename to Detectors/Upgrades/ITS3/simulation/include/ITS3Simulation/ITS3Layer.h diff --git a/Detectors/Upgrades/IT3/simulation/src/DescriptorInnerBarrelITS3.cxx b/Detectors/Upgrades/ITS3/simulation/src/DescriptorInnerBarrelITS3.cxx similarity index 98% rename from Detectors/Upgrades/IT3/simulation/src/DescriptorInnerBarrelITS3.cxx rename to Detectors/Upgrades/ITS3/simulation/src/DescriptorInnerBarrelITS3.cxx index e82b28f50574a..6e2a0eb26e409 100644 --- a/Detectors/Upgrades/IT3/simulation/src/DescriptorInnerBarrelITS3.cxx +++ b/Detectors/Upgrades/ITS3/simulation/src/DescriptorInnerBarrelITS3.cxx @@ -60,7 +60,7 @@ void DescriptorInnerBarrelITS3::configure() mNumLayers = 5; } - mSensorLayerThickness = SegmentationSuperAlpide::SensorLayerThickness; + mSensorLayerThickness = SegmentationSuperAlpide::mSensorLayerThickness; // build ITS3 upgrade detector mLayer.resize(mNumLayers); diff --git a/Detectors/Upgrades/IT3/simulation/src/DescriptorInnerBarrelITS3Param.cxx b/Detectors/Upgrades/ITS3/simulation/src/DescriptorInnerBarrelITS3Param.cxx similarity index 100% rename from Detectors/Upgrades/IT3/simulation/src/DescriptorInnerBarrelITS3Param.cxx rename to Detectors/Upgrades/ITS3/simulation/src/DescriptorInnerBarrelITS3Param.cxx diff --git a/Detectors/Upgrades/IT3/simulation/src/Digitizer.cxx b/Detectors/Upgrades/ITS3/simulation/src/Digitizer.cxx similarity index 88% rename from Detectors/Upgrades/IT3/simulation/src/Digitizer.cxx rename to Detectors/Upgrades/ITS3/simulation/src/Digitizer.cxx index 184d8e45fecae..5740ac8d66795 100644 --- a/Detectors/Upgrades/IT3/simulation/src/Digitizer.cxx +++ b/Detectors/Upgrades/ITS3/simulation/src/Digitizer.cxx @@ -18,6 +18,7 @@ #include "MathUtils/Cartesian.h" #include "SimulationDataFormat/MCTruthContainer.h" #include "DetectorsRaw/HBFUtils.h" +#include "CommonConstants/MathConstants.h" #include #include @@ -37,17 +38,23 @@ using namespace o2::its3; void Digitizer::init() { - for (int iLayer{0}; iLayer < SegmentationSuperAlpide::NLayers; ++iLayer) { - mSuperSegmentations.push_back(SegmentationSuperAlpide(iLayer)); + mLayerID.clear(); + mSuperSegmentations.clear(); + for (int iLayer{0}; iLayer < mGeometry->getNumberOfLayers() - 4; ++iLayer) { + for (int iChip{0}; iChip < mGeometry->getNumberOfChipsPerLayer(iLayer); ++iChip) { + LOGP(info, "layer: {} chip: {}", iLayer, iChip); + mSuperSegmentations.push_back(SegmentationSuperAlpide(iLayer)); + mLayerID.push_back(iLayer); + } } - const Int_t numOfChips = mGeometry->getNumberOfChips() + SegmentationSuperAlpide::NLayers; + const int numOfChips = mGeometry->getNumberOfChips(); mChips.resize(numOfChips); for (int i = numOfChips; i--;) { mChips[i].setChipIndex(i); } if (!mParams.getAlpSimResponse()) { - std::string responseFile = "$(O2_ROOT)/share/Detectors/ITSMFT/data/alpideResponseData/AlpideResponseData.root"; + std::string responseFile = "$(O2_ROOT)/share/Detectors/ITSMFT/data/AlpideResponseData/AlpideResponseData.root"; auto file = TFile::Open(responseFile.data()); mAlpSimResp = (o2::itsmft::AlpideSimResponse*)file->Get("response1"); mParams.setAlpSimResponse(mAlpSimResp); @@ -143,8 +150,8 @@ void Digitizer::fillOutputContainer(uint32_t frameLast) auto& extra = *(mExtraBuff.front().get()); for (int iChip{0}; iChip < mChips.size(); ++iChip) { auto& chip = mChips[iChip]; - if (iChip < SegmentationSuperAlpide::NLayers) { - chip.addNoise(mROFrameMin, mROFrameMin, &mParams, mSuperSegmentations[iChip].NRows, mSuperSegmentations[iChip].NCols); + if (iChip < mSuperSegmentations.size()) { // Check if this is a chip of ITS3 + chip.addNoise(mROFrameMin, mROFrameMin, &mParams, mSuperSegmentations[iChip].mNRows, mSuperSegmentations[iChip].mNCols); } else { chip.addNoise(mROFrameMin, mROFrameMin, &mParams); } @@ -225,34 +232,35 @@ void Digitizer::processHit(const o2::itsmft::Hit& hit, uint32_t& maxFr, int evID float nStepsInv = mParams.getNSimStepsInv(); int nSteps = mParams.getNSimSteps(); short detID{hit.GetDetectorID()}; - const auto& matrix = mGeometry->getMatrixL2G(detID); - bool innerBarrel{detID < SegmentationSuperAlpide::NLayers}; - auto startPos = hit.GetPosStart(); - float startPhi{std::atan2(-startPos.Y(), -startPos.X())}; - auto endPos = hit.GetPos(); - float endPhi{std::atan2(-endPos.Y(), -endPos.X())}; + const auto& matrix = mGeometry->getMatrixL2G(detID); // <<<< ????? + bool innerBarrel{detID < mSuperSegmentations.size()}; math_utils::Vector3D xyzLocS, xyzLocE; + xyzLocS = matrix ^ (hit.GetPosStart()); + xyzLocE = matrix ^ (hit.GetPos()); if (innerBarrel) { - xyzLocS = {SegmentationSuperAlpide::Radii[detID] * startPhi, 0.f, startPos.Z()}; - xyzLocE = {SegmentationSuperAlpide::Radii[detID] * endPhi, 0.f, endPos.Z()}; - } else { - xyzLocS = matrix ^ (hit.GetPosStart()); - xyzLocE = matrix ^ (hit.GetPos()); + // transform the point on the curved surface to a flat one + float xFlatE{0.f}, yFlatE{0.f}, xFlatS{0.f}, yFlatS{0.f}; + mSuperSegmentations[detID].curvedToFlat(xyzLocS.X(), xyzLocS.Y(), xFlatS, yFlatS); + mSuperSegmentations[detID].curvedToFlat(xyzLocE.X(), xyzLocE.Y(), xFlatE, yFlatE); + // update the local coordinates with the flattened ones + xyzLocS.SetXYZ(xFlatS, yFlatS, xyzLocS.Z()); + xyzLocE.SetXYZ(xFlatE, yFlatE, xyzLocE.Z()); } math_utils::Vector3D step(xyzLocE); step -= xyzLocS; step *= nStepsInv; // position increment at each step - // the electrons will injected in the middle of each step + // the electrons will be injected in the middle of each step math_utils::Vector3D stepH(step * 0.5); - xyzLocS += stepH; - xyzLocE -= stepH; + xyzLocS += stepH; // Adjust start position to the middle of the first step + xyzLocE -= stepH; // Adjust end position to the middle of the last step int rowS = -1, colS = -1, rowE = -1, colE = -1, nSkip = 0; if (innerBarrel) { // get entrance pixel row and col while (!mSuperSegmentations[detID].localToDetector(xyzLocS.X(), xyzLocS.Z(), rowS, colS)) { // guard-ring ? if (++nSkip >= nSteps) { + LOGP(info, "Start: detId {}", detID); return; // did not enter to sensitive matrix } xyzLocS += step; @@ -260,10 +268,12 @@ void Digitizer::processHit(const o2::itsmft::Hit& hit, uint32_t& maxFr, int evID // get exit pixel row and col while (!mSuperSegmentations[detID].localToDetector(xyzLocE.X(), xyzLocE.Z(), rowE, colE)) { // guard-ring ? if (++nSkip >= nSteps) { + LOGP(info, "End: detId {}", detID); return; // did not enter to sensitive matrix } xyzLocE -= step; } + // LOGP(info, "x: {}, z: {}, col: {}, row: {}, detID: {}", xyzLocS.X(), xyzLocS.Z(), colS, rowS, detID); } else { // get entrance pixel row and col while (!Segmentation::localToDetector(xyzLocS.X(), xyzLocS.Z(), rowS, colS)) { // guard-ring ? @@ -293,8 +303,8 @@ void Digitizer::processHit(const o2::itsmft::Hit& hit, uint32_t& maxFr, int evID rowS = 0; } - int maxNrows{innerBarrel ? mSuperSegmentations[detID].NRows : Segmentation::NRows}; - int maxNcols{innerBarrel ? mSuperSegmentations[detID].NCols : Segmentation::NCols}; + int maxNrows{innerBarrel ? mSuperSegmentations[detID].mNRows : Segmentation::NRows}; + int maxNcols{innerBarrel ? mSuperSegmentations[detID].mNCols : Segmentation::NCols}; if (rowE >= maxNrows) { rowE = maxNrows - 1; } @@ -348,8 +358,8 @@ void Digitizer::processHit(const o2::itsmft::Hit& hit, uint32_t& maxFr, int evID } bool flipCol, flipRow; // note that response needs coordinates along column row (locX) (locZ) then depth (locY) - float rowMax{0.5f * (innerBarrel ? mSuperSegmentations[detID].PitchRow : Segmentation::PitchRow)}; - float colMax{0.5f * (innerBarrel ? mSuperSegmentations[detID].PitchCol : Segmentation::PitchCol)}; + float rowMax{0.5f * (innerBarrel ? mSuperSegmentations[detID].mPitchRow : Segmentation::PitchRow)}; + float colMax{0.5f * (innerBarrel ? mSuperSegmentations[detID].mPitchCol : Segmentation::PitchCol)}; auto rspmat = resp->getResponse(xyzLocS.X() - cRowPix, xyzLocS.Z() - cColPix, xyzLocS.Y(), flipRow, flipCol, rowMax, colMax); xyzLocS += step; diff --git a/Detectors/Upgrades/IT3/simulation/src/ITS3Layer.cxx b/Detectors/Upgrades/ITS3/simulation/src/ITS3Layer.cxx similarity index 100% rename from Detectors/Upgrades/IT3/simulation/src/ITS3Layer.cxx rename to Detectors/Upgrades/ITS3/simulation/src/ITS3Layer.cxx diff --git a/Detectors/Upgrades/IT3/simulation/src/ITS3SimulationLinkDef.h b/Detectors/Upgrades/ITS3/simulation/src/ITS3SimulationLinkDef.h similarity index 94% rename from Detectors/Upgrades/IT3/simulation/src/ITS3SimulationLinkDef.h rename to Detectors/Upgrades/ITS3/simulation/src/ITS3SimulationLinkDef.h index 73dd65761fa1a..4febbb8f853e3 100644 --- a/Detectors/Upgrades/IT3/simulation/src/ITS3SimulationLinkDef.h +++ b/Detectors/Upgrades/ITS3/simulation/src/ITS3SimulationLinkDef.h @@ -18,6 +18,6 @@ #pragma link C++ class o2::its3::ITS3Layer + ; #pragma link C++ class o2::its3::DescriptorInnerBarrelITS3 + ; #pragma link C++ class o2::its3::DescriptorInnerBarrelITS3Param + ; -//#pragma link C++ class o2::its3::Digitizer + ; +#pragma link C++ class o2::its3::Digitizer + ; #endif diff --git a/Detectors/Upgrades/IT3/workflow/CMakeLists.txt b/Detectors/Upgrades/ITS3/workflow/CMakeLists.txt similarity index 81% rename from Detectors/Upgrades/IT3/workflow/CMakeLists.txt rename to Detectors/Upgrades/ITS3/workflow/CMakeLists.txt index 7aca1a47c7743..13566c6a33f59 100644 --- a/Detectors/Upgrades/IT3/workflow/CMakeLists.txt +++ b/Detectors/Upgrades/ITS3/workflow/CMakeLists.txt @@ -23,15 +23,16 @@ o2_add_library(ITS3Workflow # src/VertexReaderSpec.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2::SimConfig - O2::DataFormatsITS + O2::DataFormatsITSMFT O2::DataFormatsITS3 O2::SimulationDataFormat O2::ITS3Simulation - # O2::ITStracking - # O2::ITSReconstruction + O2::ITStracking + O2::ITSMFTReconstruction O2::ITS3Reconstruction - # O2::ITSMFTWorkflow + O2::ITSWorkflow O2::GPUTracking + O2::ITSBase ) # o2_add_executable(digit-writer-workflow @@ -44,10 +45,10 @@ o2_add_library(ITS3Workflow # COMPONENT_NAME its3 # PUBLIC_LINK_LIBRARIES O2::ITS3Workflow) -# o2_add_executable(reco-workflow -# SOURCES src/its3-reco-workflow.cxx -# COMPONENT_NAME its3 -# PUBLIC_LINK_LIBRARIES O2::ITS3Workflow) +o2_add_executable(reco-workflow + SOURCES src/its3-reco-workflow.cxx + COMPONENT_NAME its3 + PUBLIC_LINK_LIBRARIES O2::ITS3Workflow) # o2_add_executable(cluster-writer-workflow # SOURCES src/its-cluster-writer-workflow.cxx @@ -57,5 +58,4 @@ o2_add_library(ITS3Workflow # o2_add_executable(cluster-reader-workflow # SOURCES src/its-cluster-reader-workflow.cxx # COMPONENT_NAME its -# PUBLIC_LINK_LIBRARIES O2::ITSWorkflow) - +# PUBLIC_LINK_LIBRARIES O2::ITSWorkflow) \ No newline at end of file diff --git a/Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/ClusterReaderSpec.h b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/ClusterReaderSpec.h similarity index 100% rename from Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/ClusterReaderSpec.h rename to Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/ClusterReaderSpec.h diff --git a/Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/ClusterWriterSpec.h b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/ClusterWriterSpec.h similarity index 100% rename from Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/ClusterWriterSpec.h rename to Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/ClusterWriterSpec.h diff --git a/Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/ClusterWriterWorkflow.h b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/ClusterWriterWorkflow.h similarity index 100% rename from Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/ClusterWriterWorkflow.h rename to Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/ClusterWriterWorkflow.h diff --git a/Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/ClustererSpec.h b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/ClustererSpec.h similarity index 95% rename from Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/ClustererSpec.h rename to Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/ClustererSpec.h index c1dd586a86766..3a96921812c67 100644 --- a/Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/ClustererSpec.h +++ b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/ClustererSpec.h @@ -11,8 +11,8 @@ /// @file ClustererSpec.h -#ifndef O2_ITS_CLUSTERERDPL -#define O2_ITS_CLUSTERERDPL +#ifndef O2_ITS3_CLUSTERERDPL +#define O2_ITS3_CLUSTERERDPL #include diff --git a/Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/DigitReaderSpec.h b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/DigitReaderSpec.h similarity index 100% rename from Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/DigitReaderSpec.h rename to Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/DigitReaderSpec.h diff --git a/Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/DigitWriterSpec.h b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/DigitWriterSpec.h similarity index 100% rename from Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/DigitWriterSpec.h rename to Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/DigitWriterSpec.h diff --git a/Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/RecoWorkflow.h b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/RecoWorkflow.h similarity index 62% rename from Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/RecoWorkflow.h rename to Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/RecoWorkflow.h index 0d4dbe0564e41..eeeb2d3742df2 100644 --- a/Detectors/Upgrades/IT3/workflow/include/ITS3Workflow/RecoWorkflow.h +++ b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/RecoWorkflow.h @@ -9,8 +9,8 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. -#ifndef O2_ITS_RECOWORKFLOW_H -#define O2_ITS_RECOWORKFLOW_H +#ifndef O2_ITS3_RECOWORKFLOW_H +#define O2_ITS3_RECOWORKFLOW_H /// @file RecoWorkflow.h @@ -28,8 +28,13 @@ namespace its3 namespace reco_workflow { -framework::WorkflowSpec getWorkflow(bool useMC, const std::string& trmode, o2::gpu::GPUDataTypes::DeviceType dType = o2::gpu::GPUDataTypes::DeviceType::CPU, - bool upstreamDigits = false, bool upstreamClusters = false, bool disableRootOutput = false, bool eencode = false); +framework::WorkflowSpec getWorkflow(bool useMC, + const std::string& trmode, + o2::gpu::GPUDataTypes::DeviceType dtype, + bool upstreamDigits, + bool upstreamClusters, + bool disableRootOutput, + int useTrig); } } // namespace its3 diff --git a/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/TrackerSpec.h b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/TrackerSpec.h new file mode 100644 index 0000000000000..ffeb0b535bfed --- /dev/null +++ b/Detectors/Upgrades/ITS3/workflow/include/ITS3Workflow/TrackerSpec.h @@ -0,0 +1,76 @@ +// 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 TrackerSpec.h + +#ifndef O2_ITS3_TRACKERDPL +#define O2_ITS3_TRACKERDPL + +#include "DataFormatsParameters/GRPObject.h" +#include "DataFormatsITSMFT/TopologyDictionary.h" + +#include "Framework/DataProcessorSpec.h" +#include "Framework/Task.h" + +#include "ITStracking/TimeFrame.h" +#include "ITStracking/Tracker.h" +#include "ITStracking/TrackerTraits.h" +#include "ITStracking/Vertexer.h" +#include "ITStracking/VertexerTraits.h" + +#include "GPUO2Interface.h" +#include "GPUReconstruction.h" +#include "GPUChainITS.h" +#include "CommonUtils/StringUtils.h" +#include "TStopwatch.h" +#include "DetectorsBase/GRPGeomHelper.h" + +namespace o2 +{ +namespace its3 +{ + +class TrackerDPL : public framework::Task +{ + public: + TrackerDPL(std::shared_ptr gr, bool isMC, int trgType, const std::string& trModeS, o2::gpu::GPUDataTypes::DeviceType dType = o2::gpu::GPUDataTypes::DeviceType::CPU); + ~TrackerDPL() override = default; + void init(framework::InitContext& ic) final; + void run(framework::ProcessingContext& pc) final; + void endOfStream(framework::EndOfStreamContext& ec) final; + void finaliseCCDB(framework::ConcreteDataMatcher& matcher, void* obj) final; + void setClusterDictionary(const o2::itsmft::TopologyDictionary* d) { mDict = d; } + + private: + void updateTimeDependentParams(framework::ProcessingContext& pc); + + bool mIsMC = false; + bool mRunVertexer = true; + bool mCosmicsProcessing = false; + int mUseTriggers = 0; + std::string mMode = "sync"; + std::shared_ptr mGGCCDBRequest; + const o2::itsmft::TopologyDictionary* mDict = nullptr; + std::unique_ptr mRecChain = nullptr; + std::unique_ptr mChainITS = nullptr; + std::unique_ptr mTracker = nullptr; + std::unique_ptr mVertexer = nullptr; + TStopwatch mTimer; +}; + +/// create a processor spec +/// run ITS CA tracker +framework::DataProcessorSpec getTrackerSpec(bool useMC, int useTrig, const std::string& trModeS, o2::gpu::GPUDataTypes::DeviceType dType); + +} // namespace its3 +} // namespace o2 + +#endif /* O2_ITS_TRACKERDPL */ diff --git a/Detectors/Upgrades/IT3/workflow/src/ClusterReaderSpec.cxx b/Detectors/Upgrades/ITS3/workflow/src/ClusterReaderSpec.cxx similarity index 100% rename from Detectors/Upgrades/IT3/workflow/src/ClusterReaderSpec.cxx rename to Detectors/Upgrades/ITS3/workflow/src/ClusterReaderSpec.cxx diff --git a/Detectors/Upgrades/IT3/workflow/src/ClusterWriterSpec.cxx b/Detectors/Upgrades/ITS3/workflow/src/ClusterWriterSpec.cxx similarity index 100% rename from Detectors/Upgrades/IT3/workflow/src/ClusterWriterSpec.cxx rename to Detectors/Upgrades/ITS3/workflow/src/ClusterWriterSpec.cxx diff --git a/Detectors/Upgrades/IT3/workflow/src/ClusterWriterWorkflow.cxx b/Detectors/Upgrades/ITS3/workflow/src/ClusterWriterWorkflow.cxx similarity index 100% rename from Detectors/Upgrades/IT3/workflow/src/ClusterWriterWorkflow.cxx rename to Detectors/Upgrades/ITS3/workflow/src/ClusterWriterWorkflow.cxx diff --git a/Detectors/Upgrades/IT3/workflow/src/ClustererSpec.cxx b/Detectors/Upgrades/ITS3/workflow/src/ClustererSpec.cxx similarity index 98% rename from Detectors/Upgrades/IT3/workflow/src/ClustererSpec.cxx rename to Detectors/Upgrades/ITS3/workflow/src/ClustererSpec.cxx index 93890fa01fffe..cf0c7c38b9aff 100644 --- a/Detectors/Upgrades/IT3/workflow/src/ClustererSpec.cxx +++ b/Detectors/Upgrades/ITS3/workflow/src/ClustererSpec.cxx @@ -18,7 +18,6 @@ #include "ITS3Workflow/ClustererSpec.h" #include "DataFormatsITSMFT/Digit.h" #include "ITSMFTReconstruction/ChipMappingITS.h" -#include "ITS3Base/SegmentationSuperAlpide.h" #include "ITSMFTReconstruction/ClustererParam.h" #include "DataFormatsITS3/CompCluster.h" #include "SimulationDataFormat/MCCompLabel.h" @@ -40,7 +39,7 @@ namespace its3 void ClustererDPL::init(InitContext& ic) { mClusterer = std::make_unique(); - mClusterer->setNChips(o2::itsmft::ChipMappingITS::getNChips(o2::itsmft::ChipMappingITS::MB) + o2::itsmft::ChipMappingITS::getNChips(o2::itsmft::ChipMappingITS::OB) + SegmentationSuperAlpide::NLayers); + mClusterer->setNChips(o2::itsmft::ChipMappingITS::getNChips(o2::itsmft::ChipMappingITS::MB) + o2::itsmft::ChipMappingITS::getNChips(o2::itsmft::ChipMappingITS::OB) + 6); /// Fix MEEEE auto filenameGRP = ic.options().get("grp-file"); const auto grp = o2::parameters::GRPObject::loadFrom(filenameGRP.c_str()); diff --git a/Detectors/Upgrades/IT3/workflow/src/DigitReaderSpec.cxx b/Detectors/Upgrades/ITS3/workflow/src/DigitReaderSpec.cxx similarity index 100% rename from Detectors/Upgrades/IT3/workflow/src/DigitReaderSpec.cxx rename to Detectors/Upgrades/ITS3/workflow/src/DigitReaderSpec.cxx diff --git a/Detectors/Upgrades/IT3/workflow/src/DigitWriterSpec.cxx b/Detectors/Upgrades/ITS3/workflow/src/DigitWriterSpec.cxx similarity index 100% rename from Detectors/Upgrades/IT3/workflow/src/DigitWriterSpec.cxx rename to Detectors/Upgrades/ITS3/workflow/src/DigitWriterSpec.cxx diff --git a/Detectors/Upgrades/IT3/workflow/src/RecoWorkflow.cxx b/Detectors/Upgrades/ITS3/workflow/src/RecoWorkflow.cxx similarity index 78% rename from Detectors/Upgrades/IT3/workflow/src/RecoWorkflow.cxx rename to Detectors/Upgrades/ITS3/workflow/src/RecoWorkflow.cxx index 325687aafac48..87a067739c19f 100644 --- a/Detectors/Upgrades/IT3/workflow/src/RecoWorkflow.cxx +++ b/Detectors/Upgrades/ITS3/workflow/src/RecoWorkflow.cxx @@ -9,14 +9,11 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. -/// @file RecoWorkflow.cxx - #include "ITS3Workflow/RecoWorkflow.h" #include "ITS3Workflow/ClustererSpec.h" #include "ITS3Workflow/ClusterWriterSpec.h" -// #include "ITSWorkflow/TrackerSpec.h" -// #include "ITSWorkflow/TrackWriterSpec.h" -// #include "ITSMFTWorkflow/EntropyEncoderSpec.h" +#include "ITSWorkflow/TrackerSpec.h" +#include "ITSWorkflow/TrackWriterSpec.h" #include "ITS3Workflow/DigitReaderSpec.h" namespace o2 @@ -28,8 +25,7 @@ namespace reco_workflow { framework::WorkflowSpec getWorkflow(bool useMC, const std::string& trmode, o2::gpu::GPUDataTypes::DeviceType dtype, - bool upstreamDigits, bool upstreamClusters, bool disableRootOutput, - bool) + bool upstreamDigits, bool upstreamClusters, bool disableRootOutput, int useTrig) { framework::WorkflowSpec specs; @@ -40,17 +36,16 @@ framework::WorkflowSpec getWorkflow(bool useMC, const std::string& trmode, o2::g if (!upstreamClusters) { specs.emplace_back(o2::its3::getClustererSpec(useMC)); } + if (!disableRootOutput) { specs.emplace_back(o2::its3::getClusterWriterSpec(useMC)); } - // specs.emplace_back(o2::its::getTrackerSpec(useMC, trmode, dtype)); + + // specs.emplace_back(o2::its::getTrackerSpec(useMC, useTrig, trmode, dtype)); // if (!disableRootOutput) { // specs.emplace_back(o2::its::getTrackWriterSpec(useMC)); // } - // if (eencode) { - // specs.emplace_back(o2::itsmft::getEntropyEncoderSpec("ITS")); - // } return specs; } diff --git a/Detectors/Upgrades/ITS3/workflow/src/TrackerSpec.cxx b/Detectors/Upgrades/ITS3/workflow/src/TrackerSpec.cxx new file mode 100644 index 0000000000000..d854f136a75cc --- /dev/null +++ b/Detectors/Upgrades/ITS3/workflow/src/TrackerSpec.cxx @@ -0,0 +1,414 @@ +// 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 TrackerSpec.cxx + +#include + +#include "TGeoGlobalMagField.h" + +#include "Framework/ControlService.h" +#include "Framework/ConfigParamRegistry.h" +#include "Framework/CCDBParamSpec.h" +#include "ITS3Workflow/TrackerSpec.h" +#include "DataFormatsITS3/CompCluster.h" +#include "DataFormatsITS/TrackITS.h" +#include "SimulationDataFormat/MCCompLabel.h" +#include "SimulationDataFormat/MCTruthContainer.h" +#include "DataFormatsITSMFT/ROFRecord.h" +#include "DataFormatsITSMFT/PhysTrigger.h" + +#include "ITStracking/ROframe.h" +#include "ITStracking/IOUtils.h" +#include "ITStracking/TrackingConfigParam.h" +#include "ITSMFTBase/DPLAlpideParam.h" +#include "ITSMFTReconstruction/ClustererParam.h" + +#include "Field/MagneticField.h" +#include "DetectorsBase/GeometryManager.h" +#include "DetectorsBase/Propagator.h" +#include "ITSBase/GeometryTGeo.h" +#include "CommonDataFormat/IRFrame.h" +#include "DetectorsCommonDataFormats/DetectorNameConf.h" +#include "DataFormatsTRD/TriggerRecord.h" +#include "ITSReconstruction/FastMultEstConfig.h" +#include "ITSReconstruction/FastMultEst.h" +#include + +namespace o2 +{ +using namespace framework; +using its::FastMultEst; +using its::FastMultEstConfig; +using its::TimeFrame; +using its::Tracker; +using its::TrackingParameters; +using its::TrackITSExt; +using its::Vertexer; + +namespace its3 +{ +using Vertex = o2::dataformats::Vertex>; + +TrackerDPL::TrackerDPL(std::shared_ptr gr, bool isMC, int trgType, const std::string& trModeS, o2::gpu::GPUDataTypes::DeviceType dType) : mGGCCDBRequest(gr), mIsMC{isMC}, mUseTriggers{trgType}, mMode{trModeS}, mRecChain{o2::gpu::GPUReconstruction::CreateInstance(dType, true)} +{ + std::transform(mMode.begin(), mMode.end(), mMode.begin(), [](unsigned char c) { return std::tolower(c); }); +} + +void TrackerDPL::init(InitContext& ic) +{ + mTimer.Stop(); + mTimer.Reset(); + o2::base::GRPGeomHelper::instance().setRequest(mGGCCDBRequest); + mChainITS.reset(mRecChain->AddChain()); + mVertexer = std::make_unique(mChainITS->GetITSVertexerTraits()); + mTracker = std::make_unique(mChainITS->GetITSTrackerTraits()); + mRunVertexer = true; + mCosmicsProcessing = false; + std::vector trackParams; + + if (mMode == "async") { + + trackParams.resize(3); + trackParams[1].TrackletMinPt = 0.2f; + trackParams[1].CellDeltaTanLambdaSigma *= 2.; + trackParams[2].TrackletMinPt = 0.1f; + trackParams[2].CellDeltaTanLambdaSigma *= 4.; + trackParams[2].MinTrackLength = 4; + LOG(info) << "Initializing tracker in async. phase reconstruction with " << trackParams.size() << " passes"; + + } else if (mMode == "sync_misaligned") { + + trackParams.resize(3); + trackParams[0].PhiBins = 32; + trackParams[0].ZBins = 64; + trackParams[0].CellDeltaTanLambdaSigma *= 3.; + trackParams[0].LayerMisalignment[0] = 1.e-2; + trackParams[0].LayerMisalignment[1] = 1.e-2; + trackParams[0].LayerMisalignment[2] = 1.e-2; + trackParams[0].LayerMisalignment[3] = 3.e-2; + trackParams[0].LayerMisalignment[4] = 3.e-2; + trackParams[0].LayerMisalignment[5] = 3.e-2; + trackParams[0].LayerMisalignment[6] = 3.e-2; + trackParams[0].FitIterationMaxChi2[0] = 50.; + trackParams[0].FitIterationMaxChi2[1] = 25.; + trackParams[1] = trackParams[0]; + trackParams[2] = trackParams[0]; + trackParams[1].MinTrackLength = 6; + trackParams[2].MinTrackLength = 4; + LOG(info) << "Initializing tracker in misaligned sync. phase reconstruction with " << trackParams.size() << " passes"; + + } else if (mMode == "sync") { + trackParams.resize(1); + LOG(info) << "Initializing tracker in sync. phase reconstruction with " << trackParams.size() << " passes"; + } else if (mMode == "cosmics") { + mCosmicsProcessing = true; + mRunVertexer = false; + trackParams.resize(1); + trackParams[0].MinTrackLength = 4; + trackParams[0].CellDeltaTanLambdaSigma *= 10; + trackParams[0].PhiBins = 4; + trackParams[0].ZBins = 16; + trackParams[0].PVres = 1.e5f; + trackParams[0].LayerMisalignment[0] = 1.e-2; + trackParams[0].LayerMisalignment[1] = 1.e-2; + trackParams[0].LayerMisalignment[2] = 1.e-2; + trackParams[0].LayerMisalignment[3] = 3.e-2; + trackParams[0].LayerMisalignment[4] = 3.e-2; + trackParams[0].LayerMisalignment[5] = 3.e-2; + trackParams[0].LayerMisalignment[6] = 3.e-2; + trackParams[0].FitIterationMaxChi2[0] = 50.; + trackParams[0].FitIterationMaxChi2[1] = 25.; + trackParams[0].TrackletsPerClusterLimit = 100.; + trackParams[0].CellsPerClusterLimit = 100.; + LOG(info) << "Initializing tracker in reconstruction for cosmics with " << trackParams.size() << " passes"; + + } else { + throw std::runtime_error(fmt::format("Unsupported ITS3 tracking mode {:s} ", mMode)); + } + + for (auto& params : trackParams) { + params.CorrType = o2::base::PropagatorImpl::MatCorrType::USEMatCorrLUT; + } + mTracker->setParameters(trackParams); +} + +void TrackerDPL::run(ProcessingContext& pc) +{ + mTimer.Start(false); + updateTimeDependentParams(pc); + auto compClusters = pc.inputs().get>("compClusters"); + gsl::span patterns = pc.inputs().get>("patterns"); + gsl::span physTriggers; + std::vector fromTRD; + if (mUseTriggers == 2) { // use TRD triggers + o2::InteractionRecord ir{0, pc.services().get().firstTForbit}; + auto trdTriggers = pc.inputs().get>("phystrig"); + for (const auto& trig : trdTriggers) { + if (trig.getBCData() >= ir && trig.getNumberOfTracklets()) { + ir = trig.getBCData(); + fromTRD.emplace_back(o2::itsmft::PhysTrigger{ir, 0}); + } + } + physTriggers = gsl::span(fromTRD.data(), fromTRD.size()); + } else if (mUseTriggers == 1) { // use Phys triggers from ITS stream + physTriggers = pc.inputs().get>("phystrig"); + } + + // code further down does assignment to the rofs and the altered object is used for output + // we therefore need a copy of the vector rather than an object created directly on the input data, + // the output vector however is created directly inside the message memory thus avoiding copy by + // snapshot + auto rofsinput = pc.inputs().get>("ROframes"); + auto& rofs = pc.outputs().make>(Output{"ITS3", "ITS3TrackROF", 0, Lifetime::Timeframe}, rofsinput.begin(), rofsinput.end()); + + auto& irFrames = pc.outputs().make>(Output{"ITS3", "IRFRAMES", 0, Lifetime::Timeframe}); + + const auto& alpParams = o2::itsmft::DPLAlpideParam::Instance(); // RS: this should come from CCDB + int nBCPerTF = alpParams.roFrameLengthInBC; + + LOG(info) << "ITS3Tracker pulled " << compClusters.size() << " clusters, " << rofs.size() << " RO frames"; + + const dataformats::MCTruthContainer* labels = nullptr; + gsl::span mc2rofs; + if (mIsMC) { + labels = pc.inputs().get*>("labels").release(); + // get the array as read-only span, a snapshot is send forward + mc2rofs = pc.inputs().get>("MC2ROframes"); + LOG(info) << labels->getIndexedSize() << " MC label objects , in " << mc2rofs.size() << " MC events"; + } + + std::vector tracks; + auto& allClusIdx = pc.outputs().make>(Output{"ITS3", "TRACKCLSID", 0, Lifetime::Timeframe}); + std::vector trackLabels; + std::vector verticesLabels; + auto& allTracks = pc.outputs().make>(Output{"ITS3", "TRACKS", 0, Lifetime::Timeframe}); + std::vector allTrackLabels; + std::vector allVerticesLabels; + + auto& vertROFvec = pc.outputs().make>(Output{"ITS3", "VERTICESROF", 0, Lifetime::Timeframe}); + auto& vertices = pc.outputs().make>(Output{"ITS3", "VERTICES", 0, Lifetime::Timeframe}); + + std::uint32_t roFrame = 0; + + bool continuous = o2::base::GRPGeomHelper::instance().getGRPECS()->isDetContinuousReadOut(o2::detectors::DetID::IT3); + LOG(info) << "ITS3Tracker RO: continuous=" << continuous; + TimeFrame* timeFrame = mChainITS->GetITSTimeframe(); + mTracker->adoptTimeFrame(*timeFrame); + + mTracker->setBz(o2::base::Propagator::Instance()->getNominalBz()); + mVertexer->adoptTimeFrame(*timeFrame); + gsl::span::iterator pattIt = patterns.begin(); + + gsl::span rofspan(rofs); + timeFrame->loadROFrameDataITS3(rofspan, compClusters, pattIt, mDict, labels); + pattIt = patterns.begin(); + std::vector savedROF; + auto logger = [&](std::string s) { LOG(info) << s; }; + auto errorLogger = [&](std::string s) { LOG(error) << s; }; + + FastMultEst multEst; // mult estimator + std::vector processingMask; + int cutVertexMult{0}, cutRandomMult = int(rofs.size()) - multEst.selectROFs(rofs, compClusters, physTriggers, processingMask); + timeFrame->setMultiplicityCutMask(processingMask); + float vertexerElapsedTime{0.f}; + if (mRunVertexer) { + // Run seeding vertexer + vertexerElapsedTime = mVertexer->clustersToVertices(logger); + } + const auto& multEstConf = FastMultEstConfig::Instance(); // parameters for mult estimation and cuts + for (auto iRof{0}; iRof < rofspan.size(); ++iRof) { + std::vector vtxVecLoc; + auto& vtxROF = vertROFvec.emplace_back(rofspan[iRof]); + vtxROF.setFirstEntry(vertices.size()); + if (mRunVertexer) { + auto vtxSpan = timeFrame->getPrimaryVertices(iRof); + vtxROF.setNEntries(vtxSpan.size()); + bool selROF = vtxSpan.size() == 0; + for (auto iV{0}; iV < vtxSpan.size(); ++iV) { + auto& v = vtxSpan[iV]; + if (multEstConf.isVtxMultCutRequested() && !multEstConf.isPassingVtxMultCut(v.getNContributors())) { + continue; // skip vertex of unwanted multiplicity + } + selROF = true; + vertices.push_back(v); + if (mIsMC) { + auto vLabels = timeFrame->getPrimaryVerticesLabels(iRof)[iV]; + std::copy(vLabels.begin(), vLabels.end(), std::back_inserter(allVerticesLabels)); + } + } + if (processingMask[iRof] && !selROF) { // passed selection in clusters and not in vertex multiplicity + LOG(debug) << fmt::format("ROF {} rejected by the vertex multiplicity selection [{},{}]", + iRof, + multEstConf.cutMultVtxLow, + multEstConf.cutMultVtxHigh); + processingMask[iRof] = selROF; + cutVertexMult++; + } + } else { // cosmics + vtxVecLoc.emplace_back(Vertex()); + vtxVecLoc.back().setNContributors(1); + vtxROF.setNEntries(vtxVecLoc.size()); + for (auto& v : vtxVecLoc) { + vertices.push_back(v); + } + timeFrame->addPrimaryVertices(vtxVecLoc); + } + } + LOG(info) << fmt::format(" - rejected {}/{} ROFs: random/mult.sel:{} (seed {}), vtx.sel:{}", cutRandomMult + cutVertexMult, rofspan.size(), cutRandomMult, multEst.lastRandomSeed, cutVertexMult); + LOG(info) << fmt::format(" - Vertex seeding total elapsed time: {} ms for {} vertices found in {} ROFs", vertexerElapsedTime, timeFrame->getPrimaryVerticesNum(), rofspan.size()); + LOG(info) << fmt::format(" - Beam position computed for the TF: {}, {}", timeFrame->getBeamX(), timeFrame->getBeamY()); + + if (mCosmicsProcessing && compClusters.size() > 1500 * rofspan.size()) { + LOG(error) << "Cosmics processing was requested with an average detector occupancy exceeding 1.e-7, skipping TF processing."; + } else { + + timeFrame->setMultiplicityCutMask(processingMask); + // Run CA tracker + mTracker->clustersToTracks(logger, errorLogger); + if (timeFrame->hasBogusClusters()) { + LOG(warning) << fmt::format(" - The processed timeframe had {} clusters with wild z coordinates, check the dictionaries", timeFrame->hasBogusClusters()); + } + + for (unsigned int iROF{0}; iROF < rofs.size(); ++iROF) { + auto& rof{rofs[iROF]}; + tracks = timeFrame->getTracks(iROF); + trackLabels = timeFrame->getTracksLabel(iROF); + auto number{tracks.size()}; + auto first{allTracks.size()}; + int offset = -rof.getFirstEntry(); // cluster entry!!! + rof.setFirstEntry(first); + rof.setNEntries(number); + + if (processingMask[iROF]) { + irFrames.emplace_back(rof.getBCData(), rof.getBCData() + nBCPerTF - 1).info = tracks.size(); + } + + std::copy(trackLabels.begin(), trackLabels.end(), std::back_inserter(allTrackLabels)); + // Some conversions that needs to be moved in the tracker internals + for (unsigned int iTrk{0}; iTrk < tracks.size(); ++iTrk) { + auto& trc{tracks[iTrk]}; + trc.setFirstClusterEntry(allClusIdx.size()); // before adding tracks, create final cluster indices + int ncl = trc.getNumberOfClusters(), nclf = 0; + for (int ic = TrackITSExt::MaxClusters; ic--;) { // track internally keeps in->out cluster indices, but we want to store the references as out->in!!! + auto clid = trc.getClusterIndex(ic); + if (clid >= 0) { + allClusIdx.push_back(clid); + nclf++; + } + } + assert(ncl == nclf); + allTracks.emplace_back(trc); + } + } + LOGP(info, "ITS3Tracker pushed {} tracks and {} vertices", allTracks.size(), vertices.size()); + if (mIsMC) { + LOGP(info, "ITS3Tracker pushed {} track labels", allTrackLabels.size()); + LOGP(info, "ITS3Tracker pushed {} vertex labels", allVerticesLabels.size()); + + pc.outputs().snapshot(Output{"ITS3", "TRACKSMCTR", 0, Lifetime::Timeframe}, allTrackLabels); + pc.outputs().snapshot(Output{"ITS3", "VERTICESMCTR", 0, Lifetime::Timeframe}, allVerticesLabels); + pc.outputs().snapshot(Output{"ITS3", "ITSTrackMC2ROF", 0, Lifetime::Timeframe}, mc2rofs); + } + } + mTimer.Stop(); +} + +///_______________________________________ +void TrackerDPL::updateTimeDependentParams(ProcessingContext& pc) +{ + o2::base::GRPGeomHelper::instance().checkUpdates(pc); + static bool initOnceDone = false; + if (!initOnceDone) { // this params need to be queried only once + initOnceDone = true; + pc.inputs().get("cldict"); // just to trigger the finaliseCCDB + pc.inputs().get*>("alppar"); + GeometryTGeo* geom = GeometryTGeo::Instance(); + geom->fillMatrixCache(o2::math_utils::bit2Mask(o2::math_utils::TransformType::T2L, o2::math_utils::TransformType::T2GRot, o2::math_utils::TransformType::T2G)); + mVertexer->getGlobalConfiguration(); + mTracker->getGlobalConfiguration(); + } +} + +///_______________________________________ +void TrackerDPL::finaliseCCDB(ConcreteDataMatcher& matcher, void* obj) +{ + if (o2::base::GRPGeomHelper::instance().finaliseCCDB(matcher, obj)) { + return; + } + if (matcher == ConcreteDataMatcher("ITS", "CLUSDICT", 0)) { + LOG(info) << "cluster dictionary updated"; + setClusterDictionary((const o2::itsmft::TopologyDictionary*)obj); + return; + } + // Note: strictly speaking, for Configurable params we don't need finaliseCCDB check, the singletons are updated at the CCDB fetcher level + if (matcher == ConcreteDataMatcher("ITS", "ALPIDEPARAM", 0)) { + LOG(info) << "Alpide param updated"; + const auto& par = o2::itsmft::DPLAlpideParam::Instance(); + par.printKeyValues(); + return; + } +} + +void TrackerDPL::endOfStream(EndOfStreamContext& ec) +{ + LOGF(info, "ITS3 CA-Tracker total timing: Cpu: %.3e Real: %.3e s in %d slots", + mTimer.CpuTime(), mTimer.RealTime(), mTimer.Counter() - 1); +} + +DataProcessorSpec getTrackerSpec(bool useMC, int trgType, const std::string& trModeS, o2::gpu::GPUDataTypes::DeviceType dType) +{ + std::vector inputs; + inputs.emplace_back("compClusters", "ITS3", "COMPCLUSTERS", 0, Lifetime::Timeframe); + inputs.emplace_back("patterns", "ITS3", "PATTERNS", 0, Lifetime::Timeframe); + inputs.emplace_back("ROframes", "ITS3", "CLUSTERSROF", 0, Lifetime::Timeframe); + if (trgType == 1) { + inputs.emplace_back("phystrig", "ITS3", "PHYSTRIG", 0, Lifetime::Timeframe); + } else if (trgType == 2) { + inputs.emplace_back("phystrig", "TRD", "TRKTRGRD", 0, Lifetime::Timeframe); + } + inputs.emplace_back("cldict", "ITS", "CLUSDICT", 0, Lifetime::Condition, ccdbParamSpec("ITS/Calib/ClusterDictionary")); + inputs.emplace_back("alppar", "ITS", "ALPIDEPARAM", 0, Lifetime::Condition, ccdbParamSpec("ITS/Config/AlpideParam")); + auto ggRequest = std::make_shared(false, // orbitResetTime + true, // GRPECS=true + false, // GRPLHCIF + true, // GRPMagField + true, // askMatLUT + o2::base::GRPGeomRequest::Aligned, // geometry + inputs, + true); + std::vector outputs; + outputs.emplace_back("ITS3", "TRACKS", 0, Lifetime::Timeframe); + outputs.emplace_back("ITS3", "TRACKCLSID", 0, Lifetime::Timeframe); + outputs.emplace_back("ITS3", "ITSTrackROF", 0, Lifetime::Timeframe); + outputs.emplace_back("ITS3", "VERTICES", 0, Lifetime::Timeframe); + outputs.emplace_back("ITS3", "VERTICESROF", 0, Lifetime::Timeframe); + outputs.emplace_back("ITS3", "IRFRAMES", 0, Lifetime::Timeframe); + + if (useMC) { + inputs.emplace_back("labels", "ITS3", "CLUSTERSMCTR", 0, Lifetime::Timeframe); + inputs.emplace_back("MC2ROframes", "ITS3", "CLUSTERSMC2ROF", 0, Lifetime::Timeframe); + outputs.emplace_back("ITS3", "VERTICESMCTR", 0, Lifetime::Timeframe); + outputs.emplace_back("ITS3", "TRACKSMCTR", 0, Lifetime::Timeframe); + outputs.emplace_back("ITS3", "ITSTrackMC2ROF", 0, Lifetime::Timeframe); + outputs.emplace_back("ITS3", "VERTICES", 0, Lifetime::Timeframe); + } + + return DataProcessorSpec{ + "its3-tracker", + inputs, + outputs, + AlgorithmSpec{adaptFromTask(ggRequest, useMC, trgType, trModeS, dType)}, + Options{}}; +} + +} // namespace its3 +} // namespace o2 diff --git a/Detectors/Upgrades/IT3/workflow/src/digit-reader-workflow.cxx b/Detectors/Upgrades/ITS3/workflow/src/digit-reader-workflow.cxx similarity index 100% rename from Detectors/Upgrades/IT3/workflow/src/digit-reader-workflow.cxx rename to Detectors/Upgrades/ITS3/workflow/src/digit-reader-workflow.cxx diff --git a/Detectors/Upgrades/IT3/workflow/src/digit-writer-workflow.cxx b/Detectors/Upgrades/ITS3/workflow/src/digit-writer-workflow.cxx similarity index 100% rename from Detectors/Upgrades/IT3/workflow/src/digit-writer-workflow.cxx rename to Detectors/Upgrades/ITS3/workflow/src/digit-writer-workflow.cxx diff --git a/Detectors/Upgrades/IT3/workflow/src/its3-reco-workflow.cxx b/Detectors/Upgrades/ITS3/workflow/src/its3-reco-workflow.cxx similarity index 82% rename from Detectors/Upgrades/IT3/workflow/src/its3-reco-workflow.cxx rename to Detectors/Upgrades/ITS3/workflow/src/its3-reco-workflow.cxx index b202cf67b81f1..5e8c44808e8a0 100644 --- a/Detectors/Upgrades/IT3/workflow/src/its3-reco-workflow.cxx +++ b/Detectors/Upgrades/ITS3/workflow/src/its3-reco-workflow.cxx @@ -44,7 +44,9 @@ void customize(std::vector& workflowOptions) {"clusters-from-upstream", o2::framework::VariantType::Bool, false, {"clusters will be provided from upstream, skip clusterizer"}}, {"disable-root-output", o2::framework::VariantType::Bool, false, {"do not write output root files"}}, {"disable-mc", o2::framework::VariantType::Bool, false, {"disable MC propagation even if available"}}, + {"select-with-triggers", o2::framework::VariantType::String, "none", {"use triggers to prescale processed ROFs: phys, trd, none"}}, {"tracking-mode", o2::framework::VariantType::String, "sync", {"sync,async,cosmics"}}, + {"disable-tracking", o2::framework::VariantType::Bool, false, {"disable tracking step"}}, {"entropy-encoding", o2::framework::VariantType::Bool, false, {"produce entropy encoded data"}}, {"configKeyValues", VariantType::String, "", {"Semicolon separated key=value strings"}}, {"gpuDevice", o2::framework::VariantType::Int, 1, {"use gpu device: CPU=1,CUDA=2,HIP=3 (default: CPU)"}}}; @@ -62,15 +64,28 @@ WorkflowSpec defineDataProcessing(ConfigContext const& configcontext) // Update the (declared) parameters if changed from the command line auto useMC = !configcontext.options().get("disable-mc"); auto trmode = configcontext.options().get("tracking-mode"); + auto selTrig = configcontext.options().get("select-with-triggers"); auto gpuDevice = static_cast(configcontext.options().get("gpuDevice")); auto extDigits = configcontext.options().get("digits-from-upstream"); auto extClusters = configcontext.options().get("clusters-from-upstream"); auto disableRootOutput = configcontext.options().get("disable-root-output"); - auto eencode = configcontext.options().get("entropy-encoding"); + if (configcontext.options().get("disable-tracking")) { + trmode = ""; + } + std::transform(trmode.begin(), trmode.end(), trmode.begin(), [](unsigned char c) { return std::tolower(c); }); o2::conf::ConfigurableParam::updateFromString(configcontext.options().get("configKeyValues")); - - auto wf = o2::its3::reco_workflow::getWorkflow(useMC, trmode, gpuDevice, extDigits, extClusters, disableRootOutput, eencode); + int trType = 0; + if (!selTrig.empty() && selTrig != "none") { + if (selTrig == "phys") { + trType = 1; + } else if (selTrig == "trd") { + trType = 2; + } else { + LOG(fatal) << "Unknown trigger type requested for events prescaling: " << selTrig; + } + } + auto wf = o2::its3::reco_workflow::getWorkflow(useMC, trmode, gpuDevice, extDigits, extClusters, disableRootOutput, trType); // configure dpl timer to inject correct firstTForbit: start from the 1st orbit of TF containing 1st sampled orbit o2::raw::HBFUtilsInitializer hbfIni(configcontext, wf); diff --git a/Steer/DigitizerWorkflow/CMakeLists.txt b/Steer/DigitizerWorkflow/CMakeLists.txt index 6f9a8d314b7ad..7ac1701c9764e 100644 --- a/Steer/DigitizerWorkflow/CMakeLists.txt +++ b/Steer/DigitizerWorkflow/CMakeLists.txt @@ -19,7 +19,7 @@ o2_add_executable(digitizer-workflow src/GRPUpdaterSpec.cxx src/HMPIDDigitizerSpec.cxx src/ITSMFTDigitizerSpec.cxx - # src/ITS3DigitizerSpec.cxx + src/ITS3DigitizerSpec.cxx src/MCHDigitizerSpec.cxx src/MIDDigitizerSpec.cxx src/PHOSDigitizerSpec.cxx @@ -64,8 +64,8 @@ o2_add_executable(digitizer-workflow O2::ZDCSimulation O2::ZDCWorkflow O2::DetectorsRaw - # O2::ITS3Simulation - # O2::ITS3Workflow + O2::ITS3Simulation + O2::ITS3Workflow ) else() o2_add_executable(digitizer-workflow diff --git a/Steer/DigitizerWorkflow/src/ITS3DigitizerSpec.cxx b/Steer/DigitizerWorkflow/src/ITS3DigitizerSpec.cxx index 5a8af269ff3a9..9f03bd6733729 100644 --- a/Steer/DigitizerWorkflow/src/ITS3DigitizerSpec.cxx +++ b/Steer/DigitizerWorkflow/src/ITS3DigitizerSpec.cxx @@ -27,7 +27,7 @@ #include "ITS3Simulation/Digitizer.h" #include "ITSMFTSimulation/DPLDigitizerParam.h" #include "ITSMFTBase/DPLAlpideParam.h" -#include "ITS3Base/GeometryTGeo.h" +#include "ITSBase/GeometryTGeo.h" #include #include #include @@ -81,7 +81,7 @@ class ITS3DPLDigitizerTask : BaseDPLDigitizer << " RO mode"; // configure digitizer - o2::its3::GeometryTGeo* geom = o2::its3::GeometryTGeo::Instance(); + o2::its::GeometryTGeo* geom = o2::its::GeometryTGeo::Instance(); geom->fillMatrixCache(o2::math_utils::bit2Mask(o2::math_utils::TransformType::L2G)); // make sure L2G matrices are loaded mDigitizer.setGeometry(geom); diff --git a/Steer/DigitizerWorkflow/src/SimpleDigitizerWorkflow.cxx b/Steer/DigitizerWorkflow/src/SimpleDigitizerWorkflow.cxx index d44649bb7141a..28cf109ca483a 100644 --- a/Steer/DigitizerWorkflow/src/SimpleDigitizerWorkflow.cxx +++ b/Steer/DigitizerWorkflow/src/SimpleDigitizerWorkflow.cxx @@ -40,11 +40,11 @@ #include "ITSMFTDigitizerSpec.h" #include "ITSMFTWorkflow/DigitWriterSpec.h" -// #ifdef ENABLE_UPGRADES -// // for ITS3 -// #include "ITS3DigitizerSpec.h" -// #include "ITS3Workflow/DigitWriterSpec.h" -// #endif +#ifdef ENABLE_UPGRADES +// for ITS3 +#include "ITS3DigitizerSpec.h" +#include "ITS3Workflow/DigitWriterSpec.h" +#endif // for TOF #include "TOFDigitizerSpec.h" @@ -590,16 +590,16 @@ WorkflowSpec defineDataProcessing(ConfigContext const& configcontext) writerSpecs.emplace_back(o2::itsmft::getITSDigitWriterSpec(mctruth)); } - // #ifdef ENABLE_UPGRADES - // // the ITS3 part - // if (isEnabled(o2::detectors::DetID::IT3)) { - // detList.emplace_back(o2::detectors::DetID::IT3); - // // connect the ITS digitization - // specs.emplace_back(o2::its3::getITS3DigitizerSpec(fanoutsize++, mctruth)); - // // // connect ITS digit writer - // specs.emplace_back(o2::its3::getITS3DigitWriterSpec(mctruth)); - // } - // #endif +#ifdef ENABLE_UPGRADES + // the ITS3 part + if (isEnabled(o2::detectors::DetID::IT3)) { + detList.emplace_back(o2::detectors::DetID::IT3); + // connect the ITS digitization + specs.emplace_back(o2::its3::getITS3DigitizerSpec(fanoutsize++, mctruth)); + // // connect ITS digit writer + specs.emplace_back(o2::its3::getITS3DigitWriterSpec(mctruth)); + } +#endif // the MFT part if (isEnabled(o2::detectors::DetID::MFT)) {