Using TMath::Abs instead of fabs
[u/mrichter/AliRoot.git] / HLT / misc / AliL3TransBit.cxx
1 //$Id$
2
3 // Original author: Marian Ivanov, GSI Darmstadt for AliROOT
4 // migrated to L3 by Anders Vestbo <mailto:vestbo@fi.uib.no>
5 //*-- Copyright & copy ASV
6
7 #include "AliL3StandardIncludes.h"
8
9 #include "AliL3TransBit.h"
10
11 #if GCCVERSION == 3
12 using namespace std;
13 #endif
14
15 //_____________________________________________________________
16 // AliL3Transbit (taken from the offline AliROOT code)
17 //
18 // Time Projection Chamber ADC bit compresion lookup table
19 //
20 //  Conversion equation:
21 //    For AliTransBit_v1
22 //    dy/dx= Int_t(1+x/fX0)
23 //
24 //    For AliTransBit_v2
25 //    y  =(2**bit0) ln(1+x/fX0) / ln(1+(2**bit1-1)/fX0)
26 //
27 //    where x0 is calculated in function GetOptimumX0()
28 //
29 //  Example session in aliroot:
30 //    Int_t b0=10;  // original number of bits
31 //    Int_t b1=8;   // compressed
32 //
33 //    AliTransBit_v1 trans;
34 //    Int_t x0=TMath::Nint(TMath::Exp(b0*TMath::Log(2)));
35 //    Int_t x1=TMath::Nint(TMath::Exp(b1*TMath::Log(2)));
36 //    trans.SetBits(b0,b1);
37 //    trans.FindOptimumX0();
38 //    trans.Update();
39 //    cout<<trans.Get0to1(x0-2)<<"\n";
40 //    cout<<trans.Get1to0(x1-2)<<"\n";
41 //
42 //    // to produce table
43 //    for( Int_t i=0;i<x1;i++) cout<<i<<"\t"<<trans.Get1to0(i)<<"\n"; > table1to0.txt
44 //    for( Int_t i=0;i<x0;i++) cout<<i<<"\t"<<trans.Get0to1(i)<<"\n"; > table0to1.txt
45 //
46 //    for( Int_t i=0;i<x1-1;i++) cout<<i<<"\t"<<trans.Get1to0(i+1)-trans.Get1to0(i)<<"\n"; > tabled.txt
47 //
48
49 ClassImp(AliL3TransBit)
50 ClassImp(AliL3TransBit_v1)
51 ClassImp(AliL3TransBit_v2)
52
53 AliL3TransBit::AliL3TransBit()
54 {
55   fTable0 = 0;
56   fTable1 = 0;
57   fBit0   = 10;
58   fBit1   = 8;
59   fX0     = 0;
60 }
61
62 AliL3TransBit::~AliL3TransBit() 
63 {
64   //default destructor
65   if (fTable0!=0) delete [] fTable0;
66   if (fTable1!=0) delete [] fTable1;
67 }
68
69 Double_t AliL3TransBit_v1::FindOptimumX0()
70 {
71   //find x0 for which derivation at xder1 is equal 1
72   Int_t x0=(Int_t)rint(exp(fBit0*log(2.)));
73   Int_t x1=(Int_t)rint(exp(fBit1*log(2.)));
74
75   fX0 = ((x1-2)*(x1-2)/2.)/(x0-x1-1);  //starting fX0
76   Int_t digit=0;
77   Int_t j;
78   Double_t weight=1;
79   for(j=0;(j<50)&&(digit!=(x0-1));j++){  
80     Int_t olddigit=digit;
81     digit=0;
82     for (Int_t  i=0; i<x1-1;i++) digit+= Int_t(1+Double_t(i)/fX0);      
83     fX0*= (1.+weight*Double_t(digit)/Double_t(x0-1))/(1.+weight);    
84     if ( ((olddigit-(x0-1)) * (digit-(x0-1))) <0 ) 
85       weight*=0.25;
86     // cout<<j<<"\t"<<fX0<<"\t"<<digit<<"\n";
87   }
88   //  cout<<"End of iteration "<<j<<"\n";
89   //cout<<"digit..."<<digit<<"\n";  
90   return fX0;
91 }
92
93 void AliL3TransBit_v1::Update()
94 {
95   //construct lookup tables for loosy compression from 
96   if (fX0<1) fX0 = FindOptimumX0();  
97   Int_t x0=(Int_t)rint(exp(fBit0*log(2.)));
98   Int_t x1=(Int_t)rint(exp(fBit1*log(2.)));
99   
100   //fTable0 - conversion from bit0 coding to bit1 coding
101   if (fTable0!=0) delete fTable0;
102   fTable0= new Int_t[x0];
103   //fTable1 - conversion table from bit1 to bit0 coding
104   if (fTable1!=0) delete [] fTable1;
105   fTable1= new Int_t[x1];
106
107   Int_t digit=0;
108   for (Int_t  i=0; i<x1-1;i++){    
109     Int_t ddig=Int_t(1+Double_t(i)/fX0);
110     for (Int_t j=0;j<ddig;j++) if ((digit+j)<x0) fTable0[digit+j] = i;    
111     fTable1[i]= digit+ddig/2;
112     digit+= ddig;
113   }
114
115   for (Int_t i=digit;i<x0-1;i++) fTable0[i]=x1-2;
116   fTable1[x1-1]=x0-1; //overflow 
117   fTable0[x0-1]=x1-1; //overflow    
118   return;
119 }
120
121 Double_t AliL3TransBit_v2::FindOptimumX0()
122 {
123   //find x0 for which derivation at xder1 is equal 1
124   const Float_t xder1=1;
125   const Float_t dx=0.1;
126
127   Float_t x0=exp(fBit0*log(2.));
128   Float_t x1=exp(fBit1*log(2.));
129   Float_t deriv = 0;
130   Float_t x;
131   for (x=x1;( (x>1)&&(deriv<1)) ;x-=dx)
132     {
133       deriv = (x1-1)/( log(1.+x0/x) *x *(1+xder1/x));
134     }
135   x+=dx/2.;
136   fX0 = x;
137   return fX0;
138 }
139
140 void AliL3TransBit_v2::Update()
141 {
142   //construct lookup tables for loosy compresion from 
143   if (fX0<1) fX0 = FindOptimumX0();  
144   //Float_t x0=(Int_t)rint(exp(fBit0*log(2.)));
145   //Float_t x1=(Int_t)rint(exp(fBit1*log(2.)));
146   Int_t x0=(Int_t)rint(exp(fBit0*log(2.)));
147   Int_t x1=(Int_t)rint(exp(fBit1*log(2.)));
148   
149   //fTable0 - conversion from bit0 coding to bit1 coding
150   if (fTable0!=0) delete fTable0;
151   fTable0= new Int_t[x0];
152
153   //fTable1 - conversion table from bit1 to bit0 coding
154   if (fTable1!=0) delete [] fTable1;
155   fTable1= new Int_t[x1];
156
157   Int_t i;
158
159   for ( i=0; i<x0;i++)
160       fTable0[i] =(Int_t)rint((x1-0.501)*log(1.+Float_t(i)/fX0)/
161                                  log(1.+(x0-1)/fX0));
162
163   Int_t old0=-1;
164   Int_t old1=-1;
165   Int_t new1;
166   for ( i=0; i<x0;i++){
167       new1 = fTable0[i];
168       if (new1!=old1){
169         if (old1>=0) fTable1[old1] =(old0+i)/2;
170         old0 = i;
171         old1 = new1;
172       }
173   }
174   fTable1[old1]=(Int_t)rint((Float_t)(old0+x0)/2);
175   
176   return;
177 }