]> git.uio.no Git - u/mrichter/AliRoot.git/blob - HLT/TPCLib/tracking-ca/AliHLTTPCCATrackerComponent.cxx
799e32a298ca7f8b3ee9239ad00b833d365e1aa4
[u/mrichter/AliRoot.git] / HLT / TPCLib / tracking-ca / AliHLTTPCCATrackerComponent.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 <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
21 ///////////////////////////////////////////////////////////////////////////////
22 //                                                                           //
23 // a TPC tracker processing component for the HLT based on CA by Ivan Kisel  //
24 //                                                                           //
25 ///////////////////////////////////////////////////////////////////////////////
26
27 #if __GNUC__>= 3
28 using namespace std;
29 #endif
30
31 #include "AliHLTTPCCATrackerComponent.h"
32 #include "AliHLTTPCTransform.h"
33 #include "AliHLTTPCCATrackerFramework.h"
34 #include "AliHLTTPCCAOutTrack.h"
35 #include "AliHLTTPCCAParam.h"
36 #include "AliHLTTPCCATrackConvertor.h"
37 #include "AliHLTArray.h"
38
39 #include "AliHLTTPCSpacePointData.h"
40 #include "AliHLTTPCClusterDataFormat.h"
41 #include "AliHLTTPCCACompressedInputData.h"
42 #include "AliHLTTPCTransform.h"
43 #include "AliHLTTPCTrackSegmentData.h"
44 #include "AliHLTTPCTrackArray.h"
45 #include "AliHLTTPCTrackletDataFormat.h"
46 #include "AliHLTTPCDefinitions.h"
47 #include "AliExternalTrackParam.h"
48 #include "TStopwatch.h"
49 #include "TMath.h"
50 #include "AliCDBEntry.h"
51 #include "AliCDBManager.h"
52 #include "TObjString.h"
53 #include "TObjArray.h"
54 #include "AliHLTTPCCASliceOutput.h"
55 #include "AliHLTTPCCAClusterData.h"
56
57 const AliHLTComponentDataType AliHLTTPCCADefinitions::fgkTrackletsDataType = AliHLTComponentDataTypeInitializer( "CATRACKL", kAliHLTDataOriginTPC );
58
59 /** ROOT macro for the implementation of ROOT specific class methods */
60 ClassImp( AliHLTTPCCATrackerComponent )
61
62 AliHLTTPCCATrackerComponent::AliHLTTPCCATrackerComponent()
63     :
64     fTracker( NULL ),
65     fSolenoidBz( 0 ),
66     fMinNTrackClusters( 0 ),
67     fClusterZCut( 500. ),
68     fNeighboursSearchArea( 0 ), 
69     fClusterErrorCorrectionY(0), 
70     fClusterErrorCorrectionZ(0),
71     fFullTime( 0 ),
72     fRecoTime( 0 ),
73     fNEvents( 0 ),
74     fOutputTRAKSEGS( 0 )
75 {
76   // see header file for class documentation
77   // or
78   // refer to README to build package
79   // or
80   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
81 }
82
83 AliHLTTPCCATrackerComponent::AliHLTTPCCATrackerComponent( const AliHLTTPCCATrackerComponent& )
84     :
85     AliHLTProcessor(),
86     fTracker( NULL ),
87     fSolenoidBz( 0 ),
88     fMinNTrackClusters( 30 ),
89     fClusterZCut( 500. ),
90     fNeighboursSearchArea(0),
91     fClusterErrorCorrectionY(0), 
92     fClusterErrorCorrectionZ(0),
93     fFullTime( 0 ),
94     fRecoTime( 0 ),
95     fNEvents( 0 ),
96     fOutputTRAKSEGS( 0 )
97 {
98   // see header file for class documentation
99   HLTFatal( "copy constructor untested" );
100 }
101
102 AliHLTTPCCATrackerComponent& AliHLTTPCCATrackerComponent::operator=( const AliHLTTPCCATrackerComponent& )
103 {
104   // see header file for class documentation
105   HLTFatal( "assignment operator untested" );
106   return *this;
107 }
108
109 AliHLTTPCCATrackerComponent::~AliHLTTPCCATrackerComponent()
110 {
111   // see header file for class documentation
112   if (fTracker) delete fTracker;
113 }
114
115 //
116 // Public functions to implement AliHLTComponent's interface.
117 // These functions are required for the registration process
118 //
119
120 const char* AliHLTTPCCATrackerComponent::GetComponentID()
121 {
122   // see header file for class documentation
123   return "TPCCATracker";
124 }
125
126 void AliHLTTPCCATrackerComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list )
127 {
128   // see header file for class documentation
129   list.clear();
130   list.push_back( AliHLTTPCDefinitions::fgkClustersDataType );
131   list.push_back( AliHLTTPCCADefinitions::fgkCompressedInputDataType );
132 }
133
134 AliHLTComponentDataType AliHLTTPCCATrackerComponent::GetOutputDataType()
135 {
136   // see header file for class documentation
137   if ( fOutputTRAKSEGS ) return AliHLTTPCDefinitions::fgkTrackSegmentsDataType;
138   else return AliHLTTPCCADefinitions::fgkTrackletsDataType;
139 }
140
141 void AliHLTTPCCATrackerComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
142 {
143   // define guess for the output data size
144   constBase = 200;       // minimum size
145   inputMultiplier = 3.; // size relative to input
146 }
147
148 AliHLTComponent* AliHLTTPCCATrackerComponent::Spawn()
149 {
150   // see header file for class documentation
151   return new AliHLTTPCCATrackerComponent;
152 }
153
154 void AliHLTTPCCATrackerComponent::SetDefaultConfiguration()
155 {
156   // Set default configuration for the CA tracker component
157   // Some parameters can be later overwritten from the OCDB
158
159   fSolenoidBz = -5.00668;
160   fMinNTrackClusters = 0;
161   fClusterZCut = 500.;
162   fNeighboursSearchArea = 0;
163   fClusterErrorCorrectionY = 0;
164   fClusterErrorCorrectionZ = 0;
165   fOutputTRAKSEGS = 0;
166   fFullTime = 0;
167   fRecoTime = 0;
168   fNEvents = 0;
169 }
170
171 int AliHLTTPCCATrackerComponent::ReadConfigurationString(  const char* arguments )
172 {
173   // Set configuration parameters for the CA tracker component from the string
174
175   int iResult = 0;
176   if ( !arguments ) return iResult;
177
178   TString allArgs = arguments;
179   TString argument;
180   int bMissingParam = 0;
181
182   TObjArray* pTokens = allArgs.Tokenize( " " );
183
184   int nArgs =  pTokens ? pTokens->GetEntries() : 0;
185
186   for ( int i = 0; i < nArgs; i++ ) {
187     argument = ( ( TObjString* )pTokens->At( i ) )->GetString();
188     if ( argument.IsNull() ) continue;
189
190     if ( argument.CompareTo( "-solenoidBz" ) == 0 ) {
191       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
192       fSolenoidBz = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
193       HLTInfo( "Magnetic Field set to: %f", fSolenoidBz );
194       continue;
195     }
196
197     if ( argument.CompareTo( "-minNClustersOnTrack" ) == 0 ) {
198       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
199       fMinNTrackClusters = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
200       HLTInfo( "minNClustersOnTrack set to: %d", fMinNTrackClusters );
201       continue;
202     }
203
204     if ( argument.CompareTo( "-clusterZCut" ) == 0 ) {
205       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
206       fClusterZCut = TMath::Abs( ( ( TObjString* )pTokens->At( i ) )->GetString().Atof() );
207       HLTInfo( "ClusterZCut set to: %f", fClusterZCut );
208       continue;
209     }
210  
211    if ( argument.CompareTo( "-neighboursSearchArea" ) == 0 ) {
212       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
213       fNeighboursSearchArea = TMath::Abs( ( ( TObjString* )pTokens->At( i ) )->GetString().Atof() );
214       HLTInfo( "NeighboursSearchArea set to: %f", fNeighboursSearchArea );
215       continue;
216     }
217
218    if ( argument.CompareTo( "-errorCorrectionY" ) == 0 ) {
219      if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
220      fClusterErrorCorrectionY = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
221      HLTInfo( "Cluster Y error correction factor set to: %f", fClusterErrorCorrectionY );
222      continue;
223    }
224    
225    if ( argument.CompareTo( "-errorCorrectionZ" ) == 0 ) {
226      if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
227      fClusterErrorCorrectionZ = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
228      HLTInfo( "Cluster Z error correction factor set to: %f", fClusterErrorCorrectionZ );
229      continue;
230    }
231
232    if ( argument.CompareTo( "-outputTRAKSEGS" ) == 0 ) {
233       fOutputTRAKSEGS = 1;
234       HLTInfo( "The special output type \"TRAKSEGS\" is set" );
235       continue;
236     }
237
238     HLTError( "Unknown option \"%s\"", argument.Data() );
239     iResult = -EINVAL;
240   }
241   delete pTokens;
242
243   if ( bMissingParam ) {
244     HLTError( "Specifier missed for parameter \"%s\"", argument.Data() );
245     iResult = -EINVAL;
246   }
247
248   return iResult;
249 }
250
251
252 int AliHLTTPCCATrackerComponent::ReadCDBEntry( const char* cdbEntry, const char* chainId )
253 {
254   // see header file for class documentation
255
256   const char* defaultNotify = "";
257
258   if ( !cdbEntry ) {
259     cdbEntry = "HLT/ConfigTPC/TPCCATracker";
260     defaultNotify = " (default)";
261     chainId = 0;
262   }
263
264   HLTInfo( "configure from entry \"%s\"%s, chain id %s", cdbEntry, defaultNotify, ( chainId != NULL && chainId[0] != 0 ) ? chainId : "<none>" );
265   AliCDBEntry *pEntry = AliCDBManager::Instance()->Get( cdbEntry );//,GetRunNo());
266
267   if ( !pEntry ) {
268     HLTError( "cannot fetch object \"%s\" from CDB", cdbEntry );
269     return -EINVAL;
270   }
271
272   TObjString* pString = dynamic_cast<TObjString*>( pEntry->GetObject() );
273
274   if ( !pString ) {
275     HLTError( "configuration object \"%s\" has wrong type, required TObjString", cdbEntry );
276     return -EINVAL;
277   }
278
279   HLTInfo( "received configuration object string: \"%s\"", pString->GetString().Data() );
280
281   return  ReadConfigurationString( pString->GetString().Data() );
282 }
283
284
285 int AliHLTTPCCATrackerComponent::Configure( const char* cdbEntry, const char* chainId, const char *commandLine )
286 {
287   // Configure the component
288   // There are few levels of configuration,
289   // parameters which are set on one step can be overwritten on the next step
290
291   //* read hard-coded values
292
293   SetDefaultConfiguration();
294
295   //* read the default CDB entry
296
297   int iResult1 = ReadCDBEntry( NULL, chainId );
298
299   //* read magnetic field
300
301   int iResult2 = ReadCDBEntry( kAliHLTCDBSolenoidBz, chainId );
302
303   //* read the actual CDB entry if required
304
305   int iResult3 = ( cdbEntry ) ? ReadCDBEntry( cdbEntry, chainId ) : 0;
306
307   //* read extra parameters from input (if they are)
308
309   int iResult4 = 0;
310
311   if ( commandLine && commandLine[0] != '\0' ) {
312     HLTInfo( "received configuration string from HLT framework: \"%s\"", commandLine );
313     iResult4 = ReadConfigurationString( commandLine );
314   }
315
316   // Initialise the tracker here
317
318   return iResult1 ? iResult1 : ( iResult2 ? iResult2 : ( iResult3 ? iResult3 : iResult4 ) );
319 }
320
321
322
323 int AliHLTTPCCATrackerComponent::DoInit( int argc, const char** argv )
324 {
325   // Configure the CA tracker component
326
327   if ( fTracker ) return EINPROGRESS;
328
329
330   //fTracker = new AliHLTTPCCATrackerFramework();
331   //Do not initialize the TrackerFramework here since the CUDA framework is thread local and DoInit is called from different thread than DoEvent
332
333   TString arguments = "";
334   for ( int i = 0; i < argc; i++ ) {
335     if ( !arguments.IsNull() ) arguments += " ";
336     arguments += argv[i];
337   }
338
339   return Configure( NULL, NULL, arguments.Data() );
340 }
341
342
343 int AliHLTTPCCATrackerComponent::DoDeinit()
344 {
345   // see header file for class documentation
346   if (fTracker) delete fTracker;
347   fTracker = NULL;
348   return 0;
349 }
350
351
352
353 int AliHLTTPCCATrackerComponent::Reconfigure( const char* cdbEntry, const char* chainId )
354 {
355   // Reconfigure the component from OCDB .
356
357   return Configure( cdbEntry, chainId, NULL );
358 }
359
360
361 bool AliHLTTPCCATrackerComponent::CompareClusters( AliHLTTPCSpacePointData *a, AliHLTTPCSpacePointData *b )
362 {
363   //* Comparison function for sort clusters
364
365   if ( a->fPadRow < b->fPadRow ) return 1;
366   if ( a->fPadRow > b->fPadRow ) return 0;
367   return ( a->fZ < b->fZ );
368 }
369
370
371
372 int AliHLTTPCCATrackerComponent::DoEvent
373 (
374   const AliHLTComponentEventData& evtData,
375   const AliHLTComponentBlockData* blocks,
376   AliHLTComponentTriggerData& /*trigData*/,
377   AliHLTUInt8_t* outputPtr,
378   AliHLTUInt32_t& size,
379   vector<AliHLTComponentBlockData>& outputBlocks )
380 {
381   //* process event
382
383   AliHLTUInt32_t maxBufferSize = size;
384   size = 0; // output size
385
386   if ( GetFirstInputBlock( kAliHLTDataTypeSOR ) || GetFirstInputBlock( kAliHLTDataTypeEOR ) ) {
387     return 0;
388   }
389
390   TStopwatch timer;
391
392   // Event reconstruction in one TPC slice with CA Tracker
393
394   //Logging( kHLTLogWarning, "HLT::TPCCATracker::DoEvent", "DoEvent", "CA::DoEvent()" );
395   if ( evtData.fBlockCnt <= 0 ) {
396     HLTWarning( "no blocks in event" );
397     return 0;
398   }
399
400   const AliHLTComponentBlockData* iter = NULL;
401   unsigned long ndx;
402
403   // Determine the slice number
404
405   //Find min and max slice number with now slice missing in between (for default output)
406   int minslice = -1, maxslice = -1;
407   int slice = -1;
408   {
409     std::vector<int> slices;
410     std::vector<int>::iterator slIter;
411     std::vector<unsigned> sliceCnts;
412     std::vector<unsigned>::iterator slCntIter;
413
414     for ( ndx = 0; ndx < evtData.fBlockCnt; ndx++ ) {
415       iter = blocks + ndx;
416       if ( iter->fDataType != AliHLTTPCDefinitions::fgkClustersDataType
417            && iter->fDataType != AliHLTTPCCADefinitions::fgkCompressedInputDataType
418            ) continue;
419
420       slice = AliHLTTPCDefinitions::GetMinSliceNr( *iter );
421           if (slice < minslice || minslice == -1) minslice = slice;
422           if (slice > maxslice) maxslice = slice;
423
424       bool found = 0;
425       slCntIter = sliceCnts.begin();
426       for ( slIter = slices.begin(); slIter != slices.end(); slIter++, slCntIter++ ) {
427         if ( *slIter == slice ) {
428           found = kTRUE;
429           break;
430         }
431       }
432       if ( !found ) {
433         slices.push_back( slice );
434         sliceCnts.push_back( 1 );
435       } else *slCntIter++;
436     }
437
438           if ( slices.size() == 0 ) {
439                 HLTWarning( "no slices found in event" );
440                 return 0;
441           }
442
443
444     // Determine slice number to really use. (for obsolete output)
445     if ( slices.size() > 1 ) {
446                 Logging( fOutputTRAKSEGS ? kHLTLogError : kHLTLogDebug, "HLT::TPCSliceTracker::DoEvent", "Multiple slices found in event",
447                "Multiple slice numbers found in event 0x%08lX (%lu). Determining maximum occuring slice number...",
448                evtData.fEventID, evtData.fEventID );
449       unsigned maxCntSlice = 0;
450       slCntIter = sliceCnts.begin();
451       for ( slIter = slices.begin(); slIter != slices.end(); slIter++, slCntIter++ ) {
452         Logging( fOutputTRAKSEGS ? kHLTLogError : kHLTLogDebug, "HLT::TPCSliceTracker::DoEvent", "Multiple slices found in event",
453                  "Slice %lu found %lu times.", *slIter, *slCntIter );
454         if ( fOutputTRAKSEGS && maxCntSlice < *slCntIter ) {
455           maxCntSlice = *slCntIter;
456           slice = *slIter;
457         }
458       }
459       if (fOutputTRAKSEGS)
460       {
461             Logging( kHLTLogError, "HLT::TPCSliceTracker::DoEvent", "Multiple slices found in event",
462                "Using slice %lu.", slice );
463         }
464     } else if ( slices.size() > 0 ) {
465       slice = *( slices.begin() );
466     }
467
468
469           if (fOutputTRAKSEGS)
470           {
471                   minslice = maxslice = slice;
472           }
473           else
474           {
475                   for (int islice = minslice;islice <= maxslice;islice++)
476                   {
477                           bool found = false;
478                           for(slIter = slices.begin(); slIter != slices.end();slIter++)
479                           {
480                                   if (*slIter == islice)
481                                   {
482                                           found = true;
483                                           break;
484                                   }
485                           }
486                           if (!found)
487                           {
488                                   maxslice = islice - 1;
489                                   break;
490                           }
491                   }
492          }
493   }
494
495   if ( !fTracker ) fTracker = new AliHLTTPCCATrackerFramework;
496
497   int slicecount = maxslice + 1 - minslice;
498   if (slicecount > fTracker->MaxSliceCount())
499   {
500         maxslice = minslice + (slicecount = fTracker->MaxSliceCount());
501   }
502   int nClustersTotalSum = 0;
503   AliHLTTPCCAClusterData* clusterData = new AliHLTTPCCAClusterData[slicecount];
504
505
506   // min and max patch numbers and row numbers
507   int* slicerow = new int[slicecount * 2];
508   int* sliceminPatch = new int[slicecount];
509   int* slicemaxPatch = new int[slicecount];
510   memset(slicerow, 0, slicecount * 2 * sizeof(int));
511   for (int i = 0;i < slicecount;i++)
512   {
513           sliceminPatch[i] = 100;
514           slicemaxPatch[i] = -1;
515   }
516
517   //Prepare everything for all slices
518
519   for (int islice = 0;islice < slicecount;islice++)
520   {
521           slice = minslice + islice;
522
523           // Initialize the tracker slice
524           {
525                 int iSec = slice;
526                 float inRmin = 83.65;
527                 //    float inRmax = 133.3;
528                 //    float outRmin = 133.5;
529                 float outRmax = 247.7;
530                 float plusZmin = 0.0529937;
531                 float plusZmax = 249.778;
532                 float minusZmin = -249.645;
533                 float minusZmax = -0.0799937;
534                 float dalpha = 0.349066;
535                 float alpha = 0.174533 + dalpha * iSec;
536
537                 bool zPlus = ( iSec < 18 );
538                 float zMin =  zPlus ? plusZmin : minusZmin;
539                 float zMax =  zPlus ? plusZmax : minusZmax;
540                 //TPCZmin = -249.645, ZMax = 249.778
541                 //    float rMin =  inRmin;
542                 //    float rMax =  outRmax;
543                 int nRows = AliHLTTPCTransform::GetNRows();
544
545                 float padPitch = 0.4;
546                 float sigmaZ = 0.228808;
547
548                 float *rowX = new float [nRows];
549                 for ( int irow = 0; irow < nRows; irow++ ) {
550                   rowX[irow] = AliHLTTPCTransform::Row2X( irow );
551                 }
552
553                 AliHLTTPCCAParam param;
554
555                 param.Initialize( iSec, nRows, rowX, alpha, dalpha,
556                                                   inRmin, outRmax, zMin, zMax, padPitch, sigmaZ, fSolenoidBz );
557                 param.SetHitPickUpFactor( 2 );
558                 if( fNeighboursSearchArea>0 ) param.SetNeighboursSearchArea( fNeighboursSearchArea );
559                 if( fClusterErrorCorrectionY>1.e-4 ) param.SetClusterError2CorrectionY( fClusterErrorCorrectionY*fClusterErrorCorrectionY );
560                 if( fClusterErrorCorrectionZ>1.e-4 ) param.SetClusterError2CorrectionZ( fClusterErrorCorrectionZ*fClusterErrorCorrectionZ );
561                 param.Update();
562                 fTracker->InitializeSliceParam( slice, param );
563                 delete[] rowX;
564           }
565
566           // total n Hits
567           int nClustersTotal = 0;
568
569           // sort patches
570           std::vector<unsigned long> patchIndices;
571
572           for ( ndx = 0; ndx < evtData.fBlockCnt; ndx++ ) {
573                 iter = blocks + ndx;
574                 if ( slice != AliHLTTPCDefinitions::GetMinSliceNr( *iter ) ) continue;
575                 if ( iter->fDataType == AliHLTTPCDefinitions::fgkClustersDataType ){
576                   AliHLTTPCClusterData* inPtrSP = ( AliHLTTPCClusterData* )( iter->fPtr );
577                   nClustersTotal += inPtrSP->fSpacePointCnt;
578                 } 
579                 else if ( iter->fDataType == AliHLTTPCCADefinitions::fgkCompressedInputDataType){
580                   const AliHLTUInt8_t * inPtr =  (const AliHLTUInt8_t *)iter->fPtr;
581                   while( inPtr< ((const AliHLTUInt8_t *) iter->fPtr) + iter->fSize ){
582                     AliHLTTPCCACompressedClusterRow *row = (AliHLTTPCCACompressedClusterRow*)inPtr;
583                     nClustersTotal+= row->fNClusters;     
584                     inPtr = (const AliHLTUInt8_t *)(row->fClusters+row->fNClusters);
585                   }
586                 }
587                 else continue;
588
589                 int patch = AliHLTTPCDefinitions::GetMinPatchNr( *iter );
590                 if ( sliceminPatch[islice] > patch ) {
591                   sliceminPatch[islice] = patch;
592                   slicerow[2 * islice + 0] = AliHLTTPCTransform::GetFirstRow( patch );
593                 }
594                 if ( slicemaxPatch[islice] < patch ) {
595                   slicemaxPatch[islice] = patch;
596                   slicerow[2 * islice + 1] = AliHLTTPCTransform::GetLastRow( patch );
597                 }
598                 std::vector<unsigned long>::iterator pIter = patchIndices.begin();
599                 while ( pIter != patchIndices.end() && AliHLTTPCDefinitions::GetMinPatchNr( blocks[*pIter] ) < patch ) {
600                   pIter++;
601                 }
602                 patchIndices.insert( pIter, ndx );
603           }
604
605
606           // pass event to CA Tracker
607
608
609           Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Reading hits",
610                            "Total %d hits to read for slice %d", nClustersTotal, slice );
611
612
613           clusterData[islice].StartReading( slice, nClustersTotal );
614
615           for ( std::vector<unsigned long>::iterator pIter = patchIndices.begin(); pIter != patchIndices.end(); pIter++ ) {
616                 ndx = *pIter;
617                 iter = blocks + ndx;
618                 int patch = AliHLTTPCDefinitions::GetMinPatchNr( *iter );
619                 int nPatchClust = 0;
620                 if ( iter->fDataType == AliHLTTPCDefinitions::fgkClustersDataType ){
621                   AliHLTTPCClusterData* inPtrSP = ( AliHLTTPCClusterData* )( iter->fPtr );
622                   nPatchClust = inPtrSP->fSpacePointCnt;
623                   for ( unsigned int i = 0; i < inPtrSP->fSpacePointCnt; i++ ) {
624                     AliHLTTPCSpacePointData *c = &( inPtrSP->fSpacePoints[i] );
625                     if ( CAMath::Abs( c->fZ ) > fClusterZCut ) continue;
626                     if ( c->fPadRow > 159 ) {
627                       HLTError( "Wrong TPC cluster with row number %d received", c->fPadRow );
628                       continue;
629                     }
630                     clusterData[islice].ReadCluster( c->fID, c->fPadRow, c->fX, c->fY, c->fZ, c->fCharge );
631                   }           
632                 } 
633                 else if ( iter->fDataType == AliHLTTPCCADefinitions::fgkCompressedInputDataType){
634                   const AliHLTUInt8_t * inPtr = (const AliHLTUInt8_t *)iter->fPtr;
635                   nPatchClust=0;
636                   while( inPtr< ((const AliHLTUInt8_t *)iter->fPtr) + iter->fSize ){
637                     AliHLTTPCCACompressedClusterRow *row = (AliHLTTPCCACompressedClusterRow*)inPtr;
638                     UInt_t id = row->fSlicePatchRowID;
639                     UInt_t jslice = id>>10;       
640                     UInt_t jpatch = (id>>6) & 0x7;
641                     UInt_t jrow   =  id     & 0x3F;     
642                     jrow+= AliHLTTPCTransform::GetFirstRow( jpatch );
643                     Double_t rowX = AliHLTTPCTransform::Row2X( jrow );
644                     //cout<<"Read row: s "<<jslice<<" p "<<jpatch<<" r "<<jrow<<" x "<<row->fX<<" nclu "<<row->fNClusters<<" :"<<endl;
645                     if( jrow > 159 ) {
646                       HLTError( "Wrong TPC cluster with row number %d received", jrow );
647                       continue;
648                     }
649                     for ( unsigned int i = 0; i < row->fNClusters; i++ ) {
650                       AliHLTTPCCACompressedCluster *c = &( row->fClusters[i] );
651                       
652                       UInt_t ix0 = c->fP0 >>24;
653                       UInt_t ix1 = c->fP1 >>24;
654                       Double_t x = (ix1<<8) + ix0;
655                       Double_t y = c->fP0 & 0x00FFFFFF;
656                       Double_t z = c->fP1 & 0x00FFFFFF;
657                       x = (x - 32768.)*1.e-4 + rowX;
658                       y = (y - 8388608.)*1.e-4;
659                       z = (z - 8388608.)*1.e-4;
660                       
661                       UInt_t cluId = nPatchClust + ((jslice&0x7f)<<25)+((jpatch&0x7)<<22);
662                       //cout<<"clu "<<i<<": "<<x<<" "<<y<<" "<<z<<" "<<cluId<<endl;
663                       if ( CAMath::Abs( z ) <= fClusterZCut ){
664                         clusterData[islice].ReadCluster( cluId, jrow, x, y, z, 0 );
665                       }
666                       nPatchClust++;              
667                     }
668                     inPtr = (const AliHLTUInt8_t *)(row->fClusters+row->fNClusters);
669                   }
670                 }
671                 Logging( kHLTLogInfo, "HLT::TPCCATracker::DoEvent", "Reading hits",
672                          "Read %d hits for slice %d - patch %d", nPatchClust, slice, patch );
673           }
674
675           clusterData[islice].FinishReading();
676           nClustersTotalSum += nClustersTotal;
677   }
678
679   //Prepare Output
680   AliHLTTPCCASliceOutput::outputControlStruct outputControl;
681   //Set tracker output so tracker does not have to output both formats!
682   outputControl.fObsoleteOutput = fOutputTRAKSEGS;
683   outputControl.fDefaultOutput = !fOutputTRAKSEGS;
684
685   //For new output we can write directly to output buffer
686   outputControl.fOutputPtr = fOutputTRAKSEGS ? NULL : (char*) outputPtr;
687   outputControl.fOutputMaxSize = maxBufferSize;
688
689   AliHLTTPCCASliceOutput** sliceOutput = new AliHLTTPCCASliceOutput*[slicecount];
690   memset(sliceOutput, 0, slicecount * sizeof(AliHLTTPCCASliceOutput*));
691
692   // reconstruct the event
693   TStopwatch timerReco;
694   fTracker->SetOutputControl(&outputControl);
695   fTracker->ProcessSlices(minslice, slicecount, clusterData, sliceOutput);
696   timerReco.Stop();
697   
698   int ret = 0;
699   unsigned int mySize = 0;
700   int ntracks = 0;
701   int error = 0;
702
703   for (int islice = 0;islice < slicecount;islice++)
704   {
705           slice = minslice + islice;
706
707           if (sliceOutput[islice])
708           {
709                   // write reconstructed tracks
710
711                   if ( fOutputTRAKSEGS ) {
712
713                   Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Reconstruct",
714                                    "%d tracks found for slice %d", sliceOutput[islice]->NOutTracks(), slice );
715
716                         ntracks = sliceOutput[islice]->NOutTracks();
717
718                         AliHLTTPCTrackletData* outPtr = ( AliHLTTPCTrackletData* )( outputPtr );
719
720                         AliHLTTPCTrackSegmentData* currOutTracklet = outPtr->fTracklets;
721
722                         mySize =   ( ( AliHLTUInt8_t * )currOutTracklet ) -  ( ( AliHLTUInt8_t * )outputPtr );
723
724                         outPtr->fTrackletCnt = 0;
725
726                         for ( int itr = 0; itr < ntracks; itr++ ) {
727
728                           AliHLTTPCCAOutTrack &t = sliceOutput[islice]->OutTracks()[itr];
729
730                           //Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Wrtite output","track %d with %d hits", itr, t.NHits());
731
732                           if ( t.NHits() < fMinNTrackClusters ) continue;
733
734                           // calculate output track size
735
736                           unsigned int dSize = sizeof( AliHLTTPCTrackSegmentData ) + t.NHits() * sizeof( unsigned int );
737
738                           if ( mySize + dSize > maxBufferSize ) {
739                                 HLTWarning( "Output buffer size exceed (buffer size %d, current size %d), %d tracks are not stored", maxBufferSize, mySize, ntracks - itr + 1 );
740                                 ret = -ENOSPC;
741                                 error = 1;
742                                 break;
743                           }
744
745                           // convert CA track parameters to HLT Track Segment
746
747                           int iFirstRow = 1000;
748                           int iLastRow = -1;
749                           int iFirstHit = sliceOutput[islice]->OutTrackHit(t.FirstHitRef());
750                           int iLastHit = iFirstHit;
751                           for ( int ih = 0; ih < t.NHits(); ih++ ) {
752                                 int hitID = sliceOutput[islice]->OutTrackHit(t.FirstHitRef() + ih);
753                                 int iRow = clusterData[islice].RowNumber( hitID );
754                                 if ( iRow < iFirstRow ) {  iFirstRow = iRow; iFirstHit = hitID; }
755                                 if ( iRow > iLastRow ) { iLastRow = iRow; iLastHit = hitID; }
756                           }
757
758                           AliHLTTPCCATrackParam par = t.StartPoint();
759
760                           par.TransportToX( clusterData[islice].X( iFirstHit ), .99 );
761
762                           AliExternalTrackParam tp;
763                           AliHLTTPCCATrackConvertor::GetExtParam( par, tp, 0 );
764
765                           currOutTracklet->fX = tp.GetX();
766                           currOutTracklet->fY = tp.GetY();
767                           currOutTracklet->fZ = tp.GetZ();
768                           currOutTracklet->fCharge = ( int ) tp.GetSign();
769                           currOutTracklet->fPt = TMath::Abs( tp.GetSignedPt() );
770                           float snp =  tp.GetSnp() ;
771                           if ( snp > .999 ) snp = .999;
772                           if ( snp < -.999 ) snp = -.999;
773                           currOutTracklet->fPsi = TMath::ASin( snp );
774                           currOutTracklet->fTgl = tp.GetTgl();
775
776                           currOutTracklet->fY0err = tp.GetSigmaY2();
777                           currOutTracklet->fZ0err = tp.GetSigmaZ2();
778                           float h = -currOutTracklet->fPt * currOutTracklet->fPt;
779                           currOutTracklet->fPterr = h * h * tp.GetSigma1Pt2();
780                           h = 1. / TMath::Sqrt( 1 - snp * snp );
781                           currOutTracklet->fPsierr = h * h * tp.GetSigmaSnp2();
782                           currOutTracklet->fTglerr = tp.GetSigmaTgl2();
783
784                           if ( par.TransportToX( clusterData[islice].X( iLastHit ), .99 ) ) {
785                                 currOutTracklet->fLastX = par.GetX();
786                                 currOutTracklet->fLastY = par.GetY();
787                                 currOutTracklet->fLastZ = par.GetZ();
788                           } else {
789                                 currOutTracklet->fLastX = clusterData[islice].X( iLastHit );
790                                 currOutTracklet->fLastY = clusterData[islice].Y( iLastHit );
791                                 currOutTracklet->fLastZ = clusterData[islice].Z( iLastHit );
792                           }
793                           //if( currOutTracklet->fLastX<10. ) {
794                           //HLTError("CA last point: hitxyz=%f,%f,%f, track=%f,%f,%f, tracklet=%f,%f,%f, nhits=%d",clusterData[islice].X( iLastHit ),clusterData[islice].Y( iLastHit],clusterData[islice].Z( iLastHit],
795                           //par.GetX(), par.GetY(),par.GetZ(),currOutTracklet->fLastX,currOutTracklet->fLastY ,currOutTracklet->fLastZ, t.NHits());
796                           //}
797                 #ifdef INCLUDE_TPC_HOUGH
798                 #ifdef ROWHOUGHPARAMS
799                           currOutTracklet->fTrackID = 0;
800                           currOutTracklet->fRowRange1 = clusterData[islice].RowNumber( iFirstHit );
801                           currOutTracklet->fRowRange2 = clusterData[islice].RowNumber( iLastHit );
802                           currOutTracklet->fSector = slice;
803                           currOutTracklet->fPID = 211;
804                 #endif
805                 #endif // INCLUDE_TPC_HOUGH
806
807
808                           currOutTracklet->fNPoints = t.NHits();
809
810                           for ( int i = 0; i < t.NHits(); i++ ) {
811                                 currOutTracklet->fPointIDs[i] = clusterData[islice].Id( sliceOutput[islice]->OutTrackHit(t.FirstHitRef()+i) );
812                           }
813
814                           currOutTracklet = ( AliHLTTPCTrackSegmentData* )( ( Byte_t * )currOutTracklet + dSize );
815                           mySize += dSize;
816                           outPtr->fTrackletCnt++;
817                         }
818
819                   } else { // default output type
820
821                   Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Reconstruct",
822                                    "%d tracks found for slice %d", sliceOutput[islice]->NTracks(), slice );
823
824                           mySize += sliceOutput[islice]->OutputMemorySize();
825                           ntracks += sliceOutput[islice]->NTracks();
826                   }
827           }
828           else
829           {
830                   HLTWarning( "Output buffer size exceed (buffer size %d, current size %d), tracks are not stored", maxBufferSize, mySize );
831                   mySize = 0;
832                   ret = -ENOSPC;
833                   ntracks = 0;
834                   error = 1;
835                   break;
836           }
837   }
838
839   size = 0;
840   if (error == 0)
841   {
842           for (int islice = 0;islice < slicecount;islice++)
843           {
844                   slice = minslice + islice;
845                   if (!fOutputTRAKSEGS) mySize = sliceOutput[islice]->OutputMemorySize();
846                   if (mySize > 0)
847                   {
848                         AliHLTComponentBlockData bd;
849                         FillBlockData( bd );
850                         bd.fOffset = fOutputTRAKSEGS ? 0 : ((char*) sliceOutput[islice] - (char*) outputPtr);
851                         bd.fSize = mySize;
852                         bd.fSpecification = AliHLTTPCDefinitions::EncodeDataSpecification( slice, slice, sliceminPatch[islice], slicemaxPatch[islice] );
853                         bd.fDataType = GetOutputDataType();
854                         outputBlocks.push_back( bd );
855                         size += mySize;
856                   }
857           }
858   }
859
860   //No longer needed
861   delete[] clusterData;
862   //These are only temporary pointers to the output and no longer needed
863   delete[] sliceOutput;
864
865   timer.Stop();
866
867   fFullTime += timer.RealTime();
868   fRecoTime += timerReco.RealTime();
869   fNEvents++;
870
871   // Set log level to "Warning" for on-line system monitoring
872   int hz = ( int ) ( fFullTime > 1.e-10 ? fNEvents / fFullTime : 100000 );
873   int hz1 = ( int ) ( fRecoTime > 1.e-10 ? fNEvents / fRecoTime : 100000 );
874   //Min and Max Patch are taken for first slice processed...
875   HLTInfo( "CATracker slices %d-%d: output %d tracks;  input %d clusters, patches %d..%d, rows %d..%d; time: full %d / reco %d Hz",
876            minslice, maxslice, ntracks, nClustersTotalSum, sliceminPatch[0], slicemaxPatch[0], slicerow[0], slicerow[1], hz, hz1 );
877
878   return ret;
879 }
880
881