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