Line data Source code
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 : // for The ALICE HLT Project. *
7 : // *
8 : // Permission to use, copy, modify and distribute this software and its *
9 : // documentation strictly for non-commercial purposes is hereby granted *
10 : // without fee, provided that the above copyright notice appears in all *
11 : // copies and that both the copyright notice and this permission notice *
12 : // appear in the supporting documentation. The authors make no claims *
13 : // about the suitability of this software for any purpose. It is *
14 : // provided "as is" without express or implied warranty. *
15 : // *
16 : //***************************************************************************
17 :
18 : /// @file AliHLTTPCdEdxComponent.cxx
19 : /// @author Sergey Gorbunov <sergey.gorbunov@kip.uni-heidelberg.de>
20 : /// @date June 2009
21 : /// @brief An ITS tracker processing component for the HLT
22 :
23 :
24 : /////////////////////////////////////////////////////
25 : // //
26 : // dEdx calculation component for the HLT TPC //
27 : // //
28 : /////////////////////////////////////////////////////
29 :
30 : #include "AliHLTTPCdEdxComponent.h"
31 : #include "TStopwatch.h"
32 : #include "TMath.h"
33 : #include "AliCDBEntry.h"
34 : #include "AliCDBManager.h"
35 : #include "TObjString.h"
36 : #include "TObjArray.h"
37 : #include "AliHLTDataTypes.h"
38 : #include "AliHLTExternalTrackParam.h"
39 : #include "AliHLTGlobalBarrelTrack.h"
40 : #include "AliHLTTPCGeometry.h"
41 : #include "AliHLTTPCSpacePointData.h"
42 : #include "AliHLTTPCClusterDataFormat.h"
43 : #include "AliHLTTPCDefinitions.h"
44 : #include "AliTPCseed.h"
45 : #include "AliTPCclusterMI.h"
46 : #include "TGeoGlobalMagField.h"
47 : #include "AliTPCcalibDB.h"
48 : #include "AliTPCRecoParam.h"
49 : #include "AliTPCTransform.h"
50 :
51 : /** ROOT macro for the implementation of ROOT specific class methods */
52 6 : ClassImp( AliHLTTPCdEdxComponent )
53 3 : AliHLTTPCdEdxComponent::AliHLTTPCdEdxComponent()
54 : :
55 3 : fSolenoidBz( 0 ),
56 3 : fStatTime( 0 ),
57 3 : fStatNEvents( 0 )
58 15 : {
59 : // see header file for class documentation
60 : // or
61 : // refer to README to build package
62 : // or
63 : // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
64 1302 : for( int i=0; i<fkNPatches; i++ ){
65 648 : fPatchClusters[i] = 0;
66 648 : fNPatchClusters[i] = 0;
67 : }
68 6 : }
69 :
70 : AliHLTTPCdEdxComponent::AliHLTTPCdEdxComponent( const AliHLTTPCdEdxComponent& )
71 : :
72 0 : AliHLTProcessor(),
73 0 : fSolenoidBz( 0 ),
74 0 : fStatTime( 0 ),
75 0 : fStatNEvents( 0 )
76 0 : {
77 : // dummy
78 0 : for( int i=0; i<fkNPatches; i++ ){
79 0 : fPatchClusters[i] = 0;
80 0 : fNPatchClusters[i] = 0;
81 : }
82 0 : HLTFatal( "copy constructor untested" );
83 0 : }
84 :
85 : AliHLTTPCdEdxComponent& AliHLTTPCdEdxComponent::operator=( const AliHLTTPCdEdxComponent& )
86 : {
87 : // see header file for class documentation
88 0 : HLTFatal( "assignment operator untested" );
89 0 : for( int i=0; i<fkNPatches; i++ ){
90 0 : fPatchClusters[i] = 0;
91 0 : fNPatchClusters[i] = 0;
92 : }
93 0 : return *this;
94 : }
95 :
96 : AliHLTTPCdEdxComponent::~AliHLTTPCdEdxComponent()
97 18 : {
98 : // see header file for class documentation
99 1302 : for( int i=0; i<fkNPatches; i++ ){
100 648 : delete[] fPatchClusters[i];
101 : }
102 9 : }
103 :
104 : //
105 : // Public functions to implement AliHLTComponent's interface.
106 : // These functions are required for the registration process
107 : //
108 :
109 : const char* AliHLTTPCdEdxComponent::GetComponentID()
110 : {
111 : // see header file for class documentation
112 528 : return "TPCdEdx";
113 : }
114 :
115 : void AliHLTTPCdEdxComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list )
116 : {
117 : // see header file for class documentation
118 0 : list.clear();
119 0 : list.push_back( kAliHLTDataTypeTrack|kAliHLTDataOriginTPC );
120 0 : list.push_back( AliHLTTPCDefinitions::fgkClustersDataType );
121 0 : }
122 :
123 : AliHLTComponentDataType AliHLTTPCdEdxComponent::GetOutputDataType()
124 : {
125 : // see header file for class documentation
126 0 : return kAliHLTDataTypedEdx|kAliHLTDataOriginTPC;
127 : }
128 :
129 : void AliHLTTPCdEdxComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier )
130 : {
131 : // define guess for the output data size
132 0 : constBase = 200; // minimum size
133 0 : inputMultiplier = 0.1; // size relative to input
134 0 : }
135 :
136 : AliHLTComponent* AliHLTTPCdEdxComponent::Spawn()
137 : {
138 : // see header file for class documentation
139 0 : return new AliHLTTPCdEdxComponent;
140 0 : }
141 :
142 : void AliHLTTPCdEdxComponent::SetDefaultConfiguration()
143 : {
144 : // Set default configuration for the CA tracker component
145 : // Some parameters can be later overwritten from the OCDB
146 :
147 0 : fSolenoidBz = -5.00668;
148 0 : fStatTime = 0;
149 0 : fStatNEvents = 0;
150 0 : }
151 :
152 : int AliHLTTPCdEdxComponent::ReadConfigurationString( const char* arguments )
153 : {
154 : // Set configuration parameters for the CA tracker component from the string
155 :
156 : int iResult = 0;
157 0 : if ( !arguments ) return iResult;
158 :
159 0 : TString allArgs = arguments;
160 0 : TString argument;
161 : int bMissingParam = 0;
162 :
163 0 : TObjArray* pTokens = allArgs.Tokenize( " " );
164 :
165 0 : int nArgs = pTokens ? pTokens->GetEntries() : 0;
166 :
167 0 : for ( int i = 0; i < nArgs; i++ ) {
168 0 : argument = ( ( TObjString* )pTokens->At( i ) )->GetString();
169 0 : if ( argument.IsNull() ) continue;
170 :
171 0 : if ( argument.CompareTo( "-solenoidBz" ) == 0 ) {
172 0 : if ( ( bMissingParam = ( ++i >= pTokens->GetEntries() ) ) ) break;
173 0 : HLTWarning("argument -solenoidBz is deprecated, magnetic field set up globally (%f)", GetBz());
174 : continue;
175 : }
176 :
177 0 : HLTError( "Unknown option \"%s\"", argument.Data() );
178 : iResult = -EINVAL;
179 0 : }
180 0 : delete pTokens;
181 :
182 0 : if ( bMissingParam ) {
183 0 : HLTError( "Specifier missed for parameter \"%s\"", argument.Data() );
184 : iResult = -EINVAL;
185 0 : }
186 :
187 : return iResult;
188 0 : }
189 :
190 :
191 : int AliHLTTPCdEdxComponent::ReadCDBEntry( const char* cdbEntry, const char* chainId )
192 : {
193 : // see header file for class documentation
194 :
195 : const char* defaultNotify = "";
196 :
197 0 : if ( !cdbEntry ) {
198 0 : return 0;// need to add the HLT/ConfigTPC/TPCdEdx directory to cdb SG!!!
199 : //cdbEntry = "HLT/ConfigTPC/TPCdEdx";
200 : //defaultNotify = " (default)";
201 : //chainId = 0;
202 : }
203 :
204 0 : HLTInfo( "configure from entry \"%s\"%s, chain id %s", cdbEntry, defaultNotify, ( chainId != NULL && chainId[0] != 0 ) ? chainId : "<none>" );
205 0 : AliCDBEntry *pEntry = AliCDBManager::Instance()->Get( cdbEntry );//,GetRunNo());
206 :
207 0 : if ( !pEntry ) {
208 0 : HLTError( "cannot fetch object \"%s\" from CDB", cdbEntry );
209 0 : return -EINVAL;
210 : }
211 :
212 0 : TObjString* pString = dynamic_cast<TObjString*>( pEntry->GetObject() );
213 :
214 0 : if ( !pString ) {
215 0 : HLTError( "configuration object \"%s\" has wrong type, required TObjString", cdbEntry );
216 0 : return -EINVAL;
217 : }
218 :
219 0 : HLTInfo( "received configuration object string: \"%s\"", pString->GetString().Data() );
220 :
221 0 : return ReadConfigurationString( pString->GetString().Data() );
222 0 : }
223 :
224 :
225 : int AliHLTTPCdEdxComponent::Configure( const char* cdbEntry, const char* chainId, const char *commandLine )
226 : {
227 : // Configure the component
228 : // There are few levels of configuration,
229 : // parameters which are set on one step can be overwritten on the next step
230 :
231 : //* read hard-coded values
232 :
233 0 : SetDefaultConfiguration();
234 :
235 : //* read the default CDB entry
236 :
237 0 : int iResult1 = ReadCDBEntry( NULL, chainId );
238 :
239 : //* read magnetic field
240 :
241 0 : fSolenoidBz=GetBz();
242 :
243 : //* read the actual CDB entry if required
244 :
245 0 : int iResult2 = ( cdbEntry ) ? ReadCDBEntry( cdbEntry, chainId ) : 0;
246 :
247 : //* read extra parameters from input (if they are)
248 :
249 : int iResult3 = 0;
250 :
251 0 : if ( commandLine && commandLine[0] != '\0' ) {
252 0 : HLTInfo( "received configuration string from HLT framework: \"%s\"", commandLine );
253 0 : iResult3 = ReadConfigurationString( commandLine );
254 0 : }
255 :
256 : // Initialise the tracker here
257 :
258 0 : return iResult1 ? iResult1 : ( iResult2 ? iResult2 : iResult3 );
259 : }
260 :
261 :
262 :
263 : int AliHLTTPCdEdxComponent::DoInit( int argc, const char** argv )
264 : {
265 : // Configure the component
266 :
267 0 : TString arguments = "";
268 0 : for ( int i = 0; i < argc; i++ ) {
269 0 : if ( !arguments.IsNull() ) arguments += " ";
270 0 : arguments += argv[i];
271 : }
272 :
273 0 : int ret = Configure( NULL, NULL, arguments.Data() );
274 :
275 : // Check field
276 0 : if (!TGeoGlobalMagField::Instance()) {
277 0 : HLTError("magnetic field not initialized, please set up TGeoGlobalMagField and AliMagF");
278 0 : return -ENODEV;
279 : }
280 :
281 0 : AliTPCcalibDB::Instance()->SetRun(GetRunNo());
282 0 : AliTPCcalibDB::Instance()->UpdateRunInformations(GetRunNo());
283 :
284 0 : AliTPCcalibDB::Instance()->GetTransform()->SetCurrentRun(GetRunNo());
285 0 : AliTPCcalibDB::Instance()->GetTransform()->SetCurrentTimeStamp( GetTimeStamp() );
286 0 : AliTPCcalibDB::Instance()->GetTransform()->SetCurrentRecoParam(AliTPCRecoParam::GetHLTParam());
287 :
288 : //AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform();
289 : //if( transform ){
290 : //AliTPCRecoParam *reco = transform->GetCurrentRecoParam();
291 : //if( ! reco ){
292 : //reco = new AliTPCRecoParam;
293 : //reco->SetUseTotalCharge(0);
294 : //transform->SetCurrentRecoParam( reco );
295 : //}
296 :
297 0 : return ret;
298 0 : }
299 :
300 :
301 : int AliHLTTPCdEdxComponent::DoDeinit()
302 : {
303 : // see header file for class documentation
304 0 : return 0;
305 : }
306 :
307 :
308 :
309 : int AliHLTTPCdEdxComponent::Reconfigure( const char* cdbEntry, const char* chainId )
310 : {
311 : // Reconfigure the component from OCDB .
312 :
313 0 : return Configure( cdbEntry, chainId, NULL );
314 : }
315 :
316 :
317 :
318 : int AliHLTTPCdEdxComponent::DoEvent
319 : (
320 : const AliHLTComponentEventData& evtData,
321 : const AliHLTComponentBlockData* blocks,
322 : AliHLTComponentTriggerData& /*trigData*/,
323 : AliHLTUInt8_t* outputPtr,
324 : AliHLTUInt32_t& size,
325 : vector<AliHLTComponentBlockData>& outputBlocks )
326 : {
327 : //* process the event
328 :
329 0 : AliHLTUInt32_t maxBufferSize = size;
330 0 : size = 0; // output size
331 :
332 0 : if (!IsDataEvent()) return 0;
333 :
334 0 : if ( evtData.fBlockCnt <= 0 ) {
335 0 : HLTWarning( "no blocks in event" );
336 0 : return 0;
337 : }
338 0 : if ( !outputPtr ) {
339 0 : return -ENOSPC;
340 : }
341 :
342 : int iResult=0;
343 :
344 0 : TStopwatch timer;
345 :
346 : // Initialise the data pointers
347 :
348 0 : for( int i=0; i<fkNPatches; i++ ){
349 0 : delete[] fPatchClusters[i];
350 0 : fPatchClusters[i] = 0;
351 0 : fNPatchClusters[i] = 0;
352 : }
353 :
354 0 : int nBlocks = (int)evtData.fBlockCnt;
355 :
356 : int nInputClusters = 0;
357 : int nInputTracks = 0;
358 :
359 : // first read all the clusters
360 :
361 0 : for (int ndx=0; ndx<nBlocks && iResult>=0; ndx++) {
362 0 : const AliHLTComponentBlockData* iter = blocks+ndx;
363 0 : if ( iter->fDataType != AliHLTTPCDefinitions::fgkClustersDataType ) continue;
364 0 : Int_t slice=AliHLTTPCDefinitions::GetMinSliceNr(iter->fSpecification);
365 0 : Int_t patch=AliHLTTPCDefinitions::GetMinPatchNr(iter->fSpecification);
366 0 : Int_t slicepatch=slice*6+patch;
367 0 : if( slicepatch >= fkNPatches ){
368 0 : HLTWarning("Wrong header of TPC cluster data, slice %d, patch %d",
369 : slice, patch );
370 0 : continue;
371 : }
372 0 : AliHLTTPCClusterData* inPtrSP = ( AliHLTTPCClusterData* )( iter->fPtr );
373 0 : nInputClusters += inPtrSP->fSpacePointCnt;
374 :
375 0 : delete[] fPatchClusters[slicepatch];
376 0 : fPatchClusters[slicepatch] = new AliTPCclusterMI[inPtrSP->fSpacePointCnt];
377 0 : fNPatchClusters[slicepatch] = inPtrSP->fSpacePointCnt;
378 :
379 : // create off-line clusters out of the HLT clusters
380 : // todo: check which cluster information is really needed for the dEdx
381 :
382 0 : for ( unsigned int i = 0; i < inPtrSP->fSpacePointCnt; i++ ) {
383 0 : AliHLTTPCSpacePointData *chlt = &( inPtrSP->fSpacePoints[i] );
384 0 : AliTPCclusterMI *c = fPatchClusters[slicepatch]+i;
385 0 : c->SetX(chlt->fX);
386 0 : c->SetY(chlt->fY);
387 0 : c->SetZ(chlt->fZ);
388 0 : c->SetSigmaY2(chlt->fSigmaY2);
389 0 : c->SetSigmaYZ( 0 );
390 0 : c->SetSigmaZ2(chlt->fSigmaZ2);
391 0 : c->SetQ( chlt->fCharge );
392 0 : c->SetMax( chlt->fQMax );
393 0 : Int_t sector, row;
394 0 : Float_t padtime[3]={0,chlt->fY,chlt->fZ};
395 0 : AliHLTTPCGeometry::Slice2Sector(slice,chlt->fPadRow, sector, row);
396 0 : AliHLTTPCGeometry::Local2Raw( padtime, sector, row);
397 0 : c->SetDetector( sector );
398 0 : c->SetRow( row );
399 0 : c->SetPad( (Int_t) padtime[1] );
400 0 : c->SetTimeBin( (Int_t) padtime[2] );
401 0 : }
402 0 : }
403 :
404 :
405 : // loop over the input tracks: calculate dEdx and write output
406 :
407 : unsigned int outSize = 0;
408 0 : AliHLTFloat32_t *outPtr = ( AliHLTFloat32_t * )( outputPtr );
409 :
410 0 : for (const AliHLTComponentBlockData* pBlock=GetFirstInputBlock(kAliHLTDataTypeTrack|kAliHLTDataOriginTPC);
411 0 : pBlock!=NULL; pBlock=GetNextInputBlock()) {
412 :
413 0 : AliHLTComponentBlockData outBlock;
414 0 : FillBlockData( outBlock );
415 0 : outBlock.fOffset = outSize;
416 0 : outBlock.fSize = 0;
417 0 : outBlock.fDataType = kAliHLTDataTypedEdx|kAliHLTDataOriginTPC;
418 0 : outBlock.fSpecification = pBlock->fSpecification;
419 :
420 0 : AliHLTTracksData* dataPtr = ( AliHLTTracksData* ) pBlock->fPtr;
421 0 : int nTracks = dataPtr->fCount;
422 0 : AliHLTExternalTrackParam* currTrack = dataPtr->fTracklets;
423 0 : nInputTracks+=nTracks;
424 :
425 0 : for( int itr=0;
426 0 : itr<nTracks && ( (AliHLTUInt8_t *)currTrack < ((AliHLTUInt8_t *) pBlock->fPtr)+pBlock->fSize);
427 0 : itr++ ){
428 :
429 : // create an off-line track
430 0 : AliHLTGlobalBarrelTrack gb(*currTrack);
431 0 : AliTPCseed tTPC;
432 0 : tTPC.Set( gb.GetX(), gb.GetAlpha(),
433 0 : gb.GetParameter(), gb.GetCovariance() );
434 :
435 : // set the clusters
436 :
437 0 : for( UInt_t ic=0; ic<currTrack->fNPoints; ic++){
438 0 : UInt_t id = currTrack->fPointIDs[ic];
439 0 : int iSlice = AliHLTTPCSpacePointData::GetSlice(id);
440 0 : int iPatch = AliHLTTPCSpacePointData::GetPatch(id);
441 0 : int iCluster = AliHLTTPCSpacePointData::GetNumber(id);
442 0 : if( iSlice<0 || iSlice>36 || iPatch<0 || iPatch>5 ){
443 0 : HLTError("Corrupted TPC cluster Id: slice %d, patch %d, cluster %d",
444 : iSlice, iPatch,iCluster );
445 0 : continue;
446 : }
447 0 : AliTPCclusterMI *patchClusters = fPatchClusters[iSlice*6 + iPatch];
448 0 : if( !patchClusters ){
449 0 : HLTError("Clusters are missed for slice %d, patch %d",
450 : iSlice, iPatch );
451 0 : continue;
452 : }
453 0 : if( iCluster >= fNPatchClusters[iSlice*6 + iPatch] ){
454 0 : HLTError("TPC slice %d, patch %d: ClusterID==%d >= N Cluaters==%d ",
455 : iSlice, iPatch,iCluster, fNPatchClusters[iSlice*6 + iPatch] );
456 0 : continue;
457 : }
458 0 : AliTPCclusterMI *c = &(patchClusters[iCluster]);
459 :
460 0 : int sec = c->GetDetector();
461 0 : int row = c->GetRow();
462 0 : if ( sec >= 36 ) {
463 0 : row = row + AliHLTTPCGeometry::GetNRowLow();
464 0 : }
465 0 : tTPC.SetClusterPointer( row, c );
466 :
467 0 : AliTPCTrackerPoints::Point &point = *( (AliTPCTrackerPoints::Point*)tTPC.GetTrackPoint( row ) );
468 0 : tTPC.Propagate( TMath::DegToRad()*(sec%18*20.+10.), c->GetX(), fSolenoidBz );
469 0 : Double_t angle2 = tTPC.GetSnp()*tTPC.GetSnp();
470 0 : angle2 = (angle2<1) ?TMath::Sqrt(angle2/(1-angle2)) :10.;
471 0 : point.SetAngleY( angle2 );
472 0 : point.SetAngleZ( tTPC.GetTgl() );
473 0 : }
474 :
475 : // Cook dEdx
476 :
477 0 : if( outSize+3*sizeof( AliHLTFloat32_t ) > maxBufferSize ){
478 0 : HLTWarning( "Output buffer size %d exceed", maxBufferSize );
479 : iResult = -ENOSPC;
480 0 : break;
481 : }
482 :
483 : //Old method
484 : /*
485 : tTPC.CookdEdx( 0.02, 0.6 );
486 : outPtr[0] = tTPC.GetdEdx();
487 : outPtr[1] = tTPC.GetSDEDX(0);
488 : outPtr[2] = tTPC.GetNCDEDX(0);
489 : */
490 :
491 : //New method
492 0 : outPtr[0] = tTPC.CookdEdxAnalytical(0.02,0.6,1,0,159,0);
493 0 : outPtr[1] = 0.;
494 0 : outPtr[2] = 159;
495 :
496 0 : outPtr+=3;
497 0 : outSize+=3*sizeof( AliHLTFloat32_t );
498 0 : outBlock.fSize+=3*sizeof( AliHLTFloat32_t );
499 :
500 0 : unsigned int step = sizeof( AliHLTExternalTrackParam ) + currTrack->fNPoints * sizeof( unsigned int );
501 0 : currTrack = ( AliHLTExternalTrackParam* )( (( Byte_t * )currTrack) + step );
502 0 : }
503 0 : if( iResult<0 ) break;
504 0 : outputBlocks.push_back( outBlock );
505 0 : size = outSize;
506 0 : }
507 :
508 0 : timer.Stop();
509 0 : fStatTime += timer.RealTime();
510 0 : fStatNEvents++;
511 :
512 : // Set log level to "Warning" for on-line system monitoring
513 0 : int hz = ( int ) ( fStatTime > 1.e-10 ? fStatNEvents / fStatTime : 100000 );
514 :
515 0 : HLTInfo( "TPCdEdx : %d clusters and %d tracks processed; average time %d Hz",
516 : nInputClusters, nInputTracks, hz );
517 :
518 : return iResult;
519 0 : }
|
