SPD vertexer component created which reconstructs the XYZ primary vertex with SPD...
[u/mrichter/AliRoot.git] / HLT / ITS / AliHLTITSVertexerSPDComponent.cxx
1 // $Id: AliHLTITSVertexerSPDComponent.cxx 32659 2009-06-02 16:08:40Z sgorbuno $
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 <sergey.gorbunov@kip.uni-heidelberg.de> *
7 //                  Ivan Kisel <kisel@kip.uni-heidelberg.de>                *
8 //                  for The ALICE HLT Project.                              *
9 //                                                                          *
10 // Permission to use, copy, modify and distribute this software and its     *
11 // documentation strictly for non-commercial purposes is hereby granted     *
12 // without fee, provided that the above copyright notice appears in all     *
13 // copies and that both the copyright notice and this permission notice     *
14 // appear in the supporting documentation. The authors make no claims       *
15 // about the suitability of this software for any purpose. It is            *
16 // provided "as is" without express or implied warranty.                    *
17 //                                                                          *
18 //***************************************************************************
19
20 ///  @file   AliHLTITSVertexerSPDComponent.cxx
21 ///  @author Sergey Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de>
22 ///  @date   June 2009
23 ///  @brief  An ITS tracker processing component for the HLT
24
25
26 /////////////////////////////////////////////////////
27 //                                                 //
28 // a ITS tracker processing component for the HLT  //
29 //                                                 //
30 /////////////////////////////////////////////////////
31
32 #if __GNUC__>= 3
33 using namespace std;
34 #endif
35
36 #include "AliHLTITSVertexerSPDComponent.h"
37 #include "AliHLTArray.h"
38 #include "AliExternalTrackParam.h"
39 #include "TStopwatch.h"
40 #include "TMath.h"
41 #include "AliCDBEntry.h"
42 #include "AliCDBManager.h"
43 #include "AliGeomManager.h"
44 #include "TObjString.h"
45 #include "TObjArray.h"
46 #include "AliITStrackerHLT.h"
47 #include "AliHLTITSSpacePointData.h"
48 #include "AliHLTITSClusterDataFormat.h"
49 #include "AliHLTDataTypes.h"
50 #include "AliHLTExternalTrackParam.h"
51 #include "AliHLTGlobalBarrelTrack.h"
52 #include "AliGeomManager.h"
53 #include "TH1F.h"
54 #include "TH2F.h"
55 #include "AliESDVertex.h"
56
57 /** ROOT macro for the implementation of ROOT specific class methods */
58 ClassImp( AliHLTITSVertexerSPDComponent )
59 AliHLTITSVertexerSPDComponent::AliHLTITSVertexerSPDComponent()
60     :
61     fSolenoidBz( 0 ),
62     fProduceHistos(1),
63     fAutoCalibration(1000),
64     fFullTime( 0 ),
65     fRecoTime( 0 ),
66     fNEvents( 0 ),    
67     fHistoVertexXY(0),
68     fHistoVertexX(0),
69     fHistoVertexY(0),
70   fHistoVertexZ(0)
71
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   fDefRunVtx[0] = 0;
80   fDefRunVtx[1] = 0;
81   fDefRunVtx[2] = 0;
82 }
83
84 AliHLTITSVertexerSPDComponent::AliHLTITSVertexerSPDComponent( const AliHLTITSVertexerSPDComponent& )
85     :
86     AliHLTProcessor(),
87     fSolenoidBz( 0 ),
88     fProduceHistos(1),
89     fAutoCalibration(1000),
90     fFullTime( 0 ),
91     fRecoTime( 0 ),
92     fNEvents( 0 ),    
93     fHistoVertexXY(0),
94     fHistoVertexX(0),
95     fHistoVertexY(0),
96   fHistoVertexZ(0)
97 {
98   // see header file for class documentation
99   HLTFatal( "copy constructor untested" );
100 }
101
102 AliHLTITSVertexerSPDComponent& AliHLTITSVertexerSPDComponent::operator=( const AliHLTITSVertexerSPDComponent& )
103 {
104   // see header file for class documentation
105   HLTFatal( "assignment operator untested" );
106   return *this;
107 }
108
109 AliHLTITSVertexerSPDComponent::~AliHLTITSVertexerSPDComponent()
110 {
111   // see header file for class documentation
112   
113   delete fHistoVertexXY;
114   delete fHistoVertexX;
115   delete fHistoVertexY;
116   delete fHistoVertexZ;
117 }
118
119 //
120 // Public functions to implement AliHLTComponent's interface.
121 // These functions are required for the registration process
122 //
123
124 const char* AliHLTITSVertexerSPDComponent::GetComponentID()
125 {
126   // see header file for class documentation
127   return "ITSVertexerSPD";
128 }
129
130 void AliHLTITSVertexerSPDComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list )
131 {
132   // see header file for class documentation
133   list.clear();
134   list.push_back( kAliHLTDataTypeClusters|kAliHLTDataOriginITSSPD );
135 }
136
137 AliHLTComponentDataType AliHLTITSVertexerSPDComponent::GetOutputDataType()
138 {
139   // see header file for class documentation  
140   return kAliHLTMultipleDataType;
141 }
142
143 int AliHLTITSVertexerSPDComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList)
144
145 {
146   // see header file for class documentation
147   tgtList.clear();
148   tgtList.push_back( kAliHLTDataTypeESDVertex|kAliHLTDataOriginITS);
149   tgtList.push_back(kAliHLTDataTypeHistogram);
150   return tgtList.size();
151 }
152
153 void AliHLTITSVertexerSPDComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
154 {
155   // define guess for the output data size
156   constBase = 10000;       // minimum size
157   inputMultiplier = 0.5; // size relative to input
158 }
159
160 AliHLTComponent* AliHLTITSVertexerSPDComponent::Spawn()
161 {
162   // see header file for class documentation
163   return new AliHLTITSVertexerSPDComponent;
164 }
165
166 void AliHLTITSVertexerSPDComponent::SetDefaultConfiguration()
167 {
168   // Set default configuration for the CA tracker component
169   // Some parameters can be later overwritten from the OCDB
170
171   fSolenoidBz = -5.00668;
172   fProduceHistos=1;
173   fAutoCalibration = 1000;
174   fDefRunVtx[0] = 0;
175   fDefRunVtx[1] = 0;
176   fDefRunVtx[2] = 0;
177   fFullTime = 0;
178   fRecoTime = 0;
179   fNEvents = 0;
180 }
181
182 int AliHLTITSVertexerSPDComponent::ReadConfigurationString(  const char* arguments )
183 {
184   // Set configuration parameters for the CA tracker component from the string
185
186   int iResult = 0;
187   if ( !arguments ) return iResult;
188
189   TString allArgs = arguments;
190   TString argument;
191   int bMissingParam = 0;
192
193   TObjArray* pTokens = allArgs.Tokenize( " " );
194
195   int nArgs =  pTokens ? pTokens->GetEntries() : 0;
196
197   for ( int i = 0; i < nArgs; i++ ) {
198     argument = ( ( TObjString* )pTokens->At( i ) )->GetString();
199     if ( argument.IsNull() ) continue;
200
201     if ( argument.CompareTo( "-solenoidBz" ) == 0 ) {
202       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
203       fSolenoidBz = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
204       HLTInfo( "Magnetic Field set to: %f", fSolenoidBz );
205       continue;
206     }
207
208     if ( argument.CompareTo( "-produceHistos" ) == 0 ) {
209       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
210       fProduceHistos = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
211       HLTInfo( "fProduceHistos set to: %d", fProduceHistos );
212       continue;
213     }
214
215     if ( argument.CompareTo( "-runVertex" ) == 0 ) {
216       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
217       fDefRunVtx[0] = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
218       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
219       fDefRunVtx[1] = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
220       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
221       fDefRunVtx[2] = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
222       HLTInfo( "Default run vertex is set to (%f,%f,%f)",fDefRunVtx[0],
223                fDefRunVtx[1],fDefRunVtx[2] );
224       for( int i=0; i<3; i++){
225         fRunVtx[i] = fDefRunVtx[i];
226         fRunVtxNew[i] = fDefRunVtx[i];
227       }
228       fRunVtx[3] = 0.;
229       fRunVtxNew[3] = 0.;
230       continue;
231     }
232
233     if ( argument.CompareTo( "-beamDiamondCalibration" ) == 0 ) {
234       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
235       fAutoCalibration = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
236       HLTInfo( "N events for recalibration of the run vertex is set to: %d", fAutoCalibration );
237       continue;
238     }
239     
240
241     HLTError( "Unknown option \"%s\"", argument.Data() );
242     iResult = -EINVAL;
243   }
244   delete pTokens;
245
246   if ( bMissingParam ) {
247     HLTError( "Specifier missed for parameter \"%s\"", argument.Data() );
248     iResult = -EINVAL;
249   }
250
251   return iResult;
252 }
253
254
255 int AliHLTITSVertexerSPDComponent::ReadCDBEntry( const char* cdbEntry, const char* chainId )
256 {
257   // see header file for class documentation
258
259   const char* defaultNotify = "";
260
261   if ( !cdbEntry ) {
262     return 0;// need to add the HLT/ConfigITS/ITSTracker directory to cdb SG!!!
263     cdbEntry = "HLT/ConfigITS/ITSTracker";
264     defaultNotify = " (default)";
265     chainId = 0;
266   }
267
268   HLTInfo( "configure from entry \"%s\"%s, chain id %s", cdbEntry, defaultNotify, ( chainId != NULL && chainId[0] != 0 ) ? chainId : "<none>" );
269   AliCDBEntry *pEntry = AliCDBManager::Instance()->Get( cdbEntry );//,GetRunNo());
270
271   if ( !pEntry ) {
272     HLTError( "cannot fetch object \"%s\" from CDB", cdbEntry );
273     return -EINVAL;
274   }
275
276   TObjString* pString = dynamic_cast<TObjString*>( pEntry->GetObject() );
277
278   if ( !pString ) {
279     HLTError( "configuration object \"%s\" has wrong type, required TObjString", cdbEntry );
280     return -EINVAL;
281   }
282
283   HLTInfo( "received configuration object string: \"%s\"", pString->GetString().Data() );
284
285   return  ReadConfigurationString( pString->GetString().Data() );
286 }
287
288
289 int AliHLTITSVertexerSPDComponent::Configure( const char* cdbEntry, const char* chainId, const char *commandLine )
290 {
291   // Configure the component
292   // There are few levels of configuration,
293   // parameters which are set on one step can be overwritten on the next step
294
295   //* read hard-coded values
296
297   SetDefaultConfiguration();
298
299   //* read the default CDB entry
300
301   int iResult1 = ReadCDBEntry( NULL, chainId );
302
303   //* read magnetic field
304
305   int iResult2 = ReadCDBEntry( kAliHLTCDBSolenoidBz, chainId );
306
307   //* read the actual CDB entry if required
308
309   int iResult3 = ( cdbEntry ) ? ReadCDBEntry( cdbEntry, chainId ) : 0;
310
311   //* read extra parameters from input (if they are)
312
313   int iResult4 = 0;
314
315   if ( commandLine && commandLine[0] != '\0' ) {
316     HLTInfo( "received configuration string from HLT framework: \"%s\"", commandLine );
317     iResult4 = ReadConfigurationString( commandLine );
318   }
319
320   // Initialise the tracker here
321
322   return iResult1 ? iResult1 : ( iResult2 ? iResult2 : ( iResult3 ? iResult3 : iResult4 ) );
323 }
324
325
326
327 int AliHLTITSVertexerSPDComponent::DoInit( int argc, const char** argv )
328 {
329   // Configure the ITS tracker component
330   
331   fProduceHistos = 1;
332   fAutoCalibration = 1000;
333
334   if(AliGeomManager::GetGeometry()==NULL){
335     AliGeomManager::LoadGeometry("");
336   }
337   AliGeomManager::ApplyAlignObjsFromCDB("ITS");
338   
339   TString arguments = "";
340   for ( int i = 0; i < argc; i++ ) {
341     if ( !arguments.IsNull() ) arguments += " ";
342     arguments += argv[i];
343   }
344
345   int ret = Configure( NULL, NULL, arguments.Data() );
346
347   if( fProduceHistos ){
348     fHistoVertexXY = new TH2F("hITSvertexXY", "ITSSPD vertex in XY", 100,-2,2,100,-2,2);
349     fHistoVertexX = new TH1F("hITSvertexX", "ITSSPD vertex X", 100,-2,2);
350     fHistoVertexY = new TH1F("hITSvertexY", "ITSSPD vertex Y", 100,-2,2);
351     fHistoVertexZ = new TH1F("hITSvertexZ", "ITSSPD vertex Z", 100,-15,15);
352   }
353
354   for( int i=0; i<3; i++){
355     fRunVtx[i] = fDefRunVtx[i];
356     fRunVtxNew[i] = fDefRunVtx[i];
357   }
358   fRunVtx[3] = 0.;
359   fRunVtxNew[3] = 0.;
360
361   return ret;
362 }
363
364
365 int AliHLTITSVertexerSPDComponent::DoDeinit()
366 {
367   // see header file for class documentation
368
369   delete fHistoVertexXY;
370   delete fHistoVertexX;
371   delete fHistoVertexY;
372   delete fHistoVertexZ;
373
374   fHistoVertexXY = 0;
375   fHistoVertexX = 0;
376   fHistoVertexY = 0;
377   fHistoVertexZ = 0;
378   fAutoCalibration = 1000;
379   fProduceHistos = 1;
380   return 0;
381 }
382
383
384
385 int AliHLTITSVertexerSPDComponent::Reconfigure( const char* cdbEntry, const char* chainId )
386 {
387   // Reconfigure the component from OCDB .
388
389   return Configure( cdbEntry, chainId, NULL );
390 }
391
392
393 int AliHLTITSVertexerSPDComponent::DoEvent
394 (
395   const AliHLTComponentEventData& evtData,
396   const AliHLTComponentBlockData* blocks,
397   AliHLTComponentTriggerData& /*trigData*/,
398   AliHLTUInt8_t* /*outputPtr*/,
399   AliHLTUInt32_t& size,
400   vector<AliHLTComponentBlockData>& /*outputBlocks*/ )
401 {
402   //* process event
403
404   AliHLTUInt32_t maxBufferSize = size;
405   size = 0; // output size
406
407   if (!IsDataEvent()) return 0;
408
409   if ( evtData.fBlockCnt <= 0 ) {
410     HLTWarning( "no blocks in event" );
411     return 0;
412   }
413
414
415   TStopwatch timer;
416
417   // Event reconstruction in ITS
418
419   int iResult=0;
420
421   int nBlocks = evtData.fBlockCnt;
422   int nClustersTotal = 0;
423  
424
425   const int kNPhiBins = 20;
426   const double kDPhi = TMath::TwoPi() / kNPhiBins;
427   double vtxX = fRunVtx[0], vtxY = fRunVtx[1], vtxZ = fRunVtx[2];
428
429
430   std::vector<AliHLTITSVertexerSPDComponent::AliHLTITSVZCluster> clusters[2][ kNPhiBins ];
431
432   for (int ndx=0; ndx<nBlocks && iResult>=0; ndx++) {
433
434     const AliHLTComponentBlockData* iter = blocks+ndx;
435  
436     // Read ITS SPD clusters
437
438     if ( iter->fDataType == (kAliHLTDataTypeClusters|kAliHLTDataOriginITSSPD) ){
439
440       AliHLTITSClusterData *inPtr=reinterpret_cast<AliHLTITSClusterData*>( iter->fPtr );
441       int nClusters = inPtr->fSpacePointCnt;
442       nClustersTotal+=nClusters;
443       for( int icl=0; icl<nClusters; icl++ ){   
444         AliHLTITSSpacePointData &d = inPtr->fSpacePoints[icl];
445         if( d.fLayer>1 ) continue;// SPD only;
446         if( d.fLayer<0 ) continue;// SPD only;
447         Int_t lab[4] = { d.fTracks[0], d.fTracks[1], d.fTracks[2], d.fIndex };
448         Int_t info[3] = { d.fNy, d.fNz, d.fLayer };
449         Float_t hit[6] = { d.fY, d.fZ, d.fSigmaY2, d.fSigmaZ2, d.fQ, d.fSigmaYZ };
450         if( d.fLayer==4 ) hit[5] = -hit[5];
451
452
453         AliITSRecPoint p( lab, hit, info );
454
455         if (!p.Misalign()){
456           HLTWarning("Can not misalign an SPD cluster");
457         }
458         Float_t xyz[3];      
459         p.GetGlobalXYZ(xyz);      
460
461         // get phi bin
462         double phi = TMath::ATan2( xyz[1]-vtxY, xyz[0]-vtxX );
463         if( phi<0 ) phi+=TMath::TwoPi();
464         int iphi = (int ) phi/kDPhi;
465         if( iphi<0 ) iphi = 0;
466         if( iphi>=kNPhiBins ) iphi = kNPhiBins-1;
467         AliHLTITSVZCluster c;
468         c.fX = xyz[0];
469         c.fY = xyz[1];
470         c.fZ = xyz[2];
471         clusters[d.fLayer][iphi].push_back( c );        
472      }   
473     }    
474   }// end read input blocks
475   
476
477   // Reconstruct the vertex
478
479   TStopwatch timerReco;
480     
481
482   const double kZBinMin = -10;
483   const int kNZBins = 40;
484   const double kZStep = 1./kNZBins;
485
486   double histX[kNZBins];
487   double histY[kNZBins];
488   double histZ[kNZBins];
489   double histW[kNZBins];
490   int histN[kNZBins];
491   
492   // fit few times with decreasing of the radius cut
493
494   double vtxDR2[2] = { 2.*2., .2*.2};
495   double  maxW = 0;  
496   int bestBin = -1;
497   
498   for( int iter = 0; iter<2; iter++ ){
499     
500     for( int i=0; i<kNZBins; i++ ){
501       histX[i] = 0;
502       histY[i] = 0;
503       histZ[i] = 0;
504       histW[i] = 0;
505       histN[i] = 0;
506     }
507  
508     for( int iPhi=0; iPhi<kNPhiBins; iPhi++ ){
509       int nCluUp = clusters[1][iPhi].size();    
510       int nCluDn = clusters[0][iPhi].size();    
511       for( int icUp=0; icUp<nCluUp; icUp++ ){
512         AliHLTITSVZCluster &cu = clusters[1][iPhi][icUp];
513         double x0 = cu.fX - vtxX;
514         double y0 = cu.fY - vtxY;
515         double z0 = cu.fZ - vtxZ;
516         double bestR2 = 1.e10;
517         int bestDn=-1, bestBin =0;
518         double bestV[3]={0,0,0};
519         for( int icDn=0; icDn<nCluDn; icDn++ ){
520           AliHLTITSVZCluster &cd = clusters[0][iPhi][icDn];
521           double x1 = cd.fX - vtxX;
522           double y1 = cd.fY - vtxY;
523           double z1 = cd.fZ - vtxZ;
524           double dx = x1 - x0;
525           double dy = y1 - y0;
526           double l2 = 1./(dx*dx + dy*dy);
527           double a = x1*y0 - x0*y1;
528           double r2 = a*a*l2; 
529           if( r2>vtxDR2[iter] ) continue;
530           double xv = -dy*a*l2;
531           double yv =  dx*a*l2;
532           double zv = ( (x1*z0-x0*z1)*dx + (y1*z0-y0*z1)*dy )*l2;
533           int zbin = (int)((zv - kZBinMin)*kZStep);
534           if( zbin<0 || zbin>=kNZBins ) continue;
535           if( r2<bestR2 ){
536             bestR2 = r2;
537             bestDn = icDn;
538             bestV[0] = xv;
539             bestV[1] = yv;
540             bestV[2] = zv;
541             bestBin = zbin;
542           }
543         }
544         if( bestDn>=0 ){
545           double w = 1./(1.+bestR2);
546           histX[bestBin]+=bestV[0]*w;
547           histY[bestBin]+=bestV[1]*w;
548           histZ[bestBin]+=bestV[2]*w;
549           histW[bestBin]+=w;
550           histN[bestBin]+=1;
551         }
552       }
553     }
554     
555     maxW = 0;  
556     bestBin = -1;
557     for( int i=0; i<kNZBins; i++ ){
558       if( histW[i]>maxW ){
559         bestBin = i;
560         maxW = histW[i];
561       }
562     }
563     if( bestBin<0 || histN[bestBin] <3 ){
564       bestBin = -1;
565       break;
566     }
567     vtxX +=histX[bestBin]/maxW;
568     vtxY +=histY[bestBin]/maxW;
569     vtxZ +=histZ[bestBin]/maxW;    
570   }
571
572   if( bestBin>=0 ){
573     double nv = 1;
574     double nNew = fRunVtxNew[3] + nv;
575     double v[3] = {vtxX, vtxY, vtxZ};
576     for( int i=0; i<3; i++){
577       fRunVtxNew[i] = ( fRunVtxNew[i]*fRunVtxNew[3] + v[i]*nv )/nNew;
578     }
579     fRunVtxNew[3] = nNew;
580   }  
581   
582   if( fAutoCalibration>0 && fRunVtxNew[3] >= fAutoCalibration ){
583     for( int i=0; i<4; i++ ){
584       fRunVtx[i] = fRunVtxNew[i];
585       fRunVtxNew[i] = 0;
586     }
587     //cout<<"ITSVertexerSPD: set run vtx to "<<fRunVtx[0]<<" "<<fRunVtx[1]<<" "<<fRunVtx[2]<<endl;  
588     if( fRunVtx[0]>3. ) fRunVtx[0] = 3;
589     if( fRunVtx[0]<-3. ) fRunVtx[0] = -3;
590     if( fRunVtx[1]>3. ) fRunVtx[1] = 3;
591     if( fRunVtx[1]<-3. ) fRunVtx[1] = -3;
592     if( fRunVtx[2]>30. ) fRunVtx[2] = 30;
593     if( fRunVtx[2]<30. ) fRunVtx[2] = -30;
594   }
595
596   timerReco.Stop();
597   
598   // Fill output 
599
600   if( bestBin>=0 ){
601     double pos[3] = {vtxX, vtxY, vtxZ};
602     double s = 400.E-4;
603     double cov[6] = {s*s,0,s*s,0,0,s*s};
604     AliESDVertex v(pos, cov, 0, histN[bestBin]);
605     PushBack( (TObject*) &v, kAliHLTDataTypeESDVertex|kAliHLTDataOriginITS,0 );
606
607     //cout<<"ITSVertexerSPD: vtx found "<<vtxX<<" "<<vtxY<<" "<<vtxZ<<endl;
608
609     if( fHistoVertexXY ) fHistoVertexXY->Fill( vtxX, vtxY );
610     if( fHistoVertexX ) fHistoVertexX->Fill( vtxX );
611     if( fHistoVertexY ) fHistoVertexY->Fill( vtxY );
612     if( fHistoVertexZ ) fHistoVertexZ->Fill( vtxZ );
613   }
614
615   if( fHistoVertexXY ) PushBack( (TObject*) fHistoVertexXY, kAliHLTDataTypeHistogram,0);
616   if( fHistoVertexX ) PushBack( (TObject*) fHistoVertexX, kAliHLTDataTypeHistogram,0);
617   if( fHistoVertexY ) PushBack( (TObject*) fHistoVertexY, kAliHLTDataTypeHistogram,0);
618   if( fHistoVertexZ ) PushBack( (TObject*) fHistoVertexZ, kAliHLTDataTypeHistogram,0);
619
620
621   timer.Stop();
622   fFullTime += timer.RealTime();
623   fRecoTime += timerReco.RealTime();
624   fNEvents++;
625
626   // Set log level to "Warning" for on-line system monitoring
627   int hz = ( int ) ( fFullTime > 1.e-10 ? fNEvents / fFullTime : 100000 );
628   int hz1 = ( int ) ( fRecoTime > 1.e-10 ? fNEvents / fRecoTime : 100000 );
629   HLTInfo( "ITS Z Vertexer: input %d clusters; time: full %d / reco %d Hz",
630            nClustersTotal, hz, hz1 );
631
632   return iResult;
633 }