track alpha angle normalized to [-Pi,+Pi]
[u/mrichter/AliRoot.git] / HLT / TPCLib / tracking-ca / AliHLTTPCCAGlobalMergerComponent.cxx
1 // **************************************************************************
2 // This file is property of and copyright by the ALICE HLT Project          *
3 // ALICE Experiment at CERN, All rights reserved.                           *
4 //                                                                          *
5 // Primary Authors: Sergey Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de> *
6 //                  Ivan Kisel <kisel@kip.uni-heidelberg.de>                *
7 //                  Matthias Kretz <kretz@kde.org>                          *
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 /** @file   AliHLTTPCCAGlobalMergerComponent.cxx
22     @author Matthias Kretz
23     @date
24     @brief  HLT TPC CA global merger component.
25 */
26
27 #if __GNUC__>= 3
28 using namespace std;
29 #endif
30
31 #include "AliHLTTPCCAGlobalMergerComponent.h"
32 #include "AliHLTTPCCAMergerOutput.h"
33 #include "AliHLTTPCTransform.h"
34 #include "AliHLTTPCCAMerger.h"
35 #include "AliHLTTPCVertex.h"
36 #include "AliHLTTPCVertexData.h"
37 #include "AliHLTTPCTrackSegmentData.h"
38 #include "AliHLTTPCTrack.h"
39 #include "AliHLTTPCTrackArray.h"
40 #include "AliHLTTPCTrackletDataFormat.h"
41 #include "AliHLTTPCCADef.h"
42 #include "AliHLTTPCDefinitions.h"
43 #include "AliHLTTPCCATrackConvertor.h"
44 #include "AliHLTTPCCASliceOutput.h"
45
46 #include "AliCDBEntry.h"
47 #include "AliCDBManager.h"
48 #include "TObjString.h"
49 #include "TObjArray.h"
50 #include "AliHLTExternalTrackParam.h"
51
52 #include <climits>
53 #include <cstdlib>
54 #include <cerrno>
55
56
57 // ROOT macro for the implementation of ROOT specific class methods
58 ClassImp( AliHLTTPCCAGlobalMergerComponent )
59
60
61 AliHLTTPCCAGlobalMergerComponent::AliHLTTPCCAGlobalMergerComponent()
62   : fGlobalMerger( 0 ), fSolenoidBz( 0 ), fClusterErrorCorrectionY(0), fClusterErrorCorrectionZ(0)
63 {
64   // see header file for class documentation
65 }
66
67 // Public functions to implement AliHLTComponent's interface.
68 // These functions are required for the registration process
69
70 const char *AliHLTTPCCAGlobalMergerComponent::GetComponentID()
71 {
72   // see header file for class documentation
73   return "TPCCAGlobalMerger";
74 }
75
76 void AliHLTTPCCAGlobalMergerComponent::GetInputDataTypes( AliHLTComponentDataTypeList &list )
77 {
78   // see header file for class documentation
79   list.clear();
80   list.push_back( AliHLTTPCCADefinitions::fgkTrackletsDataType );
81   //list.push_back( AliHLTTPCDefinitions::fgkTrackSegmentsDataType );
82   //list.push_back( AliHLTTPCDefinitions::fgkVertexDataType );
83 }
84
85 AliHLTComponentDataType AliHLTTPCCAGlobalMergerComponent::GetOutputDataType()
86 {
87   // see header file for class documentation
88   //return AliHLTTPCDefinitions::fgkTracksDataType; // old
89   return kAliHLTDataTypeTrack|kAliHLTDataOriginTPC;
90 }
91
92 void AliHLTTPCCAGlobalMergerComponent::GetOutputDataSize( unsigned long &constBase, double &inputMultiplier )
93 {
94   // see header file for class documentation
95   // XXX TODO: Find more realistic values.
96   constBase = 0;
97   inputMultiplier = 1.0;
98 }
99
100 AliHLTComponent *AliHLTTPCCAGlobalMergerComponent::Spawn()
101 {
102   // see header file for class documentation
103   return new AliHLTTPCCAGlobalMergerComponent;
104 }
105
106
107
108
109 void AliHLTTPCCAGlobalMergerComponent::SetDefaultConfiguration()
110 {
111   // Set default configuration for the CA merger component
112   // Some parameters can be later overwritten from the OCDB
113
114   fSolenoidBz = -5.00668;
115   fClusterErrorCorrectionY = 0;
116   fClusterErrorCorrectionZ = 1.1;
117 }
118
119 int AliHLTTPCCAGlobalMergerComponent::ReadConfigurationString(  const char* arguments )
120 {
121   // Set configuration parameters for the CA merger component from the string
122
123   int iResult = 0;
124   if ( !arguments ) return iResult;
125
126   TString allArgs = arguments;
127   TString argument;
128   int bMissingParam = 0;
129
130   TObjArray* pTokens = allArgs.Tokenize( " " );
131
132   int nArgs =  pTokens ? pTokens->GetEntries() : 0;
133
134   for ( int i = 0; i < nArgs; i++ ) {
135     argument = ( ( TObjString* )pTokens->At( i ) )->GetString();
136     if ( argument.IsNull() ) continue;
137
138     if ( argument.CompareTo( "-solenoidBz" ) == 0 ) {
139       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
140       fSolenoidBz = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
141       HLTInfo( "Magnetic Field set to: %f", fSolenoidBz );
142       continue;
143     }
144
145     if ( argument.CompareTo( "-errorCorrectionY" ) == 0 ) {
146       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
147       fClusterErrorCorrectionY = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
148       HLTInfo( "Cluster Y error correction factor set to: %f", fClusterErrorCorrectionY );
149       continue;
150     }
151
152    if ( argument.CompareTo( "-errorCorrectionZ" ) == 0 ) {
153       if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
154       fClusterErrorCorrectionZ = ( ( TObjString* )pTokens->At( i ) )->GetString().Atof();
155       HLTInfo( "Cluster Z error correction factor set to: %f", fClusterErrorCorrectionZ );
156       continue;
157     }
158
159     HLTError( "Unknown option \"%s\"", argument.Data() );
160     iResult = -EINVAL;
161   }
162   delete pTokens;
163
164   if ( bMissingParam ) {
165     HLTError( "Specifier missed for parameter \"%s\"", argument.Data() );
166     iResult = -EINVAL;
167   }
168
169   return iResult;
170 }
171
172
173 int AliHLTTPCCAGlobalMergerComponent::ReadCDBEntry( const char* cdbEntry, const char* chainId )
174 {
175   // see header file for class documentation
176
177   const char* defaultNotify = "";
178
179   if ( !cdbEntry ) {
180     cdbEntry = "HLT/ConfigTPC/TPCCAGlobalMerger";
181     defaultNotify = " (default)";
182     chainId = 0;
183   }
184
185   HLTInfo( "configure from entry \"%s\"%s, chain id %s", cdbEntry, defaultNotify, ( chainId != NULL && chainId[0] != 0 ) ? chainId : "<none>" );
186   AliCDBEntry *pEntry = AliCDBManager::Instance()->Get( cdbEntry );//,GetRunNo());
187
188   if ( !pEntry ) {
189     HLTError( "cannot fetch object \"%s\" from CDB", cdbEntry );
190     return -EINVAL;
191   }
192
193   TObjString* pString = dynamic_cast<TObjString*>( pEntry->GetObject() );
194
195   if ( !pString ) {
196     HLTError( "configuration object \"%s\" has wrong type, required TObjString", cdbEntry );
197     return -EINVAL;
198   }
199
200   HLTInfo( "received configuration object string: \"%s\"", pString->GetString().Data() );
201
202   return  ReadConfigurationString( pString->GetString().Data() );
203 }
204
205
206
207 int AliHLTTPCCAGlobalMergerComponent::Configure( const char* cdbEntry, const char* chainId, const char *commandLine )
208 {
209   // Configure the component
210   // There are few levels of configuration,
211   // parameters which are set on one step can be overwritten on the next step
212
213   //* read hard-coded values
214
215   SetDefaultConfiguration();
216
217   //* read the default CDB entry
218
219   int iResult1 = ReadCDBEntry( NULL, chainId );
220
221   //* read magnetic field
222
223   int iResult2 = ReadCDBEntry( kAliHLTCDBSolenoidBz, chainId );
224
225   //* read the actual CDB entry if required
226
227   int iResult3 = ( cdbEntry ) ? ReadCDBEntry( cdbEntry, chainId ) : 0;
228
229   //* read extra parameters from input (if they are)
230
231   int iResult4 = 0;
232
233   if ( commandLine && commandLine[0] != '\0' ) {
234     HLTInfo( "received configuration string from HLT framework: \"%s\"", commandLine );
235     iResult4 = ReadConfigurationString( commandLine );
236   }
237
238
239   // Initialize the merger
240
241   AliHLTTPCCAParam param;
242
243   {
244     // get gemetry
245     int iSec = 0;
246     float inRmin = 83.65;
247     float outRmax = 247.7;
248     float plusZmin = 0.0529937;
249     float plusZmax = 249.778;
250     float minusZmin = -249.645;
251     float minusZmax = -0.0799937;
252     float dalpha = 0.349066;
253     float alpha = 0.174533 + dalpha * iSec;
254     bool zPlus = ( iSec < 18 );
255     float zMin =  zPlus ? plusZmin : minusZmin;
256     float zMax =  zPlus ? plusZmax : minusZmax;
257     int nRows = AliHLTTPCTransform::GetNRows();
258     float padPitch = 0.4;
259     float sigmaZ = 0.228808;
260     float *rowX = new float [nRows];
261     for ( int irow = 0; irow < nRows; irow++ ) {
262       rowX[irow] = AliHLTTPCTransform::Row2X( irow );
263     }
264
265     param.Initialize( iSec, nRows, rowX, alpha, dalpha,
266                       inRmin, outRmax, zMin, zMax, padPitch, sigmaZ, fSolenoidBz );
267
268     if( fClusterErrorCorrectionY>1.e-4 ) param.SetClusterError2CorrectionY( fClusterErrorCorrectionY*fClusterErrorCorrectionY );
269     if( fClusterErrorCorrectionZ>1.e-4 ) param.SetClusterError2CorrectionZ( fClusterErrorCorrectionZ*fClusterErrorCorrectionZ );
270     param.Update();
271
272     delete[] rowX;
273   }
274
275
276   fGlobalMerger->SetSliceParam( param );
277
278   return iResult1 ? iResult1 : ( iResult2 ? iResult2 : ( iResult3 ? iResult3 : iResult4 ) );
279 }
280
281
282
283
284 int AliHLTTPCCAGlobalMergerComponent::DoInit( int argc, const char** argv )
285 {
286   // see header file for class documentation
287
288   if ( fGlobalMerger ) {
289     return EINPROGRESS;
290   }
291
292   fGlobalMerger = new AliHLTTPCCAMerger();
293
294   TString arguments = "";
295   for ( int i = 0; i < argc; i++ ) {
296     if ( !arguments.IsNull() ) arguments += " ";
297     arguments += argv[i];
298   }
299
300   return Configure( NULL, NULL, arguments.Data()  );
301 }
302
303 int AliHLTTPCCAGlobalMergerComponent::Reconfigure( const char* cdbEntry, const char* chainId )
304 {
305   // Reconfigure the component from OCDB
306
307   return Configure( cdbEntry, chainId, NULL );
308 }
309
310
311
312 int AliHLTTPCCAGlobalMergerComponent::DoDeinit()
313 {
314   // see header file for class documentation
315   delete fGlobalMerger;
316   fGlobalMerger = 0;
317
318   return 0;
319 }
320
321 int AliHLTTPCCAGlobalMergerComponent::DoEvent( const AliHLTComponentEventData &evtData,
322     const AliHLTComponentBlockData *blocks, AliHLTComponentTriggerData &/*trigData*/,
323     AliHLTUInt8_t *outputPtr, AliHLTUInt32_t &size, AliHLTComponentBlockDataList &outputBlocks )
324 {
325   // see header file for class documentation
326   int iResult = 0;
327   unsigned int maxBufferSize = size;
328
329   size = 0;
330
331   if ( !outputPtr ) {
332     return -ENOSPC;
333   }
334   if ( !IsDataEvent() ) {
335     return 0;
336   }
337
338   fGlobalMerger->Clear();
339
340   const AliHLTComponentBlockData *const blocksEnd = blocks + evtData.fBlockCnt;
341   for ( const AliHLTComponentBlockData *block = blocks; block < blocksEnd; ++block ) {
342     if ( block->fDataType != AliHLTTPCCADefinitions::fgkTrackletsDataType ) {
343       continue;
344     }
345
346     int slice = AliHLTTPCDefinitions::GetMinSliceNr( *block );
347     if ( slice < 0 || slice >= AliHLTTPCTransform::GetNSlice() ) {
348       HLTError( "invalid slice number %d extracted from specification 0x%08lx,  skipping block of type %s",
349                 slice, block->fSpecification, DataType2Text( block->fDataType ).c_str() );
350       // just remember the error, if there are other valid blocks ignore the error, return code otherwise
351       iResult = -EBADF;
352       continue;
353     }
354
355     if ( slice != AliHLTTPCDefinitions::GetMaxSliceNr( *block ) ) {
356       // the code was not written for/ never used with multiple slices in one data block/ specification
357       HLTWarning( "specification 0x%08lx indicates multiple slices in data block %s: never used before, please audit the code",
358                   block->fSpecification, DataType2Text( block->fDataType ).c_str() );
359     }
360     AliHLTTPCCASliceOutput *sliceOut =  reinterpret_cast<AliHLTTPCCASliceOutput *>( block->fPtr );
361     sliceOut->SetPointers();
362     fGlobalMerger->SetSliceData( slice, sliceOut );
363
364         /*char filename[256];
365         sprintf(filename, "debug%d.out", slice);
366         FILE* fp = fopen(filename, "w+b");
367         if (fp == NULL) printf("Error!!!\n");
368         fwrite(sliceOut, 1, sliceOut->EstimateSize(sliceOut->NTracks(), sliceOut->NTrackClusters()), fp);
369         fclose(fp);*/
370   }
371   fGlobalMerger->Reconstruct();
372
373   const AliHLTTPCCAMergerOutput *mergerOutput = fGlobalMerger->Output();
374
375
376   // Fill output tracks
377
378   unsigned int mySize = 0;
379   {
380     AliHLTTracksData* outPtr = ( AliHLTTracksData* )( outputPtr );
381
382     AliHLTExternalTrackParam* currOutTrack = outPtr->fTracklets;
383
384     mySize =   ( ( AliHLTUInt8_t * )currOutTrack ) -  ( ( AliHLTUInt8_t * )outputPtr );
385
386     outPtr->fCount = 0;
387    
388     int nTracks = mergerOutput->NTracks();
389
390     for ( int itr = 0; itr < nTracks; itr++ ) {
391
392       // convert AliHLTTPCCAMergedTrack to AliHLTTrack
393
394       const AliHLTTPCCAMergedTrack &track = mergerOutput->Track( itr );
395
396       unsigned int dSize = sizeof( AliHLTExternalTrackParam ) + track.NClusters() * sizeof( unsigned int );
397
398       if ( mySize + dSize > maxBufferSize ) {
399         HLTWarning( "Output buffer size exceed (buffer size %d, current size %d), %d tracks are not stored", maxBufferSize, mySize, nTracks - itr + 1 );
400         iResult = -ENOSPC;
401         break;
402       }
403
404       // first convert to AliExternalTrackParam
405
406       AliExternalTrackParam tp, tpEnd;
407       AliHLTTPCCATrackConvertor::GetExtParam( track.InnerParam(), tp,  track.InnerAlpha() );
408       AliHLTTPCCATrackConvertor::GetExtParam( track.OuterParam(), tpEnd, track.OuterAlpha() );
409       
410       // normalize the angle to +-Pi
411               
412       currOutTrack->fAlpha = tp.GetAlpha() - CAMath::Nint(tp.GetAlpha()/CAMath::TwoPi())*CAMath::TwoPi();      
413       currOutTrack->fX = tp.GetX();
414       currOutTrack->fY = tp.GetY();
415       currOutTrack->fZ = tp.GetZ();      
416       {
417         float sinA = TMath::Sin( track.OuterAlpha() - track.InnerAlpha());
418         float cosA = TMath::Cos( track.OuterAlpha() - track.InnerAlpha());
419         currOutTrack->fLastX = tpEnd.GetX()*cosA - tpEnd.GetY()*sinA;
420         currOutTrack->fLastY = tpEnd.GetX()*sinA + tpEnd.GetY()*cosA;
421         currOutTrack->fLastZ = tpEnd.GetZ();
422       }
423       currOutTrack->fq1Pt = tp.GetSigned1Pt();
424       currOutTrack->fSinPsi = tp.GetSnp();
425       currOutTrack->fTgl = tp.GetTgl();
426       for( int i=0; i<15; i++ ) currOutTrack->fC[i] = tp.GetCovariance()[i];
427       currOutTrack->fTrackID = itr;
428       currOutTrack->fFlags = 0;
429       currOutTrack->fNPoints = track.NClusters();    
430       for ( int i = 0; i < track.NClusters(); i++ ) currOutTrack->fPointIDs[i] = mergerOutput->ClusterId( track.FirstClusterRef() + i );
431
432       currOutTrack = ( AliHLTExternalTrackParam* )( (( Byte_t * )currOutTrack) + dSize );
433       mySize += dSize;
434       outPtr->fCount++;
435     }
436   
437
438     AliHLTComponentBlockData resultData;
439     FillBlockData( resultData );
440     resultData.fOffset = 0;
441     resultData.fSize = mySize;
442     resultData.fDataType = kAliHLTDataTypeTrack|kAliHLTDataOriginTPC;
443     resultData.fSpecification = AliHLTTPCDefinitions::EncodeDataSpecification( 0, 35, 0, 5 );
444     outputBlocks.push_back( resultData );
445     size = resultData.fSize;
446   }
447
448   HLTInfo( "CAGlobalMerger:: output %d tracks", mergerOutput->NTracks() );
449
450
451   /* old format
452   {
453     // check if there was enough space in the output buffer
454
455     int nTracks = mergerOutput->NTracks();
456
457     AliHLTTPCTrackArray array( nTracks );
458
459     int nClusters = 0;
460     for ( int itr = 0; itr < nTracks; itr++ ) {
461
462       // convert AliHLTTPCCAMergedTrack to AliHLTTPCTrack
463
464       const AliHLTTPCCAMergedTrack &track = mergerOutput->Track( itr );
465
466       AliHLTTPCTrack out;
467
468       // first convert to AliExternalTrackParam
469
470       AliExternalTrackParam tp, tpEnd;
471       AliHLTTPCCATrackConvertor::GetExtParam( track.InnerParam(), tp, 0 );
472       AliHLTTPCCATrackConvertor::GetExtParam( track.OuterParam(), tpEnd, 0 );
473
474       // set parameters, with rotation to global coordinates
475
476       out.SetCharge( ( int ) tp.GetSign() );
477       out.SetPt( TMath::Abs( tp.GetSignedPt() ) );
478       out.SetPsi( fmod( TMath::ASin( tp.GetSnp() ) + track.InnerAlpha() , 2*TMath::Pi() ) );
479       out.SetTgl( tp.GetTgl() );
480       {
481         float sinA = TMath::Sin( track.InnerAlpha() );
482         float cosA = TMath::Cos( track.InnerAlpha() );
483
484         out.SetFirstPoint( tp.GetX()*cosA - tp.GetY()*sinA,
485                            tp.GetX()*sinA + tp.GetY()*cosA,
486                            tp.GetZ() );
487       }
488
489       {
490         float sinA = TMath::Sin( track.OuterAlpha() );
491         float cosA = TMath::Cos( track.OuterAlpha() );
492
493         out.SetLastPoint( tpEnd.GetX()*cosA - tpEnd.GetY()*sinA,
494                           tpEnd.GetX()*sinA + tpEnd.GetY()*cosA,
495                           tpEnd.GetZ() );
496       }
497
498       // set parameter errors w/o rotation, as it is done in AliHLTTPCTrackArray
499
500       out.SetY0err( tp.GetSigmaY2() );
501       out.SetZ0err( tp.GetSigmaZ2() );
502       out.SetPterr( tp.GetSigma1Pt2() );
503       out.SetPsierr( tp.GetSigmaSnp2() );
504       out.SetTglerr( tp.GetSigmaTgl2() );
505
506       // set cluster ID's
507
508       unsigned int hitID[1000];
509       for ( int i = 0; i < track.NClusters(); i++ ) hitID[i] = mergerOutput->ClusterId( track.FirstClusterRef() + i );
510
511       out.SetNHits( track.NClusters() );
512       out.SetHits( track.NClusters(), hitID );
513
514       out.SetSector( -1 );
515       out.CalculateHelix();
516       if ( !out.CheckConsistency() )  *( array.NextTrack() ) = out;
517       nClusters += track.NClusters();
518     }
519
520
521     if ( sizeof( AliHLTTPCTrackletData ) + nTracks*sizeof( AliHLTTPCTrackSegmentData ) + nClusters*sizeof( unsigned int )
522          > maxBufferSize ) {
523       iResult = -ENOSPC;
524     } else {
525       AliHLTTPCTrackletData *outPtr = ( AliHLTTPCTrackletData* )( outputPtr );
526       unsigned int nOutTracks = 0;
527       mySize = array.WriteTracks( nOutTracks, outPtr->fTracklets );
528       mySize += sizeof( AliHLTTPCTrackletData );
529       outPtr->fTrackletCnt = nOutTracks;
530     }
531   }
532
533   AliHLTComponentBlockData resultData;
534   FillBlockData( resultData );
535   resultData.fOffset = 0;
536   resultData.fSize = mySize;
537   resultData.fSpecification = AliHLTTPCDefinitions::EncodeDataSpecification( 0, 35, 0, 5 );
538   outputBlocks.push_back( resultData );
539   size = resultData.fSize;
540
541   HLTInfo( "CAGlobalMerger:: output %d tracks", mergerOutput->NTracks() );
542   */
543   return iResult;
544 }
545