]>
Commit | Line | Data |
---|---|---|
4d5ee3db | 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: S.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 | /** @file AliHLTGlobalVertexerComponent.cxx | |
18 | @author Sergey Gorbunov | |
19 | @brief Component for reconstruct primary vertex and V0's | |
20 | */ | |
21 | ||
4d5ee3db | 22 | #include "AliHLTGlobalVertexerComponent.h" |
d9386025 | 23 | #include "AliHLTDataTypes.h" |
24 | #include "AliHLTExternalTrackParam.h" | |
25 | #include "AliHLTGlobalBarrelTrack.h" | |
4d5ee3db | 26 | #include "AliCDBEntry.h" |
27 | #include "AliCDBManager.h" | |
28 | #include <TFile.h> | |
29 | #include <TString.h> | |
30 | #include "TObjString.h" | |
31 | #include "TObjArray.h" | |
4d5ee3db | 32 | #include "AliESDEvent.h" |
33 | #include "AliESDtrack.h" | |
34 | #include "AliESDVertex.h" | |
35 | #include "AliESDv0.h" | |
36 | #include "AliHLTMessage.h" | |
4d5ee3db | 37 | #include "TMath.h" |
38 | #include "AliKFParticle.h" | |
39 | #include "AliKFVertex.h" | |
9222a93a | 40 | #include "TStopwatch.h" |
4d5ee3db | 41 | |
a7f38894 | 42 | using namespace std; |
43 | ||
4d5ee3db | 44 | /** ROOT macro for the implementation of ROOT specific class methods */ |
45 | ClassImp(AliHLTGlobalVertexerComponent) | |
46 | ||
47 | AliHLTGlobalVertexerComponent::AliHLTGlobalVertexerComponent() | |
48 | : | |
d9386025 | 49 | fNTracks(0), |
4d5ee3db | 50 | fTrackInfos(0), |
57a4102f | 51 | fPrimaryVtx(), |
608cfbda | 52 | fFitTracksToVertex(1), |
4d5ee3db | 53 | fConstrainedTrackDeviation(4.), |
54 | fV0DaughterPrimDeviation( 2.5 ), | |
55 | fV0PrimDeviation( 3.5 ), | |
56 | fV0Chi(3.5), | |
9222a93a | 57 | fV0DecayLengthInSigmas(3.), |
d9386025 | 58 | fV0TimeLimit(10.e-3), |
57a4102f | 59 | fBenchmark("GlobalVertexer") |
4d5ee3db | 60 | { |
61 | // see header file for class documentation | |
62 | // or | |
63 | // refer to README to build package | |
64 | // or | |
65 | // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt | |
66 | ||
67 | } | |
68 | ||
69 | AliHLTGlobalVertexerComponent::~AliHLTGlobalVertexerComponent() | |
70 | { | |
71 | // see header file for class documentation | |
72 | ||
dc546924 | 73 | if( fTrackInfos ) delete[] fTrackInfos; |
4d5ee3db | 74 | } |
75 | ||
76 | // Public functions to implement AliHLTComponent's interface. | |
77 | // These functions are required for the registration process | |
78 | ||
79 | const char* AliHLTGlobalVertexerComponent::GetComponentID() | |
80 | { | |
81 | // see header file for class documentation | |
82 | ||
83 | return "GlobalVertexer"; | |
84 | } | |
85 | ||
86 | void AliHLTGlobalVertexerComponent::GetInputDataTypes(AliHLTComponentDataTypeList& list) | |
87 | { | |
88 | // see header file for class documentation | |
89 | list.clear(); | |
d9386025 | 90 | list.push_back( kAliHLTDataTypeESDObject ); |
91 | list.push_back( kAliHLTDataTypeTrack|kAliHLTDataOriginITS ); | |
92 | list.push_back( kAliHLTDataTypeTrack|kAliHLTDataOriginTPC ); | |
4d5ee3db | 93 | } |
94 | ||
95 | AliHLTComponentDataType AliHLTGlobalVertexerComponent::GetOutputDataType() | |
96 | { | |
97 | // see header file for class documentation | |
98 | return kAliHLTMultipleDataType; | |
99 | } | |
100 | ||
101 | int AliHLTGlobalVertexerComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList) | |
102 | ||
103 | { | |
104 | // see header file for class documentation | |
105 | tgtList.clear(); | |
d9386025 | 106 | tgtList.push_back( kAliHLTDataTypeGlobalVertexer|kAliHLTDataOriginOut); |
107 | tgtList.push_back( kAliHLTDataTypeESDVertex|kAliHLTDataOriginOut); | |
108 | tgtList.push_back( kAliHLTDataTypeESDObject|kAliHLTDataOriginOut); | |
4d5ee3db | 109 | return tgtList.size(); |
110 | } | |
111 | ||
112 | void AliHLTGlobalVertexerComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier ) | |
113 | { | |
114 | // see header file for class documentation | |
115 | // XXX TODO: Find more realistic values. | |
116 | constBase = 80000; | |
117 | inputMultiplier = 2.; | |
118 | } | |
119 | ||
120 | AliHLTComponent* AliHLTGlobalVertexerComponent::Spawn() | |
121 | { | |
122 | // see header file for class documentation | |
123 | return new AliHLTGlobalVertexerComponent; | |
124 | } | |
125 | ||
126 | int AliHLTGlobalVertexerComponent::DoInit( int argc, const char** argv ) | |
127 | { | |
128 | // init | |
9222a93a | 129 | |
57a4102f | 130 | fBenchmark.Reset(); |
131 | fBenchmark.SetTimer(0,"total"); | |
132 | fBenchmark.SetTimer(1,"convert"); | |
133 | fBenchmark.SetTimer(2,"vprim"); | |
134 | fBenchmark.SetTimer(3,"v0"); | |
d9386025 | 135 | fV0TimeLimit = 10.e-3; |
4d5ee3db | 136 | |
4d5ee3db | 137 | int iResult=0; |
138 | TString configuration=""; | |
139 | TString argument=""; | |
140 | for (int i=0; i<argc && iResult>=0; i++) { | |
141 | argument=argv[i]; | |
142 | if (!configuration.IsNull()) configuration+=" "; | |
143 | configuration+=argument; | |
144 | } | |
145 | ||
146 | if (!configuration.IsNull()) { | |
147 | iResult=Configure(configuration.Data()); | |
148 | } | |
149 | ||
150 | return iResult; | |
151 | } | |
152 | ||
153 | int AliHLTGlobalVertexerComponent::DoDeinit() | |
154 | { | |
155 | // see header file for class documentation | |
156 | ||
dc546924 | 157 | if( fTrackInfos ) delete[] fTrackInfos; |
4d5ee3db | 158 | |
159 | fTrackInfos = 0; | |
608cfbda | 160 | fFitTracksToVertex = 1; |
4d5ee3db | 161 | fConstrainedTrackDeviation = 4.; |
162 | fV0DaughterPrimDeviation = 2.5 ; | |
163 | fV0PrimDeviation =3.5; | |
164 | fV0Chi = 3.5; | |
165 | fV0DecayLengthInSigmas = 3.; | |
d9386025 | 166 | fV0TimeLimit = 10.e-3; |
9222a93a | 167 | |
4d5ee3db | 168 | return 0; |
169 | } | |
170 | ||
d9386025 | 171 | int AliHLTGlobalVertexerComponent::DoEvent( const AliHLTComponentEventData& /*evtData*/, |
172 | AliHLTComponentTriggerData& /*trigData*/ ) | |
4d5ee3db | 173 | { |
9222a93a | 174 | //cout<<"AliHLTGlobalVertexerComponent::DoEvent called"<<endl; |
9be6e605 | 175 | |
4d5ee3db | 176 | if ( GetFirstInputBlock( kAliHLTDataTypeSOR ) || GetFirstInputBlock( kAliHLTDataTypeEOR ) ) |
177 | return 0; | |
178 | ||
dc546924 | 179 | int iResult = 0; |
57a4102f | 180 | //cout<<"SG: GlobalVertexer:: DoEvent called"<<endl; |
181 | fBenchmark.StartNewEvent(); | |
182 | fBenchmark.Start(0); | |
183 | ||
d9386025 | 184 | vector< AliExternalTrackParam > tracks; |
185 | vector< int > trackId; | |
186 | vector< pair<int,int> > v0s; | |
4d5ee3db | 187 | |
d9386025 | 188 | AliESDEvent *event = 0; |
4d5ee3db | 189 | |
57a4102f | 190 | for( const AliHLTComponentBlockData *i= GetFirstInputBlock(kAliHLTDataTypeESDObject); i!=NULL; i=GetNextInputBlock() ){ |
191 | fBenchmark.AddInput(i->fSize); | |
192 | } | |
193 | ||
194 | ||
d9386025 | 195 | for ( const TObject *iter = GetFirstInputObject(kAliHLTDataTypeESDObject); iter != NULL; iter = GetNextInputObject() ) { |
196 | event = dynamic_cast<AliESDEvent*>(const_cast<TObject*>( iter ) ); | |
dc546924 | 197 | if( !event ) continue; |
57a4102f | 198 | |
4d5ee3db | 199 | event->GetStdContent(); |
d9386025 | 200 | Int_t nESDTracks=event->GetNumberOfTracks(); |
201 | ||
202 | for (Int_t iTr=0; iTr<nESDTracks; iTr++){ | |
203 | AliESDtrack *pTrack = event->GetTrack(iTr); | |
204 | if( !pTrack ) continue; | |
205 | if( pTrack->GetKinkIndex(0)>0) continue; | |
206 | if( !( pTrack->GetStatus()&AliESDtrack::kTPCin ) ) continue; | |
207 | tracks.push_back(*pTrack); | |
208 | trackId.push_back(iTr); | |
209 | } | |
210 | break; | |
211 | } | |
212 | ||
213 | if( tracks.size()==0 ){ | |
214 | ||
215 | for (const AliHLTComponentBlockData* pBlock=GetFirstInputBlock(kAliHLTDataTypeTrack|kAliHLTDataOriginITS); | |
216 | pBlock!=NULL; pBlock=GetNextInputBlock()) { | |
57a4102f | 217 | |
218 | fBenchmark.AddInput(pBlock->fSize); | |
219 | ||
d9386025 | 220 | AliHLTTracksData* dataPtr = reinterpret_cast<AliHLTTracksData*>( pBlock->fPtr ); |
221 | int nTracks = dataPtr->fCount; | |
d9386025 | 222 | AliHLTExternalTrackParam* currOutTrack = dataPtr->fTracklets; |
223 | for( int itr=0; itr<nTracks; itr++ ){ | |
224 | AliHLTGlobalBarrelTrack t(*currOutTrack); | |
225 | tracks.push_back( t ); | |
9be6e605 | 226 | //!SW trackId.push_back( currOutTrack->fTrackID ); |
227 | trackId.push_back( itr ); | |
d9386025 | 228 | unsigned int dSize = sizeof( AliHLTExternalTrackParam ) + currOutTrack->fNPoints * sizeof( unsigned int ); |
229 | currOutTrack = ( AliHLTExternalTrackParam* )( (( Byte_t * )currOutTrack) + dSize ); | |
230 | } | |
231 | } | |
232 | } | |
4d5ee3db | 233 | |
d9386025 | 234 | if( tracks.size()==0 ){ |
235 | for (const AliHLTComponentBlockData* pBlock=GetFirstInputBlock(kAliHLTDataTypeTrack|kAliHLTDataOriginTPC); | |
236 | pBlock!=NULL; pBlock=GetNextInputBlock()) { | |
57a4102f | 237 | |
238 | fBenchmark.AddInput(pBlock->fSize); | |
239 | ||
d9386025 | 240 | AliHLTTracksData* dataPtr = reinterpret_cast<AliHLTTracksData*>( pBlock->fPtr ); |
9be6e605 | 241 | int nTracks = dataPtr->fCount; |
d9386025 | 242 | AliHLTExternalTrackParam* currOutTrack = dataPtr->fTracklets; |
243 | for( int itr=0; itr<nTracks; itr++ ){ | |
244 | AliHLTGlobalBarrelTrack t(*currOutTrack); | |
245 | tracks.push_back( t ); | |
9be6e605 | 246 | //!SW trackId.push_back( currOutTrack->fTrackID ); |
247 | trackId.push_back( itr ); | |
d9386025 | 248 | unsigned int dSize = sizeof( AliHLTExternalTrackParam ) + currOutTrack->fNPoints * sizeof( unsigned int ); |
249 | currOutTrack = ( AliHLTExternalTrackParam* )( (( Byte_t * )currOutTrack) + dSize ); | |
250 | } | |
251 | } | |
252 | } | |
253 | ||
254 | ||
255 | // primary vertex & V0's | |
608cfbda | 256 | |
d9386025 | 257 | AliKFParticle::SetField( GetBz() ); |
258 | ||
57a4102f | 259 | fBenchmark.Start(1); |
260 | ||
d9386025 | 261 | { //* Fill fTrackInfo array |
262 | ||
d9386025 | 263 | if( fTrackInfos ) delete[] fTrackInfos; |
264 | fTrackInfos = 0; | |
265 | fNTracks=tracks.size(); | |
266 | fTrackInfos = new AliESDTrackInfo[ fNTracks ]; | |
267 | for (Int_t iTr=0; iTr<fNTracks; iTr++){ | |
268 | AliESDTrackInfo &info = fTrackInfos[iTr]; | |
269 | info.fOK = 1; | |
270 | info.fPrimUsedFlag = 0; | |
271 | info.fParticle = AliKFParticle( tracks[iTr], 211 ); | |
57a4102f | 272 | for( int i=0; i<8; i++ ) if( !finite(info.fParticle.GetParameter(i)) ) info.fOK = 0; |
273 | for( int i=0; i<36; i++ ) if( !finite(info.fParticle.GetCovariance(i)) ) info.fOK = 0; | |
d9386025 | 274 | } |
d9386025 | 275 | } |
57a4102f | 276 | |
277 | fBenchmark.Stop(1); | |
278 | fBenchmark.Start(2); | |
d9386025 | 279 | FindPrimaryVertex(); |
57a4102f | 280 | fBenchmark.Stop(2); |
281 | fBenchmark.Start(3); | |
d9386025 | 282 | FindV0s( v0s ); |
57a4102f | 283 | fBenchmark.Stop(3); |
d9386025 | 284 | |
285 | int *buf = new int[sizeof(AliHLTGlobalVertexerData)/sizeof(int)+1 + fNTracks + 2*v0s.size()]; | |
286 | AliHLTGlobalVertexerData *data = reinterpret_cast<AliHLTGlobalVertexerData*>(buf); | |
287 | ||
288 | if( data) { // fill the output structure | |
289 | ||
290 | data->fFitTracksToVertex = fFitTracksToVertex; | |
291 | for( int i=0; i<3; i++ ) data->fPrimP[i] = fPrimaryVtx.Parameters()[i]; | |
292 | for( int i=0; i<6; i++ ) data->fPrimC[i] = fPrimaryVtx.CovarianceMatrix()[i]; | |
293 | data->fPrimChi2 = fPrimaryVtx.GetChi2(); | |
294 | data->fPrimNContributors = fPrimaryVtx.GetNContributors(); | |
295 | data->fNPrimTracks = 0; | |
296 | for( Int_t i = 0; i<fNTracks; i++ ){ | |
297 | if( !fTrackInfos[i].fPrimUsedFlag ) continue; | |
298 | if( fTrackInfos[i].fPrimDeviation > fConstrainedTrackDeviation ) continue; | |
299 | data->fTrackIndices[ (data->fNPrimTracks)++ ] = trackId[i]; | |
300 | } | |
301 | int *listV0 = data->fTrackIndices + data->fNPrimTracks; | |
302 | data->fNV0s = v0s.size(); | |
303 | for( int i=0; i<data->fNV0s; i++ ){ | |
304 | listV0[2*i] = trackId[v0s[i].first]; | |
305 | listV0[2*i+1] = trackId[v0s[i].second]; | |
306 | } | |
4d5ee3db | 307 | |
d9386025 | 308 | unsigned int blockSize = sizeof(AliHLTGlobalVertexerData) + (data->fNPrimTracks + 2*data->fNV0s)*sizeof(int); |
309 | ||
310 | iResult = PushBack( reinterpret_cast<void*>(data), blockSize, kAliHLTDataTypeGlobalVertexer|kAliHLTDataOriginOut ); | |
57a4102f | 311 | fBenchmark.AddOutput(blockSize); |
d9386025 | 312 | } |
313 | ||
314 | ||
315 | // output the vertex if found | |
316 | { | |
317 | if( iResult==0 && data && data->fPrimNContributors >=3 ){ | |
318 | AliESDVertex vESD( data->fPrimP, data->fPrimC, data->fPrimChi2, data->fPrimNContributors ); | |
319 | iResult = PushBack( (TObject*) &vESD, kAliHLTDataTypeESDVertex|kAliHLTDataOriginOut,0 ); | |
57a4102f | 320 | fBenchmark.AddOutput(GetLastObjectSize()); |
d9386025 | 321 | } |
322 | } | |
323 | ||
324 | // output the ESD event | |
325 | if( iResult==0 && event && data ){ | |
326 | FillESD( event, data ); | |
dc546924 | 327 | iResult = PushBack( event, kAliHLTDataTypeESDObject|kAliHLTDataOriginOut, 0); |
57a4102f | 328 | fBenchmark.AddOutput(GetLastObjectSize()); |
4d5ee3db | 329 | } |
d9386025 | 330 | |
331 | delete[] buf; | |
332 | ||
57a4102f | 333 | fBenchmark.Stop(0); |
334 | HLTInfo(fBenchmark.GetStatistics()); | |
dc546924 | 335 | return iResult; |
4d5ee3db | 336 | } |
337 | ||
d9386025 | 338 | |
4d5ee3db | 339 | int AliHLTGlobalVertexerComponent::Configure(const char* arguments) |
340 | { | |
341 | ||
342 | int iResult=0; | |
343 | if (!arguments) return iResult; | |
344 | ||
345 | TString allArgs=arguments; | |
346 | TString argument; | |
347 | ||
348 | TObjArray* pTokens=allArgs.Tokenize(" "); | |
349 | int bMissingParam=0; | |
350 | ||
351 | if (pTokens) { | |
352 | for (int i=0; i<pTokens->GetEntries() && iResult>=0; i++) { | |
353 | argument=((TObjString*)pTokens->At(i))->GetString(); | |
354 | if (argument.IsNull()) continue; | |
355 | ||
356 | if (argument.CompareTo("-fitTracksToVertex")==0) { | |
357 | if ((bMissingParam=(++i>=pTokens->GetEntries()))) break; | |
358 | HLTInfo("fitTracksToVertex is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data()); | |
608cfbda | 359 | fFitTracksToVertex=((TObjString*)pTokens->At(i))->GetString().Atoi(); |
4d5ee3db | 360 | continue; |
361 | } | |
4d5ee3db | 362 | else if (argument.CompareTo("-constrainedTrackDeviation")==0) { |
363 | if ((bMissingParam=(++i>=pTokens->GetEntries()))) break; | |
364 | HLTInfo("constrainedTrackDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data()); | |
365 | fConstrainedTrackDeviation=((TObjString*)pTokens->At(i))->GetString().Atof(); | |
366 | continue; | |
367 | } | |
368 | else if (argument.CompareTo("-v0DaughterPrimDeviation")==0) { | |
369 | if ((bMissingParam=(++i>=pTokens->GetEntries()))) break; | |
370 | HLTInfo("v0DaughterPrimDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data()); | |
371 | fV0DaughterPrimDeviation=((TObjString*)pTokens->At(i))->GetString().Atof(); | |
372 | continue; | |
373 | } | |
374 | else if (argument.CompareTo("-v0PrimDeviation")==0) { | |
375 | if ((bMissingParam=(++i>=pTokens->GetEntries()))) break; | |
376 | HLTInfo("v0PrimDeviation is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data()); | |
377 | fV0PrimDeviation=((TObjString*)pTokens->At(i))->GetString().Atof(); | |
378 | continue; | |
379 | } | |
380 | else if (argument.CompareTo("-v0Chi")==0) { | |
381 | if ((bMissingParam=(++i>=pTokens->GetEntries()))) break; | |
382 | HLTInfo("v0Chi is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data()); | |
383 | fV0Chi=((TObjString*)pTokens->At(i))->GetString().Atof(); | |
384 | continue; | |
385 | } | |
386 | else if (argument.CompareTo("-v0DecayLengthInSigmas")==0) { | |
387 | if ((bMissingParam=(++i>=pTokens->GetEntries()))) break; | |
388 | HLTInfo("v0DecayLengthInSigmas is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data()); | |
389 | fV0DecayLengthInSigmas=((TObjString*)pTokens->At(i))->GetString().Atof(); | |
390 | continue; | |
391 | } | |
9222a93a | 392 | else if (argument.CompareTo("-v0MinEventRate")==0) { |
393 | if ((bMissingParam=(++i>=pTokens->GetEntries()))) break; | |
394 | HLTInfo("Minimum event rate for V0 finder is set set to: %s", ((TObjString*)pTokens->At(i))->GetString().Data()); | |
395 | Double_t rate = ((TObjString*)pTokens->At(i))->GetString().Atof(); | |
396 | fV0TimeLimit = (rate >0 ) ?1./rate :60; // 1 minute maximum time | |
397 | continue; | |
398 | } | |
4d5ee3db | 399 | else { |
400 | HLTError("unknown argument %s", argument.Data()); | |
401 | iResult=-EINVAL; | |
402 | break; | |
403 | } | |
404 | } | |
405 | delete pTokens; | |
406 | } | |
407 | if (bMissingParam) { | |
408 | HLTError("missing parameter for argument %s", argument.Data()); | |
409 | iResult=-EINVAL; | |
410 | } | |
411 | ||
412 | return iResult; | |
413 | } | |
414 | ||
412559b1 | 415 | int AliHLTGlobalVertexerComponent::Reconfigure(const char* /*cdbEntry*/, const char* /*chainId*/) |
4d5ee3db | 416 | { |
417 | // see header file for class documentation | |
9222a93a | 418 | |
419 | return 0; // no CDB path is set so far | |
5af75aae | 420 | /* |
9222a93a | 421 | int iResult=0; |
4d5ee3db | 422 | const char* path="HLT/ConfigTPC/KryptonHistoComponent"; |
423 | const char* defaultNotify=""; | |
424 | if (cdbEntry) { | |
425 | path=cdbEntry; | |
426 | defaultNotify=" (default)"; | |
427 | } | |
428 | if (path) { | |
429 | HLTInfo("reconfigure from entry %s%s, chain id %s", path, defaultNotify,(chainId!=NULL && chainId[0]!=0)?chainId:"<none>"); | |
412559b1 | 430 | AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(path); |
4d5ee3db | 431 | if (pEntry) { |
432 | TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject()); | |
433 | if (pString) { | |
434 | HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data()); | |
435 | iResult=Configure(pString->GetString().Data()); | |
436 | } else { | |
437 | HLTError("configuration object \"%s\" has wrong type, required TObjString", path); | |
438 | } | |
439 | } else { | |
440 | HLTError("can not fetch object \"%s\" from CDB", path); | |
441 | } | |
442 | } | |
443 | ||
444 | return iResult; | |
5af75aae | 445 | */ |
4d5ee3db | 446 | } |
447 | ||
448 | ||
2e56583f | 449 | struct AliHLTGlobalVertexerDeviation |
450 | { | |
451 | int fI; // track index | |
452 | float fD; // deviation from vertex | |
453 | ||
454 | bool operator<(const AliHLTGlobalVertexerDeviation &a) const { return fD<a.fD; } | |
455 | }; | |
456 | ||
4d5ee3db | 457 | |
d9386025 | 458 | void AliHLTGlobalVertexerComponent::FindPrimaryVertex() |
4d5ee3db | 459 | { |
460 | //* Find event primary vertex | |
461 | ||
4d5ee3db | 462 | fPrimaryVtx.Initialize(); |
d9386025 | 463 | //fPrimaryVtx.SetBeamConstraint(fESD->GetDiamondX(),fESD->GetDiamondY(),0, |
464 | //TMath::Sqrt(fESD->GetSigma2DiamondX()),TMath::Sqrt(fESD->GetSigma2DiamondY()),5.3); | |
900f84fe | 465 | |
466 | // select rough region (in sigmas) in which the vertex could be found, all tracks outside these limits are rejected | |
467 | // from the primary vertex finding | |
d9386025 | 468 | fPrimaryVtx.SetBeamConstraint( 0, 0, 0, 3., 3., 5.3 ); |
2e56583f | 469 | |
470 | const AliKFParticle **vSelected = new const AliKFParticle*[fNTracks]; //* Selected particles for the vertex fit | |
471 | AliHLTGlobalVertexerDeviation *dev = new AliHLTGlobalVertexerDeviation[fNTracks]; | |
472 | ||
4d5ee3db | 473 | Int_t nSelected = 0; |
2e56583f | 474 | |
d9386025 | 475 | for( Int_t i = 0; i<fNTracks; i++){ |
4d5ee3db | 476 | if(!fTrackInfos[i].fOK ) continue; |
477 | //if( fESD->GetTrack(i)->GetTPCNcls()<60 ) continue; | |
478 | const AliKFParticle &p = fTrackInfos[i].fParticle; | |
479 | Double_t chi = p.GetDeviationFromVertex( fPrimaryVtx ); | |
480 | if( chi > fConstrainedTrackDeviation ) continue; | |
2e56583f | 481 | dev[nSelected].fI = i; |
482 | dev[nSelected].fD = chi; | |
4d5ee3db | 483 | vSelected[nSelected] = &(fTrackInfos[i].fParticle); |
4d5ee3db | 484 | nSelected++; |
2e56583f | 485 | } |
486 | ||
2e56583f | 487 | // fit |
488 | ||
2e56583f | 489 | while( nSelected>2 ){ |
490 | ||
491 | //* Primary vertex finder with rejection of outliers | |
492 | ||
493 | for( Int_t i = 0; i<nSelected; i++){ | |
494 | vSelected[i] = &(fTrackInfos[dev[i].fI].fParticle); | |
495 | } | |
44fe6239 | 496 | |
497 | double xv = fPrimaryVtx.GetX(); | |
498 | double yv = fPrimaryVtx.GetY(); | |
900f84fe | 499 | double zv = fPrimaryVtx.GetZ(); // values from previous iteration of calculations |
2e56583f | 500 | fPrimaryVtx.Initialize(); |
501 | fPrimaryVtx.SetBeamConstraint( 0, 0, 0, 3., 3., 5.3 ); | |
44fe6239 | 502 | fPrimaryVtx.SetVtxGuess( xv, yv, zv ); |
2e56583f | 503 | |
900f84fe | 504 | fPrimaryVtx.Construct( vSelected, nSelected, 0, -1, 1 ); // refilled for every iteration |
44fe6239 | 505 | |
2e56583f | 506 | for( Int_t it=0; it<nSelected; it++ ){ |
507 | const AliKFParticle &p = fTrackInfos[dev[it].fI].fParticle; | |
508 | if( nSelected <= 20 ){ | |
900f84fe | 509 | AliKFVertex tmp = fPrimaryVtx - p; // exclude the current track from the sample and recalculate the vertex |
2e56583f | 510 | dev[it].fD = p.GetDeviationFromVertex( tmp ); |
511 | } else { | |
512 | dev[it].fD = p.GetDeviationFromVertex( fPrimaryVtx ); | |
513 | } | |
514 | } | |
900f84fe | 515 | sort(dev,dev+nSelected); // sort tracks with increasing chi2 (used for rejection) |
2e56583f | 516 | |
900f84fe | 517 | int nRemove = (int) ( 0.3*nSelected ); //remove 30% of the tracks (done for performance, only if there are more than 20 tracks) |
518 | if( nSelected - nRemove <=20 ) nRemove = 1; // removal based on the chi2 of every track | |
2e56583f | 519 | int firstRemove = nSelected - nRemove; |
520 | while( firstRemove<nSelected ){ | |
521 | if( dev[firstRemove].fD >= fConstrainedTrackDeviation ) break; | |
522 | firstRemove++; | |
523 | } | |
524 | if( firstRemove>=nSelected ) break; | |
2e56583f | 525 | nSelected = firstRemove; |
526 | } | |
527 | ||
2e56583f | 528 | for( Int_t i = 0; i<fNTracks; i++){ |
529 | fTrackInfos[i].fPrimUsedFlag = 0; | |
530 | } | |
531 | ||
900f84fe | 532 | if( nSelected < 3 ){ // no vertex for fewer than 3 contributors |
2e56583f | 533 | fPrimaryVtx.NDF() = -3; |
534 | fPrimaryVtx.Chi2() = 0; | |
535 | nSelected = 0; | |
536 | } | |
537 | ||
4d5ee3db | 538 | for( Int_t i = 0; i<nSelected; i++){ |
44fe6239 | 539 | AliESDTrackInfo &info = fTrackInfos[dev[i].fI]; |
540 | info.fPrimUsedFlag = 1; | |
541 | info.fPrimDeviation = dev[i].fD; | |
4d5ee3db | 542 | } |
543 | ||
d9386025 | 544 | for( Int_t i = 0; i<fNTracks; i++ ){ |
4d5ee3db | 545 | AliESDTrackInfo &info = fTrackInfos[i]; |
44fe6239 | 546 | if( info.fPrimUsedFlag ) continue; |
4d5ee3db | 547 | info.fPrimDeviation = info.fParticle.GetDeviationFromVertex( fPrimaryVtx ); |
548 | } | |
4d5ee3db | 549 | |
4d5ee3db | 550 | delete[] vSelected; |
2e56583f | 551 | delete[] dev; |
4d5ee3db | 552 | } |
553 | ||
554 | ||
2e56583f | 555 | |
f0d12ea5 | 556 | void AliHLTGlobalVertexerComponent::FindV0s( vector<pair<int,int> > &v0s ) |
4d5ee3db | 557 | { |
558 | //* V0 finder | |
559 | ||
4d5ee3db | 560 | AliKFVertex &primVtx = fPrimaryVtx; |
561 | if( primVtx.GetNContributors()<3 ) return; | |
562 | ||
9222a93a | 563 | TStopwatch timer; |
564 | Int_t statN = 0; | |
565 | Bool_t run = 1; | |
4d5ee3db | 566 | |
d9386025 | 567 | for( Int_t iTr = 0; iTr<fNTracks && run; iTr++ ){ //* first daughter |
4d5ee3db | 568 | |
569 | AliESDTrackInfo &info = fTrackInfos[iTr]; | |
570 | if( !info.fOK ) continue; | |
571 | if( info.fParticle.GetQ() >0 ) continue; | |
572 | if( info.fPrimDeviation < fV0DaughterPrimDeviation ) continue; | |
573 | ||
d9386025 | 574 | for( Int_t jTr = 0; jTr<fNTracks; jTr++ ){ //* second daughter |
9222a93a | 575 | |
576 | ||
4d5ee3db | 577 | AliESDTrackInfo &jnfo = fTrackInfos[jTr]; |
578 | if( !jnfo.fOK ) continue; | |
579 | if( jnfo.fParticle.GetQ() < 0 ) continue; | |
580 | if( jnfo.fPrimDeviation < fV0DaughterPrimDeviation ) continue; | |
9222a93a | 581 | |
582 | // check the time once a while... | |
583 | ||
584 | if( (++statN)%100 ==0 ){ | |
585 | if( timer.RealTime()>= fV0TimeLimit ){ run = 0; break; } | |
57a4102f | 586 | timer.Continue(); |
9222a93a | 587 | } |
588 | ||
589 | //* check if the particles fit | |
590 | ||
591 | if( info.fParticle.GetDeviationFromParticle(jnfo.fParticle) > fV0Chi ) continue; | |
592 | ||
4d5ee3db | 593 | //* construct V0 mother |
594 | ||
595 | AliKFParticle v0( info.fParticle, jnfo.fParticle ); | |
596 | ||
597 | //* check V0 Chi^2 | |
598 | ||
599 | if( v0.GetChi2()<0 || v0.GetChi2() > fV0Chi*fV0Chi*v0.GetNDF() ) continue; | |
600 | ||
601 | //* subtruct daughters from primary vertex | |
602 | ||
603 | AliKFVertex primVtxCopy = primVtx; | |
604 | ||
605 | if( info.fPrimUsedFlag ){ | |
606 | if( primVtxCopy.GetNContributors()<=2 ) continue; | |
607 | primVtxCopy -= info.fParticle; | |
608 | } | |
609 | if( jnfo.fPrimUsedFlag ){ | |
610 | if( primVtxCopy.GetNContributors()<=2 ) continue; | |
611 | primVtxCopy -= jnfo.fParticle; | |
612 | } | |
613 | //* Check v0 Chi^2 deviation from primary vertex | |
614 | ||
615 | if( v0.GetDeviationFromVertex( primVtxCopy ) > fV0PrimDeviation ) continue; | |
616 | ||
617 | //* Add V0 to primary vertex to improve the primary vertex resolution | |
618 | ||
619 | primVtxCopy += v0; | |
620 | ||
621 | //* Set production vertex for V0 | |
622 | ||
623 | v0.SetProductionVertex( primVtxCopy ); | |
624 | ||
625 | //* Get V0 decay length with estimated error | |
626 | ||
627 | Double_t length, sigmaLength; | |
628 | if( v0.GetDecayLength( length, sigmaLength ) ) continue; | |
629 | ||
630 | //* Reject V0 if it decays too close[sigma] to the primary vertex | |
631 | ||
632 | if( length < fV0DecayLengthInSigmas*sigmaLength ) continue; | |
4d5ee3db | 633 | |
d9386025 | 634 | //* keep v0 |
635 | ||
636 | v0s.push_back(pair<int,int>(iTr,jTr)); | |
4d5ee3db | 637 | } |
638 | } | |
4d5ee3db | 639 | } |
640 | ||
d9386025 | 641 | |
642 | ||
643 | ||
644 | void AliHLTGlobalVertexerComponent::FillESD( AliESDEvent *event, AliHLTGlobalVertexerData *data | |
645 | ) | |
646 | { | |
647 | //* put output of a vertexer to the esd event | |
648 | ||
649 | Int_t nESDTracks = event->GetNumberOfTracks(); | |
650 | ||
651 | const int *listPrim = data->fTrackIndices; | |
652 | const int *listV0 = data->fTrackIndices + data->fNPrimTracks; | |
653 | ||
654 | std::map<int,int> mapId; | |
655 | bool *constrainedToVtx = new bool[nESDTracks]; | |
656 | ||
657 | for( int i=0; i<nESDTracks; i++ ){ | |
658 | constrainedToVtx[i] = 0; | |
659 | if( !event->GetTrack(i) ) continue; | |
660 | mapId[ event->GetTrack(i)->GetID() ] = i; | |
661 | } | |
662 | ||
663 | if( data->fPrimNContributors >=3 ){ | |
664 | ||
665 | AliESDVertex vESD( data->fPrimP, data->fPrimC, data->fPrimChi2, data->fPrimNContributors ); | |
39a4205e | 666 | event->SetPrimaryVertexTPC( &vESD ); |
d9386025 | 667 | event->SetPrimaryVertexTracks( &vESD ); |
668 | ||
669 | // relate tracks to the primary vertex | |
670 | ||
671 | if( data->fFitTracksToVertex ){ | |
672 | for( Int_t i = 0; i<data->fNPrimTracks; i++ ){ | |
673 | Int_t id = listPrim[ i ]; | |
674 | map<int,int>::iterator it = mapId.find(id); | |
675 | if( it==mapId.end() ) continue; | |
676 | Int_t itr = it->second; | |
677 | event->GetTrack(itr)->RelateToVertex( &vESD, event->GetMagneticField(),100. ); | |
678 | constrainedToVtx[ itr ] = 1; | |
679 | } | |
680 | } | |
681 | } | |
682 | ||
683 | //* add ESD v0s and relate tracks to v0s | |
684 | ||
685 | ||
686 | for( int i=0; i<data->fNV0s; i++ ){ | |
687 | ||
688 | Int_t id1 = listV0[ 2*i ]; | |
689 | Int_t id2 = listV0[ 2*i + 1]; | |
690 | map<int,int>::iterator it = mapId.find(id1); | |
691 | if( it==mapId.end() ) continue; | |
692 | Int_t iTr = it->second; | |
693 | it = mapId.find(id2); | |
694 | if( it==mapId.end() ) continue; | |
695 | Int_t jTr = it->second; | |
696 | ||
697 | AliESDv0 v0( *event->GetTrack( iTr ), iTr, *event->GetTrack( jTr ), jTr ); | |
698 | event->AddV0( &v0 ); | |
699 | ||
700 | // relate the tracks to the vertex | |
701 | ||
702 | if( data->fFitTracksToVertex ){ | |
703 | if( constrainedToVtx[iTr] || constrainedToVtx[jTr] ) continue; | |
704 | double pos[3]; | |
705 | double sigma[3] = {.1,.1,.1}; | |
706 | v0.XvYvZv(pos); | |
707 | AliESDVertex vESD(pos, sigma); | |
708 | event->GetTrack(iTr)->RelateToVertex( &vESD, event->GetMagneticField(),100. ); | |
709 | event->GetTrack(jTr)->RelateToVertex( &vESD, event->GetMagneticField(),100. ); | |
710 | constrainedToVtx[iTr] = 1; | |
711 | constrainedToVtx[jTr] = 1; | |
712 | } | |
713 | } | |
714 | ||
715 | delete[] constrainedToVtx; | |
716 | } |