]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCDigitizer.cxx
Update TPCCEda to write output file in parts (to avoid too big files produced in...
[u/mrichter/AliRoot.git] / TPC / AliTPCDigitizer.cxx
index 7d54fa91f333744cf2d473cb806fad9b20c4d50d..ed62f3033bc403eb67d424ea2ab9cb82a797606e 100644 (file)
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
+/* $Id$ */
+
+/*
+  Class for creating of the sumable digits and digits from MC data
+  //
+  The input :  ideal signals (Hits->Diffusion->Attachment -Ideal signal)
+  The output:  raw digits
+
+  Effect implemented:
+  1. Pad by pad gain map
+  2. Noise map
+  3. The dead channels identified  - zerro noise for corresponding pads
+     In this case the outpu equal zerro
+*/
+
 
 
+
+#include <stdlib.h>
 #include <TTree.h> 
 #include <TObjArray.h>
 #include <TFile.h>
 #include <TDirectory.h>
-#include <iostream.h>
+#include <Riostream.h>
+#include <TParameter.h>
 
 #include "AliTPCDigitizer.h"
 
 #include "AliTPC.h"
 #include "AliTPCParam.h"
+#include "AliTPCParamSR.h" 
 #include "AliRun.h"
 #include "AliPDG.h"
 #include "AliRunDigitizer.h"
 #include "AliSimDigits.h"
+#include "AliLog.h"
 
+#include "AliTPCcalibDB.h"
+#include "AliTPCCalPad.h"
+#include "AliTPCCalROC.h"
 
 ClassImp(AliTPCDigitizer)
 
 //___________________________________________
-AliTPCDigitizer::AliTPCDigitizer() :AliDigitizer()
+  AliTPCDigitizer::AliTPCDigitizer() :AliDigitizer(),fDebug(0)
 {
+  //
 // Default ctor - don't use it
+//
+  
 }
 
 //___________________________________________
 AliTPCDigitizer::AliTPCDigitizer(AliRunDigitizer* manager) 
-    :AliDigitizer(manager)
+  :AliDigitizer(manager),fDebug(0)
 {
+  //
 // ctor which should be used
-  if (GetDebug()>2) 
-    cerr<<"AliTPCDigitizer::AliTPCDigitizer"
-       <<"(AliRunDigitizer* manager) was processed"<<endl;
+//  
+  AliDebug(2,"(AliRunDigitizer* manager) was processed");
 }
 
 //------------------------------------------------------------------------
@@ -68,19 +95,292 @@ Bool_t AliTPCDigitizer::Init()
 
 //------------------------------------------------------------------------
 void AliTPCDigitizer::Exec(Option_t* option)
+{
+  ExecFast(option);  
+}
+//------------------------------------------------------------------------
+void AliTPCDigitizer::ExecFast(Option_t* option)
 {
   
   // merge input tree's with summable digits
   //output stored in TreeTPCD
-
+  char s[100]; 
+  char ss[100];
   TString optionString = option;
-  if (optionString.Data() == "deb") {
+  if (!strcmp(optionString.Data(),"deb")) {
     cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
     fDebug = 3;
   }
   //get detector and geometry
+
+
+  AliRunLoader *rl, *orl;
+  AliLoader *gime, *ogime;
+  
+  if (gAlice == 0x0)
+   {
+     Warning("ExecFast","gAlice is NULL. Loading from input 0");
+     rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
+     if (rl == 0x0)
+      {
+        Error("ExecFast","Can not find Run Loader for input 0. Can not proceed.");
+        return;
+      }
+     rl->LoadgAlice();
+     rl->GetAliRun();
+   }
   AliTPC *pTPC  = (AliTPC *) gAlice->GetModule("TPC");
   AliTPCParam * param = pTPC->GetParam();
+  
+  sprintf(s,param->GetTitle());
+  sprintf(ss,"75x40_100x60");
+  if(strcmp(s,ss)==0){
+    printf("2 pad-length geom hits with 3 pad-lenght geom digits...\n");
+    delete param;
+    param=new AliTPCParamSR();
+  }
+  else{
+   sprintf(ss,"75x40_100x60_150x60");
+   if(strcmp(s,ss)!=0) {
+     printf("No TPC parameters found...\n");
+     exit(2); 
+   }
+  }
+  
+  pTPC->GenerNoise(500000); //create table with noise
+  //
+  Int_t nInputs = fManager->GetNinputs();
+  Int_t * masks = new Int_t[nInputs];
+  for (Int_t i=0; i<nInputs;i++)
+    masks[i]= fManager->GetMask(i);
+  Short_t **pdig= new Short_t*[nInputs];   //pointers to the expanded digits array
+  Int_t **ptr=  new Int_t*[nInputs];       //pointers to the expanded tracks array
+  Bool_t *active=  new Bool_t[nInputs];    //flag for active input segments
+  Char_t phname[100];
+  
+  //create digits array for given sectors
+  // make indexes
+  //
+  //create branch's in TPC treeD
+  orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
+  ogime = orl->GetLoader("TPCLoader");
+  TTree * tree  = ogime->TreeD();
+  AliSimDigits * digrow = new AliSimDigits;  
+
+  if (tree == 0x0)
+   {
+     ogime->MakeTree("D");
+     tree  = ogime->TreeD();
+   }
+  tree->Branch("Segment","AliSimDigits",&digrow);
+  //  
+  AliSimDigits ** digarr = new AliSimDigits*[nInputs]; 
+  for (Int_t i1=0;i1<nInputs; i1++)
+    {
+      digarr[i1]=0;
+     //    intree[i1]
+      rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i1));
+      gime = rl->GetLoader("TPCLoader");
+      gime->LoadSDigits("read");
+      TTree * treear =  gime->TreeS();
+     
+      if (!treear) 
+       {
+        cerr<<"AliTPCDigitizer: Input tree with SDigits not found in"
+            <<" input "<< i1<<endl;
+        return;
+       }
+
+      sprintf(phname,"lhcphase%d",i1);
+      TParameter<float> *ph = (TParameter<float>*)treear->GetUserInfo()
+                              ->FindObject("lhcphase0");
+      if(!ph){
+        cerr<<"AliTPCDigitizer: LHC phase  not found in"
+            <<" input "<< i1<<endl;
+        return;
+      }
+      tree->GetUserInfo()->Add(new TParameter<float>(phname,ph->GetVal()));
+             //
+      if (treear->GetIndex()==0)  
+       treear->BuildIndex("fSegmentID","fSegmentID");
+      treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
+    }
+
+
+
+
+  //
+
+  param->SetZeroSup(2);
+
+  Int_t zerosup = param->GetZeroSup(); 
+  AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor(); 
+  AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise(); 
+  //
+  //Loop over segments of the TPC
+    
+  for (Int_t segmentID=0; segmentID<param->GetNRowsTotal(); segmentID++) 
+   {
+    Int_t sec, row;
+    if (!param->AdjustSectorRow(segmentID,sec,row)) 
+     {
+      cerr<<"AliTPC warning: invalid segment ID ! "<<segmentID<<endl;
+      continue;
+     }
+    AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec);  // pad gains per given sector
+    AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(sec);  // noise per given sector
+    digrow->SetID(segmentID);
+
+    Int_t nrows = 0;
+    Int_t ncols = 0;
+
+    Bool_t digitize = kFALSE;
+    for (Int_t i=0;i<nInputs; i++) 
+     { 
+
+      rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
+      gime = rl->GetLoader("TPCLoader");
+      
+      if (gime->TreeS()->GetEntryWithIndex(segmentID,segmentID) >= 0) {
+       digarr[i]->ExpandBuffer();
+       digarr[i]->ExpandTrackBuffer();
+        nrows = digarr[i]->GetNRows();
+        ncols = digarr[i]->GetNCols();
+       active[i] = kTRUE;
+       if (!fRegionOfInterest || (i == 0)) digitize = kTRUE;
+      } else {
+       active[i] = kFALSE;
+      }
+      if (fRegionOfInterest && !digitize) break;
+     }   
+    if (!digitize) continue;
+
+    digrow->Allocate(nrows,ncols);
+    digrow->AllocateTrack(3);
+
+    Float_t q=0;
+    Int_t label[1000]; //stack for 300 events 
+    Int_t labptr = 0;
+
+    Int_t nElems = nrows*ncols;     
+    for (Int_t i=0;i<nInputs; i++)
+     if (active[i]) { 
+       pdig[i] = digarr[i]->GetDigits();
+       ptr[i]  = digarr[i]->GetTracks();
+      }
+     
+    Short_t *pdig1= digrow->GetDigits();
+    Int_t   *ptr1= digrow->GetTracks() ;
+
+    
+
+    for (Int_t elem=0;elem<nElems; elem++)
+     {    
+
+       q=0;
+       labptr=0;
+       // looop over digits 
+        for (Int_t i=0;i<nInputs; i++) if (active[i]) 
+         { 
+          //          q  += digarr[i]->GetDigitFast(rows,col);
+            q  += *(pdig[i]);
+         
+           for (Int_t tr=0;tr<3;tr++) 
+            {
+             //             Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
+             Int_t lab = ptr[i][tr*nElems];
+             if ( (lab > 1) && *(pdig[i])>zerosup) 
+              {
+                label[labptr]=lab+masks[i];
+                labptr++;
+              }          
+            }
+           pdig[i]++;
+           ptr[i]++;
+         }
+        q/=16.;  //conversion factor
+       Float_t gain = gainROC->GetValue(row,elem/nrows);  // get gain for given - pad-row pad
+       //if (gain<0.5){
+         //printf("problem\n");
+       //}
+       q*= gain;
+       Float_t noisePad = noiseROC->GetValue(row,elem/nrows);
+        //       Float_t noise  = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());  
+        Float_t noise  = pTPC->GetNoise();
+        q+=noise*noisePad;
+        q=TMath::Nint(q);
+        if (q > zerosup)
+         { 
+          if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
+          //digrow->SetDigitFast((Short_t)q,rows,col);  
+          *pdig1 =Short_t(q);
+          for (Int_t tr=0;tr<3;tr++)
+           {
+            if (tr<labptr) 
+             ptr1[tr*nElems] = label[tr];
+           }
+          }
+        pdig1++;
+        ptr1++;
+     }
+    //
+    //  glitch filter
+    //
+    digrow->GlitchFilter();
+    //
+    digrow->CompresBuffer(1,zerosup);
+    digrow->CompresTrackBuffer(1);
+    tree->Fill();
+    if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";  
+   } //for (Int_t n=0; n<param->GetNRowsTotal(); n++) 
+  
+
+  orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
+  ogime = orl->GetLoader("TPCLoader");
+  ogime->WriteDigits("OVERWRITE");
+  
+  //fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
+  
+  delete digrow;     
+  for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
+  delete []masks;
+  delete []pdig;
+  delete []ptr;
+  delete []active;
+  delete []digarr;  
+}
+
+
+
+//------------------------------------------------------------------------
+void AliTPCDigitizer::ExecSave(Option_t* option)
+{
+  
+  // merge input tree's with summable digits
+  //output stored in TreeTPCD
+
+  TString optionString = option;
+  if (!strcmp(optionString.Data(),"deb")) {
+    cout<<"AliTPCDigitizer::Exec: called with option deb "<<endl;
+    fDebug = 3;
+  }
+  //get detector and geometry 
+  AliRunLoader *rl, *orl;
+  AliLoader *gime, *ogime;
+
+  
+  orl = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
+  ogime = orl->GetLoader("TPCLoader");
+  
+  rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
+  gime = rl->GetLoader("TPCLoader");
+  
+  rl->LoadgAlice();
+  AliRun* alirun = rl->GetAliRun();
+  
+  AliTPC *pTPC  = (AliTPC *) alirun->GetModule("TPC");
+  AliTPCParam * param = pTPC->GetParam();
   pTPC->GenerNoise(500000); //create teble with noise
   printf("noise %f \n",  param->GetNoise()*param->GetNoiseNormFac());
   //
@@ -88,52 +388,70 @@ void AliTPCDigitizer::Exec(Option_t* option)
   Int_t * masks = new Int_t[nInputs];
   for (Int_t i=0; i<nInputs;i++)
     masks[i]= fManager->GetMask(i);
-  Short_t **pdig= new Short_t*[nInputs];   //pointers to the expanded digits array
-  Int_t **ptr=  new Int_t*[nInputs];       //pointers to teh expanded tracks array
 
-  //create digits array for given sectors
   AliSimDigits ** digarr = new AliSimDigits*[nInputs]; 
-  for (Int_t i1=0;i1<nInputs; i1++){
+  for (Int_t i1=0;i1<nInputs; i1++)
+   {
     digarr[i1]=0;
     //    intree[i1]
-    TTree * treear =  fManager->GetInputTreeTPCS(i1);
+    rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i1));
+    gime = rl->GetLoader("TPCLoader");
+
+    TTree * treear =  gime->TreeS();
+    TBranch * br = treear->GetBranch("fSegmentID");
+    if (br) br->GetFile()->cd();
     if (!treear) {      
       cerr<<" TPC -  not existing input = \n"<<i1<<" ";      
-    }
+    } 
     treear->GetBranch("Segment")->SetAddress(&digarr[i1]);
   }
