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 //                                                                           //
  19 // class for MUON reconstruction                                             //
  20 //                                                                           //
  21 ///////////////////////////////////////////////////////////////////////////////
  22 #include <TParticle.h>
  23 #include <TArrayF.h>
  24
  25 #include "AliRunLoader.h"
  26 #include "AliHeader.h"
  27 #include "AliGenEventHeader.h"
  28 #include "AliESD.h"
  29 #include "AliMUONReconstructor.h"
  30
  31 #include "AliMUONData.h"
  32 #include "AliMUONEventReconstructor.h"
  33 #include "AliMUONClusterReconstructor.h"
  34 #include "AliMUONClusterFinderVS.h"
  35 #include "AliMUONTrack.h"
  36 #include "AliMUONTrackParam.h"
  37 #include "AliMUONTriggerTrack.h"
  38 #include "AliESDMuonTrack.h"
  39 #include "AliMUONRawData.h"
  40
  41 #include "AliRawReader.h"
  42
  43 ClassImp(AliMUONReconstructor)
  44 //_____________________________________________________________________________
  45 AliMUONReconstructor::AliMUONReconstructor()
  46   : AliReconstructor()
  47 {
  48 }
  49 //_____________________________________________________________________________
  50 AliMUONReconstructor::~AliMUONReconstructor()
  51 {
  52 }
  53 //_____________________________________________________________________________
  54 void AliMUONReconstructor::Reconstruct(AliRunLoader* runLoader) const
  55 {
  56 //  AliLoader
  57   AliLoader* loader = runLoader->GetLoader("MUONLoader");
  58   Int_t nEvents = runLoader->GetNumberOfEvents();
  59
  60 // used local container for each method
  61 // passing fLoader as argument, could be avoided ???
  62   AliMUONEventReconstructor* recoEvent = new AliMUONEventReconstructor(loader);
  63   AliMUONData* dataEvent = recoEvent->GetMUONData();
  64
  65   AliMUONClusterReconstructor* recoCluster = new AliMUONClusterReconstructor(loader);
  66   AliMUONData* dataCluster = recoCluster->GetMUONData();
  67   AliMUONClusterFinderVS *recModel = recoCluster->GetRecoModel();
  68   recModel->SetGhostChi2Cut(10);
  69
  70   loader->LoadDigits("READ");
  71   loader->LoadRecPoints("RECREATE");
  72   loader->LoadTracks("RECREATE");
  73
  74   //   Loop over events
  75   for(Int_t ievent = 0; ievent < nEvents; ievent++) {
  76     printf("Event %d\n",ievent);
  77     runLoader->GetEvent(ievent);
  78
  79     //----------------------- digit2cluster & Trigger2Trigger -------------------
  80     if (!loader->TreeR()) loader->MakeRecPointsContainer();
  81
  82     // tracking branch
  83     dataCluster->MakeBranch("RC");
  84     dataCluster->SetTreeAddress("D,RC");
  85     recoCluster->Digits2Clusters();
  86     dataCluster->Fill("RC");
  87
  88     // trigger branch
  89     dataCluster->MakeBranch("TC");
  90     dataCluster->SetTreeAddress("TC");
  91     recoCluster->Trigger2Trigger();
  92     dataCluster->Fill("TC");
  93
  94     loader->WriteRecPoints("OVERWRITE");
  95
  96     //---------------------------- Track & TriggerTrack ---------------------
  97     if (!loader->TreeT()) loader->MakeTracksContainer();
  98
  99     // trigger branch
 100     dataEvent->MakeBranch("RL"); //trigger track
 101     dataEvent->SetTreeAddress("RL");
 102     recoEvent->EventReconstructTrigger();
 103     dataEvent->Fill("RL");
 104
 105     // tracking branch
 106     dataEvent->MakeBranch("RT"); //track
 107     dataEvent->SetTreeAddress("RT");
 108     recoEvent->EventReconstruct();
 109     dataEvent->Fill("RT");
 110
 111     loader->WriteTracks("OVERWRITE");
 112
 113     //--------------------------- Resetting branches -----------------------
 114     dataCluster->ResetDigits();
 115     dataCluster->ResetRawClusters();
 116     dataCluster->ResetTrigger();
 117
 118     dataEvent->ResetRawClusters();
 119     dataEvent->ResetTrigger();
 120     dataEvent->ResetRecTracks();
 121     dataEvent->ResetRecTriggerTracks();
 122
 123   }
 124   loader->UnloadDigits();
 125   loader->UnloadRecPoints();
 126   loader->UnloadTracks();
 127
 128   delete recoCluster;
 129   delete recoEvent;
 130 }
 131
 132 //_____________________________________________________________________________
 133 void AliMUONReconstructor::Reconstruct(AliRunLoader* runLoader, AliRawReader* rawReader) const
 134 {
 135 //  AliLoader
 136   AliLoader* loader = runLoader->GetLoader("MUONLoader");
 137
 138 // used local container for each method
 139 // passing fLoader as argument, could be avoided ???
 140   AliMUONEventReconstructor* recoEvent = new AliMUONEventReconstructor(loader);
 141   AliMUONData* dataEvent = recoEvent->GetMUONData();
 142
 143   AliMUONRawData* rawData = new AliMUONRawData(loader);
 144   AliMUONData* dataCluster = rawData->GetMUONData();
 145
 146   AliMUONClusterReconstructor* recoCluster = new AliMUONClusterReconstructor(loader, dataCluster);
 147   AliMUONClusterFinderVS *recModel = recoCluster->GetRecoModel();
 148   recModel->SetGhostChi2Cut(10);
 149
 150   loader->LoadRecPoints("RECREATE");
 151   loader->LoadTracks("RECREATE");
 152
 153   //   Loop over events
 154   Int_t iEvent = 0;
 155
 156   while (rawReader->NextEvent()) {
 157     printf("Event %d\n",iEvent);
 158     runLoader->GetEvent(iEvent++);
 159     rawData->Raw2Digits(rawReader);
 160
 161     //----------------------- digit2cluster & Trigger2Trigger -------------------
 162     if (!loader->TreeR()) loader->MakeRecPointsContainer();
 163
 164     // tracking branch
 165     dataCluster->MakeBranch("RC");
 166     dataCluster->SetTreeAddress("D,RC");
 167     recoCluster->Digits2Clusters();
 168     dataCluster->Fill("RC");
 169
 170     // trigger branch
 171     dataCluster->MakeBranch("TC");
 172     dataCluster->SetTreeAddress("TC");
 173     recoCluster->Trigger2Trigger();
 174     dataCluster->Fill("TC");
 175
 176     loader->WriteRecPoints("OVERWRITE");
 177
 178     //---------------------------- Track & TriggerTrack ---------------------
 179     if (!loader->TreeT()) loader->MakeTracksContainer();
 180
 181     // trigger branch
 182     dataEvent->MakeBranch("RL"); //trigger track
 183     dataEvent->SetTreeAddress("RL");
 184     recoEvent->EventReconstructTrigger();
 185     dataEvent->Fill("RL");
 186
 187     // tracking branch
 188     dataEvent->MakeBranch("RT"); //track
 189     dataEvent->SetTreeAddress("RT");
 190     recoEvent->EventReconstruct();
 191     dataEvent->Fill("RT");
 192
 193     loader->WriteTracks("OVERWRITE");
 194
 195     //--------------------------- Resetting branches -----------------------
 196     dataCluster->ResetDigits();
 197     dataCluster->ResetRawClusters();
 198     dataCluster->ResetTrigger();
 199
 200     dataEvent->ResetRawClusters();
 201     dataEvent->ResetTrigger();
 202     dataEvent->ResetRecTracks();
 203     dataEvent->ResetRecTriggerTracks();
 204
 205   }
 206   loader->UnloadRecPoints();
 207   loader->UnloadTracks();
 208
 209   delete recoCluster;
 210   delete recoEvent;
 211 }
 212
 213 //_____________________________________________________________________________
 214 void AliMUONReconstructor::FillESD(AliRunLoader* runLoader, AliESD* esd) const
 215 {
 216   TClonesArray* recTracksArray = 0;
 217   TClonesArray* recTrigTracksArray = 0;
 218
 219   AliLoader* loader = runLoader->GetLoader("MUONLoader");
 220   loader->LoadTracks("READ");
 221   AliMUONData* muonData = new AliMUONData(loader,"MUON","MUON");
 222
 223    // declaration
 224   Int_t iEvent, nPart;
 225   Int_t nTrackHits, nPrimary;
 226   Double_t fitFmin;
 227   TArrayF vertex(3);
 228
 229   Double_t bendingSlope, nonBendingSlope, inverseBendingMomentum;
 230   Double_t xRec, yRec, zRec, chi2MatchTrigger;
 231   Bool_t matchTrigger;
 232
 233   // setting pointer for tracks, triggertracks & trackparam at vertex
 234   AliMUONTrack* recTrack = 0;
 235   AliMUONTrackParam* trackParam = 0;
 236   AliMUONTriggerTrack* recTriggerTrack = 0;
 237   TParticle* particle = new TParticle();
 238   AliGenEventHeader* header = 0;
 239   iEvent = runLoader->GetEventNumber();
 240   runLoader->GetEvent(iEvent);
 241
 242   // vertex calculation (maybe it exists already somewhere else)
 243   vertex[0] = vertex[1] = vertex[2] = 0.;
 244   nPrimary = 0;
 245   if ( (header = runLoader->GetHeader()->GenEventHeader()) ) {
 246     header->PrimaryVertex(vertex);
 247   } else {
 248     runLoader->LoadKinematics("READ");
 249     runLoader->TreeK()->GetBranch("Particles")->SetAddress(&particle);
 250     nPart = (Int_t)runLoader->TreeK()->GetEntries();
 251     for(Int_t iPart = 0; iPart < nPart; iPart++) {
 252       runLoader->TreeK()->GetEvent(iPart);
 253       if (particle->GetFirstMother() == -1) {
 254         vertex[0] += particle->Vx();
 255         vertex[1] += particle->Vy();
 256         vertex[2] += particle->Vz();
 257         nPrimary++;
 258       }
 259       if (nPrimary) {
 260         vertex[0] /= (double)nPrimary;
 261         vertex[1] /= (double)nPrimary;
 262         vertex[2] /= (double)nPrimary;
 263       }
 264     }
 265   }
 266   // setting ESD MUON class
 267   AliESDMuonTrack* theESDTrack = new  AliESDMuonTrack() ;
 268
 269   //-------------------- trigger tracks-------------
 270   Long_t trigPat = 0;
 271   muonData->SetTreeAddress("RL");
 272   muonData->GetRecTriggerTracks();
 273   recTrigTracksArray = muonData->RecTriggerTracks();
 274
 275   // ready global trigger pattern from first track
 276   if (recTrigTracksArray)
 277     recTriggerTrack = (AliMUONTriggerTrack*) recTrigTracksArray->First();
 278   if (recTriggerTrack) trigPat = recTriggerTrack->GetGTPattern();
 279
 280   //printf(">>> Event %d Number of Recconstructed tracks %d \n",iEvent, nrectracks);
 281
 282   // -------------------- tracks-------------
 283   muonData->SetTreeAddress("RT");
 284   muonData->GetRecTracks();
 285   recTracksArray = muonData->RecTracks();
 286
 287   Int_t nRecTracks = 0;
 288   if (recTracksArray)
 289     nRecTracks = (Int_t) recTracksArray->GetEntriesFast(); //
 290
 291   // loop over tracks
 292   for (Int_t iRecTracks = 0; iRecTracks <  nRecTracks;  iRecTracks++) {
 293
 294     // reading info from tracks
 295     recTrack = (AliMUONTrack*) recTracksArray->At(iRecTracks);
 296
 297     trackParam = (AliMUONTrackParam*) (recTrack->GetTrackParamAtHit())->First();
 298     trackParam->ExtrapToVertex(vertex[0],vertex[1],vertex[2]);
 299
 300     bendingSlope            = trackParam->GetBendingSlope();
 301     nonBendingSlope         = trackParam->GetNonBendingSlope();
 302     inverseBendingMomentum = trackParam->GetInverseBendingMomentum();
 303     xRec  = trackParam->GetNonBendingCoor();
 304     yRec  = trackParam->GetBendingCoor();
 305     zRec  = trackParam->GetZ();
 306
 307     nTrackHits       = recTrack->GetNTrackHits();
 308     fitFmin          = recTrack->GetFitFMin();
 309     matchTrigger     = recTrack->GetMatchTrigger();
 310     chi2MatchTrigger = recTrack->GetChi2MatchTrigger();
 311
 312     // setting data member of ESD MUON
 313     theESDTrack->SetInverseBendingMomentum(inverseBendingMomentum);
 314     theESDTrack->SetThetaX(TMath::ATan(nonBendingSlope));
 315     theESDTrack->SetThetaY(TMath::ATan(bendingSlope));
 316     theESDTrack->SetZ(zRec);
 317     theESDTrack->SetBendingCoor(yRec); // calculate vertex at ESD or Tracking level ?
 318     theESDTrack->SetNonBendingCoor(xRec);
 319     theESDTrack->SetChi2(fitFmin);
 320     theESDTrack->SetNHit(nTrackHits);
 321     theESDTrack->SetMatchTrigger(matchTrigger);
 322     theESDTrack->SetChi2MatchTrigger(chi2MatchTrigger);
 323
 324     // storing ESD MUON Track into ESD Event
 325     if (nRecTracks != 0)
 326       esd->AddMuonTrack(theESDTrack);
 327   } // end loop tracks
 328
 329   // add global trigger pattern
 330   if (nRecTracks != 0)
 331     esd->SetTrigger(trigPat);
 332
 333   // reset muondata
 334   muonData->ResetRecTracks();
 335   muonData->ResetRecTriggerTracks();
 336
 337   //} // end loop on event
 338   loader->UnloadTracks();
 339   if (!header)
 340     runLoader->UnloadKinematics();
 341   delete theESDTrack;
 342   delete muonData;
 343   // delete particle;
 344 }