Update to the current version in the Bergen CVS. Most important
[u/mrichter/AliRoot.git] / HLT / comp / AliL3Compress.cxx
1 // @(#) $Id$
2
3 // Author: Anders Vestbo <mailto:vestbo$fi.uib.no>
4 //*-- Copyright &copy ALICE HLT Group
5
6 #include "AliL3StandardIncludes.h"
7
8 #include "bitio.h"
9 #include "AliL3RootTypes.h"
10 #include "AliL3Models.h"
11 #include "AliL3DigitData.h"
12 #include "AliL3Logging.h"
13 #include "AliL3TrackArray.h"
14 #include "AliL3ModelTrack.h"
15 #include "AliL3Transform.h"
16 #include "AliL3MemHandler.h"
17 #include "AliL3DataCompressor.h"
18 #include "AliL3SpacePointData.h"
19
20 #if 0
21 #ifdef use_root
22 #include <TH1.h>
23 #include <TH2.h>
24 #include <TRandom.h>
25 #endif
26 #ifdef use_aliroot
27 #include "AliL3FileHandler.h"
28 #endif
29 #endif
30
31 #include "AliL3Compress.h"
32
33 #if GCCVERSION == 3
34 using namespace std;
35 #endif
36
37 //_____________________________________________________________
38 //
39 //  AliL3Compress
40 //
41 // Class for compressing and uncompressing data.
42
43 ClassImp(AliL3Compress)
44
45 AliL3Compress::AliL3Compress()
46 {
47   fTracks=0;
48   fSlice =0;
49   fPatch=0;
50   fWriteShape=kFALSE;
51   fEvent=-1;
52 }
53
54 AliL3Compress::AliL3Compress(Int_t slice,Int_t patch,Char_t *path,Bool_t writeshape,Int_t event)
55 {
56   fEvent=event;
57   fSlice=slice;
58   fPatch=patch;
59   fTracks=0;
60   sprintf(fPath,"%s",path);
61   fWriteShape=writeshape;
62 }
63
64 AliL3Compress::~AliL3Compress()
65 {
66   if(fTracks)
67     delete fTracks;
68 }
69
70 Bool_t AliL3Compress::WriteFile(AliL3TrackArray *tracks,Char_t *filename)
71 {
72   Char_t fname[1024];
73   if(filename)
74     sprintf(fname,"%s/comp/%s",fPath,filename);
75   else if(fEvent<0)
76     sprintf(fname,"%s/comp/tracks_m_%d_%d.raw",fPath,fSlice,fPatch);
77   else
78     sprintf(fname,"%s/comp/tracks_m_%d_%d_%d.raw",fPath,fEvent,fSlice,fPatch);
79
80   FILE *file = fopen(fname,"w");
81   if(!file)
82     {
83       cerr<<"AliL3Compress::WriteFile : Error opening file "<<fname<<endl;
84       return kFALSE;
85     }
86   Short_t ntracks = tracks->GetNTracks();
87     
88   Int_t count=0;
89   AliL3ClusterModel *clusters=0;
90   AliL3TrackModel *model=0;
91   for(Int_t i=0; i<ntracks; i++)
92     {
93       AliL3ModelTrack *track = (AliL3ModelTrack*)tracks->GetCheckedTrack(i);
94       if(!track) continue;
95
96       //Do not save useless tracks or clusters:
97       //if(track->GetNPresentClusters() == 0)
98       //continue;
99       
100       track->FillModel();
101       model = track->GetModel();
102       //if(model->fNClusters==0) continue;
103       clusters = track->GetClusters();
104       if(fwrite(model,sizeof(AliL3TrackModel),1,file)!=1) break;
105       //if(fwrite(clusters,model->fNClusters*sizeof(AliL3ClusterModel),1,file)!=1) break;
106       if(fwrite(clusters,AliL3Transform::GetNRows(fPatch)*sizeof(AliL3ClusterModel),1,file)!=1) break;
107       count++;
108       
109     }
110   fclose(file);
111   return kTRUE;
112 }
113
114 Bool_t AliL3Compress::ReadFile(Char_t which,Char_t *filename)
115 {
116   //Read the trackfile.
117
118   Char_t fname[1024];
119   if(filename)
120     sprintf(fname,"%s/comp/%s",fPath,filename);
121   else
122     {
123       if(which == 'm')
124         {
125           if(fEvent<0)
126             sprintf(fname,"%s/comp/tracks_m_%d_%d.raw",fPath,fSlice,fPatch);
127           else
128             sprintf(fname,"%s/comp/tracks_m_%d_%d_%d.raw",fPath,fEvent,fSlice,fPatch);
129         }
130       else if(which == 'u')
131         {
132           if(fEvent<0)
133             sprintf(fname,"%s/comp/tracks_u_%d_%d.raw",fPath,fSlice,fPatch);
134           else
135             sprintf(fname,"%s/comp/tracks_u_%d_%d_%d.raw",fPath,fEvent,fSlice,fPatch);
136         }
137       else
138         {
139           cerr<<"AliL3Compress::ReadFile() : Wrong option"<<endl;
140           return kFALSE;
141         }
142     }
143
144   FILE *file = fopen(fname,"r");
145   if(!file)
146     {
147       cerr<<"AliL3Compress::ReadFile : Cannot open file "<<fname<<endl;
148       return kFALSE;
149     }
150
151   if(fTracks)
152     delete fTracks;
153   fTracks = new AliL3TrackArray("AliL3ModelTrack");
154   
155   while(!feof(file))
156     {
157       AliL3ModelTrack *track = (AliL3ModelTrack*)fTracks->NextTrack();
158       track->Init(fSlice,fPatch);
159       AliL3TrackModel *model = track->GetModel();
160       AliL3ClusterModel *clusters = track->GetClusters();
161       if(fread(model,sizeof(AliL3TrackModel),1,file)!=1) break;
162       //if(fread(clusters,model->fNClusters*sizeof(AliL3ClusterModel),1,file)!=1) break;
163       if(fread(clusters,AliL3Transform::GetNRows(fPatch)*sizeof(AliL3ClusterModel),1,file)!=1) break;
164       track->FillTrack();
165     }
166
167   fTracks->RemoveLast();
168   fclose(file);
169   return kTRUE;
170 }
171
172 Bool_t AliL3Compress::CompressFile()
173 {
174   Char_t fname[100];
175   if(fEvent<0)
176     sprintf(fname,"%s/comp/tracks_c_%d_%d.raw",fPath,fSlice,fPatch);
177   else
178     sprintf(fname,"%s/comp/tracks_c_%d_%d_%d.raw",fPath,fEvent,fSlice,fPatch);
179   BIT_FILE *output = OpenOutputBitFile(fname);
180   
181   if(fEvent<0)
182     sprintf(fname,"%s/comp/tracks_m_%d_%d.raw",fPath,fSlice,fPatch);
183   else
184     sprintf(fname,"%s/comp/tracks_m_%d_%d_%d.raw",fPath,fEvent,fSlice,fPatch);
185   
186   FILE *input = fopen(fname,"r");
187   if(!input)
188     {
189       cerr<<"AliL3Compress::CompressFile() : Error opening file: "<<fname<<endl;
190       return kFALSE;
191     }
192
193   AliL3TrackModel track;
194   AliL3ClusterModel cluster;
195   Int_t temp;
196   Int_t power;
197   
198   Int_t timeo,pado,chargeo,padshapeo,timeshapeo;
199   timeo=pado=chargeo=padshapeo=timeshapeo=0;
200   while(!feof(input))
201     {
202       if(fread(&track,sizeof(AliL3TrackModel),1,input)!=1) break;
203       
204       if(output->mask != 0x80) //Write the current byte to file.
205         {
206           //cerr<<"\nAliL3Compress::CompressFile() : Writing overhead bits!!!"<<endl;
207           if(putc(output->rack,output->file )!=output->rack)
208             cerr<<"AliL3Compress::ComressFile : Error writing to bitfile"<<endl;
209           output->mask=0x80;
210           output->rack=0;
211         }
212       
213       //Write track parameters:
214       fwrite(&track,sizeof(AliL3TrackModel),1,output->file);
215       
216       Int_t origslice=-1,slice,clustercount=0;
217       for(Int_t i=0; i<AliL3Transform::GetNRows(fPatch); i++)
218         {
219           if(fread(&cluster,sizeof(AliL3ClusterModel),1,input)!=1) break;
220           
221           //Write empty flag:
222           temp = (Int_t)cluster.fPresent;
223           OutputBit(output,temp);
224           if(!temp) continue;
225           
226           if(cluster.fSlice<0 || cluster.fSlice>35)
227             {
228               cerr<<"AliL3DataCompress::CompressFile : Fucked up slice number :"<<cluster.fSlice<<endl;
229               exit(5);
230             }
231           
232           //Write slice number of first point
233           if(clustercount==0)
234             {
235               origslice = cluster.fSlice;
236               OutputBits(output,origslice,6); //Need 6 bits to encode slice number
237             }
238           else
239             {
240               slice = cluster.fSlice;
241               if(slice == origslice)
242                 OutputBit(output,0);  //No change of slice
243               else
244                 {
245                   OutputBit(output,1);
246                   OutputBits(output,slice,6);
247                   origslice=slice;
248                 }
249             }
250           
251           //Write time information:
252           temp = (Int_t)rint(cluster.fDTime);
253           if(temp<0)
254             OutputBit(output,0);
255           else
256             OutputBit(output,1);
257           power = 1<<(AliL3DataCompressor::GetNTimeBits()-1);
258           if(abs(temp)>=power)
259             {
260               timeo++;
261               temp=power - 1;
262             }
263           OutputBits(output,abs(temp),(AliL3DataCompressor::GetNTimeBits()-1));
264           
265           //Write pad information:
266           temp = (Int_t)rint(cluster.fDPad);
267           if(temp<0)
268             OutputBit(output,0);
269           else
270             OutputBit(output,1);
271           power = 1<<(AliL3DataCompressor::GetNPadBits()-1);
272           if(abs(temp)>=power)
273             {
274               pado++;
275               temp=power - 1;
276             }
277           OutputBits(output,abs(temp),(AliL3DataCompressor::GetNPadBits()-1));
278           
279           //Write charge information:
280           temp = (Int_t)cluster.fDCharge;
281           power = 1<<(AliL3DataCompressor::GetNChargeBits());
282           if(abs(temp)>=power)
283             {
284               chargeo++;
285               temp=power - 1;
286             }
287           OutputBits(output,abs(temp),(AliL3DataCompressor::GetNChargeBits()));
288           
289           if(fWriteShape)
290             {
291               //Write shape information:
292               temp = (Int_t)cluster.fDSigmaY2;
293               if(temp<0)
294                 OutputBit(output,0);
295               else
296                 OutputBit(output,1);
297               power = 1<<(AliL3DataCompressor::GetNPadShapeBits()-1);
298               if(abs(temp) >= power)
299                 {
300                   padshapeo++;
301                   temp = power - 1;
302                 }
303               OutputBits(output,abs(temp),(AliL3DataCompressor::GetNPadShapeBits()-1));
304               
305               temp = (Int_t)cluster.fDSigmaZ2;
306               if(temp<0)
307                 OutputBit(output,0);
308               else
309                 OutputBit(output,1);
310               power = 1<<(AliL3DataCompressor::GetNTimeShapeBits()-1);
311               if(abs(temp) >= power)
312                 {
313                   timeshapeo++;
314                   temp=power - 1;
315                 }
316               OutputBits(output,abs(temp),(AliL3DataCompressor::GetNTimeShapeBits()-1));
317             }
318           
319           clustercount++;
320         }
321     }
322   
323   fclose(input);
324   CloseOutputBitFile(output);
325   if(pado || timeo || chargeo || padshapeo || timeshapeo)
326     {
327       cout<<endl<<"Saturations: "<<endl
328           <<"Pad "<<pado<<endl
329           <<"Time "<<timeo<<endl
330           <<"Charge "<<chargeo<<endl
331           <<"Padshape "<<padshapeo<<endl
332           <<"Timeshape "<<timeshapeo<<endl<<endl;
333     }
334   return kTRUE;
335 }
336
337 Bool_t AliL3Compress::ExpandFile()
338 {
339   Char_t fname[100];
340   if(fEvent<0)
341     sprintf(fname,"%s/comp/tracks_c_%d_%d.raw",fPath,fSlice,fPatch);
342   else
343     sprintf(fname,"%s/comp/tracks_c_%d_%d_%d.raw",fPath,fEvent,fSlice,fPatch);
344   BIT_FILE *input = OpenInputBitFile(fname);
345   
346   if(fEvent<0)
347     sprintf(fname,"%s/comp/tracks_u_%d_%d.raw",fPath,fSlice,fPatch);
348   else
349     sprintf(fname,"%s/comp/tracks_u_%d_%d_%d.raw",fPath,fEvent,fSlice,fPatch);
350   FILE *output = fopen(fname,"w");
351   if(!output)
352     {
353       cerr<<"AliL3Compress::ExpandFile() : Error opening file: "<<fname<<endl;
354       return kFALSE;
355     }
356
357   AliL3TrackModel trackmodel;
358   AliL3ClusterModel *clusters=0;
359   Int_t count=0;
360   
361   clusters = new AliL3ClusterModel[(AliL3Transform::GetNRows(fPatch))];
362   while(!feof(input->file))
363     {
364       input->mask=0x80;//make sure we read a new byte from file.
365       
366       //Read and write track:
367       if(fread(&trackmodel,sizeof(AliL3TrackModel),1,input->file)!=1) break;
368       fwrite(&trackmodel,sizeof(AliL3TrackModel),1,output);
369
370       memset(clusters,0,AliL3Transform::GetNRows(fPatch)*sizeof(AliL3ClusterModel));
371       Int_t origslice=-1,clustercount=0;
372       for(Int_t i=0; i<AliL3Transform::GetNRows(fPatch); i++)
373         {
374           Int_t temp,sign;
375           
376           //Read empty flag:
377           temp = InputBit(input);
378           if(!temp) 
379             {
380               clusters[i].fPresent=kFALSE;
381               continue;
382             }
383           clusters[i].fPresent=kTRUE;
384           
385           //Read slice information
386           if(clustercount==0)
387             {
388               temp = InputBits(input,6);
389               clusters[i].fSlice = temp;
390               origslice = temp;
391             }
392           else
393             {
394               temp = InputBit(input);
395               if(!temp)//no change
396                 clusters[i].fSlice = origslice;
397               else
398                 {
399                   temp = InputBits(input,6);//read new slice
400                   clusters[i].fSlice = temp;
401                   origslice = temp;//store new slice
402                 }
403             }
404           
405           //Read time information:
406           sign=InputBit(input);
407           temp = InputBits(input,(AliL3DataCompressor::GetNTimeBits()-1));
408           if(!sign)
409             temp*=-1;
410           clusters[i].fDTime = temp;
411           
412           //Read pad information:
413           sign=InputBit(input);
414           temp = InputBits(input,(AliL3DataCompressor::GetNPadBits()-1));
415           if(!sign)
416             temp*=-1;
417           clusters[i].fDPad = temp;
418           
419           //Read charge information:
420           temp=InputBits(input,(AliL3DataCompressor::GetNChargeBits()));
421           clusters[i].fDCharge = temp;
422           
423           if(fWriteShape)
424             {
425               //Read shape information:
426               sign = InputBit(input);
427               temp = InputBits(input,(AliL3DataCompressor::GetNPadShapeBits()-1));
428               if(!sign)
429                 temp*=-1;
430               clusters[i].fDSigmaY2 = temp;
431               
432               sign = InputBit(input);
433               temp = InputBits(input,(AliL3DataCompressor::GetNTimeShapeBits()-1));
434               if(!sign)
435                 temp*=-1;
436               clusters[i].fDSigmaZ2 = temp;
437             }
438           clustercount++;
439         }
440       count++;
441       //fwrite(clusters,(trackmodel.fNClusters)*sizeof(AliL3ClusterModel),1,output);
442       fwrite(clusters,AliL3Transform::GetNRows(fPatch)*sizeof(AliL3ClusterModel),1,output);
443     }
444   
445   delete [] clusters;
446   fclose(output);
447   CloseInputBitFile(input);
448   return kTRUE;
449 }
450
451 void AliL3Compress::CompressRemaining(AliL3SpacePointData *clusters[36][6],UInt_t nclusters[36][6])
452 {
453   //Write the remaining clusters in a compressed format.
454
455   Char_t filename[1024];
456   Int_t nrows = AliL3Transform::GetNRows();
457   Int_t *npoints = new Int_t[nrows];
458   for(Int_t slice=0; slice<=35; slice++)
459     {
460       for(Int_t patch=0; patch < 1; patch++)
461         {
462           sprintf(filename,"%s/comp/remains_%d_%d_%d.raw",fPath,fEvent,slice,-1);
463           BIT_FILE *output = OpenOutputBitFile(filename);
464           if(!output)
465             {
466               cerr<<"AliL3Compress::CompressRemaining : Cannot open file "<<filename<<endl;
467               exit(5);
468             }
469           AliL3SpacePointData *cl = clusters[slice][patch];
470           memset(npoints,0,nrows*sizeof(Int_t));
471           
472           UInt_t i;
473           for(i=0; i<nclusters[slice][patch]; i++)
474             {
475               if(cl[i].fCharge == 0) continue; //has been used
476               npoints[cl[i].fPadRow]++;
477             }
478           Int_t rowspresent=0;
479           for(Int_t j=0; j<nrows; j++)
480             {
481               if(!npoints[j]) continue;
482               rowspresent++;
483             }
484           
485           //Write number of padrows with clusters
486           OutputBits(output,rowspresent,8);
487           
488           Int_t last_padrow=-1;
489           for(i=0; i<nclusters[slice][patch]; i++)
490             {
491               if(cl[i].fCharge == 0) continue; //has been used
492               Int_t padrow = cl[i].fPadRow;
493               if(padrow != last_padrow)
494                 {
495                   OutputBits(output,padrow,8);//Write padrow #
496                   if(npoints[padrow] >= 1<<10)
497                     {
498                       cerr<<"AliL3Compress::CompressRemaining : Too many remaining clusters "<<npoints[padrow]<<endl;
499                       exit(5);
500                     }
501                   OutputBits(output,npoints[padrow],10);//Write number of clusters on this padrow
502                   last_padrow = padrow;
503                 }
504               
505               Float_t xyz[3] = {cl[i].fX,cl[i].fY,cl[i].fZ};
506               Int_t sector,row,buff;
507               AliL3Transform::Slice2Sector(slice,padrow,sector,row);
508               AliL3Transform::Global2Raw(xyz,sector,row);
509               
510               Float_t padw = cl[i].fSigmaY2 / AliL3Transform::GetPadPitchWidth(AliL3Transform::GetPatch(padrow));
511               Float_t timew = cl[i].fSigmaZ2 / AliL3Transform::GetZWidth();
512               
513               //Write pad
514               buff = (Int_t)rint(xyz[1]*10);
515               if(buff<0)
516                 {
517                   cerr<<"AliL3Compress:CompressRemaining : Wrong pad value "<<buff<<endl;
518                   exit(5);
519                 }
520               OutputBits(output,buff,11);
521
522               //Write time
523               buff = (Int_t)rint(xyz[2]*10);
524               if(buff<0)
525                 {
526                   cerr<<"AliL3Compress:CompressRemaining : Wrong time value "<<buff<<endl;
527                   exit(5);
528                 }
529               OutputBits(output,buff,13);
530
531               //Write widths
532               buff = (Int_t)rint(padw*100);
533               OutputBits(output,buff,8);
534               buff = (Int_t)rint(timew*100);
535               OutputBits(output,buff,8);
536               
537               //Write charge 
538               buff = cl[i].fCharge;
539               if(buff >= 1<<14)
540                 buff = (1<<14)-1;
541               OutputBits(output,buff,14);
542             }
543           
544           CloseOutputBitFile(output);
545         }
546       
547     }
548   delete [] npoints;
549 }
550
551 void AliL3Compress::ExpandRemaining(TempCluster **clusters,Int_t *ncl,const Int_t maxpoints)
552 {
553   //Expand the remaining clusters stored using function CompressRemaining
554   
555   Char_t filename[1024];
556   Int_t buff;
557   for(Int_t slice=0; slice<=35; slice++)
558     {
559       for(Int_t p=0; p<1; p++)
560         {
561           sprintf(filename,"%s/comp/remains_%d_%d_%d.raw",fPath,fEvent,slice,-1);
562           BIT_FILE *input = OpenInputBitFile(filename);
563           
564           //Read number of padrows 
565           buff = InputBits(input,8);
566           Int_t nrows = buff;
567           
568           for(Int_t i=0; i<nrows; i++)
569             {
570               //Read padrow
571               buff = InputBits(input,8);
572               Int_t padrow = buff;
573               
574               //Read nclusters;
575               buff = InputBits(input,10);
576               Int_t npoints = buff;
577               
578               for(Int_t i=0; i<npoints; i++)
579                 {
580                   clusters[slice][ncl[slice]].padrow = padrow;
581
582                   //Read pad
583                   buff = InputBits(input,11);
584                   clusters[slice][ncl[slice]].pad = (Float_t)buff/10.;
585                   
586                   //Read time
587                   buff = InputBits(input,13);
588                   clusters[slice][ncl[slice]].time = (Float_t)buff/10.;
589                   
590                   //Read widths 
591                   buff = InputBits(input,8);
592                   clusters[slice][ncl[slice]].sigmaY2 = (Float_t)buff/100.;
593                   buff = InputBits(input,8);
594                   clusters[slice][ncl[slice]].sigmaZ2 = (Float_t)buff/100.;
595                   
596                   //Read charge
597                   buff = InputBits(input,14);
598                   clusters[slice][ncl[slice]].charge = buff;
599                   //cout<<"padrow "<<padrow<<" pad "<<clusters[slice][ncl[slice]].pad<<" time "<<clusters[slice][ncl[slice]].time<<" charge "<<clusters[slice][ncl[slice]].charge<<" widths "<<clusters[slice][ncl[slice]].sigmaY2<<" "<<clusters[slice][ncl[slice]].sigmaZ2<<endl;
600                   ncl[slice]++;
601                 }
602               
603             }
604           CloseInputBitFile(input);
605         }
606     }
607 }
608
609 void AliL3Compress::PrintCompRatio(ofstream *outfile)
610 {
611   AliL3MemHandler *mem = new AliL3MemHandler();
612   Char_t fname[1024];
613   UInt_t remain_size=0,digit_size=0;
614   for(Int_t i=0; i<36; i++)
615     {
616       if(fEvent<0)
617         sprintf(fname,"%s/comp/remains_%d_%d.raw",fPath,i,-1);
618       else
619         sprintf(fname,"%s/comp/remains_%d_%d_%d.raw",fPath,fEvent,i,-1);
620       mem->SetBinaryInput(fname);
621       remain_size += mem->GetFileSize();
622       mem->CloseBinaryInput();
623
624       sprintf(fname,"%s/binaries/digits_c8_%d_%d_%d.raw",fPath,fEvent,i,-1);
625       mem->SetBinaryInput(fname);
626       digit_size += mem->GetFileSize();
627       mem->CloseBinaryInput();
628     }
629   
630   
631   if(fEvent<0)
632     sprintf(fname,"%s/comp/tracks_c_%d_%d.raw",fPath,fSlice,fPatch);
633   else
634     sprintf(fname,"%s/comp/tracks_c_%d_%d_%d.raw",fPath,fEvent,fSlice,fPatch);
635
636   mem->SetBinaryInput(fname);
637   UInt_t compress_size = mem->GetFileSize();
638   mem->CloseBinaryInput();
639   
640   if(digit_size==0)
641     {
642       cerr<<"AliL3Compress::PrintCompRatio : Zero digit size, not able to obtain comp. ratios!"<<endl;
643       return;
644     }
645   
646   Float_t compratio = (Float_t)(compress_size + remain_size)/(Float_t)digit_size;
647   if(outfile)
648     {
649       ofstream &out = *outfile;
650       out<<compress_size<<' '<<remain_size<<' '<<digit_size<<endl;
651     }
652
653   cout<<"=========================================="<<endl;
654   cout<<"Original digits size : "<<digit_size/1000<<" kByte ( 100 % )"<<endl;
655   cout<<"Compressed file size : "<<compress_size/1000<<" kByte ( "<<(Float_t)compress_size*100/(Float_t)digit_size<<" % )"<<endl;
656   cout<<"Remainig file size   : "<<remain_size/1000<<" kByte ( "<<(Float_t)remain_size*100/(Float_t)digit_size<<" % )"<<endl;
657   cout<<"---------------------- "<<endl;
658   cout<<"Compression ratio    : "<<compratio*100<<" %"<<endl;
659   cout<<"=========================================="<<endl;
660 }