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