Project
Loading...
Searching...
No Matches
SACCCDBHelper.cxx
Go to the documentation of this file.
1// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
2// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
3// All rights not expressly granted are reserved.
4//
5// This software is distributed under the terms of the GNU General Public
6// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
7//
8// In applying this license CERN does not waive the privileges and immunities
9// granted to it by virtue of its status as an Intergovernmental Organization
10// or submit itself to any jurisdiction.
11
17
18template <typename DataT>
20{
21 return mSACDelta ? (mSACDelta->mSACDelta[side].getNIDCs()) / (GEMSTACKSPERSECTOR * SECTORSPERSIDE) : 0;
22}
23
24template <typename DataT>
26{
27 return mSACOne ? mSACOne->mSACOne[side].getNIDCs() : 0;
28}
29
30template <typename DataT>
31float o2::tpc::SACCCDBHelper<DataT>::getSACZeroVal(const unsigned int sector, const unsigned int stack) const
32{
33 return !mSACZero ? -1 : mSACZero->getValueIDCZero(Sector(sector).side(), SACFactorization::getStackInSide(sector, stack));
34}
35
36template <typename DataT>
37float o2::tpc::SACCCDBHelper<DataT>::getSACDeltaVal(const unsigned int sector, const unsigned int stack, unsigned int integrationInterval) const
38{
39 return !mSACDelta ? -1 : mSACDelta->getValue(Sector(sector).side(), SACFactorization::getSACDeltaIndex(SACFactorization::getStackInSide(sector, stack), integrationInterval));
40}
41
42template <typename DataT>
43float o2::tpc::SACCCDBHelper<DataT>::getSACOneVal(const o2::tpc::Side side, const unsigned int integrationInterval) const
44{
45 return !mSACOne ? -1 : mSACOne->getValue(side, integrationInterval);
46}
47
48template <typename DataT>
49float o2::tpc::SACCCDBHelper<DataT>::getSACVal(const unsigned int sector, const unsigned int stack, unsigned int integrationInterval) const
50{
51 return (getSACDeltaVal(sector, stack, integrationInterval) + 1.f) * getSACZeroVal(sector, stack) * getSACOneVal(Sector(sector).side(), integrationInterval);
53
54template <typename DataT>
55void o2::tpc::SACCCDBHelper<DataT>::drawSACZeroHelper(const bool type, const o2::tpc::Sector sector, const std::string filename, const float minZ, const float maxZ) const
56{
57 std::function<float(const unsigned int, const unsigned int)> SACFunc = [this](const unsigned int sector, const unsigned int stack) {
58 return this->getSACZeroVal(sector, stack);
59 };
62 drawFun.mSACFunc = SACFunc;
63 const std::string zAxisTitle = SACDrawHelper::getZAxisTitle(SACType::IDCZero);
64 type ? SACDrawHelper::drawSide(drawFun, sector.side(), zAxisTitle, filename, minZ, maxZ) : SACDrawHelper::drawSector(drawFun, sector, zAxisTitle, filename, minZ, maxZ);
65}
67template <typename DataT>
68void o2::tpc::SACCCDBHelper<DataT>::drawSACDeltaHelper(const bool type, const Sector sector, const unsigned int integrationInterval, const std::string filename, const float minZ, const float maxZ) const
69{
70 std::function<float(const unsigned int, const unsigned int)> SACFunc = [this, integrationInterval](const unsigned int sector, const unsigned int stack) {
71 return this->getSACDeltaVal(sector, stack, integrationInterval);
72 };
73
75 drawFun.mSACFunc = SACFunc;
76 const std::string zAxisTitle = SACDrawHelper::getZAxisTitle(SACType::IDCDelta);
77 type ? SACDrawHelper::drawSide(drawFun, sector.side(), zAxisTitle, filename, minZ, maxZ) : SACDrawHelper::drawSector(drawFun, sector, zAxisTitle, filename, minZ, maxZ);
78}
79
80template <typename DataT>
81void o2::tpc::SACCCDBHelper<DataT>::drawSACHelper(const bool type, const Sector sector, const unsigned int integrationInterval, const std::string filename, const float minZ, const float maxZ) const
82{
83 std::function<float(const unsigned int, const unsigned int)> SACFunc = [this, integrationInterval](const unsigned int sector, const unsigned int stack) {
84 return this->getSACVal(sector, stack, integrationInterval);
85 };
86
87 SACDrawHelper::SACDraw drawFun;
88 drawFun.mSACFunc = SACFunc;
89 const std::string zAxisTitle = SACDrawHelper::getZAxisTitle(SACType::IDC);
90 type ? SACDrawHelper::drawSide(drawFun, sector.side(), zAxisTitle, filename, minZ, maxZ) : SACDrawHelper::drawSector(drawFun, sector, zAxisTitle, filename, minZ, maxZ);
91}
92
93template <typename DataT>
94void o2::tpc::SACCCDBHelper<DataT>::dumpToTree(const char* outFileName) const
95{
96 o2::utils::TreeStreamRedirector pcstream(outFileName, "RECREATE");
97 pcstream.GetFile()->cd();
98
99 const auto intervals = std::min(getNIntegrationIntervalsSACDelta(Side::A), getNIntegrationIntervalsSACDelta(Side::C));
100 for (unsigned int integrationInterval = 0; integrationInterval < intervals; ++integrationInterval) {
101 std::vector<int32_t> vSACs(GEMSTACKS);
102 std::vector<float> vSACsZero(GEMSTACKS);
103 std::vector<float> vSACsDelta(GEMSTACKS);
104 std::vector<unsigned int> vStack(GEMSTACKS);
105
106 unsigned int index = 0;
107 for (unsigned int sector = 0; sector < Mapper::NSECTORS; ++sector) {
108 for (unsigned int stack = 0; stack < GEMSTACKSPERSECTOR; ++stack) {
109 vSACs[index] = getSACVal(sector, stack, integrationInterval);
110 vSACsZero[index] = getSACZeroVal(sector, stack);
111 vSACsDelta[index] = getSACDeltaVal(sector, stack, integrationInterval);
112 vStack[index++] = SACFactorization::getStack(sector, stack);
113 }
114 }
115 float sacOneA = getSACOneVal(Side::A, integrationInterval);
116 float sacOneC = getSACOneVal(Side::C, integrationInterval);
117
118 pcstream << "tree"
119 << "integrationInterval=" << integrationInterval
120 << "SAC.=" << vSACs
121 << "SAC0.=" << vSACsZero
122 << "SAC1A=" << sacOneA
123 << "SAC1C=" << sacOneC
124 << "SACDelta.=" << vSACsDelta
125 << "stack.=" << vStack
126 << "\n";
127 }
128 pcstream.Close();
129}
130
131template <typename DataT>
133{
134 o2::utils::TreeStreamRedirector pcstream("fourierCoeff.root", "RECREATE");
135 pcstream.GetFile()->cd();
136
137 for (int iside = 0; iside < SIDES; ++iside) {
138 const Side side = (iside == 0) ? Side::A : Side::C;
139
140 if (!mFourierCoeff) {
141 continue;
142 }
143
144 const int nTFs = mFourierCoeff->mCoeff[side].getNCoefficients() / mFourierCoeff->mCoeff[side].getNCoefficientsPerTF();
145 for (int iTF = 0; iTF < nTFs; ++iTF) {
146 std::vector<float> coeff;
147 std::vector<int> ind;
148 int coeffPerTF = mFourierCoeff->mCoeff[side].getNCoefficientsPerTF();
149 for (int i = 0; i < coeffPerTF; ++i) {
150 const int index = mFourierCoeff->mCoeff[side].getIndex(iTF, i);
151 coeff.emplace_back((mFourierCoeff->mCoeff[side])(index));
152 ind.emplace_back(i);
153 }
154
155 pcstream << "tree"
156 << "iTF=" << iTF
157 << "index=" << ind
158 << "coeffPerTF=" << coeffPerTF
159 << "coeff.=" << coeff
160 << "side=" << iside
161 << "\n";
162 }
163 }
164 pcstream.Close();
165}
166
167template class o2::tpc::SACCCDBHelper<float>;
int32_t i
This file provides the structs for storing the factorized IDC values and fourier coefficients to be s...
uint32_t side
Definition RawData.h:0
uint32_t stack
Definition RawData.h:1
helper class for accessing SACs from CCDB
helper class for drawing SACs per sector/side
TPC factorization of SACs.
void dumpToTree(const char *outFileName="SACCCDBTree.root") const
unsigned int getNIntegrationIntervalsSACDelta(const o2::tpc::Side side) const
float getSACVal(const unsigned int sector, const unsigned int stack, unsigned int integrationInterval) const
float getSACOneVal(const o2::tpc::Side side, const unsigned int integrationInterval) const
void dumpToFourierCoeffToTree(const char *outFileName="FourierCCDBTree.root") const
float getSACDeltaVal(const unsigned int sector, const unsigned int stack, unsigned int integrationInterval) const
unsigned int getNIntegrationIntervalsSACOne(const o2::tpc::Side side) const
float getSACZeroVal(const unsigned int sector, const unsigned int stack) const
Side side() const
Definition Sector.h:96
GLuint index
Definition glcorearb.h:781
GLint GLint GLsizei GLint GLenum GLenum type
Definition glcorearb.h:275
GLfloat GLfloat GLfloat GLfloat GLfloat GLfloat maxZ
Definition glcorearb.h:2910
GLfloat GLfloat minZ
Definition glcorearb.h:2910
constexpr unsigned char SECTORSPERSIDE
Definition Defs.h:40
constexpr unsigned short GEMSTACKS
Definition Defs.h:60
constexpr unsigned char SIDES
Definition Defs.h:41
Side
TPC readout sidE.
Definition Defs.h:35
constexpr unsigned short GEMSTACKSPERSECTOR
Definition Defs.h:57
std::string filename()
std::function< float(const unsigned int, const unsigned int)> mSACFunc
function returning the value which will be drawn for sector, stack