AliRoot Core  v5-06-30 (35d6c57)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliTPCAltro.C
Go to the documentation of this file.
1 #if !defined(__CINT__)
2 #include <Riostream.h>
3 #include <TFile.h>
4 #include <TTree.h>
5 #include "AliTPCParamSR.h"
6 #include "AliTPCDigitsArray.h"
7 #include "AliSimDigits.h"
8 #include "AliTPCBuffer160.h"
9 #endif
10 
12 
13 int AliTPCAltro(Int_t eth=0){
14 
17 
18  Int_t offset=1;
19  /*
20  NB the amplitude values stored in the ALTRO file are shifted by offset
21  because the range for each word goes from 0 to 1023, now due to zero suppression
22  values lower that the threshold never appear.
23  */
24 
25  Int_t PSecNumber=-1;
26  Int_t PRowNumber=-1;
27  Int_t PPadNumber=-1;
28  Int_t PTimeBin=-1;
29  Int_t BunchLength=0;
30  //AliTPCBuffer160 is used in write mode to generate AltroFormat.dat file
31  AliTPCBuffer160 Buffer("AltroFormat.dat",1);
32  ULong_t Count=0;
33  Int_t nwords=0;
34  Int_t numPackets=0;
35 
36  const char * inFile_new = "galice.root";
37  AliRunLoader *rl = AliRunLoader::Open(inFile_new,"Event","read");
38 
39  Int_t nevents=rl->GetNumberOfEvents();
40  cout<<"Number of Events:"<<nevents<<endl;
41  Int_t choice=0;
42  do{
43  cout<<"Insert the event number: ";
44  cin>>choice;
45  }while (choice<=0 || choice>nevents);
46  rl->GetEvent(choice-1);
47  AliLoader *tpcloader=rl->GetLoader("TPCLoader");
48  tpcloader->LoadDigits();
49  TTree *digitsTree=tpcloader->TreeD();
50 
51  AliSimDigits digrows, *dummy=&digrows;
52  digitsTree->GetBranch("Segment")->SetAddress(&dummy);
53  Stat_t nrows = digitsTree->GetEntries();
54  cout<<"Number of entries (rows):"<<nrows<<endl;
55  // get the TPC parameters
56  rl->CdGAFile();
57  AliTPCParamSR* param = AliTPC::LoadTPCParam(gFile);
58  if (!param)
59  cout<<"No TPC parameter"<<endl;
61  digarr->Setup(param);
62  digarr->ConnectTree(digitsTree);
63  //ofstream ftxt("Data.txt");
64  for (Int_t n=0; n<nrows; n++) {
65  Int_t sec,row;
66  AliSimDigits *digrow=(AliSimDigits*)digarr->LoadEntry(n);
67  param->AdjustSectorRow(digrow->GetID(),sec,row);
68 
69  //cout<<"Sector:"<<sec<<" Row:"<<row<<endl;
70  digrow->First();
71  do{
72  Short_t dig=digrow->CurrentDigit(); //adc
73  Int_t time=digrow->CurrentRow(); //time
74  Int_t pad =digrow->CurrentColumn(); // pad
75  //cout<<"dig:"<<dig<<" time:"<<time<<" pad:"<<pad<<endl;
76  if(dig>eth){
77  Count++;
78  //ftxt<<"Sec: "<<sec<<" Row: "<<row<<" Pad:"<<pad<<" Time: "<<time<<" ADC:"<<dig<<endl;
79  //cout<<"Sec: "<<sec<<" Row: "<<row<<" Pad:"<<pad<<" Time: "<<time<<" ADC:"<<dig<<endl;
80  if (PPadNumber==-1){
81  PSecNumber=sec;
82  PRowNumber=row;
83  PPadNumber=pad;
84  // PAmplitude=dig;
85  PTimeBin=time;
86  BunchLength=1;
87  Buffer.FillBuffer(dig-offset);
88  nwords++;
89  }//end if
90  else{
91  if ( (time==(PTimeBin+1)) &&
92  (PPadNumber==pad) &&
93  (PRowNumber==row) &&
94  (PSecNumber==sec)){
95  BunchLength++;
96  }//end if
97  else{
98  Buffer.FillBuffer(PTimeBin);
99  Buffer.FillBuffer(BunchLength+2);
100  nwords+=2;
101  if ((PPadNumber!=pad)||(PRowNumber!=row)||(PSecNumber!=sec)){
102  //Trailer is formatted and inserted!!
103  Buffer.WriteTrailer(nwords,PPadNumber,PRowNumber,PSecNumber);
104  numPackets++;
105  nwords=0;
106  }//end if
107 
108  BunchLength=1;
109  PPadNumber=pad;
110  PRowNumber=row;
111  PSecNumber=sec;
112  }//end else
113  PTimeBin=time;
114  Buffer.FillBuffer(dig-offset);
115  nwords++;
116  }//end else
117  }//end if
118  } while (digrow->Next());
119  }//end for
120 
121  Buffer.FillBuffer(PTimeBin);
122  Buffer.FillBuffer(BunchLength+2);
123  nwords+=2;
124  Buffer.WriteTrailer(nwords,PPadNumber,PRowNumber,PSecNumber);
125 
126  numPackets++;
127  cout<<"There are "<<Count<<" Digits\n";
128  cout<<"Packets "<<numPackets<<"\n";
129  //ftxt.close();
130  return 0;
131 }//end macro
Int_t CurrentDigit()
Definition: AliDigits.h:42
Manager and of geomety classes for set: TPC.
Definition: AliTPCParamSR.h:15
Int_t CurrentRow()
Definition: AliDigits.h:40
Time Projection Chamber clusters objects.
virtual Bool_t Next()
Definition: AliDigits.cxx:280
virtual Bool_t First()
Definition: AliDigits.cxx:271
int AliTPCAltro(Int_t eth=0)
Definition: AliTPCAltro.C:13
virtual Bool_t ConnectTree(const char *treeName)
Int_t GetID()
Definition: AliSegmentID.h:19
virtual AliSegmentID * LoadEntry(Int_t index)
Bool_t AdjustSectorRow(Int_t index, Int_t &sector, Int_t &row) const
Bool_t Setup(AliDetectorParam *param)
Int_t CurrentColumn()
Definition: AliDigits.h:41