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