Project
Loading...
Searching...
No Matches
DataContainer3D.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
16#include "TPCBase/Mapper.h"
17#include "Framework/Logger.h"
18#include "TFile.h"
19#include "ROOT/RDataFrame.hxx"
20#include "TStopwatch.h"
21#include "TTree.h"
22
23#include <memory>
24#include <iomanip>
25#include <algorithm>
26
27using namespace o2::tpc;
28
29template <typename DataT>
30template <typename DataTOut>
31int DataContainer3D<DataT>::writeToFile(TFile& outf, const char* name) const
32{
33 if (outf.IsZombie()) {
34 LOGP(error, "Failed to write to file: {}", outf.GetName());
35 return -1;
36 }
37
38 DataContainer3D<DataTOut> containerTmp(mZVertices, mRVertices, mPhiVertices);
39 containerTmp.getData() = std::vector<DataTOut>(mData.begin(), mData.end());
40
41 outf.WriteObjectAny(&containerTmp, DataContainer3D<DataTOut>::Class(), name);
42 return 0;
43}
44
45template <typename DataT>
46int DataContainer3D<DataT>::writeToFile(std::string_view file, std::string_view option, std::string_view name, const int nthreads) const
47{
48 // max number of floats per Entry
49 const size_t maxvalues = sizeof(float) * 1024 * 1024;
50
51 // total number of values to be stored
52 const size_t nsize = getNDataPoints();
53
54 // calculate number of entries in the tree and restrict if the number of values per threads exceeds max size
55 size_t entries = ((nsize / nthreads) > maxvalues) ? (nsize / maxvalues) : nthreads;
56
57 if (entries > nsize) {
58 entries = nsize;
59 }
60
61 // calculate numbers to store per entry
62 const size_t values_per_entry = nsize / entries;
63
64 // in case of remainder add additonal entry
65 const size_t values_lastEntry = nsize % entries;
66 if (values_lastEntry) {
67 entries += 1;
68 }
69
70 // in case EnableImplicitMT was already called with different number of threads, perform reset
71 if (ROOT::IsImplicitMTEnabled() && (ROOT::GetThreadPoolSize() != nthreads)) {
72 ROOT::DisableImplicitMT();
73 }
74 ROOT::EnableImplicitMT(nthreads);
75
76 // define dataframe which will be stored in the TTree
77 ROOT::RDataFrame dFrame(entries);
78
79 // define function which is used to fill the data frame
80 auto dfStore = dFrame.DefineSlotEntry(name, [&data = std::as_const(mData), entries, values_per_entry](unsigned int, ULong64_t entry) { return DataContainer3D<DataT>::getDataSlice(data, entries, values_per_entry, entry); });
81 dfStore = dfStore.Define("nz", [mZVertices = mZVertices]() { return mZVertices; });
82 dfStore = dfStore.Define("nr", [mRVertices = mRVertices]() { return mRVertices; });
83 dfStore = dfStore.Define("nphi", [mPhiVertices = mPhiVertices]() { return mPhiVertices; });
84
85 // define options of TFile
86 ROOT::RDF::RSnapshotOptions opt;
87 opt.fMode = option;
88 opt.fOverwriteIfExists = true; // overwrite if already exists
89
90 TStopwatch timer;
91 // note: first call has some overhead (~2s)
92 dfStore.Snapshot(name, file, {name.data(), "nz", "nr", "nphi"}, opt);
93 timer.Print("u");
94 return 0;
95}
96
97template <typename DataT>
98bool DataContainer3D<DataT>::initFromFile(std::string_view file, std::string_view name, const int nthreads)
99{
100 // in case EnableImplicitMT was already called with different number of threads, perform reset
101 if (ROOT::IsImplicitMTEnabled() && (ROOT::GetThreadPoolSize() != nthreads)) {
102 ROOT::DisableImplicitMT();
103 }
104 ROOT::EnableImplicitMT(nthreads);
105
106 // compare first the meta data (is the number of vertices the same)
107 // define data frame from imput file
108 ROOT::RDataFrame dFrame(name, file);
109
110 // compare vertices
111 auto comp = [mZVertices = mZVertices, mRVertices = mRVertices, mPhiVertices = mPhiVertices](const unsigned short nz, const unsigned short nr, const unsigned short nphi) {
112 if ((nz == mZVertices) && (nr == mRVertices) && (nphi == mPhiVertices)) {
113 return false;
114 }
115 return true;
116 };
117
118 auto count = dFrame.Filter(comp, {"nz", "nr", "nphi"}).Count();
119 if (*count != 0) {
120 LOGP(error, "Data from input file has different number of vertices! Found {} same vertices", *count);
121 return false;
122 }
123
124 // define lambda function which is used to copy the data
125 auto readData = [&mData = mData](const std::pair<long, std::vector<float>>& data) {
126 std::copy(data.second.begin(), data.second.end(), mData.begin() + data.first);
127 };
128
129 LOGP(info, "Reading {} from file {}", name, file);
130
131 // fill data from RDataFrame
132 TStopwatch timer;
133 dFrame.Foreach(readData, {name.data()});
134 timer.Print("u");
135 return true;
136}
137
139template <typename DataT>
140template <typename DataTIn>
141bool DataContainer3D<DataT>::initFromFile(TFile& inpf, const char* name)
142{
143 if (inpf.IsZombie()) {
144 LOGP(error, "Failed to read from file: {}", inpf.GetName());
145 return false;
146 }
147 DataContainer3D<DataTIn>* dataCont{nullptr};
148 dataCont = reinterpret_cast<DataContainer3D<DataTIn>*>(inpf.GetObjectChecked(name, DataContainer3D<DataTIn>::Class()));
150 if (!dataCont) {
151 LOGP(error, "Failed to load {} from {}", name, inpf.GetName());
152 return false;
153 }
154
155 if (mZVertices != dataCont->getNZ() || mRVertices != dataCont->getNR() || mPhiVertices != dataCont->getNPhi()) {
156 LOGP(error, "Data from input file has different definition of vertices!");
157 LOGP(error, "set vertices before creating the sc object to: SpaceCharge<>::setGrid({}, {}, {})", dataCont->getNZ(), dataCont->getNR(), dataCont->getNPhi());
158 delete dataCont;
159 return false;
160 }
161
162 mData = std::vector<DataT>(dataCont->getData().begin(), dataCont->getData().end());
163 delete dataCont;
164 return true;
165}
166
167template <typename DataT>
169{
170 if (inpf.IsZombie()) {
171 LOGP(error, "Failed to read from file {}", inpf.GetName());
172 return nullptr;
173 }
174 DataContainer3D<DataT>* dataCont{nullptr};
176 dataCont = reinterpret_cast<DataContainer3D<DataT>*>(inpf.GetObjectChecked(name, DataContainer3D<DataT>::Class()));
177 if (!dataCont) {
178 LOGP(error, "Failed to load {} from {}", name, inpf.GetName());
179 return nullptr;
180 }
181 return dataCont;
182}
183
184template <typename DataT>
186{
187 std::stringstream stream;
188 stream.precision(3);
189 auto&& w = std::setw(9);
190 stream << std::endl;
192 for (unsigned int iz = 0; iz < mPhiVertices; ++iz) {
193 stream << "z layer: " << iz << "\n";
194 // print top x row
195 stream << "⎡" << w << (*this)(0, 0, iz);
196 for (unsigned int ix = 1; ix < mZVertices; ++ix) {
197 stream << ", " << w << (*this)(ix, 0, iz);
198 }
199 stream << " ⎤ \n";
200
201 for (unsigned int iy = 1; iy < mRVertices - 1; ++iy) {
202 stream << "⎢" << w << (*this)(0, iy, iz);
203 for (unsigned int ix = 1; ix < mZVertices; ++ix) {
204 stream << ", " << w << (*this)(ix, iy, iz);
205 }
206 stream << " ⎥ \n";
207 }
208
209 stream << "⎣" << w << (*this)(0, mRVertices - 1, iz);
210 for (unsigned int ix = 1; ix < mZVertices; ++ix) {
211 stream << ", " << w << (*this)(ix, mRVertices - 1, iz);
212 }
213 stream << " ⎦ \n \n";
214 }
215 LOGP(info, "{} \n \n", stream.str());
216}
217
218template <typename DataT>
219auto DataContainer3D<DataT>::getDataSlice(const std::vector<DataT>& data, size_t entries, const size_t values_per_entry, ULong64_t entry)
220{
221 const long indStart = entry * values_per_entry;
222 if (entry < (entries - 1)) {
223 return std::pair(indStart, std::vector<float>(data.begin() + indStart, data.begin() + indStart + values_per_entry));
224 } else if (entry == (entries - 1)) {
225 // last entry might have different number of values. just copy the rest...
226 return std::pair(indStart, std::vector<float>(data.begin() + indStart, data.end()));
227 }
228 return std::pair(indStart, std::vector<float>());
229};
230
231template <typename DataT>
233{
234 std::transform(mData.begin(), mData.end(), mData.begin(), [value = value](auto& val) { return val * value; });
235 return *this;
236}
237
238template <typename DataT>
240{
241 std::transform(mData.begin(), mData.end(), other.mData.begin(), mData.begin(), std::plus<>());
242 return *this;
243}
244
245template <typename DataT>
247{
248 std::transform(mData.begin(), mData.end(), other.mData.begin(), mData.begin(), std::minus<>());
249 return *this;
250}
251
252template <typename DataT>
254{
255 std::transform(mData.begin(), mData.end(), other.mData.begin(), mData.begin(), std::multiplies<>());
256 return *this;
257}
258
259template <typename DataT>
260size_t DataContainer3D<DataT>::getIndexZ(size_t index, const int nz, const int nr, const int nphi)
261{
262 const size_t iphi = index / (nz * nr);
263 index -= (iphi * nz * nr);
264 const size_t iz = index % nz;
265 return iz;
266}
267
268template <typename DataT>
269size_t DataContainer3D<DataT>::getIndexR(size_t index, const int nz, const int nr, const int nphi)
270{
271 const size_t iphi = index / (nz * nr);
272 index -= (iphi * nz * nr);
273 const size_t ir = index / nz;
274 return ir;
275}
276
277template <typename DataT>
278size_t DataContainer3D<DataT>::getIndexPhi(size_t index, const int nz, const int nr, const int nphi)
279{
280 return index / (nz * nr);
281}
282
283template <typename DataT>
285{
286 tree->SetAlias("ir", "o2::tpc::DataContainer3D<float>::getIndexR(first + Iteration$, nz, nr, nphi)");
287 tree->SetAlias("iz", "o2::tpc::DataContainer3D<float>::getIndexZ(first + Iteration$, nz, nr, nphi)");
288 tree->SetAlias("iphi", "o2::tpc::DataContainer3D<float>::getIndexPhi(first + Iteration$, nz, nr, nphi)");
289 tree->SetAlias("r", "o2::tpc::GridProperties<float>::getRMin() + o2::tpc::GridProperties<float>::getGridSpacingR(nr) * ir");
290 tree->SetAlias("z", "o2::tpc::GridProperties<float>::getZMin() + o2::tpc::GridProperties<float>::getGridSpacingZ(nz) * iz");
291 tree->SetAlias("phi", "o2::tpc::GridProperties<float>::getPhiMin() + o2::tpc::GridProperties<float>::getGridSpacingPhi(nphi) * iphi");
292}
293
294template <typename DataT>
296{
297 // actuall stored value
298 tree->SetAlias("val", "_0");
299
300 // some meta data
301 tree->SetAlias("iz", "_1");
302 tree->SetAlias("ir", "_2");
303 tree->SetAlias("iphi", "_3");
304 tree->SetAlias("z", "_4");
305 tree->SetAlias("r", "_5");
306 tree->SetAlias("phi", "_6");
307 tree->SetAlias("lpos", "_7");
308 tree->SetAlias("lx", "lpos.fCoordinates.fX");
309 tree->SetAlias("ly", "lpos.fCoordinates.fY");
310 tree->SetAlias("index", "_8");
311}
312
313template <typename DataT>
314void DataContainer3D<DataT>::setGrid(unsigned short nZ, unsigned short nR, unsigned short nPhi, const bool resize)
315{
316 mZVertices = nZ;
317 mRVertices = nR;
318 mPhiVertices = nPhi;
319 if (resize) {
320 mData.resize(nZ * nR * static_cast<size_t>(nPhi));
321 }
322}
323
324template <typename DataT>
325void DataContainer3D<DataT>::dumpSlice(std::string_view treename, std::string_view fileIn, std::string_view fileOut, std::string_view option, std::pair<unsigned short, unsigned short> rangeiR, std::pair<unsigned short, unsigned short> rangeiZ, std::pair<unsigned short, unsigned short> rangeiPhi, const int nthreads)
326{
327 if (ROOT::IsImplicitMTEnabled() && (ROOT::GetThreadPoolSize() != nthreads)) {
328 ROOT::DisableImplicitMT();
329 }
330 ROOT::EnableImplicitMT(nthreads);
331 ROOT::RDataFrame dFrame(treename, fileIn);
332
333 auto df = dFrame.Define("slice", [rangeiZ, rangeiR, rangeiPhi](const std::pair<long, std::vector<float>>& values, unsigned short nz, unsigned short nr, unsigned short nphi) {
334 const bool simOneSectorOnly = MGParameters::normalizeGridToOneSector;
335 std::vector<size_t> ir;
336 std::vector<size_t> iphi;
337 std::vector<size_t> iz;
338 std::vector<float> r;
339 std::vector<float> phi;
340 std::vector<float> z;
341 std::vector<float> vals;
342 std::vector<size_t> globalIdx;
343 std::vector<LocalPosition3D> lPos;
344 const auto nvalues = values.second.size();
345 ir.reserve(nvalues);
346 iphi.reserve(nvalues);
347 iz.reserve(nvalues);
348 r.reserve(nvalues);
349 phi.reserve(nvalues);
350 z.reserve(nvalues);
351 vals.reserve(nvalues);
352 lPos.reserve(nvalues);
353 globalIdx.reserve(nvalues);
354 for (size_t i = 0; i < nvalues; ++i) {
355 const size_t idx = values.first + i;
356 const auto iZTmp = o2::tpc::DataContainer3D<float>::getIndexZ(idx, nz, nr, nphi);
357 if ((rangeiZ.first < rangeiZ.second) && ((iZTmp < rangeiZ.first) || (iZTmp > rangeiZ.second))) {
358 continue;
359 }
360
361 const auto iRTmp = o2::tpc::DataContainer3D<float>::getIndexR(idx, nz, nr, nphi);
362 if ((rangeiR.first < rangeiR.second) && ((iRTmp < rangeiR.first) || (iRTmp > rangeiR.second))) {
363 continue;
364 }
365
366 const auto iPhiTmp = o2::tpc::DataContainer3D<float>::getIndexPhi(idx, nz, nr, nphi);
367 if ((rangeiPhi.first < rangeiPhi.second) && ((iPhiTmp < rangeiPhi.first) || (iPhiTmp > rangeiPhi.second))) {
368 continue;
369 }
370
373 const float phiTmp = o2::tpc::GridProperties<float>::getPhiMin() + o2::tpc::GridProperties<float>::getGridSpacingPhi(nphi) / (simOneSectorOnly ? SECTORSPERSIDE : 1) * iPhiTmp;
374
375 const float x = rTmp * std::cos(phiTmp);
376 const float y = rTmp * std::sin(phiTmp);
377 const LocalPosition3D pos(x, y, zTmp);
378 unsigned char secNum = simOneSectorOnly ? 0 : std::floor(phiTmp / SECPHIWIDTH);
379 Sector sector(secNum + (pos.Z() < 0) * SECTORSPERSIDE);
380 LocalPosition3D lPosTmp = Mapper::GlobalToLocal(pos, sector);
381
382 lPos.emplace_back(lPosTmp);
383 ir.emplace_back(iRTmp);
384 iphi.emplace_back(iPhiTmp);
385 iz.emplace_back(iZTmp);
386 r.emplace_back(rTmp);
387 phi.emplace_back(phiTmp);
388 z.emplace_back(zTmp);
389 vals.emplace_back(values.second[i]);
390 globalIdx.emplace_back(idx);
391 }
392 return std::make_tuple(vals, iz, ir, iphi, z, r, phi, lPos, globalIdx);
393 },
394 {treename.data(), "nz", "nr", "nphi"});
395
396 // define options of TFile
397 ROOT::RDF::RSnapshotOptions opt;
398 opt.fMode = option;
399 df.Snapshot(treename, fileOut, {"slice"}, opt);
400}
401
402template <typename DataT>
404{
405 TriCubicInterpolator<DataT> interpolator(*this, grid);
406 return interpolator(z, r, phi);
407}
408
409template <typename DataT>
410void DataContainer3D<DataT>::dumpInterpolation(std::string_view treename, std::string_view fileIn, std::string_view fileOut, std::string_view option, std::pair<float, float> rangeR, std::pair<float, float> rangeZ, std::pair<float, float> rangePhi, const int nR, const int nZ, const int nPhi, const int nthreads)
411{
412 if (ROOT::IsImplicitMTEnabled() && (ROOT::GetThreadPoolSize() != nthreads)) {
413 ROOT::DisableImplicitMT();
414 }
415 ROOT::EnableImplicitMT(nthreads);
416 ROOT::RDataFrame dFrame(nPhi);
417
418 // get vertices for input TTree which is needed to define the grid for interpolation
419 unsigned short nr, nz, nphi;
420 if (!getVertices(treename, fileIn, nr, nz, nphi)) {
421 return;
422 }
423
424 // load data from input TTree
426 data.setGrid(nz, nr, nphi, true);
427 data.initFromFile(fileIn, treename, nthreads);
428
429 // define grid for interpolation
430 using GridProp = GridProperties<DataT>;
431 const RegularGrid3D<DataT> mGrid3D(GridProp::ZMIN, GridProp::RMIN, GridProp::PHIMIN, GridProp::getGridSpacingZ(nz), GridProp::getGridSpacingR(nr), o2::tpc::GridProperties<float>::getGridSpacingPhi(nphi) / (MGParameters::normalizeGridToOneSector ? SECTORSPERSIDE : 1), ParamSpaceCharge{nr, nz, nphi});
432
433 auto interpolate = [&mGrid3D = std::as_const(mGrid3D), &data = std::as_const(data), rangeR, rangeZ, rangePhi, nR, nZ, nPhi](unsigned int, ULong64_t iPhi) {
434 const bool simOneSectorOnly = MGParameters::normalizeGridToOneSector;
435 std::vector<size_t> ir;
436 std::vector<size_t> iphi;
437 std::vector<size_t> iz;
438 std::vector<float> r;
439 std::vector<float> phi;
440 std::vector<float> z;
441 std::vector<float> vals;
442 std::vector<size_t> globalIdx;
443 std::vector<LocalPosition3D> lPos;
444 const auto nvalues = nR * nZ;
445 ir.reserve(nvalues);
446 iphi.reserve(nvalues);
447 iz.reserve(nvalues);
448 r.reserve(nvalues);
449 phi.reserve(nvalues);
450 z.reserve(nvalues);
451 vals.reserve(nvalues);
452 lPos.reserve(nvalues);
453 globalIdx.reserve(nvalues);
454
455 const float rSpacing = (rangeR.second - rangeR.first) / (nR - 1);
456 const float zSpacing = (rangeZ.second - rangeZ.first) / (nZ - 1);
457 const float phiSpacing = (rangePhi.second - rangePhi.first) / (nPhi - 1);
458 const DataT phiPos = rangePhi.first + iPhi * phiSpacing;
459 // loop over grid and interpolate values
460 for (int iR = 0; iR < nR; ++iR) {
461 const DataT rPos = rangeR.first + iR * rSpacing;
462 for (int iZ = 0; iZ < nZ; ++iZ) {
463 const size_t idx = (iZ + nZ * (iR + iPhi * nR)); // unique index to Build index with other friend TTrees
464 const DataT zPos = rangeZ.first + iZ * zSpacing;
465 ir.emplace_back(iR);
466 iphi.emplace_back(iPhi);
467 iz.emplace_back(iZ);
468 r.emplace_back(rPos);
469 phi.emplace_back(phiPos);
470 z.emplace_back(zPos);
471 vals.emplace_back(data.interpolate(zPos, rPos, phiPos, mGrid3D)); // interpolated values
472 globalIdx.emplace_back(idx);
473 const float x = rPos * std::cos(phiPos);
474 const float y = rPos * std::sin(phiPos);
475 const LocalPosition3D pos(x, y, zPos);
476 unsigned char secNum = simOneSectorOnly ? 0 : std::floor(phiPos / SECPHIWIDTH);
477 Sector sector(secNum + (pos.Z() < 0) * SECTORSPERSIDE);
478 LocalPosition3D lPosTmp = Mapper::GlobalToLocal(pos, sector);
479 lPos.emplace_back(lPosTmp);
480 }
481 }
482 return std::make_tuple(vals, iz, ir, iphi, z, r, phi, lPos, globalIdx);
483 };
484
485 // define RDataFrame entry
486 auto dfStore = dFrame.DefineSlotEntry(treename, interpolate);
487
488 // define options of TFile
489 ROOT::RDF::RSnapshotOptions opt;
490 opt.fMode = option;
491
492 TStopwatch timer;
493 // note: first call has some overhead (~2s)
494 dfStore.Snapshot(treename, fileOut, {treename.data()}, opt);
495 timer.Print("u");
496}
497
498template <typename DataT>
499bool DataContainer3D<DataT>::getVertices(std::string_view treename, std::string_view fileIn, unsigned short& nR, unsigned short& nZ, unsigned short& nPhi)
500{
501 TFile fTmp(fileIn.data(), "READ");
502 TTree* tree = (TTree*)fTmp.Get(treename.data());
503 if (!tree) {
504 LOGP(warning, "Tree {} not found in input file {}", treename, fileIn);
505 return false;
506 }
507 tree->SetBranchAddress("nz", &nZ);
508 tree->SetBranchAddress("nr", &nR);
509 tree->SetBranchAddress("nphi", &nPhi);
510 tree->GetEntry(0);
511 delete tree;
512 return true;
513}
514
517
518// deprecated functions (to be removed...)
519template int o2::tpc::DataContainer3D<float>::writeToFile<float>(TFile&, const char*) const;
520template int o2::tpc::DataContainer3D<float>::writeToFile<double>(TFile&, const char*) const;
521template int o2::tpc::DataContainer3D<double>::writeToFile<float>(TFile&, const char*) const;
522template int o2::tpc::DataContainer3D<double>::writeToFile<double>(TFile&, const char*) const;
523template bool o2::tpc::DataContainer3D<float>::initFromFile<float>(TFile&, const char*);
524template bool o2::tpc::DataContainer3D<float>::initFromFile<double>(TFile&, const char*);
525template bool o2::tpc::DataContainer3D<double>::initFromFile<float>(TFile&, const char*);
526template bool o2::tpc::DataContainer3D<double>::initFromFile<double>(TFile&, const char*);
This class provides a simple method to store values on a large 3-Dim grid with ROOT io functionality.
int32_t i
uint16_t pos
Definition RawData.h:3
Definition of RegularGrid3D class.
Definition of TriCubic class.
static LocalPosition3D GlobalToLocal(const GlobalPosition3D &pos, const double alpha)
Definition Mapper.h:468
GLint GLenum GLint x
Definition glcorearb.h:403
GLint GLsizei count
Definition glcorearb.h:399
GLuint entry
Definition glcorearb.h:5735
GLuint index
Definition glcorearb.h:781
GLuint const GLchar * name
Definition glcorearb.h:781
GLsizei const GLfloat * value
Definition glcorearb.h:819
GLenum GLsizei GLsizei GLint * values
Definition glcorearb.h:1576
GLboolean * data
Definition glcorearb.h:298
GLuint GLfloat * val
Definition glcorearb.h:1582
GLboolean r
Definition glcorearb.h:1233
GLuint GLuint stream
Definition glcorearb.h:1806
GLubyte GLubyte GLubyte GLubyte w
Definition glcorearb.h:852
GLdouble GLdouble GLdouble z
Definition glcorearb.h:843
Global TPC definitions and constants.
Definition SimTraits.h:167
constexpr double SECPHIWIDTH
Definition Defs.h:45
constexpr unsigned char SECTORSPERSIDE
Definition Defs.h:40
void readData(o2::tpc::GBTFrameContainer &container, std::vector< std::ofstream * > &outfiles, int &run, int &done)
const auto & getData() const
static size_t getIndexR(size_t index, const int nz, const int nr, const int nphi)
static void dumpSlice(std::string_view treename, std::string_view fileIn, std::string_view fileOut, std::string_view option, std::pair< unsigned short, unsigned short > rangeiR, std::pair< unsigned short, unsigned short > rangeiZ, std::pair< unsigned short, unsigned short > rangeiPhi, const int nthreads=1)
void print() const
print the matrix
static void dumpInterpolation(std::string_view treename, std::string_view fileIn, std::string_view fileOut, std::string_view option, std::pair< float, float > rangeR, std::pair< float, float > rangeZ, std::pair< float, float > rangePhi, const int nR, const int nZ, const int nPhi, const int nthreads=1)
bool initFromFile(TFile &inpf, const char *name="data")
set values from file
int writeToFile(TFile &outf, const char *name="data") const
void setGrid(unsigned short nZ, unsigned short nR, unsigned short nPhi, const bool resize)
set the grid points
static void setAliasesForDump(TTree *tree)
static void setAliases(TTree *tree)
DataT interpolate(const DataT z, const DataT r, const DataT phi, const o2::tpc::RegularGrid3D< DataT > &grid) const
DataContainer3D< DataT > & operator*=(const DataT value)
operator overload
static size_t getIndexZ(size_t index, const int nz, const int nr, const int nphi)
static size_t getIndexPhi(size_t index, const int nz, const int nr, const int nphi)
static bool getVertices(std::string_view treename, std::string_view fileIn, unsigned short &nR, unsigned short &nZ, unsigned short &nPhi)
static DataContainer3D< DataT > * loadFromFile(TFile &inpf, const char *name="data")
get pointer to object from file (deprecated!)
DataContainer3D< DataT > & operator-=(const DataContainer3D< DataT > &other)
DataContainer3D< DataT > & operator+=(const DataContainer3D< DataT > &other)
static constexpr DataT getRMin()
static constexpr DataT getZMin()
static constexpr DataT getGridSpacingR(const unsigned int nR)
static constexpr DataT getPhiMin()
static constexpr DataT getGridSpacingZ(const unsigned int nZ)
static constexpr DataT getGridSpacingPhi(const unsigned int nPhi)
static bool normalizeGridToOneSector
the grid in phi direction is squashed from 2 Pi to (2 Pi / SECTORSPERSIDE). This can used to get the ...
VectorOfTObjectPtrs other
o2::InteractionRecord ir(0, 0)
std::unique_ptr< TTree > tree((TTree *) flIn.Get(std::string(o2::base::NameConf::CTFTREENAME).c_str()))
pedsdata resize(norbits)