ALIROOT-5433 Transition to CDHv3 in HLT
[u/mrichter/AliRoot.git] / HLT / TPCLib / HWCFemulator / AliHLTTPCHWCFEmulatorComponent.cxx
CommitLineData
33861fe0 1// $Id$
c012881c 2//****************************************************************************
3//* This file is property of and copyright by the ALICE HLT Project *
4//* ALICE Experiment at CERN, All rights reserved. *
5//* *
6//* Primary Authors: Sergey Gorbunov, Torsten Alt *
7//* Developers: Sergey Gorbunov <sergey.gorbunov@fias.uni-frankfurt.de> *
8//* Torsten Alt <talt@cern.ch> *
9//* for The ALICE HLT Project. *
10//* *
11//* Permission to use, copy, modify and distribute this software and its *
12//* documentation strictly for non-commercial purposes is hereby granted *
13//* without fee, provided that the above copyright notice appears in all *
14//* copies and that both the copyright notice and this permission notice *
15//* appear in the supporting documentation. The authors make no claims *
16//* about the suitability of this software for any purpose. It is *
17//* provided "as is" without express or implied warranty. *
18//****************************************************************************
19
20// @file AliHLTTPCHWCFEmulatorComponent.cxx
21// @author Sergey Gorbunov <sergey.gorbunov@fias.uni-frankfurt.de>
22// @author Torsten Alt <talt@cern.ch>
23// @brief HLT Component interface for for FPGA ClusterFinder Emulator for TPC
24// @brief ( see AliHLTTPCHWCFEmulator class )
25// @note
26
27
c012881c 28#include "AliHLTTPCHWCFEmulatorComponent.h"
29
30#include "AliHLTTPCDefinitions.h"
31#include "AliHLTTPCHWCFDataTypes.h"
32#include "AliHLTTPCClusterMCData.h"
33861fe0 33#include "AliHLTTPCHWCFData.h"
c012881c 34
35#include "AliGRPObject.h"
36#include "AliCDBEntry.h"
37#include "AliCDBManager.h"
16e6f752 38#include "AliHLTCDHWrapper.h"
c012881c 39#include <cstdlib>
40#include <cerrno>
e0f315d9 41#include <memory>
c012881c 42#include "TString.h"
43#include "TObjString.h"
44#include "TObjArray.h"
45
46
47#include <sys/time.h>
48#include "TFile.h"
49
d5cf9283 50using namespace std;
51
c012881c 52AliHLTTPCHWCFEmulatorComponent::AliHLTTPCHWCFEmulatorComponent()
53 :
54 AliHLTProcessor(),
b2d72465 55 fDoDeconvTime(0),
56 fDoDeconvPad(0),
57 fDoMC(0),
c012881c 58 fDoFlowControl(0),
b2d72465 59 fDoSinglePadSuppression(0),
c012881c 60 fBypassMerger(0),
61 fClusterLowerLimit(0),
62 fSingleSeqLimit(0),
b2d72465 63 fMergerDistance(0),
64 fUseTimeBinWindow(0),
65 fUseTimeFollow(0),
ef19b824 66 fChargeFluctuation(0),
c012881c 67 fDebug(0),
68 fCFSupport(),
69 fCFEmulator(),
70 fBenchmark("TPCHWClusterFinderEmulator")
71{
72 // see header file for class documentation
73 // or
74 // refer to README to build package
75 // or
76 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
77}
78
79
80AliHLTTPCHWCFEmulatorComponent::AliHLTTPCHWCFEmulatorComponent(const AliHLTTPCHWCFEmulatorComponent&)
81 :
82 AliHLTProcessor(),
b2d72465 83 fDoDeconvTime(0),
84 fDoDeconvPad(0),
85 fDoMC(0),
c012881c 86 fDoFlowControl(0),
b2d72465 87 fDoSinglePadSuppression(0),
c012881c 88 fBypassMerger(0),
89 fClusterLowerLimit(0),
90 fSingleSeqLimit(0),
b2d72465 91 fMergerDistance(0),
92 fUseTimeBinWindow(0),
93 fUseTimeFollow(0),
ef19b824 94 fChargeFluctuation(0),
c012881c 95 fDebug(0),
96 fCFSupport(),
97 fCFEmulator(),
98 fBenchmark("TPCHWClusterFinderEmulator")
99{
100 // dummy
101}
102
103AliHLTTPCHWCFEmulatorComponent& AliHLTTPCHWCFEmulatorComponent::operator=(const AliHLTTPCHWCFEmulatorComponent&)
104{
105 // dummy
106 return *this;
107}
108
109AliHLTTPCHWCFEmulatorComponent::~AliHLTTPCHWCFEmulatorComponent()
110{
111 // see header file for class documentation
112}
113
114// Public functions to implement AliHLTComponent's interface.
115// These functions are required for the registration process
116
117const char* AliHLTTPCHWCFEmulatorComponent::GetComponentID()
118{
119 // see header file for class documentation
120 return "TPCHWClusterFinderEmulator";
121}
122
123void AliHLTTPCHWCFEmulatorComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list)
124{
125 // see header file for class documentation
126 list.clear();
127 list.push_back( AliHLTTPCDefinitions::fgkUnpackedRawDataType );
128 list.push_back( kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTPC );
129}
130
131AliHLTComponentDataType AliHLTTPCHWCFEmulatorComponent::GetOutputDataType()
132{
133 // see header file for class documentation
134 return kAliHLTMultipleDataType;
135}
136
137int AliHLTTPCHWCFEmulatorComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList)
138{
139 // see header file for class documentation
140 tgtList.clear();
141 tgtList.push_back(AliHLTTPCDefinitions::fgkHWClustersDataType | kAliHLTDataOriginTPC );
142 tgtList.push_back(AliHLTTPCDefinitions::fgkAliHLTDataTypeClusterMCInfo | kAliHLTDataOriginTPC );
143 return tgtList.size();
144}
145
146
147void AliHLTTPCHWCFEmulatorComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
148{
149 // see header file for class documentation
150 // XXX TODO: Find more realistic values.
151 constBase = 0;
38a1b995 152 inputMultiplier = (6 * 0.7);
c012881c 153}
154
155
156AliHLTComponent* AliHLTTPCHWCFEmulatorComponent::Spawn()
157{
158 // see header file for class documentation
159 return new AliHLTTPCHWCFEmulatorComponent();
160}
161
162void AliHLTTPCHWCFEmulatorComponent::GetOCDBObjectDescription( TMap* const targetMap){
163// Get a list of OCDB object description needed for the particular component
164
165 if (!targetMap) return;
166
167 // OCDB entries for component arguments
168 targetMap->Add(new TObjString("HLT/ConfigTPC/TPCHWClusterFinder"), new TObjString("component arguments, empty at the moment"));
169}
170
171
172int AliHLTTPCHWCFEmulatorComponent::DoInit( int argc, const char** argv )
173{
174 // see header file for class documentation
175
176 TString arguments = "";
177 for ( int i = 0; i < argc; i++ ) {
178 if ( !arguments.IsNull() ) arguments += " ";
179 arguments += argv[i];
180 }
181
182 return Configure( NULL, NULL, arguments.Data() );
183}
184
185int AliHLTTPCHWCFEmulatorComponent::Reconfigure( const char* cdbEntry, const char* chainId )
186{
187 // Reconfigure the component from OCDB
188
189 return Configure( cdbEntry, chainId, NULL );
190}
191
192int AliHLTTPCHWCFEmulatorComponent::ScanConfigurationArgument(int argc, const char** argv)
193{
194 // see header file for class documentation
195 TString arguments = "";
196 for ( int i = 0; i < argc; i++ ) {
197 if ( !arguments.IsNull() ) arguments += " ";
198 arguments += argv[i];
199 }
200 return ReadConfigurationString(arguments);
201}
202
203
204
205void AliHLTTPCHWCFEmulatorComponent::SetDefaultConfiguration()
206{
207 // Set default configuration for the FPGA ClusterFinder Emulator component
208 // Some parameters can be later overwritten from the OCDB
209
b2d72465 210 fDoDeconvTime = 1;
211 fDoDeconvPad = 1;
212 fDoMC = 0;
c012881c 213 fDoFlowControl = 0;
b2d72465 214 fDoSinglePadSuppression = 0;
c012881c 215 fBypassMerger = 0;
b2d72465 216 fClusterLowerLimit = 10;
c012881c 217 fSingleSeqLimit = 0;
b2d72465 218 fMergerDistance = 4;
25080052 219 fUseTimeBinWindow = 1;
4d67c3e3 220 fUseTimeFollow = 1;
ef19b824 221 fChargeFluctuation = 0;
b7ea7fe9 222 fDebug = 0;
c012881c 223 fBenchmark.Reset();
224 fBenchmark.SetTimer(0,"total");
225 fBenchmark.SetTimer(1,"reco");
226}
227
228int AliHLTTPCHWCFEmulatorComponent::ReadConfigurationString( const char* arguments )
229{
230 // Set configuration parameters for the FPGA ClusterFinder Emulator component
231 // from the string
232
233 int iResult = 0;
234 if ( !arguments ) return iResult;
6dcdce44 235 //cout<<"["<<arguments<<"]"<<endl;
c012881c 236 TString allArgs = arguments;
237 TString argument;
238 int bMissingParam = 0;
239
240 TObjArray* pTokens = allArgs.Tokenize( " " );
241
242 int nArgs = pTokens ? pTokens->GetEntries() : 0;
243
244 for ( int i = 0; i < nArgs; i++ ) {
245 argument = ( ( TObjString* )pTokens->At( i ) )->GetString();
246 if ( argument.IsNull() ) continue;
247
248 if ( argument.CompareTo( "-deconvolute-time" ) == 0 ) {
249 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
250 fDoDeconvTime = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
8366e77e 251 fUseTimeBinWindow = fDoDeconvTime;
252 HLTInfo( "Time deconvolution and using of TimeBin window are set to: %d", fDoDeconvTime );
c012881c 253 continue;
254 }
255
256 if ( argument.CompareTo( "-deconvolute-pad" ) == 0 ) {
257 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
258 fDoDeconvPad = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
259 HLTInfo( "Pad deconvolution is set to: %d", fDoDeconvPad );
260 continue;
261 }
262
263 if ( argument.CompareTo( "-deconvolute" ) == 0 ) {
264 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
265 fDoDeconvTime = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
8366e77e 266 fUseTimeBinWindow = fDoDeconvTime;
c012881c 267 fDoDeconvPad = fDoDeconvTime;
8366e77e 268 HLTInfo( "Time and pad deconvolution and using of TimeBin window are set to: %d", fDoDeconvPad );
c012881c 269 continue;
270 }
271
272 if ( argument.CompareTo( "-do-mc" ) == 0 ) {
273 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
274 fDoMC = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
275 HLTInfo( "MC processing is set to: %d", fDoMC );
276 continue;
277 }
278
279 if ( argument.CompareTo( "-flow-control" ) == 0 ) {
280 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
281 fDoFlowControl = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
282 HLTInfo( "Flow control is set to: %d", fDoFlowControl );
283 continue;
284 }
285
286 if ( argument.CompareTo( "-single-pad-suppression" ) == 0 ) {
287 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
288 fDoSinglePadSuppression = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
289 HLTInfo( "Single pad suppression is set to: %d", fDoSinglePadSuppression );
290 continue;
291 }
292
293 if ( argument.CompareTo( "-bypass-merger" ) == 0 ) {
294 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
295 fBypassMerger = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
296 HLTInfo( "Bypassing merger is set to: %d", fBypassMerger );
297 continue;
298 }
299
300 if ( argument.CompareTo( "-cluster-lower-limit" ) == 0 ) {
301 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
302 fClusterLowerLimit = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
303 HLTInfo( "Cluster lower limit is set to: %d", fClusterLowerLimit );
304 continue;
305 }
306
307 if ( argument.CompareTo( "-single-sequence-limit" ) == 0 ) {
308 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
309 fSingleSeqLimit = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
310 HLTInfo( "Single sequence limit is set to: %d", fSingleSeqLimit );
311 continue;
312 }
313
eb62889b 314 if ( argument.CompareTo( "-merger-distance" ) == 0 ) {
315 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
316 fMergerDistance = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
317 HLTInfo( "Merger distance is set to: %d", fMergerDistance );
318 continue;
319 }
ee4d1bf0 320
25080052 321 if ( argument.CompareTo( "-use-timebin-window" ) == 0 ) {
ee4d1bf0 322 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
25080052 323 fUseTimeBinWindow = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
8366e77e 324 fDoDeconvTime = fUseTimeBinWindow;
325 HLTInfo( "Using TimeBin window and Time deconvolution are set to: %d", fUseTimeBinWindow );
ee4d1bf0 326 continue;
327 }
328
ef19b824 329 if ( argument.CompareTo( "-charge-fluctuation" ) == 0 ) {
330 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
331 fChargeFluctuation = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
332 HLTInfo( "Charge fluctuation is set to: %d", fChargeFluctuation );
333 continue;
334 }
4d67c3e3 335
336 if ( argument.CompareTo( "-use-time-follow" ) == 0 ) {
337 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
338 fUseTimeFollow = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
339 HLTInfo( "Use of time follow algorithm is set to: %d", fUseTimeFollow );
340 continue;
341 }
342
c012881c 343 if ( argument.CompareTo( "-debug-level" ) == 0 ) {
344 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
345 fDebug = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
346 HLTInfo( "Debug level is set to: %d", fDebug );
347 continue;
348 }
349
350 HLTError( "Unknown option \"%s\"", argument.Data() );
351 iResult = -EINVAL;
352 }
353 delete pTokens;
354
355 if ( bMissingParam ) {
356 HLTError( "Specifier missed for parameter \"%s\"", argument.Data() );
357 iResult = -EINVAL;
358 }
359
360 return iResult;
361}
362
363
364int AliHLTTPCHWCFEmulatorComponent::ReadCDBEntry( const char* cdbEntry, const char* chainId )
365{
366 // Read configuration from OCDB
367
368 const char* defaultNotify = "";
369
370 if ( !cdbEntry ){
371 cdbEntry = "HLT/ConfigTPC/TPCHWClusterFinder";
372 defaultNotify = " (default)";
373 chainId = 0;
374 }
375
376 HLTInfo( "configure from entry \"%s\"%s, chain id %s", cdbEntry, defaultNotify, ( chainId != NULL && chainId[0] != 0 ) ? chainId : "<none>" );
377 AliCDBEntry *pEntry = AliCDBManager::Instance()->Get( cdbEntry );//,GetRunNo());
378
379 if ( !pEntry ) {
380 HLTError( "cannot fetch object \"%s\" from CDB", cdbEntry );
381 return -EINVAL;
382 }
383
384 TObjString* pString = dynamic_cast<TObjString*>( pEntry->GetObject() );
385
386 if ( !pString ) {
387 HLTError( "configuration object \"%s\" has wrong type, required TObjString", cdbEntry );
388 return -EINVAL;
389 }
390
391 HLTInfo( "received configuration object string: \"%s\"", pString->GetString().Data() );
392
393 return ReadConfigurationString( pString->GetString().Data() );
394}
395
396
397int AliHLTTPCHWCFEmulatorComponent::Configure( const char* cdbEntry, const char* chainId, const char *commandLine )
398{
399 // Configure the component
400 // There are few levels of configuration,
401 // parameters which are set on one step can be overwritten on the next step
402
403 //* read hard-coded values
404
405 SetDefaultConfiguration();
406
407 //* read the default CDB entry
408
409 int iResult1 = ReadCDBEntry( NULL, chainId );
410
411 //* read the actual CDB entry if required
412
413 int iResult2 = ( cdbEntry ) ? ReadCDBEntry( cdbEntry, chainId ) : 0;
414
415 //* read extra parameters from input (if they are)
416
417 int iResult3 = 0;
418
419 if ( commandLine && commandLine[0] != '\0' ) {
420 HLTInfo( "received configuration string from HLT framework: \"%s\"", commandLine );
421 iResult3 = ReadConfigurationString( commandLine );
422 }
b7ea7fe9 423
424 if( fDebug>1 ) fCFEmulator.SetDebugLevel( fDebug );
425 else fCFEmulator.SetDebugLevel(0);
c012881c 426
427 return iResult1 ? iResult1 : ( iResult2 ? iResult2 : iResult3 );
428}
429
430
431int AliHLTTPCHWCFEmulatorComponent::DoDeinit()
432{
433 // see header file for class documentation
434 return 0;
435}
436
437
438int AliHLTTPCHWCFEmulatorComponent::DoEvent( const AliHLTComponentEventData& evtData,
439 const AliHLTComponentBlockData* blocks,
440 AliHLTComponentTriggerData& /*trigData*/, AliHLTUInt8_t* outputPtr,
441 AliHLTUInt32_t& size,
442 vector<AliHLTComponentBlockData>& outputBlocks )
443{
444 // see header file for class documentation
445
446 int iResult=0;
447 AliHLTUInt32_t maxSize = size;
448 size = 0;
449
450 if(!IsDataEvent()){
451 return 0;
452 }
453
454 fBenchmark.StartNewEvent();
455 fBenchmark.Start(0);
456
eb62889b 457 AliHLTUInt32_t configWord1=0, configWord2=0;
458 AliHLTTPCHWCFEmulator::CreateConfiguration
4d67c3e3 459 ( fDoDeconvTime, fDoDeconvPad, fDoFlowControl, fDoSinglePadSuppression, fBypassMerger, fClusterLowerLimit, fSingleSeqLimit, fMergerDistance, fUseTimeBinWindow, fChargeFluctuation, fUseTimeFollow, configWord1, configWord2 );
c012881c 460
619825aa 461 AliHLTUInt8_t *outBlock = NULL;
462 AliHLTTPCClusterMCLabel *allocOutMC = NULL;
463
c012881c 464 for ( unsigned long ndx = 0; ndx < evtData.fBlockCnt; ndx++ )
465 {
466 const AliHLTComponentBlockData* iter = blocks+ndx;
467
468 if ( iter->fDataType != (kAliHLTDataTypeDDLRaw | kAliHLTDataOriginTPC)
469 && iter->fDataType != AliHLTTPCDefinitions::fgkUnpackedRawDataType ) continue;
470
471 int slice = AliHLTTPCDefinitions::GetMinSliceNr( *iter );
472 int patch = AliHLTTPCDefinitions::GetMinPatchNr( *iter );
473
c012881c 474 fBenchmark.AddInput(iter->fSize);
475
476 if (!iter->fPtr) continue;
477
478 // create input block for the HW cluster finder
479
480 const AliHLTUInt32_t *rawEvent=0;
481 const AliHLTTPCClusterMCLabel *mcLabels = 0;
482 AliHLTUInt32_t rawEventSize32 = 0;
483 AliHLTUInt32_t nMCLabels = 0;
484
485 if( fCFSupport.CreateRawEvent( iter, rawEvent, rawEventSize32, mcLabels, nMCLabels )<0 ) continue;
486 if( !fDoMC ){
487 mcLabels = 0;
488 nMCLabels = 0;
489 }
490
16e6f752 491 AliHLTCDHWrapper header(iter->fPtr);
c012881c 492 // book memory for the output
493
494 AliHLTUInt32_t maxNClusters = rawEventSize32 + 1; // N 32-bit words in input
33861fe0 495 AliHLTUInt32_t clustersSize32 = maxNClusters*AliHLTTPCHWCFData::fgkAliHLTTPCHWClusterSize;
c012881c 496 AliHLTUInt32_t nOutputMC = maxNClusters;
497
16e6f752 498 AliHLTUInt32_t headerSize = header.GetHeaderSize();
e0f315d9 499 AliHLTUInt32_t outBlockSize=headerSize+clustersSize32*sizeof(AliHLTUInt32_t);
619825aa 500
501 if( outBlock ) delete[] outBlock;
502 if( allocOutMC ) delete[] allocOutMC;
503 outBlock = new AliHLTUInt8_t[outBlockSize];
504 allocOutMC = new AliHLTTPCClusterMCLabel[nOutputMC+1];
505
506 if (!outBlock || !allocOutMC) {
507 iResult=-ENOMEM;
508 break;
509 }
e0f315d9 510
3afd5d44 511 memset(outBlock, 0, outBlockSize*sizeof(AliHLTUInt8_t));
512 memset(allocOutMC, 0, (nOutputMC+1)*sizeof(AliHLTTPCClusterMCLabel));
513 AliHLTTPCClusterMCData *outMC = reinterpret_cast<AliHLTTPCClusterMCData *>(allocOutMC);
c012881c 514
515 // fill CDH header here, since the HW clusterfinder does not receive it
16e6f752 516 memcpy(outBlock, iter->fPtr, headerSize);
517 memset(outBlock,0xFF,4);
518
519 //AliRawDataHeader *cdhHeader = reinterpret_cast<AliRawDataHeader*>(iter->fPtr);
520 //AliRawDataHeader *outCDHHeader = reinterpret_cast<AliRawDataHeader*>(outBlock);
521 //*outCDHHeader = *cdhHeader;
522 //outCDHHeader->fSize = 0xFFFFFFFF;
c012881c 523
3afd5d44 524 AliHLTUInt32_t *outClusters = reinterpret_cast<AliHLTUInt32_t*> (outBlock + headerSize);
c012881c 525
526 fBenchmark.Start(1);
527 fCFEmulator.Init
eb62889b 528 ( fCFSupport.GetMapping(slice,patch), configWord1, configWord2 );
c012881c 529
530 int err = fCFEmulator.FindClusters( rawEvent, rawEventSize32,
531 outClusters, clustersSize32,
532 mcLabels, nMCLabels,
533 outMC );
534 fBenchmark.Stop(1);
535 if( err==-1 ){ HLTWarning("NULL input pointer (warning %d)",err);}
536 else if( err==-2 ){ HLTWarning("No space left in the output buffer (warning %d)",err); }
537 else if( err<0 ){ HLTWarning("HWCF emulator finished with error code %d",err); }
538 if( err<0 ){
c012881c 539 continue;
540 }
541
542 if( fDebug ){
33861fe0 543 int elsize=AliHLTTPCHWCFData::fgkAliHLTTPCHWClusterSize;
c012881c 544 printf("\nHWCF Emulator: output clusters for slice%d patch %d:\n",slice,patch);
33861fe0 545 for( AliHLTUInt32_t i=0; i<clustersSize32; i+=elsize ){
c012881c 546 AliHLTUInt32_t *c = outClusters+i;
547 AliHLTUInt32_t flag = (c[0]>>30);
548 if( flag == 0x3){ //beginning of a cluster
549 int padRow = (c[0]>>24)&0x3f;
b02b7690 550 int q = c[1];
551 double p = *((AliHLTFloat32_t*)&c[2]);
552 double t = *((AliHLTFloat32_t*)&c[3]);
553 AliHLTFloat32_t p2 = *((AliHLTFloat32_t*)&c[4]);
554 AliHLTFloat32_t t2 = *((AliHLTFloat32_t*)&c[5]);
c012881c 555 printf("N: %3d R: %3d C: %4d P: %7.4f T: %8.4f DP: %6.4f DT: %6.4f\n",
33861fe0 556 i/elsize+1, padRow, q, p, t, sqrt(fabs(p2-p*p)), sqrt(fabs(t2-t*t)));
c012881c 557
558 if( outMC && outMC->fCount>0 ){
559 printf(" MC: (%3d,%6.1f) (%3d,%6.1f) (%3d,%6.1f)\n",
33861fe0 560 outMC->fLabels[i/elsize].fClusterID[0].fMCID,outMC->fLabels[i/elsize].fClusterID[0].fWeight,
561 outMC->fLabels[i/elsize].fClusterID[1].fMCID,outMC->fLabels[i/elsize].fClusterID[1].fWeight,
562 outMC->fLabels[i/elsize].fClusterID[2].fMCID,outMC->fLabels[i/elsize].fClusterID[2].fWeight
c012881c 563 );
564 }
565 }
566 }
567 }
568
569
570 AliHLTUInt32_t outSize = headerSize + clustersSize32*sizeof(AliHLTUInt32_t);
571
572 if( size + outSize <= maxSize ){
573
3afd5d44 574 memcpy( outputPtr, outBlock, outSize );
c012881c 575
576 AliHLTComponentBlockData bd;
577 FillBlockData( bd );
578 bd.fOffset = size;
579 bd.fSize = outSize;
580 bd.fSpecification = iter->fSpecification;
581 bd.fDataType = AliHLTTPCDefinitions::fgkHWClustersDataType | kAliHLTDataOriginTPC;
582 outputBlocks.push_back( bd );
583 fBenchmark.AddOutput(bd.fSize);
584 size+= bd.fSize;
585 outputPtr+=bd.fSize;
586 } else {
587 HLTWarning( "Output buffer (%db) is too small, required %db", maxSize, size+outSize);
588 iResult=-ENOSPC;
ba9bc1e2 589 break;
c012881c 590 }
591
592 if( fDoMC && outMC && outMC->fCount>0 ){
593 int s = sizeof(AliHLTTPCClusterMCData) + outMC->fCount*sizeof(AliHLTTPCClusterMCLabel);
594 if( size + s <= maxSize ){
595 memcpy( outputPtr, outMC, s );
596 AliHLTComponentBlockData bdMCInfo;
597 FillBlockData( bdMCInfo );
598 bdMCInfo.fOffset = size;
599 bdMCInfo.fSize = s;
600 bdMCInfo.fSpecification = iter->fSpecification;
601 bdMCInfo.fDataType = AliHLTTPCDefinitions::fgkAliHLTDataTypeClusterMCInfo | kAliHLTDataOriginTPC;
602 outputBlocks.push_back( bdMCInfo );
603 fBenchmark.AddOutput(bdMCInfo.fSize);
604 size+=bdMCInfo.fSize;
605 outputPtr+=bdMCInfo.fSize;
606 } else {
607 HLTWarning( "Output buffer (%db) is too small, required %db", maxSize, size+s);
3afd5d44 608 iResult=-ENOSPC;
ba9bc1e2 609 break;
c012881c 610 }
611 }
c012881c 612 }
619825aa 613
614 if( outBlock ) delete[] outBlock;
615 if( allocOutMC ) delete[] allocOutMC;
616
c012881c 617 fBenchmark.Stop(0);
618 HLTInfo(fBenchmark.GetStatistics());
ba9bc1e2 619 fCFSupport.ReleaseEventMemory();
c012881c 620 return iResult;
621}
622
623