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