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