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