]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/UPGRADE/AliITSUTrackerGlo.cxx
A skeleton of a comparison macro for pileup studies
[u/mrichter/AliRoot.git] / ITS / UPGRADE / AliITSUTrackerGlo.cxx
CommitLineData
32d38de2 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/*
17 TO CHECK
18 GetBz usage
19
20 */
21
22//-------------------------------------------------------------------------
23// Implementation of the ITS Upgrade tracker mother class.
24//-------------------------------------------------------------------------
25#include <TTree.h>
26#include <Riostream.h>
27#include <TMath.h>
76390254 28#include <TFile.h>
32d38de2 29
30#include "AliITSUTrackerGlo.h"
31#include "AliESDEvent.h"
32#include "AliESDtrack.h"
33#include "AliITSURecoDet.h"
34#include "AliITSURecoSens.h"
35#include "AliITSUReconstructor.h"
36#include "AliITSReconstructor.h"
37#include "AliITSUSeed.h"
c61e50c3 38#include "AliITSUClusterPix.h"
cb50e082 39#include "AliITSUGeomTGeo.h"
bdb92766 40#include "AliCodeTimer.h"
9cdcba2c 41#include "AliRefArray.h"
32d38de2 42using namespace AliITSUAux;
43using namespace TMath;
44
c8d1f258 45//----------------- tmp stuff -----------------
c8d1f258 46
32d38de2 47ClassImp(AliITSUTrackerGlo)
e2d2481c 48
49const Double_t AliITSUTrackerGlo::fgkToler = 1e-6;// tolerance for layer on-surface check
50
51
32d38de2 52//_________________________________________________________________________
53AliITSUTrackerGlo::AliITSUTrackerGlo(AliITSUReconstructor* rec)
54: fReconstructor(rec)
55 ,fITS(0)
e1ef49ad 56 ,fMatLUT(0)
0091e9f0 57 ,fCurrESDtrack(0)
0ddbf657 58 ,fCurrESDtrMClb(kDummyLabel)
32d38de2 59 ,fCurrMass(kPionMass)
0ddbf657 60 ,fCountProlongationTrials(0)
61 ,fCountITSin(0)
62 ,fCountITSout(0)
63 ,fCountITSrefit(0)
d99994b2 64 ,fNTracksESD(0)
716ccba7 65 ,fHypStore(100)
6cd80116 66 ,fLayerMaxCandidates(1000)
67 ,fLayerCandidates(0)
38997c3c 68 ,fNBranchesAdded(0)
69 ,fNCandidatesAdded(0)
716ccba7 70 ,fCurrHyp(0)
6cd80116 71 ,fWorkHyp(0)
c61e50c3 72 ,fSeedsPool("AliITSUSeed",0)
b8b59e05 73 ,fFreeSeedsID(0)
9cdcba2c 74 ,fESDIndex(0)
b8b59e05 75 ,fNFreeSeeds(0)
76 ,fLastSeedID(0)
6cd80116 77 ,fNLrActive(0)
42c3d4bd 78 ,fDefTrackConds(0)
79 ,fCurrTrackCond(0)
c03e4f8a 80 ,fCurrActLrID(-1)
81 ,fCurrLayer(0)
781c79ac 82 ,fTrackPhaseID(-1)
83 ,fCurrPassID(-1)
e1ef49ad 84 ,fUseMatLUT(kFALSE)
08419930 85#ifdef _ITSU_TUNING_MODE_
86 ,fCHistoArrCorr(0)
87 ,fCHistoArrFake(0)
53870004 88#endif
32d38de2 89{
90 // Default constructor
91 if (rec) Init(rec);
92}
93
94//_________________________________________________________________________
95AliITSUTrackerGlo::~AliITSUTrackerGlo()
96{
97 // Default destructor
98 //
9cdcba2c 99 delete[] fLayerCandidates;
100 if (fWorkHyp) fWorkHyp->SetTPCSeed(0); // this hypothesis does not own the seed
101 delete fWorkHyp;
e1ef49ad 102 delete fMatLUT;
32d38de2 103 //
08419930 104#ifdef _ITSU_TUNING_MODE_
105 if (fCHistoArrCorr || fCHistoArrFake) {
76390254 106 TFile* ctrOut = TFile::Open("itsuControlHistos.root","recreate");
107 ctrOut->cd();
108 AliInfo("Storing control histos");
53870004 109 // ctrOut->WriteObject(fCHistoArr,"controlH","kSingleKey");
08419930 110 if (fCHistoArrCorr) {
111 fCHistoArrCorr->Write();
112 delete fCHistoArrCorr;
113 }
114 if (fCHistoArrFake) {
115 fCHistoArrFake->Write();
116 delete fCHistoArrFake;
117 }
76390254 118 ctrOut->Close();
119 delete ctrOut;
08419930 120 fCHistoArrCorr = 0;
121 fCHistoArrFake = 0;
76390254 122 }
123#endif
124 //
32d38de2 125}
126
127//_________________________________________________________________________
128void AliITSUTrackerGlo::Init(AliITSUReconstructor* rec)
129{
130 // init with external reconstructor
c61e50c3 131 //
3d4dc3e2 132 fITS = rec->GetITSInterface();
6cd80116 133 fNLrActive = fITS->GetNLayersActive();
6cd80116 134 fWorkHyp = new AliITSUTrackHyp(fNLrActive);
e1ef49ad 135 // create material lookup table
136 const int kNTest = 1000;
137 const double kStepsPerCM=5;
138 fMatLUT = new AliITSUMatLUT(fITS->GetRMin(),fITS->GetRMax(),Nint(kStepsPerCM*(fITS->GetRMax()-fITS->GetRMin())));
139 double zmn = 1e6;
140 for (int ilr=fITS->GetNLayers();ilr--;) {
141 AliITSURecoLayer* lr = fITS->GetLayer(ilr);
142 if (zmn>Abs(lr->GetZMin())) zmn = Abs(lr->GetZMin());
143 if (zmn>Abs(lr->GetZMax())) zmn = Abs(lr->GetZMax());
144 }
145 fMatLUT->FillData(kNTest,-zmn,zmn);
c61e50c3 146 //
6cd80116 147 if (fLayerMaxCandidates<1) fLayerMaxCandidates = 1000;
148 fLayerCandidates = new AliITSUSeed*[fLayerMaxCandidates];
c61e50c3 149 fSeedsPool.ExpandCreateFast(1000); // RS TOCHECK
b8b59e05 150 fFreeSeedsID.Set(1000);
9cdcba2c 151 fESDIndex.Set(1000);
c61e50c3 152 //
32d38de2 153}
154
42c3d4bd 155//_________________________________________________________________________
156void AliITSUTrackerGlo::CreateDefaultTrackCond()
157{
158 // creates default tracking conditions to be used when recoparam does not provide them
6cd80116 159
42c3d4bd 160 AliITSUTrackCond* cond = new AliITSUTrackCond();
161 //
6cd80116 162 cond->SetNLayers(fNLrActive);
163 cond->AddNewCondition(fNLrActive);
70f61d86 164 cond->AddGroupPattern( 0xffff, fNLrActive ); // require all layers hit
ee54014a 165 cond->Init();
42c3d4bd 166 //
167 fDefTrackConds.AddLast(cond);
168 //
169 AliInfo("Created conditions: ");
170 for (int i=0;i<fDefTrackConds.GetEntriesFast();i++) fDefTrackConds[i]->Print();
171 //
172}
173
174
32d38de2 175//_________________________________________________________________________
176Int_t AliITSUTrackerGlo::Clusters2Tracks(AliESDEvent *esdEv)
177{
178 //
bdb92766 179 AliCodeTimerAuto("",0);
76390254 180 SetTrackingPhase(kClus2Tracks);
181 //
08419930 182#ifdef _ITSU_TUNING_MODE_
781c79ac 183 if (!fCHistoArrCorr) fCHistoArrCorr = BookControlHistos("Corr");
184 if (!fCHistoArrFake) fCHistoArrFake = BookControlHistos("Fake");
76390254 185#endif
c03e4f8a 186 static int evID = 0;
9cdcba2c 187 static TArrayF esdTrPt(fESDIndex.GetSize());
bdb92766 188 //
42c3d4bd 189 TObjArray *trackConds = 0;
32d38de2 190 //
0ddbf657 191 fCountProlongationTrials = 0;
192 fCountITSin = 0;
193 fCountITSout = 0;
194 fCountITSrefit = 0;
195 //
b8b59e05 196 ResetSeedsPool();
d99994b2 197 fNTracksESD = esdEv->GetNumberOfTracks();
198 AliInfo(Form("Will try to find prolongations for %d tracks",fNTracksESD));
42c3d4bd 199 int nTrackCond = AliITSUReconstructor::GetRecoParam()->GetNTrackingConditions();
e1ef49ad 200 fUseMatLUT = AliITSUReconstructor::GetRecoParam()->GetUseMatLUT(fTrackPhaseID);
42c3d4bd 201 if (nTrackCond<1) {
202 if (!fDefTrackConds.GetEntriesFast()) {
203 AliInfo("No tracking conditions found in recoparams, creating default one requesting all layers hit");
204 CreateDefaultTrackCond();
205 }
206 trackConds = &fDefTrackConds;
207 nTrackCond = trackConds->GetEntriesFast();
208 }
209 else trackConds = AliITSUReconstructor::GetRecoParam()->GetTrackingConditions();
210 //
15b02d69 211 static Bool_t first = kTRUE;
212 if (first) {
213 AliITSUReconstructor::GetRecoParam()->Print();
214 first = kFALSE;
215 }
216 fHypStore.Delete();
d99994b2 217 if (fHypStore.GetSize()<fNTracksESD) fHypStore.Expand(fNTracksESD+100);
15b02d69 218 //
219 fITS->ProcessClusters();
220 //
08419930 221#ifdef _ITSU_TUNING_MODE_
53870004 222 FlagSplitClusters(); // tmp RS
bdb92766 223#endif
224 //
225 // the tracks will be reconstructed in decreasing pt order, sort them
d99994b2 226 if (fESDIndex.GetSize()<fNTracksESD) {
227 fESDIndex.Set(fNTracksESD+200);
228 esdTrPt.Set(fNTracksESD+200);
bdb92766 229 }
9cdcba2c 230 int* esdTrackIndex = fESDIndex.GetArray();
bdb92766 231 float* trPt = esdTrPt.GetArray();
891efa00 232 for (int itr=fNTracksESD;itr--;) {
233 AliESDtrack* esdTr = esdEv->GetTrack(itr);
0e40ce66 234 esdTr->SetStatus(AliESDtrack::kITSupg); // Flag all tracks as participating in the ITSup reco
891efa00 235 trPt[itr] = esdTr->Pt();
236 }
d99994b2 237 Sort(fNTracksESD,trPt,esdTrackIndex,kTRUE);
53870004 238 //
42c3d4bd 239 for (int icnd=0;icnd<nTrackCond;icnd++) {
781c79ac 240 fCurrPassID = icnd;
42c3d4bd 241 fCurrTrackCond = (AliITSUTrackCond*)trackConds->UncheckedAt(icnd);
ee54014a 242 if (!fCurrTrackCond->IsInitDone()) fCurrTrackCond->Init();
42c3d4bd 243 // select ESD tracks to propagate
d99994b2 244 for (int itr=0;itr<fNTracksESD;itr++) {
bdb92766 245 int trID = esdTrackIndex[itr];
246 fCurrESDtrack = esdEv->GetTrack(trID);
0ddbf657 247 fCurrESDtrMClb = fCurrESDtrack->GetLabel();
248 //
0e40ce66 249
250 if (!NeedToProlong(fCurrESDtrack,trID)) continue; // are we interested in this track?
c03e4f8a 251 /*
9cdcba2c 252 // if specific tracks should be checked, create a global array int wtc[] = {evITS*10000+trID, ...};
c03e4f8a 253 Bool_t dbg = kFALSE;
254 int nwtc = sizeof(wtc)/sizeof(int);
255
256 for (int k=0;k<nwtc;k++) {
bdb92766 257 if (wtc[k]==evID*10000+trID) {
c03e4f8a 258 dbg = kTRUE;
c03e4f8a 259 printf("At esdTr: %d MC: %d\n",wtc[k],fCurrESDtrMClb);
c03e4f8a 260 break;
261 }
262 }
263 AliLog::SetClassDebugLevel("AliITSUTrackerGlo",dbg ? 10:0);
264 */
e1ef49ad 265#ifdef _ITSU_DEBUG_
1d26da6d 266 AliDebug(1,Form("Processing track %d(esd%d) | M=%.3f Pt=%.3f | MCLabel: %d",itr,trID,fCurrESDtrack->GetMassForTracking(kTRUE),fCurrESDtrack->Pt(),fCurrESDtrMClb));//RS
e1ef49ad 267#endif
bdb92766 268 FindTrack(fCurrESDtrack, trID);
42c3d4bd 269 }
270 //
e1ef49ad 271#ifdef _ITSU_DEBUG_
38997c3c 272 if (AliDebugLevelClass()>+2) {
15b02d69 273 AliInfo(Form("SeedsPool: %d, BookedUpTo: %d, free: %d",fSeedsPool.GetSize(),fSeedsPool.GetEntriesFast(),fNFreeSeeds));
274 fHypStore.Print();
275 }
e1ef49ad 276#endif
42c3d4bd 277 FinalizeHypotheses();
891efa00 278#ifdef _ITSU_TUNING_MODE_
d99994b2 279 CheckClusterUsage(); //!!RS
891efa00 280#endif
c03e4f8a 281 //
282 //AliLog::SetClassDebugLevel("AliITSUTrackerGlo",0); // in case wtc array was defined, uncomment this
42c3d4bd 283 }
716ccba7 284 //
d99994b2 285 AliInfo(Form("%d ITSin for %d tried TPC seeds out of %d ESD tracks\n",fCountITSin,fCountProlongationTrials,fNTracksESD));
c03e4f8a 286 evID++;
32d38de2 287 return 0;
288}
289
290//_________________________________________________________________________
c51c3124 291Int_t AliITSUTrackerGlo::PropagateBack(AliESDEvent *esdEv)
32d38de2 292{
293 //
c51c3124 294 // Do outward fits in ITS
bdb92766 295 AliCodeTimerAuto("",0);
c51c3124 296 //
e7d83d38 297 SetTrackingPhase(kPropBack);
d99994b2 298 fNTracksESD = esdEv->GetNumberOfTracks();
e1ef49ad 299 fUseMatLUT = AliITSUReconstructor::GetRecoParam()->GetUseMatLUT(fTrackPhaseID);
300#ifdef _ITSU_DEBUG_
d99994b2 301 AliDebug(1,Form("Will propagate back %d tracks",fNTracksESD));
e1ef49ad 302#endif
c51c3124 303 //
304 double bz0 = GetBz();
305 Double_t xyzTrk[3],xyzVtx[3]={GetX(),GetY(),GetZ()};
9cdcba2c 306 AliITSUTrackHyp dummyTr;
c51c3124 307 const double kWatchStep=10.; // for larger steps watch arc vs segment difference
d4da4017 308 Double_t times[AliPID::kSPECIESC];
c51c3124 309 //
d99994b2 310 for (int itr=0;itr<fNTracksESD;itr++) {
0ddbf657 311 fCurrESDtrack = esdEv->GetTrack(itr);
312 fCurrESDtrMClb = fCurrESDtrack->GetLabel();
c51c3124 313 // Start time integral and add distance from current position to vertex
0ddbf657 314 if (fCurrESDtrack->IsOn(AliESDtrack::kITSout)) continue; // already done
c64acbbe 315 if (!fCurrESDtrack->IsOn(AliESDtrack::kTPCin)) continue; // skip ITS s.a.
c51c3124 316 //
0ddbf657 317 fCurrESDtrack->GetXYZ(xyzTrk);
c51c3124 318 Double_t dst = 0.; // set initial track lenght, tof
319 {
320 double dxs = xyzTrk[0] - xyzVtx[0];
321 double dys = xyzTrk[1] - xyzVtx[1];
7ce7f445 322 // double dzs = xyzTrk[2] - xyzVtx[2];
c51c3124 323 // RS: for large segment steps d use approximation of cicrular arc s by
324 // s = 2R*asin(d/2R) = d/p asin(p) \approx d/p (p + 1/6 p^3) = d (1+1/6 p^2)
325 // where R is the track radius, p = d/2R = 2C*d (C is the abs curvature)
326 // Hence s^2/d^2 = (1+1/6 p^2)^2
327 dst = dxs*dxs + dys*dys;
328 if (dst > kWatchStep*kWatchStep) { // correct circular part for arc/segment factor
0ddbf657 329 double crv = Abs(fCurrESDtrack->GetC(bz0));
c51c3124 330 double fcarc = 1.+crv*crv*dst/6.;
331 dst *= fcarc*fcarc;
332 }
7ce7f445 333 // RS: temporary hack since we don't have SPD vertex:
334 // dst += dzs*dzs;
335 dst *= 1+fCurrESDtrack->GetTgl()*fCurrESDtrack->GetTgl();
c51c3124 336 dst = Sqrt(dst);
337 }
338 //
0ddbf657 339 fCurrESDtrack->SetStatus(AliESDtrack::kTIME);
c51c3124 340 //
0ddbf657 341 if (!fCurrESDtrack->IsOn(AliESDtrack::kITSin)) { // case of tracks w/o ITS prolongation: just set the time integration
342 dummyTr.AliExternalTrackParam::operator=(*fCurrESDtrack);
c51c3124 343 dummyTr.StartTimeIntegral();
e2d2481c 344 dummyTr.AddTimeStep(dst);
fc9b31a7 345 dummyTr.GetIntegratedTimes(times,AliPID::kSPECIESC);
0ddbf657 346 fCurrESDtrack->SetIntegratedTimes(times);
347 fCurrESDtrack->SetIntegratedLength(dummyTr.GetIntegratedLength());
c51c3124 348 continue;
349 }
350 //
9cdcba2c 351 fCurrHyp = GetTrackHyp(itr);
70cb7fe4 352 fCurrMass = fCurrHyp->GetMass();
9cdcba2c 353 fCurrHyp->StartTimeIntegral();
354 fCurrHyp->AddTimeStep(dst);
9cdcba2c 355 fCurrHyp->ResetCovariance(10000);
5e3587be 356 int nclFit = 0;
357 double chi2 = RefitTrack(fCurrHyp,fITS->GetRMax(),nclFit);
9cdcba2c 358 if (chi2>0) { // propagate to exit from the ITS/TPC screen
5e3587be 359 int ndf = nclFit*2-5;
70cb7fe4 360 if (ndf>0) chi2 /= ndf;
9cdcba2c 361 fCurrHyp->SetChi2(chi2);
362 UpdateESDTrack(fCurrHyp,AliESDtrack::kITSout);
68a0f687 363 }
e7d83d38 364 else {
e1ef49ad 365#ifdef _ITSU_DEBUG_
0ddbf657 366 AliDebug(2,Form("Refit Failed for track %d | ESDtrack#%d (MClb:%d)",itr,fCurrESDtrack->GetID(),fCurrESDtrMClb));
e1ef49ad 367#endif
9cdcba2c 368 //fCurrHyp->AliExternalTrackParam::Print();
369 //fCurrHyp->GetWinner()->Print();
e7d83d38 370 }
c51c3124 371 //
372 }
32d38de2 373 //
0ddbf657 374 AliInfo(Form("%d ITSout in %d ITSin tracks for %d tried TPC seeds out of %d ESD tracks\n",
d99994b2 375 fCountITSout,fCountITSin,fCountProlongationTrials,fNTracksESD));
0ddbf657 376 //
32d38de2 377 return 0;
378}
379
380//_________________________________________________________________________
e7d83d38 381Int_t AliITSUTrackerGlo::RefitInward(AliESDEvent *esdEv)
32d38de2 382{
383 //
e7d83d38 384 // refit the tracks inward, using current cov.matrix
bdb92766 385 AliCodeTimerAuto("",0);
e7d83d38 386 //
387 SetTrackingPhase(kRefitInw);
e1ef49ad 388 fUseMatLUT = AliITSUReconstructor::GetRecoParam()->GetUseMatLUT(fTrackPhaseID);
d99994b2 389 fNTracksESD = esdEv->GetNumberOfTracks();
b515ad7e 390 // AliLog::SetClassDebugLevel("AliITSUTrackerGlo",10);
391
d99994b2 392 AliDebug(1,Form("Will refit inward %d tracks",fNTracksESD));
c7e63626 393 // Bool_t uselogMS = AliExternalTrackParam::GetUseLogTermMS();
394 // AliExternalTrackParam::SetUseLogTermMS(kTRUE);
e7d83d38 395 //
d99994b2 396 for (int itr=0;itr<fNTracksESD;itr++) {
0ddbf657 397 fCurrESDtrack = esdEv->GetTrack(itr);
398 fCurrESDtrMClb = fCurrESDtrack->GetLabel();
e7d83d38 399 // Start time integral and add distance from current position to vertex
0ddbf657 400 UInt_t trStat = fCurrESDtrack->GetStatus();
e7d83d38 401 if ( !(trStat & AliESDtrack::kITSout) ) continue;
402 if ( trStat & AliESDtrack::kITSrefit ) continue; // already done
c64acbbe 403 if ( !(trStat & AliESDtrack::kTPCin) ) continue; // skip ITS s.a.
e7d83d38 404 if ( (trStat & AliESDtrack::kTPCout) && !(trStat & AliESDtrack::kTPCrefit) ) continue;
405 //
7ce7f445 406 fCurrHyp = GetTrackHyp(itr);
70cb7fe4 407 fCurrMass = fCurrHyp->GetMass();
7ce7f445 408 *fCurrHyp = *fCurrESDtrack; // fetch current ESDtrack kinematics
70cb7fe4 409 //
5e3587be 410 int nclFit = 0;
eb18cdc6 411 // fCurrHyp->ResetCovariance(1000); // in case we want to do SA refit
5e3587be 412 double chi2 = RefitTrack(fCurrHyp,fITS->GetRMin(),nclFit);
70cb7fe4 413 if (chi2>0) { // propagate up to inside radius of the beam pipe
5e3587be 414 int ndf = nclFit*2-5;
415 if (ndf>0) chi2 /= ndf;
9cdcba2c 416 fCurrHyp->SetChi2(chi2);
417 UpdateESDTrack(fCurrHyp,AliESDtrack::kITSrefit);
e7d83d38 418 }
419 else {
e1ef49ad 420#ifdef _ITSU_DEBUG_
0ddbf657 421 AliDebug(2,Form("Refit Failed for track %d |ESDtrack#%d (MClb:%d)",itr,fCurrESDtrack->GetID(),fCurrESDtrMClb));
e1ef49ad 422#endif
e7d83d38 423 }
424 }
c7e63626 425 // AliExternalTrackParam::SetUseLogTermMS(uselogMS);
32d38de2 426 //
0ddbf657 427 AliInfo(Form("%d ITSrefit in %d ITSout in %d ITSin tracks for %d tried TPC seeds out of %d ESD tracks\n",
d99994b2 428 fCountITSrefit,fCountITSout,fCountITSin,fCountProlongationTrials,fNTracksESD));
0ddbf657 429 //
b515ad7e 430 // AliLog::SetClassDebugLevel("AliITSUTrackerGlo",0);
32d38de2 431 return 0;
432}
433
434//_________________________________________________________________________
435Int_t AliITSUTrackerGlo::LoadClusters(TTree * treeRP)
436{
437 // read from tree (if pointer provided) or directly from the ITS reco interface
438 //
439 return fReconstructor->LoadClusters(treeRP);
440}
441
442//_________________________________________________________________________
443void AliITSUTrackerGlo::UnloadClusters()
444{
445 //
9473f43a 446 // Remove clusters from the memory
32d38de2 447 //
9473f43a 448 AliITSURecoDet *det=fReconstructor->GetITSInterface();
449 Int_t nlayers=det->GetNLayersActive();
450 for (Int_t i=0; i<nlayers; i++) {
451 TClonesArray *clusters=*(det->GetLayerActive(i)->GetClustersAddress());
452 clusters->Delete();
453 }
32d38de2 454}
9473f43a 455
32d38de2 456//_________________________________________________________________________
457AliCluster * AliITSUTrackerGlo::GetCluster(Int_t /*index*/) const
458{
459 //
460 // To be implemented
461 //
462
463 Info("GetCluster","To be implemented");
464 return 0x0;
465}
466
467//_________________________________________________________________________
0e40ce66 468Bool_t AliITSUTrackerGlo::NeedToProlong(AliESDtrack* esdTr, Int_t esdID)
32d38de2 469{
470 // do we need to match this track to ITS?
471 //
472 static double bz = GetBz();
0e40ce66 473 //
474 // is it already prolonged
475 if (esdTr->IsOn(AliESDtrack::kITSin)) return kFALSE;
476 AliITSUTrackHyp* hyp = GetTrackHyp(esdID); // in case there is unfinished hypothesis from prev.pass, clean it
477 if (hyp) {
478 CleanHypothesis(hyp);
479 if (hyp->GetSkip()) return kFALSE; // need to skip this
480 return kTRUE;
481 }
482 //
32d38de2 483 if (!esdTr->IsOn(AliESDtrack::kTPCin) ||
484 esdTr->IsOn(AliESDtrack::kTPCout) ||
485 esdTr->IsOn(AliESDtrack::kITSin) ||
486 esdTr->GetKinkIndex(0)>0) return kFALSE;
487 //
c61e50c3 488 if (esdTr->Pt()<AliITSUReconstructor::GetRecoParam()->GetMinPtForProlongation()) return kFALSE;
32d38de2 489 //
490 float dtz[2];
491 esdTr->GetDZ(GetX(),GetY(),GetZ(),bz,dtz);
492 // if track is not V0 candidata but has large offset wrt IP, reject it. RS TOCHECK
c61e50c3 493 if ( !(esdTr->GetV0Index(0)>0 && dtz[0]>AliITSUReconstructor::GetRecoParam()->GetMaxDforV0dghtrForProlongation())
494 && (Abs(dtz[0])>AliITSUReconstructor::GetRecoParam()->GetMaxDForProlongation() ||
495 Abs(dtz[1])>AliITSUReconstructor::GetRecoParam()->GetMaxDZForProlongation())) return kFALSE;
32d38de2 496 //
cb50e082 497 // if (esdTr->Pt()<3) return kFALSE;//RS
32d38de2 498 return kTRUE;
499}
500
501//_________________________________________________________________________
716ccba7 502void AliITSUTrackerGlo::FindTrack(AliESDtrack* esdTr, Int_t esdID)
32d38de2 503{
504 // find prolongaion candidates finding for single seed
505 //
3e4e3c23 506 AliITSUSeed seedUC; // copy of the seed from the upper layer
507 AliITSUSeed seedT; // transient seed between the seedUC and new prolongation hypothesis
508 //
d99994b2 509 // AliLog::SetClassDebugLevel("AliITSUTrackerGlo",10);
6cd80116 510 AliITSUTrackHyp* hypTr = InitHypothesis(esdTr,esdID); // initialize prolongations hypotheses tree
9cdcba2c 511 if (!hypTr || hypTr->GetSkip()) return;
512 //
d99994b2 513 double bz = GetBz();
6cd80116 514 fCurrHyp = fWorkHyp;
515 fCurrHyp->InitFrom(hypTr);
32d38de2 516 //
68b7631d 517 AliITSURecoSens *hitSens[AliITSURecoLayer::kMaxSensMatching];
f8832015 518 //
d99994b2 519 int ilaUp = fNLrActive; // previous active layer really checked (some may be excluded!)
32d38de2 520 //
d99994b2 521 for (int ila=fNLrActive;ila--;ilaUp=fCurrActLrID) {
522 if (fCurrTrackCond->IsLayerExcluded(ila)) continue;
c03e4f8a 523 fCurrActLrID = ila;
524 fCurrLayer = fITS->GetLayerActive(ila);
525 Bool_t noClSharing = fCurrTrackCond->GetClSharing(ila)==0;
3e4e3c23 526 //
527 // for the outermost layer the seed is created from the ESD track
6cd80116 528 int nSeedsUp = (ilaUp==fNLrActive) ? 1 : fCurrHyp->GetNSeeds(ilaUp);
38997c3c 529 int maxNBranches = fCurrTrackCond->GetMaxBranches(ila);
530 int maxNCandidates = fCurrTrackCond->GetMaxCandidates(ila);
3e4e3c23 531 //
32d38de2 532 for (int isd=0;isd<nSeedsUp;isd++) {
3e4e3c23 533 AliITSUSeed* seedU;
6cd80116 534 if (ilaUp==fNLrActive) {
3e4e3c23 535 seedU = 0;
70cb7fe4 536 seedUC.InitFromSeed(fCurrHyp->GetTPCSeed()); // init with TPC seed from ref.R
3e4e3c23 537 }
538 else {
539 seedU = fCurrHyp->GetSeed(ilaUp,isd); // seed on prev.active layer to prolong
c03e4f8a 540 if (seedU->IsKilled()) continue;
3e4e3c23 541 seedUC = *seedU; // its copy will be prolonged
b8b59e05 542 seedUC.SetParent(seedU);
3e4e3c23 543 }
943e1898 544 seedUC.ResetFMatrix(); // reset the matrix for propagation to next layer
32d38de2 545 // go till next active layer
e1ef49ad 546#ifdef _ITSU_DEBUG_
0ddbf657 547 AliDebug(2,Form("working on Lr:%d Seed:%d of %d for esdID=%d (MClb:%d) | pT=%.3f",ila,isd,nSeedsUp,esdID,fCurrESDtrMClb,seedUC.Pt()));
e1ef49ad 548#endif
70cb7fe4 549 if (!TransportToLayer(&seedUC, fITS->GetLrIDActive(ilaUp), fITS->GetLrIDActive(ila)) ) { // external seed already prolonged
32d38de2 550 //
e1ef49ad 551#ifdef _ITSU_DEBUG_
0ddbf657 552 AliDebug(2,Form("Transport failed | esdID=%d (MClb:%d)",esdID,fCurrESDtrMClb));
e1ef49ad 553#endif
32d38de2 554 // Check if the seed satisfies to track definition
b8b59e05 555 if (NeedToKill(&seedUC,kTransportFailed) && seedU) KillSeed(seedU,kTRUE);
32d38de2 556 continue; // RS TODO: decide what to do with tracks stopped on higher layers w/o killing
557 }
f8832015 558 if (!GetRoadWidth(&seedUC, ila)) { // failed to find road width on the layer
b8b59e05 559 if (NeedToKill(&seedUC,kRWCheckFailed) && seedU) KillSeed(seedU,kTRUE);
32d38de2 560 continue;
561 }
5785a9d8 562 /*
563 //RS toremove
564 int mcquest = -1;
c03e4f8a 565 if (!seedUC.ContainsFake() && AliDebugLevelClass()>2) {
5785a9d8 566 mcquest = fCurrESDtrMClb;
567 seedUC.Print("etp");
568 printf("FSParams: "); for (int ip=0;ip<kNTrImpData;ip++) printf("%+e ",fTrImpData[ip]); printf("\n");
569 }
570 //
c03e4f8a 571 int nsens = fCurrLayer->FindSensors(&fTrImpData[kTrPhi0], hitSens, mcquest); // find detectors which may be hit by the track
5785a9d8 572 */
c03e4f8a 573 int nsens = fCurrLayer->FindSensors(&fTrImpData[kTrPhi0], hitSens); // find detectors which may be hit by the track
e1ef49ad 574#ifdef _ITSU_DEBUG_
0ddbf657 575 AliDebug(2,Form("Will check %d sensors on lr:%d | esdID=%d (MClb:%d)",nsens,ila,esdID,fCurrESDtrMClb));
e1ef49ad 576#endif
32d38de2 577 //
38997c3c 578 seedUC.SetLr(ila);
579 //
32d38de2 580 for (int isn=nsens;isn--;) {
c61e50c3 581 AliITSURecoSens* sens = hitSens[isn];
38997c3c 582 int ncl = sens->GetNClusters();
583 if (!ncl) continue;
584 seedT = seedUC;
f8832015 585 //
c03e4f8a 586 // We need to propagate the seed to sensor on fCurrLayer staying in the frame of the sensor from prev.layer,
943e1898 587 // since the transport matrix should be defined in this frame.
588 double xs; // X in the TF of current seed, corresponding to intersection with sensor plane
53870004 589 if (!seedT.GetTrackingXAtXAlpha(sens->GetXTF(),sens->GetPhiTF(),bz, xs)) {
e1ef49ad 590#ifdef _ITSU_DEBUG_
53870004 591 if (AliDebugLevelClass()>2) {
592 printf("Failed on GetTrackingXAtXAlpha: X=%.4f alp=%.3f\n",sens->GetXTF(),sens->GetPhiTF());
593 seedT.Print("etp");
594 }
e1ef49ad 595#endif
53870004 596 continue;
597 }
598 if (xs<seedT.GetX()) {
599 if (!PropagateSeed(&seedT,xs,fCurrMass)) continue;
600 }
601 else { // some low precision tracks may hit the sensor plane outside of the layer radius
e1ef49ad 602#ifdef _ITSU_DEBUG_
53870004 603 if (AliDebugLevelClass()>2) {
604 if (!seedT.ContainsFake()) {
605 printf("WRONG PROP on L%d, sens%d of %d: %.4f > %.4f\n",ila,isn,nsens,xs,seedT.GetX());
606 sens->Print();
607 seedT.Print("etp");
608 }
609 }
e1ef49ad 610#endif
53870004 611 if (!seedT.PropagateParamOnlyTo(xs,bz)) continue;
612 }
8b16dbae 613 // if (!seedT.PropagateToX(xs,bz)) continue;
716ccba7 614 // if (!seedT.Rotate(sens->GetPhiTF())) continue;
615 if (!seedT.RotateToAlpha(sens->GetPhiTF())) continue;
943e1898 616 //
f8832015 617 int clID0 = sens->GetFirstClusterId();
cb50e082 618 for (int icl=ncl;icl--;) { // don't change the order, clusters are sorted
c03e4f8a 619 //AliLog::SetClassDebugLevel("AliITSUTrackerGlo",10);
620 int clID = clID0 + icl;
621 if (noClSharing && fCurrLayer->GetCluster(clID)->IsClusterUsed()) continue;
f8832015 622 int res = CheckCluster(&seedT,ila,clID0+icl);
c03e4f8a 623 //AliLog::SetClassDebugLevel("AliITSUTrackerGlo", 0);
f8832015 624 //
625 if (res==kStopSearchOnSensor) break; // stop looking on this sensor
626 if (res==kClusterNotMatching) continue; // cluster does not match
627 // cluster is matching and it was added to the hypotheses tree
c61e50c3 628 }
32d38de2 629 }
38997c3c 630 // cluster search is done. Do we need to have a version of this seed skipping current layer
631 if (!NeedToKill(&seedUC,kMissingCluster)) {
b2d3bc3a 632 AliITSUSeed* seedSkp = NewSeedFromPool(&seedUC);
097e4cdf 633 if (nsens) { // to do: make penalty to account for probability to miss the cluster for good reason
634 double penalty = -fCurrTrackCond->GetMissPenalty(ila);
635 seedSkp->SetChi2Cl(penalty);
636 }
08419930 637 if (seedSkp->GetChi2GloNrm()>fCurrTrackCond->GetMaxChi2GloNrm(ila)) {
638 MarkSeedFree(seedSkp);
639 }
640 else AddSeedBranch(seedSkp);
716ccba7 641 }
38997c3c 642 // transfer the new branches of the seed to the hypothesis container
643 if (fNBranchesAdded) ValidateAllowedBranches(maxNBranches);
644 //
32d38de2 645 }
38997c3c 646 if (fNCandidatesAdded) ValidateAllowedCandidates(ila,maxNCandidates);
647 // ((TObjArray*)fCurrHyp->GetLayerSeeds(ila))->Sort();
e1ef49ad 648#ifdef _ITSU_DEBUG_
38997c3c 649 if (AliDebugLevelClass()>2) { //RS
15b02d69 650 printf(">>> All hypotheses on lr %d: \n",ila);
651 for (int ih=0;ih<fCurrHyp->GetNSeeds(ila);ih++) {printf(" #%3d ",ih); fCurrHyp->GetSeed(ila,ih)->Print();}
652 }
e1ef49ad 653#endif
15b02d69 654 //
655 /*
716ccba7 656 for (int ih=0;ih<fCurrHyp->GetNSeeds(ila);ih++) {
716ccba7 657 if (ila!=0) continue;
658 double vecL[5] = {0};
659 double matL[15] = {0};
660 AliITSUSeed* sp = fCurrHyp->GetSeed(ila,ih);
661 while(sp->GetParent()) {
662 sp->Smooth(vecL,matL);
6cd80116 663 if (sp->GetLayerID()>=fNLrActive-1) break;
716ccba7 664 sp = (AliITSUSeed*)sp->GetParent();
665 }
666 */
32d38de2 667 }
668 //
6cd80116 669 SaveReducedHypothesesTree(hypTr);
e1ef49ad 670#ifdef _ITSU_DEBUG_
6cd80116 671 if (AliDebugLevelClass()>1) {
672 printf("\nSaved hypotheses for esdTrack %d (MCLab:%d)\n",esdID,fCurrESDtrMClb);
673 hypTr->Print("l");
674 }
e1ef49ad 675#endif
6cd80116 676 fCurrHyp = 0;
d99994b2 677 // AliLog::SetClassDebugLevel("AliITSUTrackerGlo",0);
716ccba7 678 //
32d38de2 679}
680
681//_________________________________________________________________________
6cd80116 682AliITSUTrackHyp* AliITSUTrackerGlo::InitHypothesis(AliESDtrack *esdTr, Int_t esdID)
32d38de2 683{
c61e50c3 684 // init prolongaion candidates finding for single seed
6cd80116 685 AliITSUTrackHyp* hyp = GetTrackHyp(esdID);
686 if (hyp) return hyp;
716ccba7 687 //
0ddbf657 688 fCountProlongationTrials++;
689 //
1d26da6d 690 fCurrMass = esdTr->GetMassForTracking();
c61e50c3 691 //
6cd80116 692 hyp = new AliITSUTrackHyp(fNLrActive);
693 hyp->SetESDTrack(esdTr);
694 hyp->SetUniqueID(esdID);
695 hyp->SetMass(fCurrMass);
9cdcba2c 696 hyp->SetTPCSeed( new AliExternalTrackParam(*esdTr) );
6cd80116 697 SetTrackHyp(hyp,esdID);
e1ef49ad 698 Bool_t lut = fUseMatLUT; // for propagation from TPC use TGeo material budget always
699 fUseMatLUT = kFALSE;
70cb7fe4 700 if (!TransportToLayer(hyp->GetTPCSeed(),fITS->GetNLayers(), fITS->GetLrIDActive(fNLrActive-1), fITS->GetRITSTPCRef())) hyp->SetSkip(); // propagate to outer R of ITS
e1ef49ad 701 fUseMatLUT = lut;
3e4e3c23 702 //
6cd80116 703 return hyp;
32d38de2 704 // TO DO
705}
706
707//_________________________________________________________________________
70cb7fe4 708Bool_t AliITSUTrackerGlo::TransportToLayer(AliITSUSeed* seed, Int_t lFrom, Int_t lTo, Double_t rLim)
32d38de2 709{
70cb7fe4 710 // transport seed from layerFrom to the entrance (along the direction of the track) of layerTo or to rLim (if>0), wathever is closer
32d38de2 711 //
8b16dbae 712 //
713 if (lTo==lFrom) AliFatal(Form("was called with lFrom=%d lTo=%d",lFrom,lTo));
714 //
32d38de2 715 int dir = lTo > lFrom ? 1:-1;
716 AliITSURecoLayer* lrFr = fITS->GetLayer(lFrom); // this can be 0 when extrapolation from TPC to ITS is requested
717 Bool_t checkFirst = kTRUE;
70cb7fe4 718 Bool_t limReached = kFALSE;
32d38de2 719 while(lFrom!=lTo) {
32d38de2 720 if (lrFr) {
70cb7fe4 721 if (!GoToExitFromLayer(seed,lrFr,dir,checkFirst)) return kFALSE; // go till the end of current layer
e2d2481c 722 checkFirst = kFALSE;
32d38de2 723 }
724 AliITSURecoLayer* lrTo = fITS->GetLayer( (lFrom+=dir) );
725 if (!lrTo) AliFatal(Form("Layer %d does not exist",lFrom));
726 //
727 // go the entrance of the layer, assuming no materials in between
8b16dbae 728 double xToGo = lrTo->GetR(-dir);
70cb7fe4 729 if (rLim>0) {
730 if (dir>0) {
731 if (rLim<xToGo) {xToGo = rLim; limReached = kTRUE;}
732 }
733 else {
734 if (rLim>xToGo) {xToGo = rLim; limReached = kTRUE;}
735 }
736 }
0ddbf657 737 // double xts = xToGo;
b2d3bc3a 738 if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) {
0ddbf657 739 //printf("FailHere1: %f %f %d\n",xts,xToGo,dir);
740 //seed->Print("etp");
70cb7fe4 741 //
b2d3bc3a 742 return kFALSE;
743 }
e1ef49ad 744#ifdef _ITSU_DEBUG_
b515ad7e 745 AliDebug(2,Form("go in dir=%d to R=%.4f(X:%.4f)",dir,lrTo->GetR(-dir), xToGo));
e1ef49ad 746#endif
b2d3bc3a 747 if (!PropagateSeed(seed,xToGo,fCurrMass,100, kFALSE )) {
0ddbf657 748 //printf("FailHere2: %f %f %d\n",xts,xToGo,dir);
749 //seed->Print("etp");
b2d3bc3a 750 return kFALSE;
751 }
32d38de2 752 lrFr = lrTo;
70cb7fe4 753 if (limReached) break;
32d38de2 754 }
755 return kTRUE;
756 //
757}
758
3e4e3c23 759//_________________________________________________________________________
70cb7fe4 760Bool_t AliITSUTrackerGlo::TransportToLayer(AliExternalTrackParam* seed, Int_t lFrom, Int_t lTo, Double_t rLim)
3e4e3c23 761{
70cb7fe4 762 // transport track from layerFrom to the entrance of layerTo or to rLim (if>0), wathever is closer
3e4e3c23 763 //
8b16dbae 764 if (lTo==lFrom) AliFatal(Form("was called with lFrom=%d lTo=%d",lFrom,lTo));
765 //
3e4e3c23 766 int dir = lTo > lFrom ? 1:-1;
767 AliITSURecoLayer* lrFr = fITS->GetLayer(lFrom); // this can be 0 when extrapolation from TPC to ITS is requested
768 Bool_t checkFirst = kTRUE;
70cb7fe4 769 Bool_t limReached = kFALSE;
3e4e3c23 770 while(lFrom!=lTo) {
3e4e3c23 771 if (lrFr) {
e2d2481c 772 if (!GoToExitFromLayer(seed,lrFr,dir,checkFirst)) return kFALSE; // go till the end of current layer
773 checkFirst = kFALSE;
3e4e3c23 774 }
775 AliITSURecoLayer* lrTo = fITS->GetLayer( (lFrom+=dir) );
776 if (!lrTo) AliFatal(Form("Layer %d does not exist",lFrom));
777 //
0ddbf657 778 // go the entrance of the layer, assuming no materials in between
779 double xToGo = lrTo->GetR(-dir);
70cb7fe4 780 if (rLim>0) {
781 if (dir>0) {
782 if (rLim<xToGo) {xToGo = rLim; limReached = kTRUE;}
783 }
784 else {
785 if (rLim>xToGo) {xToGo = rLim; limReached = kTRUE;}
786 }
787 }
0ddbf657 788 // double xts = xToGo;
789 if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) {
790 // printf("FailHere1: %f %f %d\n",xts,xToGo,dir);
791 // seed->Print("etp");
792 return kFALSE;
793 }
e1ef49ad 794#ifdef _ITSU_DEBUG_
b515ad7e 795 AliDebug(2,Form("go in dir=%d to R=%.4f(X:%.4f)",dir,lrTo->GetR(-dir), xToGo));
e1ef49ad 796#endif
b515ad7e 797 if (!PropagateSeed(seed,xToGo,fCurrMass,100, kFALSE )) {
798 //printf("FailHere2: %f %f %d\n",xts,xToGo,dir);
799 //seed->Print("etp");
800 return kFALSE;
801 }
802 lrFr = lrTo;
70cb7fe4 803 if (limReached) break;
b515ad7e 804 }
805 return kTRUE;
806 //
807}
808
809//_________________________________________________________________________
810Bool_t AliITSUTrackerGlo::TransportToLayerX(AliExternalTrackParam* seed, Int_t lFrom, Int_t lTo, Double_t xStop)
811{
812 // transport track from layerFrom to the entrance of layerTo but do not pass control parameter X
813 //
814 if (lTo==lFrom) AliFatal(Form("was called with lFrom=%d lTo=%d",lFrom,lTo));
815 //
816 int dir = lTo > lFrom ? 1:-1;
817 AliITSURecoLayer* lrFr = fITS->GetLayer(lFrom); // this can be 0 when extrapolation from TPC to ITS is requested
818 Bool_t checkFirst = kTRUE;
819 while(lFrom!=lTo) {
820 if (lrFr) {
821 if (!GoToExitFromLayer(seed,lrFr,dir,checkFirst)) return kFALSE; // go till the end of current layer
822 checkFirst = kFALSE;
823 }
824 AliITSURecoLayer* lrTo = fITS->GetLayer( (lFrom+=dir) );
825 if (!lrTo) AliFatal(Form("Layer %d does not exist",lFrom));
826 //
827 // go the entrance of the layer, assuming no materials in between
828 double xToGo = lrTo->GetR(-dir); // R of the entrance to layer
829 //
830 // double xts = xToGo;
831 if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) {
832 // printf("FailHere1: %f %f %d\n",xts,xToGo,dir);
833 // seed->Print("etp");
834 return kFALSE;
835 }
836 if ( (dir>0&&xToGo>xStop) || (dir<0&&xToGo<xStop) ) xToGo = xStop;
837 //
e1ef49ad 838#ifdef _ITSU_DEBUG_
b515ad7e 839 AliDebug(2,Form("go in dir=%d to R=%.4f(X:%.4f)",dir,lrTo->GetR(-dir), xToGo));
e1ef49ad 840#endif
0ddbf657 841 if (!PropagateSeed(seed,xToGo,fCurrMass,100, kFALSE )) {
842 //printf("FailHere2: %f %f %d\n",xts,xToGo,dir);
843 //seed->Print("etp");
844 return kFALSE;
845 }
3e4e3c23 846 lrFr = lrTo;
847 }
848 return kTRUE;
849 //
850}
851
e2d2481c 852//_________________________________________________________________________
853Bool_t AliITSUTrackerGlo::GoToExitFromLayer(AliITSUSeed* seed, AliITSURecoLayer* lr, Int_t dir, Bool_t check)
854{
855 // go to the exit from lr in direction dir, applying material corrections in steps specific for this layer
856 // If check is requested, do this only provided the track has not exited the layer already
857 double xToGo = lr->GetR(dir);
858 if (check) { // do we need to track till the surface of the current layer ?
859 double curR2 = seed->GetX()*seed->GetX() + seed->GetY()*seed->GetY(); // current radius
e1ef49ad 860#ifdef _ITSU_DEBUG_
15b02d69 861 AliDebug(3,Form(" dir:%d Cur: %e Tgt: %e",dir,Sqrt(curR2),xToGo));
e1ef49ad 862#endif
0ddbf657 863 if (dir>0) { if (curR2-xToGo*xToGo>-fgkToler) return kTRUE; } // on the surface or outside of the layer
864 else if (dir<0) { if (xToGo*xToGo-curR2>-fgkToler) return kTRUE; } // on the surface or outside of the layer
e2d2481c 865 }
866 if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
867 // go via layer to its boundary, applying material correction.
868 if (!PropagateSeed(seed,xToGo,fCurrMass, lr->GetMaxStep())) return kFALSE;
e7d83d38 869 //
e2d2481c 870 return kTRUE;
871 //
872}
873
874//_________________________________________________________________________
875Bool_t AliITSUTrackerGlo::GoToExitFromLayer(AliExternalTrackParam* seed, AliITSURecoLayer* lr, Int_t dir, Bool_t check)
876{
877 // go to the exit from lr in direction dir, applying material corrections in steps specific for this layer
878 // If check is requested, do this only provided the track has not exited the layer already
879 double xToGo = lr->GetR(dir);
880 if (check) { // do we need to track till the surface of the current layer ?
881 double curR2 = seed->GetX()*seed->GetX() + seed->GetY()*seed->GetY(); // current radius
0ddbf657 882 if (dir>0) { if (curR2-xToGo*xToGo>-fgkToler) return kTRUE; } // on the surface or outside of the layer
883 else if (dir<0) { if (xToGo*xToGo-curR2>-fgkToler) return kTRUE; } // on the surface or outside of the layer
e2d2481c 884 }
e1ef49ad 885#ifdef _ITSU_DEBUG_
b515ad7e 886 AliDebug(2,Form(" dir=%d : from R=%.4f -> R=%.4f",dir,Sqrt(seed->GetX()*seed->GetX() + seed->GetY()*seed->GetY()), xToGo));
e1ef49ad 887#endif
e2d2481c 888 if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
889 // go via layer to its boundary, applying material correction.
890 if (!PropagateSeed(seed,xToGo,fCurrMass, lr->GetMaxStep())) return kFALSE;
e7d83d38 891 //
e2d2481c 892 return kTRUE;
893 //
894}
895
896
897//_________________________________________________________________________
898Bool_t AliITSUTrackerGlo::GoToEntranceToLayer(AliITSUSeed* seed, AliITSURecoLayer* lr, Int_t dir, Bool_t check)
899{
900 // go to the entrance of lr in direction dir, w/o applying material corrections.
901 // If check is requested, do this only provided the track did not reach the layer already
902 double xToGo = lr->GetR(-dir);
903 if (check) { // do we need to track till the surface of the current layer ?
904 double curR2 = seed->GetX()*seed->GetX() + seed->GetY()*seed->GetY(); // current radius
0ddbf657 905 if (dir>0) { if (curR2-xToGo*xToGo>-fgkToler) return kTRUE; } // already passed
906 else if (dir<0) { if (xToGo*xToGo-curR2>-fgkToler) return kTRUE; } // already passed
e2d2481c 907 }
908 if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
909 // go via layer to its boundary, applying material correction.
910 if (!PropagateSeed(seed,xToGo,fCurrMass, 100, kFALSE)) return kFALSE;
911 return kTRUE;
912 //
913}
914
915//_________________________________________________________________________
916Bool_t AliITSUTrackerGlo::GoToEntranceToLayer(AliExternalTrackParam* seed, AliITSURecoLayer* lr, Int_t dir, Bool_t check)
917{
918 // go to the entrance of lr in direction dir, w/o applying material corrections.
919 // If check is requested, do this only provided the track did not reach the layer already
920 double xToGo = lr->GetR(-dir);
921 if (check) { // do we need to track till the surface of the current layer ?
922 double curR2 = seed->GetX()*seed->GetX() + seed->GetY()*seed->GetY(); // current radius
0ddbf657 923 if (dir>0) { if (curR2-xToGo*xToGo>-fgkToler) return kTRUE; } // already passed
924 else if (dir<0) { if (xToGo*xToGo-curR2>-fgkToler) return kTRUE; } // already passed
e2d2481c 925 }
926 if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
927 // go via layer to its boundary, applying material correction.
928 if (!PropagateSeed(seed,xToGo,fCurrMass, 100, kFALSE)) return kFALSE;
929 return kTRUE;
930 //
931}
932
32d38de2 933//_________________________________________________________________________
c61e50c3 934Bool_t AliITSUTrackerGlo::GetRoadWidth(AliITSUSeed* seed, int ilrA)
32d38de2 935{
936 // calculate road width in terms of phi and z for the track which MUST be on the external radius of the layer
937 // as well as some aux info
938 double bz = GetBz();
939 AliITSURecoLayer* lrA = fITS->GetLayerActive(ilrA);
c61e50c3 940 seed->GetXYZ(&fTrImpData[kTrXIn]); // lab position at the entrance from above
53870004 941 static AliExternalTrackParam sc; // seed copy for manipulations
943e1898 942 sc = *seed;
32d38de2 943 //
c61e50c3 944 fTrImpData[kTrPhiIn] = ATan2(fTrImpData[kTrYIn],fTrImpData[kTrXIn]);
943e1898 945 if (!sc.Rotate(fTrImpData[kTrPhiIn])) return kFALSE; // go to the frame of the entry point into the layer
68b7631d 946 BringTo02Pi(fTrImpData[kTrPhiIn]);
c61e50c3 947 double dr = lrA->GetDR(); // approximate X dist at the inner radius
943e1898 948 if (!sc.GetXYZAt(sc.GetX()-dr, bz, fTrImpData + kTrXOut)) {
32d38de2 949 // special case: track does not reach inner radius, might be tangential
943e1898 950 double r = sc.GetD(0,0,bz);
32d38de2 951 double x;
bdb92766 952 if (!sc.GetXatLabR(r,x,bz,-1)) return kFALSE;
943e1898 953 dr = Abs(sc.GetX() - x);
bdb92766 954 if (!sc.GetXYZAt(x, bz, fTrImpData + kTrXOut)) return kFALSE;
32d38de2 955 }
956 //
c61e50c3 957 fTrImpData[kTrPhiOut] = ATan2(fTrImpData[kTrYOut],fTrImpData[kTrXOut]);
68b7631d 958 BringTo02Pi(fTrImpData[kTrPhiOut]);
943e1898 959 double sgy = sc.GetSigmaY2() + dr*dr*sc.GetSigmaSnp2() + AliITSUReconstructor::GetRecoParam()->GetSigmaY2(ilrA);
960 double sgz = sc.GetSigmaZ2() + dr*dr*sc.GetSigmaTgl2() + AliITSUReconstructor::GetRecoParam()->GetSigmaZ2(ilrA);
ee54014a 961 sgy = Sqrt(sgy)*fCurrTrackCond->GetNSigmaRoadY(ilrA);
962 sgz = Sqrt(sgz)*fCurrTrackCond->GetNSigmaRoadZ(ilrA);
68b7631d 963 double phi0 = MeanPhiSmall(fTrImpData[kTrPhiOut],fTrImpData[kTrPhiIn]);
964 double dphi0 = DeltaPhiSmall(fTrImpData[kTrPhiOut],fTrImpData[kTrPhiIn]);
965 //
5785a9d8 966 fTrImpData[kTrPhi0] = phi0;
173b3073 967 fTrImpData[kTrZ0] = 0.5*(fTrImpData[kTrZOut]+fTrImpData[kTrZIn]);
31e2f5d4 968 dphi0 += sgy/lrA->GetR();
969 fTrImpData[kTrDPhi] = dphi0<PiOver2() ? dphi0 : PiOver2();
173b3073 970 fTrImpData[kTrDZ] = 0.5*Abs(fTrImpData[kTrZOut]-fTrImpData[kTrZIn]) + sgz;
32d38de2 971 //
972 return kTRUE;
973}
974
b8b59e05 975//________________________________________
976void AliITSUTrackerGlo::ResetSeedsPool()
977{
978 // mark all seeds in the pool as unused
0ddbf657 979 AliDebug(-1,Form("CurrentSize: %d, BookedUpTo: %d, free: %d",fSeedsPool.GetSize(),fSeedsPool.GetEntriesFast(),fNFreeSeeds));
b8b59e05 980 fNFreeSeeds = 0;
981 fSeedsPool.Clear(); // seeds don't allocate memory
982}
983
984
985//________________________________________
986void AliITSUTrackerGlo::MarkSeedFree(AliITSUSeed *sd)
32d38de2 987{
b8b59e05 988 // account that this seed is "deleted"
989 int id = sd->GetPoolID();
990 if (id<0) {
991 AliError(Form("Freeing of seed %p NOT from the pool is requested",sd));
992 return;
993 }
994 // AliInfo(Form("%d %p",id, seed));
995 fSeedsPool.RemoveAt(id);
996 if (fFreeSeedsID.GetSize()<=fNFreeSeeds) fFreeSeedsID.Set( 2*fNFreeSeeds + 100 );
997 fFreeSeedsID.GetArray()[fNFreeSeeds++] = id;
32d38de2 998}
f8832015 999
1000//_________________________________________________________________________
1001Int_t AliITSUTrackerGlo::CheckCluster(AliITSUSeed* track, Int_t lr, Int_t clID)
1002{
1003 // Check if the cluster (in tracking frame!) is matching to track.
1004 // The track must be already propagated to sensor tracking frame.
1005 // Returns: kStopSearchOnSensor if the search on given sensor should be stopped,
1006 // kClusterMatching if the cluster is matching
1007 // kClusterMatching otherwise
1008 //
f8832015 1009 const double kTolerX = 5e-4;
c03e4f8a 1010 // AliCluster *cl = fITS->GetLayerActive(lr)->GetCluster(clID);
1011 AliCluster *cl = fCurrLayer->GetCluster(clID);
f8832015 1012 //
c8d1f258 1013 Bool_t goodCl = kFALSE;
0ddbf657 1014 int currLabel = Abs(fCurrESDtrMClb);
716ccba7 1015 //
8b16dbae 1016 if (cl->GetLabel(0)>=0) {
1017 for (int i=0;i<3;i++) if (cl->GetLabel(i)>=0 && cl->GetLabel(i)==currLabel) {goodCl = kTRUE; break;}
1018 }
c8d1f258 1019 //
f8832015 1020 if (TMath::Abs(cl->GetX())>kTolerX) { // if due to the misalingment X is large, propagate track only
c8d1f258 1021 if (!track->PropagateParamOnlyTo(track->GetX()+cl->GetX(),GetBz())) {
08419930 1022 //
e1ef49ad 1023#ifdef _ITSU_DEBUG_
08419930 1024 if (AliDebugLevelClass()>2 && goodCl && !track->ContainsFake()) {
53870004 1025 AliDebug(2,Form("Lost good cl on L:%d failed propagation. |ESDtrack#%d (MClb:%d)",lr,fCurrESDtrack->GetID(),fCurrESDtrMClb));
1026 track->Print("etp");
15b02d69 1027 cl->Print();
1028 }
e1ef49ad 1029#endif
c8d1f258 1030 return kStopSearchOnSensor; // propagation failed, seedT is intact
1031 }
f8832015 1032 }
3dd9c283 1033 double dy = cl->GetY()-track->GetY();
c8d1f258 1034 double dz = cl->GetZ()-track->GetZ();
0091e9f0 1035 //
f8832015 1036 double dy2 = dy*dy;
1037 double tol2 = (track->GetSigmaY2() + AliITSUReconstructor::GetRecoParam()->GetSigmaY2(lr))*
ee54014a 1038 fCurrTrackCond->GetNSigmaRoadY(lr)*fCurrTrackCond->GetNSigmaRoadY(lr); // RS TOOPTIMIZE
f8832015 1039 if (dy2>tol2) { // the clusters are sorted in Z(col) then in Y(row).
08419930 1040 //
e1ef49ad 1041#ifdef _ITSU_DEBUG_
08419930 1042 if (AliDebugLevelClass()>2 && goodCl && !track->ContainsFake()) {
53870004 1043 AliDebug(2,Form("Lost good cl: dy2=%e > tol2=%e |ESDtrack#%d (MClb:%d)",dy2,tol2,fCurrESDtrack->GetID(),fCurrESDtrMClb));
1044 track->Print("etp");
15b02d69 1045 cl->Print();
53870004 1046 PrintSeedClusters(track);
15b02d69 1047 }
e1ef49ad 1048#endif
cb50e082 1049 // clusters are sorted in increasing Y and the loop goes from last (highers Y) to 1st.
1050 // if the track has too large y for this cluster (dy<0) it will be so for all other clusters of the sensor
1051 if (dy<0) return kStopSearchOnSensor; // No chance that other cluster of this sensor will match (all Y's will be even larger)
f8832015 1052 else return kClusterNotMatching; // Other clusters may match
1053 }
c8d1f258 1054 double dz2 = dz*dz;
f8832015 1055 tol2 = (track->GetSigmaZ2() + AliITSUReconstructor::GetRecoParam()->GetSigmaZ2(lr))*
ee54014a 1056 fCurrTrackCond->GetNSigmaRoadZ(lr)*fCurrTrackCond->GetNSigmaRoadZ(lr); // RS TOOPTIMIZE
c8d1f258 1057 if (dz2>tol2) {
e1ef49ad 1058#ifdef _ITSU_DEBUG_
08419930 1059 if (AliDebugLevelClass()>2 && goodCl && !track->ContainsFake()) {
53870004 1060 AliDebug(2,Form("Lost good cl on L:%d : dz2=%e > tol2=%e |ESDtrack#%d (MClb:%d)",lr,dz2,tol2,fCurrESDtrack->GetID(),fCurrESDtrMClb));
1061 track->Print("etp");
15b02d69 1062 cl->Print();
53870004 1063 PrintSeedClusters(track);
15b02d69 1064 }
e1ef49ad 1065#endif
c8d1f258 1066 return kClusterNotMatching; // Other clusters may match
1067 }
f8832015 1068 //
1069 // check chi2
1070 Double_t p[2]={cl->GetY(), cl->GetZ()};
1071 Double_t cov[3]={cl->GetSigmaY2(), cl->GetSigmaYZ(), cl->GetSigmaZ2()};
1072 double chi2 = track->GetPredictedChi2(p,cov);
76390254 1073 //
ee54014a 1074 if (chi2>fCurrTrackCond->GetMaxTr2ClChi2(lr)) {
e1ef49ad 1075#ifdef _ITSU_DEBUG_
08419930 1076 if (AliDebugLevelClass()>2 && goodCl && !track->ContainsFake()) {
53870004 1077 AliDebug(2,Form("Lost good cl on L:%d : Chi2=%e > Chi2Max=%e |dy: %+.3e dz: %+.3e |ESDtrack#%d (MClb:%d)\n",
ee54014a 1078 lr,chi2,fCurrTrackCond->GetMaxTr2ClChi2(lr),dy,dz,fCurrESDtrack->GetID(),fCurrESDtrMClb));
3dd9c283 1079 track->Print("etp");
1080 cl->Print("");
53870004 1081 PrintSeedClusters(track);
3dd9c283 1082 }
e1ef49ad 1083#endif
c8d1f258 1084 return kClusterNotMatching;
1085 }
f8832015 1086 //
1087 track = NewSeedFromPool(track); // input track will be reused, use its clone for updates
943e1898 1088 if (!track->Update()) {
e1ef49ad 1089#ifdef _ITSU_DEBUG_
08419930 1090 if (AliDebugLevelClass()>2 && goodCl && !track->ContainsFake()) {
53870004 1091 AliDebug(2,Form("Lost good cluster on L:%d : failed to update",lr));
38997c3c 1092 track->Print("etp");
1093 cl->Print("");
53870004 1094 PrintSeedClusters(track);
38997c3c 1095 }
e1ef49ad 1096#endif
b8b59e05 1097 MarkSeedFree(track);
c8d1f258 1098 return kClusterNotMatching;
1099 }
f8832015 1100 track->SetChi2Cl(chi2);
1101 track->SetLrClusterID(lr,clID);
08419930 1102 //
1103 if (track->GetChi2GloNrm()>fCurrTrackCond->GetMaxChi2GloNrm(lr)) {
e1ef49ad 1104#ifdef _ITSU_DEBUG_
08419930 1105 if (AliDebugLevelClass()>2 && goodCl && !track->ContainsFake()) {
1106 AliDebug(2,Form("Lost good cl on L:%d : Chi2Glo=%e > Chi2Max=%e |dy: %+.3e dz: %+.3e |ESDtrack#%d (MClb:%d)\n",
1107 lr,track->GetChi2GloNrm(),fCurrTrackCond->GetMaxChi2GloNrm(lr),dy,dz,fCurrESDtrack->GetID(),fCurrESDtrMClb));
1108 track->Print("etp");
1109 cl->Print("");
1110 PrintSeedClusters(track);
1111 }
e1ef49ad 1112#endif
08419930 1113 MarkSeedFree(track);
1114 return kClusterNotMatching;
1115 }
c03e4f8a 1116 // cl->IncreaseClusterUsage(); // do this only for winners
f8832015 1117 //
716ccba7 1118 track->SetFake(!goodCl);
1119 //
e1ef49ad 1120#ifdef _ITSU_DEBUG_
0ddbf657 1121 AliDebug(2,Form("AddCl(%d) Cl%d lr:%d: dY:%+8.4f dZ:%+8.4f (MC: %5d %5d %5d) |Chi2=%f(%c)| ",
716ccba7 1122 goodCl,clID,lr,dy,dz2,cl->GetLabel(0),cl->GetLabel(1),cl->GetLabel(2), chi2, track->IsFake() ? '-':'+'));
e1ef49ad 1123#endif
c8d1f258 1124 //
38997c3c 1125 AddSeedBranch(track);
f8832015 1126 //
1127 return kClusterMatching;
1128}
c8d1f258 1129
1130//_________________________________________________________________________
1131Bool_t AliITSUTrackerGlo::NeedToKill(AliITSUSeed *seed, Int_t flag)
1132{
1133 // check if the seed should not be discarded
1134 const UShort_t kMask = 0xffff;
1135 if (flag==kMissingCluster) {
1136 int lastChecked = seed->GetLayerID();
1137 UShort_t patt = seed->GetHitsPattern();
097e4cdf 1138 if (lastChecked) patt |= (~(kMask<<lastChecked))&fCurrTrackCond->GetAllowedLayers(); // not all layers were checked, complete unchecked ones by potential hits
42c3d4bd 1139 Bool_t seedOK = fCurrTrackCond->CheckPattern(patt);
c8d1f258 1140 return !seedOK;
1141 }
1142 return kTRUE;
1143}
44785f3e 1144
1145//______________________________________________________________________________
1146Bool_t AliITSUTrackerGlo::PropagateSeed(AliITSUSeed *seed, Double_t xToGo, Double_t mass, Double_t maxStep, Bool_t matCorr)
1147{
1148 // propagate seed to given x applying material correction if requested
1149 const Double_t kEpsilon = 1e-5;
1150 Double_t xpos = seed->GetX();
1151 Int_t dir = (xpos<xToGo) ? 1:-1;
e1ef49ad 1152 Double_t xyz0[3],xyz1[3];
44785f3e 1153 //
8b16dbae 1154 Bool_t updTime = dir>0 && seed->IsStartedTimeIntegral();
f497470c 1155 if (matCorr || updTime) seed->GetXYZ(xyz1); //starting global position
44785f3e 1156 while ( (xToGo-xpos)*dir > kEpsilon){
1157 Double_t step = dir*TMath::Min(TMath::Abs(xToGo-xpos), maxStep);
1158 Double_t x = xpos+step;
1159 Double_t bz=GetBz(); // getting the local Bz
0ddbf657 1160 if (!seed->PropagateToX(x,bz)) return kFALSE;
8b16dbae 1161 double ds = 0;
1162 if (matCorr || updTime) {
44785f3e 1163 xyz0[0]=xyz1[0]; // global pos at the beginning of step
1164 xyz0[1]=xyz1[1];
1165 xyz0[2]=xyz1[2];
1166 seed->GetXYZ(xyz1); // // global pos at the end of step
8b16dbae 1167 if (matCorr) {
e1ef49ad 1168 Double_t xrho,xx0;
1169 ds = GetMaterialBudget(xyz0,xyz1,xx0,xrho);
8b16dbae 1170 if (dir>0) xrho = -xrho; // outward should be negative
b2d3bc3a 1171 if (!seed->ApplyMaterialCorrection(xx0,xrho,mass,kFALSE)) {AliInfo("Fail2"); seed->Print("etp"); return kFALSE;}
8b16dbae 1172 }
1173 else { // matCorr is not requested but time integral is
1174 double d0 = xyz1[0]-xyz0[0];
1175 double d1 = xyz1[1]-xyz0[1];
1176 double d2 = xyz1[2]-xyz0[2];
1177 ds = TMath::Sqrt(d0*d0+d1*d1+d2*d2);
1178 }
44785f3e 1179 }
8b16dbae 1180 if (updTime) seed->AddTimeStep(ds);
44785f3e 1181 xpos = seed->GetX();
1182 }
1183 return kTRUE;
1184}
716ccba7 1185
3e4e3c23 1186//______________________________________________________________________________
1187Bool_t AliITSUTrackerGlo::PropagateSeed(AliExternalTrackParam *seed, Double_t xToGo, Double_t mass, Double_t maxStep, Bool_t matCorr)
1188{
1189 // propagate seed to given x applying material correction if requested
1190 const Double_t kEpsilon = 1e-5;
1191 Double_t xpos = seed->GetX();
1192 Int_t dir = (xpos<xToGo) ? 1:-1;
e1ef49ad 1193 Double_t xyz0[3],xyz1[3];
3e4e3c23 1194 //
e1ef49ad 1195#ifdef _ITSU_DEBUG_
b515ad7e 1196 if (AliDebugLevelClass()>1) {
1197 AliDebug(2,Form("Before Propagate to X=%f with M=%.3f MaxStep=%.4f MatCorr=%d",xToGo,mass,maxStep,matCorr));
1198 seed->AliExternalTrackParam::Print();
1199 }
e1ef49ad 1200#endif
3e4e3c23 1201 Bool_t updTime = dir>0 && seed->IsStartedTimeIntegral();
1202 if (matCorr || updTime) seed->GetXYZ(xyz1); //starting global position
1203 while ( (xToGo-xpos)*dir > kEpsilon){
1204 Double_t step = dir*TMath::Min(TMath::Abs(xToGo-xpos), maxStep);
1205 Double_t x = xpos+step;
1206 Double_t bz=GetBz(); // getting the local Bz
1207 if (!seed->PropagateTo(x,bz)) return kFALSE;
1208 double ds = 0;
1209 if (matCorr || updTime) {
1210 xyz0[0]=xyz1[0]; // global pos at the beginning of step
1211 xyz0[1]=xyz1[1];
1212 xyz0[2]=xyz1[2];
1213 seed->GetXYZ(xyz1); // // global pos at the end of step
1214 //
1215 if (matCorr) {
e1ef49ad 1216 Double_t xrho,xx0;
1217 ds = GetMaterialBudget(xyz0,xyz1,xx0,xrho);
3e4e3c23 1218 if (dir>0) xrho = -xrho; // outward should be negative
1219 if (!seed->CorrectForMeanMaterial(xx0,xrho,mass)) return kFALSE;
3e4e3c23 1220 }
1221 else { // matCorr is not requested but time integral is
1222 double d0 = xyz1[0]-xyz0[0];
1223 double d1 = xyz1[1]-xyz0[1];
1224 double d2 = xyz1[2]-xyz0[2];
1225 ds = TMath::Sqrt(d0*d0+d1*d1+d2*d2);
1226 }
1227 }
1228 if (updTime) seed->AddTimeStep(ds);
1229 //
1230 xpos = seed->GetX();
1231 }
e1ef49ad 1232#ifdef _ITSU_DEBUG_
b515ad7e 1233 if (AliDebugLevelClass()>1) {
1234 AliDebug(2,Form("After Propagate to X=%f",xToGo));
1235 seed->AliExternalTrackParam::Print();
1236 }
e1ef49ad 1237#endif
3e4e3c23 1238 return kTRUE;
1239}
1240
9cdcba2c 1241//______________________________________________________________________________
1242Bool_t AliITSUTrackerGlo::FinalizeHypothesis(AliITSUTrackHyp* hyp)
1243{
1244 // finalize single TPC track hypothesis
1245 if (!hyp || hyp->GetESDTrack()->IsOn(AliESDtrack::kITSin)) return kFALSE;
1246 AliITSUSeed* winner = 0;
1247 fCurrHyp = hyp;
70cb7fe4 1248 fCurrMass = hyp->GetMass();
d99994b2 1249 if (!(winner=hyp->DefineWinner(fCurrTrackCond->GetActiveLrInner()))) return kFALSE;
9cdcba2c 1250 CookMCLabel(hyp);
1251 //
08419930 1252#ifdef _ITSU_TUNING_MODE_ // fill tuning histos
1253 TObjArray* dest = hyp->GetLabel()>=0 ? fCHistoArrCorr : fCHistoArrFake;
1254 if (dest) {
1255 AliITSUSeed* sd = winner;
1256 double htrPt = hyp->Pt();
1257 do {
1258 int clID,lrID;
1259 if ( (clID=sd->GetLrCluster(lrID))<0 ) continue;
781c79ac 1260 ((TH2*)dest->At( GetHistoID(lrID,kHChi2Nrm ,fCurrPassID) ))->Fill(htrPt,sd->GetChi2GloNrm());
1261 ((TH2*)dest->At( GetHistoID(lrID,kHBestInBranch,fCurrPassID) ))->Fill(htrPt,sd->GetOrdBranch());
1262 ((TH2*)dest->At( GetHistoID(lrID,kHBestInCand ,fCurrPassID) ))->Fill(htrPt,sd->GetOrdCand());
08419930 1263 //
1264 if (dest==fCHistoArrFake && !sd->IsFake()) continue; // for the fake seeds fill only fake clusters part
1265 //
781c79ac 1266 ((TH2*)dest->At( GetHistoID(lrID,kHResY ,fCurrPassID) ))->Fill(htrPt,sd->GetResidY());
1267 ((TH2*)dest->At( GetHistoID(lrID,kHResZ ,fCurrPassID) ))->Fill(htrPt,sd->GetResidZ());
1268 ((TH2*)dest->At( GetHistoID(lrID,kHResYP ,fCurrPassID) ))->Fill(htrPt,sd->GetPullY());
1269 ((TH2*)dest->At( GetHistoID(lrID,kHResZP ,fCurrPassID) ))->Fill(htrPt,sd->GetPullZ());
1270 ((TH2*)dest->At( GetHistoID(lrID,kHChi2Cl ,fCurrPassID) ))->Fill(htrPt,sd->GetChi2Cl());
08419930 1271 //
1272 } while((sd=(AliITSUSeed*)sd->GetParent()));
1273 //
781c79ac 1274 ((TH2*)dest->At( GetHistoID(-1,kHChiMatch,fCurrPassID) ))->Fill(htrPt,winner->GetChi2ITSTPC());
1275 ((TH2*)dest->At( GetHistoID(-1,kHChiITSSA,fCurrPassID) ))->Fill(htrPt,winner->GetChi2ITSSA());
9cdcba2c 1276 }
9cdcba2c 1277 //
08419930 1278#endif
1279 //
1280 CheckClusterSharingConflicts(hyp);
1281 return hyp->GetWinner(); // winner might change of disappear after resolving conflicts
1282 //
9cdcba2c 1283}
1284
3e4e3c23 1285//______________________________________________________________________________
1286void AliITSUTrackerGlo::FinalizeHypotheses()
1287{
1288 // select winner for each hypothesis, remove cl. sharing conflicts
9cdcba2c 1289 AliCodeTimerAuto("",0);
3e4e3c23 1290 //
d99994b2 1291 int* esdTrackIndex = fESDIndex.GetArray();
1292 for (int ih=0;ih<fNTracksESD;ih++) FinalizeHypothesis(GetTrackHyp(esdTrackIndex[ih])); // finlize in decreasing pt order
c03e4f8a 1293 //
d99994b2 1294 for (int ih=fNTracksESD;ih--;) UpdateESDTrack(GetTrackHyp(esdTrackIndex[ih]),AliESDtrack::kITSin);
08419930 1295 //
1296}
bdb92766 1297
08419930 1298//______________________________________________________________________________
1299void AliITSUTrackerGlo::CheckClusterSharingConflicts(AliITSUTrackHyp* hyp)
1300{
1301 // remove eventual cluster sharing conflicts
1302 AliITSUSeed* winner = 0;
1303 if (!(winner=hyp->GetWinner())) return;
1304 UShort_t idH = (UShort_t)hyp->GetUniqueID()+1;
1305 int lrID,clID;
d99994b2 1306 int inLr = fCurrTrackCond->GetActiveLrInner();
08419930 1307 AliITSUSeed* winSD=winner;
1308 do {
1309 if ( (clID=winSD->GetLrCluster(lrID))<0 ) continue;
1310 AliITSUClusterPix* cl = (AliITSUClusterPix*)fITS->GetLayerActive(lrID)->GetCluster(clID);
1311 Int_t refID = cl->GetRecoInfo(); // was it already referred by some track (with refID-1 if >0)
1312 if (!refID) {cl->SetRecoInfo(idH); continue;} // if not, refer from track to cluster
1313 if (refID==idH) continue; // ignore reference to itself
1314 // the cluster is already used by other track, need to resolve the conflict
1315 AliITSUTrackHyp* hypC = GetTrackHyp(refID-1); // competitor
1316 // printf("Ref to %d (%p) from %d Cl %d %d\n",refID-1,hypC, idH-1,clID,lrID);
1317 AliITSUSeed* winnerC = hypC->GetWinner();
1318 if (!winnerC) {
1319 printf("Missing winner, ERROR\n");
1320 continue;
1321 }
6091ad4e 1322 //
e1ef49ad 1323#ifdef _ITSU_DEBUG_
08419930 1324 if (AliDebugLevelClass()>1) { //------------------------ DEBUG -------------------
1325 AliInfo(Form("Shared cl#%4d on Lr:%d: Hyp%3d/Q:%6.3f vs Hyp%3d/Q:%6.3f",
1326 clID,lrID,idH-1,winner->GetQualityVar(),refID-1,winnerC->GetQualityVar()));
1327 // dump winner of hyp
1328 printf("#%4d Pt:%.3f Chi:%6.3f/%6.3f/%6.3f Ncl:%d MCits%+5d MCtpc:%+5d |",
1329 idH-1,winner->Pt(),winner->GetChi2GloNrm(),winner->GetChi2ITSTPC(),winner->GetChi2ITSSA(),
1330 winner->GetNLayersHit(),hyp->GetITSLabel(),hyp->GetESDTrack()->GetTPCLabel());
1331 int prevL=-1;
1332 AliITSUSeed* sd = winner;
1333 do {
1334 int lrs;
1335 int clIDt = sd->GetLrCluster(lrs);
1336 if (clIDt<0) continue;
1337 while( lrs>++prevL ) printf("%4s ","----");
1338 printf("%4d (%5.1f)",clIDt,sd->GetChi2Cl());
1339 } while ((sd=(AliITSUSeed*)sd->GetParent()));
1340 printf("|\n");
1341 // dump winner of hypC
1342 printf("#%4d Pt:%.3f Chi:%6.3f/%6.3f/%6.3f Ncl:%d MCits%+5d MCtpc:%+5d |",
1343 refID-1,winnerC->Pt(),winnerC->GetChi2GloNrm(),winnerC->GetChi2ITSTPC(),winnerC->GetChi2ITSSA(),
1344 winnerC->GetNLayersHit(),hypC->GetITSLabel(),hypC->GetESDTrack()->GetTPCLabel());
1345 prevL=-1;
1346 sd = winnerC;
6091ad4e 1347 do {
08419930 1348 int lrs;
1349 int clIDt = sd->GetLrCluster(lrs);
1350 if (clIDt<0) continue;
1351 while( lrs>++prevL ) printf("%4s ","----");
1352 printf("%4d (%5.1f)",clIDt,sd->GetChi2Cl());
1353 } while ((sd=(AliITSUSeed*)sd->GetParent()));
1354 printf("|\n");
1355 } //-------------------------------------------------- DEBUG -------------------
e1ef49ad 1356#endif
08419930 1357 //
1358 if (winner->GetQualityVar()<winnerC->GetQualityVar()) { // current tracks is better then competitor track
1359 FlagSeedClusters(winnerC,kFALSE,refID); // unflag cluster usage by loser
1360 cl->SetRecoInfo(idH);
1361 //
1362 //if ( hypC->GetLabel()>=0) AliInfo(Form("KILLING CORRECT: %d",refID-1));
1363 winnerC->Kill();
d99994b2 1364 if (hypC->DefineWinner(inLr)) { // find new winner instead of suppressed candidate
08419930 1365 CookMCLabel(hypC);
1366 CheckClusterSharingConflicts(hypC); // and check its sharing conflicts again
1367 if (winner->IsKilled()) break; // the current winner might have been killed during check of new winner of hypC!
1368 }
6091ad4e 1369 }
08419930 1370 else { // competitor hypC is better than the hyp
1371 FlagSeedClusters(winner,kFALSE,idH); // unflag cluster usage by loser
1372 winner->Kill();
1373 //if ( hyp->GetLabel()>=0) AliInfo(Form("KILLING CORRECT: %d",idH-1));
d99994b2 1374 if (hyp->DefineWinner(inLr)) {
08419930 1375 CookMCLabel(hyp);
1376 CheckClusterSharingConflicts(hyp);
9cdcba2c 1377 }
08419930 1378 break;
9cdcba2c 1379 }
08419930 1380 //
1381 } while ((winSD=(AliITSUSeed*)winSD->GetParent()));
1382 //
1383 return;
3e4e3c23 1384}
c51c3124 1385
1386//______________________________________________________________________________
68a0f687 1387void AliITSUTrackerGlo::UpdateESDTrack(AliITSUTrackHyp* hyp,Int_t flag)
c51c3124 1388{
1389 // update ESD track with current best hypothesis
08419930 1390 if (!hyp) return;
c51c3124 1391 AliESDtrack* esdTr = hyp->GetESDTrack();
d99994b2 1392 if (!esdTr || esdTr->IsOn(flag)) return;
c51c3124 1393 AliITSUSeed* win = hyp->GetWinner();
08419930 1394 if (!win || win->IsKilled()) return;
70cb7fe4 1395 double chiSav;
1396 //
68a0f687 1397 switch (flag) {
1398 case AliESDtrack::kITSin:
1399 esdTr->UpdateTrackParams(hyp,flag); // update kinematics
08419930 1400 fCountITSin++;
dc443a34 1401 // set fakes cluster info
1402 {
1403 UShort_t clfake = 0;
c7e63626 1404 Int_t clSplit = 0;
dc443a34 1405 AliITSUSeed* sd = win;
c7e63626 1406 int ip=0;
dc443a34 1407 do {
c7e63626 1408 int lr, clID = sd->GetLrCluster(lr);
1409 if (sd->IsFake()) clfake |= 0x1<<lr;
1410 if (clID>=0) {
1411 esdTr->SetITSModuleIndex(ip++, sd->GetLrClusterID());
1412 AliITSUClusterPix *cl = (AliITSUClusterPix*)fITS->GetLayerActive(lr)->GetCluster(clID);
1413#ifdef _ITSU_TUNING_MODE_
1414 if (cl->IsSplit()) clSplit |= 0x1<<lr;
1415#endif
1416 }
dc443a34 1417 } while ((sd=(AliITSUSeed*)sd->GetParent()));
0e40ce66 1418 //
1419 // RS: Temporary set special flag for tracks from the afterburner
1420 if (fCurrPassID>0) clfake |= 0x1<<7;
1421 //
dc443a34 1422 esdTr->SetITSSharedMap(clfake);
c7e63626 1423 esdTr->SetITSModuleIndex(10,clSplit);
dc443a34 1424 }
097e4cdf 1425 // TEMPORARY: store iteration id
1426 esdTr->SetITSModuleIndex(11,fCurrPassID);
68a0f687 1427 break;
1428 //
1429 case AliESDtrack::kITSout:
70cb7fe4 1430 // here the stored chi2 will correspond to backward ITS-SA fit
68a0f687 1431 esdTr->UpdateTrackParams(hyp,flag); // update kinematics
08419930 1432 fCountITSout++;
68a0f687 1433 // TODO: avoid setting friend
1434 break;
1435 //
1436 case AliESDtrack::kITSrefit:
70cb7fe4 1437
1438 // at the moment fill as a chi2 the TPC-ITS matching chi2
1439 chiSav = hyp->GetChi2();
1440 hyp->SetChi2(win->GetChi2ITSTPC());
68a0f687 1441 esdTr->UpdateTrackParams(hyp,flag); // update kinematics
70cb7fe4 1442 hyp->SetChi2(chiSav);
08419930 1443 fCountITSrefit++;
68a0f687 1444 // TODO: avoid setting cluster info
1445 break;
1446 default:
1447 AliFatal(Form("Unknown flag %d",flag));
1448 }
c51c3124 1449 //
bdb92766 1450 esdTr->SetITSLabel(hyp->GetITSLabel());
852af72e 1451 // transfer chip indices
c51c3124 1452 // TODO
1453}
1454
1455//______________________________________________________________________________
7ce7f445 1456Double_t AliITSUTrackerGlo::RefitTrack(AliITSUTrackHyp* trc, Double_t rDest, Int_t &nclFit, Int_t stopCond)
c51c3124 1457{
9cdcba2c 1458 // refit track till radius rDest. The cluster,mass info is taken from fCurrHyp (and its winner seed)
1459 // if stopCond<0 : propagate till last cluster then stop
1460 // if stopCond==0: propagate till last cluster then try to go till limiting rDest, don't mind if fail
1461 // if stopCond>0 : rDest must be reached
8b16dbae 1462 //
c51c3124 1463 double rCurr = Sqrt(trc->GetX()*trc->GetX() + trc->GetY()*trc->GetY());
1464 int dir,lrStart,lrStop;
1465 //
8b16dbae 1466 dir = rCurr<rDest ? 1 : -1;
1467 lrStart = fITS->FindFirstLayerID(rCurr,dir);
70cb7fe4 1468 lrStop = fITS->FindLastLayerID(rDest,dir); // lr id before which we have to stop
1469 //
e1ef49ad 1470#ifdef _ITSU_DEBUG_
15b02d69 1471 if (AliDebugLevelClass()>2) {
1472 printf("Refit %d: Lr:%d (%f) -> Lr:%d (%f)\n",dir,lrStart,rCurr, lrStop,rDest);
9cdcba2c 1473 printf("Before refit: "); trc->Print();
15b02d69 1474 }
e1ef49ad 1475#endif
8b16dbae 1476 if (lrStop<0 || lrStart<0) AliFatal(Form("Failed to find start(%d) or last(%d) layers. Track from %.3f to %.3f",lrStart,lrStop,rCurr,rDest));
1477 //
9cdcba2c 1478 Int_t clInfo[2*AliITSUAux::kMaxLayers];
1479 Int_t nCl = fCurrHyp->FetchClusterInfo(clInfo);
1480 fCurrMass = fCurrHyp->GetMass();
7ce7f445 1481 AliITSUTrackHyp tmpTr(*(AliKalmanTrack*)trc);
9cdcba2c 1482 double chi2 = 0;
5e3587be 1483 int iclLr[2],nclLr;
1484 nclFit = 0;
c51c3124 1485 //
c7e63626 1486 int lrStop1 = lrStop+dir;
70cb7fe4 1487 for (int ilr=lrStart;ilr!=lrStop1;ilr+=dir) {
c51c3124 1488 AliITSURecoLayer* lr = fITS->GetLayer(ilr);
1489 if ( dir*(rCurr-lr->GetR(dir))>0) continue; // this layer is already passed
8b16dbae 1490 int ilrA2,ilrA = lr->GetActiveID();
1491 // passive layer or active w/o hits will be traversed on the way to next cluster
9cdcba2c 1492 if (!lr->IsActive() || clInfo[ilrA2=(ilrA<<1)]<0) continue;
8b16dbae 1493 //
8b16dbae 1494 // select the order in which possible 2 clusters (in case of the overlap) will be traversed and fitted
1495 nclLr=0;
1496 if (dir>0) { // clusters are stored in increasing radius order
9cdcba2c 1497 iclLr[nclLr++]=clInfo[ilrA2++];
1498 if (clInfo[ilrA2]>=0) iclLr[nclLr++]=clInfo[ilrA2];
8b16dbae 1499 }
1500 else {
9cdcba2c 1501 if ( clInfo[ilrA2+1]>=0 ) iclLr[nclLr++]=clInfo[ilrA2+1];
1502 iclLr[nclLr++]=clInfo[ilrA2];
8b16dbae 1503 }
1504 //
b515ad7e 1505 Bool_t transportedToLayer = kFALSE;
8b16dbae 1506 for (int icl=0;icl<nclLr;icl++) {
1507 AliITSUClusterPix* clus = (AliITSUClusterPix*)lr->GetCluster(iclLr[icl]);
1508 AliITSURecoSens* sens = lr->GetSensorFromID(clus->GetVolumeId());
e7d83d38 1509 if (!tmpTr.Rotate(sens->GetPhiTF())) {
e1ef49ad 1510#ifdef _ITSU_DEBUG_
fb0eae48 1511 AliESDtrack* trESD = fCurrHyp->GetESDTrack();
9cdcba2c 1512 AliDebug(2,Form("Failed on rotate to %f | ESDtrack#%d (MClb:%d)",sens->GetPhiTF(),trESD->GetID(),trESD->GetTPCLabel()));
e1ef49ad 1513#endif
9cdcba2c 1514 return -1;
e7d83d38 1515 }
b515ad7e 1516 //
1517 double xClus = sens->GetXTF()+clus->GetX();
1518 if (!transportedToLayer) {
1519 if (ilr!=lrStart && !TransportToLayerX(&tmpTr,lrStart,ilr,xClus)) {
e1ef49ad 1520#ifdef _ITSU_DEBUG_
846c381e 1521 AliESDtrack* trESD = fCurrHyp->GetESDTrack();
9cdcba2c 1522 AliDebug(2,Form("Failed to transport %d -> %d | ESDtrack#%d (MClb:%d)\n",lrStart,ilr,trESD->GetID(),trESD->GetTPCLabel()));
e1ef49ad 1523#endif
9cdcba2c 1524 return -1; // go to the entrance to the layer
b515ad7e 1525 }
1526 lrStart = ilr;
1527 transportedToLayer = kTRUE;
1528 }
1529 //
e1ef49ad 1530#ifdef _ITSU_DEBUG_
b515ad7e 1531 if (AliDebugLevelClass()>1) {
1532 AliDebug(2,Form("Propagate to cl:%d on lr %d Need to go %.4f -> %.4f",icl,ilrA,tmpTr.GetX(),xClus));
b515ad7e 1533 }
e1ef49ad 1534#endif
b515ad7e 1535 //
1536 if (!PropagateSeed(&tmpTr,xClus,fCurrMass)) {
e1ef49ad 1537#ifdef _ITSU_DEBUG_
846c381e 1538 AliESDtrack* trESD = fCurrHyp->GetESDTrack();
9cdcba2c 1539 AliDebug(2,Form("Failed on propagate to %f (dir=%d) | ESDtrack#%d (MClb:%d)",xClus,dir,trESD->GetID(),trESD->GetTPCLabel()));
e1ef49ad 1540#endif
9cdcba2c 1541 return -1;
e7d83d38 1542 }
76390254 1543 //
9cdcba2c 1544 Double_t p[2]={clus->GetY(), clus->GetZ()};
1545 Double_t cov[3]={clus->GetSigmaY2(), clus->GetSigmaYZ(), clus->GetSigmaZ2()};
1546 double chi2cl = tmpTr.GetPredictedChi2(p,cov);
1547 chi2 += chi2cl;
1548 //
08419930 1549#ifdef _ITSU_TUNING_MODE_
1550 TObjArray* dest = trc->GetLabel()>=0 ? fCHistoArrCorr : fCHistoArrFake;
781c79ac 1551 if (dest && fTrackPhaseID>kClus2Tracks) {
08419930 1552 //
1553 double htrPt = tmpTr.Pt();
9cdcba2c 1554 double dy = p[0]-tmpTr.GetY();
1555 double dz = p[1]-tmpTr.GetZ();
781c79ac 1556 ((TH2*)dest->At( GetHistoID(ilrA,kHResY,0,fTrackPhaseID) ))->Fill(htrPt,dy);
1557 ((TH2*)dest->At( GetHistoID(ilrA,kHResZ,0,fTrackPhaseID) ))->Fill(htrPt,dz);
76390254 1558 double errY = tmpTr.GetSigmaY2();
1559 double errZ = tmpTr.GetSigmaZ2();
781c79ac 1560 if (errY>0) ((TH2*)dest->At( GetHistoID(ilrA,kHResYP,0,fTrackPhaseID) ))->Fill(htrPt,dy/Sqrt(errY));
1561 if (errZ>0) ((TH2*)dest->At( GetHistoID(ilrA,kHResZP,0,fTrackPhaseID) ))->Fill(htrPt,dz/Sqrt(errZ));
1562 ((TH2*)dest->At( GetHistoID(ilrA,kHChi2Cl,0,fTrackPhaseID) ))->Fill(htrPt,chi2cl);
76390254 1563 }
1564#endif
9cdcba2c 1565 //
1566 if ( !tmpTr.Update(p,cov) ) {
e1ef49ad 1567#ifdef _ITSU_DEBUG_
846c381e 1568 AliESDtrack* trESD = fCurrHyp->GetESDTrack();
9cdcba2c 1569 AliDebug(2,Form("Failed on Update | ESDtrack#%d (MClb:%d)",trESD->GetID(),trESD->GetTPCLabel()));
e1ef49ad 1570#endif
9cdcba2c 1571 return -1;
e7d83d38 1572 }
e1ef49ad 1573#ifdef _ITSU_DEBUG_
b515ad7e 1574 if (AliDebugLevelClass()>1) {
9cdcba2c 1575 printf("AfterRefit: "); tmpTr.Print();
1576 }
e1ef49ad 1577#endif
5e3587be 1578 if (++nclFit==nCl && stopCond<0) {
7ce7f445 1579 *trc = (AliKalmanTrack&)tmpTr;
9cdcba2c 1580 return chi2; // it was requested to not propagate after last update
b515ad7e 1581 }
8b16dbae 1582 }
c51c3124 1583 //
1584 }
8b16dbae 1585 // All clusters were succesfully fitted. Even if the track does not reach rDest, this is enough to validate it.
1586 // Still, try to go as close as possible to rDest.
1587 //
70cb7fe4 1588 if (lrStart!=lrStop) {
e7d83d38 1589 if (!TransportToLayer(&tmpTr,lrStart,lrStop)) {
e1ef49ad 1590#ifdef _ITSU_DEBUG_
2947434b 1591 AliDebug(1,Form("Failed on TransportToLayer %d->%d | ESDtrack#%d (MClb:%d), X=%f",lrStart,lrStop,fCurrESDtrack->GetID(),fCurrESDtrMClb,tmpTr.GetX()));
e1ef49ad 1592#endif
2947434b 1593 return (stopCond>0) ? -chi2 : chi2; // rDest was obligatory
e7d83d38 1594 }
1595 if (!GoToExitFromLayer(&tmpTr,fITS->GetLayer(lrStop),dir)) {
e1ef49ad 1596#ifdef _ITSU_DEBUG_
2947434b 1597 AliDebug(1,Form("Failed on GoToExitFromLayer %d | ESDtrack#%d (MClb:%d), X=%f",lrStop,fCurrESDtrack->GetID(),fCurrESDtrMClb,tmpTr.GetX()));
e1ef49ad 1598#endif
2947434b 1599 return (stopCond>0) ? -chi2 : chi2; // rDest was obligatory
9cdcba2c 1600 }
70cb7fe4 1601 }
1602 // go to the destination radius. Note that here we don't select direction to avoid precision problems
1603 if (!tmpTr.GetXatLabR(rDest,rDest,GetBz(),0) || !PropagateSeed(&tmpTr,rDest,fCurrMass, 100, kFALSE)) {
2947434b 1604 return (stopCond>0) ? -chi2 : chi2; // rDest was obligatory
70cb7fe4 1605 }
7ce7f445 1606 *trc = (AliKalmanTrack&)tmpTr;
e1ef49ad 1607#ifdef _ITSU_DEBUG_
15b02d69 1608 if (AliDebugLevelClass()>2) {
9cdcba2c 1609 printf("After refit (now at lr %d): ",lrStop); trc->Print();
15b02d69 1610 }
e1ef49ad 1611#endif
9cdcba2c 1612 return chi2;
c51c3124 1613}
dd2117a2 1614
1615//__________________________________________________________________
1616void AliITSUTrackerGlo::CookMCLabel(AliITSUTrackHyp* hyp)
1617{
1618 // build MC label
1619 //
1620 const int kMaxLbPerCl = 3;
1621 int lbID[kMaxLayers*6],lbStat[kMaxLayers*6];
1622 Int_t lr,nLab=0,nCl=0;
1623 AliITSUSeed *seed = hyp->GetWinner();
1624 while(seed) {
1625 int clID = seed->GetLrCluster(lr);
1626 if (clID>=0) {
1627 AliCluster *cl = fITS->GetLayerActive(lr)->GetCluster(clID);
1628 nCl++;
1629 for (int imc=0;imc<kMaxLbPerCl;imc++) { // labels within single cluster
1630 int trLb = cl->GetLabel(imc);
704d9398 1631 if (trLb<0) break;
dd2117a2 1632 // search this mc track in already accounted ones
1633 int iLab;
1634 for (iLab=0;iLab<nLab;iLab++) if (lbID[iLab]==trLb) break;
1635 if (iLab<nLab) lbStat[iLab]++;
1636 else {
1637 lbID[nLab] = trLb;
1638 lbStat[nLab++] = 1;
1639 }
1640 } // loop over given cluster's labels
1641 }
1642 seed = (AliITSUSeed*)seed->GetParent();
1643 } // loop over clusters
1644 //
5785a9d8 1645 AliESDtrack* esdTr = hyp->GetESDTrack();
1646 int tpcLab = esdTr ? Abs(esdTr->GetTPCLabel()) : -kDummyLabel;
41343470 1647 if (nCl && nLab) {
dd2117a2 1648 int maxLab=0,nTPCok=0;
dd2117a2 1649 for (int ilb=nLab;ilb--;) {
1650 int st = lbStat[ilb];
1651 if (lbStat[maxLab]<st) maxLab=ilb;
1652 if (tpcLab==lbID[ilb]) nTPCok += st;
1653 }
1654 hyp->SetFakeRatio(1.-float(nTPCok)/nCl);
1655 hyp->SetLabel( nTPCok==nCl ? tpcLab : -tpcLab);
bdb92766 1656 hyp->SetITSLabel( lbStat[maxLab]==nCl ? lbID[maxLab] : -lbID[maxLab]); // winner label
dd2117a2 1657 return;
1658 }
1659 //
1660 hyp->SetFakeRatio(-1.);
5785a9d8 1661 hyp->SetLabel( -tpcLab );
dd2117a2 1662 hyp->SetITSLabel( kDummyLabel );
1663}
38997c3c 1664
1665//__________________________________________________________________
1666Bool_t AliITSUTrackerGlo::AddSeedBranch(AliITSUSeed* seed)
1667{
1668 // add new prolongation branch for the seed to temporary array. The seeds are sorted according to Compare f-n
6cd80116 1669 if (fNCandidatesAdded+fNBranchesAdded>=fLayerMaxCandidates ) {
1670 AliInfo(Form("Number of candidates at layer %d for current seed reached %d, increasing buffer",fCurrActLrID,fLayerMaxCandidates));
1671 fLayerMaxCandidates = 2*(fLayerMaxCandidates+1);
1672 AliITSUSeed** tmpArr = fLayerCandidates;
1673 fLayerCandidates = new AliITSUSeed*[fLayerMaxCandidates];
1674 memcpy(fLayerCandidates,tmpArr,(fNCandidatesAdded+fNBranchesAdded)*sizeof(AliITSUSeed*));
026bab91 1675 delete[] tmpArr; // delete only array, not objects
38997c3c 1676 }
1677 AliITSUSeed** branches = &fLayerCandidates[fNCandidatesAdded]; // note: fNCandidatesAdded is incremented after adding all branches of current seed
1678 int slot=fNBranchesAdded++;
1679 for (int slotF=slot;slotF--;) { // slotF is always slot-1
1680 AliITSUSeed* si = branches[slotF];
1681 if (si->Compare(seed)<0) break; // found the last seed with better quality
1682 // otherwise, shift the worse seed to the next slot
1683 branches[slot] = si;
1684 slot = slotF; // slot should be slotF+1
1685 }
1686 // if needed, move worse seeds
1687 branches[slot] = seed;
1688 return kTRUE;
1689 //
1690}
1691
1692//__________________________________________________________________
1693void AliITSUTrackerGlo::ValidateAllowedBranches(Int_t acceptMax)
1694{
1695 // keep allowed number of branches for current seed and disable extras
70cb7fe4 1696 AliCodeTimerAuto("",0);
38997c3c 1697 int nb = Min(fNBranchesAdded,acceptMax);
1698 // if (nb<fNBranchesAdded) printf("ValidateAllowedBranches: %d of %d (%d) on lr %d\n",nb,fNBranchesAdded,fNCandidatesAdded,ilr);
1699 // disable unused branches
1700 AliITSUSeed** branches = &fLayerCandidates[fNCandidatesAdded];
08419930 1701 //
1702#ifdef _ITSU_TUNING_MODE_
1703 for (int ib=0;ib<fNBranchesAdded;ib++) branches[ib]->SetOrdBranch(ib);
c03e4f8a 1704#endif
bdb92766 1705 //
08419930 1706 for (int ib=nb;ib<fNBranchesAdded;ib++) MarkSeedFree(branches[ib]);
38997c3c 1707 fNCandidatesAdded += nb; // update total candidates counter
1708 fNBranchesAdded = 0; // reset branches counter
1709 //
1710}
1711
1712//__________________________________________________________________
1713void AliITSUTrackerGlo::ValidateAllowedCandidates(Int_t ilr, Int_t acceptMax)
1714{
1715 // transfer allowed number of branches to hypothesis container
1716 //
9cdcba2c 1717 AliCodeTimerAuto("",0);
38997c3c 1718 // sort candidates in increasing order of chi2
6cd80116 1719 static int lastSize = 0;
1720 static int *index = 0;
1721 static float *chi2 = 0;
1722 if (fLayerMaxCandidates>lastSize) {
1723 lastSize = fLayerMaxCandidates;
1724 delete[] index;
1725 delete[] chi2;
1726 index = new int[lastSize];
1727 chi2 = new float[lastSize];
1728 }
38997c3c 1729 for (int i=fNCandidatesAdded;i--;) chi2[i] = fLayerCandidates[i]->GetChi2GloNrm();
1730 Sort(fNCandidatesAdded,chi2,index,kFALSE);
1731 //
9cdcba2c 1732 int nacc=0,nb=0;
097e4cdf 1733 if (ilr>fCurrTrackCond->GetActiveLrInner()) { // just take 1st acceptMax candidates
9cdcba2c 1734 nb = Min(fNCandidatesAdded,acceptMax);
1735 for (int ib=0;ib<nb;ib++) AddProlongationHypothesis(fLayerCandidates[index[ib]],ilr);
1736 }
1737 else { // for innermost layer test ITS SA backward chi2 and matching to TPC
1738 AliITSUSeed* wn0 = fCurrHyp->GetWinner(); // in principle, should be NULL
1739 for (nb=0;nb<fNCandidatesAdded;nb++) {
1740 AliITSUSeed* sd = fLayerCandidates[index[nb]];
1741 fCurrHyp->SetWinner(sd);
1742 if (!CheckBackwardMatching(sd)) MarkSeedFree(sd); // discard bad backward tracks
1743 else {
1744 AddProlongationHypothesis(sd,ilr);
1745 if (++nacc==acceptMax) {nb++; break;} // the rest will be discarded
1746 }
1747 }
1748 fCurrHyp->SetWinner(wn0); // restore original winner (NULL?)
1749 }
38997c3c 1750 // disable unused candidates
08419930 1751 for (int ib=nb;ib<fNCandidatesAdded;ib++) MarkSeedFree(fLayerCandidates[index[ib]]);
38997c3c 1752 fNCandidatesAdded = 0; // reset candidates counter
1753 //
1754}
1755
9cdcba2c 1756//__________________________________________________________________
1757Bool_t AliITSUTrackerGlo::CheckBackwardMatching(AliITSUSeed* seed)
1758{
1759 // check seed backward propagation chi2 and matching to TPC
1760 double bz0 = GetBz();
70cb7fe4 1761 double rDest = fITS->GetRITSTPCRef(); // reference radius for comparison
7ce7f445 1762 AliITSUTrackHyp trback;
1763 trback.AliExternalTrackParam::operator=(*seed);
9cdcba2c 1764 trback.ResetCovariance(10000);
5e3587be 1765 int nclFit = 0;
1766 double chi2sa = RefitTrack(&trback,rDest,nclFit,1);
9cdcba2c 1767 if (chi2sa<0) return kFALSE;
5e3587be 1768 int ndf = nclFit*2-5;
9cdcba2c 1769 if (ndf>0) chi2sa /= ndf;
5e3587be 1770 if (chi2sa>fCurrTrackCond->GetMaxITSSAChi2(nclFit)) return kFALSE;
9cdcba2c 1771 //
1772 // relate to TPC track at outer layer
1773 AliExternalTrackParam* tpcSeed = fCurrHyp->GetTPCSeed();
1774 if (!trback.Rotate(tpcSeed->GetAlpha()) || !trback.PropagateParamOnlyTo(tpcSeed->GetX(),bz0)) {
e1ef49ad 1775#ifdef _ITSU_DEBUG_
9cdcba2c 1776 if (AliDebugLevelClass()>+1 && !seed->ContainsFake()) {
1777 AliInfo(Form("Failed to propagate seed to TPC track @ X:%.3f Alpha:%.3f",tpcSeed->GetX(),tpcSeed->GetAlpha()));
1778 seed->Print("etp");
1779 trback.Print();
1780 }
e1ef49ad 1781#endif
9cdcba2c 1782 return kFALSE;
1783 }
1784 double chi2Match = trback.GetPredictedChi2(tpcSeed)/5;
1785 if (chi2Match>fCurrTrackCond->GetMaxITSTPCMatchChi2()) return kFALSE;
1786 //
1787 seed->SetChi2ITSSA(chi2sa);
1788 seed->SetChi2ITSTPC(chi2Match);
1789 return kTRUE;
1790}
1791
6cd80116 1792//__________________________________________________________________
1793void AliITSUTrackerGlo::SaveReducedHypothesesTree(AliITSUTrackHyp* dest)
1794{
1795 // remove those hypothesis seeds which dont lead to candidates at final layer
1796 //
1797 // 1st, flag the seeds to save
0e40ce66 1798 int lr0 = fCurrTrackCond->GetActiveLrInner();
6cd80116 1799 for (int isd=0;isd<fCurrHyp->GetNSeeds(lr0);isd++) {
1800 AliITSUSeed* seed = fCurrHyp->RemoveSeed(lr0,isd);
1801 if (!seed) continue;
1802 seed->FlagTree(AliITSUSeed::kSave);
1803 dest->AddSeed(seed,lr0);
1804 }
0e40ce66 1805 for (int ilr=0;ilr<fNLrActive;ilr++) {
6cd80116 1806 int nsd = fCurrHyp->GetNSeeds(ilr);
1807 for (int isd=0;isd<nsd;isd++) {
1808 AliITSUSeed* seed = fCurrHyp->RemoveSeed(ilr,isd);
1809 if (!seed) continue; // already discarded or saved
1810 if (seed->IsSaved()) dest->AddSeed(seed,ilr);
1811 else MarkSeedFree(seed);
1812 }
1813 }
1814 //
1815 // AliInfo(Form("SeedsPool: %d, BookedUpTo: %d, free: %d",fSeedsPool.GetSize(),fSeedsPool.GetEntriesFast(),fNFreeSeeds));
1816}
1817
0e40ce66 1818//__________________________________________________________________
1819void AliITSUTrackerGlo::CleanHypothesis(AliITSUTrackHyp* hyp)
1820{
1821 // clean hypothesis
1822 hyp->SetWinner(0);
1823 for (int ilr=0;ilr<fNLrActive;ilr++) {
1824 int nsd = hyp->GetNSeeds(ilr);
1825 for (int isd=0;isd<nsd;isd++) {
1826 AliITSUSeed* seed = hyp->RemoveSeed(ilr,isd);
1827 if (!seed) continue; // already discarded or saved
1828 MarkSeedFree(seed);
1829 }
1830 }
1831}
1832
53870004 1833//__________________________________________________________________
1834void AliITSUTrackerGlo::FlagSplitClusters()
1835{
1836 // set special bit on split clusters using MC info
6cd80116 1837 for (int ilr=fNLrActive;ilr--;) {
53870004 1838 int nsplit=0;
1839 AliITSURecoLayer* lr = fITS->GetLayerActive(ilr);
1840 for (int isn=lr->GetNSensors();isn--;) {
1841 AliITSURecoSens* sens = lr->GetSensor(isn);
1842 int nCl = sens->GetNClusters();
1843 if (!nCl) continue;
1844 int cl0 = sens->GetFirstClusterId();
1845 for (int icl=nCl;icl--;) {
1846 AliITSUClusterPix *cl = (AliITSUClusterPix*)lr->GetCluster(cl0+icl);
1847 for (int icl1=icl;icl1--;) {
1848 AliITSUClusterPix *cl1 = (AliITSUClusterPix*)lr->GetCluster(cl0+icl1);
1849 if (cl->HasCommonTrack(cl1)) {
1850 if (!cl->IsSplit()) nsplit++;
1851 if (!cl1->IsSplit()) nsplit++;
1852 cl->SetSplit();
1853 cl1->SetSplit();
1854 }
1855 }
1856 }
1857 }
1858 AliInfo(Form("%4d out of %4d clusters are split",nsplit,lr->GetNClusters()));
1859 }
1860 //
1861}
1862
1863//__________________________________________________________________
1864Bool_t AliITSUTrackerGlo::ContainsSplitCluster(const AliITSUSeed* seed, Int_t maxSize)
1865{
1866 // check if the seed contains split cluster with size < maxSize
1867 int lrID,clID;
1868 if ( (clID=seed->GetLrCluster(lrID))>=0 ) {
1869 AliITSUClusterPix* cl = (AliITSUClusterPix*)fITS->GetLayerActive(lrID)->GetCluster(clID);
1870 if (cl->IsSplit() && cl->GetNPix()<maxSize ) return kTRUE;
1871 }
1872 return seed->GetParent() ? ContainsSplitCluster((AliITSUSeed*)seed->GetParent(),maxSize) : kFALSE;
1873 //
1874}
1875
1876//__________________________________________________________________
1877void AliITSUTrackerGlo::PrintSeedClusters(const AliITSUSeed* seed, Option_t* option)
1878{
1879 // print seeds clusters
1880 int lrID,clID;
1881 if ( (clID=seed->GetLrCluster(lrID))>=0 ) {
1882 AliITSUClusterPix* cl = (AliITSUClusterPix*)fITS->GetLayerActive(lrID)->GetCluster(clID);
1883 cl->Print(option);
1884 }
1885 if (seed->GetParent()) PrintSeedClusters((AliITSUSeed*)seed->GetParent(), option);
1886 //
1887}
1888
c03e4f8a 1889//__________________________________________________________________
08419930 1890void AliITSUTrackerGlo::FlagSeedClusters(const AliITSUSeed* seed, Bool_t flg, UShort_t ref)
c03e4f8a 1891{
1892 // mark used clusters
1893 int lrID,clID;
1894 while (seed) {
08419930 1895 if ( (clID=seed->GetLrCluster(lrID))>=0 ) {
1896 AliITSUClusterPix* cl = (AliITSUClusterPix*)fITS->GetLayerActive(lrID)->GetCluster(clID);
1897 if (ref) { // do we need to set or delete cluster->track ref?
1898 if (flg) {
1899 if (!cl->GetRecoInfo()) cl->SetRecoInfo(ref); // set ref only if cluster already does not refer to other track, inc.counter
1900 }
1901 else {
1902 if (cl->GetRecoInfo()==ref) cl->SetRecoInfo(0); // unset reference only if it refers to ref, decrease counter
1903 }
1904 }
1905 }
c03e4f8a 1906 seed = (AliITSUSeed*)seed->GetParent();
1907 }
1908 //
1909}
1910
d99994b2 1911//__________________________________________________________________
1912void AliITSUTrackerGlo::CheckClusterUsage()
1913{
1914 // check cluster usage info
eb18cdc6 1915 printf("ClusterUsage at pass %d\n",fCurrPassID);
d99994b2 1916 for (int ilr=0;ilr<fNLrActive;ilr++) {
1917 AliITSURecoLayer* lr = fITS->GetLayerActive(ilr);
1918 int ncl = lr->GetNClusters();
1919 int nclUsed = 0;
1920 int nclUs0CntShr1 = 0;
1921 int nclUs1Cnt0 = 0;
1922 int nclUs1Shr1 = 0;
1923 int nclUseL = 0;
1924 for (int icl=0;icl<ncl;icl++) {
1925 AliITSUClusterPix *cl = (AliITSUClusterPix*)lr->GetCluster(icl);
1926 int usc = cl->GetClUsage();
1927 if (usc<1) {
1928 if (cl->IsClusterUsed() || cl->IsClusterShared()) {
1929 nclUs0CntShr1++;
1930 printf("Lr%d Cl%4d: UseCnt=0 but UseFlg=%d UseShr=%d\n",ilr,icl,cl->IsClusterUsed(),cl->IsClusterShared());
1931 }
1932 continue;
1933 }
1934 nclUsed++;
1935 if (!cl->IsClusterUsed()) {
1936 nclUs1Cnt0++;
1937 printf("Lr%d Cl%4d: UseCnt=%d but UseFlg=%d UseShr=%d\n",ilr,icl,usc, cl->IsClusterUsed(),cl->IsClusterShared());
1938 }
1939 if (cl->IsClusterShared()) {
1940 nclUs1Shr1++;
1941 printf("Lr%d Cl%4d: UseCnt=%d but UseFlg=%d UseShr=%d\n",ilr,icl,usc, cl->IsClusterUsed(),cl->IsClusterShared());
1942 }
1943 if (usc>1) {
1944 nclUseL++;
1945 printf("Lr%d Cl%4d: UseCnt=%d!! but UseFlg=%d UseShr=%d\n",ilr,icl,usc, cl->IsClusterUsed(),cl->IsClusterShared());
1946 }
1947 }
1948 printf("ClUsage Lr%d: Ncl=%5d Nusd=%5d (%5.2f) Use0CS1:%4d Use1C0:%4d Use1S1:%4d UseL:%d\n",
1949 ilr,ncl,nclUsed,ncl? float(nclUsed)/ncl : 0, nclUs0CntShr1,nclUs1Cnt0,nclUs1Shr1,nclUseL);
1950 }
1951}
1952
1953
08419930 1954#ifdef _ITSU_TUNING_MODE_
53870004 1955//__________________________________________________________________
781c79ac 1956TObjArray* AliITSUTrackerGlo::BookControlHistos(const char* pref)
53870004 1957{
1958 // book special control histos
08419930 1959 TString prefS = pref;
1960 prefS.ToLower();
781c79ac 1961 TObjArray* dest = new TObjArray;
08419930 1962 dest->SetOwner(kTRUE);
1963 //
bdb92766 1964 const int kNResDef=7;
1965 const double kResDef[kNResDef]={0.05,0.05,0.3, 0.05,1,0.5,1.5};
bdb92766 1966 const double ptMax=10;
1967 const double plMax=10;
1968 const double chiMax=100;
1969 const int nptbins=50;
1970 const int nresbins=400;
1971 const int nplbins=50;
1972 const int nchbins=200;
1973 const int maxBr = 15;
1974 const int maxCand = 200;
1975 TString ttl;
781c79ac 1976 for (int phase=0;phase<kNTrackingPhases;phase++) {
1977 for (int pass=0;pass<AliITSUReconstructor::GetRecoParam()->GetNTrackingConditions();pass++) {
bdb92766 1978 //
781c79ac 1979 for (int ilr=0;ilr<fNLrActive;ilr++) {
1980 //
1981 // ----------------- These are histos to be filled in Cluster2Tracks of each pass.
1982 // PropagateBack and RefitInward will be stored among the histos of 1st pass
1983 if (pass>0 && phase!=kClus2Tracks) continue;
1984 //
1985 double mxdf = ilr>=kNResDef ? kResDef[kNResDef-1] : kResDef[ilr];
1986 ttl = Form("Pass%d_S%d_residY%d_%s",pass,phase,ilr,pref);
1987 TH2F* hdy = new TH2F(ttl.Data(),ttl.Data(),nptbins,0,ptMax,nresbins,-mxdf,mxdf);
1988 dest->AddAtAndExpand(hdy,GetHistoID(ilr,kHResY,pass,phase));
1989 hdy->SetDirectory(0);
1990 //
1991 ttl = Form("Pass%d_S%d_residYPull%d_%s",pass,phase,ilr,pref);
1992 TH2F* hdyp = new TH2F(ttl.Data(),ttl.Data(),nptbins,0,ptMax,nplbins,-plMax,plMax);
1993 dest->AddAtAndExpand(hdyp,GetHistoID(ilr,kHResYP,pass,phase));
1994 hdyp->SetDirectory(0);
1995 //
1996 ttl = Form("Pass%d_S%d_residZ%d_%s",pass,phase,ilr,pref);
1997 TH2F* hdz = new TH2F(ttl.Data(),ttl.Data(),nptbins,0,ptMax,nresbins,-mxdf,mxdf);
1998 dest->AddAtAndExpand(hdz,GetHistoID(ilr,kHResZ,pass,phase));
1999 hdz->SetDirectory(0);
2000 //
2001 ttl = Form("Pass%d_S%d_residZPull%d_%s",pass,phase,ilr,pref);
2002 TH2F* hdzp = new TH2F(ttl.Data(),ttl.Data(),nptbins,0,ptMax,nplbins,-plMax,plMax);
2003 hdzp->SetDirectory(0);
2004 dest->AddAtAndExpand(hdzp,GetHistoID(ilr,kHResZP,pass,phase));
2005 //
2006 ttl = Form("Pass%d_S%d_chi2Cl%d_%s",pass,phase,ilr,pref);
2007 TH2F* hchi = new TH2F(ttl.Data(),ttl.Data(),nptbins,0,ptMax, nchbins,0.,chiMax);
2008 hchi->SetDirectory(0);
2009 dest->AddAtAndExpand(hchi,GetHistoID(ilr,kHChi2Cl,pass,phase));
2010 //
2011 // ------------------- These histos are filled for Clusters2Tracks only
2012 if (phase!=kClus2Tracks) continue;
2013 //
2014 ttl = Form("Pass%d_S%d_chi2Nrm%d_%s",pass,phase,ilr,pref);
2015 TH2F* hchiN = new TH2F(ttl.Data(),ttl.Data(),nptbins,0,ptMax, nchbins,0.,chiMax);
2016 hchiN->SetDirectory(0);
2017 dest->AddAtAndExpand(hchiN,GetHistoID(ilr,kHChi2Nrm,pass,phase));
2018 //
2019 ttl = Form("Pass%d_S%d_bestInBranch%d_%s",pass,phase,ilr,pref);
bdb92766 2020 TH2* hnbr = new TH2F(ttl.Data(),ttl.Data(),nptbins,0,ptMax, maxBr,-0.5,maxBr-0.5);
2021 hnbr->SetDirectory(0);
781c79ac 2022 dest->AddAtAndExpand(hnbr,GetHistoID(ilr,kHBestInBranch,pass,phase));
c03e4f8a 2023 //
781c79ac 2024 ttl = Form("Pass%d_S%d_bestInCands%d_%s",pass,phase,ilr,pref);
bdb92766 2025 TH2* hncn = new TH2F(ttl.Data(),ttl.Data(),nptbins,0,ptMax, maxCand,-0.5,maxCand-0.5);
2026 hncn->SetDirectory(0);
781c79ac 2027 dest->AddAtAndExpand(hncn,GetHistoID(ilr,kHBestInCand,pass,phase));
c03e4f8a 2028 //
781c79ac 2029 } // loop over layers
2030 //
2031 // these are special histos, filled not per layer but in the end of track fit in Clusters2Tracks in EVERY pass
2032 //
2033 if (phase==kClus2Tracks) {
2034 TH2* hchiMatch = 0;
2035 ttl = Form("Pass%d_Chi2Match_%s",pass,pref);
2036 hchiMatch = new TH2F(ttl.Data(),ttl.Data(),nptbins,0,ptMax, nchbins,0.,chiMax);
2037 hchiMatch->SetDirectory(0);
2038 dest->AddAtAndExpand(hchiMatch,GetHistoID(-1,kHChiMatch,pass,phase));
2039 //
2040 TH2* hchiSA = 0;
2041 ttl = Form("Pass%d_Chi2ITSSA_%s",pass,pref);
2042 hchiSA = new TH2F(ttl.Data(),ttl.Data(),nptbins,0,ptMax, nchbins/2,0.,chiMax/2);
2043 hchiSA->SetDirectory(0);
2044 dest->AddAtAndExpand(hchiSA,GetHistoID(-1,kHChiITSSA,pass,phase));
2045 //
53870004 2046 }
781c79ac 2047 } // loop over tracking passes
2048 }// loop over tracking phases
2049 //
2050 return dest;
53870004 2051}
781c79ac 2052
2053//__________________________________________________________________
2054Int_t AliITSUTrackerGlo::GetHistoID(Int_t lr, Int_t hid, Int_t pass, Int_t phase)
2055{
2056 // get id for the requested histo
2057 if (lr<0) lr=-1;
2058 lr++;
2059 return pass*kHistosPass + phase*kHistosPhase + lr*kMaxHID + hid;
2060 //
2061}
2062
53870004 2063#endif