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