HLT/PHOS/AliHLTPHOSDigitMaker.cxx

   1 /**************************************************************************
   2  * This file is property of and copyright by the ALICE HLT Project        *
   3  * All rights reserved.                                                   *
   4  *                                                                        *
   5  * Primary Authors: Oystein Djuvsland                                     *
   6  s*                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15  /**
  16  * @file   AliHLTPHOSClusterizer.cxx
  17  * @author Oystein Djuvsland
  18  * @date
  19  * @brief  Digit maker for PHOS HLT
  20  */
  21
  22
  23
  24
  25 // see header file for class documentation
  26 // or
  27 // refer to README to build package
  28 // or
  29 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
  30
  31 #include "AliHLTPHOSDigitMaker.h"
  32 #include "AliHLTPHOSDigit.h"
  33 #include "AliHLTPHOSConstants.h"
  34 #include "AliHLTPHOSBaseline.h"
  35 #include "AliHLTPHOSMapper.h"
  36 #include "TTree.h"
  37 #include "TBranch.h"
  38 #include "TClonesArray.h"
  39 #include "TFile.h"
  40 #include "TH2F.h"
  41
  42 #include "AliHLTPHOSValidCellDataStruct.h"
  43 #include "AliHLTPHOSRcuCellEnergyDataStruct.h"
  44 #include "AliHLTPHOSChannelDataStruct.h"
  45 #include "AliHLTPHOSChannelDataHeaderStruct.h"
  46 #include "AliHLTPHOSDigitDataStruct.h"
  47 #include "AliHLTPHOSDigitContainerDataStruct.h"
  48 #include "AliHLTPHOSSharedMemoryInterfacev2.h" // added by PTH
  49
  50
  51 ClassImp(AliHLTPHOSDigitMaker);
  52
  53 using namespace PhosHLTConst;
  54
  55 AliHLTPHOSDigitMaker::AliHLTPHOSDigitMaker() :
  56   AliHLTPHOSBase(),
  57   fCellDataPtr(0),
  58   fDigitContainerStructPtr(0),
  59   fDigitArrayPtr(0),
  60   fDigitStructPtr(0),
  61   fDigitCount(0),
  62   fShmPtr(0),
  63   fOrdered(true)
  64 {
  65   // See header file for documentation
  66
  67   fShmPtr = new AliHLTPHOSSharedMemoryInterfacev2();
  68
  69   for(int x = 0; x < N_XCOLUMNS_MOD; x++)
  70     {
  71       for(int z = 0; z < N_ZROWS_MOD; z++)
  72         {
  73           fHighGainFactors[x][z] = 0.005;
  74           fLowGainFactors[x][z] = 0.08;
  75           fDigitThresholds[x][z][HIGH_GAIN] = 2;
  76           fDigitThresholds[x][z][LOW_GAIN] = 2;
  77           fBadChannelMask[x][z][HIGH_GAIN] = 1;
  78           fBadChannelMask[x][z][LOW_GAIN] = 1;
  79         }
  80     }
  81
  82     for(int x = 0; x < N_XCOLUMNS_RCU; x++)
  83       {
  84         for(int z = 0; z < N_ZROWS_RCU; z++)
  85           {
  86             for(int gain = 0; gain < N_GAINS; gain++)
  87               {
  88                 fEnergyArray[x][z][gain] = 0;
  89               }
  90           }
  91       }
  92
  93
  94 }
  95
  96 AliHLTPHOSDigitMaker::~AliHLTPHOSDigitMaker()
  97 {
  98   //See header file for documentation
  99 }
 100
 101 Int_t
 102 AliHLTPHOSDigitMaker::MakeDigits(AliHLTPHOSChannelDataHeaderStruct* channelDataHeader)
 103 {
 104   //See header file for documentation
 105
 106   Int_t j = 0;
 107   Float_t tmpEnergy = 0;
 108
 109 //   Int_t xMod = -1;
 110 //   Int_t zMod = -1;
 111
 112   AliHLTPHOSMapper mapper;
 113   UShort_t coord1[4];
 114   UShort_t coord2[4];
 115
 116
 117   AliHLTPHOSChannelDataStruct* currentchannel = 0;
 118   AliHLTPHOSChannelDataStruct* currentchannelLG = 0;
 119   AliHLTPHOSChannelDataStruct* tmpchannel = 0;
 120
 121   Reset();
 122
 123   fShmPtr->SetMemory(channelDataHeader);
 124   currentchannel = fShmPtr->NextChannel();
 125
 126 //   while(currentchannel != 0)
 127 //     {
 128 //       fEnergyArray[currentchannel->fX][currentchannel->fZ][currentchannel->fGain] = currentchannel->fEnergy;
 129 //       currentchannel = fShmPtr->NextChannel();
 130 //     }
 131
 132
 133   while(currentchannel != 0)
 134     {
 135       AliHLTPHOSMapper::GetChannelCoord(currentchannel->fChannelID, coord1);
 136       if(fOrdered)
 137         {
 138           tmpchannel = currentchannel;
 139           if(coord1[2] == HIGH_GAIN)
 140             {
 141
 142               if(currentchannel->fEnergy < MAX_BIN_VALUE)
 143                 {
 144                   fDigitStructPtr = &(fDigitContainerStructPtr->fDigitDataStruct[j+fDigitCount]);
 145                   fDigitStructPtr->fX = coord1[0];
 146                   fDigitStructPtr->fZ = coord1[1];
 147                   fDigitStructPtr->fAmplitude = currentchannel->fEnergy;
 148                   fDigitStructPtr->fEnergy = currentchannel->fEnergy * fHighGainFactors[coord1[0]][coord1[1]];
 149                   //TODO: fix time //CRAP
 150                   fDigitStructPtr->fTime = currentchannel->fTime * 0.0000001;
 151                   fDigitStructPtr->fCrazyness = currentchannel->fCrazyness;
 152                   fDigitStructPtr->fModule = coord1[3];
 153                   j++;
 154                   currentchannel = fShmPtr->NextChannel();
 155                   AliHLTPHOSMapper::GetChannelCoord(currentchannel->fChannelID, coord2);
 156                   if(coord1[0] == coord2[0] && coord1[1] == coord2[1])
 157                     {
 158                       fShmPtr->NextChannel();
 159                     }
 160                 }
 161               else
 162                 {
 163                   currentchannel = fShmPtr->NextChannel();
 164                   AliHLTPHOSMapper::GetChannelCoord(currentchannel->fChannelID, coord2);
 165                   if(coord2[2] == LOW_GAIN)
 166                     {
 167                       fDigitStructPtr = &(fDigitContainerStructPtr->fDigitDataStruct[j+fDigitCount]);
 168                       fDigitStructPtr->fX = coord2[0];
 169                       fDigitStructPtr->fZ = coord2[1];
 170                       fDigitStructPtr->fAmplitude = currentchannel->fEnergy;
 171                       fDigitStructPtr->fEnergy = currentchannel->fEnergy * fHighGainFactors[coord2[0]][coord2[1]];
 172                       //TODO: fix time //CRAP
 173                       fDigitStructPtr->fTime = currentchannel->fTime * 0.0000001;
 174                       fDigitStructPtr->fCrazyness = currentchannel->fCrazyness;
 175                       fDigitStructPtr->fModule = coord2[3];
 176                       j++;
 177                       currentchannel = fShmPtr->NextChannel();
 178                     }
 179                   else
 180                     {
 181                       fDigitStructPtr = &(fDigitContainerStructPtr->fDigitDataStruct[j+fDigitCount]);
 182                       fDigitStructPtr->fX = coord1[0];
 183                       fDigitStructPtr->fZ = coord1[1];
 184                       fDigitStructPtr->fAmplitude = tmpchannel->fEnergy;
 185                       fDigitStructPtr->fEnergy = tmpchannel->fEnergy * fHighGainFactors[coord1[0]][coord1[1]];
 186                       //TODO: fix time //CRAP
 187                       fDigitStructPtr->fTime = tmpchannel->fTime * 0.0000001;
 188                       fDigitStructPtr->fCrazyness = 100;
 189                       fDigitStructPtr->fModule = coord1[3];
 190                       j++;
 191                     }
 192
 193                 }
 194             }
 195           else
 196             {
 197               fDigitStructPtr = &(fDigitContainerStructPtr->fDigitDataStruct[j+fDigitCount]);
 198               fDigitStructPtr->fX = coord1[0];
 199               fDigitStructPtr->fZ = coord1[1];
 200               fDigitStructPtr->fAmplitude = tmpchannel->fEnergy;
 201               fDigitStructPtr->fEnergy = tmpchannel->fEnergy * fHighGainFactors[coord1[0]][coord1[1]];
 202               //TODO: fix time //CRAP
 203               fDigitStructPtr->fTime = tmpchannel->fTime * 0.0000001;
 204               fDigitStructPtr->fCrazyness = 100;
 205               fDigitStructPtr->fModule = coord1[3];
 206               j++;
 207               currentchannel = fShmPtr->NextChannel();
 208             }
 209
 210         }
 211       else
 212         {
 213           if(coord1[2] == LOW_GAIN)
 214             {
 215               currentchannelLG = currentchannel;
 216               currentchannel = fShmPtr->NextChannel();
 217               AliHLTPHOSMapper::GetChannelCoord(currentchannel->fChannelID, coord2);
 218
 219               if(coord1[0] == coord2[0] && coord1[1] == coord2[1])
 220                 {
 221                   if(currentchannel->fEnergy < MAX_BIN_VALUE)
 222                     {
 223                       fDigitStructPtr = &(fDigitContainerStructPtr->fDigitDataStruct[j+fDigitCount]);
 224                       fDigitStructPtr->fX = coord2[0];
 225                       fDigitStructPtr->fZ = coord2[1];
 226                       fDigitStructPtr->fAmplitude = currentchannel->fEnergy;
 227                       fDigitStructPtr->fEnergy = currentchannel->fEnergy * fHighGainFactors[coord2[0]][coord2[1]];
 228                       //TODO: fix time //CRAP
 229                       fDigitStructPtr->fTime = currentchannel->fTime * 0.0000001;
 230                       fDigitStructPtr->fCrazyness = currentchannel->fCrazyness;
 231                       fDigitStructPtr->fModule = coord2[3];
 232                       j++;
 233                       currentchannel = fShmPtr->NextChannel();
 234                     }
 235                   else
 236                     {
 237                       fDigitStructPtr = &(fDigitContainerStructPtr->fDigitDataStruct[j+fDigitCount]);
 238                       fDigitStructPtr->fX = coord1[0];
 239                       fDigitStructPtr->fZ = coord1[1];
 240                       fDigitStructPtr->fAmplitude = currentchannelLG->fEnergy;
 241                       fDigitStructPtr->fEnergy = currentchannelLG->fEnergy * fHighGainFactors[coord1[0]][coord1[1]];
 242                       //TODO: fix time //CRAP
 243                       fDigitStructPtr->fTime = currentchannelLG->fTime * 0.0000001;
 244                       fDigitStructPtr->fCrazyness = currentchannel->fCrazyness;
 245                       fDigitStructPtr->fModule = coord1[3];
 246                       j++;
 247                       currentchannel = fShmPtr->NextChannel();
 248                     }
 249                 }
 250               else
 251                 {
 252                   fDigitStructPtr = &(fDigitContainerStructPtr->fDigitDataStruct[j+fDigitCount]);
 253                   fDigitStructPtr->fX = coord1[0];
 254                   fDigitStructPtr->fZ = coord1[1];
 255                   fDigitStructPtr->fAmplitude = currentchannelLG->fEnergy;
 256                   fDigitStructPtr->fEnergy = currentchannelLG->fEnergy * fHighGainFactors[coord1[0]][coord1[1]];
 257                   //TODO: fix time //CRAP
 258                   fDigitStructPtr->fTime = currentchannelLG->fTime * 0.0000001;
 259                   fDigitStructPtr->fCrazyness = currentchannel->fCrazyness;
 260                   fDigitStructPtr->fModule = coord1[3];
 261                   j++;
 262                 }
 263             }
 264           else
 265             {
 266               fDigitStructPtr = &(fDigitContainerStructPtr->fDigitDataStruct[j+fDigitCount]);
 267               fDigitStructPtr->fX = coord1[0];
 268               fDigitStructPtr->fZ = coord1[1];
 269               fDigitStructPtr->fAmplitude = currentchannel->fEnergy;
 270               fDigitStructPtr->fEnergy = currentchannel->fEnergy * fHighGainFactors[coord1[0]][coord1[1]];
 271               //TODO: fix time //CRAP
 272               fDigitStructPtr->fTime = currentchannel->fTime * 0.0000001;
 273               fDigitStructPtr->fCrazyness = currentchannel->fCrazyness;
 274               fDigitStructPtr->fModule = coord2[3];
 275               j++;
 276               fShmPtr->NextChannel();
 277             }
 278         }
 279     }
 280
 281
 282
 283
 284 //   for(int x = 0; x < N_XCOLUMNS_RCU; x++)
 285 //     {
 286 //       for(int z = 0; z < N_ZROWS_RCU; z++)
 287 //      {
 288 //        xMod = x + rcuData->fRcuX * N_XCOLUMNS_RCU;
 289 //        zMod = z + rcuData->fRcuZ * N_ZROWS_RCU;
 290
 291 //        if(fEnergyArray[x][z][HIGH_GAIN] > fDigitThresholds[xMod][zMod][HIGH_GAIN] && fEnergyArray[x][z][HIGH_GAIN] < MAX_BIN_VALUE && fBadChannelMask[xMod][zMod][HIGH_GAIN])
 292 //          {
 293 //            fDigitStructPtr = &(fDigitContainerStructPtr->fDigitDataStruct[j+fDigitCount]);
 294 //            fDigitStructPtr->fX = xMod;
 295 //            fDigitStructPtr->fZ = zMod;
 296 //            fDigitStructPtr->fAmplitude = fEnergyArray[x][z][HIGH_GAIN];
 297 //            fDigitStructPtr->fEnergy = fEnergyArray[x][z][HIGH_GAIN] * fHighGainFactors[xMod][zMod];
 298 //            //TODO: fix time //CRAP
 299 //            fDigitStructPtr->fTime = fCellDataPtr->fTime * 0.0000001;
 300 //            fDigitStructPtr->fCrazyness = fCellDataPtr->fCrazyness;
 301 //            fDigitStructPtr->fModule = rcuData->fModuleID;
 302 //            j++;
 303 //          }
 304 //        else if(fEnergyArray[x][z][LOW_GAIN] >= MAX_BIN_VALUE && fBadChannelMask[xMod][zMod][LOW_GAIN])
 305 //          {
 306 //            if(fEnergyArray[x][z][LOW_GAIN] > fDigitThresholds[xMod][zMod][LOW_GAIN])
 307 //              {
 308 //                fDigitStructPtr = &(fDigitContainerStructPtr->fDigitDataStruct[j+fDigitCount]);
 309 //                fDigitStructPtr->fX = xMod;
 310 //                fDigitStructPtr->fZ = zMod;
 311 //                fDigitStructPtr->fAmplitude = fEnergyArray[x][z][LOW_GAIN];
 312 //                fDigitStructPtr->fEnergy = fEnergyArray[x][z][LOW_GAIN] * fLowGainFactors[xMod][zMod];
 313 //                //TODO: fix time //CRAP
 314 //                fDigitStructPtr->fTime = fCellDataPtr->fTime  * 0.0000001;;
 315 //                fDigitStructPtr->fCrazyness = fCellDataPtr->fCrazyness;
 316 //                fDigitStructPtr->fModule = rcuData->fModuleID;
 317 //                j++;
 318 //              }
 319 //          }
 320 //      }
 321 //     }
 322
 323   fDigitCount += j;
 324   fDigitContainerStructPtr->fNDigits = fDigitCount;
 325   return fDigitCount;
 326 }
 327
 328
 329 void
 330 AliHLTPHOSDigitMaker::Reset()
 331 {
 332   //See header file for documentation
 333
 334   for(int x = 0; x < N_XCOLUMNS_RCU; x++)
 335     {
 336       for(int z = 0; z < N_ZROWS_RCU; z++)
 337         {
 338           for(int gain = 0; gain < N_GAINS; gain++)
 339             {
 340               fEnergyArray[x][z][gain] = 0;
 341             }
 342         }
 343     }
 344
 345   fDigitCount = 0;
 346 }
 347
 348 void
 349 AliHLTPHOSDigitMaker::SetGlobalHighGainFactor(Float_t factor)
 350 {
 351   //See header file for documentation
 352   for(int x = 0; x < N_XCOLUMNS_MOD; x++)
 353     {
 354       for(int z = 0; z < N_ZROWS_MOD; z++)
 355         {
 356           fHighGainFactors[x][z] = factor;
 357         }
 358     }
 359 }
 360
 361 void
 362 AliHLTPHOSDigitMaker::SetGlobalLowGainFactor(Float_t factor)
 363 {
 364   //See header file for documentation
 365   for(int x = 0; x < N_XCOLUMNS_MOD; x++)
 366     {
 367       for(int z = 0; z < N_ZROWS_MOD; z++)
 368         {
 369           fLowGainFactors[x][z] = factor;
 370         }
 371     }
 372 }
 373
 374 void
 375 AliHLTPHOSDigitMaker::SetDigitThresholds(const char* filepath, Int_t nSigmas)
 376 {
 377   //See header file for documentation
 378   TFile *histFile = new TFile(filepath);
 379
 380   TH2F *lgHist = (TH2F*)histFile->Get("RMSLGMapHist");
 381   TH2F *hgHist = (TH2F*)histFile->Get("RMSHGMapHist");
 382
 383   for(int x = 0; x < N_XCOLUMNS_MOD; x++)
 384     {
 385       for(int z = 0; z < N_ZROWS_MOD; z++)
 386         {
 387           fDigitThresholds[x][z][LOW_GAIN] = lgHist->GetBinContent(x, z) * nSigmas;
 388           fDigitThresholds[x][z][HIGH_GAIN] = hgHist->GetBinContent(x, z) * nSigmas;
 389         }
 390     }
 391 }
 392
 393 void
 394 AliHLTPHOSDigitMaker::SetDigitThresholds(const Float_t thresholdHG, const Float_t thresholdLG)
 395 {
 396   //See header file for documentation
 397   for(int x = 0; x < N_XCOLUMNS_MOD; x++)
 398     {
 399       for(int z = 0; z < N_ZROWS_MOD; z++)
 400         {
 401           fDigitThresholds[x][z][LOW_GAIN] = thresholdHG;
 402           fDigitThresholds[x][z][HIGH_GAIN] = thresholdLG;
 403         }
 404     }
 405 }
 406
 407
 408 void
 409 AliHLTPHOSDigitMaker::SetBadChannelMask(TH2F* badChannelHGHist, TH2F* badChannelLGHist, Float_t qCut)
 410 {
 411  for(int x = 0; x < N_XCOLUMNS_MOD; x++)
 412     {
 413       for(int z = 0; z < N_ZROWS_MOD; z++)
 414         {
 415           if(badChannelHGHist->GetBinContent(x, z) < qCut && badChannelHGHist->GetBinContent(x, z) > 0)
 416             {
 417               fBadChannelMask[x][z][HIGH_GAIN] = 1;
 418             }
 419           else
 420             {
 421               fBadChannelMask[x][z][HIGH_GAIN] = 0;
 422             }
 423           if(badChannelLGHist->GetBinContent(x, z) < qCut && badChannelLGHist->GetBinContent(x, z) > 0)
 424             {
 425               fBadChannelMask[x][z][LOW_GAIN] = 0;
 426             }
 427           else
 428             {
 429               fBadChannelMask[x][z][LOW_GAIN] = 0;
 430             }
 431         }
 432     }
 433 }