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->BeginEvent(event)) {
 304       AliError(Form("fails to read MC objects for event %d", event));
 305       return 0x0;
 306     }
 307     MakeTrackRefs();
 308     fMCEventHandler->FinishEvent();
 309     return fTrackRefStore;
 310   }
 311 }
 312
 313 //_____________________________________________________________________________
 314 AliMUONVTriggerTrackStore* AliMUONRecoCheck::TriggerableTracks(Int_t event)
 315 {
 316   /// Return a trigger track store containing the triggerable track references (converted into
 317   /// AliMUONTriggerTrack objects) for a given event
 318
 319   if (!fESDEventOwner) {
 320     if (fRecoTriggerRefStore == 0x0) MakeTriggerableTracks();
 321     return fRecoTriggerRefStore;
 322   }
 323
 324   if (event != fCurrentEvent) {
 325     ResetStores();
 326     fCurrentEvent = event;
 327   }
 328
 329   if (fRecoTriggerRefStore != 0x0) return fRecoTriggerRefStore;
 330   else {
 331     if (!fMCEventHandler->BeginEvent(event)) {
 332       AliError(Form("fails to read MC objects for event %d", event));
 333       return 0x0;
 334     }
 335     MakeTriggerableTracks();
 336     fMCEventHandler->FinishEvent();
 337     return fRecoTriggerRefStore;
 338   }
 339 }
 340
 341
 342 //_____________________________________________________________________________
 343 AliMUONVTrackStore* AliMUONRecoCheck::ReconstructibleTracks(Int_t event, UInt_t requestedStationMask,
 344                                                             Bool_t request2ChInSameSt45, Bool_t hitInFrontOfPad)
 345 {
 346   /// Return a track store containing the reconstructible tracks for a given event,
 347   /// according to the mask of requested stations and the minimum number of chambers hit in stations 4 & 5.
 348
 349   if (!fESDEventOwner) {
 350     if (fRecoTrackRefStore == 0x0) {
 351       if (TrackRefs(event) == 0x0) return 0x0;
 352       MakeReconstructibleTracks(requestedStationMask, request2ChInSameSt45, hitInFrontOfPad);
 353     }
 354     return fRecoTrackRefStore;
 355   }
 356
 357   if (event != fCurrentEvent) {
 358     ResetStores();
 359     fCurrentEvent = event;
 360   }
 361
 362   if (fRecoTrackRefStore != 0x0) return fRecoTrackRefStore;
 363   else {
 364     if (TrackRefs(event) == 0x0) return 0x0;
 365     MakeReconstructibleTracks(requestedStationMask, request2ChInSameSt45, hitInFrontOfPad);
 366     return fRecoTrackRefStore;
 367   }
 368 }
 369
 370
 371 //_____________________________________________________________________________
 372 void AliMUONRecoCheck::MakeReconstructedTracks(Bool_t refit)
 373 {
 374   /// Make reconstructed tracks
 375   if (!(fRecoTrackStore = AliMUONESDInterface::NewTrackStore())) return;
 376
 377   // loop over all reconstructed tracks and add them to the store (skip ghosts)
 378   Int_t nTracks = (Int_t) fESDEvent->GetNumberOfMuonTracks();
 379   for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
 380     AliESDMuonTrack* esdTrack = fESDEvent->GetMuonTrack(iTrack);
 381     if (esdTrack->ContainTrackerData()) AliMUONESDInterface::Add(*esdTrack, *fRecoTrackStore, refit);
 382   }
 383
 384 }
 385
 386
 387 //_____________________________________________________________________________
 388 void AliMUONRecoCheck::MakeTriggeredTracks()
 389 {
 390   /// Make reconstructed trigger tracks
 391   if (!(fRecoTriggerTrackStore = AliMUONESDInterface::NewTriggerTrackStore())) return;
 392
 393   AliMUONVTriggerStore* tmpTriggerStore = AliMUONESDInterface::NewTriggerStore();
 394   if ( ! tmpTriggerStore ) return;
 395
 396   // loop over all reconstructed tracks and add them to the store (include ghosts)
 397   Int_t nTracks = (Int_t) fESDEvent->GetNumberOfMuonTracks();
 398   for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
 399     AliESDMuonTrack* esdTrack = fESDEvent->GetMuonTrack(iTrack);
 400     if (esdTrack->ContainTriggerData()) AliMUONESDInterface::Add(*esdTrack, *tmpTriggerStore);
 401   }
 402
 403   if ( InitCircuit() ) {
 404
 405     AliMUONVTrackReconstructor* tracker = AliMUONESDInterface::GetTracker();
 406     tracker->EventReconstructTrigger(*fTriggerCircuit, *tmpTriggerStore, *fRecoTriggerTrackStore);
 407   }
 408
 409   delete tmpTriggerStore;
 410 }
 411
 412 //_____________________________________________________________________________
 413 Bool_t AliMUONRecoCheck::TriggerToTrack(const AliMUONLocalTrigger& locTrg, AliMUONTriggerTrack& triggerTrack)
 414 {
 415   /// Make trigger track from local trigger info
 416   if ( ! InitCircuit() ) return kFALSE;
 417   AliMUONVTrackReconstructor* tracker = AliMUONESDInterface::GetTracker();
 418   return tracker->TriggerToTrack(*fTriggerCircuit, locTrg, triggerTrack);
 419 }
 420
 421
 422 //_____________________________________________________________________________
 423 void AliMUONRecoCheck::MakeTrackRefs()
 424 {
 425   /// Make reconstructible tracks
 426   AliMUONVTrackStore *tmpTrackRefStore = AliMUONESDInterface::NewTrackStore();
 427   if (!tmpTrackRefStore) return;
 428
 429   Double_t x, y, z, pX, pY, pZ, bendingSlope, nonBendingSlope, inverseBendingMomentum;
 430   TParticle* particle;
 431   TClonesArray* trackRefs;
 432   Int_t nTrackRef = fMCEventHandler->MCEvent()->GetNumberOfTracks();
 433   AliMUONVCluster* hit = AliMUONESDInterface::NewCluster();
 434
 435   // loop over simulated tracks
 436   for (Int_t iTrackRef  = 0; iTrackRef < nTrackRef; ++iTrackRef) {
 437     Int_t nHits = fMCEventHandler->GetParticleAndTR(iTrackRef, particle, trackRefs);
 438
 439     // skip empty trackRefs
 440     if (nHits < 1) continue;
 441
 442     // get the particle charge for further calculation
 443     TParticlePDG* ppdg = particle->GetPDG();
 444     Int_t charge = ppdg != NULL ? (Int_t)(ppdg->Charge()/3.0) : 0;
 445
 446     AliMUONTrack track;
 447
 448     // loop over simulated track hits
 449     for (Int_t iHit = 0; iHit < nHits; ++iHit) {
 450       AliTrackReference* trackReference = static_cast<AliTrackReference*>(trackRefs->UncheckedAt(iHit));
 451
 452       // skip trackRefs not in MUON
 453       if (trackReference->DetectorId() != AliTrackReference::kMUON) continue;
 454
 455       // Get track parameters of current hit
 456       x = trackReference->X();
 457       y = trackReference->Y();
 458       z = trackReference->Z();
 459       pX = trackReference->Px();
 460       pY = trackReference->Py();
 461       pZ = trackReference->Pz();
 462
 463       // check chamberId of current trackReference
 464       Int_t detElemId = trackReference->UserId();
 465       Int_t chamberId = detElemId / 100 - 1;
 466       if (chamberId < 0 || chamberId >= AliMUONConstants::NTrackingCh()) continue;
 467
 468       // set hit parameters
 469       hit->SetUniqueID(AliMUONVCluster::BuildUniqueID(chamberId, detElemId, iHit));
 470       hit->SetXYZ(x,y,z);
 471       hit->SetErrXY(0.,0.);
 472
 473       // compute track parameters at hit
 474       bendingSlope = 0;
 475       nonBendingSlope = 0;
 476       inverseBendingMomentum = 0;
 477       if (TMath::Abs(pZ) > 0) {
 478         bendingSlope = pY/pZ;
 479         nonBendingSlope = pX/pZ;
 480       }
 481       Double_t pYZ = TMath::Sqrt(pY*pY+pZ*pZ);
 482       if (pYZ >0) inverseBendingMomentum = 1/pYZ;
 483       inverseBendingMomentum *= charge;
 484
 485       // set track parameters at hit
 486       AliMUONTrackParam trackParam;
 487       trackParam.SetNonBendingCoor(x);
 488       trackParam.SetBendingCoor(y);
 489       trackParam.SetZ(z);
 490       trackParam.SetBendingSlope(bendingSlope);
 491       trackParam.SetNonBendingSlope(nonBendingSlope);
 492       trackParam.SetInverseBendingMomentum(inverseBendingMomentum);
 493
 494       // add track parameters at current hit to the track
 495       track.AddTrackParamAtCluster(trackParam, *hit, kTRUE);
 496     }
 497
 498     // if none of the track hits was in MUON, goto the next track
 499     if (track.GetNClusters() < 1) continue;
 500
 501     // get track parameters at particle's vertex
 502     x = particle->Vx();
 503     y = particle->Vy();
 504     z = particle->Vz();
 505     pX = particle->Px();
 506     pY = particle->Py();
 507     pZ = particle->Pz();
 508
 509     // compute rest of track parameters at particle's vertex
 510     bendingSlope = 0;
 511     nonBendingSlope = 0;
 512     inverseBendingMomentum = 0;
 513     if (TMath::Abs(pZ) > 0) {
 514       bendingSlope = pY/pZ;
 515       nonBendingSlope = pX/pZ;
 516     }
 517     Double_t pYZ = TMath::Sqrt(pY*pY+pZ*pZ);
 518     if (pYZ >0) inverseBendingMomentum = 1/pYZ;
 519     inverseBendingMomentum *= charge;
 520
 521     // set track parameters at particle's vertex
 522     AliMUONTrackParam trackParamAtVertex;
 523     trackParamAtVertex.SetNonBendingCoor(x);
 524     trackParamAtVertex.SetBendingCoor(y);
 525     trackParamAtVertex.SetZ(z);
 526     trackParamAtVertex.SetBendingSlope(bendingSlope);
 527     trackParamAtVertex.SetNonBendingSlope(nonBendingSlope);
 528     trackParamAtVertex.SetInverseBendingMomentum(inverseBendingMomentum);
 529
 530     // add track parameters at vertex
 531     track.SetTrackParamAtVertex(&trackParamAtVertex);
 532
 533     // store the track
 534     track.SetUniqueID(iTrackRef);
 535     tmpTrackRefStore->Add(track);
 536   }
 537
 538   CleanMuonTrackRef(tmpTrackRefStore);
 539
 540   delete hit;
 541   delete tmpTrackRefStore;
 542 }
 543
 544 //_____________________________________________________________________________
 545 void AliMUONRecoCheck::MakeTriggerableTracks()
 546 {
 547   /// Make triggerable tracks
 548  if (!(fRecoTriggerRefStore = AliMUONESDInterface::NewTriggerTrackStore()))
 549    return;
 550
 551   Double_t x, y, z, slopeX, slopeY, pZ, xLoc, yLoc, zLoc;
 552   TParticle* particle;
 553   TClonesArray* trackRefs;
 554   Int_t nTrackRef = fMCEventHandler->MCEvent()->GetNumberOfTracks();
 555
 556   // loop over simulated tracks
 557   for (Int_t iTrackRef  = 0; iTrackRef < nTrackRef; ++iTrackRef) {
 558     Int_t nHits = fMCEventHandler->GetParticleAndTR(iTrackRef, particle, trackRefs);
 559
 560     // skip empty trackRefs
 561     if (nHits < 1) continue;
 562
 563     AliMUONTriggerTrack track;
 564     AliMUONDigitStoreV2S digitStore;
 565     Int_t hitsOnTrigger = 0;
 566     Int_t currCh = -1;
 567
 568     // loop over simulated track hits
 569     for (Int_t iHit = 0; iHit < nHits; ++iHit) {
 570       AliTrackReference* trackReference = static_cast<AliTrackReference*>(trackRefs->UncheckedAt(iHit));
 571
 572       // skip trackRefs not in MUON
 573       if (trackReference->DetectorId() != AliTrackReference::kMUON) continue;
 574
 575       // check chamberId of current trackReference
 576       Int_t detElemId = trackReference->UserId();
 577       Int_t chamberId = detElemId / 100 - 1;
 578       if (chamberId < AliMUONConstants::NTrackingCh() || chamberId >= AliMUONConstants::NCh() ) continue;
 579
 580       // Get track parameters of current hit
 581       x = trackReference->X();
 582       y = trackReference->Y();
 583       z = trackReference->Z();
 584
 585       if ( InitTriggerResponse() ) {
 586         fGeometryTransformer->Global2Local(detElemId, x, y, z, xLoc, yLoc, zLoc);
 587
 588         Int_t nboard = 0;
 589         for ( Int_t cath = AliMp::kCath0; cath <= AliMp::kCath1; ++cath )
 590         {
 591           const AliMpVSegmentation* seg
 592           = AliMpSegmentation::Instance()
 593           ->GetMpSegmentation(detElemId,AliMp::GetCathodType(cath));
 594
 595           AliMpPad pad = seg->PadByPosition(xLoc,yLoc,kFALSE);
 596           Int_t ix = pad.GetIx();
 597           Int_t iy = pad.GetIy();
 598
 599           if ( !pad.IsValid() ) continue;
 600
 601           if ( cath == AliMp::kCath0 ) nboard = pad.GetLocalBoardId(0);
 602
 603           AliMUONVDigit* digit = digitStore.Add(detElemId, nboard, pad.GetLocalBoardChannel(0),
 604                                                 cath, AliMUONVDigitStore::kDeny);
 605           if (digit) {
 606             digit->SetPadXY(ix,iy);
 607             digit->SetCharge(1.);
 608           }
 609         }
 610       }
 611
 612       if ( hitsOnTrigger == 0 ) {
 613         pZ = trackReference->Pz();
 614         slopeX = ( pZ == 0. ) ? 99999. : trackReference->Px() / pZ;
 615         slopeY = ( pZ == 0. ) ? 99999. : trackReference->Py() / pZ;
 616
 617         track.SetX11(x);
 618         track.SetY11(y);
 619         track.SetZ11(z);
 620         track.SetSlopeX(slopeX);
 621         track.SetSlopeY(slopeY);
 622       }
 623
 624       if ( currCh != chamberId ) {
 625         hitsOnTrigger++;
 626         currCh = chamberId;
 627       }
 628
 629     } // loop on hits
 630
 631     if ( hitsOnTrigger < 3 ) continue;
 632
 633     // Check if the track passes the trigger algorithm
 634     if ( InitTriggerResponse() ) {
 635       AliMUONTriggerStoreV1 triggerStore;
 636       fTriggerElectronics->Digits2Trigger(digitStore,triggerStore);
 637
 638       TIter next(triggerStore.CreateIterator());
 639       AliMUONLocalTrigger* locTrg(0x0);
 640
 641       Int_t ptCutLevel = 0;
 642
 643       while ( ( locTrg = static_cast<AliMUONLocalTrigger*>(next()) ) )
 644       {
 645         if ( locTrg->IsTrigX() && locTrg->IsTrigY() )
 646         {
 647           ptCutLevel = TMath::Max(ptCutLevel, 1);
 648           if ( locTrg->LoHpt() ) ptCutLevel = TMath::Max(ptCutLevel, 3);
 649           else if ( locTrg->LoLpt() ) ptCutLevel = TMath::Max(ptCutLevel, 2);
 650         } // board is fired
 651       } // end of loop on Local Trigger
 652       track.SetPtCutLevel(ptCutLevel);
 653     }
 654     track.SetUniqueID(iTrackRef);
 655     fRecoTriggerRefStore->Add(track);
 656   }
 657 }
 658
 659
 660 //_____________________________________________________________________________
 661 void AliMUONRecoCheck::CleanMuonTrackRef(const AliMUONVTrackStore *tmpTrackRefStore)
 662 {
 663   /// Re-calculate hits parameters because two AliTrackReferences are recorded for
 664   /// each chamber (one when particle is entering + one when particle is leaving
 665   /// the sensitive volume)
 666   if (!(fTrackRefStore = AliMUONESDInterface::NewTrackStore())) return;
 667
 668   Double_t maxGasGap = 1.; // cm
 669   Double_t x, y, z, pX, pY, pZ, x1, y1, z1, pX1, pY1, pZ1, z2;
 670   Double_t bendingSlope,nonBendingSlope,inverseBendingMomentum;
 671   AliMUONVCluster* hit = AliMUONESDInterface::NewCluster();
 672
 673   // create iterator
 674   TIter next(tmpTrackRefStore->CreateIterator());
 675
 676   // loop over tmpTrackRef
 677   AliMUONTrack* track;
 678   while ( ( track = static_cast<AliMUONTrack*>(next()) ) ) {
 679
 680     AliMUONTrack newTrack;
 681
 682     // loop over tmpTrackRef's hits
 683     Int_t iHit1 = 0;
 684     Int_t nTrackHits = track->GetNClusters();
 685     while (iHit1 < nTrackHits) {
 686       AliMUONTrackParam *trackParam1 = (AliMUONTrackParam*) track->GetTrackParamAtCluster()->UncheckedAt(iHit1);
 687
 688       // get track parameters at hit1
 689       x1  = trackParam1->GetNonBendingCoor();
 690       y1  = trackParam1->GetBendingCoor();
 691       z1  = trackParam1->GetZ();
 692       pX1 = trackParam1->Px();
 693       pY1 = trackParam1->Py();
 694       pZ1 = trackParam1->Pz();
 695
 696       // prepare new track parameters
 697       x  = x1;
 698       y  = y1;
 699       z  = z1;
 700       pX = pX1;
 701       pY = pY1;
 702       pZ = pZ1;
 703
 704       // loop over next tmpTrackRef's hits
 705       Int_t nCombinedHits = 1;
 706       for (Int_t iHit2 = iHit1+1; iHit2 < nTrackHits; iHit2++) {
 707         AliMUONTrackParam *trackParam2 = (AliMUONTrackParam*) track->GetTrackParamAtCluster()->UncheckedAt(iHit2);
 708
 709         // get z position of hit2
 710         z2 = trackParam2->GetZ();
 711
 712         // complete new track parameters if hit2 is on the same detection element
 713         if ( TMath::Abs(z2-z1) < maxGasGap ) {
 714           x  += trackParam2->GetNonBendingCoor();
 715           y  += trackParam2->GetBendingCoor();
 716           z  += z2;
 717           pX += trackParam2->Px();
 718           pY += trackParam2->Py();
 719           pZ += trackParam2->Pz();
 720           nCombinedHits++;
 721           iHit1 = iHit2;
 722         }
 723
 724       }
 725
 726       // finalize new track parameters
 727       x  /= (Double_t)nCombinedHits;
 728       y  /= (Double_t)nCombinedHits;
 729       z  /= (Double_t)nCombinedHits;
 730       pX /= (Double_t)nCombinedHits;
 731       pY /= (Double_t)nCombinedHits;
 732       pZ /= (Double_t)nCombinedHits;
 733       bendingSlope = 0;
 734       nonBendingSlope = 0;
 735       inverseBendingMomentum = 0;
 736       if (TMath::Abs(pZ) > 0) {
 737         bendingSlope = pY/pZ;
 738         nonBendingSlope = pX/pZ;
 739       }
 740       Double_t pYZ = TMath::Sqrt(pY*pY+pZ*pZ);
 741       if (pYZ >0) inverseBendingMomentum = 1/pYZ;
 742       inverseBendingMomentum *= trackParam1->GetCharge();
 743
 744       // set hit parameters
 745       hit->SetUniqueID(trackParam1->GetClusterPtr()->GetUniqueID());
 746       hit->SetXYZ(x,y,z);
 747       hit->SetErrXY(0.,0.);
 748
 749       // set new track parameters at new hit
 750       AliMUONTrackParam trackParam;
 751       trackParam.SetNonBendingCoor(x);
 752       trackParam.SetBendingCoor(y);
 753       trackParam.SetZ(z);
 754       trackParam.SetBendingSlope(bendingSlope);
 755       trackParam.SetNonBendingSlope(nonBendingSlope);
 756       trackParam.SetInverseBendingMomentum(inverseBendingMomentum);
 757
 758       // add track parameters at current hit to the track
 759       newTrack.AddTrackParamAtCluster(trackParam, *hit, kTRUE);
 760
 761       iHit1++;
 762     }
 763
 764     newTrack.SetUniqueID(track->GetUniqueID());
 765     newTrack.SetTrackParamAtVertex(track->GetTrackParamAtVertex());
 766     fTrackRefStore->Add(newTrack);
 767
 768   }
 769
 770   delete hit;
 771 }
 772
 773 //_____________________________________________________________________________
 774 void AliMUONRecoCheck::MakeReconstructibleTracks(UInt_t requestedStationMask, Bool_t request2ChInSameSt45,
 775                                                  Bool_t hitInFrontOfPad)
 776 {
 777   /// Isolate the reconstructible tracks
 778   if (!(fRecoTrackRefStore = AliMUONESDInterface::NewTrackStore())) return;
 779
 780   // need geometry transformer and segmentation to check if trackrefs are located in front of pad(s)
 781   if (hitInFrontOfPad) {
 782     if (!InitGeometryTransformer()) return;
 783     if (!AliMpSegmentation::Instance(kFALSE) && !AliMUONCDB::LoadMapping(kTRUE)) return;
 784   }
 785
 786   // create iterator on trackRef
 787   TIter next(fTrackRefStore->CreateIterator());
 788
 789   // loop over trackRef and add reconstructible tracks to fRecoTrackRefStore
 790   AliMUONTrack* track;
 791   while ( ( track = static_cast<AliMUONTrack*>(next()) ) ) {
 792
 793     // check if the track is valid as it is
 794     if (track->IsValid(requestedStationMask, request2ChInSameSt45)) {
 795
 796       if (hitInFrontOfPad) {
 797
 798         AliMUONTrack trackTmp(*track);
 799
 800         // loop over clusters and remove those which are not in front of pad(s)
 801         AliMUONTrackParam *param = static_cast<AliMUONTrackParam*>(trackTmp.GetTrackParamAtCluster()->First());
 802         while (param) {
 803           AliMUONTrackParam *nextParam = static_cast<AliMUONTrackParam*>(trackTmp.GetTrackParamAtCluster()->After(param));
 804           if (!IsHitInFrontOfPad(param)) trackTmp.RemoveTrackParamAtCluster(param);
 805           param = nextParam;
 806         }
 807
 808         // add the track if it is still valid
 809         if (trackTmp.IsValid(requestedStationMask, request2ChInSameSt45)) fRecoTrackRefStore->Add(trackTmp);
 810
 811       } else {
 812
 813         // add the track
 814         fRecoTrackRefStore->Add(*track);
 815
 816       }
 817
 818     }
 819
 820   }
 821
 822 }
 823
 824 //_____________________________________________________________________________
 825 AliMUONTrack* AliMUONRecoCheck::FindCompatibleTrack(AliMUONTrack &track, AliMUONVTrackStore &trackStore,
 826                                                     Int_t &nMatchClusters, Bool_t useLabel, Double_t sigmaCut)
 827 {
 828   /// Return the track from the store matched with the given track (or 0x0) and the number of matched clusters.
 829   /// Matching is done by using the MC label of by comparing cluster/TrackRef positions according to the flag "useLabel".
 830   /// WARNING: Who match who matters since the matching algorithm uses the *fraction* of matched clusters of the given track
 831
 832   AliMUONTrack *matchedTrack = 0x0;
 833   nMatchClusters = 0;
 834
 835   if (useLabel) { // by using the MC label
 836
 837     // get the corresponding simulated track if any
 838     Int_t label = track.GetMCLabel();
 839     matchedTrack = (AliMUONTrack*) trackStore.FindObject(label);
 840
 841     // get the fraction of matched clusters
 842     if (matchedTrack) {
 843       Int_t nClusters = track.GetNClusters();
 844       for (Int_t iCl = 0; iCl < nClusters; iCl++)
 845         if (((AliMUONTrackParam*) track.GetTrackParamAtCluster()->UncheckedAt(iCl))->GetClusterPtr()->GetMCLabel() == label)
 846           nMatchClusters++;
 847     }
 848
 849   } else { // by comparing cluster/TrackRef positions
 850
 851     // look for the corresponding simulated track if any
 852     TIter next(trackStore.CreateIterator());
 853     AliMUONTrack* track2;
 854     while ( ( track2 = static_cast<AliMUONTrack*>(next()) ) ) {
 855
 856       // check compatibility
 857       Int_t n = 0;
 858       if (track.Match(*track2, sigmaCut, n)) {
 859         matchedTrack = track2;
 860         nMatchClusters = n;
 861         break;
 862       }
 863
 864     }
 865
 866   }
 867
 868   return matchedTrack;
 869
 870 }
 871
 872
 873 //_____________________________________________________________________________
 874 AliMUONTriggerTrack* AliMUONRecoCheck::FindCompatibleTrack(AliMUONTriggerTrack &track, const AliMUONVTriggerTrackStore &triggerTrackStore,
 875                                                            Double_t sigmaCut)
 876 {
 877   /// Return the trigger track from the store matched with the given track (or 0x0).
 878   /// Matching is done by comparing cluster/TrackRef positions.
 879
 880   AliMUONTriggerTrack *matchedTrack = 0x0;
 881
 882   // look for the corresponding simulated track if any
 883   TIter next(triggerTrackStore.CreateIterator());
 884   AliMUONTriggerTrack* track2;
 885   while ( ( track2 = static_cast<AliMUONTriggerTrack*>(next()) ) ) {
 886     // check compatibility
 887     if (track.Match(*track2, sigmaCut)) {
 888       matchedTrack = track2;
 889       break;
 890     }
 891   }
 892
 893   return matchedTrack;
 894
 895 }
 896
 897
 898 //____________________________________________________________________________
 899 Bool_t AliMUONRecoCheck::InitTriggerResponse()
 900 {
 901   /// Initialize trigger electronics
 902   /// for building of triggerable tracks from MC
 903
 904   if ( fTriggerElectronics ) return kTRUE;
 905
 906   if ( ! InitGeometryTransformer() ) return kFALSE;
 907
 908   if ( ! AliMpDDLStore::Instance(false) && ! AliMpCDB::LoadDDLStore()) return kFALSE;
 909
 910   if ( ! InitCalibrationData() ) return kFALSE;
 911
 912   fTriggerElectronics = new AliMUONTriggerElectronics(fCalibrationData);
 913
 914   return kTRUE;
 915 }
 916
 917
 918 //____________________________________________________________________________
 919 Bool_t AliMUONRecoCheck::InitCalibrationData()
 920 {
 921   /// Initialize calibration data
 922   if ( ! fCalibrationData ) {
 923     if ( !AliMUONCDB::CheckOCDB() ) return kFALSE;
 924     fCalibrationData = new AliMUONCalibrationData(AliCDBManager::Instance()->GetRun());
 925   }
 926   return kTRUE;
 927 }
 928
 929
 930 //____________________________________________________________________________
 931 Bool_t AliMUONRecoCheck::InitGeometryTransformer()
 932 {
 933   /// Return geometry transformer
 934   /// (create it if necessary)
 935   if ( ! fGeometryTransformer ) {
 936
 937     if ( !AliMUONCDB::CheckOCDB() ) return kFALSE;
 938
 939     if (!AliGeomManager::GetGeometry()) {
 940       AliGeomManager::LoadGeometry();
 941       if (!AliGeomManager::GetGeometry()) return kFALSE;
 942       if (!AliGeomManager::ApplyAlignObjsFromCDB("MUON")) return kFALSE;
 943     }
 944
 945     fGeometryTransformer = new AliMUONGeometryTransformer();
 946     fGeometryTransformer->LoadGeometryData();
 947   }
 948
 949   return kTRUE;
 950 }
 951
 952 //________________________________________________________________________
 953 Bool_t AliMUONRecoCheck::IsHitInFrontOfPad(AliMUONTrackParam *param) const
 954 {
 955   /// Return kTRUE if the hit is located in front of at least 1 pad
 956
 957   Int_t deId = param->GetClusterPtr()->GetDetElemId();
 958
 959   Double_t xL, yL, zL;
 960   fGeometryTransformer->Global2Local(deId, param->GetNonBendingCoor(), param->GetBendingCoor(), param->GetZ(), xL, yL, zL);
 961
 962   const AliMpVSegmentation* seg1 = AliMpSegmentation::Instance()->GetMpSegmentation(deId, AliMp::kCath0);
 963   const AliMpVSegmentation* seg2 = AliMpSegmentation::Instance()->GetMpSegmentation(deId, AliMp::kCath1);
 964   if (!seg1 || !seg2) return kFALSE;
 965
 966   AliMpPad pad1 = seg1->PadByPosition(xL, yL, kFALSE);
 967   AliMpPad pad2 = seg2->PadByPosition(xL, yL, kFALSE);
 968
 969   return (pad1.IsValid() || pad2.IsValid());
 970 }
 971