VZERO/AliVZERODigitizer.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   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 /* $Id$ */
  17
  18 ///_________________________________________________________________________
  19 ///
  20 /// This class constructs Digits out of Hits
  21 ///
  22 ///
  23
  24 // --- Standard library ---
  25
  26 // --- ROOT system ---
  27 #include <TMath.h>
  28 #include <TTree.h>
  29 #include <TMap.h>
  30 #include <TGeoManager.h>
  31 #include <TGeoPhysicalNode.h>
  32 #include <AliGeomManager.h>
  33 #include <TRandom.h>
  34
  35 // --- AliRoot header files ---
  36 #include "AliVZEROConst.h"
  37 #include "AliRun.h"
  38 #include "AliVZERO.h"
  39 #include "AliVZEROhit.h"
  40 #include "AliRunLoader.h"
  41 #include "AliLoader.h"
  42 #include "AliGRPObject.h"
  43 #include "AliRunDigitizer.h"
  44 #include "AliCDBManager.h"
  45 #include "AliCDBStorage.h"
  46 #include "AliCDBEntry.h"
  47 #include "AliVZEROCalibData.h"
  48
  49 #include "AliVZEROdigit.h"
  50 #include "AliVZERODigitizer.h"
  51
  52 ClassImp(AliVZERODigitizer)
  53
  54  AliVZERODigitizer::AliVZERODigitizer()
  55                    :AliDigitizer(),
  56                     fCalibData(GetCalibData()),
  57                     fPhotoCathodeEfficiency(0.18),
  58                     fPMVoltage(768.0),
  59                     fPMGain(TMath::Power((fPMVoltage / 112.5) ,7.04277)),
  60                     fNdigits(0),
  61                     fDigits(0),
  62                     fCollisionMode(0),
  63                     fBeamEnergy(0.)
  64
  65 {
  66   // default constructor
  67
  68 //    fNdigits = 0;
  69 //    fDigits  = 0;
  70 //
  71 //    fPhotoCathodeEfficiency =   0.18;
  72 //    fPMVoltage              =  768.0;
  73 //    fPMGain = TMath::Power((fPMVoltage / 112.5) ,7.04277);
  74
  75 //   fCalibData = GetCalibData();
  76 }
  77
  78 //____________________________________________________________________________
  79   AliVZERODigitizer::AliVZERODigitizer(AliRunDigitizer* manager)
  80                     :AliDigitizer(manager),
  81                      fCalibData(GetCalibData()),
  82                      fPhotoCathodeEfficiency(0.18),
  83                      fPMVoltage(768.0),
  84                      fPMGain(TMath::Power((fPMVoltage / 112.5) ,7.04277)),
  85                      fNdigits(0),
  86                      fDigits(0),
  87                      fCollisionMode(0),
  88                      fBeamEnergy(0.)
  89
  90 {
  91   // constructor
  92
  93 //   fNdigits = 0;
  94 //   fDigits  = 0;
  95 //
  96 //   fPhotoCathodeEfficiency =   0.18;
  97 //   fPMVoltage              =  768.0;
  98 //   fPMGain = TMath::Power( (fPMVoltage / 112.5) ,7.04277 );
  99
 100 //  fCalibData = GetCalibData();
 101
 102 }
 103
 104 //____________________________________________________________________________
 105   AliVZERODigitizer::~AliVZERODigitizer()
 106 {
 107   // destructor
 108
 109   if (fDigits) {
 110     fDigits->Delete();
 111     delete fDigits;
 112     fDigits=0;
 113   }
 114 }
 115
 116 //_____________________________________________________________________________
 117 Bool_t AliVZERODigitizer::Init()
 118 {
 119   // Initialises the digitizer
 120
 121   // Initialises the Digit array
 122   fDigits = new TClonesArray ("AliVZEROdigit", 1000);
 123
 124   //  TGeoHMatrix *im = AliGeomManager::GetMatrix("VZERO/V0C");
 125   //  im->Print();
 126
 127   GetCollisionMode();
 128   return kTRUE;
 129 }
 130
 131 //____________________________________________________________________________
 132 void AliVZERODigitizer::Exec(Option_t* /*option*/)
 133 {
 134   // Creates digits from hits
 135
 136   Int_t       map[80];    // 48 values on V0C + 32 on V0A
 137 //  Int_t       pmNumber[80];
 138   Int_t       adc[64];    // 32 PMs on V0C + 32 PMs on V0A
 139   Float_t     time[80], time_ref[80], time2[64];
 140   Float_t     adc_gain[80];
 141   Float_t     adc_pedestal[64],adc_sigma[64];
 142   fNdigits     =    0;
 143   Float_t     pmGain_smeared[64];
 144   Float_t     cPM[80];
 145
 146   // Smearing of the PM tubes intrinsic characteristics
 147
 148   for(Int_t ii=0; ii<64; ii++) {
 149       pmGain_smeared[ii] = gRandom->Gaus(fPMGain, fPMGain/5); }
 150
 151   // Retrieval of ADC gain values and pedestal information from local CDB
 152   // I use only the first 64th values of the calibration array in CDB
 153   // as I have no beam burst structure - odd or even beam burst number
 154
 155   // Reminder : We have 16 scintillating cells mounted on 8 PMs
 156   // on Ring 3 and Ring 4 in V0C -  added to produce  ADC outputs
 157   // on these rings...
 158
 159   for(Int_t i=0; i<16; i++) {
 160             adc_gain[i] = fCalibData->GetGain(i);
 161             cPM[i]      = fPhotoCathodeEfficiency * pmGain_smeared[i];}
 162
 163   for(Int_t j=16; j<48; j=j+2) {
 164             Int_t i=(j+17)/2;
 165             adc_gain[j]   = fCalibData->GetGain(i);
 166             adc_gain[j+1] = fCalibData->GetGain(i);
 167             cPM[j]        = fPhotoCathodeEfficiency * pmGain_smeared[i];
 168             cPM[j+1]      = fPhotoCathodeEfficiency * pmGain_smeared[i]; }
 169
 170   for(Int_t i=48; i<80; i++){
 171             adc_gain[i] = fCalibData->GetGain(i-16);
 172             cPM[i]      = fPhotoCathodeEfficiency * pmGain_smeared[i-16];};
 173
 174   for(Int_t  i=0; i<64; i++){ adc_pedestal[i] = fCalibData->GetPedestal(i);
 175                               adc_sigma[i]    = fCalibData->GetSigma(i); };
 176
 177 //  for(Int_t i=0; i<64; i++) { printf(" i = %d pedestal = %f sigma = %f \n\n",
 178 //                                       i, adc_pedestal[i], adc_sigma[i] );}
 179
 180   AliRunLoader* outRunLoader =
 181     AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
 182   if (!outRunLoader) {
 183     Error("Exec", "Can not get output Run Loader");
 184     return;}
 185
 186   AliLoader* outLoader = outRunLoader->GetLoader("VZEROLoader");
 187   if (!outLoader) {
 188     Error("Exec", "Can not get output VZERO Loader");
 189     return;}
 190
 191   const char* mode = "update";
 192   if(outRunLoader->GetEventNumber() == 0) mode = "recreate";
 193   outLoader->LoadDigits(mode);
 194   if (!outLoader->TreeD()) outLoader->MakeTree("D");
 195   outLoader->MakeDigitsContainer();
 196   TTree* treeD  = outLoader->TreeD();
 197   Int_t bufsize = 16000;
 198   treeD->Branch("VZERODigit", &fDigits, bufsize);
 199
 200   for (Int_t iInput = 0; iInput < fManager->GetNinputs(); iInput++) {
 201      AliRunLoader* runLoader =
 202      AliRunLoader::GetRunLoader(fManager->GetInputFolderName(iInput));
 203      AliLoader* loader = runLoader->GetLoader("VZEROLoader");
 204      if (!loader) {
 205        Error("Exec", "Can not get VZERO Loader for input %d", iInput);
 206        continue;}
 207
 208      if (!runLoader->GetAliRun()) runLoader->LoadgAlice();
 209
 210      AliVZERO* vzero = (AliVZERO*) runLoader->GetAliRun()->GetDetector("VZERO");
 211      if (!vzero) {
 212        Error("Exec", "No VZERO detector for input %d", iInput);
 213        continue;}
 214
 215      loader->LoadHits();
 216      TTree* treeH = loader->TreeH();
 217      if (!treeH) {
 218        Error("Exec", "Cannot get TreeH for input %d", iInput);
 219        continue; }
 220
 221      for(Int_t i=0; i<80; i++) {map[i] = 0; time[i] = 0.0;}
 222
 223      TClonesArray* hits = vzero->Hits();
 224
 225 //  Now makes Digits from hits
 226
 227      Int_t nTracks = (Int_t) treeH->GetEntries();
 228      for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
 229          for (Int_t i=0; i<80; i++) {time_ref[i] = 999999.0;}
 230          vzero->ResetHits();
 231          treeH->GetEvent(iTrack);
 232          Int_t nHits = hits->GetEntriesFast();
 233          for (Int_t iHit = 0; iHit < nHits; iHit++) {
 234              AliVZEROhit* hit = (AliVZEROhit *)hits->UncheckedAt(iHit);
 235              Int_t nPhot = hit->Nphot();
 236              Int_t cell  = hit->Cell();
 237              map[cell] += nPhot;
 238              Float_t dt_scintillator = gRandom->Gaus(0,0.7);
 239              Float_t t = dt_scintillator + 1e9*hit->Tof();
 240              if (t > 0.0) {
 241                 if(t < time_ref[cell]) time_ref[cell] = t;
 242                 time[cell] = TMath::Min(t,time_ref[cell]); }
 243          }           // hit   loop
 244      }             // track loop
 245
 246      loader->UnloadHits();
 247
 248   }               // input loop
 249
 250 // Now builds the scintillator cell response (80 cells i.e. 80 responses)
 251
 252    for (Int_t i=0; i<80; i++) {
 253         Float_t q1 = Float_t ( map[i] )* cPM[i] * kQe;
 254         Float_t noise = gRandom->Gaus(10.5,3.22);
 255         Float_t pmResponse  =  q1/kC*TMath::Power(ktheta/kthau,1/(1-ktheta/kthau))
 256         + noise*1e-3;
 257         if(fCollisionMode >0) adc_gain[i] = adc_gain[i]/70.0; // reduce dynamics in Ion Collision Mode
 258         map[i] = Int_t( pmResponse * adc_gain[i]);
 259         Float_t MIP = 1.0/fCalibData->GetMIPperADC(GetPMNumber(i));
 260         if(fCollisionMode >0) MIP=2.0;
 261 //      printf("cell = %d,  ADC = %d, TDC = %f \n",i,map[i], time[i]*10.0 );
 262         if(map[i] > (int(( MIP/2 ) + 0.5)) )
 263                   {map[i] = Int_t(gRandom->Gaus(map[i], (int(( MIP/6 ) + 0.5)) ));}
 264    }
 265
 266 // Now transforms 80 cell responses into 64 photomultiplier responses
 267 // Also adds the ADC pedestals taken out of the calibration data base
 268
 269    for (Int_t j=0; j<16; j++){
 270         adc[j]  = static_cast<Int_t>(map [j] + gRandom->Gaus(adc_pedestal[j], adc_sigma[j]));
 271         time2[j]= time[j];}
 272
 273    for (Int_t j=48; j<80; j++){
 274         adc[j-16]  = static_cast<Int_t>(map [j]
 275                                         + gRandom->Gaus(adc_pedestal[j-16],adc_sigma[j-16]));
 276         time2[j-16]= time[j]; }
 277
 278    for (Int_t j=0; j<16; j++){
 279         adc[16+j] = static_cast<Int_t>(map [16+2*j]+ map [16+2*j+1]
 280                                        + gRandom->Gaus(adc_pedestal[16+j], adc_sigma[16+j]));
 281         Float_t min_time = TMath::Min(time [16+2*j],time [16+2*j+1]);
 282         time2[16+j] = min_time;
 283         if(min_time==0.0){time2[16+j]=TMath::Max(time[16+2*j],time[16+2*j+1]);}
 284    }
 285
 286
 287 // Now add digits to the digit Tree
 288
 289    for (Int_t i=0; i<64; i++) {
 290       if(adc[i] > 0) {
 291 //           printf(" Event, cell, adc, tof = %d %d %d %f\n",
 292 //                    outRunLoader->GetEventNumber(),i, map[i], time2[i]*10.0);
 293 //           multiply by 10 to have 100 ps per channel :
 294              AddDigit(i, adc[i], Int_t((time2[i]*10.0) +0.5)) ;}
 295    }
 296
 297   treeD->Fill();
 298   outLoader->WriteDigits("OVERWRITE");
 299   outLoader->UnloadDigits();
 300   ResetDigit();
 301 }
 302
 303 //____________________________________________________________________________
 304 void AliVZERODigitizer::AddDigit(Int_t PMnumber, Int_t adc, Int_t time)
 305  {
 306
 307 // Adds Digit
 308
 309   TClonesArray &ldigits = *fDigits;
 310   new(ldigits[fNdigits++]) AliVZEROdigit(PMnumber,adc,time);
 311 }
 312 //____________________________________________________________________________
 313 void AliVZERODigitizer::ResetDigit()
 314 {
 315
 316 // Clears Digits
 317
 318   fNdigits = 0;
 319   if (fDigits) fDigits->Delete();
 320 }
 321
 322 //____________________________________________________________________________
 323 void AliVZERODigitizer::GetCollisionMode()
 324 {
 325 // Retrieves the collision mode from GRP data
 326
 327 // Initialization of the GRP entry
 328
 329    Int_t run = AliCDBManager::Instance()->GetRun();
 330
 331 //   printf("\n ++++++ Run Number retrieved as %d \n",run);
 332
 333   AliCDBEntry*  entry = AliCDBManager::Instance()->Get("GRP/GRP/Data",run);
 334   AliGRPObject* grpData = 0x0;
 335
 336   if(entry){
 337     TMap* m = dynamic_cast<TMap*>(entry->GetObject());  // old GRP entry
 338     if(m){
 339        m->Print();
 340        grpData = new AliGRPObject();
 341        grpData->ReadValuesFromMap(m);
 342     }
 343     else{
 344        grpData = dynamic_cast<AliGRPObject*>(entry->GetObject());  // new GRP entry
 345        entry->SetOwner(0);
 346     }
 347     AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
 348   }
 349
 350   if(!grpData) AliError("No GRP entry found in OCDB!");
 351
 352 // Retrieval of collision mode
 353
 354   TString beamType = grpData->GetBeamType();
 355   if(beamType==AliGRPObject::GetInvalidString()){
 356      AliError("GRP/GRP/Data entry:  missing value for the beam type !");
 357      AliError("\t VZERO cannot retrieve beam type\n");
 358      return;
 359   }
 360
 361    if( (beamType.CompareTo("P-P") ==0)  || (beamType.CompareTo("p-p") ==0) ){
 362        fCollisionMode=0;
 363   }
 364    else if( (beamType.CompareTo("Pb-Pb") ==0)  || (beamType.CompareTo("A-A") ==0) ){
 365        fCollisionMode=1;
 366    }
 367
 368   fBeamEnergy = grpData->GetBeamEnergy();
 369   if(fBeamEnergy==AliGRPObject::GetInvalidFloat()) {
 370      AliError("GRP/GRP/Data entry:  missing value for the beam energy ! Using 0");
 371      fBeamEnergy = 0.;
 372   }
 373
 374 //     printf("\n ++++++ Beam type and collision mode retrieved as %s %d @ %1.3f GeV ++++++\n\n",beamType.Data(), fCollisionMode, fBeamEnergy);
 375
 376 }
 377
 378 //____________________________________________________________________________
 379 AliVZEROCalibData* AliVZERODigitizer::GetCalibData() const
 380
 381 {
 382   AliCDBManager *man = AliCDBManager::Instance();
 383
 384   AliCDBEntry *entry=0;
 385
 386   entry = man->Get("VZERO/Calib/Data");
 387
 388 //   if(!entry){
 389 //     AliWarning("Load of calibration data from default storage failed!");
 390 //     AliWarning("Calibration data will be loaded from local storage ($ALICE_ROOT)");
 391 //     Int_t runNumber = man->GetRun();
 392 //     entry = man->GetStorage("local://$ALICE_ROOT/OCDB")
 393 //       ->Get("VZERO/Calib/Data",runNumber);
 394 //
 395 //   }
 396
 397   // Retrieval of data in directory VZERO/Calib/Data:
 398
 399
 400   AliVZEROCalibData *calibdata = 0;
 401
 402   if (entry) calibdata = (AliVZEROCalibData*) entry->GetObject();
 403   if (!calibdata)  AliFatal("No calibration data from calibration database !");
 404
 405   return calibdata;
 406
 407 }
 408
 409 //____________________________________________________________________________
 410 Int_t AliVZERODigitizer::GetPMNumber(Int_t cell) const
 411
 412 {
 413
 414   Int_t pmNumber[80] = { 0,  1,  2,  3,  4,  5,  6,  7,
 415                           8,  9, 10, 11, 12, 13, 14, 15,
 416                          16, 16, 17, 17, 18, 18, 19, 19, 20, 20, 21, 21, 22, 22, 23, 23,
 417                          24, 24, 25, 25, 26, 26, 27, 27, 28, 28, 29, 29, 30, 30, 31, 31,
 418                          32, 33, 34, 35, 36, 37, 38, 39,
 419                          40, 41, 42, 43, 44, 45, 46, 47,
 420                          48, 49, 50, 51, 52, 53, 54, 55,
 421                          56, 57, 58, 59, 60, 61, 62, 63 };
 422
 423   return pmNumber[cell];
 424 }
 425
 426