2 // **************************************************************************
3 // This file is property of and copyright by the ALICE HLT Project *
4 // ALICE Experiment at CERN, All rights reserved. *
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. *
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. *
18 //***************************************************************************
21 ///////////////////////////////////////////////////////////////////////////////
23 // a TPC tracker processing component for the HLT based on CA by Ivan Kisel //
25 ///////////////////////////////////////////////////////////////////////////////
31 #include "AliHLTTPCCATrackerComponent.h"
32 #include "AliHLTTPCTransform.h"
33 #include "AliHLTTPCCATracker.h"
34 #include "AliHLTTPCCAOutTrack.h"
35 #include "AliHLTTPCCAParam.h"
36 #include "AliHLTTPCCATrackConvertor.h"
37 #include "AliHLTArray.h"
39 #include "AliHLTTPCSpacePointData.h"
40 #include "AliHLTTPCClusterDataFormat.h"
41 #include "AliHLTTPCTransform.h"
42 #include "AliHLTTPCTrackSegmentData.h"
43 #include "AliHLTTPCTrackArray.h"
44 #include "AliHLTTPCTrackletDataFormat.h"
45 #include "AliHLTTPCDefinitions.h"
46 #include "AliExternalTrackParam.h"
47 #include "TStopwatch.h"
49 #include "AliCDBEntry.h"
50 #include "AliCDBManager.h"
51 #include "TObjString.h"
52 #include "TObjArray.h"
53 #include "AliHLTTPCCASliceOutput.h"
54 #include "AliHLTTPCCAClusterData.h"
56 const AliHLTComponentDataType AliHLTTPCCADefinitions::fgkTrackletsDataType = AliHLTComponentDataTypeInitializer( "CATRACKL", kAliHLTDataOriginTPC );
58 /** ROOT macro for the implementation of ROOT specific class methods */
59 ClassImp( AliHLTTPCCATrackerComponent )
61 AliHLTTPCCATrackerComponent::AliHLTTPCCATrackerComponent()
65 fMinNTrackClusters( 0 ),
72 // see header file for class documentation
74 // refer to README to build package
76 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
79 AliHLTTPCCATrackerComponent::AliHLTTPCCATrackerComponent( const AliHLTTPCCATrackerComponent& )
84 fMinNTrackClusters( 30 ),
91 // see header file for class documentation
92 HLTFatal( "copy constructor untested" );
95 AliHLTTPCCATrackerComponent& AliHLTTPCCATrackerComponent::operator=( const AliHLTTPCCATrackerComponent& )
97 // see header file for class documentation
98 HLTFatal( "assignment operator untested" );
102 AliHLTTPCCATrackerComponent::~AliHLTTPCCATrackerComponent()
104 // see header file for class documentation
109 // Public functions to implement AliHLTComponent's interface.
110 // These functions are required for the registration process
113 const char* AliHLTTPCCATrackerComponent::GetComponentID()
115 // see header file for class documentation
116 return "TPCCATracker";
119 void AliHLTTPCCATrackerComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list )
121 // see header file for class documentation
123 list.push_back( AliHLTTPCDefinitions::fgkClustersDataType );
126 AliHLTComponentDataType AliHLTTPCCATrackerComponent::GetOutputDataType()
128 // see header file for class documentation
129 if ( fNewOutputType ) return AliHLTTPCCADefinitions::fgkTrackletsDataType;
130 else return AliHLTTPCDefinitions::fgkTrackSegmentsDataType;
133 void AliHLTTPCCATrackerComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
135 // define guess for the output data size
136 constBase = 200; // minimum size
137 inputMultiplier = 0.5; // size relative to input
140 AliHLTComponent* AliHLTTPCCATrackerComponent::Spawn()
142 // see header file for class documentation
143 return new AliHLTTPCCATrackerComponent;
146 int AliHLTTPCCATrackerComponent::DoInit( int argc, const char** argv )
148 // Configure the CA tracker component
151 fMinNTrackClusters = 0;
158 if ( fTracker ) return EINPROGRESS;
159 fTracker = new AliHLTTPCCATracker();
163 TString arguments = "";
164 for ( int i = 0; i < argc; i++ ) {
165 TString argument = argv[i];
166 if ( !arguments.IsNull() ) arguments += " ";
167 arguments += argument;
169 if ( !arguments.IsNull() ) {
170 iResult = Configure( arguments.Data() );
172 iResult = Reconfigure( NULL, NULL );
178 int AliHLTTPCCATrackerComponent::DoDeinit()
180 // see header file for class documentation
181 if ( fTracker ) delete fTracker;
186 int AliHLTTPCCATrackerComponent::Reconfigure( const char* /*cdbEntry*/, const char* /*chainId*/ )
188 // see header file for class documentation
190 HLTInfo( "TODO: dummy Reconfigure() method" );
197 const char* path="HLT/ConfigTPC/CATrackerComponent";
198 const char* defaultNotify="";
201 defaultNotify=" (default)";
204 HLTInfo("reconfigure from entry %s%s, chain id %s", path, defaultNotify,(chainId!=NULL && chainId[0]!=0)?chainId:"<none>");
205 AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(path);//,GetRunNo());
207 TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
209 HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data());
210 iResult=Configure(pString->GetString().Data());
212 HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
215 HLTError("cannot fetch object \"%s\" from CDB", path);
219 const char* pathBField=kAliHLTCDBSolenoidBz;
222 HLTInfo("reconfigure B-Field from entry %s, chain id %s", path,(chainId!=NULL && chainId[0]!=0)?chainId:"<none>");
223 AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(pathBField);//,GetRunNo());
225 TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
227 HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data());
228 iResult=Configure(pString->GetString().Data());
230 HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
233 HLTError("cannot fetch object \"%s\" from CDB", path);
241 bool AliHLTTPCCATrackerComponent::CompareClusters( AliHLTTPCSpacePointData *a, AliHLTTPCSpacePointData *b )
243 //* Comparison function for sorting clusters
244 if ( a->fPadRow < b->fPadRow ) return 1;
245 if ( a->fPadRow > b->fPadRow ) return 0;
246 return ( a->fZ < b->fZ );
250 int AliHLTTPCCATrackerComponent::Configure( const char* arguments )
255 if ( !arguments ) return iResult;
257 TString allArgs = arguments;
259 int bMissingParam = 0;
261 TObjArray* pTokens = allArgs.Tokenize( " " );
263 int nArgs = pTokens ? pTokens->GetEntries() : 0;
265 for ( int i = 0; i < nArgs; i++ ) {
266 argument = ( ( TObjString* )pTokens->At( i ) )->GetString();
267 if ( argument.IsNull() ) {
268 } else if ( argument.CompareTo( "-solenoidBz" ) == 0 ) {
269 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
270 fSolenoidBz = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
271 HLTInfo( "Magnetic Field set to: %f", fSolenoidBz );
272 } else if ( argument.CompareTo( "-minNClustersOnTrack" ) == 0 ||
273 argument.CompareTo( "-minNTrackClusters" ) == 0 ) {
274 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
275 fMinNTrackClusters = ( ( TObjString* )pTokens->At( i ) )->GetString().Atoi();
276 HLTInfo( "minNClustersOnTrack set to: %d", fMinNTrackClusters );
277 } else if ( argument.CompareTo( "-clusterZCut" ) == 0 ) {
278 if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
279 fClusterZCut = TMath::Abs( ( ( TObjString* )pTokens->At( i ) )->GetString().Atof() );
280 HLTInfo( "ClusterZCut set to: %f", fClusterZCut );
281 } else if ( argument.CompareTo( "-newOutputType" ) == 0 ) {
283 HLTInfo( "NewOutputType is set" );
285 HLTError( "Unknown option %s ", argument.Data() );
291 if ( bMissingParam ) {
292 HLTError( "Specifier missed for %s", argument.Data() );
302 int AliHLTTPCCATrackerComponent::DoEvent
304 const AliHLTComponentEventData& evtData,
305 const AliHLTComponentBlockData* blocks,
306 AliHLTComponentTriggerData& /*trigData*/,
307 AliHLTUInt8_t* outputPtr,
308 AliHLTUInt32_t& size,
309 vector<AliHLTComponentBlockData>& outputBlocks )
312 AliHLTUInt32_t maxBufferSize = size;
313 size = 0; // output size
315 if ( GetFirstInputBlock( kAliHLTDataTypeSOR ) || GetFirstInputBlock( kAliHLTDataTypeEOR ) ) {
321 // Event reconstruction in one TPC slice with CA Tracker
323 //Logging( kHLTLogWarning, "HLT::TPCCATracker::DoEvent", "DoEvent", "CA::DoEvent()" );
324 if ( evtData.fBlockCnt <= 0 ) {
325 HLTWarning( "no blocks in event" );
329 const AliHLTComponentBlockData* iter = NULL;
331 AliHLTTPCClusterData* inPtrSP;
333 // Determine the slice number
337 std::vector<int> slices;
338 std::vector<int>::iterator slIter;
339 std::vector<unsigned> sliceCnts;
340 std::vector<unsigned>::iterator slCntIter;
342 for ( ndx = 0; ndx < evtData.fBlockCnt; ndx++ ) {
344 if ( iter->fDataType != AliHLTTPCDefinitions::fgkClustersDataType ) continue;
346 slice = AliHLTTPCDefinitions::GetMinSliceNr( *iter );
349 slCntIter = sliceCnts.begin();
350 for ( slIter = slices.begin(); slIter != slices.end(); slIter++, slCntIter++ ) {
351 if ( *slIter == slice ) {
357 slices.push_back( slice );
358 sliceCnts.push_back( 1 );
363 // Determine slice number to really use.
364 if ( slices.size() > 1 ) {
365 Logging( kHLTLogError, "HLT::TPCSliceTracker::DoEvent", "Multiple slices found in event",
366 "Multiple slice numbers found in event 0x%08lX (%lu). Determining maximum occuring slice number...",
367 evtData.fEventID, evtData.fEventID );
368 unsigned maxCntSlice = 0;
369 slCntIter = sliceCnts.begin();
370 for ( slIter = slices.begin(); slIter != slices.end(); slIter++, slCntIter++ ) {
371 Logging( kHLTLogError, "HLT::TPCSliceTracker::DoEvent", "Multiple slices found in event",
372 "Slice %lu found %lu times.", *slIter, *slCntIter );
373 if ( maxCntSlice < *slCntIter ) {
374 maxCntSlice = *slCntIter;
378 Logging( kHLTLogError, "HLT::TPCSliceTracker::DoEvent", "Multiple slices found in event",
379 "Using slice %lu.", slice );
380 } else if ( slices.size() > 0 ) {
381 slice = *( slices.begin() );
386 HLTWarning( "no slices found in event" );
391 // Initialize the tracker
395 if ( !fTracker ) fTracker = new AliHLTTPCCATracker;
397 float inRmin = 83.65;
398 // float inRmax = 133.3;
399 // float outRmin = 133.5;
400 float outRmax = 247.7;
401 float plusZmin = 0.0529937;
402 float plusZmax = 249.778;
403 float minusZmin = -249.645;
404 float minusZmax = -0.0799937;
405 float dalpha = 0.349066;
406 float alpha = 0.174533 + dalpha * iSec;
408 bool zPlus = ( iSec < 18 );
409 float zMin = zPlus ? plusZmin : minusZmin;
410 float zMax = zPlus ? plusZmax : minusZmax;
411 //TPCZmin = -249.645, ZMax = 249.778
412 // float rMin = inRmin;
413 // float rMax = outRmax;
414 int nRows = AliHLTTPCTransform::GetNRows();
416 float padPitch = 0.4;
417 float sigmaZ = 0.228808;
419 float *rowX = new float [nRows];
420 for ( int irow = 0; irow < nRows; irow++ ) {
421 rowX[irow] = AliHLTTPCTransform::Row2X( irow );
424 const double kCLight = 0.000299792458;
426 AliHLTTPCCAParam param;
428 param.Initialize( iSec, nRows, rowX, alpha, dalpha,
429 inRmin, outRmax, zMin, zMax, padPitch, sigmaZ, fSolenoidBz*kCLight );
430 param.SetHitPickUpFactor( 2 );
432 fTracker->Initialize( param );
437 // min and max patch numbers and row numbers
440 int minPatch = 100, maxPatch = -1;
444 int nClustersTotal = 0;
448 std::vector<unsigned long> patchIndices;
450 for ( ndx = 0; ndx < evtData.fBlockCnt; ndx++ ) {
452 if ( iter->fDataType != AliHLTTPCDefinitions::fgkClustersDataType ) continue;
453 if ( slice != AliHLTTPCDefinitions::GetMinSliceNr( *iter ) ) continue;
454 inPtrSP = ( AliHLTTPCClusterData* )( iter->fPtr );
455 nClustersTotal += inPtrSP->fSpacePointCnt;
456 int patch = AliHLTTPCDefinitions::GetMinPatchNr( *iter );
457 if ( minPatch > patch ) {
459 row[0] = AliHLTTPCTransform::GetFirstRow( patch );
461 if ( maxPatch < patch ) {
463 row[1] = AliHLTTPCTransform::GetLastRow( patch );
465 std::vector<unsigned long>::iterator pIter = patchIndices.begin();
466 while ( pIter != patchIndices.end() && AliHLTTPCDefinitions::GetMinPatchNr( blocks[*pIter] ) < patch ) {
469 patchIndices.insert( pIter, ndx );
473 // pass event to CA Tracker
476 Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Reading hits",
477 "Total %d hits to read for slice %d", nClustersTotal, slice );
480 AliHLTTPCCAClusterData clusterData;
481 clusterData.StartReading( slice, nClustersTotal );
483 for ( std::vector<unsigned long>::iterator pIter = patchIndices.begin(); pIter != patchIndices.end(); pIter++ ) {
487 int patch = AliHLTTPCDefinitions::GetMinPatchNr( *iter );
488 inPtrSP = ( AliHLTTPCClusterData* )( iter->fPtr );
490 Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Reading hits",
491 "Reading %d hits for slice %d - patch %d", inPtrSP->fSpacePointCnt, slice, patch );
493 for ( unsigned int i = 0; i < inPtrSP->fSpacePointCnt; i++ ) {
494 AliHLTTPCSpacePointData *c = &( inPtrSP->fSpacePoints[i] );
495 if ( CAMath::Abs( c->fZ ) > fClusterZCut ) continue;
496 clusterData.ReadCluster( c->fID, c->fPadRow, c->fX, c->fY, c->fZ, c->fCharge );
500 clusterData.FinishReading();
502 // reconstruct the event
504 TStopwatch timerReco;
506 fTracker->ReadEvent( &clusterData );
508 fTracker->Reconstruct();
514 Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Reconstruct",
515 "%d tracks found for slice %d", fTracker->NOutTracks(), slice );
518 // write reconstructed tracks
520 unsigned int mySize = 0;
521 int ntracks = *fTracker->NOutTracks();
524 if ( !fNewOutputType ) {
526 AliHLTTPCTrackletData* outPtr = ( AliHLTTPCTrackletData* )( outputPtr );
528 AliHLTTPCTrackSegmentData* currOutTracklet = outPtr->fTracklets;
530 mySize = ( ( AliHLTUInt8_t * )currOutTracklet ) - ( ( AliHLTUInt8_t * )outputPtr );
532 outPtr->fTrackletCnt = 0;
534 for ( int itr = 0; itr < ntracks; itr++ ) {
536 AliHLTTPCCAOutTrack &t = fTracker->OutTracks()[itr];
538 //Logging( kHLTLogDebug, "HLT::TPCCATracker::DoEvent", "Wrtite output","track %d with %d hits", itr, t.NHits());
540 if ( t.NHits() < fMinNTrackClusters ) continue;
542 // calculate output track size
544 unsigned int dSize = sizeof( AliHLTTPCTrackSegmentData ) + t.NHits() * sizeof( unsigned int );
546 if ( mySize + dSize > maxBufferSize ) {
547 HLTWarning( "Output buffer size exceed (buffer size %d, current size %d), %d tracks are not stored", maxBufferSize, mySize, ntracks - itr + 1 );
552 // convert CA track parameters to HLT Track Segment
554 int iFirstRow = 1000;
556 int iFirstHit = fTracker->OutTrackHits()[t.FirstHitRef()];
557 int iLastHit = iFirstHit;
558 for ( int ih = 0; ih < t.NHits(); ih++ ) {
559 int hitID = fTracker->OutTrackHits()[t.FirstHitRef() + ih ];
560 int iRow = clusterData.RowNumber( hitID );
561 if ( iRow < iFirstRow ) { iFirstRow = iRow; iFirstHit = hitID; }
562 if ( iRow > iLastRow ) { iLastRow = iRow; iLastHit = hitID; }
565 AliHLTTPCCATrackParam par = t.StartPoint();
567 par.TransportToX( clusterData.X( iFirstHit ), .99 );
569 AliExternalTrackParam tp;
570 AliHLTTPCCATrackConvertor::GetExtParam( par, tp, 0 );
572 currOutTracklet->fX = tp.GetX();
573 currOutTracklet->fY = tp.GetY();
574 currOutTracklet->fZ = tp.GetZ();
575 currOutTracklet->fCharge = ( int ) tp.GetSign();
576 currOutTracklet->fPt = TMath::Abs( tp.GetSignedPt() );
577 float snp = tp.GetSnp() ;
578 if ( snp > .999 ) snp = .999;
579 if ( snp < -.999 ) snp = -.999;
580 currOutTracklet->fPsi = TMath::ASin( snp );
581 currOutTracklet->fTgl = tp.GetTgl();
583 currOutTracklet->fY0err = tp.GetSigmaY2();
584 currOutTracklet->fZ0err = tp.GetSigmaZ2();
585 float h = -currOutTracklet->fPt * currOutTracklet->fPt;
586 currOutTracklet->fPterr = h * h * tp.GetSigma1Pt2();
587 h = 1. / TMath::Sqrt( 1 - snp * snp );
588 currOutTracklet->fPsierr = h * h * tp.GetSigmaSnp2();
589 currOutTracklet->fTglerr = tp.GetSigmaTgl2();
591 if ( par.TransportToX( clusterData.X( iLastHit ), .99 ) ) {
592 currOutTracklet->fLastX = par.GetX();
593 currOutTracklet->fLastY = par.GetY();
594 currOutTracklet->fLastZ = par.GetZ();
596 currOutTracklet->fLastX = clusterData.X( iLastHit );
597 currOutTracklet->fLastY = clusterData.Y( iLastHit );
598 currOutTracklet->fLastZ = clusterData.Z( iLastHit );
600 //if( currOutTracklet->fLastX<10. ) {
601 //HLTError("CA last point: hitxyz=%f,%f,%f, track=%f,%f,%f, tracklet=%f,%f,%f, nhits=%d",clusterData.X( iLastHit ),clusterData.Y( iLastHit],clusterData.Z( iLastHit],
602 //par.GetX(), par.GetY(),par.GetZ(),currOutTracklet->fLastX,currOutTracklet->fLastY ,currOutTracklet->fLastZ, t.NHits());
604 #ifdef INCLUDE_TPC_HOUGH
605 #ifdef ROWHOUGHPARAMS
606 currOutTracklet->fTrackID = 0;
607 currOutTracklet->fRowRange1 = clusterData.RowNumber( iFirstHit );
608 currOutTracklet->fRowRange2 = clusterData.RowNumber( iLastHit );
609 currOutTracklet->fSector = slice;
610 currOutTracklet->fPID = 211;
612 #endif // INCLUDE_TPC_HOUGH
615 currOutTracklet->fNPoints = t.NHits();
617 for ( int i = 0; i < t.NHits(); i++ ) {
618 currOutTracklet->fPointIDs[i] = clusterData.Id( fTracker->OutTrackHits()[t.FirstHitRef()+i] );
621 currOutTracklet = ( AliHLTTPCTrackSegmentData* )( ( Byte_t * )currOutTracklet + dSize );
623 outPtr->fTrackletCnt++;
625 } else { // new output type
627 mySize = fTracker->Output()->EstimateSize( fTracker->Output()->NTracks(),
628 fTracker->Output()->NTrackClusters() );
629 if ( mySize <= maxBufferSize ) {
630 const AliHLTUInt8_t* outputevent = reinterpret_cast<const AliHLTUInt8_t*>( fTracker->Output() );
631 for ( unsigned int i = 0; i < mySize; i++ ) outputPtr[i] = outputevent[i];
633 HLTWarning( "Output buffer size exceed (buffer size %d, current size %d), tracks are not stored", maxBufferSize, mySize );
640 AliHLTComponentBlockData bd;
644 bd.fSpecification = AliHLTTPCDefinitions::EncodeDataSpecification( slice, slice, minPatch, maxPatch );
645 outputBlocks.push_back( bd );
651 fFullTime += timer.RealTime();
652 fRecoTime += timerReco.RealTime();
655 // Set log level to "Warning" for on-line system monitoring
656 int hz = ( int ) ( fFullTime > 1.e-10 ? fNEvents / fFullTime : 100000 );
657 int hz1 = ( int ) ( fRecoTime > 1.e-10 ? fNEvents / fRecoTime : 100000 );
658 HLTInfo( "CATracker slice %d: output %d tracks; input %d clusters, patches %d..%d, rows %d..%d; reco time %d/%d Hz",
659 slice, ntracks, nClustersTotal, minPatch, maxPatch, row[0], row[1], hz, hz1 );