Project
Loading...
Searching...
No Matches
Digit.cxx
Go to the documentation of this file.
1// Copyright 2020-2022 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
18/* ------ HISTORY ---------
19 10/03/2021 / complete review
20 05/11/2021 Add and review for the Cluster class
21*/
22
23#include <iostream>
24#include <TRandom.h>
26#include "HMPIDBase/Geo.h"
27#include "HMPIDBase/Param.h"
29
31
32namespace o2
33{
34namespace hmpid
35{
36
37// ============= Digit Class implementation =======
43Digit::Digit(int pad, uint16_t charge)
44{
45 mQ = charge > 0x0FFF ? 0x0FFF : charge;
46 mCh = a2C(pad);
47 mPh = a2P(pad);
48 mX = a2X(pad);
49 mY = a2Y(pad);
50}
51
60Digit::Digit(int chamber, int photo, int x, int y, uint16_t charge)
61{
62 mQ = charge > 0x0FFF ? 0x0FFF : charge;
63 mCh = chamber;
64 mPh = photo;
65 mX = x;
66 mY = y;
67}
68
77Digit::Digit(uint16_t charge, int equipment, int column, int dilogic, int channel)
78{
79 mQ = charge > 0x0FFF ? 0x0FFF : charge;
80 pad2Photo(equipment2Pad(equipment, column, dilogic, channel), &mCh, &mPh, &mX, &mY);
81}
82
90Digit::Digit(uint16_t charge, int module, int x, int y)
91{
92 mQ = charge > 0x0FFF ? 0x0FFF : charge;
93 pad2Photo(absolute2Pad(module, x, y), &mCh, &mPh, &mX, &mY);
94}
95
96// Digit ASCCI format Dump := [Chamber,PhotoCathod,X,Y]@(Orbit,BunchCrossing)=Charge
97std::ostream& operator<<(std::ostream& os, const o2::hmpid::Digit& d)
98{
99 os << "[" << (int)d.mCh << "," << (int)d.mPh << "," << (int)d.mX << "," << (int)d.mY << "]=" << d.mQ;
100 return os;
101};
102
103// ----- Coordinate Conversion -----
104
111uint32_t Digit::equipment2Pad(int Equi, int Colu, int Dilo, int Chan)
112{
113 // Check the input data
114 if (Equi < 0 || Equi >= Geo::MAXEQUIPMENTS || Colu < 0 || Colu >= Geo::N_COLUMNS ||
115 Dilo < 0 || Dilo >= Geo::N_DILOGICS || Chan < 0 || Chan >= Geo::N_CHANNELS) {
116 return -1;
117 }
118
119 int chan2y[6] = {3, 2, 4, 1, 5, 0}; // y coordinate translation for a channel address (index position) for even chamber
120
121 bool isEven = (Equi % Geo::EQUIPMENTSPERMODULE) == 0 ? true : false; // Calculate the odd/even of geometry
122 int ch = Equi / Geo::EQUIPMENTSPERMODULE; // The Module
123
124 // Calculate the x,y photo cathode relative coords For Odd equipment
125 int pc = (Colu / Geo::N_COLXSEGMENT) * 2 + 1; // col [0..23] -> [1,3,5]
126 int px = Geo::MAXXPHOTO - ((Dilo * Geo::DILOPADSROWS) + (Chan / Geo::DILOPADSCOLS));
127 int py = (Colu % Geo::DILOPADSROWS) * Geo::DILOPADSCOLS + chan2y[Chan % Geo::DILOPADSCOLS];
128 if (isEven) {
129 pc = 5 - pc;
130 py = Geo::MAXYPHOTO - py;
131 }
132 return abs(ch, pc, px, py); // Pack the coords into the PadID word
133}
134
141void Digit::pad2Equipment(uint32_t pad, int* Equi, int* Colu, int* Dilo, int* Chan)
142{
143 uint8_t ch, ph, px, py;
144 int y2chan[6] = {5, 3, 1, 0, 2, 4};
145
146 pad2Photo(pad, &ch, &ph, &px, &py); // Unpak the pad ID in the photo cathode coords
147
148 bool isEven = (ph % 2) == 0 ? true : false;
149 int eq = ch * Geo::EQUIPMENTSPERMODULE + 1;
150 px = Geo::MAXXPHOTO - px; // revert the X coord
151 if (isEven) {
152 eq--; // Correct the equipment number
153 py = Geo::MAXYPHOTO - py; // revert the Y coord
154 ph = 5 - ph; // revert the photo cathode index [0,2,4] -> [5,3,1]
155 }
156 *Dilo = px / Geo::DILOPADSROWS; // Calculate the Dilogic x [0..79] -> dil [0..9]
157 *Colu = ((ph / 2) * Geo::N_COLXSEGMENT) + (py / Geo::DILOPADSCOLS); // calculate the column (ph [1,3,5], y [0..47]) -> col [0..23]
158 *Chan = ((px % Geo::DILOPADSROWS) * Geo::DILOPADSCOLS) + y2chan[py % Geo::DILOPADSCOLS];
159 *Equi = eq;
160 return;
161}
162
171void Digit::absolute2Equipment(int Module, int x, int y, int* Equi, int* Colu, int* Dilo, int* Chan)
172{
173 uint32_t pad = absolute2Pad(Module, x, y);
174 pad2Equipment(pad, Equi, Colu, Dilo, Chan);
175 return;
176}
177
186void Digit::equipment2Absolute(int Equi, int Colu, int Dilo, int Chan, int* Module, int* x, int* y)
187{
188 uint32_t pad = equipment2Pad(Equi, Colu, Dilo, Chan);
189 pad2Absolute(pad, Module, x, y);
190 return;
191}
192
198uint32_t Digit::absolute2Pad(int Module, int x, int y)
199{
200 int ph = (y / Geo::N_PHOTOCATODSY) * 2 + ((x >= Geo::HALFXROWS) ? 1 : 0);
201 int px = x % Geo::HALFXROWS;
202 int py = y % Geo::N_PHOTOCATODSY;
203 return abs(Module, ph, px, py);
204}
205
211void Digit::pad2Absolute(uint32_t pad, int* Module, int* x, int* y)
212{
213 *Module = a2C(pad);
214 int ph = a2P(pad);
215 int px = a2X(pad);
216 int py = a2Y(pad);
217 *x = px + ((ph % 2 == 1) ? Geo::HALFXROWS : 0);
218 *y = ((ph >> 1) * Geo::N_PHOTOCATODSY) + py;
219 return;
220}
221
228void Digit::pad2Photo(uint32_t pad, uint8_t* chamber, uint8_t* photo, uint8_t* x, uint8_t* y)
229{
230 *chamber = a2C(pad);
231 *photo = a2P(pad);
232 *x = a2X(pad);
233 *y = a2Y(pad);
234 return;
235}
236
245void Digit::getPadAndTotalCharge(HitType const& hit, int& chamber, int& pc, int& px, int& py, float& totalcharge)
246{
247 double localX;
248 double localY;
249 chamber = hit.GetDetectorID();
250 double tmp[3] = {hit.GetX(), hit.GetY(), hit.GetZ()};
251 Param::instance()->mars2Lors(chamber, tmp, localX, localY);
252 Param::lors2Pad(localX, localY, pc, px, py);
253
254 totalcharge = Digit::qdcTot(hit.GetEnergyLoss(), hit.GetTime(), pc, px, py, localX, localY);
255 return;
256}
257
264{
265 double localX;
266 double localY;
267
268 const auto chamber = hit.GetDetectorID(); // chamber number is in detID
269 double tmp[3] = {hit.GetX(), hit.GetY(), hit.GetZ()};
270 // converting chamber id and hit coordiates to local coordinates
271 Param::instance()->mars2Lors(chamber, tmp, localX, localY);
272 // calculate charge fraction in given pad
273 return Digit::intMathieson(localX, localY, somepad);
274}
275
286Double_t Digit::qdcTot(Double_t e, Double_t time, Int_t pc, Int_t px, Int_t py, Double_t& localX, Double_t& localY)
287{
288 //
289 // Arguments: e-
290 // Returns:
291 double Q = 0;
292 if (time > 1.2e-6) {
293 Q = 0;
294 }
295 if (py < 0) {
296 return 0;
297 } else {
298 double y = Param::lorsY(pc, py);
299 localY = ((y - localY) > 0) ? y - 0.2 : y + 0.2; // shift to the nearest anod wire
300
301 double x = (localX > 66.6) ? localX - 66.6 : localX; // sagita is for PC (0-64) and not for chamber
302 double qdcEle = 34.06311 + 0.2337070 * x + 5.807476e-3 * x * x - 2.956471e-04 * x * x * x + 2.310001e-06 * x * x * x * x; // reparametrised from DiMauro
303
304 int iNele = int((e / 26e-9) * 0.8);
305 if (iNele < 1) {
306 iNele = 1; // number of electrons created by hit, if photon e=0 implies iNele=1
307 }
308 for (Int_t i = 1; i <= iNele; i++) {
309 double rnd = gRandom->Rndm();
310 if (rnd == 0) {
311 rnd = 1e-12; // 1e-12 is a protection against 0 from rndm
312 }
313 Q -= qdcEle * TMath::Log(rnd);
314 }
315 }
316 return Q;
317}
318
326double Digit::intPartMathiX(double x, int pad)
327{
328 double shift1 = -lorsX(pad) + o2::hmpid::Param::sizeHalfPadX();
329 double shift2 = -lorsX(pad) - o2::hmpid::Param::sizeHalfPadX();
330
333
334 return o2::hmpid::Param::k4x() * (atan(ux2) - atan(ux1));
335}
336
344double Digit::intPartMathiY(double y, int pad)
345{
346 double shift1 = -lorsY(pad) + o2::hmpid::Param::sizeHalfPadY();
347 double shift2 = -lorsY(pad) - o2::hmpid::Param::sizeHalfPadY();
348
351
352 return o2::hmpid::Param::k4y() * (atan(uy2) - atan(uy1));
353}
354
363double Digit::intMathieson(double localX, double localY, int pad)
364{
365 return 4. * intPartMathiX(localX, pad) * intPartMathiY(localY, pad);
366}
367
368// Mathieson function.
369// This is the answer to electrostatic problem of charge distrubution in MWPC described elsewhere. (NIM A370(1988)602-603)
370// Arguments: x- position of the center of Mathieson distribution
371// Returns: value of the Mathieson function
372double Digit::mathiesonX(double x)
373{
374 double lambda = x / o2::hmpid::Param::pitchAnodeCathode();
375 double tanh_v = tanh(o2::hmpid::Param::k2x() * lambda);
376 double a = 1 - tanh_v * tanh_v;
377 double b = 1 + o2::hmpid::Param::sqrtK3x() * o2::hmpid::Param::sqrtK3x() * tanh_v * tanh_v;
378 double mathi = o2::hmpid::Param::k1x() * a / b;
379 return mathi;
380}
381
382// Mathieson function.
383// This is the answer to electrostatic problem of charge distrubution in MWPC described elsewhere. (NIM A370(1988)602-603)
384// Arguments: x- position of the center of Mathieson distribution
385// Returns: value of the Mathieson function
386double Digit::mathiesonY(double y)
387{
388 double lambda = y / o2::hmpid::Param::pitchAnodeCathode();
389 double tanh_v = tanh(o2::hmpid::Param::k2y() * lambda);
390 double a = 1 - tanh_v * tanh_v;
391 double b = 1 + o2::hmpid::Param::sqrtK3y() * o2::hmpid::Param::sqrtK3y() * tanh_v * tanh_v;
392 double mathi = o2::hmpid::Param::k1y() * a / b;
393 return mathi;
394}
395
396/*
397// ---- Time conversion functions ----
398
405Double_t Digit::OrbitBcToTimeNs(uint32_t Orbit, uint16_t BC)
406{
407 return (BC * o2::constants::lhc::LHCBunchSpacingNS + Orbit * o2::constants::lhc::LHCOrbitNS);
408}
409
415uint32_t Digit::TimeNsToOrbit(Double_t TimeNs)
416{
417 return (uint32_t)(TimeNs / o2::constants::lhc::LHCOrbitNS);
418}
419
425uint16_t Digit::TimeNsToBc(Double_t TimeNs)
426{
427 return (uint16_t)(std::fmod(TimeNs, o2::constants::lhc::LHCOrbitNS) / o2::constants::lhc::LHCBunchSpacingNS);
428}
429
436void Digit::TimeNsToOrbitBc(double TimeNs, uint32_t& Orbit, uint16_t& Bc)
437{
438 Orbit = TimeNsToOrbit(TimeNs);
439 Bc = TimeNsToBc(TimeNs);
440 return;
441}
442
443
444// ---- Functions to manage Digit vectors ----
445
452bool Digit::eventEquipPadsComp(Digit& d1, Digit& d2)
453{
454 uint64_t t1, t2;
455 t1 = d1.getTriggerID();
456 t2 = d2.getTriggerID();
457 if (t1 < t2) {
458 return true;
459 }
460 if (t2 < t1) {
461 return false;
462 }
463 if (d1.getPadID() < d2.getPadID()) {
464 return true;
465 }
466 return false;
467};
468
469
475std::vector<o2::hmpid::Digit>* Digit::extractDigitsPerEvent(std::vector<o2::hmpid::Digit>& Digits, uint64_t EventID)
476{
477 std::vector<o2::hmpid::Digit>* subVector = new std::vector<o2::hmpid::Digit>();
478 for (const auto& digit : Digits) {
479 if (digit.getTriggerID() == EventID) {
480 subVector->push_back(digit);
481 }
482 }
483 return (subVector);
484};
485
490std::vector<uint64_t>* Digit::extractEvents(std::vector<o2::hmpid::Digit>& Digits)
491{
492 std::vector<uint64_t>* eventIds = new std::vector<uint64_t>();
493 for (const auto& digit : Digits) {
494 if (find(eventIds->begin(), eventIds->end(), digit.getTriggerID()) == eventIds->end()) {
495 eventIds->push_back(digit.getTriggerID());
496 }
497 }
498 return (eventIds);
499};
500*/
501
502} // namespace hmpid
503} // namespace o2
ClassImp(o2::hmpid::Digit)
int16_t charge
Definition RawEventData.h:5
int16_t time
Definition RawEventData.h:4
int32_t i
Header to collect LHC related constants.
V GetEnergyLoss() const
Definition BaseHits.h:103
E GetTime() const
Definition BaseHits.h:71
T GetZ() const
Definition BaseHits.h:66
T GetY() const
Definition BaseHits.h:65
T GetX() const
Definition BaseHits.h:64
short GetDetectorID() const
Definition BaseHits.h:73
HMPID Digit declaration.
Definition Digit.h:36
static uint32_t absolute2Pad(int Module, int x, int y)
Definition Digit.cxx:198
static int a2C(uint32_t pad)
Definition Digit.h:41
uint8_t mCh
Definition Digit.h:127
static void absolute2Equipment(int Module, int x, int y, int *Equi, int *Colu, int *Dilo, int *Chan)
Definition Digit.cxx:171
static int a2X(uint32_t pad)
Definition Digit.h:43
static double mathiesonY(double x)
Definition Digit.cxx:386
static int a2Y(uint32_t pad)
Definition Digit.h:44
uint8_t mX
Definition Digit.h:129
static void pad2Absolute(uint32_t pad, int *Module, int *x, int *y)
Definition Digit.cxx:211
static float lorsY(int pad)
Definition Digit.h:142
static int a2P(uint32_t pad)
Definition Digit.h:42
uint8_t mPh
Definition Digit.h:128
static double intPartMathiX(double x, int pad)
Definition Digit.cxx:326
static uint32_t equipment2Pad(int Equi, int Colu, int Dilo, int Chan)
Definition Digit.cxx:111
static void pad2Equipment(uint32_t pad, int *Equi, int *Colu, int *Dilo, int *Chan)
Definition Digit.cxx:141
uint16_t mQ
Definition Digit.h:126
static void pad2Photo(uint32_t pad, uint8_t *chamber, uint8_t *photo, uint8_t *x, uint8_t *y)
Definition Digit.cxx:228
static double qdcTot(double e, double time, int pc, int px, int py, double &localX, double &localY)
Definition Digit.cxx:286
static float lorsX(int pad)
Definition Digit.h:141
static double intPartMathiY(double y, int pad)
Definition Digit.cxx:344
static double mathiesonX(double x)
Definition Digit.cxx:372
static float getFractionalContributionForPad(o2::hmpid::HitType const &hit, int somepad)
Definition Digit.cxx:263
uint8_t mY
Definition Digit.h:130
static double intMathieson(double localX, double localY, int pad)
Definition Digit.cxx:363
static void getPadAndTotalCharge(o2::hmpid::HitType const &hit, int &chamber, int &pc, int &px, int &py, float &totalcharge)
Definition Digit.cxx:245
static uint32_t abs(int ch, int pc, int x, int y)
Definition Digit.h:39
static void equipment2Absolute(int Equi, int Colu, int Dilo, int Chan, int *Module, int *x, int *y)
Definition Digit.cxx:186
static constexpr int N_PHOTOCATODSY
Definition Geo.h:105
static constexpr int N_COLXSEGMENT
Definition Geo.h:81
static constexpr int N_COLUMNS
Definition Geo.h:82
static constexpr int N_DILOGICS
Definition Geo.h:83
static constexpr int N_CHANNELS
Definition Geo.h:84
static constexpr int DILOPADSROWS
Definition Geo.h:96
static constexpr int DILOPADSCOLS
Definition Geo.h:95
static constexpr int MAXYPHOTO
Definition Geo.h:107
static constexpr int MAXXPHOTO
Definition Geo.h:106
static constexpr int EQUIPMENTSPERMODULE
Definition Geo.h:98
static constexpr int HALFXROWS
Definition Geo.h:93
static constexpr int MAXEQUIPMENTS
Definition Geo.h:79
static float sizeHalfPadX()
Definition Param.h:72
void mars2Lors(Int_t c, double *m, double &x, double &y) const
Definition Param.h:199
static double k4y()
Definition Param.h:254
static double pitchAnodeCathode()
Definition Param.h:246
static float sizeHalfPadY()
Definition Param.h:73
static void lors2Pad(float x, float y, Int_t &pc, Int_t &px, Int_t &py)
Definition Param.cxx:488
static float lorsY(Int_t pc, Int_t pady)
Definition Param.h:88
static double k1y()
Definition Param.h:253
static double k2y()
Definition Param.h:252
static Param * instance()
Definition Param.cxx:438
static double k2x()
Definition Param.h:248
static double k1x()
Definition Param.h:249
static double sqrtK3y()
Definition Param.h:251
static double k4x()
Definition Param.h:250
static double sqrtK3x()
Definition Param.h:247
GLint GLenum GLint x
Definition glcorearb.h:403
GLboolean GLboolean GLboolean b
Definition glcorearb.h:1233
GLint y
Definition glcorearb.h:270
GLboolean GLboolean GLboolean GLboolean a
Definition glcorearb.h:1233
a couple of static helper functions to create timestamp values for CCDB queries or override obsolete ...
std::ostream & operator<<(std::ostream &stream, o2::InteractionRecord const &ir)