MUON/AliMUONReconstructor.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 /* $Id$ */
  16
  17 //-----------------------------
  18 // Class AliMUONReconstructor
  19 //-----------------------------
  20 // Class for the
  21 // MUON track reconstruction
  22
  23 #include "AliMUONReconstructor.h"
  24
  25 #include "AliMUONConstants.h"
  26 #include "AliMUONCalibrationData.h"
  27 #include "AliMUONClusterFinderAZ.h"
  28 #include "AliMUONClusterReconstructor.h"
  29 #include "AliMUONData.h"
  30 #include "AliMUONDigitCalibrator.h"
  31 #include "AliMUONEventRecoCombi.h"
  32 #include "AliMUONDigitMaker.h"
  33 #include "AliMUONTrack.h"
  34 #include "AliMUONTrackParam.h"
  35 #include "AliMUONTrackExtrap.h"
  36 #include "AliMUONTracker.h"
  37 #include "AliMUONVTrackReconstructor.h"
  38 #include "AliMUONTrackReconstructor.h"
  39 #include "AliMUONTrackReconstructorK.h"
  40 #include "AliMUONTriggerTrack.h"
  41 #include "AliMUONTriggerCircuit.h"
  42 #include "AliMUONTriggerCrateStore.h"
  43 #include "AliMUONSegFactory.h"
  44 #include "AliMUONSegmentation.h"
  45 #include "AliMUONPreClusterFinder.h"
  46 #include "AliMUONClusterFinderCOG.h"
  47 #include "AliMUONClusterFinderSimpleFit.h"
  48 #include "AliMUONClusterFinderMLEM.h"
  49
  50 #include "AliESD.h"
  51 #include "AliESDMuonTrack.h"
  52 #include "AliLog.h"
  53 #include "AliRawReader.h"
  54 #include "AliRunLoader.h"
  55 #include "AliCDBManager.h"
  56
  57 #include "TTask.h"
  58 #include "TStopwatch.h"
  59
  60 /// \cond CLASSIMP
  61 ClassImp(AliMUONReconstructor)
  62 /// \endcond
  63
  64 //_____________________________________________________________________________
  65 AliMUONReconstructor::AliMUONReconstructor()
  66   : AliReconstructor(),
  67     fDigitMaker(new AliMUONDigitMaker()),
  68     fCalibrationData(0x0),
  69     fCrateManager(new AliMUONTriggerCrateStore()),
  70     fTriggerCircuit(new TClonesArray("AliMUONTriggerCircuit", 234)),
  71     fTransformer(new AliMUONGeometryTransformer(kTRUE)),
  72     fSegmentation(0x0),
  73     fMUONData(new AliMUONData(0x0,"MUON","MUON"))
  74 {
  75 /// Default constructor
  76
  77     AliDebug(1,"");
  78     // Crate manager
  79     fCrateManager->ReadFromFile();
  80
  81     // set to digit maker
  82     fDigitMaker->SetCrateManager(fCrateManager);
  83
  84     // transformater
  85     fTransformer->ReadGeometryData("volpath.dat", "geometry.root");
  86
  87     // create segmentation and pass it to EventRecoCombi
  88     AliMUONSegFactory factory(fTransformer);
  89     fSegmentation = factory.CreateSegmentation();
  90     AliMUONEventRecoCombi::Instance(fSegmentation);
  91
  92     // trigger circuit
  93     for (Int_t i = 0; i < AliMUONConstants::NTriggerCircuit(); i++)  {
  94       AliMUONTriggerCircuit* c = new AliMUONTriggerCircuit();
  95       c->SetTransformer(fTransformer);
  96       c->Init(i,*fCrateManager);
  97       TClonesArray& circuit = *fTriggerCircuit;
  98       new(circuit[circuit.GetEntriesFast()])AliMUONTriggerCircuit(*c);
  99       delete c;
 100     }
 101
 102
 103 }
 104
 105 //_____________________________________________________________________________
 106 AliMUONReconstructor::~AliMUONReconstructor()
 107 {
 108 /// Destructor
 109
 110   AliDebug(1,"");
 111   delete fCalibrationData;
 112   delete fDigitMaker;
 113   delete fCrateManager;
 114   delete fTriggerCircuit;
 115   delete fTransformer;
 116   delete fSegmentation;
 117   delete fMUONData;
 118 }
 119
 120 //_____________________________________________________________________________
 121 TTask*
 122 AliMUONReconstructor::GetCalibrationTask() const
 123 {
 124 /// Create the calibration task(s).
 125
 126   //const AliRun* run = fRunLoader->GetAliRun();
 127   //Int_t runNumber = run->GetRunNumber();
 128   Int_t runNumber = AliCDBManager::Instance()->GetRun();
 129   AliInfo("Calibration will occur.");
 130
 131   fCalibrationData = new AliMUONCalibrationData(runNumber);
 132   if ( !fCalibrationData->IsValid() )
 133     {
 134       AliError("Could not retrieve calibrations !");
 135       delete fCalibrationData;
 136       fCalibrationData = 0x0;
 137       return 0x0;
 138     }
 139   TTask* calibration = new TTask("MUONCalibrator","MUON Digit calibrator");
 140
 141   TString opt(GetOption());
 142   opt.ToUpper();
 143   Bool_t statusMap(kTRUE);
 144
 145   if ( strstr(opt,"NOSTATUSMAP") )
 146   {
 147     AliWarning("Disconnecting status map : SHOULD BE USED FOR DEBUG ONLY. NOT FOR PRODUCTION !!!");
 148     statusMap = kFALSE;
 149   }
 150   calibration->Add(new AliMUONDigitCalibrator(fMUONData,fCalibrationData,statusMap));
 151   return calibration;
 152 }
 153
 154 //_____________________________________________________________________________
 155 AliMUONClusterReconstructor*
 156 AliMUONReconstructor::CreateClusterReconstructor() const
 157 {
 158 /// Create cluster reconstructor
 159
 160   AliMUONVClusterFinder* clusterFinder(0x0);
 161
 162   TString opt(GetOption());
 163   opt.ToUpper();
 164
 165   if ( strstr(opt,"PRECLUSTER") )
 166   {
 167     clusterFinder = new AliMUONPreClusterFinder;
 168   }
 169   else if ( strstr(opt,"COG") )
 170   {
 171     clusterFinder = new AliMUONClusterFinderCOG;
 172   }
 173   else if ( strstr(opt,"SIMPLEFIT") )
 174   {
 175     clusterFinder = new AliMUONClusterFinderSimpleFit;
 176   }
 177   else if ( strstr(opt,"MLEM:DRAW") )
 178   {
 179     clusterFinder = new AliMUONClusterFinderMLEM(kTRUE);
 180   }
 181   else if ( strstr(opt,"MLEM") )
 182   {
 183     clusterFinder = new AliMUONClusterFinderMLEM(kFALSE);
 184   }
 185
 186   if ( clusterFinder)
 187   {
 188     AliInfo(Form("Will use %s for clusterizing",clusterFinder->ClassName()));
 189   }
 190
 191   AliMUONClusterReconstructor* clusterReco =
 192     new AliMUONClusterReconstructor(fMUONData,clusterFinder,fTransformer);
 193   return clusterReco;
 194 }
 195
 196 //_____________________________________________________________________________
 197 void AliMUONReconstructor::Reconstruct(AliRunLoader* runLoader) const
 198 {
 199 /// Reconstruct
 200 /// \todo add more
 201
 202   AliLoader* loader = runLoader->GetLoader("MUONLoader");
 203   Int_t nEvents     = runLoader->GetNumberOfEvents();
 204   Int_t evtNumber   = runLoader->GetEventNumber();
 205
 206   fMUONData->SetLoader(loader);
 207
 208 // passing loader as argument.
 209   AliMUONVTrackReconstructor* recoEvent;
 210   if (strstr(GetOption(),"Original")) recoEvent = new AliMUONTrackReconstructor(fMUONData);
 211   else if (strstr(GetOption(),"Combi")) recoEvent = new AliMUONTrackReconstructorK(fMUONData,"Combi");
 212   else recoEvent = new AliMUONTrackReconstructorK(fMUONData,"Kalman");
 213
 214   recoEvent->SetTriggerCircuit(fTriggerCircuit);
 215
 216   AliMUONClusterReconstructor* recoCluster = CreateClusterReconstructor();
 217
 218   AliMUONClusterFinderVS *recModel = recoCluster->GetRecoModel();
 219
 220   if (!strstr(GetOption(),"VS")) {
 221     recModel = (AliMUONClusterFinderVS*) new AliMUONClusterFinderAZ();
 222     recoCluster->SetRecoModel(recModel);
 223   }
 224   recModel->SetGhostChi2Cut(10);
 225   recModel->SetEventNumber(evtNumber);
 226
 227   loader->LoadDigits("READ");
 228   loader->LoadRecPoints("RECREATE");
 229   loader->LoadTracks("RECREATE");
 230
 231   TTask* calibration = GetCalibrationTask();
 232
 233   Int_t chBeg = (strstr(GetOption(),"Combi") ? 6 : 0);
 234
 235   //   Loop over events
 236   for(Int_t ievent = 0; ievent < nEvents; ievent++) {
 237
 238     AliDebug(1,Form("Event %d",ievent));
 239
 240     runLoader->GetEvent(ievent);
 241
 242     //----------------------- digit2cluster & Trigger2Trigger -------------------
 243     if (!loader->TreeR()) loader->MakeRecPointsContainer();
 244
 245     // tracking branch
 246     if (!strstr(GetOption(),"Combi")) {
 247       fMUONData->MakeBranch("RC");
 248       fMUONData->SetTreeAddress("D,RC");
 249     } else {
 250       fMUONData->SetTreeAddress("D");
 251       fMUONData->SetTreeAddress("RCC");
 252     }
 253     // Important for avoiding a memory leak when reading digits ( to be investigated more in detail)
 254     // In any case the reading of GLT is needed for the Trigger2Tigger method below
 255     fMUONData->SetTreeAddress("GLT");
 256
 257     fMUONData->GetDigits();
 258
 259     if ( calibration )
 260     {
 261       calibration->ExecuteTask();
 262     }
 263
 264     recoCluster->Digits2Clusters(chBeg);
 265
 266     if (strstr(GetOption(),"Combi")) {
 267       // Combined cluster / track finder
 268       AliMUONEventRecoCombi::Instance()->FillEvent(fMUONData, (AliMUONClusterFinderAZ*)recModel);
 269       ((AliMUONClusterFinderAZ*) recModel)->SetReco(2);
 270     }
 271     else fMUONData->Fill("RC");
 272
 273     // trigger branch
 274     fMUONData->MakeBranch("TC");
 275     fMUONData->SetTreeAddress("TC");
 276     recoCluster->Trigger2Trigger();
 277     fMUONData->Fill("TC");
 278
 279     //AZ loader->WriteRecPoints("OVERWRITE");
 280
 281     //---------------------------- Track & TriggerTrack ---------------------
 282     if (!loader->TreeT()) loader->MakeTracksContainer();
 283
 284     // trigger branch
 285     fMUONData->MakeBranch("RL"); //trigger track
 286     fMUONData->SetTreeAddress("RL");
 287     recoEvent->EventReconstructTrigger();
 288     fMUONData->Fill("RL");
 289
 290     // tracking branch
 291     fMUONData->MakeBranch("RT"); //track
 292     fMUONData->SetTreeAddress("RT");
 293     recoEvent->EventReconstruct();
 294     fMUONData->Fill("RT");
 295
 296     loader->WriteTracks("OVERWRITE");
 297
 298     if (strstr(GetOption(),"Combi")) {
 299       // Combined cluster / track
 300       ((AliMUONClusterFinderAZ*) recModel)->SetReco(1);
 301       fMUONData->MakeBranch("RC");
 302       fMUONData->SetTreeAddress("RC");
 303       AliMUONEventRecoCombi::Instance()->FillRecP(fMUONData, (AliMUONTrackReconstructorK*)recoEvent);
 304       fMUONData->Fill("RC");
 305     }
 306     loader->WriteRecPoints("OVERWRITE");
 307
 308     //--------------------------- Resetting branches -----------------------
 309     fMUONData->ResetDigits();
 310     fMUONData->ResetRawClusters();
 311     fMUONData->ResetTrigger();
 312     fMUONData->ResetRecTracks();
 313     fMUONData->ResetRecTriggerTracks();
 314
 315   }
 316   loader->UnloadDigits();
 317   loader->UnloadRecPoints();
 318   loader->UnloadTracks();
 319
 320   delete recoCluster;
 321   delete recoEvent;
 322   delete calibration;
 323 }
 324
 325 //_____________________________________________________________________________
 326 void AliMUONReconstructor::Reconstruct(AliRunLoader* runLoader,
 327                                        AliRawReader* rawReader) const
 328 {
 329 /// Recontruct
 330 /// \todo add more
 331
 332   //  AliLoader
 333   AliLoader* loader = runLoader->GetLoader("MUONLoader");
 334   Int_t evtNumber   = runLoader->GetEventNumber();
 335
 336   fMUONData->SetLoader(loader);
 337
 338   // passing loader as argument.
 339   fDigitMaker->SetMUONData(fMUONData);
 340
 341   // disable trigger rawdata reading
 342   if (strstr(GetOption(),"TriggerDisable"))
 343       fDigitMaker->DisableTrigger();
 344
 345   AliMUONClusterReconstructor* recoCluster = CreateClusterReconstructor();
 346
 347   AliMUONClusterFinderVS *recModel = recoCluster->GetRecoModel();
 348
 349   if (!strstr(GetOption(),"VS"))
 350   {
 351     recModel = (AliMUONClusterFinderVS*) new AliMUONClusterFinderAZ();
 352     recoCluster->SetRecoModel(recModel);
 353   }
 354   recModel->SetGhostChi2Cut(10);
 355   recModel->SetEventNumber(evtNumber);
 356
 357   TTask* calibration = GetCalibrationTask();
 358
 359   loader->LoadRecPoints("RECREATE");
 360
 361   // Digits are not stored on disk and created on flight from rawdata.
 362   // In order to write digits on disk the following line should be uncommented
 363   // loader->LoadDigits("RECREATE");
 364
 365   //   Loop over events
 366   Int_t iEvent = 0;
 367
 368   TStopwatch totalTimer;
 369   TStopwatch rawTimer;
 370   TStopwatch calibTimer;
 371   TStopwatch clusterTimer;
 372
 373   rawTimer.Start(kTRUE); rawTimer.Stop();
 374   calibTimer.Start(kTRUE); calibTimer.Stop();
 375   clusterTimer.Start(kTRUE); clusterTimer.Stop();
 376
 377   totalTimer.Start(kTRUE);
 378
 379   while (rawReader->NextEvent())
 380   {
 381     AliDebug(1,Form("Event %d",iEvent));
 382
 383     runLoader->GetEvent(iEvent++);
 384
 385     //----------------------- raw2digits & raw2trigger-------------------
 386     //  if (!loader->TreeD())
 387     //  {
 388     //    AliDebug(1,Form("Making Digit Container for event %d",iEvent));
 389     //    loader->MakeDigitsContainer();
 390     //  }
 391     //  Digits are not stored on disk and created on flight from rawdata.
 392     //  In order to write digits on disk the following lines should be uncommented
 393     //  fMUONData->MakeBranch("D,GLT");
 394     //  fMUONData->SetTreeAddress("D,GLT");
 395     fMUONData->SetDataContainer("D, GLT");
 396     rawTimer.Start(kFALSE);
 397     fDigitMaker->Raw2Digits(rawReader);
 398     rawTimer.Stop();
 399
 400     if ( calibration )
 401     {
 402       calibTimer.Start(kFALSE);
 403       calibration->ExecuteTask();
 404       calibTimer.Stop();
 405     }
 406     // Digits are not stored on disk and created on flight from rawdata.
 407     // In order to write digits on disk the following lines should be uncommented
 408     // fMUONData->Fill("D,GLT");
 409     // loader->WriteDigits("OVERWRITE");
 410     //----------------------- digit2cluster & Trigger2Trigger -------------------
 411     clusterTimer.Start(kFALSE);
 412
 413     if (!loader->TreeR()) loader->MakeRecPointsContainer();
 414
 415     // tracking branch
 416     fMUONData->MakeBranch("RC");
 417     fMUONData->SetTreeAddress("RC");
 418     recoCluster->Digits2Clusters();
 419     fMUONData->Fill("RC");
 420
 421     // trigger branch
 422     fMUONData->MakeBranch("TC");
 423     fMUONData->SetTreeAddress("TC");
 424     fMUONData->Fill("TC");
 425
 426     loader->WriteRecPoints("OVERWRITE");
 427
 428     clusterTimer.Stop();
 429
 430
 431     //--------------------------- Resetting branches -----------------------
 432     fMUONData->ResetDigits();
 433     fMUONData->ResetRawClusters();
 434     fMUONData->ResetTrigger();
 435   }
 436
 437   totalTimer.Stop();
 438
 439   loader->UnloadRecPoints();
 440   loader->UnloadDigits();
 441
 442   delete recoCluster;
 443
 444   AliInfo(Form("Execution time for converting RAW data to digits in MUON : R:%.2fs C:%.2fs",
 445                rawTimer.RealTime(),rawTimer.CpuTime()));
 446   AliInfo(Form("Execution time for calibrating MUON : R:%.2fs C:%.2fs",
 447                calibTimer.RealTime(),calibTimer.CpuTime()));
 448   AliInfo(Form("Execution time for clusterizing MUON : R:%.2fs C:%.2fs",
 449                clusterTimer.RealTime(),clusterTimer.CpuTime()));
 450   AliInfo(Form("Total Execution time for Reconstruct(from raw) MUON : R:%.2fs C:%.2fs",
 451                totalTimer.RealTime(),totalTimer.CpuTime()));
 452 }
 453
 454 //_____________________________________________________________________________
 455 void AliMUONReconstructor::FillESD(AliRunLoader* runLoader, AliESD* esd) const
 456 {
 457 /// Fill ESD
 458 /// \todo add more
 459
 460   TClonesArray* recTracksArray = 0;
 461   TClonesArray* recTrigTracksArray = 0;
 462
 463   AliLoader* loader = runLoader->GetLoader("MUONLoader");
 464   loader->LoadTracks("READ");
 465   fMUONData->SetLoader(loader);
 466
 467    // declaration
 468   Int_t iEvent;// nPart;
 469   Int_t nTrackHits;// nPrimary;
 470   Double_t fitFmin;
 471
 472   Double_t bendingSlope, nonBendingSlope, inverseBendingMomentum;
 473   Double_t xRec, yRec, zRec, chi2MatchTrigger;
 474   Bool_t matchTrigger;
 475
 476   // setting pointer for tracks, triggertracks & trackparam at vertex
 477   AliMUONTrack* recTrack = 0;
 478   AliMUONTrackParam* trackParam = 0;
 479   AliMUONTriggerTrack* recTriggerTrack = 0;
 480
 481   iEvent = runLoader->GetEventNumber();
 482   runLoader->GetEvent(iEvent);
 483
 484   // Get vertex
 485   Double_t vertex[3] = {0};
 486   const AliESDVertex *esdVert = esd->GetVertex();
 487   if (esdVert->GetNContributors()) {
 488     esdVert->GetXYZ(vertex);
 489     AliDebug(1, "find vertex\n");
 490   }
 491   // setting ESD MUON class
 492   AliESDMuonTrack* theESDTrack = new  AliESDMuonTrack() ;
 493
 494   //-------------------- trigger tracks-------------
 495   Long_t trigPat = 0;
 496   fMUONData->SetTreeAddress("RL");
 497   fMUONData->GetRecTriggerTracks();
 498   recTrigTracksArray = fMUONData->RecTriggerTracks();
 499
 500   // ready global trigger pattern from first track
 501   if (recTrigTracksArray)
 502     recTriggerTrack = (AliMUONTriggerTrack*) recTrigTracksArray->First();
 503   if (recTriggerTrack) trigPat = recTriggerTrack->GetGTPattern();
 504
 505   //printf(">>> Event %d Number of Recconstructed tracks %d \n",iEvent, nrectracks);
 506
 507   // -------------------- tracks-------------
 508   fMUONData->SetTreeAddress("RT");
 509   fMUONData->GetRecTracks();
 510   recTracksArray = fMUONData->RecTracks();
 511
 512   Int_t nRecTracks = 0;
 513   if (recTracksArray)
 514     nRecTracks = (Int_t) recTracksArray->GetEntriesFast(); //
 515
 516   // loop over tracks
 517   for (Int_t iRecTracks = 0; iRecTracks <  nRecTracks;  iRecTracks++) {
 518
 519     // reading info from tracks
 520     recTrack = (AliMUONTrack*) recTracksArray->At(iRecTracks);
 521
 522     trackParam = (AliMUONTrackParam*) (recTrack->GetTrackParamAtHit())->First();
 523
 524     if (esdVert->GetNContributors())
 525       AliMUONTrackExtrap::ExtrapToVertex(trackParam, vertex[0],vertex[1],vertex[2]);
 526
 527     bendingSlope            = trackParam->GetBendingSlope();
 528     nonBendingSlope         = trackParam->GetNonBendingSlope();
 529     inverseBendingMomentum  = trackParam->GetInverseBendingMomentum();
 530     xRec  = trackParam->GetNonBendingCoor();
 531     yRec  = trackParam->GetBendingCoor();
 532     zRec  = trackParam->GetZ();
 533
 534     nTrackHits       = recTrack->GetNTrackHits();
 535     fitFmin          = recTrack->GetFitFMin();
 536     matchTrigger     = recTrack->GetMatchTrigger();
 537     chi2MatchTrigger = recTrack->GetChi2MatchTrigger();
 538
 539     // setting data member of ESD MUON
 540     theESDTrack->SetInverseBendingMomentum(inverseBendingMomentum);
 541     theESDTrack->SetThetaX(TMath::ATan(nonBendingSlope));
 542     theESDTrack->SetThetaY(TMath::ATan(bendingSlope));
 543     theESDTrack->SetZ(zRec);
 544     theESDTrack->SetBendingCoor(yRec); // calculate vertex at ESD or Tracking level ?
 545     theESDTrack->SetNonBendingCoor(xRec);
 546     theESDTrack->SetChi2(fitFmin);
 547     theESDTrack->SetNHit(nTrackHits);
 548     theESDTrack->SetMatchTrigger(matchTrigger);
 549     theESDTrack->SetChi2MatchTrigger(chi2MatchTrigger);
 550
 551     // storing ESD MUON Track into ESD Event
 552     if (nRecTracks != 0)
 553       esd->AddMuonTrack(theESDTrack);
 554   } // end loop tracks
 555
 556   // reset muondata
 557   fMUONData->ResetRecTracks();
 558   fMUONData->ResetRecTriggerTracks();
 559
 560   //} // end loop on event
 561   loader->UnloadTracks();
 562
 563   delete theESDTrack;
 564 }
 565
 566 //_____________________________________________________________________________
 567 void AliMUONReconstructor::FillESD(AliRunLoader* runLoader, AliRawReader* /*rawReader*/, AliESD* esd) const
 568 {
 569 /// Fill ESD
 570 /// \todo add more
 571
 572   // don't need rawReader ???
 573   FillESD(runLoader, esd);
 574 }
 575
 576 //_____________________________________________________________________________
 577 AliTracker* AliMUONReconstructor::CreateTracker(AliRunLoader* runLoader) const
 578 {
 579   /// create tracker for MUON
 580   /// go into the tracking framework finally (Ch.F)
 581
 582   AliLoader* loader = runLoader->GetLoader("MUONLoader");
 583
 584   fMUONData->SetLoader(loader);
 585
 586   AliMUONTracker* tracker = new AliMUONTracker();
 587   tracker->SetMUONData(fMUONData);
 588   tracker->SetTriggerCircuit(fTriggerCircuit);
 589   tracker->SetOption(GetOption());
 590
 591   return tracker;
 592
 593 }