]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - MFT/AliMFTTrackerMU.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
Analysis code updated
[u/mrichter/AliRoot.git] / MFT / AliMFTTrackerMU.cxx
... / ...
CommitLineData
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//
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"
32#include "AliESDMuonGlobalTrack.h"
33#include "AliMFTTrackerMU.h"
34#include "TMath.h"
35#include "AliRun.h"
36#include "AliRunLoader.h"
37#include "AliHeader.h"
38#include "AliGenEventHeader.h"
39#include "TArrayF.h"
40#include "AliMFT.h"
41#include "AliMUONTrackExtrap.h"
42#include "AliMUONTrack.h"
43#include "AliMUONESDInterface.h"
44#include "AliMuonForwardTrack.h"
45#include "AliMUONConstants.h"
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),
60 fSigmaClusterCut(2),
61 fScaleSigmaClusterCut(1.),
62 fNMaxMissingMFTClusters(0),
63 fGlobalTrackingDiverged(kFALSE),
64 fCandidateTracks(0),
65 fMUONTrack(0),
66 fCurrentTrack(0),
67 fFinalBestCandidate(0),
68 fXExtrapVertex(0),
69 fYExtrapVertex(0),
70 fZExtrapVertex(0),
71 fXExtrapVertexError(0),
72 fYExtrapVertexError(0),
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
85 AliInfo(Form("fMFT = %p, fSegmentation = %p", fMFT, fSegmentation));
86
87 for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
88 fMFTClusterArray[iPlane] = new TClonesArray("AliMFTCluster");
89 fMFTClusterArrayFront[iPlane] = new TClonesArray("AliMFTCluster");
90 fMFTClusterArrayBack[iPlane] = new TClonesArray("AliMFTCluster");
91 fMFTClusterArray[iPlane] -> SetOwner(kTRUE);
92 fMFTClusterArrayFront[iPlane] -> SetOwner(kTRUE);
93 fMFTClusterArrayBack[iPlane] -> SetOwner(kTRUE);
94 fMinResearchRadiusAtPlane[iPlane] = 0.;
95 }
96
97 fCandidateTracks = new TClonesArray("AliMuonForwardTrack",50000);
98
99}
100
101//====================================================================================================================================================
102
103AliMFTTrackerMU::~AliMFTTrackerMU() {
104
105 // destructor
106
107 for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
108 fMFTClusterArray[iPlane] -> Delete();
109 delete fMFTClusterArray[iPlane];
110 delete fMFTClusterArrayFront[iPlane];
111 delete fMFTClusterArrayBack[iPlane];
112 }
113
114 delete fCandidateTracks;
115
116}
117
118//====================================================================================================================================================
119
120Int_t AliMFTTrackerMU::LoadClusters(TTree *cTree) {
121
122 //--------------------------------------------------------------------
123 // This function loads the MFT clusters
124 //--------------------------------------------------------------------
125
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
131 if (!cTree->GetEvent()) return kFALSE;
132 for (Int_t iPlane=0; iPlane<fNPlanesMFT; iPlane++) {
133 AliInfo(Form("plane %02d: nClusters = %d\n", iPlane, fMFTClusterArray[iPlane]->GetEntries()));
134 }
135 SeparateFrontBackClusters();
136
137 return 0;
138
139}
140
141//====================================================================================================================================================
142
143void AliMFTTrackerMU::UnloadClusters() {
144
145 //--------------------------------------------------------------------
146 // This function unloads MFT clusters
147 //--------------------------------------------------------------------
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
166 // Tree of AliMuonForwardTrack objects. Created outside the ESD framework for cross-check purposes
167
168 TFile *outputFileMuonGlobalTracks = new TFile("MuonGlobalTracks.root", "update");
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
177 fXExtrapVertex = 0;
178 fYExtrapVertex = 0;
179 fZExtrapVertex = 0;
180 fXExtrapVertexError = AliMFTConstants::fXVertexTolerance;
181 fYExtrapVertexError = AliMFTConstants::fYVertexTolerance;
182
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
185 const AliESDVertex* esdVert = fESD->GetVertex();
186 if (esdVert->GetNContributors() > 0 || !strcmp(esdVert->GetTitle(),"vertexer: smearMC")) {
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));
194 }
195 else GetVertexFromMC();
196
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
208 AliInfo("****************************************************************************************");
209 AliInfo(Form("*************************** MUON TRACK %3d/%d ***************************************", iTrack, nTracksMUON));
210 AliInfo("****************************************************************************************");
211
212 fCandidateTracks -> Delete();
213
214 fNPlanesMFTAnalyzed = 0;
215
216 const AliESDMuonTrack *esdTrack = event->GetMuonTrack(iTrack);
217 if (fMUONTrack) delete fMUONTrack;
218 fMUONTrack = new AliMUONTrack();
219
220 AliMUONESDInterface::ESDToMUON(*esdTrack, *fMUONTrack, kFALSE);
221
222 if (!fMUONTrack->GetTrackParamAtCluster()->First()) {
223 AliInfo("Skipping track, no parameters available!!!");
224 iTrack++;
225 continue;
226 }
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
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
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);
255
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 }
275
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();
288 AliInfo(Form("nFinalTracks = %d", nFinalTracks));
289
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
320 new ((*muonForwardTracks)[muonForwardTracks->GetEntries()]) AliMuonForwardTrack(*newTrack);
321
322 //----------------------- Save the information to the AliESDMuonGlobalTrack object
323
324 newTrack -> EvalKinem(AliMFTConstants::fZEvalKinem);
325
326 AliDebug(2,"Creating a new Muon Global Track");
327 AliESDMuonGlobalTrack *myESDTrack = event->NewMuonGlobalTrack();
328 myESDTrack -> SetPxPyPz(newTrack->Px(), newTrack->Py(), newTrack->Pz());
329
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());
337 myESDTrack -> SetXYAtVertex(newTrack->GetOffsetX(0., AliMFTConstants::fZEvalKinem), newTrack->GetOffsetY(0., AliMFTConstants::fZEvalKinem));
338 myESDTrack -> SetRAtAbsorberEnd(newTrack->GetRAtAbsorberEnd());
339 myESDTrack -> SetCovariances(newTrack->GetTrackParamAtMFTCluster(0)->GetCovariances());
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());
345
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
354 //---------------------------------------------------------------------------------
355
356 }
357
358 fFinalBestCandidate = NULL;
359
360 }
361
362 outputTreeMuonGlobalTracks -> Fill();
363
364 Int_t myEventID = 0;
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
371 muonForwardTracks -> Delete();
372 delete muonForwardTracks;
373
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
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;
415 if (fBransonCorrection) {
416 AliMUONTrackExtrap::ExtrapToVertex(&currentParamFront, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
417 AliMUONTrackExtrap::ExtrapToVertex(&currentParamBack, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
418 }
419 else {
420 AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&currentParamFront, fZExtrapVertex);
421 AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&currentParamBack, fZExtrapVertex);
422 }
423 AliMUONTrackExtrap::ExtrapToVertex(&currentParamForResearchFront, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
424 AliMUONTrackExtrap::ExtrapToVertex(&currentParamForResearchBack, fXExtrapVertex, fYExtrapVertex, fZExtrapVertex, fXExtrapVertexError, fYExtrapVertexError);
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
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