ALIROOT-5433 Transition to CDHv3 in HLT
[u/mrichter/AliRoot.git] / HLT / TPCLib / comp / AliHLTTPCDataCompressionComponent.cxx
1 // $Id$
2 //**************************************************************************
3 //* This file is property of and copyright by the                          * 
4 //* ALICE Experiment at CERN, All rights reserved.                         *
5 //*                                                                        *
6 //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no>        *
7 //*                  for The ALICE HLT Project.                            *
8 //*                                                                        *
9 //* Permission to use, copy, modify and distribute this software and its   *
10 //* documentation strictly for non-commercial purposes is hereby granted   *
11 //* without fee, provided that the above copyright notice appears in all   *
12 //* copies and that both the copyright notice and this permission notice   *
13 //* appear in the supporting documentation. The authors make no claims     *
14 //* about the suitability of this software for any purpose. It is          *
15 //* provided "as is" without express or implied warranty.                  *
16 //**************************************************************************
17
18 /// @file   AliHLTTPCDataCompressionComponent.cxx
19 /// @author Matthias Richter
20 /// @date   2011-08-08
21 /// @brief  TPC component for data compression
22 ///
23
24 #include "AliHLTTPCDataCompressionComponent.h"
25 #include "AliHLTTPCDefinitions.h"
26 #include "AliHLTTPCDataCompressionDescriptor.h"
27 #include "AliHLTTPCRawClustersDescriptor.h"
28 #include "AliHLTTPCTrackGeometry.h"
29 #include "AliHLTTPCSpacePointContainer.h"
30 #include "AliHLTTPCRawSpacePointContainer.h"
31 #include "AliHLTGlobalBarrelTrack.h"
32 #include "AliHLTComponentBenchmark.h"
33 #include "AliHLTDataDeflaterSimple.h"
34 #include "AliHLTDataDeflaterHuffman.h"
35 #include "AliHLTTPCTransform.h"
36 #include "AliHLTTPCClusterMCData.h"
37 #include "AliHLTTPCClusterTransformation.h"
38 #include "AliHLTErrorGuard.h"
39 #include "AliCDBManager.h"
40 #include "AliCDBPath.h"
41 #include "AliCDBId.h"
42 #include "AliCDBMetaData.h"
43 #include "AliCDBEntry.h"
44 #include "TH1F.h"
45 #include "TFile.h"
46 #include <memory>
47
48 ClassImp(AliHLTTPCDataCompressionComponent)
49
50 AliHLTTPCDataCompressionComponent::AliHLTTPCDataCompressionComponent()
51   : AliHLTProcessor()
52   , fMode(kCompressionModeNone)
53   , fDeflaterMode(kDeflaterModeNone)
54   , fVerificationMode(0)
55   , fMaxDeltaPad(AliHLTTPCDefinitions::GetMaxClusterDeltaPad())
56   , fMaxDeltaTime(AliHLTTPCDefinitions::GetMaxClusterDeltaTime())
57   , fRawInputClusters(NULL)
58   , fInputClusters(NULL)
59   , fTrackGrid(NULL)
60   , fSpacePointGrid(NULL)
61   , fpDataDeflater(NULL)
62   , fHistoCompFactor(NULL)
63   , fHistoResidualPad(NULL)
64   , fHistoResidualTime(NULL)
65   , fHistoClustersOnTracks(NULL)
66   , fHistoClusterRatio(NULL)
67   , fHistoTrackClusterRatio(NULL)
68   , fHistogramFile()
69   , fTrainingTableOutput()
70   , fpBenchmark(NULL)
71   , fpWrittenAssociatedClusterIds(NULL)
72   , fDriftTimeFactorA(1.)
73   , fDriftTimeOffsetA(0.)
74   , fDriftTimeFactorC(1.)
75   , fDriftTimeOffsetC(0.)
76   , fVerbosity(0)
77 {
78 }
79
80 AliHLTTPCDataCompressionComponent::~AliHLTTPCDataCompressionComponent()
81 {
82   /// destructor
83   if (fpWrittenAssociatedClusterIds) delete fpWrittenAssociatedClusterIds;
84 }
85
86
87 const char* AliHLTTPCDataCompressionComponent::GetComponentID()
88 {
89   /// inherited from AliHLTComponent: id of the component
90   return "TPCDataCompressor";
91 }
92
93
94 void AliHLTTPCDataCompressionComponent::GetInputDataTypes( AliHLTComponentDataTypeList& tgtList)
95 {
96   /// inherited from AliHLTComponent: list of data types in the vector reference
97   tgtList.clear();
98   tgtList.push_back(AliHLTTPCDefinitions::RawClustersDataType());
99   tgtList.push_back(AliHLTTPCDefinitions::RawClustersDescriptorDataType());
100   tgtList.push_back(AliHLTTPCDefinitions::ClustersDataType());
101   tgtList.push_back(kAliHLTDataTypeTrack|kAliHLTDataOriginTPC);
102 }
103
104 AliHLTComponentDataType AliHLTTPCDataCompressionComponent::GetOutputDataType()
105 {
106   /// inherited from AliHLTComponent: output data type of the component.
107   return kAliHLTMultipleDataType;
108 }
109
110 int AliHLTTPCDataCompressionComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList)
111 {
112   /// inherited from AliHLTComponent: multiple output data types of the component.
113   tgtList.clear();
114   tgtList.push_back(AliHLTTPCDefinitions::DataCompressionDescriptorDataType());
115   tgtList.push_back(AliHLTTPCDefinitions::RawClustersDataType());
116   tgtList.push_back(AliHLTTPCDefinitions::RemainingClustersCompressedDataType());
117   tgtList.push_back(AliHLTTPCDefinitions::RemainingClusterIdsDataType());
118   tgtList.push_back(AliHLTTPCDefinitions::ClusterTracksCompressedDataType());
119   tgtList.push_back(AliHLTTPCDefinitions::ClusterIdTracksDataType());
120   return tgtList.size();
121 }
122
123 void AliHLTTPCDataCompressionComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
124 {
125   /// inherited from AliHLTComponent: output data size estimator
126   constBase=0;
127   inputMultiplier=1.;  // there should not be more data than input
128   inputMultiplier+=.3; // slightly more data when using the old HWCF data with 20 Byte and raw clusters 22 Byte
129   if (fpWrittenAssociatedClusterIds) inputMultiplier+=.3; // space for optional cluster id array
130 }
131
132 AliHLTComponent* AliHLTTPCDataCompressionComponent::Spawn()
133 {
134   /// inherited from AliHLTComponent: spawn function.
135   return new AliHLTTPCDataCompressionComponent;
136 }
137
138 void AliHLTTPCDataCompressionComponent::GetOCDBObjectDescription(TMap* const targetMap)
139 {
140   /// Get a list of OCDB object needed for the particular component
141   if (!targetMap) return;
142
143   targetMap->Add(new TObjString("HLT/ConfigTPC/TPCDataCompressor"),
144                  new TObjString("component arguments"));
145   if (fDeflaterMode==kDeflaterModeHuffman) {
146     targetMap->Add(new TObjString("HLT/ConfigTPC/TPCDataCompressorHuffmanTables"),
147                    new TObjString("huffman tables for deflater mode 'huffman'"));
148   }
149 }
150
151 int AliHLTTPCDataCompressionComponent::DoEvent( const AliHLTComponentEventData& /*evtData*/, 
152                                                 const AliHLTComponentBlockData* /*inputBlocks*/, 
153                                                 AliHLTComponentTriggerData& /*trigData*/,
154                                                 AliHLTUInt8_t* outputPtr,
155                                                 AliHLTUInt32_t& size,
156                                                 AliHLTComponentBlockDataList& outputBlocks )
157 {
158   /// inherited from AliHLTProcessor: data processing
159   int iResult=0;
160   AliHLTUInt32_t capacity=size;
161   size=0;
162
163   if (!IsDataEvent()) return 0;
164
165   if (!fRawInputClusters) {
166     return -ENODEV;
167   }
168
169   if (GetBenchmarkInstance()) {
170     GetBenchmarkInstance()->StartNewEvent();
171     GetBenchmarkInstance()->Start(0);
172   }
173
174   // Process an event
175
176   // Loop over all input blocks in the event
177   bool bHaveMC=(GetFirstInputBlock(AliHLTTPCDefinitions::fgkAliHLTDataTypeClusterMCInfo | kAliHLTDataOriginTPC))!=NULL;
178   if ((bHaveMC || fVerificationMode>0) && fpWrittenAssociatedClusterIds==NULL) {
179     fpWrittenAssociatedClusterIds=new vector<AliHLTUInt32_t>;
180   }
181
182   const AliHLTComponentBlockData* pDesc=NULL;
183
184   AliHLTUInt8_t minSlice=0xFF, maxSlice=0xFF, minPatch=0xFF, maxPatch=0xFF;
185   AliHLTUInt32_t inputRawClusterSize=0;
186   AliHLTUInt32_t outputDataSize=0;
187   int allClusters=0;
188   int associatedClusters=0;
189   float bz=GetBz();
190
191   /// input track array
192   vector<AliHLTGlobalBarrelTrack> inputTrackArray;
193
194   if (GetBenchmarkInstance()) {
195     GetBenchmarkInstance()->Start(2);
196   }
197
198   bool isInputPresent = kFALSE;
199
200   // transformed clusters
201   // the transformed clusters have not been used yet
202   if (false) { // FIXME: condition to be adjusted
203     for (pDesc=GetFirstInputBlock(AliHLTTPCDefinitions::fgkClustersDataType);
204          pDesc!=NULL; pDesc=GetNextInputBlock()) {
205       if (GetBenchmarkInstance()) {
206         GetBenchmarkInstance()->AddInput(pDesc->fSize);
207       }
208       isInputPresent = kTRUE;
209       AliHLTUInt8_t slice = 0;
210       AliHLTUInt8_t patch = 0;
211       slice = AliHLTTPCDefinitions::GetMinSliceNr( pDesc->fSpecification );
212       patch = AliHLTTPCDefinitions::GetMinPatchNr( pDesc->fSpecification );
213       if ( minSlice==0xFF || slice<minSlice )   minSlice = slice;
214       if ( maxSlice==0xFF || slice>maxSlice )   maxSlice = slice;
215       if ( minPatch==0xFF || patch<minPatch )   minPatch = patch;
216       if ( maxPatch==0xFF || patch>maxPatch )   maxPatch = patch;
217       if (fInputClusters) {
218         fInputClusters->AddInputBlock(pDesc);
219       }
220     }
221     if (GetBenchmarkInstance()) {
222       GetBenchmarkInstance()->Stop(2);
223       GetBenchmarkInstance()->Start(3);
224     }
225   }
226
227   vector<int> trackindexmap; // stores index for every track id
228
229   // track data input
230   if (fMode==kCompressionModeV1TrackModel || fMode==kCompressionModeV2TrackModel) {
231     for (pDesc=GetFirstInputBlock(kAliHLTDataTypeTrack|kAliHLTDataOriginTPC);
232          pDesc!=NULL; pDesc=GetNextInputBlock()) {
233       if (GetBenchmarkInstance()) {
234         GetBenchmarkInstance()->AddInput(pDesc->fSize);
235       }
236       AliHLTUInt8_t slice = 0;
237       AliHLTUInt8_t patch = 0;
238       slice = AliHLTTPCDefinitions::GetMinSliceNr( pDesc->fSpecification );
239       patch = AliHLTTPCDefinitions::GetMinPatchNr( pDesc->fSpecification );
240       if ( minSlice==0xFF || slice<minSlice )   minSlice = slice;
241       if ( maxSlice==0xFF || slice>maxSlice )   maxSlice = slice;
242       if ( minPatch==0xFF || patch<minPatch )   minPatch = patch;
243       if ( maxPatch==0xFF || patch>maxPatch )   maxPatch = patch;
244       const AliHLTTracksData* pTracks=reinterpret_cast<const AliHLTTracksData*>(pDesc->fPtr);
245       if ((iResult=AliHLTGlobalBarrelTrack::ConvertTrackDataArray(pTracks, pDesc->fSize, inputTrackArray))<0) {
246         return iResult;
247       }
248       trackindexmap.resize(inputTrackArray.size(), -1);
249     }
250   }
251
252   if (GetBenchmarkInstance()) {
253     GetBenchmarkInstance()->Stop(3);
254     GetBenchmarkInstance()->Start(4);
255   }
256
257   // processing
258   int trackindex=0;
259   for (vector<AliHLTGlobalBarrelTrack>::iterator track=inputTrackArray.begin();
260        track!=inputTrackArray.end();
261        track++, trackindex++) {
262     int trackID=track->GetID();
263     if (trackID<0) {
264       // FIXME: error guard
265       HLTError("invalid track ID");
266       continue;
267     }
268     if (trackID>=(int)trackindexmap.size())
269       trackindexmap.resize(trackID+1, -1);
270     trackindexmap[trackID]=trackindex;
271
272     if (fVerbosity>0) {
273       UInt_t nofPoints=track->GetNumberOfPoints();
274       const UInt_t* points=track->GetPoints();
275       for (unsigned i=0; i<nofPoints; i++) {
276         int slice=AliHLTTPCSpacePointData::GetSlice(points[i]);
277         int partition=AliHLTTPCSpacePointData::GetPatch(points[i]);
278         int number=AliHLTTPCSpacePointData::GetNumber(points[i]);
279         HLTInfo("track %d point %d id 0x%08x slice %d partition %d number %d", track->GetID(), i, points[i], slice, partition, number);
280       }
281     }
282
283     AliHLTTPCTrackGeometry* trackpoints=new AliHLTTPCTrackGeometry;
284     if (!trackpoints) continue;
285     HLTDebug("track %d id %d:", trackindex, trackID);
286
287     // in order to avoid rounding errors the track points are
288     // calculated in exactly the same way as in the decoding
289     // Thats why the track instance can not be used directly
290     // but a new instance is created from the values in the
291     // storage format.
292     // think about moving that to some common code used by
293     // both compression and decoding
294     AliHLTExternalTrackParam param;
295     memset(&param, 0, sizeof(param));    
296     float alpha=track->GetAlpha();
297     while (alpha<0.) alpha+=TMath::TwoPi();
298     while (alpha>TMath::TwoPi()) alpha-=TMath::TwoPi();
299     AliHLTUInt8_t tSlice=AliHLTUInt8_t(9*alpha/TMath::Pi());
300     param.fAlpha   =( tSlice + 0.5 ) * TMath::Pi() / 9.0;
301     if (param.fAlpha>TMath::TwoPi()) param.fAlpha-=TMath::TwoPi();
302     param.fX       = track->GetX();
303     param.fY       = track->GetY();
304     param.fZ       = track->GetZ();
305     param.fSinPsi  = track->GetSnp();
306     param.fTgl     = track->GetTgl();
307     param.fq1Pt    = track->GetSigned1Pt();
308     AliHLTGlobalBarrelTrack ctrack(param);
309     ctrack.CalculateHelixParams(bz);
310     trackpoints->InitDriftTimeTransformation(fDriftTimeFactorA, fDriftTimeOffsetA, fDriftTimeFactorC, fDriftTimeOffsetC);
311     trackpoints->SetTrackId(trackID);
312     trackpoints->CalculateTrackPoints(ctrack);
313     trackpoints->RegisterTrackPoints(fTrackGrid);
314     track->SetTrackGeometry(trackpoints);
315   }
316
317   for (vector<AliHLTGlobalBarrelTrack>::const_iterator track=inputTrackArray.begin();
318        track!=inputTrackArray.end();
319        track++) {
320     AliHLTTrackGeometry* trackpoints=track->GetTrackGeometry();
321     if (!trackpoints) continue;
322     trackpoints->FillTrackPoints(fTrackGrid);
323   }
324   if (fVerbosity>0) {
325     fTrackGrid->Print();
326   }
327
328   if (GetBenchmarkInstance()) {
329     GetBenchmarkInstance()->Stop(4);
330     GetBenchmarkInstance()->Start(5);
331   }
332
333   // loop over raw cluster blocks, assign to tracks and write
334   // unassigned clusters
335   for (pDesc=GetFirstInputBlock(AliHLTTPCDefinitions::fgkRawClustersDataType);
336        pDesc!=NULL; pDesc=GetNextInputBlock()) {    
337     if (GetBenchmarkInstance()) {
338       GetBenchmarkInstance()->Start(1);
339       GetBenchmarkInstance()->AddInput(pDesc->fSize);
340     }
341     isInputPresent = kTRUE;
342     AliHLTUInt8_t slice = 0;
343     AliHLTUInt8_t patch = 0;
344     slice = AliHLTTPCDefinitions::GetMinSliceNr( pDesc->fSpecification );
345     patch = AliHLTTPCDefinitions::GetMinPatchNr( pDesc->fSpecification );
346     if ( minSlice==0xFF || slice<minSlice )     minSlice = slice;
347     if ( maxSlice==0xFF || slice>maxSlice )     maxSlice = slice;
348     if ( minPatch==0xFF || patch<minPatch )     minPatch = patch;
349     if ( maxPatch==0xFF || patch>maxPatch )     maxPatch = patch;
350     inputRawClusterSize+=pDesc->fSize;
351
352     // add the data and populate the index grid
353     fRawInputClusters->AddInputBlock(pDesc);
354     fRawInputClusters->PopulateAccessGrid(fSpacePointGrid, pDesc->fSpecification);
355     if (fVerbosity>0 && fSpacePointGrid->GetNumberOfSpacePoints()>0) {
356       HLTInfo("index grid slice %d partition %d", slice, patch);
357       fSpacePointGrid->Print();
358       for (AliHLTSpacePointContainer::AliHLTSpacePointPropertyGrid::iterator& cl=fSpacePointGrid->begin();
359            cl!=fSpacePointGrid->end(); cl++) {
360         AliHLTUInt32_t id=cl.Data().fId;
361         float row=fRawInputClusters->GetX(id);
362         float pad=fRawInputClusters->GetY(id);
363         float time=fRawInputClusters->GetZ(id);
364         HLTInfo("    cluster id 0x%08x: row %f  pad %f  time %f", id, row, pad, time);
365       }
366     }
367     if (GetBenchmarkInstance()) {
368       GetBenchmarkInstance()->Stop(1);
369       GetBenchmarkInstance()->Start(4);
370     }
371
372     // process the clusters per padrow and check the track grid
373     // for tracks crossing that particular padrow
374     if (GetBenchmarkInstance()) {
375       GetBenchmarkInstance()->Stop(4);
376       GetBenchmarkInstance()->Start(5);
377     }
378     allClusters+=fSpacePointGrid->GetNumberOfSpacePoints();
379     iResult=ProcessTrackClusters(&inputTrackArray[0], inputTrackArray.size(), fTrackGrid, trackindexmap, fSpacePointGrid, fRawInputClusters, slice, patch);
380     int assignedInThisPartition=0;
381     if (iResult>=0) {
382       assignedInThisPartition=iResult;
383       associatedClusters+=iResult;
384     }
385     iResult=ProcessRemainingClusters(&inputTrackArray[0], inputTrackArray.size(), fTrackGrid, trackindexmap, fSpacePointGrid, fRawInputClusters, slice, patch);
386     if (iResult>=0) {
387       if (fSpacePointGrid->GetNumberOfSpacePoints()>0) {
388         if (fVerbosity>0) HLTInfo("associated %d (%d) of %d clusters in slice %d partition %d", iResult+assignedInThisPartition, assignedInThisPartition, fSpacePointGrid->GetNumberOfSpacePoints(), slice, patch);
389       }
390       associatedClusters+=iResult;
391     }
392
393     // write all remaining clusters not yet assigned to tracks
394     // the index grid is used to write sorted in padrow
395     // FIXME: decoder index instead of data specification to be used
396     // use an external access grid to reduce allocated memory
397     // set to NULL after writing the clusters
398     const char* writeoptions="";
399     if (fpWrittenAssociatedClusterIds) {
400       writeoptions="write-cluster-ids";
401     }
402     fRawInputClusters->SetSpacePointPropertyGrid(pDesc->fSpecification, fSpacePointGrid);
403     iResult=fRawInputClusters->Write(outputPtr+size, capacity-size, outputBlocks, fpDataDeflater, writeoptions);
404     fRawInputClusters->SetSpacePointPropertyGrid(pDesc->fSpecification, NULL);
405     if (iResult>=0) {
406       size+=iResult;
407       outputDataSize+=iResult;
408       // the size of the optional cluster id array must be subtracted
409       if (fpWrittenAssociatedClusterIds && outputBlocks.size()>0 &&
410           outputBlocks.back().fDataType==AliHLTTPCDefinitions::RemainingClusterIdsDataType()) {
411         outputDataSize-=outputBlocks.back().fSize;
412       }
413       if (GetBenchmarkInstance()) GetBenchmarkInstance()->AddOutput(iResult);
414     }
415     if (GetBenchmarkInstance()) {
416       GetBenchmarkInstance()->Stop(5);
417     }
418
419     fSpacePointGrid->Clear();
420   }
421   if (fHistoClusterRatio && allClusters>0) {
422     if (fVerbosity>0) HLTInfo("associated %d of %d clusters to tracks", associatedClusters, allClusters);
423     float ratio=associatedClusters; ratio/=allClusters;
424     fHistoClusterRatio->Fill(ratio);
425   }
426
427   // output of track model clusters
428   if (iResult>=0) do {
429     AliHLTUInt32_t tracksBufferOffset=sizeof(AliHLTTPCTrackModelBlock);
430     if (capacity-size<tracksBufferOffset) {
431       iResult=-ENOSPC;
432       break;
433     }
434     if (fpWrittenAssociatedClusterIds) fpWrittenAssociatedClusterIds->clear();
435     AliHLTTPCTrackModelBlock* trackModelBlock=reinterpret_cast<AliHLTTPCTrackModelBlock*>(outputPtr+size);
436     trackModelBlock->fVersion=1;
437     trackModelBlock->fDeflaterMode=fpDataDeflater?fpDataDeflater->GetDeflaterVersion():0;
438     trackModelBlock->fTrackCount=inputTrackArray.size();
439     trackModelBlock->fClusterCount=0;
440     trackModelBlock->fGlobalParameterCnt=5;
441     tracksBufferOffset+=trackModelBlock->fGlobalParameterCnt*sizeof(trackModelBlock->fGlobalParameters);
442     if (capacity-size<tracksBufferOffset) {
443       iResult=-ENOSPC;
444       break;
445     }
446
447     AliHLTUInt32_t parameterIndex=0;
448     trackModelBlock->fGlobalParameters[parameterIndex++]=bz;
449     trackModelBlock->fGlobalParameters[parameterIndex++]=fDriftTimeFactorA;
450     trackModelBlock->fGlobalParameters[parameterIndex++]=fDriftTimeOffsetA;
451     trackModelBlock->fGlobalParameters[parameterIndex++]=fDriftTimeFactorC;
452     trackModelBlock->fGlobalParameters[parameterIndex++]=fDriftTimeOffsetC;
453     if (parameterIndex!=trackModelBlock->fGlobalParameterCnt) {
454       HLTError("internal error, size of parameter array has changed without providing all values");
455       iResult=-EFAULT;
456       break;
457     }
458
459     if (trackindexmap.size()>0) {// condition for track model compression
460     iResult=WriteTrackClusters(inputTrackArray, fRawInputClusters, fpDataDeflater, outputPtr+size+tracksBufferOffset, capacity-size-tracksBufferOffset);
461     if (iResult>=0) {
462       AliHLTComponent_BlockData bd;
463       FillBlockData(bd);
464       bd.fOffset        = size;
465       bd.fSize          = tracksBufferOffset+iResult;
466       bd.fDataType      = AliHLTTPCDefinitions::ClusterTracksCompressedDataType();
467       bd.fSpecification = AliHLTTPCDefinitions::EncodeDataSpecification(minSlice, maxSlice, minPatch, maxPatch);
468       outputBlocks.push_back(bd);
469       size += bd.fSize;
470       outputDataSize+=bd.fSize;
471       HLTBenchmark("track data block of %d tracks: size %d", inputTrackArray.size(), bd.fSize);
472
473       if (fpWrittenAssociatedClusterIds && fpWrittenAssociatedClusterIds->size()>0) {
474         AliHLTComponent::FillBlockData(bd);
475         bd.fOffset        = size;
476         bd.fSize        = fpWrittenAssociatedClusterIds->size()*sizeof(vector<AliHLTUInt32_t>::value_type);
477         if (capacity-size>bd.fSize) {
478           memcpy(outputPtr+bd.fOffset, &(*fpWrittenAssociatedClusterIds)[0], bd.fSize);
479           bd.fDataType    = AliHLTTPCDefinitions::ClusterIdTracksDataType();
480           bd.fSpecification = AliHLTTPCDefinitions::EncodeDataSpecification(minSlice, maxSlice, minPatch, maxPatch);
481           outputBlocks.push_back(bd);    
482           size += bd.fSize;
483         } else {
484           iResult=-ENOSPC;
485         }
486         
487         fpWrittenAssociatedClusterIds->clear();
488       }
489     }
490     }
491   } while (0);
492
493   fRawInputClusters->Clear();
494
495   // Write header block 
496   
497   if( isInputPresent ){
498     pDesc=GetFirstInputBlock(AliHLTTPCDefinitions::RawClustersDescriptorDataType() );
499     if( pDesc ){
500       const AliHLTTPCRawClustersDescriptor &clDesc = *reinterpret_cast<const AliHLTTPCRawClustersDescriptor*>(pDesc->fPtr);
501       if( !clDesc.CheckSize( pDesc->fSize ) ){
502         HLTError("Corrupted cluster descriptor");
503       }
504      
505       AliHLTComponent_BlockData bd;
506       FillBlockData(bd);
507       bd.fOffset        = size;
508       bd.fSize          = sizeof(AliHLTTPCDataCompressionDescriptor);
509       bd.fDataType      = AliHLTTPCDefinitions::DataCompressionDescriptorDataType();
510       if( capacity < size + bd.fSize ){
511         iResult = -ENOSPC;
512         return iResult;
513       }
514
515       AliHLTTPCDataCompressionDescriptor compDesc;
516       compDesc.SetMergedClustersFlag( clDesc.GetMergedClustersFlag() );
517       *(AliHLTTPCDataCompressionDescriptor*)(outputPtr + bd.fOffset ) = compDesc; 
518
519       outputBlocks.push_back(bd);
520       size += bd.fSize;
521       outputDataSize+=bd.fSize;
522       //HLTBenchmark("header data block of size %d", bd.fSize);
523     }
524   }
525
526   float compressionFactor=(float)inputRawClusterSize;
527   if ((outputDataSize)>0) compressionFactor/=outputDataSize;
528   else compressionFactor=0.;
529   if (fHistoCompFactor) fHistoCompFactor->Fill(compressionFactor);
530
531   if (GetBenchmarkInstance() && allClusters>0) {
532     GetBenchmarkInstance()->Stop(0);
533     if (fDeflaterMode!=kDeflaterModeHuffmanTrainer) {
534       HLTBenchmark("%s - compression factor %.2f", GetBenchmarkInstance()->GetStatistics(), compressionFactor);
535     } else {
536       HLTBenchmark("%s", GetBenchmarkInstance()->GetStatistics());
537     }
538   }
539
540   if (fInputClusters) {
541     fInputClusters->Clear();
542   }
543   if (fRawInputClusters) {
544     fRawInputClusters->Clear();
545   }
546   if (fTrackGrid) {
547     fTrackGrid->Clear();
548   }
549
550   // forward MC labels
551   for (pDesc=GetFirstInputBlock(AliHLTTPCDefinitions::fgkAliHLTDataTypeClusterMCInfo | kAliHLTDataOriginTPC);
552        pDesc!=NULL; pDesc=GetNextInputBlock()) {
553     outputBlocks.push_back(*pDesc);
554   }
555
556   return iResult;
557 }
558
559 int AliHLTTPCDataCompressionComponent::ProcessTrackClusters(AliHLTGlobalBarrelTrack* pTracks, unsigned nofTracks,
560                                                             AliHLTTrackGeometry::AliHLTTrackGrid* pTrackIndex,
561                                                             const vector<int>& trackIndexMap,
562                                                             AliHLTSpacePointContainer::AliHLTSpacePointPropertyGrid* pClusterIndex,
563                                                             AliHLTSpacePointContainer* pClusters,
564                                                             int slice, int partition) const
565 {
566   // process to assigned track clusters
567   int iResult=0;
568   int assignedClusters=0;
569   if (!pTracks || nofTracks==0) return 0;
570
571   vector<int> processedTracks(nofTracks, -1);
572   for (AliHLTTrackGeometry::AliHLTTrackGrid::iterator& trackId=pTrackIndex->begin(slice, partition, -1);
573        trackId!=pTrackIndex->end(); trackId++) {
574     if (trackId.Data()>=trackIndexMap.size()) {
575       HLTError("can not find track id %d in index map of size %d", trackId.Data(), trackIndexMap.size());
576       continue;
577     }
578     int trackindex=trackIndexMap[trackId.Data()];
579     if (trackindex<0 || trackindex>=(int)nofTracks) {
580       HLTError("invalid index %d found for track id %d", trackindex, trackId.Data());
581       continue;
582     }
583     if (processedTracks[trackindex]>0) continue;
584     processedTracks[trackindex]=1;
585     AliHLTGlobalBarrelTrack& track=pTracks[trackindex];
586     if (!track.GetTrackGeometry()) {
587       HLTError("can not find track geometry for track %d", trackId.Data());
588       continue;
589     }
590     AliHLTTPCTrackGeometry* pTrackPoints=dynamic_cast<AliHLTTPCTrackGeometry*>(track.GetTrackGeometry());
591     if (!pTrackPoints) {
592       HLTError("invalid track geometry type for track %d, expecting AliHLTTPCTrackGeometry", trackId.Data());
593       continue; 
594     }
595
596     UInt_t nofTrackPoints=track.GetNumberOfPoints();
597     const UInt_t* trackPoints=track.GetPoints();
598     for (unsigned i=0; i<nofTrackPoints; i++) {
599       const AliHLTUInt32_t& clusterId=trackPoints[i];
600       if (AliHLTTPCSpacePointData::GetSlice(clusterId)!=(unsigned)slice ||
601           AliHLTTPCSpacePointData::GetPatch(clusterId)!=(unsigned)partition) {
602         // not in the current partition;
603         continue;
604       }
605
606       int clusterrow=(int)pClusters->GetX(clusterId);
607       AliHLTUInt32_t pointId=AliHLTTPCSpacePointData::GetID(slice, partition, clusterrow);
608       AliHLTTrackGeometry::AliHLTTrackPoint* point=pTrackPoints->GetRawTrackPoint(pointId);
609       if (!point) {
610         //HLTError("can not find track point slice %d partition %d padrow %d (0x%08x) of track %d", slice, partition, clusterrow, pointId, trackId.Data());
611         continue;
612       }
613       float pad=point->GetU();
614       float time=point->GetV();
615
616       iResult=FindCellClusters(trackId.Data(), clusterrow, pad, time, pClusterIndex, pClusters, point, clusterId);
617       if (iResult>0) assignedClusters+=iResult;
618     }
619   }
620   return assignedClusters;
621 }
622
623 int AliHLTTPCDataCompressionComponent::ProcessRemainingClusters(AliHLTGlobalBarrelTrack* pTracks, unsigned nofTracks,
624                                                                 AliHLTTrackGeometry::AliHLTTrackGrid* pTrackIndex,
625                                                                 const vector<int>& trackIndexMap,
626                                                                 AliHLTSpacePointContainer::AliHLTSpacePointPropertyGrid* pClusterIndex,
627                                                                 AliHLTSpacePointContainer* pClusters,
628                                                                 int slice, int partition) const
629 {
630   // assign remaining clusters to tracks
631   int iResult=0;
632   int associatedClusters=0;
633   if (!pTracks || nofTracks==0) return 0;
634
635   for (int padrow=0; padrow<AliHLTTPCTransform::GetNRows(partition); padrow++) {
636     for (AliHLTTrackGeometry::AliHLTTrackGrid::iterator& trackId=pTrackIndex->begin(slice, partition, padrow);
637          trackId!=pTrackIndex->end(); trackId++) {
638       if (trackId.Data()>=trackIndexMap.size()) {
639         HLTError("can not find track id %d in index map of size %d", trackId.Data(), trackIndexMap.size());
640         continue;
641       }
642       int trackindex=trackIndexMap[trackId.Data()];
643       if (trackindex<0 || trackindex>=(int)nofTracks) {
644         HLTError("invalid index %d found for track id %d", trackindex, trackId.Data());
645         continue;
646       }
647       AliHLTGlobalBarrelTrack& track=pTracks[trackindex];
648       if (!track.GetTrackGeometry()) {
649         HLTError("can not find track geometry for track %d", trackId.Data());
650         continue;
651       }
652       AliHLTTPCTrackGeometry* pTrackPoints=dynamic_cast<AliHLTTPCTrackGeometry*>(track.GetTrackGeometry());
653       if (!pTrackPoints) {
654         HLTError("invalid track geometry type for track %d, expecting AliHLTTPCTrackGeometry", trackId.Data());
655         continue;       
656       }
657       AliHLTUInt32_t pointId=AliHLTTPCSpacePointData::GetID(slice, partition, padrow);
658       AliHLTTrackGeometry::AliHLTTrackPoint* point=pTrackPoints->GetRawTrackPoint(pointId);
659       if (!point) {
660         //HLTError("can not find track point slice %d partition %d padrow %d (0x%08x) of track %d", slice, partition, padrow, pointId, trackId.Data());
661         continue;
662       }
663       float pad=point->GetU();
664       float time=point->GetV();
665
666       iResult=FindCellClusters(trackId.Data(), padrow, pad, time, pClusterIndex, pClusters, point);
667       if (iResult>0) associatedClusters+=iResult;
668       if (fVerbosity>0) {
669         HLTInfo("trackpoint track %d slice %d partition %d padrow %d: %.3f \t%.3f - associated %d", track.GetID(), slice, partition, padrow, pad, time, iResult);
670       }
671     }
672   }
673   if (iResult<0) return iResult;
674   return associatedClusters;
675 }
676
677 int AliHLTTPCDataCompressionComponent::FindCellClusters(int trackId, int padrow, float pad, float time,
678                                                         AliHLTSpacePointContainer::AliHLTSpacePointPropertyGrid* pClusterIndex,
679                                                         AliHLTSpacePointContainer* pClusters,
680                                                         AliHLTTrackGeometry::AliHLTTrackPoint* pTrackPoint,
681                                                         AliHLTUInt32_t clusterId) const
682 {
683   // check index cell for entries and assign to track
684   int count=0;
685   // search a 4x4 matrix out of the 9x9 matrix around the cell addressed by
686   // pad and time
687   int rowindex=pClusterIndex->GetXIndex((float)padrow);
688   int padstartindex=pClusterIndex->GetYIndex(pad);
689   int timestartindex=pClusterIndex->GetZIndex(time);
690   int cellindex=pClusterIndex->Index(rowindex, padstartindex, timestartindex);
691   float centerpad=pClusterIndex->GetCenterY(cellindex);
692   float centertime=pClusterIndex->GetCenterZ(cellindex);
693   if ((TMath::Abs(centerpad-pad)>fMaxDeltaPad && pad>0.) ||
694       (TMath::Abs(centertime-time)>fMaxDeltaTime && time>0.)) {
695     ALIHLTERRORGUARD(20, "invalid pad center calculation, please check dimensions if dimensions of index grid match the maximum possible deviation");
696   }
697
698   int paddirection=1;
699   int timedirection=1;
700   if (centerpad>pad) paddirection=-1;
701   if (centertime>time) timedirection=-1;
702   for (int padcount=0, padindex=padstartindex; padcount<2; padcount++, padindex+=paddirection) {
703     if (padindex<0) continue;
704     if (padindex>=pClusterIndex->GetDimensionY()) break;
705     for (int timecount=0, timeindex=timestartindex; timecount<2; timecount++, timeindex+=timedirection) {
706       if (timeindex<0) continue;
707       if (timeindex>=pClusterIndex->GetDimensionZ()) break;
708       cellindex=pClusterIndex->Index(rowindex, padindex, timeindex);
709       float cellpad=pClusterIndex->GetCenterY(cellindex);
710       float celltime=pClusterIndex->GetCenterZ(cellindex);
711       for (AliHLTSpacePointContainer::AliHLTSpacePointPropertyGrid::iterator& cl=pClusterIndex->begin((float)padrow, cellpad, celltime);
712            cl!=pClusterIndex->end(); cl++) {
713         if (cl.Data().fTrackId>=0) continue;
714         if (clusterId!=(~(AliHLTUInt32_t)0) && clusterId!=cl.Data().fId) continue;
715         if (TMath::Abs(padrow-pClusters->GetX(cl.Data().fId))>=1.) {
716           HLTError("cluster 0x%08x: mismatch on padrow: trackpoint %d  cluster %f", cl.Data().fId, padrow, pClusters->GetX(cl.Data().fId));
717           continue;
718         }
719         float clusterpad=pClusters->GetY(cl.Data().fId);
720         float clustertime=pClusters->GetZ(cl.Data().fId);
721         if (TMath::Abs(clusterpad-pad)<fMaxDeltaPad &&
722             TMath::Abs(clustertime-time)<fMaxDeltaTime) {
723           // add this cluster to the track point and mark in the index grid
724           cl.Data().fTrackId=trackId;
725           pTrackPoint->AddAssociatedSpacePoint(cl.Data().fId, clusterpad-pad, clustertime-time);
726           count++;
727         }
728         if (clusterId!=(~(AliHLTUInt32_t)0)) break;
729       }
730     }
731   }
732   return count;
733 }
734
735 int AliHLTTPCDataCompressionComponent::WriteTrackClusters(const vector<AliHLTGlobalBarrelTrack>& tracks,
736                                                           AliHLTSpacePointContainer* pSpacePoints,
737                                                           AliHLTDataDeflater* pDeflater,
738                                                           AliHLTUInt8_t* outputPtr,
739                                                           AliHLTUInt32_t capacity) const
740 {
741   // write the track data block including all associated clusters
742   AliHLTUInt32_t size=0;
743   for (vector<AliHLTGlobalBarrelTrack>::const_iterator track=tracks.begin();
744        track!=tracks.end();
745        track++) {
746     if (!track->GetTrackGeometry()) {
747       HLTError("can not find track geometry for track %d", track->GetID());
748       return -EBADF;
749     }
750     AliHLTTPCTrackGeometry* pTrackPoints=dynamic_cast<AliHLTTPCTrackGeometry*>(track->GetTrackGeometry());
751     if (!pTrackPoints) {
752       HLTError("invalid track geometry type for track %d, expecting AliHLTTPCTrackGeometry", track->GetID());
753       return -EBADF;
754     }
755
756     int result=pTrackPoints->Write(*track, pSpacePoints, pDeflater, outputPtr+size, capacity-size, fpWrittenAssociatedClusterIds);
757     if (result<0) return result;
758     size+=result;
759
760     UInt_t nofTrackPoints=track->GetNumberOfPoints();
761     const UInt_t* trackPoints=track->GetPoints();
762
763     int assignedPoints=0;
764     int assignedTrackPoints=0;
765     const vector<AliHLTTrackGeometry::AliHLTTrackPoint>& rawPoints=pTrackPoints->GetRawPoints();
766     for (vector<AliHLTTrackGeometry::AliHLTTrackPoint>::const_iterator point=rawPoints.begin();
767          point!=rawPoints.end(); point++) {
768       const vector<AliHLTTrackGeometry::AliHLTTrackSpacepoint>& spacePoints=point->GetSpacepoints();
769       for (vector<AliHLTTrackGeometry::AliHLTTrackSpacepoint>::const_iterator spacePoint=spacePoints.begin();
770            spacePoint!=spacePoints.end(); spacePoint++) {
771         float dpad=spacePoint->GetResidual(0);
772         float dtime=spacePoint->GetResidual(1);
773         if (dpad>-1000 && dtime>-1000 && fHistoResidualPad && fHistoResidualTime) {
774           fHistoResidualPad->Fill(dpad);
775           fHistoResidualTime->Fill(dtime);
776         }
777         assignedPoints++;
778         for (unsigned i=0; i<nofTrackPoints; i++) {
779           if (trackPoints[i]==spacePoint->fId) {
780             assignedTrackPoints++;
781             break;
782           }
783         }
784       }
785     }
786     if (fHistoClustersOnTracks) {
787       fHistoClustersOnTracks->Fill(assignedPoints);
788     }
789     if (fHistoTrackClusterRatio && nofTrackPoints>0) {
790       float ratio=assignedTrackPoints; ratio/=nofTrackPoints;
791       fHistoTrackClusterRatio->Fill(ratio);
792     }
793   }
794   return size;
795 }
796
797 int AliHLTTPCDataCompressionComponent::DoInit( int argc, const char** argv )
798 {
799   /// inherited from AliHLTComponent: component initialisation and argument scan.
800   int iResult=0;
801
802   // component configuration
803   //Stage 1: default initialization.
804   //Default values.
805
806   //Stage 2: OCDB.
807   TString cdbPath("HLT/ConfigTPC/");
808   cdbPath += GetComponentID();
809   //
810   iResult = ConfigureFromCDBTObjString(cdbPath);
811   if (iResult < 0) 
812     return iResult;
813
814   //Stage 3: command line arguments.
815   int mode=fMode; // just a backup for the info message below
816   int deflaterMode=fDeflaterMode;
817   if (argc && (iResult = ConfigureFromArgumentString(argc, argv)) < 0)
818     return iResult;
819
820   if (mode!=fMode || deflaterMode!=fDeflaterMode) {
821     HLTInfo("configured from command line: mode %d, deflater mode %d", fMode, fDeflaterMode);
822   }
823
824   std::auto_ptr<AliHLTComponentBenchmark> benchmark(new AliHLTComponentBenchmark);
825   if (benchmark.get()) {
826     benchmark->SetTimer(0,"total");
827     benchmark->SetTimer(1,"rawclusterinput");
828     benchmark->SetTimer(2,"clusterinput");
829     benchmark->SetTimer(3,"trackinput");
830     benchmark->SetTimer(4,"processing");
831     benchmark->SetTimer(5,"output");
832   } else {
833     return -ENOMEM;
834   }
835
836   unsigned spacePointContainerMode=0;
837   if (fMode==kCompressionModeV1TrackModel || fMode==kCompressionModeV2TrackModel) {
838     // initialize map data for cluster access in the track association loop
839     spacePointContainerMode|=AliHLTTPCRawSpacePointContainer::kModeCreateMap;
840   }
841   if (fMode==kCompressionModeV2 || fMode==kCompressionModeV2TrackModel) {
842     // optimized storage format: differential pad and time storage
843     spacePointContainerMode|=AliHLTTPCRawSpacePointContainer::kModeDifferentialPadTime;
844   }
845   std::auto_ptr<AliHLTTPCRawSpacePointContainer> rawInputClusters(new AliHLTTPCRawSpacePointContainer(spacePointContainerMode));
846   std::auto_ptr<AliHLTTPCSpacePointContainer> inputClusters(new AliHLTTPCSpacePointContainer);
847
848   std::auto_ptr<TH1F> histoCompFactor(new TH1F("CompressionFactor",
849                                                "HLT TPC data compression factor",
850                                                100, 0., 10.));
851   std::auto_ptr<TH1F> histoResidualPad(new TH1F("PadResidual",
852                                                 "HLT TPC pad residual",
853                                                 100, -fMaxDeltaPad, fMaxDeltaPad));
854   std::auto_ptr<TH1F> histoResidualTime(new TH1F("TimeResidual",
855                                                  "HLT TPC time residual",
856                                                  100, -fMaxDeltaTime, fMaxDeltaTime));
857   std::auto_ptr<TH1F> histoClustersOnTracks(new TH1F("ClustersOnTracks",
858                                                  "Clusters in track model compression",
859                                                  200, 0., 600));
860   std::auto_ptr<TH1F> histoClusterRatio(new TH1F("ClusterRatio",
861                                                  "Fraction of clusters in track model compression",
862                                                  100, 0., 1.));
863   std::auto_ptr<TH1F> histoTrackClusterRatio(new TH1F("UsedTrackClusters",
864                                                  "Fraction of track clusters in track model compression",
865                                                  100, 0., 1.));
866
867   // track grid: 36 slices, each 6 partitions with max 33 rows
868   fTrackGrid=new AliHLTTrackGeometry::AliHLTTrackGrid(36, 1, 6, 1, 33, 1, 20000);
869   fSpacePointGrid=AliHLTTPCRawSpacePointContainer::AllocateIndexGrid();
870
871   if (!rawInputClusters.get() ||
872       !inputClusters.get() ||
873       !fTrackGrid ||
874       !fSpacePointGrid) {
875     if (fTrackGrid) delete fTrackGrid; fTrackGrid=NULL;
876     if (fSpacePointGrid) delete fSpacePointGrid; fSpacePointGrid=NULL;
877     return -ENOMEM;
878   }
879
880   if (fDeflaterMode>0 && (iResult=InitDeflater(fDeflaterMode))<0)
881     return iResult;
882
883   fpBenchmark=benchmark.release();
884   fRawInputClusters=rawInputClusters.release();
885   fInputClusters=inputClusters.release();
886
887   // initialize the histograms if stored at the end
888   // condition might be extended
889   if (!fHistogramFile.IsNull()) {
890     fHistoCompFactor=histoCompFactor.release();
891     fHistoResidualPad=histoResidualPad.release();
892     fHistoResidualTime=histoResidualTime.release();
893     fHistoClustersOnTracks=histoClustersOnTracks.release();
894     fHistoClusterRatio=histoClusterRatio.release();
895     fHistoTrackClusterRatio=histoTrackClusterRatio.release();
896   }
897
898   if (iResult>=0 && (iResult=InitDriftTimeTransformation())<0) return iResult;
899
900   return iResult;
901 }
902
903 int AliHLTTPCDataCompressionComponent::InitDeflater(int mode)
904 {
905   /// init the data deflater
906   int iResult=0;
907   if (mode==kDeflaterModeHuffman || mode==kDeflaterModeHuffmanTrainer) {
908     // huffman deflater
909     std::auto_ptr<AliHLTDataDeflaterHuffman> deflater(new AliHLTDataDeflaterHuffman(mode==kDeflaterModeHuffmanTrainer));
910     if (!deflater.get()) return -ENOMEM;
911
912     if (!deflater->IsTrainingMode()) {
913       TString cdbPath("HLT/ConfigTPC/");
914       cdbPath += GetComponentID();
915       cdbPath += "HuffmanTables";
916       TObject* pConf=LoadAndExtractOCDBObject(cdbPath);
917       if (!pConf) return -ENOENT;
918       if (dynamic_cast<TList*>(pConf)==NULL) {
919         HLTError("huffman table configuration object of inconsistent type");
920         return -EINVAL;
921       }
922       iResult=deflater->InitDecoders(dynamic_cast<TList*>(pConf));
923       if (iResult<0) return iResult;
924     }
925     
926     if (!fHistogramFile.IsNull())
927       deflater->EnableStatistics();
928
929     unsigned nofParameters=AliHLTTPCDefinitions::GetNumberOfClusterParameterDefinitions();
930     unsigned p=0;
931     for (; p<nofParameters; p++) {
932       const AliHLTTPCDefinitions::AliClusterParameter& parameter=AliHLTTPCDefinitions::fgkClusterParameterDefinitions[p];
933       // use the pad/time length as reference for the calculation of ratio for residuals
934       unsigned refLength=0;
935       unsigned refLengthPad=0;
936       unsigned refLengthTime=0;
937       if (parameter.fId==AliHLTTPCDefinitions::kPad)               refLengthPad=parameter.fBitLength;
938       else if (parameter.fId==AliHLTTPCDefinitions::kTime)         refLengthTime=parameter.fBitLength;
939       else if (parameter.fId==AliHLTTPCDefinitions::kResidualPad)  refLength=refLengthPad;
940       else if (parameter.fId==AliHLTTPCDefinitions::kResidualTime) refLength=refLengthTime;
941
942       if (deflater->AddParameterDefinition(parameter.fName,
943                                            parameter.fBitLength,
944                                            refLength)!=(int)parameter.fId) {
945         // for performance reason the parameter id is simply used as index in the array of
946         // definitions, the position must match the id
947         HLTFatal("mismatch between parameter id and position in array for parameter %s, rearrange definitions!", parameter.fName);
948         return -EFAULT;
949       }
950     }
951     fpDataDeflater=deflater.release();
952     return 0;
953   }
954   if (mode==kDeflaterModeSimple) {
955     std::auto_ptr<AliHLTDataDeflaterSimple> deflater(new AliHLTDataDeflaterSimple);
956     if (!deflater.get()) return -ENOMEM;
957
958     if (!fHistogramFile.IsNull())
959       deflater->EnableStatistics();
960
961     unsigned nofParameters=AliHLTTPCDefinitions::GetNumberOfClusterParameterDefinitions();
962     unsigned p=0;
963     for (; p<nofParameters; p++) {
964       const AliHLTTPCDefinitions::AliClusterParameter& parameter=AliHLTTPCDefinitions::fgkClusterParameterDefinitions[p];
965       if (deflater->AddParameterDefinition(parameter.fName,
966                                            parameter.fBitLength,
967                                            parameter.fOptional)!=(int)parameter.fId) {
968         // for performance reason the parameter id is simply used as index in the array of
969         // definitions, the position must match the id
970         HLTFatal("mismatch between parameter id and position in array for parameter %s, rearrange definitions!", parameter.fName);
971         return -EFAULT;
972       }
973     }
974     fpDataDeflater=deflater.release();
975     return 0;
976   }
977   HLTError("invalid deflater mode %d, allowed 1=simple 2=huffman", mode);
978   return -EINVAL;
979 }
980
981 int AliHLTTPCDataCompressionComponent::DoDeinit()
982 {
983   /// inherited from AliHLTComponent: component cleanup
984   int iResult=0;
985   if (fpBenchmark) delete fpBenchmark; fpBenchmark=NULL;
986   if (fRawInputClusters) delete fRawInputClusters; fRawInputClusters=NULL;
987   if (fInputClusters) delete fInputClusters; fInputClusters=NULL;
988   if (!fHistogramFile.IsNull()) {
989     TFile out(fHistogramFile, "RECREATE");
990     if (!out.IsZombie()) {
991       out.cd();
992       if (fHistoCompFactor) fHistoCompFactor->Write();
993       if (fHistoResidualPad) fHistoResidualPad->Write();
994       if (fHistoResidualTime) fHistoResidualTime->Write();
995       if (fHistoClusterRatio) fHistoClusterRatio->Write();
996       if (fHistoClustersOnTracks) fHistoClustersOnTracks->Write();
997       if (fHistoTrackClusterRatio) fHistoTrackClusterRatio->Write();
998       out.Close();
999     }
1000   }
1001   if (fHistoCompFactor) delete fHistoCompFactor;
1002   fHistoCompFactor=NULL;
1003   if (fHistoResidualPad) delete fHistoResidualPad;
1004   fHistoResidualPad=NULL;
1005   if (fHistoResidualTime) delete fHistoResidualTime;
1006   fHistoResidualTime=NULL;
1007   if (fHistoClustersOnTracks) delete fHistoClustersOnTracks;
1008   fHistoClustersOnTracks=NULL;
1009   if (fHistoClusterRatio) delete fHistoClusterRatio;
1010   fHistoClusterRatio=NULL;
1011   if (fHistoTrackClusterRatio) delete fHistoTrackClusterRatio;
1012   fHistoTrackClusterRatio=NULL;
1013
1014   if (fpDataDeflater) {
1015     if (!fHistogramFile.IsNull()) {
1016       TString filename=fHistogramFile;
1017       filename.ReplaceAll(".root", "-deflater.root");
1018       fpDataDeflater->SaveAs(filename);
1019     }
1020     if (fDeflaterMode==kDeflaterModeHuffmanTrainer) {
1021       if (fTrainingTableOutput.IsNull()) {
1022         fTrainingTableOutput=GetComponentID();
1023         fTrainingTableOutput+="-huffman.root";
1024       }
1025       // TODO: currently, the code tables are also calculated in FindObject
1026       // check if a different function is more appropriate
1027       TObject* pConf=fpDataDeflater->FindObject("DeflaterConfiguration");
1028       if (pConf) {
1029         TString cdbEntryPath("HLT/ConfigTPC/");
1030         cdbEntryPath += GetComponentID();
1031         cdbEntryPath += "HuffmanTables";
1032         AliCDBPath cdbPath(cdbEntryPath);
1033         AliCDBId cdbId(cdbPath, AliCDBManager::Instance()->GetRun(), AliCDBRunRange::Infinity(), 0, 0);
1034         AliCDBMetaData* cdbMetaData=new AliCDBMetaData;
1035         cdbMetaData->SetResponsible("ALICE HLT Matthias.Richter@cern.ch");
1036         cdbMetaData->SetComment("Huffman encoder configuration");
1037         AliCDBEntry* entry=new AliCDBEntry(pConf, cdbId, cdbMetaData, kTRUE);
1038
1039         entry->SaveAs(fTrainingTableOutput);
1040         // this is a small memory leak
1041         // seg fault in ROOT object handling if the two objects are deleted
1042         // investigate later
1043         //delete entry;
1044         //delete cdbMetaData;
1045       }
1046     }
1047     delete fpDataDeflater;
1048   }
1049   fpDataDeflater=NULL;
1050
1051
1052   if (fTrackGrid) delete fTrackGrid; fTrackGrid=NULL;
1053   if (fSpacePointGrid) delete fSpacePointGrid; fSpacePointGrid=NULL;
1054
1055   return iResult;
1056 }
1057
1058 int AliHLTTPCDataCompressionComponent::ScanConfigurationArgument(int argc, const char** argv)
1059 {
1060   /// inherited from AliHLTComponent: argument scan
1061   int iResult=0;
1062   if (argc<1) return 0;
1063   int bMissingParam=0;
1064   int i=0;
1065   TString argument=argv[i];
1066
1067   do {
1068     // -mode
1069     if (argument.CompareTo("-mode")==0) {
1070       if ((bMissingParam=(++i>=argc))) break;
1071       TString parameter=argv[i];
1072       if (parameter.IsDigit()) {
1073         fMode=parameter.Atoi();
1074         return 2;
1075       } else {
1076         HLTError("invalid parameter for argument %s, expecting number instead of %s", argument.Data(), parameter.Data());
1077         return -EPROTO;
1078       }
1079     }
1080
1081     // -deflater-mode
1082     if (argument.CompareTo("-deflater-mode")==0) {
1083       if ((bMissingParam=(++i>=argc))) break;
1084       TString parameter=argv[i];
1085       if (parameter.IsDigit()) {
1086         fDeflaterMode=parameter.Atoi();
1087         return 2;
1088       } else {
1089         HLTError("invalid parameter for argument %s, expecting number instead of %s", argument.Data(), parameter.Data());
1090         return -EPROTO;
1091       }
1092     }
1093
1094     // -histogram-file
1095     if (argument.CompareTo("-histogram-file")==0) {
1096       if ((bMissingParam=(++i>=argc))) break;
1097       fHistogramFile=argv[i++];
1098       return 2;
1099     }
1100     // -save-histogram-table
1101     if (argument.CompareTo("-save-huffman-table")==0) {
1102       if ((bMissingParam=(++i>=argc))) break;
1103       fTrainingTableOutput=argv[i++];
1104       return 2;
1105     }
1106     // -cluster-verification
1107     if (argument.CompareTo("-cluster-verification")==0) {
1108       if ((bMissingParam=(++i>=argc))) break;
1109       TString parameter=argv[i];
1110       if (parameter.IsDigit()) {
1111         fVerificationMode=parameter.Atoi();
1112         return 2;
1113       } else {
1114         HLTError("invalid parameter for argument %s, expecting number instead of %s", argument.Data(), parameter.Data());
1115         return -EPROTO;
1116       }
1117     }
1118   } while (0); // using do-while only to have break available
1119
1120   if (bMissingParam) {
1121     HLTError("missing parameter for argument %s", argument.Data());
1122     iResult=-EPROTO;
1123   }
1124
1125   return iResult;
1126 }
1127
1128 int AliHLTTPCDataCompressionComponent::InitDriftTimeTransformation()
1129 {
1130   /// calculate correction factor and offset for a linear approximation of the
1131   /// drift time transformation, separately for A and C side
1132   int iResult=0;
1133   AliHLTTPCClusterTransformation transform;
1134   if ((iResult=transform.Init( GetBz(), GetTimeStamp()))<0) {
1135     HLTError("failed to init AliHLTTPCClusterTransformation: %d", iResult);
1136     return iResult;
1137   }
1138
1139   if ((iResult=CalculateDriftTimeTransformation(transform, 0, 0, fDriftTimeFactorA, fDriftTimeOffsetA))<0) return iResult;
1140   if (fVerbosity>0) HLTInfo("drift time transformation A side: m=%f n=%f", fDriftTimeFactorA, fDriftTimeOffsetA);
1141   if ((iResult=CalculateDriftTimeTransformation(transform, 18, 0, fDriftTimeFactorC, fDriftTimeOffsetC))<0) return iResult;
1142   if (fVerbosity>0) HLTInfo("drift time transformation C side: m=%f n=%f", fDriftTimeFactorC, fDriftTimeOffsetC);
1143
1144   return 0;
1145 }
1146
1147 int AliHLTTPCDataCompressionComponent::CalculateDriftTimeTransformation(AliHLTTPCClusterTransformation& transform,
1148                                                                         int slice, int padrow,
1149                                                                         float& m, float& n) const
1150 {
1151   /// calculate correction factor and offset for a linear approximation of the
1152   /// drift time transformation by just probing the range of timebins with
1153   /// AliHLTTPCClusterTransformation
1154   const int nofSteps=100;
1155   vector<float> zvalues;
1156
1157   int nofTimebins=AliHLTTPCTransform::GetNTimeBins();
1158   int stepWidth=nofTimebins/nofSteps;
1159   int time=0;
1160   int count=0;
1161   float meanT=0.;
1162   float meanZ=0.;
1163   for (time=0; time<nofTimebins; time+=stepWidth, count++) {
1164     Float_t xyz[3];
1165     transform.Transform(slice, padrow, 0, time, xyz);
1166     zvalues.push_back(xyz[2]);
1167     meanT+=time;
1168     meanZ+=xyz[2];
1169   }
1170   if (count==0) count=1;
1171   meanT/=count;
1172   meanZ/=count;
1173   float sumTZ=.0;
1174   float sumT2=.0;
1175   time=0;
1176   for (vector<float>::const_iterator z=zvalues.begin();
1177        z!=zvalues.end(); z++, time+=stepWidth) {
1178     sumTZ+=(time-meanT)*((*z)-meanZ);
1179     sumT2+=(time-meanT)*(time-meanT);
1180   }
1181   m=sumTZ/sumT2;
1182   n=meanZ-m*meanT;
1183
1184   return 0;
1185 }