TPC/AliTPCDigitizer.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2000, 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   Class for creating of the sumable digits and digits from MC data
  20   //
  21   The input :  ideal signals (Hits->Diffusion->Attachment -Ideal signal)
  22   The output:  raw digits
  23
  24   Effect implemented:
  25   1. Pad by pad gain map
  26   2. Noise map
  27   3. The dead channels identified  - zerro noise for corresponding pads
  28      In this case the outpu equal zerro
  29
  30 */
  31
  32
  33
  34
  35 #include <stdlib.h>
  36 #include <TTree.h>
  37 #include <TObjArray.h>
  38 #include <TFile.h>
  39 #include <TDirectory.h>
  40 #include <Riostream.h>
  41
  42 #include "AliTPCDigitizer.h"
  43
  44 #include "AliTPC.h"
  45 #include "AliTPCParam.h"
  46 #include "AliTPCParamSR.h"
  47 #include "AliRun.h"
  48 #include "AliPDG.h"
  49 #include "AliRunDigitizer.h"
  50 #include "AliSimDigits.h"
  51 #include "AliLog.h"
  52
  53 #include "AliTPCcalibDB.h"
  54 #include "AliTPCCalPad.h"
  55 #include "AliTPCCalROC.h"
  56
  57 ClassImp(AliTPCDigitizer)
  58
  59 //___________________________________________
  60   AliTPCDigitizer::AliTPCDigitizer() :AliDigitizer(),fDebug(0)
  61 {
  62   //
  63 // Default ctor - don't use it
  64 //
  65
  66 }
  67
  68 //___________________________________________
  69 AliTPCDigitizer::AliTPCDigitizer(AliRunDigitizer* manager)
  70   :AliDigitizer(manager),fDebug(0)
  71 {
  72   //
  73 // ctor which should be used
  74 //
  75   AliDebug(2,"(AliRunDigitizer* manager) was processed");
  76 }
  77
  78 //------------------------------------------------------------------------
  79 AliTPCDigitizer::~AliTPCDigitizer()
  80 {
  81 // Destructor
  82 }
  83
  84
  85
  86 //------------------------------------------------------------------------
  87 Bool_t AliTPCDigitizer::Init()
  88 {
  89 // Initialization
  90
  91  return kTRUE;
  92 }
  93
  94
  95 //------------------------------------------------------------------------
  96 void AliTPCDigitizer::Exec(Option_t* option)
  97 {
  98   ExecFast(option);
  99 }
 100 //------------------------------------------------------------------------
 101 void AliTPCDigitizer::ExecFast(Option_t* option)
 102 {
 103
 104   // merge input tree's with summable digits
 105   //output stored in TreeTPCD
 106   char s[100];
 107   char ss[100];
 108   TString optionString = option;
 109   if (!strcmp(optionString.Data(),"deb")) {
 110     cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
 111     fDebug = 3;
 112   }
 113   //get detector and geometry
 114
 115
 116   AliRunLoader *rl, *orl;
 117   AliLoader *gime, *ogime;
 118
 119   if (gAlice == 0x0)
 120    {
 121      Warning("ExecFast","gAlice is NULL. Loading from input 0");
 122      rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
 123      if (rl == 0x0)
 124       {
 125         Error("ExecFast","Can not find Run Loader for input 0. Can not proceed.");
 126         return;
 127       }
 128      rl->LoadgAlice();
 129      rl->GetAliRun();
 130    }
 131   AliTPC *pTPC  = (AliTPC *) gAlice->GetModule("TPC");
 132   AliTPCParam * param = pTPC->GetParam();
 133
 134   sprintf(s,param->GetTitle());
 135   sprintf(ss,"75x40_100x60");
 136   if(strcmp(s,ss)==0){
 137     printf("2 pad-length geom hits with 3 pad-lenght geom digits...\n");
 138     delete param;
 139     param=new AliTPCParamSR();
 140   }
 141   else{
 142    sprintf(ss,"75x40_100x60_150x60");
 143    if(strcmp(s,ss)!=0) {
 144      printf("No TPC parameters found...\n");
 145      exit(2);
 146    }
 147   }
 148
 149   pTPC->GenerNoise(500000); //create teble with noise
 150   //
 151   Int_t nInputs = fManager->GetNinputs();
 152   Int_t * masks = new Int_t[nInputs];
 153   for (Int_t i=0; i<nInputs;i++)
 154     masks[i]= fManager->GetMask(i);
 155   Short_t **pdig= new Short_t*[nInputs];   //pointers to the expanded digits array
 156   Int_t **ptr=  new Int_t*[nInputs];       //pointers to teh expanded tracks array
 157   Bool_t *active=  new Bool_t[nInputs];    //flag for active input segments
 158
 159
 160   //create digits array for given sectors
 161   // make indexes
 162
 163   AliSimDigits ** digarr = new AliSimDigits*[nInputs];
 164   for (Int_t i1=0;i1<nInputs; i1++)
 165     {
 166       digarr[i1]=0;
 167      //    intree[i1]
 168       rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i1));
 169       gime = rl->GetLoader("TPCLoader");
 170       gime->LoadSDigits("read");
 171       TTree * treear =  gime->TreeS();
 172
 173       if (!treear)
 174        {
 175         cerr<<"AliTPCDigitizer: Input tree with SDigits not found in"
 176             <<" input "<< i1<<endl;
 177         return;
 178        }
 179
 180       if (treear->GetIndex()==0)
 181         treear->BuildIndex("fSegmentID","fSegmentID");
 182       treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
 183     }
 184
 185
 186   //create branch's in TPC treeD
 187   AliSimDigits * digrow = new AliSimDigits;
 188
 189   orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
 190   ogime = orl->GetLoader("TPCLoader");
 191
 192   TTree * tree  = ogime->TreeD();
 193   if (tree == 0x0)
 194    {
 195      ogime->MakeTree("D");
 196      tree  = ogime->TreeD();
 197    }
 198   tree->Branch("Segment","AliSimDigits",&digrow);
 199   //
 200
 201   param->SetZeroSup(2);
 202
 203   Int_t zerosup = param->GetZeroSup();
 204   AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor();
 205   AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
 206   //
 207   //Loop over segments of the TPC
 208
 209   for (Int_t segmentID=0; segmentID<param->GetNRowsTotal(); segmentID++)
 210    {
 211     Int_t sec, row;
 212     if (!param->AdjustSectorRow(segmentID,sec,row))
 213      {
 214       cerr<<"AliTPC warning: invalid segment ID ! "<<segmentID<<endl;
 215       continue;
 216      }
 217     AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec);  // pad gains per given sector
 218     AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(sec);  // noise per given sector
 219     digrow->SetID(segmentID);
 220
 221     Int_t nrows = 0;
 222     Int_t ncols = 0;
 223
 224     Bool_t digitize = kFALSE;
 225     for (Int_t i=0;i<nInputs; i++)
 226      {
 227
 228       rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
 229       gime = rl->GetLoader("TPCLoader");
 230
 231       if (gime->TreeS()->GetEntryWithIndex(segmentID,segmentID) >= 0) {
 232         digarr[i]->ExpandBuffer();
 233         digarr[i]->ExpandTrackBuffer();
 234         nrows = digarr[i]->GetNRows();
 235         ncols = digarr[i]->GetNCols();
 236         active[i] = kTRUE;
 237         if (!fRegionOfInterest || (i == 0)) digitize = kTRUE;
 238       } else {
 239         active[i] = kFALSE;
 240       }
 241       if (fRegionOfInterest && !digitize) break;
 242      }
 243     if (!digitize) continue;
 244
 245     digrow->Allocate(nrows,ncols);
 246     digrow->AllocateTrack(3);
 247
 248     Float_t q=0;
 249     Int_t label[1000]; //stack for 300 events
 250     Int_t labptr = 0;
 251
 252     Int_t nElems = nrows*ncols;
 253
 254     for (Int_t i=0;i<nInputs; i++)
 255      if (active[i]) {
 256        pdig[i] = digarr[i]->GetDigits();
 257        ptr[i]  = digarr[i]->GetTracks();
 258       }
 259
 260     Short_t *pdig1= digrow->GetDigits();
 261     Int_t   *ptr1= digrow->GetTracks() ;
 262
 263
 264
 265     for (Int_t elem=0;elem<nElems; elem++)
 266      {
 267
 268        q=0;
 269        labptr=0;
 270        // looop over digits
 271         for (Int_t i=0;i<nInputs; i++) if (active[i])
 272          {
 273           //          q  += digarr[i]->GetDigitFast(rows,col);
 274             q  += *(pdig[i]);
 275
 276            for (Int_t tr=0;tr<3;tr++)
 277             {
 278              //             Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
 279              Int_t lab = ptr[i][tr*nElems];
 280              if ( (lab > 1) && *(pdig[i])>zerosup)
 281               {
 282                 label[labptr]=lab+masks[i];
 283                 labptr++;
 284               }
 285             }
 286            pdig[i]++;
 287            ptr[i]++;
 288          }
 289         q/=16.;  //conversion factor
 290         Float_t gain = gainROC->GetValue(row,elem/nrows);  // get gain for given - pad-row pad
 291         //if (gain<0.5){
 292           //printf("problem\n");
 293         //}
 294         q*= gain;
 295         Float_t noisePad = noiseROC->GetValue(row,elem/nrows);
 296         //       Float_t noise  = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
 297         Float_t noise  = pTPC->GetNoise();
 298         q+=noise*noisePad;
 299         q=TMath::Nint(q);
 300         if (q > zerosup)
 301          {
 302           if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
 303           //digrow->SetDigitFast((Short_t)q,rows,col);
 304           *pdig1 =Short_t(q);
 305           for (Int_t tr=0;tr<3;tr++)
 306            {
 307             if (tr<labptr)
 308              ptr1[tr*nElems] = label[tr];
 309            }
 310           }
 311         pdig1++;
 312         ptr1++;
 313      }
 314     //
 315     //  glitch filter
 316     //
 317     digrow->GlitchFilter();
 318     //
 319     digrow->CompresBuffer(1,zerosup);
 320     digrow->CompresTrackBuffer(1);
 321     tree->Fill();
 322     if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
 323    } //for (Int_t n=0; n<param->GetNRowsTotal(); n++)
 324
 325
 326   orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
 327   ogime = orl->GetLoader("TPCLoader");
 328   ogime->WriteDigits("OVERWRITE");
 329
 330   //fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
 331
 332   delete digrow;
 333   for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
 334   delete []masks;
 335   delete []pdig;
 336   delete []ptr;
 337   delete []active;
 338   delete []digarr;
 339 }
 340
 341
 342
 343 //------------------------------------------------------------------------
 344 void AliTPCDigitizer::ExecSave(Option_t* option)
 345 {
 346
 347   // merge input tree's with summable digits
 348   //output stored in TreeTPCD
 349
 350   TString optionString = option;
 351   if (!strcmp(optionString.Data(),"deb")) {
 352     cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
 353     fDebug = 3;
 354   }
 355   //get detector and geometry
 356   AliRunLoader *rl, *orl;
 357   AliLoader *gime, *ogime;
 358
 359
 360   orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
 361   ogime = orl->GetLoader("TPCLoader");
 362
 363   rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
 364   gime = rl->GetLoader("TPCLoader");
 365
 366   rl->LoadgAlice();
 367   AliRun* alirun = rl->GetAliRun();
 368
 369   AliTPC *pTPC  = (AliTPC *) alirun->GetModule("TPC");
 370   AliTPCParam * param = pTPC->GetParam();
 371   pTPC->GenerNoise(500000); //create teble with noise
 372   printf("noise %f \n",  param->GetNoise()*param->GetNoiseNormFac());
 373   //
 374   Int_t nInputs = fManager->GetNinputs();
 375   Int_t * masks = new Int_t[nInputs];
 376   for (Int_t i=0; i<nInputs;i++)
 377     masks[i]= fManager->GetMask(i);
 378
 379   AliSimDigits ** digarr = new AliSimDigits*[nInputs];
 380   for (Int_t i1=0;i1<nInputs; i1++)
 381    {
 382     digarr[i1]=0;
 383     //    intree[i1]
 384     rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i1));
 385     gime = rl->GetLoader("TPCLoader");
 386
 387     TTree * treear =  gime->TreeS();
 388     TBranch * br = treear->GetBranch("fSegmentID");
 389     if (br) br->GetFile()->cd();
 390     if (!treear) {
 391       cerr<<" TPC -  not existing input = \n"<<i1<<" ";
 392     }
 393     treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
 394   }
 395
 396   rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
 397   gime = rl->GetLoader("TPCLoader");
 398   Stat_t nentries = gime->TreeS()->GetEntries();
 399
 400
 401   //create branch's in TPC treeD
 402   AliSimDigits * digrow = new AliSimDigits;
 403   TTree * tree  = ogime->TreeD();
 404
 405   tree->Branch("Segment","AliSimDigits",&digrow);
 406   param->SetZeroSup(2);
 407
 408   Int_t zerosup = param->GetZeroSup();
 409   //Loop over segments of the TPC
 410
 411   AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor();
 412   AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
 413   for (Int_t n=0; n<nentries; n++) {
 414     rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
 415     gime = rl->GetLoader("TPCLoader");
 416     gime->TreeS()->GetEvent(n);
 417
 418     digarr[0]->ExpandBuffer();
 419     digarr[0]->ExpandTrackBuffer();
 420
 421
 422     for (Int_t i=1;i<nInputs; i++){
 423 //      fManager->GetInputTreeTPCS(i)->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
 424       rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
 425       gime = rl->GetLoader("TPCLoader");
 426       gime->TreeS()->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
 427       digarr[i]->ExpandBuffer();
 428       digarr[i]->ExpandTrackBuffer();
 429       if ((digarr[0]->GetID()-digarr[i]->GetID())>0)
 430        printf("problem\n");
 431
 432     }
 433
 434     Int_t sec, row;
 435     if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
 436       cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
 437       continue;
 438     }
 439
 440     AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec);  // pad gains per given sector
 441     AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(sec);  // noise per given sector
 442     digrow->SetID(digarr[0]->GetID());
 443
 444     Int_t nrows = digarr[0]->GetNRows();
 445     Int_t ncols = digarr[0]->GetNCols();
 446     digrow->Allocate(nrows,ncols);
 447     digrow->AllocateTrack(3);
 448
 449     Float_t q=0;
 450     Int_t label[1000]; //stack for 300 events
 451     Int_t labptr = 0;
 452
 453
 454
 455     for (Int_t rows=0;rows<nrows; rows++){
 456       for (Int_t col=0;col<ncols; col++){
 457
 458        q=0;
 459        labptr=0;
 460        // looop over digits
 461         for (Int_t i=0;i<nInputs; i++){
 462          q  += digarr[i]->GetDigitFast(rows,col);
 463           //q  += *(pdig[i]);
 464
 465           for (Int_t tr=0;tr<3;tr++) {
 466            Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
 467            //Int_t lab = ptr[i][tr*nElems];
 468             if ( (lab > 1) ) {
 469               label[labptr]=lab+masks[i];
 470               labptr++;
 471             }
 472           }
 473          // pdig[i]++;
 474          //ptr[i]++;
 475
 476         }
 477        q/=16.;  //conversion factor
 478        //       Float_t noise  = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
 479        Float_t gain = gainROC->GetValue(row,col);
 480        q*= gain;
 481        Float_t noisePad = noiseROC->GetValue(row, col);
 482
 483        Float_t noise  = pTPC->GetNoise();
 484        q+=noise*noisePad;
 485
 486         q=TMath::Nint(q);
 487         if (q > zerosup){
 488
 489          if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
 490          digrow->SetDigitFast((Short_t)q,rows,col);
 491          // *pdig1 =Short_t(q);
 492          for (Int_t tr=0;tr<3;tr++){
 493            if (tr<labptr)
 494              ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
 495            //ptr1[tr*nElems] = label[tr];
 496            //else
 497              //           ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);
 498            //  ptr1[tr*nElems] = 1;
 499          }
 500        }
 501        //pdig1++;
 502        //ptr1++;
 503     }
 504     }
 505
 506     digrow->CompresBuffer(1,zerosup);
 507     digrow->CompresTrackBuffer(1);
 508     tree->Fill();
 509     if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
 510   }
 511 //  printf("end TPC merging - end -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
 512   //fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
 513     ogime->WriteDigits("OVERWRITE");
 514
 515     for (Int_t i=1;i<nInputs; i++)
 516      {
 517       rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
 518       gime = rl->GetLoader("TPCLoader");
 519       gime->UnloadSDigits();
 520      }
 521     ogime->UnloadDigits();
 522
 523   delete digrow;
 524   for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
 525   delete []masks;
 526   delete digarr;
 527 }