-  Stat_t nentries = fManager->GetInputTreeTPCS(0)->GetEntries();
+  
+  rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
+  gime = rl->GetLoader("TPCLoader");
+  Stat_t nentries = gime->TreeS()->GetEntries();
   
 
   //create branch's in TPC treeD
   AliSimDigits * digrow = new AliSimDigits;
-  TTree * tree  = fManager->GetTreeDTPC();
-  //if (tree->GetBranch("Segment") ) tree->GetBranch("Segment")->SetAddress(&digrow);
-  //else
-  tree->Branch("Segment","AliSimDigits",&digrow);
-  //
+  TTree * tree  = ogime->TreeD();
 
+  tree->Branch("Segment","AliSimDigits",&digrow);
   param->SetZeroSup(2);
 
   Int_t zerosup = param->GetZeroSup();
   //Loop over segments of the TPC
     
+  AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor();
+  AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise();
   for (Int_t n=0; n<nentries; n++) {
-  //    for (Int_t n=0; n<300; n++) {
-           
-    for (Int_t i=0;i<nInputs; i++){ 
-      fManager->GetInputTreeTPCS(i)->GetEvent(n);      
+    rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(0));
+    gime = rl->GetLoader("TPCLoader");
+    gime->TreeS()->GetEvent(n);
+
+    digarr[0]->ExpandBuffer();
+    digarr[0]->ExpandTrackBuffer();
+
+
+    for (Int_t i=1;i<nInputs; i++){ 
+//      fManager->GetInputTreeTPCS(i)->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());      
+      rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
+      gime = rl->GetLoader("TPCLoader");
+      gime->TreeS()->GetEntryWithIndex(digarr[0]->GetID(),digarr[0]->GetID());  
       digarr[i]->ExpandBuffer();
       digarr[i]->ExpandTrackBuffer();
+      if ((digarr[0]->GetID()-digarr[i]->GetID())>0) 
+       printf("problem\n");
+    
     }   
-    if ((digarr[0]->GetID()-digarr[1]->GetID())>0) 
-      printf("problem\n");
+    
     Int_t sec, row;
     if (!param->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
       cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
       continue;
     }
 
+    AliTPCCalROC * gainROC = gainTPC->GetCalROC(sec);  // pad gains per given sector
+    AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(sec);  // noise per given sector
     digrow->SetID(digarr[0]->GetID());
 
     Int_t nrows = digarr[0]->GetNRows();
@@ -146,70 +464,77 @@ void AliTPCDigitizer::Exec(Option_t* option)
     Int_t labptr = 0;
 
     
-    for (Int_t i=0;i<nInputs; i++){ 
-      pdig[i] = digarr[i]->GetDigits();
-      ptr[i]  = digarr[i]->GetTracks();
-    }
-    Short_t *pdig1= digrow->GetDigits();
-    Int_t   *ptr1= digrow->GetTracks() ;
 
-    Int_t nElems = nrows*ncols;      
-
-    //    for (Int_t rows=0;rows<nrows; rows++){
-    //  for (Int_t col=0;col<ncols; col++){
-    for (Int_t elem=0;elem<nElems; elem++){
-
-       q=0;
-       labptr=0;
-       // looop over digits 
+    for (Int_t rows=0;rows<nrows; rows++){
+      for (Int_t col=0;col<ncols; col++){
+    
+       q=0;
+       labptr=0;
+       // looop over digits 
         for (Int_t i=0;i<nInputs; i++){ 
-         //          q  += digarr[i]->GetDigitFast(rows,col);
-          q  += *(pdig[i]);
-         
+         q  += digarr[i]->GetDigitFast(rows,col);
+          //q  += *(pdig[i]);
+         
           for (Int_t tr=0;tr<3;tr++) {
-           //             Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
-           Int_t lab = ptr[i][tr*nElems];
-            if ( lab > 1) {
+           Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
+           //Int_t lab = ptr[i][tr*nElems];
+            if ( (lab > 1) ) {
               label[labptr]=lab+masks[i];
               labptr++;
-            }     
+            }          
           }
-         pdig[i]++;
-         ptr[i]++;
-         
+         // pdig[i]++;
+         //ptr[i]++;
+         
         }
-       q/=16.;  //conversion factor
-       //      Float_t noise  = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());  
-       Float_t noise  = pTPC->GetNoise();
-       q+=noise;
+       q/=16.;  //conversion factor
+       //       Float_t noise  = gRandom->Gaus(0,param->GetNoise()*param->GetNoiseNormFac());  
+       Float_t gain = gainROC->GetValue(row,col);
+       q*= gain;
+       Float_t noisePad = noiseROC->GetValue(row, col);
+
+       Float_t noise  = pTPC->GetNoise();
+       q+=noise*noisePad;
+
         q=TMath::Nint(q);
         if (q > zerosup){ 
-         
-         if(q > param->GetADCSat()) q = (Short_t)(param->GetADCSat());
-         //digrow->SetDigitFast((Short_t)q,rows,col);  
-         *pdig1 =Short_t(q);
-         for (Int_t tr=0;tr<3;tr++){
-           if (tr<labptr) 
-             //  ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
-             ptr1[tr*nElems] = label[tr];
-           //else
-             //            ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);          
-           //  ptr1[tr*nElems] = 1;
-         }
-       }
-       pdig1++;
-       ptr1++;
+         
+         if(q >= param->GetADCSat()) q = (Short_t)(param->GetADCSat() - 1);
+         digrow->SetDigitFast((Short_t)q,rows,col);  
+         // *pdig1 =Short_t(q);
+         for (Int_t tr=0;tr<3;tr++){
+           if (tr<labptr) 
+             ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
+           //ptr1[tr*nElems] = label[tr];
+           //else
+             //           ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);          
+           //  ptr1[tr*nElems] = 1;
+         }
+       }
+       //pdig1++;
+       //ptr1++;
+    }
     }
     
     digrow->CompresBuffer(1,zerosup);
     digrow->CompresTrackBuffer(1);
     tree->Fill();
     if (fDebug>0) cerr<<sec<<"\t"<<row<<"\n";  
-  }
-  fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
+  } 
+//  printf("end TPC merging - end -Tree %s\t%p\n",fManager->GetInputTreeH(0)->GetName(),fManager->GetInputTreeH(0)->GetListOfBranches()->At(3));
+  //fManager->GetTreeDTPC()->Write(0,TObject::kOverwrite);
+    ogime->WriteDigits("OVERWRITE");
+
+    for (Int_t i=1;i<nInputs; i++) 
+     { 
+      rl = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(i));
+      gime = rl->GetLoader("TPCLoader");
+      gime->UnloadSDigits();
+     }
+    ogime->UnloadDigits();
+    
   delete digrow;     
   for (Int_t i1=0;i1<nInputs; i1++) delete digarr[i1];
-  delete []masks;
-  delete digarr;  
+  delete [] masks;
+  delete [] digarr;  
 }
-