26#include <fairlogger/Logger.h>  
   38  LOG(info) << 
"Initializing digitizer";
 
   40  mChips.resize(mNumberOfChips); 
 
   41  for (
int i = mNumberOfChips; 
i--;) {
 
   42    mChips[
i].setChipIndex(
i);
 
   44      mChips[
i].setNoiseMap(mNoiseMap);
 
   48      mChips[
i].setDeadChanMap(mDeadChanMap);
 
   53  auto file = TFile::Open(mResponseFile.data());
 
   55    LOG(fatal) << 
"Cannot open response file " << mResponseFile;
 
   60  mChipSimRespVD = mChipSimResp;   
 
   61  mChipSimRespMLOT = mChipSimResp; 
 
   65  LOG(info) << 
" Depth max VD: " << mChipSimRespVD->getDepthMax();
 
   66  LOG(info) << 
" Depth min VD: " << mChipSimRespVD->getDepthMin();
 
   68  LOG(info) << 
" Depth max MLOT: " << mChipSimRespMLOT->getDepthMax();
 
   69  LOG(info) << 
" Depth min MLOT: " << mChipSimRespMLOT->getDepthMin();
 
   71  float thicknessVD = 0.0095; 
 
   72  float thicknessMLOT = 0.1;  
 
   76  mSimRespVDShift = -mChipSimRespVD->getDepthMax(); 
 
   79  mSimRespMLOTShift = mChipSimRespMLOT->getDepthMax() - thicknessMLOT / 2.f; 
 
   80  mSimRespOrientation = 
false;
 
   85  LOGP(info, 
"TRK Digitizer is initialised.");
 
   87  LOGP(info, 
"VD shift = {}  ; ML/OT shift = {} = {} - {}", mSimRespVDShift, mSimRespMLOTShift, mChipSimRespMLOT->getDepthMax(), thicknessMLOT / 2.f);
 
   88  LOGP(info, 
"VD pixel scale on x = {} ; z = {}", mSimRespVDScaleX, mSimRespVDScaleZ);
 
   89  LOGP(info, 
"ML/OT pixel scale on x = {} ; z = {}", mSimRespMLOTScaleX, mSimRespMLOTScaleZ);
 
   90  LOGP(info, 
"Response orientation: {}", mSimRespOrientation ? 
"flipped" : 
"normal");
 
   97  if (mGeometry->getSubDetID(chipID) == 0) { 
 
   98    return mChipSimRespVD;
 
  101  else if (mGeometry->getSubDetID(chipID) == 1) { 
 
  102    return mChipSimRespMLOT;
 
  112  LOG(info) << 
" Digitizing " << mGeometry->
getName() << 
" (ID: " << mGeometry->
getDetID()
 
  113            << 
") hits of entry " << evID << 
" from source " << srcID
 
  114            << 
" at time " << mEventTime << 
" ROFrame= " << mNewROFrame << 
")" 
  116            << 
" Min/Max ROFrames " << mROFrameMin << 
"/" << mROFrameMax;
 
  122  if (mNewROFrame > mROFrameMin) {
 
  126  int nHits = hits->size();
 
  127  std::vector<int> hitIdx(nHits);
 
  128  std::iota(std::begin(hitIdx), std::end(hitIdx), 0);
 
  130  std::sort(hitIdx.begin(), hitIdx.end(),
 
  131            [hits](
auto lhs, 
auto rhs) {
 
  132              return (*hits)[lhs].GetDetectorID() < (*hits)[rhs].GetDetectorID();
 
  134  LOG(info) << 
"Processing " << nHits << 
" hits";
 
  135  for (
int i : hitIdx) {
 
  136    processHit((*hits)[
i], mROFrameMax, evID, srcID);
 
  150  LOG(info) << 
"Setting event time ";
 
  162    if (mCollisionTimeWrtROF < 0 && nbc > 0) {
 
  168    LOG(info) << 
" NewROFrame " << mNewROFrame << 
" = " << nbc << 
"/" << mParams.
getROFrameLengthInBC() << 
" (nbc/mParams.getROFrameLengthInBC()";
 
  176  if (mNewROFrame < mROFrameMin) {
 
  177    LOG(error) << 
"New ROFrame " << mNewROFrame << 
" (" << irt << 
") precedes currently cashed " << mROFrameMin;
 
  178    throw std::runtime_error(
"deduced ROFrame precedes already processed one");
 
  181  if (mParams.
isContinuous() && mROFrameMax < mNewROFrame) {
 
  182    mROFrameMax = mNewROFrame - 1; 
 
  190  if (frameLast > mROFrameMax) {
 
  191    frameLast = mROFrameMax;
 
  194  getExtraDigBuffer(mROFrameMax);
 
  195  LOG(info) << 
"Filling " << mGeometry->
getName() << 
" digits output for RO frames " << mROFrameMin << 
":" 
  201  for (; mROFrameMin <= frameLast; mROFrameMin++) {
 
  205    auto& extra = *(mExtraBuff.front().get());
 
  206    for (
auto& chip : mChips) {
 
  207      if (chip.isDisabled()) {
 
  211      auto& 
buffer = chip.getPreDigits();
 
  215      auto itBeg = 
buffer.begin();
 
  217      ULong64_t maxKey = chip.getOrderingKey(mROFrameMin + 1, 0, 0) - 1; 
 
  218      for (; iter != 
buffer.end(); ++iter) {
 
  219        if (iter->first > maxKey) {
 
  222        auto& preDig = iter->second; 
 
  224          int digID = mDigits->size();
 
  225          mDigits->emplace_back(chip.getChipIndex(), preDig.row, preDig.col, preDig.charge);
 
  226          LOG(
debug) << 
"Adding digit ID: " << digID << 
" with chipID: " << chip.getChipIndex() << 
", row: " << preDig.row << 
", col: " << preDig.col << 
", charge: " << preDig.charge;
 
  227          mMCLabels->
addElement(digID, preDig.labelRef.label);
 
  228          auto& nextRef = preDig.labelRef; 
 
  229          while (nextRef.next >= 0) {
 
  230            nextRef = extra[nextRef.next];
 
  235      buffer.erase(itBeg, iter);
 
  245      mROFRecords->push_back(rcROF);
 
  249    mExtraBuff.emplace_back(mExtraBuff.front().release());
 
  250    mExtraBuff.pop_front();
 
  255void Digitizer::processHit(
const o2::itsmft::Hit& hit, uint32_t& maxFr, 
int evID, 
int srcID)
 
  258  int subDetID = mGeometry->getSubDetID(chipID);
 
  260  int layer = mGeometry->getLayer(chipID);
 
  261  int disk = mGeometry->getDisk(chipID);
 
  264    LOG(
debug) << 
"Skipping disk " << disk;
 
  268  LOG(
debug) << 
"Processing hit for chip " << chipID;
 
  269  auto& chip = mChips[chipID];
 
  270  if (chip.isDisabled()) {
 
  271    LOG(
debug) << 
"Skipping disabled chip " << chipID;
 
  274  float timeInROF = hit.
GetTime() * sec2ns;
 
  275  LOG(
debug) << 
"timeInROF: " << timeInROF;
 
  276  if (timeInROF > 20e3) {
 
  277    const int maxWarn = 10;
 
  278    static int warnNo = 0;
 
  279    if (warnNo < maxWarn) {
 
  280      LOG(warning) << 
"Ignoring hit with time_in_event = " << timeInROF << 
" ns" 
  281                   << ((++warnNo < maxWarn) ? 
"" : 
" (suppressing further warnings)");
 
  286    timeInROF += mCollisionTimeWrtROF;
 
  290    LOG(
debug) << 
"Ignoring hit with timeInROF = " << timeInROF;
 
  303  int nFrames = roFrameRelMax + 1 - roFrameRel;
 
  304  uint32_t roFrameMax = mNewROFrame + roFrameRelMax;
 
  305  if (roFrameMax > maxFr) {
 
  324    LOG(
debug) << 
"Called curved to flat: " << xyzLocS.x() << 
" -> " << xyFlatS.x() << 
", " << xyzLocS.y() << 
" -> " << xyFlatS.y();
 
  326    xyzLocS.SetXYZ(xyFlatS.x(), xyFlatS.y(), xyzLocS.Z());
 
  327    xyzLocE.SetXYZ(xyFlatE.x(), xyFlatE.y(), xyzLocE.Z());
 
  337  LOG(
debug) << 
"local original: startPos = " << xyzLocS << 
", endPos = " << xyzLocE << std::endl;
 
  338  xyzLocS.SetY(xyzLocS.Y());
 
  339  xyzLocE.SetY(xyzLocE.Y());
 
  341  LOG(
debug) << 
"rescaled Y: startPos = " << xyzLocS << 
", endPos = " << xyzLocE << std::endl;
 
  352  LOG(
debug) << 
"Step into the sensitive volume: " << 
step << 
".  Number of steps: " << nSteps;
 
  353  int rowS = -1, colS = -1, rowE = -1, colE = -1, nSkip = 0;
 
  358    if (++nSkip >= nSteps) {
 
  359      LOG(
debug) << 
"Did not enter to sensitive matrix, " << nSkip << 
" >= " << nSteps;
 
  367    if (++nSkip >= nSteps) {
 
  368      LOG(
debug) << 
"Did not enter to sensitive matrix, " << nSkip << 
" >= " << nSteps;
 
  379    std::swap(rowS, rowE);
 
  382    std::swap(colS, colE);
 
  400  int rowSpan = rowE - rowS + 1, colSpan = colE - colS + 1; 
 
  402  float respMatrix[rowSpan][colSpan]; 
 
  403  std::fill(&respMatrix[0][0], &respMatrix[0][0] + rowSpan * colSpan, 0.f);
 
  406  nElectrons *= nStepsInv;                                                  
 
  411  int rowPrev = -1, colPrev = -1, 
row, 
col;
 
  412  float cRowPix = 0.f, cColPix = 0.f; 
 
  421  xyzLocS.SetY(xyzLocS.Y() + ((subDetID == 0) ? mSimRespVDShift : mSimRespMLOTShift));
 
  424  for (
int iStep = nSteps; iStep--;) {
 
  427    if (
row != rowPrev || 
col != colPrev) { 
 
  434    bool flipCol = 
false, flipRow = 
false;
 
  436    float rowMax{}, colMax{};
 
  439      rowMax = 0.5f * Segmentation::PitchRowVD * mSimRespVDScaleX;
 
  440      colMax = 0.5f * Segmentation::PitchColVD * mSimRespVDScaleZ;
 
  441      rspmat = resp->
getResponse(mSimRespVDScaleX * (xyzLocS.X() - cRowPix), mSimRespVDScaleZ * (xyzLocS.Z() - cColPix), xyzLocS.Y(), flipRow, flipCol, rowMax, colMax);
 
  443      rowMax = 0.5f * Segmentation::PitchRowMLOT * mSimRespMLOTScaleX;
 
  444      colMax = 0.5f * Segmentation::PitchColMLOT * mSimRespMLOTScaleZ;
 
  445      rspmat = resp->
getResponse(mSimRespMLOTScaleX * (xyzLocS.X() - cRowPix), mSimRespMLOTScaleZ * (xyzLocS.Z() - cColPix), xyzLocS.Y(), flipRow, flipCol, rowMax, colMax);
 
  450    if (rspmat == 
nullptr) {
 
  451      LOG(
debug) << 
"Error in rspmat for step " << iStep << 
" / " << nSteps;
 
  454    LOG(
debug) << 
"rspmat valid! for step " << iStep << 
" / " << nSteps << 
", (row,col) = (" << 
row << 
"," << 
col << 
")";
 
  459      if (rowDest < 0 || rowDest >= rowSpan) {
 
  464        if (colDest < 0 || colDest >= colSpan) {
 
  467        respMatrix[rowDest][colDest] += rspmat->getValue(irow, icol, mSimRespOrientation ? !flipRow : flipRow, !flipCol);
 
  474  auto roFrameAbs = mNewROFrame + roFrameRel;
 
  475  LOG(
debug) << 
"Spanning through rows and columns; rowspan = " << rowSpan << 
" colspan = " << colSpan << 
" = " << colE << 
" - " << colS << 
" +1 " << std::endl;
 
  476  for (
int irow = rowSpan; irow--;) {          
 
  477    uint16_t rowIS = irow + rowS;              
 
  478    for (
int icol = colSpan; icol--;) {        
 
  479      float nEleResp = respMatrix[irow][icol]; 
 
  480      if (nEleResp <= 1.e-36) {
 
  483      LOG(
debug) << 
"nEleResp: value " << nEleResp << 
" for pixel " << irow << 
" " << icol << std::endl;
 
  484      int nEle = gRandom->Poisson(nElectrons * nEleResp); 
 
  485      LOG(
debug) << 
"Charge detected in the pixel: " << nEle << 
" for pixel " << irow << 
" " << icol << std::endl;
 
  488        LOG(
debug) << 
"Ignoring pixel with nEle = " << nEle << 
" < min charge to account " 
  493      uint16_t colIS = icol + colS; 
 
  494      if (mNoiseMap && mNoiseMap->
isNoisy(chipID, rowIS, colIS)) {
 
  497      if (mDeadChanMap && mDeadChanMap->
isNoisy(chipID, rowIS, colIS)) {
 
  501      registerDigits(chip, roFrameAbs, timeInROF, nFrames, rowIS, colIS, nEle, lbl);
 
  513  LOG(
debug) << 
"Registering digits for chip " << chip.
getChipIndex() << 
" at ROFrame " << roFrame
 
  514             << 
" row " << 
row << 
" col " << 
col << 
" nEle " << nEle << 
" label " << lbl;
 
  516  for (
int i = 0; 
i < nROF; 
i++) {
 
  517    uint32_t roFr = roFrame + 
i;
 
  525    if (roFr > mEventROFrameMax) {
 
  526      mEventROFrameMax = roFr;
 
  528    if (roFr < mEventROFrameMin) {
 
  529      mEventROFrameMin = roFr;
 
  535      LOG(
debug) << 
"Added digit " << 
key << 
"  " << roFr << 
"  " << 
row << 
"  " << 
col << 
"  " << nEleROF;
 
  541      ExtraDig* extra = getExtraDigBuffer(roFr);
 
  545        if ((*extra)[nxt].
label == lbl) { 
 
  549        nxt = (*extra)[nxt].next;
 
  556      extra->emplace_back(lbl);
 
Definition of the SegmentationChipclass.
Definition of the TRK digitizer.
math_utils::Point3D< T > GetPos() const
unsigned short GetDetectorID() const
static const DPLDigitizerParam< N > & Instance()
const char * getName() const
const o2::detectors::DetID & getDetID() const
const Mat3D & getMatrixL2G(int sensID) const
static int constexpr NPix
float getMaxDuration() const
float getCollectedCharge(float totalNEle, float tMin, float tMax) const
bool getResponse(float vRow, float vCol, float cDepth, AlpideRespSimMat &dest) const
void addDigit(ULong64_t key, UInt_t roframe, UShort_t row, UShort_t col, int charge, o2::MCCompLabel lbl)
o2::itsmft::PreDigit * findDigit(ULong64_t key)
static ULong64_t getOrderingKey(UInt_t roframe, UShort_t row, UShort_t col)
Get global ordering key made of readout frame, column and row.
UShort_t getChipIndex() const
int getMinChargeToAccount() const
float getROFrameLengthInv() const
virtual void print() const
const SignalShape & getSignalShape() const
float getEnergyToNElectrons() const
int getROFrameLengthInBC() const
float getStrobeDelay() const
bool isContinuous() const
int getChargeThreshold() const
float getNSimStepsInv() const
float getROFrameLength() const
float getStrobeLength() const
void setEventTime(const o2::InteractionTimeRecord &irt)
auto getChipResponse(int chipID)
bool isContinuous() const
void process(const std::vector< Hit > *hits, int evID, int srcID)
Steer conversion of hits to digits.
void fillOutputContainer(uint32_t maxFrame=0xffffffff)
Int_t getNumberOfChips() const
math_utils::Point3D< Float_t > GetPosStart() const
bool isFullChipMasked(int chip) const
bool isNoisy(int chip, int row, int col) const
const BCData & getBCData() const
void setFirstEntry(int idx)
int getFirstEntry() const
void setROFrame(ROFtype rof)
static bool localToDetector(float x, float z, int &iRow, int &iCol)
static bool detectorToLocal(L row, L col, T &xRow, T &zCol)
static constexpr float PitchColVD
static constexpr float PitchColMLOT
static constexpr float PitchRowMLOT
static constexpr float PitchRowVD
GLuint GLsizei const GLchar * label
GLenum GLuint GLint GLint layer
constexpr double LHCBunchSpacingNS
constexpr std::array< int, nLayers > nRows
int64_t differenceInBC(const InteractionRecord &other) const
void setFromLong(int64_t l)
double timeInBCNS
time in NANOSECONDS relative to orbit/bc
int next
eventual next contribution to the same pixel
o2::MCCompLabel label
hit label
PreDigitLabelRef labelRef
label and reference to the next one
IR getFirstSampledTFIR() const
get TF and HB (abs) for this IR
LOG(info)<< "Compressed in "<< sw.CpuTime()<< " s"