ALIROOT-5433 Transition to CDHv3 in HLT
[u/mrichter/AliRoot.git] / HLT / TPCLib / comp / AliHLTTPCDataCompressionComponent.cxx
CommitLineData
7ab8209c 1// $Id$
2//**************************************************************************
b960733a 3//* This file is property of and copyright by the *
7ab8209c 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"
124b5fc8 26#include "AliHLTTPCDataCompressionDescriptor.h"
3ba734d3 27#include "AliHLTTPCRawClustersDescriptor.h"
7ab8209c 28#include "AliHLTTPCTrackGeometry.h"
29#include "AliHLTTPCSpacePointContainer.h"
49bdc466 30#include "AliHLTTPCRawSpacePointContainer.h"
7ab8209c 31#include "AliHLTGlobalBarrelTrack.h"
32#include "AliHLTComponentBenchmark.h"
8b1d554b 33#include "AliHLTDataDeflaterSimple.h"
80ae4557 34#include "AliHLTDataDeflaterHuffman.h"
4f573ca4 35#include "AliHLTTPCTransform.h"
0fa10601 36#include "AliHLTTPCClusterMCData.h"
9b72add5 37#include "AliHLTTPCClusterTransformation.h"
4653ce26 38#include "AliHLTErrorGuard.h"
80ae4557 39#include "AliCDBManager.h"
40#include "AliCDBPath.h"
41#include "AliCDBId.h"
42#include "AliCDBMetaData.h"
43#include "AliCDBEntry.h"
8b1d554b 44#include "TH1F.h"
45#include "TFile.h"
46#include <memory>
7ab8209c 47
0fa10601 48ClassImp(AliHLTTPCDataCompressionComponent)
49
7ab8209c 50AliHLTTPCDataCompressionComponent::AliHLTTPCDataCompressionComponent()
51 : AliHLTProcessor()
928d50f1 52 , fMode(kCompressionModeNone)
53 , fDeflaterMode(kDeflaterModeNone)
9b72add5 54 , fVerificationMode(0)
4f573ca4 55 , fMaxDeltaPad(AliHLTTPCDefinitions::GetMaxClusterDeltaPad())
56 , fMaxDeltaTime(AliHLTTPCDefinitions::GetMaxClusterDeltaTime())
7ab8209c 57 , fRawInputClusters(NULL)
58 , fInputClusters(NULL)
0fa10601 59 , fTrackGrid(NULL)
60 , fSpacePointGrid(NULL)
8b1d554b 61 , fpDataDeflater(NULL)
62 , fHistoCompFactor(NULL)
4f573ca4 63 , fHistoResidualPad(NULL)
64 , fHistoResidualTime(NULL)
65 , fHistoClustersOnTracks(NULL)
66 , fHistoClusterRatio(NULL)
67 , fHistoTrackClusterRatio(NULL)
8b1d554b 68 , fHistogramFile()
80ae4557 69 , fTrainingTableOutput()
7ab8209c 70 , fpBenchmark(NULL)
9b72add5 71 , fpWrittenAssociatedClusterIds(NULL)
72 , fDriftTimeFactorA(1.)
73 , fDriftTimeOffsetA(0.)
74 , fDriftTimeFactorC(1.)
75 , fDriftTimeOffsetC(0.)
4f573ca4 76 , fVerbosity(0)
7ab8209c 77{
78}
79
80AliHLTTPCDataCompressionComponent::~AliHLTTPCDataCompressionComponent()
81{
82 /// destructor
9b72add5 83 if (fpWrittenAssociatedClusterIds) delete fpWrittenAssociatedClusterIds;
7ab8209c 84}
85
86
87const char* AliHLTTPCDataCompressionComponent::GetComponentID()
88{
89 /// inherited from AliHLTComponent: id of the component
90 return "TPCDataCompressor";
91}
92
93
94void AliHLTTPCDataCompressionComponent::GetInputDataTypes( AliHLTComponentDataTypeList& tgtList)
95{
96 /// inherited from AliHLTComponent: list of data types in the vector reference
97 tgtList.clear();
49bdc466 98 tgtList.push_back(AliHLTTPCDefinitions::RawClustersDataType());
3ba734d3 99 tgtList.push_back(AliHLTTPCDefinitions::RawClustersDescriptorDataType());
b960733a 100 tgtList.push_back(AliHLTTPCDefinitions::ClustersDataType());
7ab8209c 101 tgtList.push_back(kAliHLTDataTypeTrack|kAliHLTDataOriginTPC);
102}
103
104AliHLTComponentDataType AliHLTTPCDataCompressionComponent::GetOutputDataType()
105{
106 /// inherited from AliHLTComponent: output data type of the component.
107 return kAliHLTMultipleDataType;
108}
109
110int AliHLTTPCDataCompressionComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList)
111{
112 /// inherited from AliHLTComponent: multiple output data types of the component.
113 tgtList.clear();
124b5fc8 114 tgtList.push_back(AliHLTTPCDefinitions::DataCompressionDescriptorDataType());
9b72add5 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());
7ab8209c 120 return tgtList.size();
121}
122
123void AliHLTTPCDataCompressionComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
124{
125 /// inherited from AliHLTComponent: output data size estimator
126 constBase=0;
9b72add5 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
7ab8209c 130}
131
132AliHLTComponent* AliHLTTPCDataCompressionComponent::Spawn()
133{
134 /// inherited from AliHLTComponent: spawn function.
135 return new AliHLTTPCDataCompressionComponent;
136}
137
8a3426fd 138void 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"));
928d50f1 145 if (fDeflaterMode==kDeflaterModeHuffman) {
8a3426fd 146 targetMap->Add(new TObjString("HLT/ConfigTPC/TPCDataCompressorHuffmanTables"),
928d50f1 147 new TObjString("huffman tables for deflater mode 'huffman'"));
8a3426fd 148 }
149}
150
b60d3f6a 151int AliHLTTPCDataCompressionComponent::DoEvent( const AliHLTComponentEventData& /*evtData*/,
152 const AliHLTComponentBlockData* /*inputBlocks*/,
7ab8209c 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
0fa10601 165 if (!fRawInputClusters) {
166 return -ENODEV;
167 }
168
7ab8209c 169 if (GetBenchmarkInstance()) {
170 GetBenchmarkInstance()->StartNewEvent();
171 GetBenchmarkInstance()->Start(0);
172 }
173
174 // Process an event
124b5fc8 175
7ab8209c 176 // Loop over all input blocks in the event
0fa10601 177 bool bHaveMC=(GetFirstInputBlock(AliHLTTPCDefinitions::fgkAliHLTDataTypeClusterMCInfo | kAliHLTDataOriginTPC))!=NULL;
9b72add5 178 if ((bHaveMC || fVerificationMode>0) && fpWrittenAssociatedClusterIds==NULL) {
179 fpWrittenAssociatedClusterIds=new vector<AliHLTUInt32_t>;
180 }
0fa10601 181
7ab8209c 182 const AliHLTComponentBlockData* pDesc=NULL;
183
8b1d554b 184 AliHLTUInt8_t minSlice=0xFF, maxSlice=0xFF, minPatch=0xFF, maxPatch=0xFF;
185 AliHLTUInt32_t inputRawClusterSize=0;
0fa10601 186 AliHLTUInt32_t outputDataSize=0;
4f573ca4 187 int allClusters=0;
188 int associatedClusters=0;
8a3426fd 189 float bz=GetBz();
8b1d554b 190
7ab8209c 191 /// input track array
192 vector<AliHLTGlobalBarrelTrack> inputTrackArray;
193
194 if (GetBenchmarkInstance()) {
7ab8209c 195 GetBenchmarkInstance()->Start(2);
196 }
5a3ff15b 197
198 bool isInputPresent = kFALSE;
7ab8209c 199
200 // transformed clusters
928d50f1 201 // the transformed clusters have not been used yet
202 if (false) { // FIXME: condition to be adjusted
7ab8209c 203 for (pDesc=GetFirstInputBlock(AliHLTTPCDefinitions::fgkClustersDataType);
204 pDesc!=NULL; pDesc=GetNextInputBlock()) {
205 if (GetBenchmarkInstance()) {
206 GetBenchmarkInstance()->AddInput(pDesc->fSize);
207 }
5a3ff15b 208 isInputPresent = kTRUE;
7ab8209c 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
4653ce26 227 vector<int> trackindexmap; // stores index for every track id
228
7ab8209c 229 // track data input
928d50f1 230 if (fMode==kCompressionModeV1TrackModel || fMode==kCompressionModeV2TrackModel) {
7ab8209c 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 }
4653ce26 248 trackindexmap.resize(inputTrackArray.size(), -1);
7ab8209c 249 }
250 }
251
252 if (GetBenchmarkInstance()) {
253 GetBenchmarkInstance()->Stop(3);
254 GetBenchmarkInstance()->Start(4);
255 }
256
257 // processing
4653ce26 258 int trackindex=0;
0fa10601 259 for (vector<AliHLTGlobalBarrelTrack>::iterator track=inputTrackArray.begin();
7ab8209c 260 track!=inputTrackArray.end();
4653ce26 261 track++, trackindex++) {
7ab8209c 262 int trackID=track->GetID();
263 if (trackID<0) {
264 // FIXME: error guard
265 HLTError("invalid track ID");
266 continue;
267 }
4653ce26 268 if (trackID>=(int)trackindexmap.size())
269 trackindexmap.resize(trackID+1, -1);
270 trackindexmap[trackID]=trackindex;
4f573ca4 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
9b72add5 283 AliHLTTPCTrackGeometry* trackpoints=new AliHLTTPCTrackGeometry;
0fa10601 284 if (!trackpoints) continue;
4653ce26 285 HLTDebug("track %d id %d:", trackindex, trackID);
8a3426fd 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);
9b72add5 310 trackpoints->InitDriftTimeTransformation(fDriftTimeFactorA, fDriftTimeOffsetA, fDriftTimeFactorC, fDriftTimeOffsetC);
0fa10601 311 trackpoints->SetTrackId(trackID);
8a3426fd 312 trackpoints->CalculateTrackPoints(ctrack);
0fa10601 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 }
4f573ca4 324 if (fVerbosity>0) {
325 fTrackGrid->Print();
7ab8209c 326 }
327
328 if (GetBenchmarkInstance()) {
329 GetBenchmarkInstance()->Stop(4);
330 GetBenchmarkInstance()->Start(5);
331 }
332
61e66346 333 // loop over raw cluster blocks, assign to tracks and write
334 // unassigned clusters
49bdc466 335 for (pDesc=GetFirstInputBlock(AliHLTTPCDefinitions::fgkRawClustersDataType);
abee47d8 336 pDesc!=NULL; pDesc=GetNextInputBlock()) {
61e66346 337 if (GetBenchmarkInstance()) {
338 GetBenchmarkInstance()->Start(1);
339 GetBenchmarkInstance()->AddInput(pDesc->fSize);
340 }
5a3ff15b 341 isInputPresent = kTRUE;
61e66346 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;
0fa10601 350 inputRawClusterSize+=pDesc->fSize;
351
352 // add the data and populate the index grid
353 fRawInputClusters->AddInputBlock(pDesc);
354 fRawInputClusters->PopulateAccessGrid(fSpacePointGrid, pDesc->fSpecification);
4f573ca4 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 }
61e66346 367 if (GetBenchmarkInstance()) {
368 GetBenchmarkInstance()->Stop(1);
0fa10601 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);
61e66346 376 GetBenchmarkInstance()->Start(5);
377 }
4f573ca4 378 allClusters+=fSpacePointGrid->GetNumberOfSpacePoints();
1e21c76d 379 iResult=ProcessTrackClusters(&inputTrackArray[0], inputTrackArray.size(), fTrackGrid, trackindexmap, fSpacePointGrid, fRawInputClusters, slice, patch);
4f573ca4 380 int assignedInThisPartition=0;
381 if (iResult>=0) {
382 assignedInThisPartition=iResult;
383 associatedClusters+=iResult;
384 }
3ca47235 385 iResult=ProcessRemainingClusters(&inputTrackArray[0], inputTrackArray.size(), fTrackGrid, trackindexmap, fSpacePointGrid, fRawInputClusters, slice, patch);
4f573ca4 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 }
0fa10601 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
9b72add5 398 const char* writeoptions="";
399 if (fpWrittenAssociatedClusterIds) {
400 writeoptions="write-cluster-ids";
401 }
0fa10601 402 fRawInputClusters->SetSpacePointPropertyGrid(pDesc->fSpecification, fSpacePointGrid);
9b72add5 403 iResult=fRawInputClusters->Write(outputPtr+size, capacity-size, outputBlocks, fpDataDeflater, writeoptions);
0fa10601 404 fRawInputClusters->SetSpacePointPropertyGrid(pDesc->fSpecification, NULL);
7ab8209c 405 if (iResult>=0) {
406 size+=iResult;
0fa10601 407 outputDataSize+=iResult;
b60d3f6a 408 // the size of the optional cluster id array must be subtracted
409 if (fpWrittenAssociatedClusterIds && outputBlocks.size()>0 &&
4653ce26 410 outputBlocks.back().fDataType==AliHLTTPCDefinitions::RemainingClusterIdsDataType()) {
b60d3f6a 411 outputDataSize-=outputBlocks.back().fSize;
412 }
7ab8209c 413 if (GetBenchmarkInstance()) GetBenchmarkInstance()->AddOutput(iResult);
414 }
61e66346 415 if (GetBenchmarkInstance()) {
416 GetBenchmarkInstance()->Stop(5);
417 }
0fa10601 418
0fa10601 419 fSpacePointGrid->Clear();
7ab8209c 420 }
4f573ca4 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 }
7ab8209c 426
0fa10601 427 // output of track model clusters
bdad7033 428 if (iResult>=0) do {
429 AliHLTUInt32_t tracksBufferOffset=sizeof(AliHLTTPCTrackModelBlock);
430 if (capacity-size<tracksBufferOffset) {
431 iResult=-ENOSPC;
432 break;
433 }
9b72add5 434 if (fpWrittenAssociatedClusterIds) fpWrittenAssociatedClusterIds->clear();
bdad7033 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;
8a3426fd 448 trackModelBlock->fGlobalParameters[parameterIndex++]=bz;
bdad7033 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
8a3426fd 459 if (trackindexmap.size()>0) {// condition for track model compression
bdad7033 460 iResult=WriteTrackClusters(inputTrackArray, fRawInputClusters, fpDataDeflater, outputPtr+size+tracksBufferOffset, capacity-size-tracksBufferOffset);
4f573ca4 461 if (iResult>=0) {
462 AliHLTComponent_BlockData bd;
463 FillBlockData(bd);
464 bd.fOffset = size;
bdad7033 465 bd.fSize = tracksBufferOffset+iResult;
4f573ca4 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);
9b72add5 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);
9b2c9b8b 477 if (capacity-size>bd.fSize) {
bdad7033 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 }
9b72add5 486
487 fpWrittenAssociatedClusterIds->clear();
488 }
4f573ca4 489 }
b60d3f6a 490 }
bdad7033 491 } while (0);
0fa10601 492
493 fRawInputClusters->Clear();
494
5a3ff15b 495 // Write header block
3ba734d3 496
5a3ff15b 497 if( isInputPresent ){
3ba734d3 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
29e2013e 515 AliHLTTPCDataCompressionDescriptor compDesc;
3ba734d3 516 compDesc.SetMergedClustersFlag( clDesc.GetMergedClustersFlag() );
29e2013e 517 *(AliHLTTPCDataCompressionDescriptor*)(outputPtr + bd.fOffset ) = compDesc;
3ba734d3 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 }
5a3ff15b 524 }
525
8b1d554b 526 float compressionFactor=(float)inputRawClusterSize;
0fa10601 527 if ((outputDataSize)>0) compressionFactor/=outputDataSize;
8b1d554b 528 else compressionFactor=0.;
529 if (fHistoCompFactor) fHistoCompFactor->Fill(compressionFactor);
530
b960733a 531 if (GetBenchmarkInstance() && allClusters>0) {
7ab8209c 532 GetBenchmarkInstance()->Stop(0);
928d50f1 533 if (fDeflaterMode!=kDeflaterModeHuffmanTrainer) {
80ae4557 534 HLTBenchmark("%s - compression factor %.2f", GetBenchmarkInstance()->GetStatistics(), compressionFactor);
535 } else {
536 HLTBenchmark("%s", GetBenchmarkInstance()->GetStatistics());
537 }
7ab8209c 538 }
539
540 if (fInputClusters) {
541 fInputClusters->Clear();
542 }
543 if (fRawInputClusters) {
544 fRawInputClusters->Clear();
545 }
0fa10601 546 if (fTrackGrid) {
547 fTrackGrid->Clear();
8b1d554b 548 }
549
b60d3f6a 550 // forward MC labels
551 for (pDesc=GetFirstInputBlock(AliHLTTPCDefinitions::fgkAliHLTDataTypeClusterMCInfo | kAliHLTDataOriginTPC);
552 pDesc!=NULL; pDesc=GetNextInputBlock()) {
553 outputBlocks.push_back(*pDesc);
554 }
555
7ab8209c 556 return iResult;
557}
558
4f573ca4 559int AliHLTTPCDataCompressionComponent::ProcessTrackClusters(AliHLTGlobalBarrelTrack* pTracks, unsigned nofTracks,
560 AliHLTTrackGeometry::AliHLTTrackGrid* pTrackIndex,
1e21c76d 561 const vector<int>& trackIndexMap,
4f573ca4 562 AliHLTSpacePointContainer::AliHLTSpacePointPropertyGrid* pClusterIndex,
563 AliHLTSpacePointContainer* pClusters,
564 int slice, int partition) const
565{
566 // process to assigned track clusters
1e21c76d 567 int iResult=0;
4f573ca4 568 int assignedClusters=0;
569 if (!pTracks || nofTracks==0) return 0;
570
1e21c76d 571 vector<int> processedTracks(nofTracks, -1);
4f573ca4 572 for (AliHLTTrackGeometry::AliHLTTrackGrid::iterator& trackId=pTrackIndex->begin(slice, partition, -1);
573 trackId!=pTrackIndex->end(); trackId++) {
1e21c76d 574 if (trackId.Data()>=trackIndexMap.size()) {
575 HLTError("can not find track id %d in index map of size %d", trackId.Data(), trackIndexMap.size());
4f573ca4 576 continue;
577 }
1e21c76d 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());
4f573ca4 581 continue;
582 }
1e21c76d 583 if (processedTracks[trackindex]>0) continue;
584 processedTracks[trackindex]=1;
4f573ca4 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 }
4f573ca4 605
606 int clusterrow=(int)pClusters->GetX(clusterId);
4f573ca4 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();
1e21c76d 615
616 iResult=FindCellClusters(trackId.Data(), clusterrow, pad, time, pClusterIndex, pClusters, point, clusterId);
617 if (iResult>0) assignedClusters+=iResult;
4f573ca4 618 }
619 }
620 return assignedClusters;
621}
622
623int AliHLTTPCDataCompressionComponent::ProcessRemainingClusters(AliHLTGlobalBarrelTrack* pTracks, unsigned nofTracks,
624 AliHLTTrackGeometry::AliHLTTrackGrid* pTrackIndex,
3ca47235 625 const vector<int>& trackIndexMap,
4f573ca4 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++) {
3ca47235 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;
4f573ca4 641 }
3ca47235 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());
4f573ca4 645 continue;
646 }
3ca47235 647 AliHLTGlobalBarrelTrack& track=pTracks[trackindex];
4f573ca4 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
677int AliHLTTPCDataCompressionComponent::FindCellClusters(int trackId, int padrow, float pad, float time,
678 AliHLTSpacePointContainer::AliHLTSpacePointPropertyGrid* pClusterIndex,
679 AliHLTSpacePointContainer* pClusters,
1e21c76d 680 AliHLTTrackGeometry::AliHLTTrackPoint* pTrackPoint,
681 AliHLTUInt32_t clusterId) const
4f573ca4 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
4653ce26 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);
1e21c76d 693 if ((TMath::Abs(centerpad-pad)>fMaxDeltaPad && pad>0.) ||
694 (TMath::Abs(centertime-time)>fMaxDeltaTime && time>0.)) {
4653ce26 695 ALIHLTERRORGUARD(20, "invalid pad center calculation, please check dimensions if dimensions of index grid match the maximum possible deviation");
696 }
697
1e21c76d 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) {
4f573ca4 703 if (padindex<0) continue;
704 if (padindex>=pClusterIndex->GetDimensionY()) break;
1e21c76d 705 for (int timecount=0, timeindex=timestartindex; timecount<2; timecount++, timeindex+=timedirection) {
4f573ca4 706 if (timeindex<0) continue;
707 if (timeindex>=pClusterIndex->GetDimensionZ()) break;
4653ce26 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);
4f573ca4 712 cl!=pClusterIndex->end(); cl++) {
713 if (cl.Data().fTrackId>=0) continue;
1e21c76d 714 if (clusterId!=(~(AliHLTUInt32_t)0) && clusterId!=cl.Data().fId) continue;
4653ce26 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 }
4f573ca4 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 }
1e21c76d 728 if (clusterId!=(~(AliHLTUInt32_t)0)) break;
4f573ca4 729 }
730 }
731 }
732 return count;
733}
734
735int 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
9b72add5 756 int result=pTrackPoints->Write(*track, pSpacePoints, pDeflater, outputPtr+size, capacity-size, fpWrittenAssociatedClusterIds);
4f573ca4 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
7ab8209c 797int 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 //
8b1d554b 810 iResult = ConfigureFromCDBTObjString(cdbPath);
7ab8209c 811 if (iResult < 0)
812 return iResult;
813
814 //Stage 3: command line arguments.
928d50f1 815 int mode=fMode; // just a backup for the info message below
816 int deflaterMode=fDeflaterMode;
7ab8209c 817 if (argc && (iResult = ConfigureFromArgumentString(argc, argv)) < 0)
818 return iResult;
819
928d50f1 820 if (mode!=fMode || deflaterMode!=fDeflaterMode) {
821 HLTInfo("configured from command line: mode %d, deflater mode %d", fMode, fDeflaterMode);
822 }
823
8b1d554b 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;
7ab8209c 834 }
835
8a3426fd 836 unsigned spacePointContainerMode=0;
928d50f1 837 if (fMode==kCompressionModeV1TrackModel || fMode==kCompressionModeV2TrackModel) {
8a3426fd 838 // initialize map data for cluster access in the track association loop
49bdc466 839 spacePointContainerMode|=AliHLTTPCRawSpacePointContainer::kModeCreateMap;
8a3426fd 840 }
928d50f1 841 if (fMode==kCompressionModeV2 || fMode==kCompressionModeV2TrackModel) {
8a3426fd 842 // optimized storage format: differential pad and time storage
49bdc466 843 spacePointContainerMode|=AliHLTTPCRawSpacePointContainer::kModeDifferentialPadTime;
8a3426fd 844 }
49bdc466 845 std::auto_ptr<AliHLTTPCRawSpacePointContainer> rawInputClusters(new AliHLTTPCRawSpacePointContainer(spacePointContainerMode));
8b1d554b 846 std::auto_ptr<AliHLTTPCSpacePointContainer> inputClusters(new AliHLTTPCSpacePointContainer);
4f573ca4 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
0fa10601 867 // track grid: 36 slices, each 6 partitions with max 33 rows
868 fTrackGrid=new AliHLTTrackGeometry::AliHLTTrackGrid(36, 1, 6, 1, 33, 1, 20000);
49bdc466 869 fSpacePointGrid=AliHLTTPCRawSpacePointContainer::AllocateIndexGrid();
0fa10601 870
871 if (!rawInputClusters.get() ||
872 !inputClusters.get() ||
0fa10601 873 !fTrackGrid ||
874 !fSpacePointGrid) {
875 if (fTrackGrid) delete fTrackGrid; fTrackGrid=NULL;
876 if (fSpacePointGrid) delete fSpacePointGrid; fSpacePointGrid=NULL;
877 return -ENOMEM;
878 }
8b1d554b 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();
4f573ca4 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 }
7ab8209c 897
9b72add5 898 if (iResult>=0 && (iResult=InitDriftTimeTransformation())<0) return iResult;
899
7ab8209c 900 return iResult;
901}
902
8b1d554b 903int AliHLTTPCDataCompressionComponent::InitDeflater(int mode)
904{
905 /// init the data deflater
80ae4557 906 int iResult=0;
928d50f1 907 if (mode==kDeflaterModeHuffman || mode==kDeflaterModeHuffmanTrainer) {
80ae4557 908 // huffman deflater
928d50f1 909 std::auto_ptr<AliHLTDataDeflaterHuffman> deflater(new AliHLTDataDeflaterHuffman(mode==kDeflaterModeHuffmanTrainer));
80ae4557 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
5e76efa5 926 if (!fHistogramFile.IsNull())
927 deflater->EnableStatistics();
928
80ae4557 929 unsigned nofParameters=AliHLTTPCDefinitions::GetNumberOfClusterParameterDefinitions();
930 unsigned p=0;
931 for (; p<nofParameters; p++) {
932 const AliHLTTPCDefinitions::AliClusterParameter& parameter=AliHLTTPCDefinitions::fgkClusterParameterDefinitions[p];
5e76efa5 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
80ae4557 942 if (deflater->AddParameterDefinition(parameter.fName,
5e76efa5 943 parameter.fBitLength,
944 refLength)!=(int)parameter.fId) {
80ae4557 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 }
928d50f1 954 if (mode==kDeflaterModeSimple) {
8b1d554b 955 std::auto_ptr<AliHLTDataDeflaterSimple> deflater(new AliHLTDataDeflaterSimple);
956 if (!deflater.get()) return -ENOMEM;
957
5e76efa5 958 if (!fHistogramFile.IsNull())
959 deflater->EnableStatistics();
960
0fa10601 961 unsigned nofParameters=AliHLTTPCDefinitions::GetNumberOfClusterParameterDefinitions();
962 unsigned p=0;
963 for (; p<nofParameters; p++) {
71300445 964 const AliHLTTPCDefinitions::AliClusterParameter& parameter=AliHLTTPCDefinitions::fgkClusterParameterDefinitions[p];
965 if (deflater->AddParameterDefinition(parameter.fName,
966 parameter.fBitLength,
967 parameter.fOptional)!=(int)parameter.fId) {
8b1d554b 968 // for performance reason the parameter id is simply used as index in the array of
969 // definitions, the position must match the id
0fa10601 970 HLTFatal("mismatch between parameter id and position in array for parameter %s, rearrange definitions!", parameter.fName);
8b1d554b 971 return -EFAULT;
972 }
973 }
974 fpDataDeflater=deflater.release();
975 return 0;
8b1d554b 976 }
977 HLTError("invalid deflater mode %d, allowed 1=simple 2=huffman", mode);
978 return -EINVAL;
979}
980
7ab8209c 981int AliHLTTPCDataCompressionComponent::DoDeinit()
982{
983 /// inherited from AliHLTComponent: component cleanup
984 int iResult=0;
8b1d554b 985 if (fpBenchmark) delete fpBenchmark; fpBenchmark=NULL;
986 if (fRawInputClusters) delete fRawInputClusters; fRawInputClusters=NULL;
987 if (fInputClusters) delete fInputClusters; fInputClusters=NULL;
4f573ca4 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();
8b1d554b 999 }
8b1d554b 1000 }
4f573ca4 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
80ae4557 1014 if (fpDataDeflater) {
9b72add5 1015 if (!fHistogramFile.IsNull()) {
1016 TString filename=fHistogramFile;
1017 filename.ReplaceAll(".root", "-deflater.root");
1018 fpDataDeflater->SaveAs(filename);
1019 }
928d50f1 1020 if (fDeflaterMode==kDeflaterModeHuffmanTrainer) {
80ae4557 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
0fa10601 1052 if (fTrackGrid) delete fTrackGrid; fTrackGrid=NULL;
1053 if (fSpacePointGrid) delete fSpacePointGrid; fSpacePointGrid=NULL;
8b1d554b 1054
7ab8209c 1055 return iResult;
1056}
1057
1058int AliHLTTPCDataCompressionComponent::ScanConfigurationArgument(int argc, const char** argv)
1059{
1060 /// inherited from AliHLTComponent: argument scan
1061 int iResult=0;
8b1d554b 1062 if (argc<1) return 0;
7ab8209c 1063 int bMissingParam=0;
1064 int i=0;
8b1d554b 1065 TString argument=argv[i];
7ab8209c 1066
8b1d554b 1067 do {
7ab8209c 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();
8b1d554b 1074 return 2;
7ab8209c 1075 } else {
1076 HLTError("invalid parameter for argument %s, expecting number instead of %s", argument.Data(), parameter.Data());
1077 return -EPROTO;
1078 }
1079 }
8b1d554b 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 }
80ae4557 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 }
9b2c9b8b 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 }
0fa10601 1118 } while (0); // using do-while only to have break available
7ab8209c 1119
1120 if (bMissingParam) {
1121 HLTError("missing parameter for argument %s", argument.Data());
1122 iResult=-EPROTO;
1123 }
1124
7ab8209c 1125 return iResult;
1126}
0fa10601 1127
9b72add5 1128int 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
1147int 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 }
00eee21e 1170 if (count==0) count=1;
9b72add5 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}