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