1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 /// \class AliMUONRecoCheck
19 /// \brief Utility class to check reconstruction
20 /// Reconstructed tracks are compared to reference tracks.
21 /// The reference tracks are built from AliTrackReference for the
22 /// hit in chamber (0..9) and from kinematics for the vertex parameters.
25 #include "AliMUONRecoCheck.h"
26 #include "AliMUONTrack.h"
27 #include "AliMUONData.h"
28 #include "AliMUONConstants.h"
30 #include "AliRun.h" // for gAlice
31 #include "AliLoader.h"
32 #include "AliRunLoader.h"
33 #include "AliHeader.h"
35 #include "AliTrackReference.h"
38 #include <TParticle.h>
41 ClassImp(AliMUONRecoCheck)
44 //_____________________________________________________________________________
45 AliMUONRecoCheck::AliMUONRecoCheck(AliRunLoader *runloader, AliMUONData *muondata)
51 fReconstructibleTracks(0),
56 fMuonTrackRef = new TClonesArray("AliMUONTrack", 10);
59 fRunLoader = runloader;
61 AliError(Form("no run loader found " ));
68 AliError(Form("no MUONData found " ));
75 //_____________________________________________________________________________
76 AliMUONRecoCheck::~AliMUONRecoCheck()
80 fMuonTrackRef->Delete();
84 //_____________________________________________________________________________
85 void AliMUONRecoCheck::MakeTrackRef()
87 /// Make reconstructible tracks
89 AliTrackReference *trackReference;
90 AliMUONTrack *muonTrack;
91 AliMUONTrackParam *trackParam;
92 AliMUONHitForRec *hitForRec;
93 Float_t x, y, z, pX, pY, pZ, pYZ;
94 Int_t track, trackSave;
96 Float_t bendingSlope = 0;
97 Float_t nonBendingSlope = 0;
98 Float_t inverseBendingMomentum = 0;
100 TTree* treeTR = fRunLoader->TreeTR();
101 if (treeTR == NULL) return;
103 TBranch* branch = treeTR->GetBranch("MUON");
104 if (branch == NULL) return;
106 TClonesArray* trackRefs = 0;
107 branch->SetAddress(&trackRefs);
108 branch->SetAutoDelete(kTRUE);
110 Int_t nTrackRef = (Int_t)branch->GetEntries();
112 track = trackSave = -999;
114 Int_t iHitMin, iChamber, detElemId;
116 trackParam = new AliMUONTrackParam();
117 hitForRec = new AliMUONHitForRec();
118 muonTrack = new AliMUONTrack();
120 Int_t max = fRunLoader->GetHeader()->Stack()->GetNtrack();
121 for (Int_t iTrackRef = 0; iTrackRef < nTrackRef; iTrackRef++) {
123 branch->GetEntry(iTrackRef);
128 if (!trackRefs->GetEntries()) continue;
131 for (Int_t iHit = iHitMin; iHit < trackRefs->GetEntries(); iHit++) {
133 trackReference = (AliTrackReference*)trackRefs->At(iHit);
134 x = trackReference->X();
135 y = trackReference->Y();
136 z = trackReference->Z();
137 pX = trackReference->Px();
138 pY = trackReference->Py();
139 pZ = trackReference->Pz();
141 track = trackReference->GetTrack();
145 << "Track ID " << track
146 << " larger than max number of particles " << max << endl;
150 if (track != trackSave && iHit != 0) {
156 // track parameters at hit
157 trackParam->SetBendingCoor(y);
158 trackParam->SetNonBendingCoor(x);
161 if (TMath::Abs(pZ) > 0) {
162 bendingSlope = pY/pZ;
163 nonBendingSlope = pX/pZ;
165 pYZ = TMath::Sqrt(pY*pY+pZ*pZ);
166 if (pYZ >0) inverseBendingMomentum = 1/pYZ;
168 trackParam->SetBendingSlope(bendingSlope);
169 trackParam->SetNonBendingSlope(nonBendingSlope);
170 trackParam->SetInverseBendingMomentum(inverseBendingMomentum);
172 hitForRec->SetBendingCoor(y);
173 hitForRec->SetNonBendingCoor(x);
175 hitForRec->SetBendingReso2(0.0);
176 hitForRec->SetNonBendingReso2(0.0);
177 detElemId = hitForRec->GetDetElemId();
178 if (detElemId) iChamber = detElemId / 100 - 1;
179 else iChamber = AliMUONConstants::ChamberNumber(z);
180 hitForRec->SetChamberNumber(iChamber);
182 muonTrack->AddTrackParamAtHit(trackParam,hitForRec);
183 muonTrack->AddHitForRecAtHit(hitForRec);
184 muonTrack->SetTrackID(track);
187 if (iHit == trackRefs->GetEntries()-1) isNewTrack = kFALSE;
190 // track parameters at vertex
191 particle = fRunLoader->GetHeader()->Stack()->Particle(muonTrack->GetTrackID());
202 trackParam->SetBendingCoor(y);
203 trackParam->SetNonBendingCoor(x);
205 if (TMath::Abs(pZ) > 0) {
206 bendingSlope = pY/pZ;
207 nonBendingSlope = pX/pZ;
209 pYZ = TMath::Sqrt(pY*pY+pZ*pZ);
210 if (pYZ >0) inverseBendingMomentum = 1/pYZ;
211 trackParam->SetBendingSlope(bendingSlope);
212 trackParam->SetNonBendingSlope(nonBendingSlope);
213 trackParam->SetInverseBendingMomentum(inverseBendingMomentum);
215 muonTrack->SetTrackParamAtVertex(trackParam);
218 AddMuonTrackReference(muonTrack);
219 muonTrack->ResetTrackParamAtHit();
220 muonTrack->ResetHitForRecAtHit();
222 } // end while isNewTrack
227 ReconstructibleTracks();
235 //____________________________________________________________________________
236 TClonesArray* AliMUONRecoCheck::GetTrackReco()
238 /// Return TClonesArray of reconstructed tracks
240 fMUONData->ResetRecTracks();
241 fMUONData->SetTreeAddress("RT");
242 fTrackReco = fMUONData->RecTracks();
243 fMUONData->GetRecTracks();
244 fRecoTracks = fTrackReco->GetEntriesFast();
247 //_____________________________________________________________________________
248 void AliMUONRecoCheck::PrintEvent() const
253 AliMUONHitForRec *hitForRec;
254 TClonesArray * hitForRecAtHit = 0;
255 Float_t xRec,yRec,zRec;
257 Int_t nTrackRef = fMuonTrackRef->GetEntriesFast();
259 printf(" ******************************* \n");
260 printf(" nb of tracks %d \n",nTrackRef);
262 for (Int_t index = 0; index < nTrackRef; index++) {
263 track = (AliMUONTrack*)fMuonTrackRef->At(index);
264 hitForRecAtHit = track->GetHitForRecAtHit();
265 Int_t nTrackHits = hitForRecAtHit->GetEntriesFast();
266 printf(" track number %d \n",index);
267 for (Int_t iHit = 0; iHit < nTrackHits; iHit++){
268 hitForRec = (AliMUONHitForRec*) hitForRecAtHit->At(iHit);
269 xRec = hitForRec->GetNonBendingCoor();
270 yRec = hitForRec->GetBendingCoor();
271 zRec = hitForRec->GetZ();
272 printf(" x,y,z: %f , %f , %f \n",xRec,yRec,zRec);
277 //_____________________________________________________________________________
278 void AliMUONRecoCheck::ResetTracks() const
282 if (fMuonTrackRef) fMuonTrackRef->Clear();
284 //_____________________________________________________________________________
285 void AliMUONRecoCheck::CleanMuonTrackRef()
287 /// Re-calculate hits parameters because two AliTrackReferences are recorded for
288 /// each chamber (one when particle is entering + one when particle is leaving
289 /// the sensitive volume)
291 Float_t maxGasGap = 1.; // cm
292 AliMUONTrack *track, *trackNew;
293 AliMUONHitForRec *hitForRec, *hitForRec1, *hitForRec2;
294 AliMUONTrackParam *trackParam, *trackParam1, *trackParam2, *trackParamAtVertex;
295 TClonesArray * hitForRecAtHit = 0;
296 TClonesArray * trackParamAtHit = 0;
297 Float_t xRec,yRec,zRec;
298 Float_t xRec1,yRec1,zRec1;
299 Float_t xRec2,yRec2,zRec2;
300 Float_t bendingSlope,nonBendingSlope,bendingMomentum;
301 Float_t bendingSlope1,nonBendingSlope1,bendingMomentum1;
302 Float_t bendingSlope2,nonBendingSlope2,bendingMomentum2;
303 TClonesArray *newMuonTrackRef = new TClonesArray("AliMUONTrack", 10);
305 Int_t iChamber = 0, detElemId = 0;
307 Int_t nTrackHits = 0;
309 hitForRec = new AliMUONHitForRec();
310 trackParam = new AliMUONTrackParam();
311 trackNew = new AliMUONTrack();
313 Int_t nTrackRef = fMuonTrackRef->GetEntriesFast();
314 for (Int_t index = 0; index < nTrackRef; index++) {
315 track = (AliMUONTrack*)fMuonTrackRef->At(index);
316 hitForRecAtHit = track->GetHitForRecAtHit();
317 trackParamAtHit = track->GetTrackParamAtHit();
318 trackParamAtVertex = track->GetTrackParamAtVertex();
319 nTrackHits = hitForRecAtHit->GetEntriesFast();
321 while (iHit1 < nTrackHits) {
322 hitForRec1 = (AliMUONHitForRec*) hitForRecAtHit->At(iHit1);
323 trackParam1 = (AliMUONTrackParam*) trackParamAtHit->At(iHit1);
324 xRec1 = hitForRec1->GetNonBendingCoor();
325 yRec1 = hitForRec1->GetBendingCoor();
326 zRec1 = hitForRec1->GetZ();
330 bendingSlope1 = trackParam1->GetBendingSlope();
331 nonBendingSlope1 = trackParam1->GetNonBendingSlope();
332 bendingMomentum1 = 0;
333 if (TMath::Abs(trackParam1->GetInverseBendingMomentum()) > 0)
334 bendingMomentum1 = 1./trackParam1->GetInverseBendingMomentum();
335 bendingSlope = bendingSlope1;
336 nonBendingSlope = nonBendingSlope1;
337 bendingMomentum = bendingMomentum1;
339 for (Int_t iHit2 = iHit1+1; iHit2 < nTrackHits; iHit2++) {
340 hitForRec2 = (AliMUONHitForRec*) hitForRecAtHit->At(iHit2);
341 trackParam2 = (AliMUONTrackParam*) trackParamAtHit->At(iHit2);
342 xRec2 = hitForRec2->GetNonBendingCoor();
343 yRec2 = hitForRec2->GetBendingCoor();
344 zRec2 = hitForRec2->GetZ();
345 bendingSlope2 = trackParam2->GetBendingSlope();
346 nonBendingSlope2 = trackParam2->GetNonBendingSlope();
347 bendingMomentum2 = 0;
348 if (TMath::Abs(trackParam2->GetInverseBendingMomentum()) > 0)
349 bendingMomentum2 = 1./trackParam2->GetInverseBendingMomentum();
351 if ( TMath::Abs(zRec2-zRec1) < maxGasGap ) {
357 bendingSlope += bendingSlope2;
358 nonBendingSlope += nonBendingSlope2;
359 bendingMomentum += bendingMomentum2;
364 xRec /= (Float_t)nRec;
365 yRec /= (Float_t)nRec;
366 zRec /= (Float_t)nRec;
367 bendingSlope /= (Float_t)nRec;
368 nonBendingSlope /= (Float_t)nRec;
369 bendingMomentum /= (Float_t)nRec;
371 hitForRec->SetNonBendingCoor(xRec);
372 hitForRec->SetBendingCoor(yRec);
373 hitForRec->SetZ(zRec);
374 detElemId = hitForRec->GetDetElemId();
375 if (detElemId) iChamber = detElemId / 100 - 1;
376 else iChamber = AliMUONConstants::ChamberNumber(zRec);
377 hitForRec->SetChamberNumber(iChamber);
378 hitForRec->SetBendingReso2(0.0);
379 hitForRec->SetNonBendingReso2(0.0);
380 trackParam->SetNonBendingCoor(xRec);
381 trackParam->SetBendingCoor(yRec);
382 trackParam->SetZ(zRec);
383 trackParam->SetNonBendingSlope(nonBendingSlope);
384 trackParam->SetBendingSlope(bendingSlope);
385 if (TMath::Abs(bendingMomentum) > 0)
386 trackParam->SetInverseBendingMomentum(1./bendingMomentum);
388 trackNew->AddHitForRecAtHit(hitForRec);
389 trackNew->AddTrackParamAtHit(trackParam,hitForRec);
394 trackNew->SetTrackID(track->GetTrackID());
395 trackNew->SetTrackParamAtVertex(trackParamAtVertex);
396 {new ((*newMuonTrackRef)[newMuonTrackRef->GetEntriesFast()]) AliMUONTrack(*trackNew);}
397 trackNew->ResetHitForRecAtHit();
398 trackNew->ResetTrackParamAtHit();
402 fMuonTrackRef->Clear();
403 nTrackRef = newMuonTrackRef->GetEntriesFast();
404 for (Int_t index = 0; index < nTrackRef; index++) {
405 track = (AliMUONTrack*)newMuonTrackRef->At(index);
406 AddMuonTrackReference(track);
413 newMuonTrackRef->Delete();
414 delete newMuonTrackRef;
418 //_____________________________________________________________________________
419 void AliMUONRecoCheck::ReconstructibleTracks()
421 /// Calculate the number of reconstructible tracks
424 TClonesArray* hitForRecAtHit = NULL;
425 AliMUONHitForRec* hitForRec;
427 Int_t nTrackRef, nTrackHits;
428 Int_t isChamberInTrack[10];
430 Bool_t isTrackOK = kTRUE;
432 fReconstructibleTracks = 0;
434 nTrackRef = fMuonTrackRef->GetEntriesFast();
435 for (Int_t index = 0; index < nTrackRef; index++) {
436 track = (AliMUONTrack*)fMuonTrackRef->At(index);
437 hitForRecAtHit = track->GetHitForRecAtHit();
438 nTrackHits = hitForRecAtHit->GetEntriesFast();
439 for (Int_t ch = 0; ch < 10; ch++) isChamberInTrack[ch] = 0;
441 for ( Int_t iHit = 0; iHit < nTrackHits; iHit++) {
442 hitForRec = (AliMUONHitForRec*) hitForRecAtHit->At(iHit);
443 zRec = hitForRec->GetZ();
444 iChamber = hitForRec->GetChamberNumber();
445 if (iChamber < 0 || iChamber > 10) continue;
446 isChamberInTrack[iChamber] = 1;
448 // track is reconstructible if the particle is depositing a hit
449 // in the following chamber combinations:
452 if (!isChamberInTrack[0] && !isChamberInTrack[1]) isTrackOK = kFALSE;
453 if (!isChamberInTrack[2] && !isChamberInTrack[3]) isTrackOK = kFALSE;
454 if (!isChamberInTrack[4] && !isChamberInTrack[5]) isTrackOK = kFALSE;
455 Int_t nHitsInLastStations=0;
456 for (Int_t ch = 6; ch < AliMUONConstants::NTrackingCh(); ch++)
457 if (isChamberInTrack[ch]) nHitsInLastStations++;
458 if(nHitsInLastStations < 3) isTrackOK = kFALSE;
460 if (isTrackOK) fReconstructibleTracks++;
461 if (!isTrackOK) fMuonTrackRef->Remove(track); // remove non reconstructible tracks
463 fMuonTrackRef->Compress();