7bd5e0fd56ed72b71afda8ce7b7a229d277731c8
[u/mrichter/AliRoot.git] / TPC / AliTPCCompression.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2003, 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 #include <TObjArray.h>
17 #include "Riostream.h"
18 #include <math.h>
19 #include "AliTPCCompression.h"
20 #include "AliTPCBuffer160.h"
21 #include "AliTPCHuffman.h"
22
23 ClassImp(AliTPCCompression)
24 //////////////////////////////////////////////////////////////////////////////////////////////////
25 AliTPCCompression::AliTPCCompression(){
26   fDimBuffer=sizeof(ULong_t)*8;
27   fFreeBitsBuffer=fDimBuffer;
28   fReadBits=0;
29   fPos=0;
30   fBuffer=0;
31   fVerbose=0;
32   fFillWords=0;
33   return;
34 }
35 //////////////////////////////////////////////////////////////////////////////////////////////////
36 AliTPCCompression::AliTPCCompression(const AliTPCCompression &source){
37   this->fDimBuffer=source.fDimBuffer;
38   this->fFreeBitsBuffer=source.fFreeBitsBuffer;
39   this->fReadBits=source.fReadBits;
40   this->fPos=source.fPos;
41   this->fBuffer=source.fBuffer;
42   this->fVerbose=source.fVerbose;
43   this->fFillWords=source.fFillWords;
44   return;
45 }
46 //////////////////////////////////////////////////////////////////////////////////////////////////
47 AliTPCCompression&  AliTPCCompression::operator=(const AliTPCCompression &source){
48   this->fDimBuffer=source.fDimBuffer;
49   this->fFreeBitsBuffer=source.fFreeBitsBuffer;
50   this->fReadBits=source.fReadBits;
51   this->fPos=source.fPos;
52   this->fBuffer=source.fBuffer;
53   this->fVerbose=source.fVerbose;
54   this->fFillWords=source.fFillWords;
55   return *this;
56
57 //////////////////////////////////////////////////////////////////////////////////////////////////
58 void AliTPCCompression::NextTable(Int_t Val,Int_t &NextTableType,Int_t &BunchLen,Int_t &Count){
59   /*
60     Table index:
61     0==> Bunch length value     
62     1==> Time Bin value 
63     2==> 1-samples bunch
64     3==> Central samples
65     4==> Border samples
66   */  
67   switch (NextTableType){
68   case 0:{
69     BunchLen=Val-2;
70     NextTableType=1;
71     break;
72   }//end case 0
73   case 1:{
74     if (BunchLen==1)NextTableType=2;
75     else{
76       NextTableType=4;
77       Count=1;
78     }
79     break;
80   }//end case 1
81   case 2:{
82     NextTableType=0;
83     break;
84   }//end case 2
85   case 3:{
86     Count++;
87     if (Count==(BunchLen-1)){
88       NextTableType=4;
89     }
90     break;
91   }//end case 3
92   case 4:{
93     if (Count==1){
94       if (BunchLen>2)
95         NextTableType=3;
96       else
97         Count++;
98     }
99     else
100       NextTableType=0;
101     break;
102   }//end case 4
103   }//end switch
104   return;
105 }
106
107 /////////////////////////////////////////////////////////////////////////////////////////////////////
108
109 /////////////////////////////////////////////////////////////////////////////////////////////////////
110
111 Int_t AliTPCCompression::FillTables(const char* fSource,AliTPCHTable* table[],const Int_t NumTables){
112   //This method is used to compute the frequencies of the symbols in the source file
113   AliTPCBuffer160 Buff(fSource,0);
114   ULong_t CountWords=0;
115   ULong_t CountTrailer=0;
116   Int_t NumWords,PadNum,RowNum,SecNum=0;
117   Int_t Value=0;
118   ULong_t Stat[5]={0,0,0,0,0};
119   Int_t EndFill=0;
120   Int_t End=1;
121   while(Buff.ReadTrailerBackward(NumWords,PadNum,RowNum,SecNum) !=-1 ){
122     if(End){
123       EndFill=Buff.GetFillWordsNum();
124       End=0;
125     }//endif
126     CountTrailer++;
127     if (NumWords%4){
128       fFillWords+=4-NumWords%4;
129       for(Int_t j=0;j<(4-NumWords%4);j++){
130         Value=Buff.GetNextBackWord();
131       }//end for
132     }//end if
133     
134     Int_t Packet[1024];
135     Int_t TimePos[345];
136     Int_t Tp=0;
137     for(Int_t i=0;i<345;i++)TimePos[i]=0;
138     for(Int_t i=0;i<1024;i++)Packet[i]=0;
139     
140     Int_t NextTableType=0;
141     Int_t BunchLen=0;
142     Int_t Count=0;
143     for(Int_t i=0;i<NumWords;i++){
144       Value=Buff.GetNextBackWord();
145       Packet[i]=Value;
146       if(NextTableType==1){
147         TimePos[Tp]=i;
148         Tp++;
149       }
150       NextTable(Value,NextTableType,BunchLen,Count);
151     }//end for
152     //computing the Time gap between two bunches
153     Int_t temp=0;
154     Tp--;
155     Int_t PreviousTime=Packet[TimePos[Tp]];
156     for(Int_t i=Tp-1;i>=0;i--){
157       Int_t TimPos=TimePos[i];
158       Int_t BunchLen=Packet[TimPos-1]-2;
159       temp=Packet[TimPos];
160       Packet[TimPos]=Packet[TimPos]-PreviousTime-BunchLen;
161       PreviousTime=temp;
162     }
163     NextTableType=0;
164     Count=0;
165     BunchLen=0;
166     for(Int_t i=0;i<NumWords;i++){
167       Value=Packet[i];
168       table[NextTableType]->SetFrequency(Value);
169       Stat[NextTableType]++;
170       NextTable(Value,NextTableType,BunchLen,Count);
171       CountWords++;
172     }//end for
173   }//end while
174   cout<<"Number of words:       "<<CountWords<<endl;
175   cout<<"Number of trailers:    "<<CountTrailer<<endl;
176   cout<<"Number of fill words   "<<fFillWords+EndFill<<endl;
177   cout<<"Total number of words: "<<CountWords+CountTrailer*4+fFillWords<<endl;
178   //STATISTICS  
179   stat.open("Statistics");
180   stat<<"Number of words:..........................................."<<CountWords<<endl;
181   stat<<"Number of trailers (4 10 bits words in each one)..........."<<CountTrailer<<endl;
182   stat<<"Number of fill words:......................................"<<fFillWords+EndFill<<endl;
183   stat<<"Total number of words:....................................."<<CountWords+CountTrailer*4+fFillWords+EndFill<<endl;
184   stat<<"-----------------------------------------"<<endl;
185   stat<<"Number of Bunches............."<<Stat[0]<<endl;
186   stat<<"Number of Time bin............"<<Stat[1]<<endl;
187   stat<<"Number of One Samples Bunch..."<<Stat[2]<<endl;
188   stat<<"Number of Central Samples....."<<Stat[3]<<endl;
189   stat<<"Number of Border Samples......"<<Stat[4]<<endl;
190   stat<<"-----------------------------------------"<<endl;
191   ULong_t FileDimension=(ULong_t)ceil(double((CountTrailer*4+CountWords+fFillWords+EndFill)*10/8));
192   stat<<"Total file Size in bytes.."<<FileDimension<<endl;
193   Double_t Percentage=ceil((fFillWords+EndFill)*125)/FileDimension;
194   stat<<"Fill Words................"<<(ULong_t)ceil((fFillWords+EndFill)*10/8)<<" bytes   "<<Percentage<<"%"<<endl;  
195   Percentage=(Double_t)CountTrailer*500/FileDimension;
196   stat<<"Trailer..................."<<CountTrailer*5<<" bytes   "<<Percentage<<"%"<<endl;
197
198   Percentage=(Double_t)((Stat[0]+Stat[1]+Stat[2]+Stat[3]+Stat[4])) *125/FileDimension;
199   stat<<"Data......................"<<(ULong_t)ceil((Stat[0]+Stat[1]+Stat[2]+Stat[3]+Stat[4])*10/8)<<" bytes   "<<Percentage<<"%"<<endl;
200
201   Percentage=(Double_t)(Stat[0]*125)/FileDimension;
202   stat<<"Bunch....................."<<(ULong_t)ceil(Stat[0]*10/8)<<" bytes  "<<Percentage<<"%"<<endl;  //  
203   Percentage=(Double_t)(Stat[1]*125)/FileDimension;
204   stat<<"Time......................"<<(ULong_t)ceil(Stat[1]*10/8)<<" bytes  "<<Percentage<<"%"<<endl;  //  
205
206
207   Percentage=(Double_t)((Stat[2]+Stat[3]+Stat[4])) *125/FileDimension;
208   stat<<"Amplitude values.........."<<(ULong_t)ceil((Stat[2]+Stat[3]+Stat[4])*10/8)<<" bytes  "<<Percentage<<"%"<<endl;
209   Percentage=(Double_t)(Stat[2]*125)/FileDimension;
210   stat<<"     One Samples..............."<<(ULong_t)ceil(Stat[2]*10/8)<<" bytes  "<<Percentage<<"%"<<endl;  //  
211   Percentage=(Double_t)(Stat[3]*125)/FileDimension;
212   stat<<"     Central Samples..........."<<(ULong_t)ceil(Stat[3]*10/8)<<" bytes  "<<Percentage<<"%"<<endl;  //  
213   Percentage=(Double_t)(Stat[4]*125)/FileDimension;
214   stat<<"     Border Samples............"<<(ULong_t)ceil(Stat[4]*10/8)<<" bytes  "<<Percentage<<"%"<<endl;  //  
215   stat.close();
216   return 0;
217 }
218 ////////////////////////////////////////////////////////////////////////////////////////
219 Int_t AliTPCCompression::StoreTables(AliTPCHTable* table[],const Int_t NumTable){
220   char filename[15];
221   ofstream fTable;
222   for(Int_t k=0;k<NumTable;k++){
223     sprintf(filename,"Table%d.dat",k); 
224     fTable.open(filename,ios::binary);
225     Int_t dim=table[k]->Size();
226     //Table dimension is written into a file
227     fTable.write((char*)(&dim),sizeof(Int_t));
228     //One table is written into a file
229     for(Int_t i=0;i<dim;i++){
230       UChar_t CodeLen=table[k]->CodeLen()[i];
231       //      ULong_t Code=(ULong_t)table[k]->Code()[i];
232       Double_t Code=table[k]->Code()[i];
233       fTable.write((char*)(&CodeLen),sizeof(UChar_t));
234       //fTable.write((char*)(&Code),sizeof(ULong_t));
235       fTable.write((char*)(&Code),sizeof(Double_t));
236     } //end for
237     fTable.close();
238   }//end for
239   return 0;
240 }
241 ////////////////////////////////////////////////////////////////////////////////////////
242 Int_t AliTPCCompression::CreateTables(const char* fSource,const Int_t NumTables){
243     Int_t n=10;// 10 bits per symbol 
244   /*
245     Table index:
246     0==> Bunch length values     
247     1==> Time Bin values 
248     2==> 1-samples bunch
249     3==> Central samples
250     4==> Border samples
251   */
252   AliTPCHTable ** table = new AliTPCHTable*[NumTables];
253   //The table is inizialized with the rigth number of rows 
254   for(Int_t i=0;i<NumTables;i++){table[i]=new  AliTPCHTable((Int_t)(pow(2,n)));}
255   //The frequencies are calculated and the tables are filled
256   if (fVerbose)
257     cout<<"Filling tables...\n";
258   //The get the frequencies 
259   FillTables(fSource,table,NumTables);
260
261   //This part will be used in the table optimization phase
262   /*
263   for(Int_t i=0;i<NumTables;i++){
264     table[i]->CompleteTable(i);
265   }
266   */
267   if(fVerbose){
268     cout<<"Entropy of Bunch length table........."<<table[0]->GetEntropy()<<endl;
269     cout<<"Entropy of Time bin table............."<<table[1]->GetEntropy()<<endl;
270     cout<<"Entropy of one Sample bunch table....."<<table[2]->GetEntropy()<<endl;
271     cout<<"Entropy of Central Sample table......."<<table[3]->GetEntropy()<<endl;
272     cout<<"Entropy Border Samples table.........."<<table[4]->GetEntropy()<<endl;
273   }
274   stat.open("Statistics",ios::app);
275   stat<<endl;
276   stat<<"----------------- ENTROPY for castomized tables --------------------------"<<endl;
277   stat<<"Entropy of Bunch length table......."<<table[0]->GetEntropy()<<endl;
278   stat<<"Entropy of Time bin table..........."<<table[1]->GetEntropy()<<endl;
279   stat<<"Entropy of one Sample bunch table..."<<table[2]->GetEntropy()<<endl;
280   stat<<"Entropy of Central Sample table....."<<table[3]->GetEntropy()<<endl;
281   stat<<"Entropy Border Samples table........"<<table[4]->GetEntropy()<<endl;
282   stat.close();
283  
284   if (fVerbose)
285     cout<<"Tables filled \n";
286   //Tables are saved in a sequence of text file and using the macro Histo.C is it possible to get
287   //a series of histograms rappresenting the frequency distribution
288   table[0]->StoreFrequencies("BunchLenFreq.txt");
289   table[1]->StoreFrequencies("TimeFreq.txt");
290   table[2]->StoreFrequencies("Sample1Freq.txt");
291   table[3]->StoreFrequencies("SCentralFreq.txt");
292   table[4]->StoreFrequencies("SBorderFreq.txt");
293   if (fVerbose)
294     cout<<"Creating Tables..\n";
295   //One Huffman tree is created for each table starting from the frequencies of the symbols
296   for(Int_t i=0;i<NumTables;i++){
297     table[i]->BuildHTable();
298     if (fVerbose==2){
299       cout<<"Number of elements inside the table:"<<table[i]->GetWordsNumber();
300       switch(i){
301       case 0:{
302         cout<<" (Bunch Length)"<<endl;
303         break;
304       }
305       case 1:{
306         cout<<" (Time Bin)"<<endl;
307         break;
308       }
309       case 2:{
310         cout<<" (1 Samples Bunch)"<<endl;
311         break;
312       }
313       case 3:{
314         cout<<" (Central Samples)"<<endl;
315         break;
316       }
317       case 4:{
318         cout<<" (Border Samples)"<<endl;
319         break;
320       }
321       }//end switch
322       table[i]->PrintTable();
323     }
324   }
325   //The tables are saved ad binary files
326   StoreTables(table,NumTables);
327   //The tables stored in memory are deleted; 
328   for(Int_t i=0;i<NumTables;i++)delete table[i];
329   delete [] table;
330   return 0;
331 }
332 /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
333 Int_t AliTPCCompression::RetrieveTables(AliTPCHTable* table[],Int_t NumTable){
334   if (fVerbose)
335     cout<<"Retrieving tables from files \n";
336   //  ULong_t Code;
337   Double_t Code;
338   UChar_t CodeLen;
339   ifstream fTable;  
340   char filename[15];
341   //The following for loop is used to generate the Huffman trees acording to the tables
342   for(Int_t k=0;k<NumTable;k++){
343     Int_t Dim;//this variable contains the table dimension
344     sprintf(filename,"Table%d.dat",k); 
345     fTable.open(filename,ios::binary);
346     fTable.read((char*)(&Dim),sizeof(Int_t));
347     if (fVerbose)
348       cout<<"Table dimension: "<<Dim<<endl;
349     table[k]=new AliTPCHTable(Dim);
350     for(Int_t i=0;i<Dim;i++){
351       fTable.read((char*)(&CodeLen),sizeof(UChar_t));
352       table[k]->SetCodeLen(CodeLen,i);
353       //      fTable.read((char*)(&Code),sizeof(ULong_t));
354       fTable.read((char*)(&Code),sizeof(Double_t));
355       table[k]->SetCode(Mirror((ULong_t)Code,CodeLen),i);
356     }//end for 
357     fTable.close();
358   }//end for 
359   if (fVerbose)
360     cout<<"Trees generated \n";
361   //At this point the trees are been built
362   return 0;
363 }
364 ////////////////////////////////////////////////////////////////////////////////////////
365 /*                               COMPRESSION                                          */
366 ////////////////////////////////////////////////////////////////////////////////////////
367
368 void AliTPCCompression::StoreValue(ULong_t val,UChar_t len){
369   if (len<=fFreeBitsBuffer){           // val is not splitted in two buffer
370     fFreeBitsBuffer-=len;
371     fBuffer=fBuffer<<len;
372     fBuffer=fBuffer|val;    
373     if(!fFreeBitsBuffer){              // if the buffer is full it is written into a file 
374       f.write((char*)(&fBuffer),sizeof(ULong_t));       
375       fFreeBitsBuffer=fDimBuffer;
376       fBuffer=0;
377     }
378   }//end if
379   else{                               //val has to be splitted in two buffers
380     fBuffer=fBuffer<<fFreeBitsBuffer;
381     ULong_t temp;
382     temp=val;
383     temp=temp>>(len-fFreeBitsBuffer);
384     fBuffer=fBuffer|temp;
385     f.write((char*)(&fBuffer),sizeof(ULong_t));
386     fFreeBitsBuffer=fDimBuffer-(len-fFreeBitsBuffer);
387     val=val<<fFreeBitsBuffer;
388     val=val>>fFreeBitsBuffer;
389     fBuffer=val;
390   }//end else
391   return;
392 }
393 //////////////////////////////////////////////////////////////////////////////////////////////////
394 void AliTPCCompression::Flush(){
395   //The last buffen cannot be completely full
396   if(fFreeBitsBuffer<fDimBuffer){
397     fBuffer=fBuffer<<fFreeBitsBuffer;
398     f.write((char*)(&fBuffer),sizeof(ULong_t));  
399   }//end if
400   return;
401 }
402 //////////////////////////////////////////////////////////////////////////////////////////////////
403 ULong_t AliTPCCompression::Mirror(ULong_t val,UChar_t len){
404   ULong_t specular=0;
405   ULong_t Mask=0x1;
406   ULong_t bit;
407   for(Int_t i=0;i<len;i++){
408     bit=val&Mask;
409     bit=bit>>i;
410     specular=specular<<1;
411     specular=specular|bit;
412     Mask=Mask<<1;
413   }
414   return specular;
415 }
416 //////////////////////////////////////////////////////////////////////////////////////////////////
417 Int_t AliTPCCompression::CompressData(AliTPCHTable* table[],Int_t NumTable,const char* fSource,const char* fDest){
418   cout<<" COMPRESSION "<<endl;
419   cout<<"compression of the file "<<fSource<<" Output File: "<<fDest<<endl;
420   //the output file is open
421   f.open(fDest,ios::binary|ios::out);
422   //Tables are written into the output file
423   for(Int_t k=0;k<NumTable;k++){
424     Int_t dim=table[k]->Size();
425     //Table dimension is written into a file
426     f.write((char*)(&dim),sizeof(Int_t));
427     //One table is written into a file
428     for(Int_t i=0;i<dim;i++){
429       UChar_t CodeLen=table[k]->CodeLen()[i];
430       ULong_t Code=(ULong_t)table[k]->Code()[i];
431       f.write((char*)(&CodeLen),sizeof(UChar_t));
432       f.write((char*)(&Code),sizeof(ULong_t));
433     } //end for
434   }//end for
435
436   // Source file is open
437   AliTPCBuffer160 Buff(fSource,0);
438   //coded words are written into the output file
439   Int_t NumWords,PadNum,RowNum,SecNum=0;
440   ULong_t StoredWords=0;
441   Int_t Value=0;
442   ULong_t NumPacket=0;
443   while(Buff.ReadTrailerBackward(NumWords,PadNum,RowNum,SecNum) !=-1 ){
444     NumPacket++;
445     if (NumWords%4){
446       for(Int_t j=0;j<(4-NumWords%4);j++){
447         Value=Buff.GetNextBackWord();
448       }//end for
449     }//end if
450
451     Int_t Packet[1024];
452     Int_t TimePos[345];
453     Int_t Tp=0;
454     for(Int_t i=0;i<345;i++)TimePos[i]=0;
455     for(Int_t i=0;i<1024;i++)Packet[i]=0;
456
457     Int_t NextTableType=0;
458     Int_t BunchLen=0;
459     Int_t Count=0;
460     for(Int_t i=0;i<NumWords;i++){
461       Value=Buff.GetNextBackWord();
462       Packet[i]=Value;
463       if(NextTableType==1){
464         TimePos[Tp]=i;
465         Tp++;
466       }
467       NextTable(Value,NextTableType,BunchLen,Count);
468     }//end for
469     //computing the Time gap between two bunches
470     Int_t temp=0;
471     Tp--;
472     Int_t PreviousTime=Packet[TimePos[Tp]];
473     for(Int_t i=Tp-1;i>=0;i--){
474       Int_t TimPos=TimePos[i];
475       Int_t BunchLen=Packet[TimPos-1]-2;
476       temp=Packet[TimPos];
477       Packet[TimPos]=Packet[TimPos]-PreviousTime-BunchLen;
478       PreviousTime=temp;
479     }//end for
480     NextTableType=0;
481     Count=0;
482     BunchLen=0;
483     Int_t TimeBin=0;
484     //All the words for one pad are compressed and stored in the compress file
485     for(Int_t i=0;i<NumWords;i++){
486       Value=Packet[i];
487       if(NextTableType==1)TimeBin=Value;
488       if(NextTableType>1){
489         //      ULong_t val=(ULong_t)table[NextTableType]->Code()[Value];     // val is the code
490         Double_t val=table[NextTableType]->Code()[Value];     // val is the code
491         UChar_t len=table[NextTableType]->CodeLen()[Value];  // len is the length (number of bits)of val
492         StoreValue(Mirror((ULong_t)val,len),len);
493         StoredWords++;
494       }//end if
495       NextTable(Value,NextTableType,BunchLen,Count);
496       if(NextTableType==0){
497         //      ULong_t val=(ULong_t)table[1]->Code()[TimeBin];     // val is the code
498         Double_t val=table[1]->Code()[TimeBin];     // val is the code
499         UChar_t len=table[1]->CodeLen()[TimeBin];  // len is the length (number of bits)of val
500         StoreValue(Mirror((ULong_t)val,len),len);
501         //val=(ULong_t)table[NextTableType]->Code()[(BunchLen+2)];     // val is the code
502         val=table[NextTableType]->Code()[(BunchLen+2)];     // val is the code
503         len=table[NextTableType]->CodeLen()[(BunchLen+2)];  // len is the length (number of bits)of val
504         StoreValue(Mirror((ULong_t)val,len),len);
505         StoredWords+=2;
506       }
507     }//end for
508     //Trailer
509     StoreValue(NumWords,10);
510     StoreValue(PadNum,10);
511     StoreValue(RowNum,10);
512     StoreValue(SecNum,9);
513     StoreValue(1,1);
514     StoredWords+=4;
515   }//end  while
516   StoreValue(NumPacket,32);
517   cout<<"Number of strored packet: "<<NumPacket<<endl;
518   StoreValue(1,1);
519   //The last buffen cannot be completely full
520   Flush();
521   cout<<"Number of stored words: "<<StoredWords<<endl;
522   f.close();
523   return 0;
524 }
525
526 /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
527 Int_t AliTPCCompression::CompressDataOptTables(Int_t NumTable,const char* fSource,const char* fDest){
528   if (fVerbose){
529     cout<<" BackWord COMPRESSION "<<endl;
530     cout<<"compression of the file "<<fSource<<" Output File: "<<fDest<<endl;
531   }
532   //Tables are read from the files (Each codeword has been "Mirrored")
533   AliTPCHTable ** table = new AliTPCHTable*[NumTable];
534   RetrieveTables(table,NumTable);
535   //the output file is open
536   f.open(fDest,ios::binary|ios::out);
537   // Source file is open
538   AliTPCBuffer160 Buff(fSource,0);
539   //coded words are written into a file
540   Int_t NumWords,PadNum,RowNum,SecNum=0;
541   ULong_t  StoredWords=0;
542   Int_t    Value=0;
543   ULong_t  NumPacket=0;
544   Double_t Stat[5]={0.,0.,0.,0.,0.};
545   ULong_t  TrailerNumber=0;
546   Double_t  NumElem[5]={0,0,0,0,0};
547   Double_t FillWords=0.;
548   stat.open("Statistics",ios::app);
549   stat<<endl;
550   stat<<"-------------------COMPRESSION STATISTICS----------"<<endl;
551   Int_t End=1;
552   while(Buff.ReadTrailerBackward(NumWords,PadNum,RowNum,SecNum) !=-1 ){
553     if(End){
554       FillWords=Buff.GetFillWordsNum();
555       End=0;
556     }//endif
557
558     NumPacket++;
559     if (NumWords%4){
560       FillWords+=4-NumWords%4;
561       for(Int_t j=0;j<(4-NumWords%4);j++){
562         Value=Buff.GetNextBackWord();
563       }//end for
564     }//end if
565
566     Int_t Packet[1024];
567     Int_t TimePos[345];
568     Int_t Tp=0;
569     for(Int_t i=0;i<345;i++)TimePos[i]=0;
570     for(Int_t i=0;i<1024;i++)Packet[i]=0;
571
572     Int_t NextTableType=0;
573     Int_t BunchLen=0;
574     Int_t Count=0;
575     for(Int_t i=0;i<NumWords;i++){
576       Value=Buff.GetNextBackWord();
577       Packet[i]=Value;
578       if(NextTableType==1){
579         TimePos[Tp]=i;
580         Tp++;
581       }
582       NextTable(Value,NextTableType,BunchLen,Count);
583     }//end for
584     //computing the Time gap between two bunches
585     Int_t temp=0;
586     Tp--;
587     Int_t PreviousTime=Packet[TimePos[Tp]];
588     for(Int_t i=Tp-1;i>=0;i--){
589       Int_t TimPos=TimePos[i];
590       Int_t BunchLen=Packet[TimPos-1]-2;
591       temp=Packet[TimPos];
592       Packet[TimPos]=Packet[TimPos]-PreviousTime-BunchLen;
593       PreviousTime=temp;
594     }//end for
595
596     NextTableType=0;
597     Count=0;
598     BunchLen=0;
599     Int_t TimeBin=0;
600     for(Int_t i=0;i<NumWords;i++){
601       Value=Packet[i];
602       if(NextTableType==1)TimeBin=Value;
603       if(NextTableType>1){
604         //ULong_t val=(ULong_t)table[NextTableType]->Code()[Value];     // val is the code
605         Double_t val=table[NextTableType]->Code()[Value];     // val is the code
606         UChar_t len=table[NextTableType]->CodeLen()[Value];  // len is the length (number of bits)of val
607         Stat[NextTableType]+=len;
608         NumElem[NextTableType]++;
609         StoreValue((ULong_t)val,len);
610         StoredWords++;
611       }//end if
612       NextTable(Value,NextTableType,BunchLen,Count);
613       if(NextTableType==0){
614         //      ULong_t val=(ULong_t)table[1]->Code()[TimeBin];     // val is the code
615         Double_t val=table[1]->Code()[TimeBin];     // val is the code
616         UChar_t len=table[1]->CodeLen()[TimeBin];  // len is the length (number of bits)of val
617         Stat[1]+=len;
618         NumElem[1]++;
619         StoreValue((ULong_t)val,len);
620         //      val=(ULong_t)table[NextTableType]->Code()[(BunchLen+2)];     // val is the code
621         val=table[NextTableType]->Code()[(BunchLen+2)];     // val is the code
622         len=table[NextTableType]->CodeLen()[(BunchLen+2)];  // len is the length (number of bits)of val
623         StoreValue((ULong_t)val,len);
624         Stat[NextTableType]+=len;
625         NumElem[NextTableType]++;
626         StoredWords+=2;
627       }
628     }//end for
629     //Trailer
630     StoreValue(NumWords,10);
631     StoreValue(PadNum,10);
632     StoreValue(RowNum,10);
633     StoreValue(SecNum,9);
634     StoreValue(1,1);
635     StoredWords+=4;
636     TrailerNumber++;
637   }//end  while
638   StoreValue(NumPacket,32);
639   if(fVerbose)
640     cout<<"Number of strored packet: "<<NumPacket<<endl;
641   StoreValue(1,1);
642   //The last buffen cannot be completely full
643   Flush();
644   if(fVerbose)
645     cout<<"Number of stored words: "<<StoredWords<<endl;
646   f.close();
647   //Tables are deleted
648   for(Int_t i=0;i<NumTable;i++){
649     delete table[i];
650   }//end for
651   delete [] table;
652   Double_t dimension=(ULong_t)ceil((Stat[0]+Stat[1]+Stat[2]+Stat[3]+Stat[4])/8)+TrailerNumber*5;
653   stat<<"Trailer Dimension in bytes......"<<TrailerNumber*5<<endl;
654   stat<<"Data Dimension in bytes........."<<(ULong_t)ceil((Stat[0]+Stat[1]+Stat[2]+Stat[3]+Stat[4])/8)<<endl;
655   stat<<"Compressed file dimension......."<<(ULong_t)dimension<<endl;
656   /*
657   stat<<(ULong_t)TrailerNumber<<endl;
658   stat<<(ULong_t)FillWords<<endl;
659   stat<<(ULong_t)NumElem[0]<<endl;
660   stat<<(ULong_t)NumElem[1]<<endl;
661   stat<<(ULong_t)NumElem[2]<<endl;
662   stat<<(ULong_t)NumElem[3]<<endl;
663   stat<<(ULong_t)NumElem[4]<<endl;
664   */
665   FillWords=(FillWords+NumElem[0]+NumElem[1]+NumElem[2]+NumElem[3]+NumElem[4]+TrailerNumber*4)*10/8;
666   stat<<"Original file dimension........."<<(ULong_t)FillWords<<endl;
667
668   Double_t ratio=(dimension/FillWords)*100;
669   stat<<"Compression ratio (Compressed/Uncompressed)..."<<ratio<<"%"<<endl;
670   stat<<endl;
671   stat<<"Bunch length size in bytes......"<<(ULong_t)ceil(Stat[0]/8)<<" Comppression.."<<(Stat[0]/NumElem[0])*10<<"%"<<endl;
672   
673   stat<<"Time gap size in bytes.........."<<(ULong_t)ceil(Stat[1]/8)<<" Comppression.."<<(Stat[1]/NumElem[1])*10<<"%"<<endl;
674   stat<<"Amplitude values in bytes......."<<(ULong_t)ceil((Stat[2]+Stat[3]+Stat[4])/8)<<" Comppression.."<<
675     ((Stat[2]+Stat[3]+Stat[4])/(NumElem[2]+NumElem[3]+NumElem[4]))*10<<"%"<<endl;
676   stat<<"     One Samples in bytes............"<<(ULong_t)ceil(Stat[2]/8)<<" Comppression.."<<(Stat[2]/NumElem[2])*10<<"%"<<endl;
677   stat<<"     Central Samples size in bytes..."<<(ULong_t)ceil(Stat[3]/8)<<" Comppression.."<<(Stat[3]/NumElem[3])*10<<"%"<<endl;
678   stat<<"     Border Samples size in bytes...."<<(ULong_t)ceil(Stat[4]/8)<<" Comppression.."<<(Stat[4]/NumElem[4])*10<<"%"<<endl;
679   stat<<endl;
680   stat<<"Average number of bits per word"<<endl;
681   stat<<"Bunch length ......"<<Stat[0]/NumElem[0]<<endl;
682   stat<<"Time gap .........."<<Stat[1]/NumElem[1]<<endl;
683   stat<<"One Samples........"<<Stat[2]/NumElem[2]<<endl;
684   stat<<"Central Samples ..."<<Stat[3]/NumElem[3]<<endl;
685   stat<<"Border Samples....."<<Stat[4]/NumElem[4]<<endl;
686   stat.close();
687   return 0;
688 }
689
690 ////////////////////////////////////////////////////////////////////////////////////////
691
692 ////////////////////////////////////////////////////////////////////////////////////////
693 /*                               DECOMPRESSION                                        */
694 ////////////////////////////////////////////////////////////////////////////////////////
695 void AliTPCCompression::CreateTrees(AliTPCHNode *RootNode[],const Int_t NumTables){
696   //The first part of the compressed file cotains the tables
697   //The following for loop is used to generate the Huffman trees acording to the tables
698   if(fVerbose)
699     cout<<"Creating the Huffman trees \n";
700   AliTPCHNode *node=0;
701   //  ULong_t Code;
702   Double_t Code;
703   UChar_t CodeLen;
704   //loop over the numbero of tables
705   for(Int_t k=0;k<NumTables;k++){
706     RootNode[k]=new AliTPCHNode(); //RootNode is the root of the tree
707     Int_t Dim;//this variable contains the table dimension
708     f.read((char*)(&Dim),sizeof(Int_t));
709     if (fVerbose)
710       cout<<"Table dimension: "<<Dim<<endl;
711     //loop over the words of a table
712     for(Int_t i=0;i<Dim;i++){
713       f.read((char*)(&CodeLen),sizeof(UChar_t));
714       //f.read((char*)(&Code),sizeof(ULong_t));
715       f.read((char*)(&Code),sizeof(Double_t));
716       node=RootNode[k];
717       for(Int_t j=1;j<=CodeLen;j++){
718         ULong_t bit,val=0;
719         val=(ULong_t)pow(2,CodeLen-j);
720         bit=(ULong_t)Code&val; 
721         AliTPCHNode *temp=node;
722         if(bit){
723           node=node->GetRight();
724           if(!node){
725             node=new AliTPCHNode();
726             temp->SetRight(node);
727           }//end if
728         }//end if
729         else{
730           node=node->GetLeft();
731           if(!node){
732             node=new AliTPCHNode();
733             temp->SetLeft(node);
734           }//end if
735         }//end else
736       }//end for
737       if(CodeLen){
738         node->SetSymbol(i);
739         node->SetFrequency(CodeLen);
740       }//end if
741       //cout<<node->GetSymbol()<<"  "<<(Int_t)node->GetFrequency()<<endl;
742     }//end for 
743   }//end for 
744   if (fVerbose)
745     cout<<"Trees generated \n";
746   //At this point the trees are been built
747 }
748 //////////////////////////////////////////////////////////////////////////////////////////////////
749 void AliTPCCompression::CreateTreesFromFile(AliTPCHNode *RootNode[],const Int_t NumTables){
750   if(fVerbose)
751     cout<<"Creating the Huffman trees \n";
752   AliTPCHNode *node=0;
753   // ULong_t Code;
754   Double_t Code;
755   UChar_t CodeLen;
756   ifstream fTable;  
757   char filename[15];
758   //The following for loop is used to generate the Huffman trees acording to the tables
759   //loop over the tables
760   for(Int_t k=0;k<NumTables;k++){
761     RootNode[k]=new AliTPCHNode(); //RootNode is the root of the tree
762     Int_t Dim=0;//this variable contains the table dimension
763     sprintf(filename,"Table%d.dat",k); 
764     fTable.open(filename,ios::binary);
765     fTable.read((char*)(&Dim),sizeof(Int_t));
766     if (fVerbose)
767       cout<<"Table dimension: "<<Dim<<endl;
768     //loop over the words of one table
769     for(Int_t i=0;i<Dim;i++){
770       fTable.read((char*)(&CodeLen),sizeof(UChar_t));
771       //fTable.read((char*)(&Code),sizeof(ULong_t));
772       fTable.read((char*)(&Code),sizeof(Double_t));
773       node=RootNode[k];
774       for(Int_t j=1;j<=CodeLen;j++){
775         ULong_t bit,val=0;
776         val=(ULong_t)pow(2,CodeLen-j);
777         bit=(ULong_t)Code&val; 
778         AliTPCHNode *temp=node;
779         if(bit){
780           node=node->GetRight();
781           if(!node){
782             node=new AliTPCHNode();
783             temp->SetRight(node);
784           }//end if 
785         }//end if
786         else{
787           node=node->GetLeft();
788           if(!node){
789             node=new AliTPCHNode();
790             temp->SetLeft(node);
791           }//end if
792         }//end else
793       }//end for
794       if(CodeLen){
795         node->SetSymbol(i);
796         node->SetFrequency(CodeLen);
797       }//end if
798     }//end for 
799     fTable.close();
800   }//end for 
801   if (fVerbose)
802     cout<<"Trees generated \n";
803   //At this point the trees are been built
804 }
805 //////////////////////////////////////////////////////////////////////////////////////////////////
806 void AliTPCCompression::DeleteHuffmanTree(AliTPCHNode* node){
807   //This function deletes all the nodes of an Huffman tree
808   //In an Huffman tree any internal node has always two children
809   if (node){
810     DeleteHuffmanTree(node->GetLeft());
811     DeleteHuffmanTree(node->GetRight());
812     //    cout<<node->GetSymbol()<<"  "<<(Int_t)node->GetFrequency()<<endl;
813     delete node;
814   }
815 }
816 //////////////////////////////////////////////////////////////////////////////////////////////////
817 void AliTPCCompression::VisitHuffmanTree(AliTPCHNode* node){
818   //This function realizes an in order visit of a binary tree
819   if (node){
820     cout<<node->GetSymbol()<<" "<<node->GetFrequency()<<endl;
821     VisitHuffmanTree(node->GetLeft());
822     VisitHuffmanTree(node->GetRight());
823   }
824 }
825 //////////////////////////////////////////////////////////////////////////////////////////////////
826 ULong_t AliTPCCompression::ReadWord(Int_t NumberOfBit){
827   ULong_t Result=0;
828   ULong_t bit=0;
829   for (Int_t i=0;i<NumberOfBit;i++){
830     if (fReadBits==32){
831       fPos-=sizeof(ULong_t);
832       f.seekg(fPos);
833       f.read((char*)(&fBuffer),sizeof(ULong_t));
834       fReadBits=0;
835     }//end if
836     ULong_t mask=0;
837     mask=(ULong_t)pow(2,fReadBits);
838     bit=fBuffer&mask;
839     bit=bit>>fReadBits;
840     fReadBits++;
841     bit=bit<<i;
842     Result=Result|bit;
843   }//end for
844   return Result;
845 }
846 //////////////////////////////////////////////////////////////////////////////////////////////////
847 void AliTPCCompression::ReadTrailer(Int_t &WordsNumber,Int_t &PadNumber,Int_t &RowNumber,Int_t &SecNumber){
848   ReadWord(1);
849   SecNumber=ReadWord(9);
850   RowNumber=ReadWord(10);
851   PadNumber=ReadWord(10);
852   WordsNumber=ReadWord(10);
853   return;
854 }
855 //////////////////////////////////////////////////////////////////////////////////////////////////
856 ULong_t AliTPCCompression::GetDecodedWord(AliTPCHNode* root){
857   AliTPCHNode *node=root;
858   ULong_t symbol=0;
859   Bool_t decoded=0;
860   while(!decoded){
861     ULong_t bit=ReadWord(1);
862     if(bit)
863       node=node->GetRight();
864     else
865       node=node->GetLeft();
866     if (!(node->GetLeft())){
867       symbol=node->GetSymbol();
868       decoded=1;
869     }
870   }//end while
871   return symbol;
872 }
873 //////////////////////////////////////////////////////////////////////////////////////////////////
874 Int_t AliTPCCompression::DecompressData(Int_t NumTables,const char* fname,char* fDest){
875   cout<<"   DECOMPRESSION:"<<endl;
876   cout<<"Source File "<<fname<<" Destination File "<<fDest<<endl; 
877   f.open(fname,ios::binary|ios::in);
878   if(!f){cout<<"File doesn't exist\n";return -1;}
879   AliTPCHNode ** RootNode = new AliTPCHNode*[NumTables];
880   //Creation of the Huffman trees
881   CreateTrees(RootNode,NumTables);
882   //to go to the end of the file
883   f.seekg(0,ios::end);
884   //to get the file dimension in byte
885   fPos=f.tellg();
886   fPos-=sizeof(ULong_t);
887   f.seekg(fPos);
888   fReadBits=0;
889   fBuffer=0;
890   f.read((char*)(&fBuffer),sizeof(ULong_t));
891   Int_t bit=0;
892   ULong_t Mask=0x1;
893   while(!bit){
894     bit=fBuffer&Mask;
895     Mask=Mask<<1;
896     fReadBits++;
897   }
898   ULong_t PacketNumber=ReadWord(sizeof(ULong_t)*8);
899   cout<<"Number of Packect: "<<PacketNumber<<endl;
900   AliTPCBuffer160 BufferFile(fDest,1);
901   ULong_t k=0;
902   ULong_t WordsRead=0; //number of read coded words 
903   while(k<PacketNumber){
904     Int_t NumWords,PadNumber,RowNumber,SecNumber=0;
905     ReadTrailer(NumWords,PadNumber,RowNumber,SecNumber);
906     k++;
907     WordsRead+=4;
908     Int_t PreviousTime=-1;
909     Int_t Time=0;
910     Int_t NextTableType=0;
911     Int_t BunchLen=0;
912     Int_t Count=0;
913     for(Int_t i=0;i<NumWords;i++){
914       ULong_t symbol=GetDecodedWord(RootNode[NextTableType]);
915       WordsRead++;
916       //Time reconstruction
917       if (NextTableType==1){
918         if (PreviousTime!=-1){
919           PreviousTime=symbol+PreviousTime+BunchLen;
920         }
921         else PreviousTime=symbol;
922         Time=PreviousTime;
923       }
924       if(NextTableType>1)
925         BufferFile.FillBuffer(symbol);
926       NextTable(symbol,NextTableType,BunchLen,Count); 
927       if(NextTableType==0){
928         BufferFile.FillBuffer(Time);
929         BufferFile.FillBuffer(BunchLen+2);
930         BunchLen=0;
931       }
932     }//end for
933     BufferFile.WriteTrailer(NumWords,PadNumber,RowNumber,SecNumber);
934   }//end while
935   cout<<"Number of decoded words:"<<WordsRead<<endl;
936   f.close();
937   //The trees are deleted 
938   for(Int_t j=0;j<NumTables;j++){
939       DeleteHuffmanTree(RootNode[j]);
940   }//end for
941   delete [] RootNode; 
942   return 0; 
943 }
944
945 //////////////////////////////////////////////////////////////////////////////////////////////////////////////////
946
947 Int_t AliTPCCompression::DecompressDataOptTables(Int_t NumTables,const char* fname,char* fDest){
948   if(fVerbose){
949     cout<<"   DECOMPRESSION:"<<endl;
950     cout<<"Source File "<<fname<<" Destination File "<<fDest<<endl; 
951   }
952   AliTPCHNode ** RootNode = new AliTPCHNode*[NumTables];
953   //Creation of the Huffman trees
954   CreateTreesFromFile(RootNode,NumTables);
955   f.open(fname,ios::binary|ios::in);
956   if(!f){cout<<"File doesn't exist\n";return -1;}
957   //to go to the end of the file
958   f.seekg(0,ios::end);
959   //to get the file dimension in byte
960   fPos=f.tellg();
961   fPos-=sizeof(ULong_t);
962   f.seekg(fPos);
963   fReadBits=0;
964   fBuffer=0;
965   f.read((char*)(&fBuffer),sizeof(ULong_t));
966   Int_t bit=0;
967   ULong_t Mask=0x1;
968   while(!bit){
969     bit=fBuffer&Mask;
970     Mask=Mask<<1;
971     fReadBits++;
972   }
973   ULong_t PacketNumber=ReadWord(sizeof(ULong_t)*8);
974   if(fVerbose){
975     cout<<"Number of Packect: "<<PacketNumber<<endl;
976   }
977   AliTPCBuffer160 BufferFile(fDest,1);
978   ULong_t k=0;
979   ULong_t WordsRead=0; //number of read coded words 
980   while(k<PacketNumber){
981     Int_t NumWords,PadNumber,RowNumber,SecNumber=0;
982     ReadTrailer(NumWords,PadNumber,RowNumber,SecNumber);
983     k++;
984     WordsRead+=4;
985     Int_t PreviousTime=-1;
986     Int_t Time=0;
987     Int_t NextTableType=0;
988     Int_t BunchLen=0;
989     Int_t Count=0;
990     for(Int_t i=0;i<NumWords;i++){
991       ULong_t symbol=GetDecodedWord(RootNode[NextTableType]);
992       WordsRead++;
993       //Time reconstruction
994       if (NextTableType==1){
995         if (PreviousTime!=-1){
996           PreviousTime=symbol+PreviousTime+BunchLen;
997         }
998         else PreviousTime=symbol;
999         Time=PreviousTime;
1000       }
1001       if(NextTableType>1)
1002         BufferFile.FillBuffer(symbol);
1003       NextTable(symbol,NextTableType,BunchLen,Count); 
1004       if(NextTableType==0){
1005         BufferFile.FillBuffer(Time);
1006         BufferFile.FillBuffer(BunchLen+2);
1007         BunchLen=0;
1008       }
1009     }//end for
1010     BufferFile.WriteTrailer(NumWords,PadNumber,RowNumber,SecNumber);
1011   }//end while
1012   if(fVerbose){
1013     cout<<"Number of decoded words:"<<WordsRead<<endl;
1014   }
1015   f.close();
1016   //The trees are deleted 
1017   for(Int_t j=0;j<NumTables;j++){
1018       DeleteHuffmanTree(RootNode[j]);
1019   }//end for
1020   delete [] RootNode;
1021   return 0; 
1022 }
1023
1024 ///////////////////////////////////////////////////////////////////////////////////////////
1025
1026 void AliTPCCompression::ReadAltroFormat(char* fileOut,char* fileIn){
1027   ofstream ftxt(fileOut);
1028   AliTPCBuffer160 Buff(fileIn,0);
1029   Int_t NumWords,PadNum,RowNum,SecNum=0;
1030   Int_t Value=0;
1031   while(Buff.ReadTrailerBackward(NumWords,PadNum,RowNum,SecNum) !=-1 ){
1032     ftxt<<"W:"<<NumWords<<" P:"<<PadNum<<" R:"<<RowNum<<" S:"<<SecNum<<endl;
1033     if (NumWords%4){
1034       for(Int_t j=0;j<(4-NumWords%4);j++){
1035         Value=Buff.GetNextBackWord();
1036       }//end for
1037     }//end if
1038     for(Int_t i=0;i<NumWords;i++){
1039       Value=Buff.GetNextBackWord();
1040       ftxt<<Value<<endl;
1041     }//end for
1042   }//end while
1043   ftxt.close();
1044   return;
1045 }
1046
1047 //////////////////////////////////////////////////////////////////////////////////////////