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