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   //
 189   //Loop over segments of the TPC
 190
 191   for (Int_t segmentID=0; segmentID<param->GetNRowsTotal(); segmentID++)
 192    {
 193     Int_t sec, row;
 194     if (!param->AdjustSectorRow(segmentID,sec,row))
 195      {
 196       cerr<<"AliTPC warning: invalid segment ID ! "<<segmentID<<endl;
 197       continue;
 198      }
 199     AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec);  // pad gains per given sector
 200     digrow->SetID(segmentID);
 201
 202     Int_t nrows = 0;
 203     Int_t ncols = 0;
 204
 205     Bool_t digitize = kFALSE;
 206     for (Int_t i=0;i<nInputs; i++)
 207      {
 208
 209       rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
 210       gime = rl->GetLoader("TPCLoader");
 211
 212       if (gime->TreeS()->GetEntryWithIndex(segmentID,segmentID) >= 0) {
 213         digarr[i]->ExpandBuffer();
 214         digarr[i]->ExpandTrackBuffer();
 215         nrows = digarr[i]->GetNRows();
 216         ncols = digarr[i]->GetNCols();
 217         active[i] = kTRUE;
 218         if (!fRegionOfInterest || (i == 0)) digitize = kTRUE;
 219       } else {
 220         active[i] = kFALSE;
 221       }
 222       if (fRegionOfInterest && !digitize) break;
 223      }
 224     if (!digitize) continue;
 225
 226     digrow->Allocate(nrows,ncols);
 227     digrow->AllocateTrack(3);
 228
 229     Float_t q=0;
 230     Int_t label[1000]; //stack for 300 events
 231     Int_t labptr = 0;
 232
 233     Int_t nElems = nrows*ncols;
 234
 235     for (Int_t i=0;i<nInputs; i++)
 236      if (active[i]) {
 237        pdig[i] = digarr[i]->GetDigits();
 238        ptr[i]  = digarr[i]->GetTracks();
 239       }
 240
 241     Short_t *pdig1= digrow->GetDigits();
 242     Int_t   *ptr1= digrow->GetTracks() ;
 243
 244
 245
 246     for (Int_t elem=0;elem<nElems; elem++)
 247      {
 248
 249        q=0;
 250        labptr=0;
 251        // looop over digits
 252         for (Int_t i=0;i<nInputs; i++) if (active[i])
 253          {
 254           //          q  += digarr[i]->GetDigitFast(rows,col);
 255             q  += *(pdig[i]);
 256
 257            for (Int_t tr=0;tr<3;tr++)
 258             {
 259              //             Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
 260              Int_t lab = ptr[i][tr*nElems];
 261              if ( (lab > 1) && *(pdig[i])>zerosup)
 262               {
 263                 label[labptr]=lab+masks[i];
 264                 labptr++;
 265               }
 266             }
 267            pdig[i]++;
 268            ptr[i]++;
 269          }
 270         q/=16.;  //conversion factor
 271         Float_t gain = gainROC->GetValue(row,elem/nrows);  // get gain for given - pad-row pad
 272         if (gain<0.5){
 273           printf("problem\n");
 274         }
 275         q*= gain;
 276         //       Float_t noise  = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
 277         Float_t noise  = pTPC->GetNoise();
 278         q+=noise;
 279         q=TMath::Nint(q);
 280         if (q > zerosup)
 281          {
 282           if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
 283           //digrow->SetDigitFast((Short_t)q,rows,col);
 284           *pdig1 =Short_t(q);
 285           for (Int_t tr=0;tr<3;tr++)
 286            {
 287             if (tr<labptr)
 288              ptr1[tr*nElems] = label[tr];
 289            }
 290           }
 291         pdig1++;
 292         ptr1++;
 293      }
 294
 295     digrow->CompresBuffer(1,zerosup);
 296     digrow->CompresTrackBuffer(1);
 297     tree->Fill();
 298     if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
 299    } //for (Int_t n=0; n<param->GetNRowsTotal(); n++)
 300
 301
 302   orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
 303   ogime = orl->GetLoader("TPCLoader");
 304   ogime->WriteDigits("OVERWRITE");
 305
 306   //fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
 307
 308   delete digrow;
 309   for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
 310   delete []masks;
 311   delete []pdig;
 312   delete []ptr;
 313   delete []active;
 314   delete []digarr;
 315 }
 316
 317
 318
 319 //------------------------------------------------------------------------
 320 void AliTPCDigitizer::ExecSave(Option_t* option)
 321 {
 322
 323   // merge input tree's with summable digits
 324   //output stored in TreeTPCD
 325
 326   TString optionString = option;
 327   if (optionString.Data() == "deb") {
 328     cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
 329     fDebug = 3;
 330   }
 331   //get detector and geometry
 332   AliRunLoader *rl, *orl;
 333   AliLoader *gime, *ogime;
 334
 335
 336   orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
 337   ogime = orl->GetLoader("TPCLoader");
 338
 339   rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
 340   gime = rl->GetLoader("TPCLoader");
 341
 342   rl->LoadgAlice();
 343   AliRun* alirun = rl->GetAliRun();
 344
 345   AliTPC *pTPC  = (AliTPC *) alirun->GetModule("TPC");
 346   AliTPCParam * param = pTPC->GetParam();
 347   pTPC->GenerNoise(500000); //create teble with noise
 348   printf("noise %f \n",  param->GetNoise()*param->GetNoiseNormFac());
 349   //
 350   Int_t nInputs = fManager->GetNinputs();
 351   Int_t * masks = new Int_t[nInputs];
 352   for (Int_t i=0; i<nInputs;i++)
 353     masks[i]= fManager->GetMask(i);
 354
 355   AliSimDigits ** digarr = new AliSimDigits*[nInputs];
 356   for (Int_t i1=0;i1<nInputs; i1++)
 357    {
 358     digarr[i1]=0;
 359     //    intree[i1]
 360     rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i1));
 361     gime = rl->GetLoader("TPCLoader");
 362
 363     TTree * treear =  gime->TreeS();
 364     TBranch * br = treear->GetBranch("fSegmentID");
 365     if (br) br->GetFile()->cd();
 366     if (!treear) {
 367       cerr<<" TPC -  not existing input = \n"<<i1<<" ";
 368     }
 369     treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
 370   }
 371
 372   rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
 373   gime = rl->GetLoader("TPCLoader");
 374   Stat_t nentries = gime->TreeS()->GetEntries();
 375
 376
 377   //create branch's in TPC treeD
 378   AliSimDigits * digrow = new AliSimDigits;
 379   TTree * tree  = ogime->TreeD();
 380
 381   tree->Branch("Segment","AliSimDigits",&digrow);
 382   param->SetZeroSup(2);
 383
 384   Int_t zerosup = param->GetZeroSup();
 385   //Loop over segments of the TPC
 386
 387   AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetPadGainFactor();
 388   for (Int_t n=0; n<nentries; n++) {
 389     rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
 390     gime = rl->GetLoader("TPCLoader");
 391     gime->TreeS()->GetEvent(n);
 392
 393     digarr[0]->ExpandBuffer();
 394     digarr[0]->ExpandTrackBuffer();
 395
 396
 397     for (Int_t i=1;i<nInputs; i++){
 398 //      fManager->GetInputTreeTPCS(i)->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
 399       rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
 400       gime = rl->GetLoader("TPCLoader");
 401       gime->TreeS()->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());
 402       digarr[i]->ExpandBuffer();
 403       digarr[i]->ExpandTrackBuffer();
 404       if ((digarr[0]->GetID()-digarr[i]->GetID())>0)
 405        printf("problem\n");
 406
 407     }
 408
 409     Int_t sec, row;
 410     if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
 411       cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
 412       continue;
 413     }
 414
 415     AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec);  // pad gains per given sector
 416     digrow->SetID(digarr[0]->GetID());
 417
 418     Int_t nrows = digarr[0]->GetNRows();
 419     Int_t ncols = digarr[0]->GetNCols();
 420     digrow->Allocate(nrows,ncols);
 421     digrow->AllocateTrack(3);
 422
 423     Float_t q=0;
 424     Int_t label[1000]; //stack for 300 events
 425     Int_t labptr = 0;
 426
 427
 428
 429     for (Int_t rows=0;rows<nrows; rows++){
 430       for (Int_t col=0;col<ncols; col++){
 431
 432        q=0;
 433        labptr=0;
 434        // looop over digits
 435         for (Int_t i=0;i<nInputs; i++){
 436          q  += digarr[i]->GetDigitFast(rows,col);
 437           //q  += *(pdig[i]);
 438
 439           for (Int_t tr=0;tr<3;tr++) {
 440            Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
 441            //Int_t lab = ptr[i][tr*nElems];
 442             if ( (lab > 1) ) {
 443               label[labptr]=lab+masks[i];
 444               labptr++;
 445             }
 446           }
 447          // pdig[i]++;
 448          //ptr[i]++;
 449
 450         }
 451        q/=16.;  //conversion factor
 452        //       Float_t noise  = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());
 453        Float_t gain = gainROC->GetValue(row,col);
 454        q*= gain;
 455        Float_t noise  = pTPC->GetNoise();
 456        q+=noise;
 457         q=TMath::Nint(q);
 458         if (q > zerosup){
 459
 460          if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
 461          digrow->SetDigitFast((Short_t)q,rows,col);
 462          // *pdig1 =Short_t(q);
 463          for (Int_t tr=0;tr<3;tr++){
 464            if (tr<labptr)
 465              ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
 466            //ptr1[tr*nElems] = label[tr];
 467            //else
 468              //           ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);
 469            //  ptr1[tr*nElems] = 1;
 470          }
 471        }
 472        //pdig1++;
 473        //ptr1++;
 474     }
 475     }
 476
 477     digrow->CompresBuffer(1,zerosup);
 478     digrow->CompresTrackBuffer(1);
 479     tree->Fill();
 480     if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";
 481   }
 482 //  printf("end TPC merging - end -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
 483   //fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
 484     ogime->WriteDigits("OVERWRITE");
 485
 486     for (Int_t i=1;i<nInputs; i++)
 487      {
 488       rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
 489       gime = rl->GetLoader("TPCLoader");
 490       gime->UnloadSDigits();
 491      }
 492     ogime->UnloadDigits();
 493
 494   delete digrow;
 495   for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
 496   delete []masks;
 497   delete digarr;
 498 }