MUON/AliMUONReconstructor.cxx

   1 /**************************************************************************
   2
   3  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   4
   5  *                                                                        *
   6
   7  * Author: The ALICE Off-line Project.                                    *
   8
   9  * Contributors are mentioned in the code where appropriate.              *
  10
  11  *                                                                        *
  12
  13  * Permission to use, copy, modify and distribute this software and its   *
  14
  15  * documentation strictly for non-commercial purposes is hereby granted   *
  16
  17  * without fee, provided that the above copyright notice appears in all   *
  18
  19  * copies and that both the copyright notice and this permission notice   *
  20
  21  * appear in the supporting documentation. The authors make no claims     *
  22
  23  * about the suitability of this software for any purpose. It is          *
  24
  25  * provided "as is" without express or implied warranty.                  *
  26
  27  **************************************************************************/
  28
  29 /* $Id$ */
  30
  31
  32
  33 ///////////////////////////////////////////////////////////////////////////////
  34
  35 //                                                                           //
  36
  37 // class for MUON reconstruction                                             //
  38
  39 //                                                                           //
  40
  41 ///////////////////////////////////////////////////////////////////////////////
  42
  43 #include <TParticle.h>
  44
  45 #include <TArrayF.h>
  46
  47
  48
  49 #include "AliRunLoader.h"
  50
  51 #include "AliHeader.h"
  52
  53 #include "AliGenEventHeader.h"
  54
  55 #include "AliESD.h"
  56
  57 #include "AliMUONReconstructor.h"
  58
  59
  60
  61 #include "AliMUONData.h"
  62
  63 #include "AliMUONTrackReconstructor.h"
  64
  65 #include "AliMUONClusterReconstructor.h"
  66
  67 #include "AliMUONClusterFinderVS.h"
  68
  69 #include "AliMUONClusterFinderAZ.h"
  70
  71 #include "AliMUONTrack.h"
  72
  73 #include "AliMUONTrackParam.h"
  74
  75 #include "AliMUONTriggerTrack.h"
  76
  77 #include "AliESDMuonTrack.h"
  78
  79 #include "AliMUONRawData.h"
  80
  81
  82
  83 #include "AliRawReader.h"
  84
  85
  86
  87
  88
  89 ClassImp(AliMUONReconstructor)
  90
  91 //_____________________________________________________________________________
  92
  93 AliMUONReconstructor::AliMUONReconstructor()
  94
  95   : AliReconstructor()
  96
  97 {
  98
  99 }
 100
 101 //_____________________________________________________________________________
 102
 103 AliMUONReconstructor::~AliMUONReconstructor()
 104
 105 {
 106
 107 }
 108
 109 //_____________________________________________________________________________
 110
 111 void AliMUONReconstructor::Reconstruct(AliRunLoader* runLoader) const
 112
 113 {
 114
 115 //  AliLoader
 116
 117   AliLoader* loader = runLoader->GetLoader("MUONLoader");
 118
 119   Int_t nEvents = runLoader->GetNumberOfEvents();
 120
 121
 122
 123 // used local container for each method
 124
 125 // passing fLoader as argument, could be avoided ???
 126
 127   AliMUONTrackReconstructor* recoEvent = new AliMUONTrackReconstructor(loader);
 128
 129   AliMUONData* dataEvent = recoEvent->GetMUONData();
 130
 131   if (strstr(GetOption(),"Kalman")) recoEvent->SetTrackMethod(2); // Kalman
 132
 133   else recoEvent->SetTrackMethod(1); // original
 134
 135
 136
 137   AliMUONClusterReconstructor* recoCluster = new AliMUONClusterReconstructor(loader);
 138
 139   AliMUONData* dataCluster = recoCluster->GetMUONData();
 140
 141   AliMUONClusterFinderVS *recModel = recoCluster->GetRecoModel();
 142
 143   if (strstr(GetOption(),"AZ")) {
 144
 145     recModel = (AliMUONClusterFinderVS*) new AliMUONClusterFinderAZ(0,1);
 146
 147     recoCluster->SetRecoModel(recModel);
 148
 149   }
 150
 151   recModel->SetGhostChi2Cut(10);
 152
 153
 154
 155   loader->LoadDigits("READ");
 156
 157   loader->LoadRecPoints("RECREATE");
 158
 159   loader->LoadTracks("RECREATE");
 160
 161
 162
 163   //   Loop over events
 164
 165   for(Int_t ievent = 0; ievent < nEvents; ievent++) {
 166
 167     printf("Event %d\n",ievent);
 168
 169     runLoader->GetEvent(ievent);
 170
 171
 172
 173     //----------------------- digit2cluster & Trigger2Trigger -------------------
 174
 175     if (!loader->TreeR()) loader->MakeRecPointsContainer();
 176
 177
 178
 179     // tracking branch
 180
 181     dataCluster->MakeBranch("RC");
 182
 183     dataCluster->SetTreeAddress("D,RC");
 184
 185     recoCluster->Digits2Clusters();
 186
 187     dataCluster->Fill("RC");
 188
 189
 190
 191     // trigger branch
 192
 193     dataCluster->MakeBranch("TC");
 194
 195     dataCluster->SetTreeAddress("TC");
 196
 197     recoCluster->Trigger2Trigger();
 198
 199     dataCluster->Fill("TC");
 200
 201
 202
 203     loader->WriteRecPoints("OVERWRITE");
 204
 205
 206
 207     //---------------------------- Track & TriggerTrack ---------------------
 208
 209     if (!loader->TreeT()) loader->MakeTracksContainer();
 210
 211
 212
 213     // trigger branch
 214
 215     dataEvent->MakeBranch("RL"); //trigger track
 216
 217     dataEvent->SetTreeAddress("RL");
 218
 219     recoEvent->EventReconstructTrigger();
 220
 221     dataEvent->Fill("RL");
 222
 223
 224
 225     // tracking branch
 226
 227     dataEvent->MakeBranch("RT"); //track
 228
 229     dataEvent->SetTreeAddress("RT");
 230
 231     recoEvent->EventReconstruct();
 232
 233     dataEvent->Fill("RT");
 234
 235
 236
 237     loader->WriteTracks("OVERWRITE");
 238
 239
 240
 241     //--------------------------- Resetting branches -----------------------
 242
 243     dataCluster->ResetDigits();
 244
 245     dataCluster->ResetRawClusters();
 246
 247     dataCluster->ResetTrigger();
 248
 249
 250
 251     dataEvent->ResetRawClusters();
 252
 253     dataEvent->ResetTrigger();
 254
 255     dataEvent->ResetRecTracks();
 256
 257     dataEvent->ResetRecTriggerTracks();
 258
 259
 260
 261   }
 262
 263   loader->UnloadDigits();
 264
 265   loader->UnloadRecPoints();
 266
 267   loader->UnloadTracks();
 268
 269
 270
 271   delete recoCluster;
 272
 273   delete recoEvent;
 274
 275 }
 276
 277
 278
 279 //_____________________________________________________________________________
 280
 281 void AliMUONReconstructor::Reconstruct(AliRunLoader* runLoader, AliRawReader* rawReader) const
 282
 283 {
 284
 285 //  AliLoader
 286
 287   AliLoader* loader = runLoader->GetLoader("MUONLoader");
 288
 289
 290
 291 // used local container for each method
 292
 293 // passing fLoader as argument, could be avoided ???
 294
 295   AliMUONTrackReconstructor* recoEvent = new AliMUONTrackReconstructor(loader);
 296
 297   AliMUONData* dataEvent = recoEvent->GetMUONData();
 298
 299
 300
 301   AliMUONRawData* rawData = new AliMUONRawData(loader);
 302
 303   AliMUONData* dataCluster = rawData->GetMUONData();
 304
 305
 306
 307   AliMUONClusterReconstructor* recoCluster = new AliMUONClusterReconstructor(loader, dataCluster);
 308
 309   AliMUONClusterFinderVS *recModel = recoCluster->GetRecoModel();
 310
 311   recModel->SetGhostChi2Cut(10);
 312
 313
 314
 315   loader->LoadRecPoints("RECREATE");
 316
 317   loader->LoadTracks("RECREATE");
 318
 319   loader->LoadDigits("RECREATE");
 320
 321
 322
 323
 324
 325   //   Loop over events
 326
 327   Int_t iEvent = 0;
 328
 329
 330
 331   while (rawReader->NextEvent()) {
 332
 333     printf("Event %d\n",iEvent);
 334
 335     runLoader->GetEvent(iEvent++);
 336
 337
 338
 339     //----------------------- raw2digits & raw2trigger-------------------
 340
 341     if (!loader->TreeD()) loader->MakeDigitsContainer();
 342
 343
 344
 345     // tracking branch
 346
 347     dataCluster->MakeBranch("D");
 348
 349     dataCluster->SetTreeAddress("D");
 350
 351     rawData->ReadTrackerDDL(rawReader);
 352
 353     dataCluster->Fill("D");
 354
 355
 356
 357     // trigger branch
 358
 359     dataCluster->MakeBranch("GLT");
 360
 361     dataCluster->SetTreeAddress("GLT");
 362
 363     rawData->ReadTriggerDDL(rawReader);
 364
 365     dataCluster->Fill("GLT");
 366
 367
 368
 369     loader->WriteDigits("OVERWRITE");
 370
 371
 372
 373     //----------------------- digit2cluster & Trigger2Trigger -------------------
 374
 375     if (!loader->TreeR()) loader->MakeRecPointsContainer();
 376
 377
 378
 379     // tracking branch
 380
 381     dataCluster->MakeBranch("RC");
 382
 383     dataCluster->SetTreeAddress("RC");
 384
 385     recoCluster->Digits2Clusters();
 386
 387     dataCluster->Fill("RC");
 388
 389
 390
 391     // trigger branch
 392
 393     dataCluster->MakeBranch("TC");
 394
 395     dataCluster->SetTreeAddress("TC");
 396
 397     recoCluster->Trigger2Trigger();
 398
 399     dataCluster->Fill("TC");
 400
 401
 402
 403     loader->WriteRecPoints("OVERWRITE");
 404
 405
 406
 407     //---------------------------- Track & TriggerTrack ---------------------
 408
 409     if (!loader->TreeT()) loader->MakeTracksContainer();
 410
 411
 412
 413     // trigger branch
 414
 415     dataEvent->MakeBranch("RL"); //trigger track
 416
 417     dataEvent->SetTreeAddress("RL");
 418
 419     recoEvent->EventReconstructTrigger();
 420
 421     dataEvent->Fill("RL");
 422
 423
 424
 425     // tracking branch
 426
 427     dataEvent->MakeBranch("RT"); //track
 428
 429     dataEvent->SetTreeAddress("RT");
 430
 431     recoEvent->EventReconstruct();
 432
 433     dataEvent->Fill("RT");
 434
 435
 436
 437     loader->WriteTracks("OVERWRITE");
 438
 439
 440
 441     //--------------------------- Resetting branches -----------------------
 442
 443     dataCluster->ResetDigits();
 444
 445     dataCluster->ResetRawClusters();
 446
 447     dataCluster->ResetTrigger();
 448
 449
 450
 451     dataEvent->ResetRawClusters();
 452
 453     dataEvent->ResetTrigger();
 454
 455     dataEvent->ResetRecTracks();
 456
 457     dataEvent->ResetRecTriggerTracks();
 458
 459
 460
 461   }
 462
 463   loader->UnloadRecPoints();
 464
 465   loader->UnloadTracks();
 466
 467   loader->UnloadDigits();
 468
 469
 470
 471   delete recoCluster;
 472
 473   delete recoEvent;
 474
 475 }
 476
 477
 478
 479 //_____________________________________________________________________________
 480
 481 void AliMUONReconstructor::FillESD(AliRunLoader* runLoader, AliESD* esd) const
 482
 483 {
 484
 485   TClonesArray* recTracksArray = 0;
 486
 487   TClonesArray* recTrigTracksArray = 0;
 488
 489
 490
 491   AliLoader* loader = runLoader->GetLoader("MUONLoader");
 492
 493   loader->LoadTracks("READ");
 494
 495   AliMUONData* muonData = new AliMUONData(loader,"MUON","MUON");
 496
 497
 498
 499    // declaration
 500
 501   Int_t iEvent;// nPart;
 502
 503   Int_t nTrackHits;// nPrimary;
 504
 505   Double_t fitFmin;
 506
 507   TArrayF vertex(3);
 508
 509
 510
 511   Double_t bendingSlope, nonBendingSlope, inverseBendingMomentum;
 512
 513   Double_t xRec, yRec, zRec, chi2MatchTrigger;
 514
 515   Bool_t matchTrigger;
 516
 517
 518
 519   // setting pointer for tracks, triggertracks & trackparam at vertex
 520
 521   AliMUONTrack* recTrack = 0;
 522
 523   AliMUONTrackParam* trackParam = 0;
 524
 525   AliMUONTriggerTrack* recTriggerTrack = 0;
 526
 527 //   TParticle* particle = new TParticle();
 528
 529 //   AliGenEventHeader* header = 0;
 530
 531   iEvent = runLoader->GetEventNumber();
 532
 533   runLoader->GetEvent(iEvent);
 534
 535
 536
 537   // vertex calculation (maybe it exists already somewhere else)
 538
 539   vertex[0] = vertex[1] = vertex[2] = 0.;
 540
 541  //  nPrimary = 0;
 542
 543 //   if ( (header = runLoader->GetHeader()->GenEventHeader()) ) {
 544
 545 //     header->PrimaryVertex(vertex);
 546
 547 //   } else {
 548
 549 //     runLoader->LoadKinematics("READ");
 550
 551 //     runLoader->TreeK()->GetBranch("Particles")->SetAddress(&particle);
 552
 553 //     nPart = (Int_t)runLoader->TreeK()->GetEntries();
 554
 555 //     for(Int_t iPart = 0; iPart < nPart; iPart++) {
 556
 557 //       runLoader->TreeK()->GetEvent(iPart);
 558
 559 //       if (particle->GetFirstMother() == -1) {
 560
 561 //      vertex[0] += particle->Vx();
 562
 563 //      vertex[1] += particle->Vy();
 564
 565 //      vertex[2] += particle->Vz();
 566
 567 //      nPrimary++;
 568
 569 //       }
 570
 571 //       if (nPrimary) {
 572
 573 //      vertex[0] /= (double)nPrimary;
 574
 575 //      vertex[1] /= (double)nPrimary;
 576
 577 //      vertex[2] /= (double)nPrimary;
 578
 579 //       }
 580
 581 //     }
 582
 583 //   }
 584
 585   // setting ESD MUON class
 586
 587   AliESDMuonTrack* theESDTrack = new  AliESDMuonTrack() ;
 588
 589
 590
 591   //-------------------- trigger tracks-------------
 592
 593   Long_t trigPat = 0;
 594
 595   muonData->SetTreeAddress("RL");
 596
 597   muonData->GetRecTriggerTracks();
 598
 599   recTrigTracksArray = muonData->RecTriggerTracks();
 600
 601
 602
 603   // ready global trigger pattern from first track
 604
 605   if (recTrigTracksArray)
 606
 607     recTriggerTrack = (AliMUONTriggerTrack*) recTrigTracksArray->First();
 608
 609   if (recTriggerTrack) trigPat = recTriggerTrack->GetGTPattern();
 610
 611
 612
 613   //printf(">>> Event %d Number of Recconstructed tracks %d \n",iEvent, nrectracks);
 614
 615
 616
 617   // -------------------- tracks-------------
 618
 619   muonData->SetTreeAddress("RT");
 620
 621   muonData->GetRecTracks();
 622
 623   recTracksArray = muonData->RecTracks();
 624
 625
 626
 627   Int_t nRecTracks = 0;
 628
 629   if (recTracksArray)
 630
 631     nRecTracks = (Int_t) recTracksArray->GetEntriesFast(); //
 632
 633
 634
 635   // loop over tracks
 636
 637   for (Int_t iRecTracks = 0; iRecTracks <  nRecTracks;  iRecTracks++) {
 638
 639
 640
 641     // reading info from tracks
 642
 643     recTrack = (AliMUONTrack*) recTracksArray->At(iRecTracks);
 644
 645
 646
 647     trackParam = (AliMUONTrackParam*) (recTrack->GetTrackParamAtHit())->First();
 648
 649     trackParam->ExtrapToVertex(vertex[0],vertex[1],vertex[2]);
 650
 651
 652
 653     bendingSlope            = trackParam->GetBendingSlope();
 654
 655     nonBendingSlope         = trackParam->GetNonBendingSlope();
 656
 657     inverseBendingMomentum = trackParam->GetInverseBendingMomentum();
 658
 659     xRec  = trackParam->GetNonBendingCoor();
 660
 661     yRec  = trackParam->GetBendingCoor();
 662
 663     zRec  = trackParam->GetZ();
 664
 665
 666
 667     nTrackHits       = recTrack->GetNTrackHits();
 668
 669     fitFmin          = recTrack->GetFitFMin();
 670
 671     matchTrigger     = recTrack->GetMatchTrigger();
 672
 673     chi2MatchTrigger = recTrack->GetChi2MatchTrigger();
 674
 675
 676
 677     // setting data member of ESD MUON
 678
 679     theESDTrack->SetInverseBendingMomentum(inverseBendingMomentum);
 680
 681     theESDTrack->SetThetaX(TMath::ATan(nonBendingSlope));
 682
 683     theESDTrack->SetThetaY(TMath::ATan(bendingSlope));
 684
 685     theESDTrack->SetZ(zRec);
 686
 687     theESDTrack->SetBendingCoor(yRec); // calculate vertex at ESD or Tracking level ?
 688
 689     theESDTrack->SetNonBendingCoor(xRec);
 690
 691     theESDTrack->SetChi2(fitFmin);
 692
 693     theESDTrack->SetNHit(nTrackHits);
 694
 695     theESDTrack->SetMatchTrigger(matchTrigger);
 696
 697     theESDTrack->SetChi2MatchTrigger(chi2MatchTrigger);
 698
 699
 700
 701     // storing ESD MUON Track into ESD Event
 702
 703     if (nRecTracks != 0)
 704
 705       esd->AddMuonTrack(theESDTrack);
 706
 707   } // end loop tracks
 708
 709
 710
 711   // add global trigger pattern
 712
 713   if (nRecTracks != 0)
 714
 715     esd->SetTrigger(trigPat);
 716
 717
 718
 719   // reset muondata
 720
 721   muonData->ResetRecTracks();
 722
 723   muonData->ResetRecTriggerTracks();
 724
 725
 726
 727   //} // end loop on event
 728
 729   loader->UnloadTracks();
 730
 731  //  if (!header)
 732
 733 //     runLoader->UnloadKinematics();
 734
 735   delete theESDTrack;
 736
 737   delete muonData;
 738
 739   // delete particle;
 740
 741 }//_____________________________________________________________________________
 742
 743 void AliMUONReconstructor::FillESD(AliRunLoader* runLoader, AliRawReader* /*rawReader*/, AliESD* esd) const
 744
 745 {
 746
 747   // don't need rawReader ???
 748
 749   FillESD(runLoader, esd);
 750
 751 }
 752