Merge branch 'master' of https://git.cern.ch/reps/AliRoot
[u/mrichter/AliRoot.git] / MFT / AliMFTTrackerMU.cxx
CommitLineData
65273893 1/*************************************************************************
4cc511f4
AU
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//
18// MFT tracker
19//
20// Class for the creation of the "global muon tracks" built from the tracks reconstructed in the
21// muon spectrometer and the clusters of the Muon Forward Tracker
22//
23// Contact author: antonio.uras@cern.ch
24//
25//====================================================================================================================================================
26
27#include "TTree.h"
28#include "AliLog.h"
29#include "AliGeomManager.h"
30#include "AliESDEvent.h"
31#include "AliESDMuonTrack.h"
026547c6 32#include "AliESDMuonGlobalTrack.h"
4cc511f4
AU
33#include "AliMFTTrackerMU.h"
34#include "TMath.h"
35#include "AliRun.h"
65273893 36#include "AliRunLoader.h"
37#include "AliHeader.h"
38#include "AliGenEventHeader.h"
39#include "TArrayF.h"
4cc511f4
AU
40#include "AliMFT.h"
41#include "AliMUONTrackExtrap.h"
42#include "AliMUONTrack.h"
43#include "AliMUONESDInterface.h"
44#include "AliMuonForwardTrack.h"
026547c6 45#include "AliMUONConstants.h"
4cc511f4
AU
46
47ClassImp(AliMFTTrackerMU)
48
49const Double_t AliMFTTrackerMU::fRadLengthSi = AliMFTConstants::fRadLengthSi;
50
51//====================================================================================================================================================
52
53AliMFTTrackerMU::AliMFTTrackerMU() :
54 AliTracker(),
55 fESD(0),
56 fMFT(0),
57 fSegmentation(0),
58 fNPlanesMFT(0),
59 fNPlanesMFTAnalyzed(0),
e21f8bf5 60 fSigmaClusterCut(2),
4cc511f4
AU
61 fScaleSigmaClusterCut(1.),
62 fNMaxMissingMFTClusters(0),
63 fGlobalTrackingDiverged(kFALSE),
64 fCandidateTracks(0),
65 fMUONTrack(0),
66 fCurrentTrack(0),
67 fFinalBestCandidate(0),
65273893 68 fXExtrapVertex(0),
69 fYExtrapVertex(0),
70 fZExtrapVertex(0),
71 fXExtrapVertexError(0),
72 fYExtrapVertexError(0),
4cc511f4
AU
73 fBransonCorrection(kFALSE)
74{
75
76 //--------------------------------------------------------------------
77 // This is the AliMFTTrackerMU constructor
78 //--------------------------------------------------------------------
79
80 fMFT = (AliMFT*) gAlice->GetDetector("MFT");
81 fSegmentation = fMFT->GetSegmentation();
82 SetNPlanesMFT(fSegmentation->GetNPlanes());
83 AliMUONTrackExtrap::SetField(); // set the magnetic field for track extrapolations
84
e21f8bf5 85 AliInfo(Form("fMFT = %p, fSegmentation = %p", fMFT, fSegmentation));
86
4cc511f4 87 for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
e21f8bf5 88 fMFTClusterArray[iPlane] = new TClonesArray("AliMFTCluster");
4cc511f4
AU
89 fMFTClusterArrayFront[iPlane] = new TClonesArray("AliMFTCluster");
90 fMFTClusterArrayBack[iPlane] = new TClonesArray("AliMFTCluster");
e21f8bf5 91 fMFTClusterArray[iPlane] -> SetOwner(kTRUE);
92 fMFTClusterArrayFront[iPlane] -> SetOwner(kTRUE);
93 fMFTClusterArrayBack[iPlane] -> SetOwner(kTRUE);
94 fMinResearchRadiusAtPlane[iPlane] = 0.;
4cc511f4
AU
95 }
96
e21f8bf5 97 fCandidateTracks = new TClonesArray("AliMuonForwardTrack",50000);
98
4cc511f4
AU
99}
100
101//====================================================================================================================================================
102
103AliMFTTrackerMU::~AliMFTTrackerMU() {
104
105 // destructor
106
107 for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
e21f8bf5 108 fMFTClusterArray[iPlane] -> Delete();
4cc511f4
AU
109 delete fMFTClusterArray[iPlane];
110 delete fMFTClusterArrayFront[iPlane];
111 delete fMFTClusterArrayBack[iPlane];
112 }
113
e21f8bf5 114 delete fCandidateTracks;
115
4cc511f4
AU
116}
117
118//====================================================================================================================================================
119
120Int_t AliMFTTrackerMU::LoadClusters(TTree *cTree) {
121
122 //--------------------------------------------------------------------
123 // This function loads the MFT clusters
124 //--------------------------------------------------------------------
125
e21f8bf5 126 for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
127 AliDebug(1, Form("Setting Address for Branch Plane_%02d", iPlane));
128 cTree->SetBranchAddress(Form("Plane_%02d",iPlane), &fMFTClusterArray[iPlane]);
129 }
130
4cc511f4
AU
131 if (!cTree->GetEvent()) return kFALSE;
132 for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
e21f8bf5 133 AliInfo(Form("plane %02d: nClusters = %d\n", iPlane, fMFTClusterArray[iPlane]->GetEntries()));
4cc511f4
AU
134 }
135 SeparateFrontBackClusters();
136
137 return 0;
138
139}
140
141//====================================================================================================================================================
142
143void AliMFTTrackerMU::UnloadClusters() {
144
65273893 145 //--------------------------------------------------------------------
146 // This function unloads MFT clusters
147 //--------------------------------------------------------------------
4cc511f4
AU
148
149 for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
150 fMFTClusterArray[iPlane] -> Clear("C");
151 fMFTClusterArrayFront[iPlane] -> Clear("C");
152 fMFTClusterArrayBack[iPlane] -> Clear("C");
153 }
154
155}
156
157//====================================================================================================================================================
158
159Int_t AliMFTTrackerMU::Clusters2Tracks(AliESDEvent *event) {
160
161 //--------------------------------------------------------------------
162 // This functions reconstructs the Muon Forward Tracks
163 // The clusters must be already loaded !
164 //--------------------------------------------------------------------
165
e21f8bf5 166 // Tree of AliMuonForwardTrack objects. Created outside the ESD framework for cross-check purposes
167
168 TFile *outputFileMuonGlobalTracks = new TFile("MuonGlobalTracks.root", "update");
4cc511f4
AU
169 TTree *outputTreeMuonGlobalTracks = new TTree("AliMuonForwardTracks", "Tree of AliMuonForwardTracks");
170 TClonesArray *muonForwardTracks = new TClonesArray("AliMuonForwardTrack");
171 outputTreeMuonGlobalTracks -> Branch("tracks", &muonForwardTracks);
172
173 //--------------------------------------------------------------------
174
175 fESD = event;
176
65273893 177 fXExtrapVertex = 0;
178 fYExtrapVertex = 0;
179 fZExtrapVertex = 0;
180 fXExtrapVertexError = AliMFTConstants::fXVertexTolerance;
181 fYExtrapVertexError = AliMFTConstants::fYVertexTolerance;
e21f8bf5 182
65273893 183 // Imposing a fixed primary vertex, in order to keep memory of its position. Taking the primary vertex would imply the risk
184 // of loosing memory of its position when passing from ESD to AOD, due to possible refitting
e21f8bf5 185 const AliESDVertex* esdVert = fESD->GetVertex();
186 if (esdVert->GetNContributors() > 0 || !strcmp(esdVert->GetTitle(),"vertexer: smearMC")) {
65273893 187 fXExtrapVertex = esdVert->GetX();
188 fYExtrapVertex = esdVert->GetY();
189 fZExtrapVertex = esdVert->GetZ();
190 fXExtrapVertexError = TMath::Max(AliMFTConstants::fXVertexTolerance, esdVert->GetXRes());
191 fYExtrapVertexError = TMath::Max(AliMFTConstants::fYVertexTolerance, esdVert->GetYRes());
192 AliInfo(Form("Found ESD vertex from %d contributors (%f +/- %f, %f +/- %f, %f)",
193 esdVert->GetNContributors(),fXExtrapVertex,fXExtrapVertexError,fYExtrapVertex,fYExtrapVertexError,fZExtrapVertex));
e21f8bf5 194 }
65273893 195 else GetVertexFromMC();
196
4cc511f4
AU
197 //----------- Read ESD MUON tracks -------------------
198
199 Int_t nTracksMUON = event->GetNumberOfMuonTracks();
200
201 AliInfo(Form("Number of ESD MUON tracks: %d\n", nTracksMUON));
202
203 Int_t iTrack=0;
204 while (iTrack<nTracksMUON) {
205
206 fNPlanesMFTAnalyzed = 0;
207
e21f8bf5 208 AliInfo("****************************************************************************************");
209 AliInfo(Form("*************************** MUON TRACK %3d/%d ***************************************", iTrack, nTracksMUON));
210 AliInfo("****************************************************************************************");
4cc511f4
AU
211
212 fCandidateTracks -> Delete();
213
214 fNPlanesMFTAnalyzed = 0;
215
216 const AliESDMuonTrack *esdTrack = event->GetMuonTrack(iTrack);
e21f8bf5 217 if (fMUONTrack) delete fMUONTrack;
218 fMUONTrack = new AliMUONTrack();
219
f7cc8591 220 AliMUONESDInterface::ESDToMUON(*esdTrack, *fMUONTrack, kFALSE);
221
222 if (!fMUONTrack->GetTrackParamAtCluster()->First()) {
223 AliInfo("Skipping track, no parameters available!!!");
224 iTrack++;
225 continue;
226 }
4cc511f4
AU
227
228 // the track we are going to build, starting from fMUONTrack and adding the MFT clusters
229 AliMuonForwardTrack *track = new ((*fCandidateTracks)[0]) AliMuonForwardTrack();
230 track -> SetMUONTrack(new AliMUONTrack(*fMUONTrack));
231 track -> SetMCLabel(fMUONTrack->GetMCLabel());
232 track -> SetMatchTrigger(fMUONTrack->GetMatchTrigger());
233
026547c6 234 // track parameters linearly extrapolated from the first tracking station to the end of the absorber
235 AliMUONTrackParam trackParamEndOfAbsorber(*((AliMUONTrackParam*)(fMUONTrack->GetTrackParamAtCluster()->First())));
236 AliMUONTrackExtrap::ExtrapToZCov(&trackParamEndOfAbsorber, AliMUONConstants::AbsZEnd()); // absorber extends from -90 to -503 cm
237 Double_t xEndOfAbsorber = trackParamEndOfAbsorber.GetNonBendingCoor();
238 Double_t yEndOfAbsorber = trackParamEndOfAbsorber.GetBendingCoor();
239 Double_t rAbsorber = TMath::Sqrt(xEndOfAbsorber*xEndOfAbsorber + yEndOfAbsorber*yEndOfAbsorber);
240 track -> SetRAtAbsorberEnd(rAbsorber);
241
4cc511f4
AU
242 //------------------------- NOW THE CYCLE OVER THE MFT PLANES STARTS ---------------------------------------
243
244 for (Int_t iPlane=fNPlanesMFT-1; iPlane>=0; iPlane--) { /* *** do not reverse the order of this cycle!!!
245 *** this reflects the fact that the extrapolation is performed
246 *** starting from the last MFT plane back to the origin */
247
248 // --------- updating the array of candidates according to the clusters available in the i-th plane ---------
249
250 fNPlanesMFTAnalyzed++;
251
252 Int_t nCandidates = fCandidateTracks->GetEntriesFast();
253 for (Int_t iCandidate=0; iCandidate<nCandidates; iCandidate++) {
254 fCurrentTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(iCandidate);
e21f8bf5 255
4cc511f4
AU
256 // if the old track is compatible with the new cluster, the track is updated and inserted as new track in the array
257 // (several new tracks can be created for one old track)
258 if (FindClusterInPlane(iPlane) == kDiverged) {
259 fGlobalTrackingDiverged = kTRUE;
260 break;
261 }
262
263 if ((fNPlanesMFTAnalyzed-fCurrentTrack->GetNMFTClusters())>fNMaxMissingMFTClusters || fIsPlaneMandatory[iPlane]) {
264 fCandidateTracks->Remove(fCurrentTrack); // the old track is removed after the check;
265 }
266 }
267 if (fGlobalTrackingDiverged) {
268 if (fScaleSigmaClusterCut>0) fScaleSigmaClusterCut -= 0.1;
269 continue;
270 }
271
272 fCandidateTracks->Compress();
273
274 }
e21f8bf5 275
4cc511f4
AU
276 // -------------------------- END OF THE CYCLE OVER THE MFT PLANES --------------------------------------------
277
278 fGlobalTrackingDiverged = kFALSE;
279 fScaleSigmaClusterCut = 1.0;
280
281 AliDebug(1, "Finished cycle over planes");
282
283 iTrack++;
284
285 // If we have several final tracks, we must find the best candidate:
286
287 Int_t nFinalTracks = fCandidateTracks->GetEntriesFast();
e21f8bf5 288 AliInfo(Form("nFinalTracks = %d", nFinalTracks));
289
4cc511f4
AU
290 Int_t nGoodClustersBestCandidate = 0;
291 Int_t idBestCandidate = 0;
292 Double_t bestChi2 = -1.; // variable defining the best candidate
293
294 for (Int_t iFinalCandidate=0; iFinalCandidate<nFinalTracks; iFinalCandidate++) {
295
296 AliMuonForwardTrack *finalTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(iFinalCandidate);
297 Int_t nMFTClusters = finalTrack->GetNMFTClusters();
298
299 Double_t chi2 = 0;
300 for (Int_t iCluster=0; iCluster<nMFTClusters; iCluster++) {
301 AliMFTCluster *localCluster = finalTrack->GetMFTCluster(iCluster);
302 chi2 += localCluster->GetLocalChi2();
303 }
304 chi2 /= nMFTClusters;
305
306 // now comparing the tracks in order to find the best one
307
308 if (chi2<bestChi2 || bestChi2<0) {
309 bestChi2 = chi2;
310 idBestCandidate = iFinalCandidate;
311 }
312
313 }
314
315 if (nFinalTracks) {
316
317 AliMuonForwardTrack *newTrack = (AliMuonForwardTrack*) fCandidateTracks->UncheckedAt(idBestCandidate);
318 newTrack -> SetNWrongClustersMC(newTrack->GetNMFTClusters() - nGoodClustersBestCandidate);
319
e21f8bf5 320 new ((*muonForwardTracks)[muonForwardTracks->GetEntries()]) AliMuonForwardTrack(*newTrack);
4cc511f4 321
e21f8bf5 322 //----------------------- Save the information to the AliESDMuonGlobalTrack object
323
65273893 324 newTrack -> EvalKinem(AliMFTConstants::fZEvalKinem);
e21f8bf5 325
326 AliDebug(2,"Creating a new Muon Global Track");
327 AliESDMuonGlobalTrack *myESDTrack = event->NewMuonGlobalTrack();
328 myESDTrack -> SetPxPyPz(newTrack->Px(), newTrack->Py(), newTrack->Pz());
026547c6 329
f7cc8591 330 myESDTrack -> SetLabel(newTrack->GetMCLabel());
331 myESDTrack -> SetChi2OverNdf(newTrack->GetChi2OverNdf());
332 myESDTrack -> SetCharge(newTrack->GetCharge());
333 myESDTrack -> SetMatchTrigger(newTrack->GetMatchTrigger());
334 myESDTrack -> SetNMFTClusters(newTrack->GetNMFTClusters());
335 myESDTrack -> SetNWrongMFTClustersMC(newTrack->GetNWrongClustersMC());
336 myESDTrack -> SetFirstTrackingPoint(newTrack->GetMFTCluster(0)->GetX(), newTrack->GetMFTCluster(0)->GetY(), newTrack->GetMFTCluster(0)->GetZ());
65273893 337 myESDTrack -> SetXYAtVertex(newTrack->GetOffsetX(0., AliMFTConstants::fZEvalKinem), newTrack->GetOffsetY(0., AliMFTConstants::fZEvalKinem));
f7cc8591 338 myESDTrack -> SetRAtAbsorberEnd(newTrack->GetRAtAbsorberEnd());
19393921 339 // myESDTrack -> SetCovariances(newTrack->GetTrackParamAtMFTCluster(0)->GetCovariances()); // waiting for commit from Peter...
f7cc8591 340 myESDTrack -> SetChi2MatchTrigger(esdTrack->GetChi2MatchTrigger());
341 myESDTrack -> SetMuonClusterMap(esdTrack->GetMuonClusterMap());
342 myESDTrack -> SetHitsPatternInTrigCh(esdTrack->GetHitsPatternInTrigCh());
343 myESDTrack -> SetHitsPatternInTrigChTrk(esdTrack->GetHitsPatternInTrigChTrk());
344 myESDTrack -> Connected(esdTrack->IsConnected());
4cc511f4 345
65273893 346 ULong_t mftClusterPattern = 0;
347 for (Int_t iCluster=0; iCluster<newTrack->GetNMFTClusters(); iCluster++) {
348 AliMFTCluster *localCluster = newTrack->GetMFTCluster(iCluster);
349 mftClusterPattern |= 1 << localCluster->GetPlane();
350 mftClusterPattern |= IsCorrectMatch(localCluster, newTrack->GetMCLabel()) << fNMaxPlanes+localCluster->GetPlane();
351 }
352 myESDTrack -> SetMFTClusterPattern(mftClusterPattern);
353
e21f8bf5 354 //---------------------------------------------------------------------------------
4cc511f4
AU
355
356 }
357
4cc511f4
AU
358 fFinalBestCandidate = NULL;
359
360 }
361
e21f8bf5 362 outputTreeMuonGlobalTracks -> Fill();
363
4cc511f4 364 Int_t myEventID = 0;
4cc511f4
AU
365 while (outputFileMuonGlobalTracks->cd(Form("Event%d",myEventID))) myEventID++;
366 outputFileMuonGlobalTracks -> mkdir(Form("Event%d",myEventID));
367 outputFileMuonGlobalTracks -> cd(Form("Event%d",myEventID));
368 outputTreeMuonGlobalTracks -> Write();
369 outputFileMuonGlobalTracks -> Close();
370
e21f8bf5 371 muonForwardTracks -> Delete();
372 delete muonForwardTracks;
373
4cc511f4
AU
374 return 0;
375
376}
377
378//=========================================================================================================================================
379
380void AliMFTTrackerMU::SeparateFrontBackClusters() {
381
382 for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
383 fMFTClusterArrayFront[iPlane]->Delete();
384 fMFTClusterArrayBack[iPlane] ->Delete();
385 for (Int_t iCluster=0; iCluster<fMFTClusterArray[iPlane]->GetEntries(); iCluster++) {
386 AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArray[iPlane]->At(iCluster);
387 if (TMath::Abs(cluster->GetZ())<TMath::Abs(fSegmentation->GetPlane(iPlane)->GetZCenter())) {
388 new ((*fMFTClusterArrayFront[iPlane])[fMFTClusterArrayFront[iPlane]->GetEntries()]) AliMFTCluster(*cluster);
389 }
390 else {
391 new ((*fMFTClusterArrayBack[iPlane])[fMFTClusterArrayBack[iPlane]->GetEntries()]) AliMFTCluster(*cluster);
392 }
393 }
394 }
395
396}
397
398//==========================================================================================================================================
399
400Int_t AliMFTTrackerMU::FindClusterInPlane(Int_t planeId) {
401
4cc511f4
AU
402 // !!!!!!!!! coordinates and errors on the interaction vertex should be taken from the event itself (ITS) if available
403
404 // propagate track to plane #planeId (both to front and back active sensors)
405 // look for compatible clusters
406 // update TrackParam at found cluster (if any) using Kalman Filter
407
408 AliMUONTrackParam currentParamFront, currentParamBack, currentParamForResearchFront, currentParamForResearchBack;
409
410 if (planeId == fNPlanesMFT-1) { // last plane of the telecope
411 currentParamFront = (*((AliMUONTrackParam*)(fMUONTrack->GetTrackParamAtCluster()->First())));
412 currentParamBack = (*((AliMUONTrackParam*)(fMUONTrack->GetTrackParamAtCluster()->First())));
413 currentParamForResearchFront = currentParamFront;
414 currentParamForResearchBack = currentParamBack;
4cc511f4 415 if (fBransonCorrection) {
65273893 416 AliMUONTrackExtrap::ExtrapToVertex(&currentParamFront, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
417 AliMUONTrackExtrap::ExtrapToVertex(&currentParamBack, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
4cc511f4
AU
418 }
419 else {
65273893 420 AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&currentParamFront, fZExtrapVertex);
421 AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&currentParamBack, fZExtrapVertex);
4cc511f4 422 }
65273893 423 AliMUONTrackExtrap::ExtrapToVertex(&currentParamForResearchFront, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
424 AliMUONTrackExtrap::ExtrapToVertex(&currentParamForResearchBack, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
4cc511f4
AU
425 }
426 else { // MFT planes others than the last one: mult. scattering correction because of the upstream MFT planes is performed
427 currentParamFront = (*((AliMUONTrackParam*)(fCurrentTrack->GetTrackParamAtCluster()->First())));
428 currentParamBack = (*((AliMUONTrackParam*)(fCurrentTrack->GetTrackParamAtCluster()->First())));
429 currentParamForResearchFront = currentParamFront;
430 currentParamForResearchBack = currentParamBack;
431 AliMUONTrackExtrap::AddMCSEffect(&currentParamFront, (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
432 fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeFront())/fRadLengthSi,-1.);
433 AliMUONTrackExtrap::AddMCSEffect(&currentParamForResearchFront,(fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
434 fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeFront())/fRadLengthSi,-1.);
435 AliMUONTrackExtrap::AddMCSEffect(&currentParamBack, (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
436 fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeBack())/fRadLengthSi,-1.);
437 AliMUONTrackExtrap::AddMCSEffect(&currentParamForResearchBack, (fSegmentation->GetPlane(planeId+1)->GetEquivalentSilicon()+
438 fSegmentation->GetPlane(planeId)->GetEquivalentSiliconBeforeBack())/fRadLengthSi,-1.);
439 }
440 // for all planes: extrapolation to the Z of the plane
441 AliMUONTrackExtrap::ExtrapToZCov(&currentParamFront, -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveFront());
442 AliMUONTrackExtrap::ExtrapToZCov(&currentParamForResearchFront, -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveFront());
443 AliMUONTrackExtrap::ExtrapToZCov(&currentParamBack, -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveBack());
444 AliMUONTrackExtrap::ExtrapToZCov(&currentParamForResearchBack, -1.*fSegmentation->GetPlane(planeId)->GetZCenterActiveBack());
445
446 //---------------------------------------------------------------------------------------
447
448 TMatrixD covFront(5,5); covFront = currentParamForResearchFront.GetCovariances();
449 TMatrixD covBack(5,5); covBack = currentParamForResearchBack.GetCovariances();
450
451 Double_t squaredError_X_Front = covFront(0,0);
452 Double_t squaredError_Y_Front = covFront(2,2);
453 Double_t squaredError_X_Back = covBack(0,0);
454 Double_t squaredError_Y_Back = covBack(2,2);
455
456 Double_t corrFact = 1.0;
457
458 Double_t researchRadiusFront = TMath::Sqrt(squaredError_X_Front + squaredError_Y_Front);
459 Double_t researchRadiusBack = TMath::Sqrt(squaredError_X_Back + squaredError_Y_Back);
460 if (0.5*(researchRadiusFront+researchRadiusBack)<fMinResearchRadiusAtPlane[planeId]) {
461 corrFact = fMinResearchRadiusAtPlane[planeId]/(0.5*(researchRadiusFront+researchRadiusBack));
462 }
463
464 //---------------------------------------------------------------------------------------
465
466 Double_t chi2cut = 2.*fScaleSigmaClusterCut*fScaleSigmaClusterCut*fSigmaClusterCut*fSigmaClusterCut; // depends on the number of variables (here, 2)
467
468 // Analyizing the clusters: FRONT ACTIVE ELEMENTS
469
470 Int_t nClustersFront = fMFTClusterArrayFront[planeId]->GetEntries();
471 AliDebug(2, Form("There are %3d clusters in plane %02d FRONT\n", nClustersFront, planeId));
472
473 for (Int_t iCluster=0; iCluster<nClustersFront; iCluster++) {
474
475 Bool_t isGoodChi2 = kFALSE;
476
477 AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArrayFront[planeId]->At(iCluster);
478 Double_t chi2 = (1./(corrFact*corrFact)) * TryOneCluster(currentParamForResearchFront, cluster); // describes the compatibility between the track and the cluster
479 if (chi2<chi2cut) isGoodChi2 = kTRUE;
480
481 if (isGoodChi2) {
482 AliDebug(3, Form("accepting cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
483 AliMuonForwardTrack *newTrack = new ((*fCandidateTracks)[fCandidateTracks->GetEntriesFast()]) AliMuonForwardTrack(*fCurrentTrack);
484 if (fCandidateTracks->GetEntriesFast() > fMaxNCandidates) return kDiverged;
485 newTrack->AddTrackParamAtMFTCluster(currentParamFront, *cluster); // creating new track param and attaching the cluster
486 AliDebug(2, Form("After plane %02d: newTrack->GetNMFTClusters() = %d (fCurrentTrack->GetNMFTClusters() = %d)",
487 planeId, newTrack->GetNMFTClusters(), fCurrentTrack->GetNMFTClusters()));
488 newTrack->SetPlaneExists(planeId);
489 }
490 else AliDebug(3, Form("discarding cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
491
492 }
493
494 // Analyizing the clusters: BACK ACTIVE ELEMENTS
495
496 Int_t nClustersBack = fMFTClusterArrayBack[planeId]->GetEntries();
497 AliDebug(2, Form("There are %3d clusters in plane %02d BACK\n", nClustersBack, planeId));
498
499 for (Int_t iCluster=0; iCluster<nClustersBack; iCluster++) {
500
501 Bool_t isGoodChi2 = kFALSE;
502
503 AliMFTCluster *cluster = (AliMFTCluster*) fMFTClusterArrayBack[planeId]->At(iCluster);
504 Double_t chi2 = (1./(corrFact*corrFact)) * TryOneCluster(currentParamForResearchBack, cluster); // describes the compatibility between the track and the cluster
505 if (chi2<chi2cut) isGoodChi2 = kTRUE;
506
507 if (isGoodChi2) {
508 AliDebug(3,Form("accepting cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
509 AliMuonForwardTrack *newTrack = new ((*fCandidateTracks)[fCandidateTracks->GetEntriesFast()]) AliMuonForwardTrack(*fCurrentTrack);
510 if (fCandidateTracks->GetEntriesFast() > fMaxNCandidates) return kDiverged;
511 newTrack->AddTrackParamAtMFTCluster(currentParamBack, *cluster); // creating new track param and attaching the cluster
512 AliDebug(2, Form("After plane %02d: newTrack->GetNMFTClusters() = %d (fCurrentTrack->GetNMFTClusters() = %d)",
513 planeId, newTrack->GetNMFTClusters(), fCurrentTrack->GetNMFTClusters()));
514 newTrack->SetPlaneExists(planeId);
515 }
516 else AliDebug(3,Form("discarding cluster: chi2=%f (cut = %f)\n", chi2, chi2cut));
517
518 }
519
520 //---------------------------------------------------------------------------------------------
521
522 return kConverged;
523
524}
525
526//==========================================================================================================================================
527
528Double_t AliMFTTrackerMU::TryOneCluster(const AliMUONTrackParam &trackParam, AliMFTCluster *cluster) {
529
530 // Test the compatibility between the track and the cluster (using trackParam's covariance matrix):
531 // return the corresponding Chi2
532 // assume the track parameters are given at the Z of the cluster
533
534 // Set differences between trackParam and cluster in the bending and non bending directions
535 Double_t dX = cluster->GetX() - trackParam.GetNonBendingCoor();
536 Double_t dY = cluster->GetY() - trackParam.GetBendingCoor();
537 AliDebug(3,Form("dX = %f, dY = %f\n", dX, dY));
538
539 // Calculate errors and covariances
540 const TMatrixD& kParamCov = trackParam.GetCovariances();
541 Double_t sigmaX2 = kParamCov(0,0) + cluster->GetErrX2();
542 Double_t sigmaY2 = kParamCov(2,2) + cluster->GetErrY2();
543 AliDebug(3, Form("dX2 = %f, dY2 = %f\n", sigmaX2, sigmaY2));
544 Double_t covXY = kParamCov(0,2);
545 Double_t det = sigmaX2 * sigmaY2 - covXY * covXY;
546
547 // Compute chi2
548 if (det==0.) return 1.e10;
549 return (dX*dX*sigmaY2 + dY*dY*sigmaX2 - 2.*dX*dY*covXY) / det;
550
551}
552
553//=========================================================================================================================================
554
65273893 555Bool_t AliMFTTrackerMU::IsCorrectMatch(AliMFTCluster *cluster, Int_t labelMC) {
556
557 Bool_t result = kFALSE;
558
559 // check if the cluster belongs to the correct MC track
560
561 for (Int_t iTrack=0; iTrack<cluster->GetNMCTracks(); iTrack++) {
562 if (cluster->GetMCLabel(iTrack)==labelMC) {
563 result = kTRUE;
564 break;
565 }
566 }
567
568 return result;
569
570}
571
572//======================================================================================================================================
573
574void AliMFTTrackerMU::GetVertexFromMC() {
575
576 AliRunLoader *fRunLoader = AliRunLoader::Open("galice.root");
577 if (!fRunLoader) {
578 AliError("no run loader found in file galice.root");
579 return;
580 }
581
582 fRunLoader->CdGAFile();
583 fRunLoader->LoadgAlice();
584 fRunLoader->LoadHeader();
585 fRunLoader->GetEvent(gAlice->GetEvNumber());
586
587 TArrayF vtx(3);
588 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(vtx);
589
590 fXExtrapVertex = gRandom->Gaus(vtx[0], AliMFTConstants::fPrimaryVertexResX);
591 fYExtrapVertex = gRandom->Gaus(vtx[1], AliMFTConstants::fPrimaryVertexResY);
592 fZExtrapVertex = gRandom->Gaus(vtx[2], AliMFTConstants::fPrimaryVertexResZ);
593 fXExtrapVertexError = AliMFTConstants::fXVertexTolerance;
594 fYExtrapVertexError = AliMFTConstants::fYVertexTolerance;
595 AliInfo(Form("Set ESD vertex from MC (%f +/- %f, %f +/- %f, %f)",
596 fXExtrapVertex,fXExtrapVertexError,fYExtrapVertex,fYExtrapVertexError,fZExtrapVertex));
597
598}
599
600//======================================================================================================================================
601