DQM configure file
[u/mrichter/AliRoot.git] / MUON / AliMUONRecoCheck.cxx
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 /* $Id$ */
17
18 //-----------------------------------------------------------------------------
19 /// \class AliMUONRecoCheck
20 /// Utility class to check reconstruction
21 /// Reconstructed tracks are compared to reference tracks. 
22 /// The reference tracks are built from AliTrackReference for the
23 /// hit in chamber (0..9) and from kinematics for the vertex parameters.     
24 //-----------------------------------------------------------------------------
25
26 // 13 Nov 2007:
27 // Added a method to create a list of reconstructed AliMUONTrack objects from
28 // ESD data. This is necessary since the track objects that are actually created
29 // during offline reconstruction are no longer stored to disk.
30 //  - Artur Szostak <artursz@iafrica.com>
31 // 25 Jan 2008:
32 // Use the new ESDInterface to create MUON objects from ESD data
33 // - Philippe Pillot
34 // 
35
36 #include "AliMUONRecoCheck.h"
37 #include "AliMUONTrack.h"
38 #include "AliMUONVTrackStore.h"
39 #include "AliMUONVCluster.h"
40 #include "AliMUONVClusterStore.h"
41 #include "AliMUONConstants.h"
42 #include "AliMUONESDInterface.h"
43 #include "AliMUONTrackParam.h"
44 #include "AliMUONTriggerTrack.h"
45 #include "AliMUONVTriggerTrackStore.h"
46 #include "AliMUONTriggerStoreV1.h"
47 #include "AliMUONLocalTrigger.h"
48 #include "AliMCEventHandler.h"
49 #include "AliMCEvent.h"
50 #include "AliStack.h"
51 #include "AliTrackReference.h"
52 #include "AliLog.h"
53 #include "AliESDEvent.h"
54 #include "AliESDMuonTrack.h"
55 #include "AliMUONDigitStoreV2S.h"
56 #include "AliMUONVDigit.h"
57 #include "AliMpVSegmentation.h"
58 #include "AliMpSegmentation.h"
59 #include "AliMpPad.h"
60
61 #include "AliGeomManager.h"
62 #include "AliCDBManager.h"
63 #include "AliMpCDB.h"
64 #include "AliMpDDLStore.h"
65 #include "AliMUONCDB.h"
66 #include "AliMUONGeometryTransformer.h"
67 #include "AliMUONTriggerCircuit.h"
68 #include "AliMUONVTrackReconstructor.h"
69 #include "AliMUONVTriggerStore.h"
70 #include "AliMUONCalibrationData.h"
71 #include "AliMUONTriggerElectronics.h"
72
73
74 #include "TGeoManager.h"
75
76 #include <TFile.h>
77 #include <TTree.h>
78 #include <TParticle.h>
79 #include <TParticlePDG.h>
80 #include <Riostream.h>
81
82 #include "AliMUONRecoCheck.h"
83
84 /// \cond CLASSIMP
85 ClassImp(AliMUONRecoCheck)
86 /// \endcond
87
88 //_____________________________________________________________________________
89 AliMUONRecoCheck::AliMUONRecoCheck(const Char_t *esdFileName, const Char_t *pathSim)
90 : TObject(),
91 fMCEventHandler(new AliMCEventHandler()),
92 fESDEvent(new AliESDEvent()),
93 fESDTree (0x0),
94 fESDFile (0x0),
95 fCurrentEvent(0),
96 fTrackRefStore(0x0),
97 fRecoTrackRefStore(0x0),
98 fRecoTriggerRefStore(0x0),
99 fRecoTrackStore(0x0),
100 fRecoTriggerTrackStore(0x0),
101 fGeometryTransformer(0x0),
102 fTriggerCircuit(0x0),
103 fCalibrationData(0x0),
104 fTriggerElectronics(0x0),
105 fESDEventOwner(kTRUE)
106 {
107   /// Normal ctor
108   
109   // TrackRefs and Particules
110   fMCEventHandler->SetInputPath(pathSim);
111   fMCEventHandler->InitIO("");
112   
113   // ESD MUON Tracks
114   fESDFile = TFile::Open(esdFileName); // open the file
115   if (!fESDFile || !fESDFile->IsOpen()) {
116     AliError(Form("opening ESD file %s failed", esdFileName));
117     fESDFile = 0x0;
118     return;
119   }
120   fESDTree = (TTree*) fESDFile->Get("esdTree"); // get the tree
121   if (!fESDTree) {
122     AliError("no ESD tree found");
123     fESDFile->Close();
124     fESDFile = 0x0;
125     return;
126   }
127   fESDEvent->ReadFromTree(fESDTree); // link fESDEvent to the tree
128 }
129
130 //_____________________________________________________________________________
131 AliMUONRecoCheck::AliMUONRecoCheck(AliESDEvent *esdEvent, AliMCEventHandler *mcEventHandler)
132 : TObject(),
133 fMCEventHandler(0),
134 fESDEvent(0),
135 fESDTree (0x0),
136 fESDFile (0x0),
137 fCurrentEvent(0),
138 fTrackRefStore(0x0),
139 fRecoTrackRefStore(0x0),
140 fRecoTriggerRefStore(0x0),
141 fRecoTrackStore(0x0),
142 fRecoTriggerTrackStore(0x0),
143 fGeometryTransformer(0x0),
144 fTriggerCircuit(0x0),
145 fCalibrationData(0x0),
146 fTriggerElectronics(0x0),
147 fESDEventOwner(kFALSE)
148 {
149   /// Normal ctor
150   
151   // TrackRefs and Particules
152   fMCEventHandler = mcEventHandler;
153   
154   // ESD MUON Tracks
155   fESDEvent = esdEvent;
156 }
157
158 //_____________________________________________________________________________
159 AliMUONRecoCheck::~AliMUONRecoCheck()
160 {
161   /// Destructor
162   if (fESDEventOwner) {
163     delete fMCEventHandler;
164     delete fESDEvent;
165     if (fESDFile) fESDFile->Close();
166   }
167   ResetStores();
168   delete fGeometryTransformer;
169   delete fTriggerCircuit;
170   delete fTriggerElectronics;
171   delete fCalibrationData;
172 }
173
174 //_____________________________________________________________________________
175 void AliMUONRecoCheck::ResetStores()
176 {
177   /// Deletes all the store objects that have been created and resets the pointers to 0x0
178   delete fTrackRefStore;      fTrackRefStore = 0x0;
179   delete fRecoTrackRefStore;  fRecoTrackRefStore = 0x0;
180   delete fRecoTriggerRefStore;  fRecoTriggerRefStore = 0x0;
181   delete fRecoTrackStore;     fRecoTrackStore = 0x0;
182   delete fRecoTriggerTrackStore; fRecoTriggerTrackStore = 0x0;
183 }
184
185 //_____________________________________________________________________________
186 Bool_t AliMUONRecoCheck::InitCircuit()
187 {
188   
189   if ( fTriggerCircuit ) return kTRUE;
190   
191   if ( ! InitGeometryTransformer() ) return kFALSE;
192   
193   // reset tracker for local trigger to trigger track conversion
194   if ( ! AliMUONESDInterface::GetTracker() ) AliMUONESDInterface::ResetTracker();
195   
196   fTriggerCircuit = new AliMUONTriggerCircuit(fGeometryTransformer);
197   
198   return kTRUE;
199 }
200
201
202 //_____________________________________________________________________________
203 Int_t AliMUONRecoCheck::GetRunNumber()
204 {
205   /// Return the run number of the current ESD event
206   
207   if (fESDEventOwner && fRecoTrackStore == 0x0 && fRecoTriggerTrackStore == 0x0) {
208     if (!fESDTree || fESDTree->GetEvent(fCurrentEvent) <= 0) {
209       AliError(Form("fails to read ESD object for event %d: cannot get the run number",fCurrentEvent));
210       return -1;
211     }
212   }
213   
214   return fESDEvent->GetRunNumber();
215 }
216
217 //_____________________________________________________________________________
218 Int_t AliMUONRecoCheck::NumberOfEvents() const
219 {
220   /// Return the number of events
221   if (fESDEventOwner && fESDTree) return fESDTree->GetEntries();
222   return 0;
223 }
224
225 //_____________________________________________________________________________
226 AliMUONVTrackStore* AliMUONRecoCheck::ReconstructedTracks(Int_t event, Bool_t refit)
227 {
228   /// Return a track store containing the reconstructed tracks (converted into 
229   /// MUONTrack objects) for a given event.
230   /// Track parameters at each clusters are computed or not depending on the flag "refit".
231   /// If not, only the track parameters at first cluster are valid.
232   
233   if (!fESDEventOwner) {
234     if (fRecoTrackStore == 0x0) MakeReconstructedTracks(refit);
235     return fRecoTrackStore;
236   }
237
238   if (event != fCurrentEvent) {
239     ResetStores();
240     fCurrentEvent = event;
241   }
242   
243   if (fRecoTrackStore != 0x0) return fRecoTrackStore;
244   else {
245     if (!fESDTree) return 0x0;
246     if (fESDTree->GetEvent(event) <= 0) {
247       AliError(Form("fails to read ESD object for event %d", event));
248       return 0x0;
249     }
250     MakeReconstructedTracks(refit);
251     return fRecoTrackStore;
252   }
253 }
254
255
256 //_____________________________________________________________________________
257 AliMUONVTriggerTrackStore* AliMUONRecoCheck::TriggeredTracks(Int_t event)
258 {
259   /// Return a track store containing the reconstructed trigger tracks (converted into 
260   /// MUONTriggerTrack objects) for a given event.
261         
262   if (!fESDEventOwner) {
263     if (fRecoTriggerTrackStore == 0x0) MakeTriggeredTracks();
264     return fRecoTriggerTrackStore;
265   }
266         
267   if (event != fCurrentEvent) {
268     ResetStores();
269     fCurrentEvent = event;
270   }
271         
272   if (fRecoTriggerTrackStore != 0x0) return fRecoTriggerTrackStore;
273   else {
274     if (!fESDTree) return 0x0;
275     if (fESDTree->GetEvent(event) <= 0) {
276       AliError(Form("fails to read ESD object for event %d", event));
277       return 0x0;
278     }
279     MakeTriggeredTracks();
280     return fRecoTriggerTrackStore;
281   }
282 }
283
284
285 //_____________________________________________________________________________
286 AliMUONVTrackStore* AliMUONRecoCheck::TrackRefs(Int_t event)
287 {
288   /// Return a track store containing the track references (converted into 
289   /// MUONTrack objects) for a given event
290   
291   if (!fESDEventOwner) {
292     if (fTrackRefStore == 0x0) MakeTrackRefs();
293     return fTrackRefStore;
294   }
295
296   if (event != fCurrentEvent) {
297     ResetStores();
298     fCurrentEvent = event;
299   }
300   
301   if (fTrackRefStore != 0x0) return fTrackRefStore;
302   else {
303     if (!fMCEventHandler->LoadEvent(event)) {
304       AliError(Form("fails to read MC objects for event %d", event));
305       return 0x0;
306     }
307     MakeTrackRefs();
308     return fTrackRefStore;
309   }
310 }
311
312 //_____________________________________________________________________________
313 AliMUONVTriggerTrackStore* AliMUONRecoCheck::TriggerableTracks(Int_t event)
314 {
315   /// Return a trigger track store containing the triggerable track references (converted into 
316   /// AliMUONTriggerTrack objects) for a given event
317         
318   if (!fESDEventOwner) {
319     if (fRecoTriggerRefStore == 0x0) MakeTriggerableTracks();
320     return fRecoTriggerRefStore;
321   }
322         
323   if (event != fCurrentEvent) {
324     ResetStores();
325     fCurrentEvent = event;
326   }
327         
328   if (fRecoTriggerRefStore != 0x0) return fRecoTriggerRefStore;
329   else {
330     if (!fMCEventHandler->LoadEvent(event)) {
331       AliError(Form("fails to read MC objects for event %d", event));
332       return 0x0;
333     }
334     MakeTriggerableTracks();
335     return fRecoTriggerRefStore;
336   }
337 }
338
339
340 //_____________________________________________________________________________
341 AliMUONVTrackStore* AliMUONRecoCheck::ReconstructibleTracks(Int_t event, UInt_t requestedStationMask,
342                                                             Bool_t request2ChInSameSt45, Bool_t hitInFrontOfPad)
343 {
344   /// Return a track store containing the reconstructible tracks for a given event,
345   /// according to the mask of requested stations and the minimum number of chambers hit in stations 4 & 5.
346
347   if (!fESDEventOwner) {
348     if (fRecoTrackRefStore == 0x0) {
349       if (TrackRefs(event) == 0x0) return 0x0;
350       MakeReconstructibleTracks(requestedStationMask, request2ChInSameSt45, hitInFrontOfPad);
351     }
352     return fRecoTrackRefStore;
353   }
354
355   if (event != fCurrentEvent) {
356     ResetStores();
357     fCurrentEvent = event;
358   }
359   
360   if (fRecoTrackRefStore != 0x0) return fRecoTrackRefStore;
361   else {
362     if (TrackRefs(event) == 0x0) return 0x0;
363     MakeReconstructibleTracks(requestedStationMask, request2ChInSameSt45, hitInFrontOfPad);
364     return fRecoTrackRefStore;
365   }
366 }
367
368
369 //_____________________________________________________________________________
370 void AliMUONRecoCheck::MakeReconstructedTracks(Bool_t refit)
371 {
372   /// Make reconstructed tracks
373   if (!(fRecoTrackStore = AliMUONESDInterface::NewTrackStore())) return;
374   
375   // loop over all reconstructed tracks and add them to the store (skip ghosts)
376   Int_t nTracks = (Int_t) fESDEvent->GetNumberOfMuonTracks();
377   for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
378     AliESDMuonTrack* esdTrack = fESDEvent->GetMuonTrack(iTrack);
379     if (esdTrack->ContainTrackerData()) AliMUONESDInterface::Add(*esdTrack, *fRecoTrackStore, refit);
380   }
381   
382 }
383
384
385 //_____________________________________________________________________________
386 void AliMUONRecoCheck::MakeTriggeredTracks()
387 {
388   /// Make reconstructed trigger tracks
389   if (!(fRecoTriggerTrackStore = AliMUONESDInterface::NewTriggerTrackStore())) return;
390         
391   AliMUONVTriggerStore* tmpTriggerStore = AliMUONESDInterface::NewTriggerStore();
392   if ( ! tmpTriggerStore ) return;
393   
394   // loop over all reconstructed tracks and add them to the store (include ghosts)
395   Int_t nTracks = (Int_t) fESDEvent->GetNumberOfMuonTracks();
396   for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
397     AliESDMuonTrack* esdTrack = fESDEvent->GetMuonTrack(iTrack);
398     if (esdTrack->ContainTriggerData()) AliMUONESDInterface::Add(*esdTrack, *tmpTriggerStore);
399   }
400         
401   if ( InitCircuit() ) {
402   
403     AliMUONVTrackReconstructor* tracker = AliMUONESDInterface::GetTracker();
404     tracker->EventReconstructTrigger(*fTriggerCircuit, *tmpTriggerStore, *fRecoTriggerTrackStore);
405   }
406   
407   delete tmpTriggerStore;
408 }
409
410 //_____________________________________________________________________________
411 void AliMUONRecoCheck::TriggerToTrack(const AliMUONLocalTrigger& locTrg, AliMUONTriggerTrack& triggerTrack)
412 {
413   /// Make trigger track from local trigger info
414   if ( ! InitCircuit() ) return;
415   AliMUONVTrackReconstructor* tracker = AliMUONESDInterface::GetTracker();
416   tracker->TriggerToTrack(*fTriggerCircuit, locTrg, triggerTrack);
417 }
418
419
420 //_____________________________________________________________________________
421 void AliMUONRecoCheck::MakeTrackRefs()
422 {
423   /// Make reconstructible tracks
424   AliMUONVTrackStore *tmpTrackRefStore = AliMUONESDInterface::NewTrackStore();
425   if (!tmpTrackRefStore) return;
426   
427   Double_t x, y, z, pX, pY, pZ, bendingSlope, nonBendingSlope, inverseBendingMomentum;
428   TParticle* particle;
429   TClonesArray* trackRefs;
430   Int_t nTrackRef = fMCEventHandler->MCEvent()->GetNumberOfTracks();
431   AliMUONVCluster* hit = AliMUONESDInterface::NewCluster();
432   
433   // loop over simulated tracks
434   for (Int_t iTrackRef  = 0; iTrackRef < nTrackRef; ++iTrackRef) {
435     Int_t nHits = fMCEventHandler->GetParticleAndTR(iTrackRef, particle, trackRefs);
436     
437     // skip empty trackRefs
438     if (nHits < 1) continue;
439     
440     // get the particle charge for further calculation
441     TParticlePDG* ppdg = particle->GetPDG();
442     Int_t charge = ppdg != NULL ? (Int_t)(ppdg->Charge()/3.0) : 0;
443     
444     AliMUONTrack track;
445     
446     // loop over simulated track hits
447     for (Int_t iHit = 0; iHit < nHits; ++iHit) {        
448       AliTrackReference* trackReference = static_cast<AliTrackReference*>(trackRefs->UncheckedAt(iHit));
449       
450       // skip trackRefs not in MUON
451       if (trackReference->DetectorId() != AliTrackReference::kMUON) continue;
452     
453       // Get track parameters of current hit
454       x = trackReference->X();
455       y = trackReference->Y();
456       z = trackReference->Z();
457       pX = trackReference->Px();
458       pY = trackReference->Py();
459       pZ = trackReference->Pz();
460       
461       // check chamberId of current trackReference
462       Int_t detElemId = trackReference->UserId();
463       Int_t chamberId = detElemId / 100 - 1;
464       if (chamberId < 0 || chamberId >= AliMUONConstants::NTrackingCh()) continue;
465       
466       // set hit parameters
467       hit->SetUniqueID(AliMUONVCluster::BuildUniqueID(chamberId, detElemId, iHit));
468       hit->SetXYZ(x,y,z);
469       hit->SetErrXY(0.,0.);
470       
471       // compute track parameters at hit
472       bendingSlope = 0;
473       nonBendingSlope = 0;
474       inverseBendingMomentum = 0;
475       if (TMath::Abs(pZ) > 0) {
476         bendingSlope = pY/pZ;
477         nonBendingSlope = pX/pZ;
478       }
479       Double_t pYZ = TMath::Sqrt(pY*pY+pZ*pZ);
480       if (pYZ >0) inverseBendingMomentum = 1/pYZ; 
481       inverseBendingMomentum *= charge;
482       
483       // set track parameters at hit
484       AliMUONTrackParam trackParam;
485       trackParam.SetNonBendingCoor(x);
486       trackParam.SetBendingCoor(y);
487       trackParam.SetZ(z);
488       trackParam.SetBendingSlope(bendingSlope);
489       trackParam.SetNonBendingSlope(nonBendingSlope);
490       trackParam.SetInverseBendingMomentum(inverseBendingMomentum);
491       
492       // add track parameters at current hit to the track
493       track.AddTrackParamAtCluster(trackParam, *hit, kTRUE);
494     }
495     
496     // if none of the track hits was in MUON, goto the next track
497     if (track.GetNClusters() < 1) continue;
498     
499     // get track parameters at particle's vertex
500     x = particle->Vx();
501     y = particle->Vy();
502     z = particle->Vz();
503     pX = particle->Px();
504     pY = particle->Py();
505     pZ = particle->Pz();
506     
507     // compute rest of track parameters at particle's vertex
508     bendingSlope = 0;
509     nonBendingSlope = 0;
510     inverseBendingMomentum = 0;
511     if (TMath::Abs(pZ) > 0) {
512       bendingSlope = pY/pZ;
513       nonBendingSlope = pX/pZ;
514     }
515     Double_t pYZ = TMath::Sqrt(pY*pY+pZ*pZ);
516     if (pYZ >0) inverseBendingMomentum = 1/pYZ;
517     inverseBendingMomentum *= charge;
518     
519     // set track parameters at particle's vertex
520     AliMUONTrackParam trackParamAtVertex;
521     trackParamAtVertex.SetNonBendingCoor(x);
522     trackParamAtVertex.SetBendingCoor(y);
523     trackParamAtVertex.SetZ(z);
524     trackParamAtVertex.SetBendingSlope(bendingSlope);
525     trackParamAtVertex.SetNonBendingSlope(nonBendingSlope);
526     trackParamAtVertex.SetInverseBendingMomentum(inverseBendingMomentum);
527     
528     // add track parameters at vertex
529     track.SetTrackParamAtVertex(&trackParamAtVertex);
530     
531     // store the track
532     track.SetUniqueID(iTrackRef);
533     tmpTrackRefStore->Add(track);
534   }
535   
536   CleanMuonTrackRef(tmpTrackRefStore);
537   
538   delete hit;
539   delete tmpTrackRefStore;
540 }
541
542 //_____________________________________________________________________________
543 void AliMUONRecoCheck::MakeTriggerableTracks()
544 {
545   /// Make triggerable tracks
546  if (!(fRecoTriggerRefStore = AliMUONESDInterface::NewTriggerTrackStore()))
547    return;
548
549   Double_t x, y, z, slopeX, slopeY, pZ, xLoc, yLoc, zLoc;
550   TParticle* particle;
551   TClonesArray* trackRefs;
552   Int_t nTrackRef = fMCEventHandler->MCEvent()->GetNumberOfTracks();
553         
554   // loop over simulated tracks
555   for (Int_t iTrackRef  = 0; iTrackRef < nTrackRef; ++iTrackRef) {
556     Int_t nHits = fMCEventHandler->GetParticleAndTR(iTrackRef, particle, trackRefs);
557                 
558     // skip empty trackRefs
559     if (nHits < 1) continue;
560                                 
561     AliMUONTriggerTrack track;
562     AliMUONDigitStoreV2S digitStore;
563     Int_t hitsOnTrigger = 0;
564     Int_t currCh = -1;
565                 
566     // loop over simulated track hits
567     for (Int_t iHit = 0; iHit < nHits; ++iHit) {        
568       AliTrackReference* trackReference = static_cast<AliTrackReference*>(trackRefs->UncheckedAt(iHit));
569                         
570       // skip trackRefs not in MUON
571       if (trackReference->DetectorId() != AliTrackReference::kMUON) continue;
572
573       // check chamberId of current trackReference
574       Int_t detElemId = trackReference->UserId();
575       Int_t chamberId = detElemId / 100 - 1;
576       if (chamberId < AliMUONConstants::NTrackingCh() || chamberId >= AliMUONConstants::NCh() ) continue;
577
578       // Get track parameters of current hit
579       x = trackReference->X();
580       y = trackReference->Y();
581       z = trackReference->Z();
582       
583       if ( InitTriggerResponse() ) {
584         fGeometryTransformer->Global2Local(detElemId, x, y, z, xLoc, yLoc, zLoc);
585       
586         Int_t nboard = 0;
587         for ( Int_t cath = AliMp::kCath0; cath <= AliMp::kCath1; ++cath )
588         {
589           const AliMpVSegmentation* seg 
590           = AliMpSegmentation::Instance()
591           ->GetMpSegmentation(detElemId,AliMp::GetCathodType(cath));
592         
593           AliMpPad pad = seg->PadByPosition(xLoc,yLoc,kFALSE);
594           Int_t ix = pad.GetIx();
595           Int_t iy = pad.GetIy();
596         
597           if ( !pad.IsValid() ) continue;
598         
599           if ( cath == AliMp::kCath0 ) nboard = pad.GetLocalBoardId(0);
600         
601           AliMUONVDigit* digit = digitStore.Add(detElemId, nboard, pad.GetLocalBoardChannel(0),
602                                                 cath, AliMUONVDigitStore::kDeny);
603           if (digit) {
604             digit->SetPadXY(ix,iy);
605             digit->SetCharge(1.);
606           }
607         }
608       }
609     
610       if ( hitsOnTrigger == 0 ) {
611         pZ = trackReference->Pz();
612         slopeX = ( pZ == 0. ) ? 99999. : trackReference->Px() / pZ;
613         slopeY = ( pZ == 0. ) ? 99999. : trackReference->Py() / pZ;
614
615         track.SetX11(x);
616         track.SetY11(y);
617         track.SetZ11(z);
618         track.SetSlopeX(slopeX);
619         track.SetSlopeY(slopeY);
620       }
621
622       if ( currCh != chamberId ) {
623         hitsOnTrigger++;
624         currCh = chamberId;
625       }
626
627     } // loop on hits
628     
629     if ( hitsOnTrigger < 3 ) continue;
630
631     // Check if the track passes the trigger algorithm
632     if ( InitTriggerResponse() ) {
633       AliMUONTriggerStoreV1 triggerStore;
634       fTriggerElectronics->Digits2Trigger(digitStore,triggerStore);
635     
636       TIter next(triggerStore.CreateIterator());
637       AliMUONLocalTrigger* locTrg(0x0);
638     
639       Int_t ptCutLevel = 0;
640     
641       while ( ( locTrg = static_cast<AliMUONLocalTrigger*>(next()) ) )
642       {
643         if ( locTrg->IsTrigX() && locTrg->IsTrigY() ) 
644         {
645           ptCutLevel = TMath::Max(ptCutLevel, 1);
646           if ( locTrg->LoHpt() ) ptCutLevel = TMath::Max(ptCutLevel, 3);
647           else if ( locTrg->LoLpt() ) ptCutLevel = TMath::Max(ptCutLevel, 2);
648         } // board is fired 
649       } // end of loop on Local Trigger
650       track.SetPtCutLevel(ptCutLevel);
651     }
652     track.SetUniqueID(iTrackRef);
653     fRecoTriggerRefStore->Add(track);
654   }
655 }
656
657
658 //_____________________________________________________________________________
659 void AliMUONRecoCheck::CleanMuonTrackRef(const AliMUONVTrackStore *tmpTrackRefStore)
660 {
661   /// Re-calculate hits parameters because two AliTrackReferences are recorded for
662   /// each chamber (one when particle is entering + one when particle is leaving 
663   /// the sensitive volume) 
664   if (!(fTrackRefStore = AliMUONESDInterface::NewTrackStore())) return;
665   
666   Double_t maxGasGap = 1.; // cm 
667   Double_t x, y, z, pX, pY, pZ, x1, y1, z1, pX1, pY1, pZ1, z2;
668   Double_t bendingSlope,nonBendingSlope,inverseBendingMomentum;
669   AliMUONVCluster* hit = AliMUONESDInterface::NewCluster();
670   
671   // create iterator
672   TIter next(tmpTrackRefStore->CreateIterator());
673   
674   // loop over tmpTrackRef
675   AliMUONTrack* track;
676   while ( ( track = static_cast<AliMUONTrack*>(next()) ) ) {
677     
678     AliMUONTrack newTrack;
679     
680     // loop over tmpTrackRef's hits
681     Int_t iHit1 = 0;
682     Int_t nTrackHits = track->GetNClusters();
683     while (iHit1 < nTrackHits) {
684       AliMUONTrackParam *trackParam1 = (AliMUONTrackParam*) track->GetTrackParamAtCluster()->UncheckedAt(iHit1);
685       
686       // get track parameters at hit1
687       x1  = trackParam1->GetNonBendingCoor();
688       y1  = trackParam1->GetBendingCoor();
689       z1  = trackParam1->GetZ();
690       pX1 = trackParam1->Px();
691       pY1 = trackParam1->Py();
692       pZ1 = trackParam1->Pz();
693       
694       // prepare new track parameters
695       x  = x1;
696       y  = y1;
697       z  = z1;
698       pX = pX1;
699       pY = pY1;
700       pZ = pZ1;
701       
702       // loop over next tmpTrackRef's hits
703       Int_t nCombinedHits = 1;
704       for (Int_t iHit2 = iHit1+1; iHit2 < nTrackHits; iHit2++) {
705         AliMUONTrackParam *trackParam2 = (AliMUONTrackParam*) track->GetTrackParamAtCluster()->UncheckedAt(iHit2);
706         
707         // get z position of hit2
708         z2 = trackParam2->GetZ();
709         
710         // complete new track parameters if hit2 is on the same detection element
711         if ( TMath::Abs(z2-z1) < maxGasGap ) {
712           x  += trackParam2->GetNonBendingCoor();
713           y  += trackParam2->GetBendingCoor();
714           z  += z2;
715           pX += trackParam2->Px();
716           pY += trackParam2->Py();
717           pZ += trackParam2->Pz();
718           nCombinedHits++;
719           iHit1 = iHit2;
720         }
721         
722       }
723       
724       // finalize new track parameters
725       x  /= (Double_t)nCombinedHits;
726       y  /= (Double_t)nCombinedHits;
727       z  /= (Double_t)nCombinedHits;
728       pX /= (Double_t)nCombinedHits;
729       pY /= (Double_t)nCombinedHits;
730       pZ /= (Double_t)nCombinedHits;
731       bendingSlope = 0;
732       nonBendingSlope = 0;
733       inverseBendingMomentum = 0;
734       if (TMath::Abs(pZ) > 0) {
735         bendingSlope = pY/pZ;
736         nonBendingSlope = pX/pZ;
737       }
738       Double_t pYZ = TMath::Sqrt(pY*pY+pZ*pZ);
739       if (pYZ >0) inverseBendingMomentum = 1/pYZ; 
740       inverseBendingMomentum *= trackParam1->GetCharge();
741       
742       // set hit parameters
743       hit->SetUniqueID(trackParam1->GetClusterPtr()->GetUniqueID());
744       hit->SetXYZ(x,y,z);
745       hit->SetErrXY(0.,0.);
746       
747       // set new track parameters at new hit
748       AliMUONTrackParam trackParam;
749       trackParam.SetNonBendingCoor(x);
750       trackParam.SetBendingCoor(y);
751       trackParam.SetZ(z);
752       trackParam.SetBendingSlope(bendingSlope);
753       trackParam.SetNonBendingSlope(nonBendingSlope);
754       trackParam.SetInverseBendingMomentum(inverseBendingMomentum);
755       
756       // add track parameters at current hit to the track
757       newTrack.AddTrackParamAtCluster(trackParam, *hit, kTRUE);
758       
759       iHit1++;
760     }
761     
762     newTrack.SetUniqueID(track->GetUniqueID());
763     newTrack.SetTrackParamAtVertex(track->GetTrackParamAtVertex());
764     fTrackRefStore->Add(newTrack);
765     
766   }
767   
768   delete hit;
769 }
770
771 //_____________________________________________________________________________
772 void AliMUONRecoCheck::MakeReconstructibleTracks(UInt_t requestedStationMask, Bool_t request2ChInSameSt45,
773                                                  Bool_t hitInFrontOfPad)
774 {
775   /// Isolate the reconstructible tracks
776   if (!(fRecoTrackRefStore = AliMUONESDInterface::NewTrackStore())) return;
777   
778   // need geometry transformer and segmentation to check if trackrefs are located in front of pad(s)
779   if (hitInFrontOfPad) {
780     if (!InitGeometryTransformer()) return;
781     if (!AliMpSegmentation::Instance(kFALSE) && !AliMUONCDB::LoadMapping(kTRUE)) return;
782   }
783   
784   // create iterator on trackRef
785   TIter next(fTrackRefStore->CreateIterator());
786   
787   // loop over trackRef and add reconstructible tracks to fRecoTrackRefStore
788   AliMUONTrack* track;
789   while ( ( track = static_cast<AliMUONTrack*>(next()) ) ) {
790     
791     // check if the track is valid as it is
792     if (track->IsValid(requestedStationMask, request2ChInSameSt45)) {
793       
794       if (hitInFrontOfPad) {
795         
796         AliMUONTrack trackTmp(*track);
797         
798         // loop over clusters and remove those which are not in front of pad(s)
799         AliMUONTrackParam *param = static_cast<AliMUONTrackParam*>(trackTmp.GetTrackParamAtCluster()->First());
800         while (param) {
801           AliMUONTrackParam *nextParam = static_cast<AliMUONTrackParam*>(trackTmp.GetTrackParamAtCluster()->After(param));
802           if (!IsHitInFrontOfPad(param)) trackTmp.RemoveTrackParamAtCluster(param);
803           param = nextParam;
804         }
805         
806         // add the track if it is still valid
807         if (trackTmp.IsValid(requestedStationMask, request2ChInSameSt45)) fRecoTrackRefStore->Add(trackTmp);
808         
809       } else {
810         
811         // add the track
812         fRecoTrackRefStore->Add(*track);
813         
814       }
815       
816     }
817       
818   }
819
820 }
821
822 //_____________________________________________________________________________
823 AliMUONTrack* AliMUONRecoCheck::FindCompatibleTrack(AliMUONTrack &track, AliMUONVTrackStore &trackStore,
824                                                     Int_t &nMatchClusters, Bool_t useLabel, Double_t sigmaCut)
825 {
826   /// Return the track from the store matched with the given track (or 0x0) and the number of matched clusters.
827   /// Matching is done by using the MC label of by comparing cluster/TrackRef positions according to the flag "useLabel".
828   /// WARNING: Who match who matters since the matching algorithm uses the *fraction* of matched clusters of the given track
829   
830   AliMUONTrack *matchedTrack = 0x0;
831   nMatchClusters = 0;
832   
833   if (useLabel) { // by using the MC label
834     
835     // get the corresponding simulated track if any
836     Int_t label = track.GetMCLabel();
837     matchedTrack = (AliMUONTrack*) trackStore.FindObject(label);
838     
839     // get the fraction of matched clusters
840     if (matchedTrack) {
841       Int_t nClusters = track.GetNClusters();
842       for (Int_t iCl = 0; iCl < nClusters; iCl++)
843         if (((AliMUONTrackParam*) track.GetTrackParamAtCluster()->UncheckedAt(iCl))->GetClusterPtr()->GetMCLabel() == label)
844           nMatchClusters++;
845     }
846     
847   } else { // by comparing cluster/TrackRef positions
848     
849     // look for the corresponding simulated track if any
850     TIter next(trackStore.CreateIterator());
851     AliMUONTrack* track2;
852     while ( ( track2 = static_cast<AliMUONTrack*>(next()) ) ) {
853       
854       // check compatibility
855       Int_t n = 0;
856       if (track.Match(*track2, sigmaCut, n)) {
857         matchedTrack = track2;
858         nMatchClusters = n;
859         break;
860       }
861       
862     }
863     
864   }
865   
866   return matchedTrack;
867   
868 }
869
870
871 //_____________________________________________________________________________
872 AliMUONTriggerTrack* AliMUONRecoCheck::FindCompatibleTrack(AliMUONTriggerTrack &track, const AliMUONVTriggerTrackStore &triggerTrackStore,
873                                                            Double_t sigmaCut)
874 {
875   /// Return the trigger track from the store matched with the given track (or 0x0).
876   /// Matching is done by comparing cluster/TrackRef positions.
877   
878   AliMUONTriggerTrack *matchedTrack = 0x0;
879   
880   // look for the corresponding simulated track if any
881   TIter next(triggerTrackStore.CreateIterator());
882   AliMUONTriggerTrack* track2;
883   while ( ( track2 = static_cast<AliMUONTriggerTrack*>(next()) ) ) {
884     // check compatibility
885     if (track.Match(*track2, sigmaCut)) {
886       matchedTrack = track2;
887       break;
888     }
889   }
890   
891   return matchedTrack;
892   
893 }
894
895
896 //____________________________________________________________________________ 
897 Bool_t AliMUONRecoCheck::InitTriggerResponse()
898 {
899   /// Initialize trigger electronics
900   /// for building of triggerable tracks from MC
901   
902   if ( fTriggerElectronics ) return kTRUE;
903   
904   if ( ! InitGeometryTransformer() ) return kFALSE;
905   
906   if ( ! AliMpDDLStore::Instance(false) && ! AliMpCDB::LoadDDLStore()) return kFALSE;
907   
908   if ( ! InitCalibrationData() ) return kFALSE;
909   
910   fTriggerElectronics = new AliMUONTriggerElectronics(fCalibrationData);
911   
912   return kTRUE;
913 }
914
915
916 //____________________________________________________________________________ 
917 Bool_t AliMUONRecoCheck::InitCalibrationData()
918 {
919   /// Initialize calibration data
920   if ( ! fCalibrationData ) {
921     if ( !AliMUONCDB::CheckOCDB() ) return kFALSE;
922     fCalibrationData = new AliMUONCalibrationData(AliCDBManager::Instance()->GetRun());
923   }
924   return kTRUE;
925 }
926
927
928 //____________________________________________________________________________ 
929 Bool_t AliMUONRecoCheck::InitGeometryTransformer()
930 {
931   /// Return geometry transformer
932   /// (create it if necessary)
933   if ( ! fGeometryTransformer ) {
934     
935     if ( !AliMUONCDB::CheckOCDB() ) return kFALSE;
936     
937     if (!AliGeomManager::GetGeometry()) {
938       AliGeomManager::LoadGeometry();
939       if (!AliGeomManager::GetGeometry()) return kFALSE;  
940       if (!AliGeomManager::ApplyAlignObjsFromCDB("MUON")) return kFALSE;
941     }
942     
943     fGeometryTransformer = new AliMUONGeometryTransformer();
944     fGeometryTransformer->LoadGeometryData();
945   }
946   
947   return kTRUE;
948 }
949
950 //________________________________________________________________________
951 Bool_t AliMUONRecoCheck::IsHitInFrontOfPad(AliMUONTrackParam *param) const
952 {
953   /// Return kTRUE if the hit is located in front of at least 1 pad
954   
955   Int_t deId = param->GetClusterPtr()->GetDetElemId();
956   
957   Double_t xL, yL, zL;
958   fGeometryTransformer->Global2Local(deId, param->GetNonBendingCoor(), param->GetBendingCoor(), param->GetZ(), xL, yL, zL);
959   
960   const AliMpVSegmentation* seg1 = AliMpSegmentation::Instance()->GetMpSegmentation(deId, AliMp::kCath0);
961   const AliMpVSegmentation* seg2 = AliMpSegmentation::Instance()->GetMpSegmentation(deId, AliMp::kCath1);
962   if (!seg1 || !seg2) return kFALSE;
963   
964   AliMpPad pad1 = seg1->PadByPosition(xL, yL, kFALSE);
965   AliMpPad pad2 = seg2->PadByPosition(xL, yL, kFALSE);
966   
967   return (pad1.IsValid() || pad2.IsValid());
968 }
969