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