make AliVVevent abstract (=0)
[u/mrichter/AliRoot.git] / HLT / global / AliHLTVertexFinderBase.cxx
CommitLineData
f078d002 1// $Id$
2
3//**************************************************************************
4//* This file is property of and copyright by the ALICE HLT Project *
5//* ALICE Experiment at CERN, All rights reserved. *
6//* *
7//* Primary Authors: Timur Pocheptsov <Timur.Pocheptsov@cern.ch> *
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/// @file AliHLTVertexFinderBase.cxx
20/// @author Timur Pocheptsov
21/// @date 2010-12-26
22/// @brief Base class for vertex finder components
23///
24
25#include <stdexcept>
26#include <cmath>
27
28//AliHLTExternalTrackParam uses typedes from Rtypes.h
29//but does not include it.
30#include "Rtypes.h"
31
32#include "AliHLTExternalTrackParam.h"
33#include "AliHLTGlobalBarrelTrack.h"
34#include "AliHLTVertexFinderBase.h"
35#include "AliExternalTrackParam.h"
36#include "AliESDVertex.h"
37#include "AliESDEvent.h"
38#include "AliESDtrack.h"
39#include "AliKFVertex.h"
40#include "AliESDv0.h"
41
42ClassImp(AliHLTVertexFinderBase)
43
44//_________________________________________________________________
45bool AliHLTVertexFinderBase::AliHLTTrackInfo::IsKFFinite()const
46{
47 //Check KF particle.
48
49 //In C99, 'isfinite' is a macro.
50 //But still, I add using directive, in case
51 //it's a function in std.
52 using namespace std;
53
54 for (int i = 0; i < 8; ++i)
55 if (!isfinite(fParticle.GetParameter(i)))
56 return false;
57 for (int i = 0; i < 36; ++i)//Not sure, probably, 27 is enough.
58 if (!isfinite(fParticle.GetCovariance(i)))
59 return false;
60
61 return true;
62}
63
64//_________________________________________________________________
bf0b4a59 65void AliHLTVertexFinderBase::ReadESDTracks(int posPID, int negPID)
f078d002 66{
67 //Try to read input tracks from AliESDEvent.
68 //Use esd track's index as "track id" (fID).
69 //IMPORTANT: fTrackInfos is _NOT_ cleared here,
70 //must be done externally (if you need to).
71
72 const TObject* iter = GetFirstInputObject(kAliHLTDataTypeESDObject);
73 for (; iter; iter = GetNextInputObject()) {
74 //Cast away constness before dynamic_cast.
75 AliESDEvent* esd = dynamic_cast<AliESDEvent*>((TObject*)iter);
76 if (!esd)//From original code. Not sure, if this is possible at all.
77 continue;
78
bf0b4a59 79 ReadESDTracks(esd, negPID, posPID);
80 //Only one esd event is taken.
81 break;
82 }
83}
f078d002 84
bf0b4a59 85//_________________________________________________________________
86void AliHLTVertexFinderBase::ReadESDTracks(AliESDEvent* esd, int posPID, int negPID)
87{
88 //Try to read input tracks from AliESDEvent.
89 //Use esd track's index as "track id" (fID).
90 //IMPORTANT: fTrackInfos is _NOT_ cleared here,
91 //must be done externally (if you need to).
92 esd->GetStdContent();
f078d002 93
bf0b4a59 94 const int nTracks = esd->GetNumberOfTracks();
95 if (nTracks)
96 fTrackInfos.reserve(nTracks + fTrackInfos.size());
f078d002 97
bf0b4a59 98 for (int i = 0; i < nTracks; ++i) {
99 AliESDtrack* pTrack = esd->GetTrack(i);
100 //This checks of track parameters are
101 //from the original global vertexer.
102 if (!pTrack)
103 continue;
104 if (pTrack->GetKinkIndex(0) > 0)
105 continue;
106 if (!(pTrack->GetStatus() & AliESDtrack::kTPCin))
107 continue;
f078d002 108
bf0b4a59 109 //i: track id, 211: pid, false: not used in primary, 0.: prim. deviation.
110 AliHLTTrackInfo newTrackInfo(i, *pTrack, pTrack->Charge() > 0 ? posPID : negPID, false, 0.);
111 if (!newTrackInfo.IsKFFinite())
112 continue;
113
114 fTrackInfos.push_back(newTrackInfo);
f078d002 115 }
116}
117
118//________________________________________________________________________
bf0b4a59 119void AliHLTVertexFinderBase::ReadHLTTracks(const AliHLTComponentDataType& blockType,
120 int posPID, int negPID)
f078d002 121{
122 //Read HLT tracks as input.
123 //IMPORTANT: fTrackInfos is _NOT_ cleared here,
124 //must be done externally (if you need it).
125
126 const AliHLTComponentBlockData* block = GetFirstInputBlock(blockType);
127 for (; block; block = GetNextInputBlock()) {
128 const AliHLTTracksData* dataPtr = static_cast<AliHLTTracksData*>(block->fPtr);
129 const AliHLTExternalTrackParam* hltTrk = dataPtr->fTracklets;
130
131 const int nTracks = dataPtr->fCount;
132 if (nTracks)
133 fTrackInfos.reserve(fTrackInfos.size() + nTracks);
134
135 for (int i = 0; i < nTracks; ++i) {
136 //Ugly conversion from one track format to another
137 //and to AliKFParitcle then, to be changed later.
138 AliHLTGlobalBarrelTrack tmpTrk(*hltTrk);
139
bf0b4a59 140 AliHLTTrackInfo newTrackInfo(hltTrk->fTrackID, tmpTrk,
141 tmpTrk.Charge() > 0 ? posPID : negPID, false, 0.);
f078d002 142 if (!newTrackInfo.IsKFFinite())
143 continue;
144
145 fTrackInfos.push_back(newTrackInfo);
146
147 const unsigned dSize = sizeof(AliHLTExternalTrackParam)
148 + hltTrk->fNPoints * sizeof(unsigned);
149 hltTrk = (AliHLTExternalTrackParam*)((char *)hltTrk + dSize);
150 }
151 }
152}
153
154//________________________________________________________________________
155void AliHLTVertexFinderBase::FillESD(AliESDEvent* esd, AliKFVertex* vtx,
156 const void* primData, const void* v0Data)
157{
158 //Put the output of a primary finder and v0 finder to the esd event.
159 //Code was taken from the original global vertexer.
160 //Code is an absolute mess.
161
162 typedef std::map<int, int> Map_t;
163 typedef Map_t::const_iterator MapIter_t;
164
165 //Map esdId -> esdTrackIndex.
166 Map_t mapId;
167
168 double params[3];
169 for (int i = 0; i < 3; ++i)
170 params[i] = vtx->Parameters()[i];
171 double cov[6];
172 for (int i = 0; i < 6; ++i)
173 cov[i] = vtx->CovarianceMatrix()[i];
174
175 AliESDVertex vESD(params, cov, vtx->GetChi2(), vtx->GetNContributors());
176 esd->SetPrimaryVertexTPC(&vESD);
177 esd->SetPrimaryVertexTracks(&vESD);
178
179 //Relate tracks to the primary vertex
180 const PrimaryFinderBlock* prim = static_cast<const PrimaryFinderBlock*>(primData);
181 const int nESDTracks = esd->GetNumberOfTracks();
182 std::vector<bool> constrainedToVtx(nESDTracks);
183
184 if (prim->fFitTracksFlag) {
185 for (int i = 0; i < nESDTracks; ++i) {
186 if (!esd->GetTrack(i))
187 continue;
188
189 mapId[esd->GetTrack(i)->GetID()] = i;
190 }
191
192 for (int i = 0; i < prim->fNPrimaryTracks; ++i) {
193 MapIter_t it = mapId.find(prim->fPrimTrackIds[i]);
194 if (it == mapId.end())
195 continue;
196 const int itr = it->second;
197 //100. is an argument for parameter maxd in AliESDtrack - cut on impact parameter.
198 esd->GetTrack(itr)->RelateToVertex(&vESD, esd->GetMagneticField(), 100.);
199 constrainedToVtx[itr] = true;
200 }
201 }
202
203 //Add v0s.
ac3cab27 204 if (v0Data) {
f078d002 205 const int* v0s = static_cast<const int*>(v0Data);
206 const int nV0s = v0s[0];
207 ++v0s;
208 for (int i = 0; i < nV0s; ++i) {
209 MapIter_t it = mapId.find(v0s[2 * i]);
210 if (it==mapId.end())
211 continue;
212 const int iTr = it->second;
213
214 it = mapId.find(v0s[2 * i + 1]);
215 if (it == mapId.end())
216 continue;
217 const int jTr = it->second;
218
219 AliESDv0 v0(*esd->GetTrack(iTr), iTr, *esd->GetTrack(jTr), jTr);
220 esd->AddV0(&v0);
221 // relate the tracks to the vertex
222 if (prim->fFitTracksFlag) {
223 if (constrainedToVtx[iTr] || constrainedToVtx[jTr])
224 continue;
225
226 double pos[3];
227 double sigma[3] = {.1, .1, .1};
228 v0.XvYvZv(pos);
229 AliESDVertex v0ESD(pos, sigma);
230 esd->GetTrack(iTr)->RelateToVertex(&v0ESD, esd->GetMagneticField(), 100.);
231 esd->GetTrack(jTr)->RelateToVertex(&v0ESD, esd->GetMagneticField(), 100.);
232 constrainedToVtx[iTr] = true;
233 constrainedToVtx[jTr] = true;
234 }
235 }
ac3cab27 236 }
f078d002 237}
238
239namespace AliHLTUtility
240{
241
242//_______________________________________________________________________
243Parameter::Parameter()
244 : fWasSet(false),
245 fConstraint(none),
246 fBool(0),
247 fInt(0),
248 fDouble(0),
249 fCompound(0),
250 fConstraintChecker(0)
251{
252 //Default ctor.
253}
254
255//_______________________________________________________________________
256Parameter::Parameter(bool* b)
257 : fWasSet(false),
258 fConstraint(none),
259 fBool(b),
260 fInt(0),
261 fDouble(0),
262 fCompound(0),
263 fConstraintChecker(0)
264{
265 //Parameter of type bool.
266}
267
268//_______________________________________________________________________
269Parameter::Parameter(int* i, Constraint c)
270 : fWasSet(false),
271 fConstraint(c),
272 fBool(0),
273 fInt(i),
274 fDouble(0),
275 fCompound(0),
276 fConstraintChecker(0)
277{
278 //Parameter of type int.
279 SetConstraintChecker();
280}
281
282//_______________________________________________________________________
283Parameter::Parameter(double* d, Constraint c)
284 : fWasSet(false),
285 fConstraint(c),
286 fBool(0),
287 fInt(0),
288 fDouble(d),
289 fCompound(0),
290 fConstraintChecker(0)
291{
292 //Parameter of type double.
293 SetConstraintChecker();
294}
295
296//_______________________________________________________________________
297Parameter::Parameter(CompoundType *ct)
298 : fWasSet(false),
299 fConstraint(none),
300 fBool(0),
301 fInt(0),
302 fDouble(0),
303 fCompound(ct),
304 fConstraintChecker(0)
305{
306 //Parameter of more complex user-defined type.
307 //All checks and conversions must be implemented
308 //by user of CompoundType.
309}
310
311//_______________________________________________________________________
312unsigned Parameter::SetParameter(unsigned argc, const char** argv, unsigned currPos)
313{
314 //Set parameter from command line tokens.
315
316 //It's up to compound parameter to parse.
317 if (fCompound)
318 return fCompound->SetParameter(argc, argv, currPos);
319
320 //Now, int, bool or double must be set from a string.
321 //Primitive checks are done here.
322 if (currPos == argc)
323 throw std::runtime_error("value expected");
324 if (fWasSet)
325 throw std::runtime_error("parameter was set already");
326
327 const TString val(argv[currPos]);
328 if (!val.Length())
329 throw std::runtime_error("value expected");
330
331 if (fBool) {
332 if (val.Length() != 1 || (val[0] != '0' && val[0] != '1'))
333 throw std::runtime_error("expected 0 or 1 for bool parameter");
334 *fBool = val[0] - '0';
335 } else if (fInt)
336 *fInt = val.Atoi();
337 else
338 *fDouble = val.Atof();
339
340 if (fConstraintChecker)
341 fConstraintChecker(*this);
342
343 fWasSet = true;
344
345 //For double, int and bool - only one element of argv is processed.
346 return 1;
347}
348
349//_______________________________________________________________________
350void Parameter::SetConstraintChecker()
351{
352 //Set function pointer.
353 if (fConstraint == positive)
354 fConstraintChecker = CheckPositive;
355 else if (fConstraint == nonNegative)
356 fConstraintChecker = CheckNonNegative;
357 //No constraints for bool or CompoundType.
358}
359
360//Aux. functions to check constraints.
361//_______________________________________________________________________
362void Parameter::CheckPositive(const Parameter& param)
363{
364 //val > 0
365 if (param.fInt) {
366 if (*param.fInt <= 0)
367 throw std::runtime_error("integer parameter must be positive");
368 } else if (param.fDouble) {
369 if (*param.fDouble <= 0.)//Hmmmm.
370 throw std::runtime_error("double parameter must be positive");
371 }
372 //For bool or CompoundType there is no check.
373}
374
375//_______________________________________________________________________
376void Parameter::CheckNonNegative(const Parameter& param)
377{
378 //val >= 0
379 if (param.fInt) {
380 if (*param.fInt < 0)
381 throw std::runtime_error("integer parameter must be non-negative");
382 } else if (param.fDouble) {
383 if (*param.fDouble < 0.)
384 throw std::runtime_error("double parameter must be non-negative");
385 }
386 //For bool or CompoundType there is no check.
387}
388
389
390//Command line arguments parser.
391
392//_______________________________________________________________________
393void CmdLineParser::Add(const TString& cmd, bool* b)
394{
395 //Command with bool parameter.
396 if (fParameters.find(cmd) == fParameters.end())
397 fParameters[cmd] = Parameter(b);
398}
399
400//_______________________________________________________________________
401void CmdLineParser::Add(const TString& cmd, int* i, Parameter::Constraint c)
402{
403 //Command with int parameter, with constraint.
404 if (fParameters.find(cmd) == fParameters.end())
405 fParameters[cmd] = Parameter(i, c);
406}
407
408//_______________________________________________________________________
409void CmdLineParser::Add(const TString& cmd, double* d, Parameter::Constraint c)
410{
411 //Coomand with a parameter of type double, possibly with a constraint.
412 if (fParameters.find(cmd) == fParameters.end())
413 fParameters[cmd] = Parameter(d, c);
414}
415
416//_______________________________________________________________________
417void CmdLineParser::Add(const TString& cmd, CompoundType *ct)
418{
419 //Command with a parameter of compound type.
420 if (fParameters.find(cmd) == fParameters.end())
421 fParameters[cmd] = Parameter(ct);
422}
423
424//_______________________________________________________________________
425void CmdLineParser::IgnoreCommand(const TString& cmd, unsigned nParams)
426{
427 //Command and number of parameters to ignore.
428 //I do not check, if it's in map already.
429 fIgnored[cmd] = nParams;
430}
431
432//_______________________________________________________________________
433int CmdLineParser::Parse(unsigned argc, const char** argv, unsigned currPos)
434{
435 //Parse. Returns number of processed arguments.
436 //Negative number - error (as it's in HLT)
437 fError.Clear();
438
439 unsigned nProcessed = 0;
440 while (currPos < argc) {
441 //Command.
442 const TString command(argv[currPos]);
443 //Check, if this command must be ignored.
444 SkipMapIter_t skipCmd = fIgnored.find(command);
445 if (skipCmd != fIgnored.end()) {
446 //Command and its parameters must be skipped
447 nProcessed += 1 + skipCmd->second;
448 currPos += 1 + skipCmd->second;
449 continue;
450 }
451
452 //Corresponding parameter.
453 MapIter_t it = fParameters.find(command);
454 if (it == fParameters.end()) {
455 fError = "Unknown command: " + command;
456 return -1;
457 }
458
459 ++currPos;
460 ++nProcessed;
461
462 try {
463 const unsigned n = it->second.SetParameter(argc, argv, currPos);
464 nProcessed += n;
465 currPos += n;
466 } catch (const std::runtime_error& e) {
467 fError = "Command " + command + " error: " + e.what();
468 return -1;
469 }
470 }
471
472 return nProcessed;
473}
474
475}